dolfin team mailing list archive

[Anders Logg <logg@xxxxxxxxx>]

This looks a like a non-optimal solution. Are you returning the mesh
by value? That will lead to creation of a temporary and copying of the
entire refined mesh.

What if we have

  FooRefinement::refine(old_mesh, new_mesh)

or similar?

--
Anders


On Sun, Feb 07, 2010 at 10:13:15PM -0000, noreply@xxxxxxxxxxxxx wrote:
> ------------------------------------------------------------
> revno: 4461
> committer: Garth N. Wells <gnw20@xxxxxxxxx>
> branch nick: dolfin-nompi
> timestamp: Sun 2010-02-07 21:53:42 +0000
> message:
>   Work on mesh refinement interface.
>
>   Do not overwrite old mesh inside refinement functions belonging to static classes. Need to keep old mesh to peform update in adapative methods, so previously the mesh had to be copied before calling refine(...).
> modified:
>   dolfin/adaptivity/AdaptiveObjects.cpp
>   dolfin/mesh/LocalMeshRefinement.cpp
>   dolfin/mesh/LocalMeshRefinement.h
>   dolfin/mesh/RivaraRefinement.cpp
>   dolfin/mesh/RivaraRefinement.h
>   dolfin/mesh/UniformMeshRefinement.cpp
>   dolfin/mesh/UniformMeshRefinement.h
>
>

> === modified file 'dolfin/adaptivity/AdaptiveObjects.cpp'
> --- dolfin/adaptivity/AdaptiveObjects.cpp	2009-11-27 01:59:42 +0000
> +++ dolfin/adaptivity/AdaptiveObjects.cpp	2010-02-07 21:53:42 +0000
> @@ -129,21 +129,21 @@
>    // Type of refinement (should perhaps be a parameter)
>    const bool recursive_refinement = true;
>
> -  // Create new mesh (copy of old mesh)
> -  Mesh new_mesh(*mesh);
> +  // Create new mesh
> +  Mesh new_mesh;
>
>    // Refine new mesh
>    if (cell_markers)
>    {
>      if (recursive_refinement)
> -      LocalMeshRefinement::refineRecursivelyByEdgeBisection(new_mesh, *cell_markers);
> +      new_mesh = LocalMeshRefinement::refineRecursivelyByEdgeBisection(*mesh, *cell_markers);
>      else
> -      LocalMeshRefinement::refineIterativelyByEdgeBisection(new_mesh, *cell_markers);
> +      new_mesh = LocalMeshRefinement::refineIterativelyByEdgeBisection(*mesh, *cell_markers);
>    }
>    else
>    {
>      info("No cells marked for refinement, assuming uniform mesh refinement.");
> -    UniformMeshRefinement::refine(new_mesh);
> +    new_mesh = UniformMeshRefinement::refine(*mesh);
>    }
>
>    // Refined mesh may not be ordered
> @@ -154,7 +154,6 @@
>
>    // Copy data from new mesh to old
>    *mesh = new_mesh;
> -
>  }
>  //-----------------------------------------------------------------------------
>  void AdaptiveObjects::refine(FunctionSpace* function_space,
>
> === modified file 'dolfin/mesh/LocalMeshRefinement.cpp'
> --- dolfin/mesh/LocalMeshRefinement.cpp	2009-10-08 19:11:44 +0000
> +++ dolfin/mesh/LocalMeshRefinement.cpp	2010-02-07 21:53:42 +0000
> @@ -34,89 +34,74 @@
>  };
>
>  //-----------------------------------------------------------------------------
> -void LocalMeshRefinement::refineMeshByEdgeBisection(Mesh& mesh,
> +dolfin::Mesh LocalMeshRefinement::refineMeshByEdgeBisection(const Mesh& mesh,
>                                                      MeshFunction<bool>& cell_marker,
>                                                      bool refine_boundary)
>  {
>    info("Refining simplicial mesh by edge bisection.");
>
> +  // Copy old mesh
> +  Mesh old_mesh(mesh);
> +
>    // Get size of old mesh
> -  const uint num_vertices = mesh.size(0);
> -  const uint num_cells = mesh.size(mesh.topology().dim());
> +  const uint num_vertices = old_mesh.size(0);
> +  const uint num_cells = old_mesh.size(mesh.topology().dim());
>
>    // Check cell marker
>    if (cell_marker.size() != num_cells)
>      error("Wrong dimension of cell_marker");
>
>    // Generate cell - edge connectivity if not generated
> -  mesh.init(mesh.topology().dim(), 1);
> +  old_mesh.init(mesh.topology().dim(), 1);
>
>    // Generate edge - vertex connectivity if not generated
> -  mesh.init(1, 0);
> +  old_mesh.init(1, 0);
>
>    // Get cell type
> -  const CellType& cell_type = mesh.type();
> +  const CellType& cell_type = old_mesh.type();
>
>    // Init new vertices and cells
>    uint num_new_vertices = 0;
>    uint num_new_cells = 0;
>
> -  // Compute number of vertices and cells
> -  /*
> -  for (CellIterator c(mesh); !c.end(); ++c)
> -  {
> -    if((cell_marker.get(*c) == true))
> -    {
> -      //cout << "marked cell: " << endl;
> -      //cout << c->midpoint() << endl;
> -      //cout << c->index() << endl;
> -    }
> -  }
> -  */
> -
>    // Create new mesh and open for editing
>    Mesh refined_mesh;
>    MeshEditor editor;
>    editor.open(refined_mesh, cell_type.cell_type(),
> -              mesh.topology().dim(), mesh.geometry().dim());
> +              old_mesh.topology().dim(), old_mesh.geometry().dim());
>
>    // Initialize mappings
> -  std::vector<int> old2new_cell(mesh.num_cells());
> -  std::vector<int> old2new_vertex(mesh.num_vertices());
> +  std::vector<int> old2new_cell(old_mesh.num_cells());
> +  std::vector<int> old2new_vertex(old_mesh.num_vertices());
>
>    // Initialise forbidden edges
> -  MeshFunction<bool> edge_forbidden(mesh);
> -  edge_forbidden.init(1);
> -  for (EdgeIterator e(mesh); !e.end(); ++e)
> -    edge_forbidden[e->index()] = false;
> +  MeshFunction<bool> edge_forbidden(old_mesh, true);
> +  for (EdgeIterator e(old_mesh); !e.end(); ++e)
> +    edge_forbidden[*e] = false;
>
>    // If refinement of boundary is forbidden
>    if ( !refine_boundary )
>    {
> -    BoundaryMesh boundary(mesh);
> +    BoundaryMesh boundary(old_mesh);
>      for (EdgeIterator e(boundary); !e.end(); ++e)
> -      edge_forbidden[e->index()] = true;
> +      edge_forbidden[*e] = true;
>    }
>
>    // Initialise forbidden cells
> -  MeshFunction<bool> cell_forbidden(mesh);
> -  cell_forbidden.init(mesh.topology().dim());
> -  for (CellIterator c(mesh); !c.end(); ++c)
> -    cell_forbidden[c->index()] = false;
> +  MeshFunction<bool> cell_forbidden(old_mesh);
> +  cell_forbidden.init(old_mesh.topology().dim());
> +  for (CellIterator c(old_mesh); !c.end(); ++c)
> +    cell_forbidden[*c] = false;
>
>    // Initialise data for finding longest edge
>    uint longest_edge_index = 0;
>    double lmax, l;
>
>    // Compute number of vertices and cells
> -  for (CellIterator c(mesh); !c.end(); ++c)
> +  for (CellIterator c(old_mesh); !c.end(); ++c)
>    {
>      if ( cell_marker[*c] && !cell_forbidden[*c] )
>      {
> -//       cout << "marked cell: " << endl;
> -//       cout << c->midpoint() << endl;
> -//       cout << c->index() << endl;
> -
>        // Find longest edge of cell c
>        lmax = 0.0;
>        for (EdgeIterator e(*c); !e.end(); ++e)
> @@ -132,7 +117,7 @@
>          }
>        }
>
> -      Edge longest_edge(mesh,longest_edge_index);
> +      Edge longest_edge(old_mesh,longest_edge_index);
>
>        // If at least one edge should be bisected
>        if ( lmax > 0.0 )
> @@ -147,7 +132,6 @@
>              // Count new cells
>              num_new_cells++;
>              // set markers of all cell neighbors of longest edge to false
> -            //if ( cn->index() != c->index() )
>              cell_forbidden[cn->index()] = true;
>              // set all the edges of cell neighbors to forbidden
>              for (EdgeIterator en(*cn); !en.end(); ++en)
> @@ -162,18 +146,15 @@
>    editor.init_vertices(num_vertices + num_new_vertices);
>    editor.init_cells(num_cells + num_new_cells);
>
> -  //cout << "Number of cells in old mesh: " << num_cells << "; to add: " << num_new_cells << endl;
> -  //cout << "Number of vertices in old mesh: " << num_vertices << "; to add: " << num_new_vertices << endl;
> -
>    // Add old vertices
>    uint current_vertex = 0;
> -  for (VertexIterator v(mesh); !v.end(); ++v)
> +  for (VertexIterator v(old_mesh); !v.end(); ++v)
>      editor.add_vertex(current_vertex++, v->point());
>
>    // Add old unrefined cells
>    uint current_cell = 0;
>    std::vector<uint> cell_vertices(cell_type.num_entities(0));
> -  for (CellIterator c(mesh); !c.end(); ++c)
> +  for (CellIterator c(old_mesh); !c.end(); ++c)
>    {
>      //if ( (cell_marker.get(*c) == false) && (cell_forbidden.get(*c) == false) )
>      if ( !cell_forbidden[*c] )
> @@ -186,23 +167,23 @@
>    }
>
>    // Reset forbidden edges
> -  for (EdgeIterator e(mesh); !e.end(); ++e)
> +  for (EdgeIterator e(old_mesh); !e.end(); ++e)
>      edge_forbidden[e->index()] = false;
>
>    // If refinement of boundary is forbidden
>    if ( !refine_boundary )
>    {
> -    BoundaryMesh boundary(mesh);
> +    BoundaryMesh boundary(old_mesh);
>      for (EdgeIterator e(boundary); !e.end(); ++e)
>        edge_forbidden[e->index()] = true;
>    }
>
>    // Reset forbidden cells
> -  for (CellIterator c(mesh); !c.end(); ++c)
> +  for (CellIterator c(old_mesh); !c.end(); ++c)
>      cell_forbidden[c->index()] = false;
>
>    // Add new vertices and cells.
> -  for (CellIterator c(mesh); !c.end(); ++c)
> +  for (CellIterator c(old_mesh); !c.end(); ++c)
>    {
>      if ( cell_marker[*c] && !cell_forbidden[*c] )
>      {
> @@ -221,7 +202,7 @@
>          }
>        }
>
> -      Edge longest_edge(mesh,longest_edge_index);
> +      Edge longest_edge(old_mesh,longest_edge_index);
>
>        // If at least one edge should be bisected
>        if ( lmax > 0.0 )
> @@ -235,7 +216,6 @@
>           bisectEdgeOfSimplexCell(*cn, longest_edge, current_vertex, editor, current_cell);
>
>           // set markers of all cell neighbors of longest edge to false
> -         //if ( cn->index() != c->index() )
>           cell_marker[cn->index()] = false;
>           // set all edges of cell neighbors to forbidden
>           for (EdgeIterator en(*cn); !en.end(); ++en)
> @@ -245,33 +225,27 @@
>      }
>    }
>
> -  // Overwrite old mesh with refined mesh
> +  // Close editor
>    editor.close();
> -  mesh = refined_mesh;
> +
> +  return refined_mesh;
>  }
>  //-----------------------------------------------------------------------------
> -void LocalMeshRefinement::refineIterativelyByEdgeBisection(Mesh& mesh,
> +dolfin::Mesh LocalMeshRefinement::refineIterativelyByEdgeBisection(const Mesh& mesh,
>                                                   MeshFunction<bool>& cell_marker)
>  {
> +  Mesh refined_mesh(mesh);
> +  MeshFunction<uint> edges(mesh, 1, 0);
>
>    begin("Iterative mesh refinement algorithm");
> -
> -  MeshFunction<uint> edges(mesh, 1);
> -  edges = 0;
> -
> -  int i = 0;
>    bool more_iterations = true;
>    while( more_iterations )
> -  {
> -    cout << "=== Iteration : " << i++ << " ===" << endl;
> -    more_iterations = iterationOfRefinement(mesh, cell_marker, edges);
> -  }
> -
> +    more_iterations = iterationOfRefinement(refined_mesh, cell_marker, edges);
>    end();
> -
> +  return refined_mesh;
>  }
>  //-----------------------------------------------------------------------------
> -void LocalMeshRefinement::refineRecursivelyByEdgeBisection(Mesh& mesh,
> +dolfin::Mesh LocalMeshRefinement::refineRecursivelyByEdgeBisection(const Mesh& mesh,
>                                                   MeshFunction<bool>& cell_marker)
>  {
>    // Transformation maps
> @@ -279,22 +253,16 @@
>    std::vector<int> facet_map;
>
>    // Create new mesh
> -  Mesh newmesh;
> -  newmesh = mesh;
> -
> -  RivaraRefinement::refine(newmesh, cell_marker, cell_map, facet_map);
> -
> -  transformMeshData(newmesh, mesh, cell_map, facet_map);
> -
> -  // Overwrite old mesh with refined mesh
> -  mesh = newmesh;
> +  Mesh refined_mesh = RivaraRefinement::refine(mesh, cell_marker, cell_map, facet_map);
> +  transformMeshData(refined_mesh, mesh, cell_map, facet_map);
> +  return refined_mesh;
>  }
>  //-----------------------------------------------------------------------------
>  void LocalMeshRefinement::bisectEdgeOfSimplexCell(const Cell& cell,
>                                                    Edge& edge,
>                                                    uint new_vertex,
>                                                    MeshEditor& editor,
> -                                                  uint& current_cell)
> +                                                  uint current_cell)
>  {
>    // Init cell vertices
>    std::vector<uint> cell1_vertices(cell.num_entities(0));
> @@ -328,8 +296,7 @@
>  //-----------------------------------------------------------------------------
>  bool LocalMeshRefinement::iterationOfRefinement(Mesh& mesh,
>                                        MeshFunction<bool>& cell_marker,
> -                                      MeshFunction<uint>& bisected_edges
> -                                      )
> +                                      MeshFunction<uint>& bisected_edges)
>  {
>    // Map used for edge function transformation
>    std::map<std::vector<uint>, uint, cmp2> edge_map;
> @@ -355,8 +322,7 @@
>    double lmax, l;
>
>    // Temporary edge marker function
> -  MeshFunction<bool> found_edge(mesh, 1);
> -  found_edge = false;
> +  MeshFunction<bool> found_edge(mesh, 1, false);
>
>    // Calculate number of new vertices
>    for (CellIterator c(mesh); !c.end(); ++c)
> @@ -391,12 +357,12 @@
>    editor.init_vertices(num_vertices + num_new_vertices);
>    editor.init_cells(num_cells + num_new_cells);
>
> -  //Rewrite old vertices
> +  // Rewrite old vertices
>    uint current_vertex = 0;
>    for (VertexIterator v(mesh); !v.end(); ++v)
>      editor.add_vertex(current_vertex++, v->point());
>
> -  //Rewrite old cells
> +  // Rewrite old cells
>    uint current_cell = 0;
>    std::vector<uint> cell_vertices(cell_type.num_entities(0));
>    for (CellIterator c(mesh); !c.end(); ++c)
> @@ -410,7 +376,7 @@
>      }
>    }
>
> -  //Main loop
> +  // Main loop
>    for (CellIterator c(mesh); !c.end(); ++c)
>    {
>      if( cell_marker[c->index()] )
> @@ -484,13 +450,12 @@
>      ev[0]= e->entities(0)[0];
>      ev[1]= e->entities(0)[1];
>      if ( (edge_map_it = edge_map.find(ev)) == edge_map.end())
> -    {
>        bisected_edges[e->index()] = 0;
> -    }
>      else
>      {
>        bisected_edges[e->index()] = edge_map_it->second;
> -      for(CellIterator c(*e); !c.end(); ++c) {
> +      for(CellIterator c(*e); !c.end(); ++c)
> +      {
>          cell_marker[c->index()] = true;
>          next_iteration = true;
>        }
> @@ -504,11 +469,10 @@
>    return next_iteration;
>  }
>  //-----------------------------------------------------------------------------
> -void LocalMeshRefinement::transformMeshData(Mesh& newmesh, Mesh& oldmesh,
> +void LocalMeshRefinement::transformMeshData(Mesh& newmesh, const Mesh& oldmesh,
>                                              MeshFunction<uint>& cell_map,
>                                              std::vector<int>& facet_map)
>  {
> -
>    newmesh.data().clear();
>
>    // Rewrite materials
> @@ -521,7 +485,7 @@
>      info("MeshData MeshFunction \"material indicators\" transformed.");
>    }
>
> -  //Rewrite boundary indicators
> +  // Rewrite boundary indicators
>    if( oldmesh.data().array("boundary facet cells")
>        && oldmesh.data().array("boundary facet numbers")
>        && oldmesh.data().array("boundary indicators") )
> @@ -533,20 +497,16 @@
>      std::vector<dolfin::uint>* bi ;
>      bfc = oldmesh.data().array("boundary facet cells");
>      bfn = oldmesh.data().array("boundary facet numbers");
> -    bi = oldmesh.data().array("boundary indicators");
> +    bi  = oldmesh.data().array("boundary indicators");
>      dolfin::uint bi_table_size = oldmesh.num_cells()*num_ent;
>      std::vector<int> bi_table;
>      bi_table.resize(bi_table_size);
>      for(dolfin::uint i=0; i< bi_table_size; i++)
> -    {
>        bi_table[i] = -1;
> -    }
>      for(dolfin::uint i=0; i< bi->size(); i++)
> -    {
>        bi_table[ (*bfc)[i]*num_ent+(*bfn)[i] ] = (*bi)[i];
> -    }
>
> -    //Empty loop to count elements
> +    // Empty loop to count elements
>      dolfin::uint bi_size = 0;
>      for(dolfin::uint c=0; c< newmesh.num_cells(); c++)
>      {
> @@ -556,9 +516,7 @@
>          {
>            dolfin::uint table_map = cell_map[c]*num_ent + facet_map[c*num_ent+f];
>            if( bi_table[ table_map ] != -1 )
> -          {
>              bi_size++;
> -          }
>          }
>        }
>      }
> @@ -569,7 +527,7 @@
>      std::vector<dolfin::uint>* bi_new ;
>      bfc_new = newmesh.data().create_array("boundary facet cells", bi_size);
>      bfn_new = newmesh.data().create_array("boundary facet numbers", bi_size);
> -    bi_new = newmesh.data().create_array("boundary indicators", bi_size);
> +    bi_new  = newmesh.data().create_array("boundary indicators", bi_size);
>
>      // Main transformation loop
>      dolfin::uint number_bi = 0;
> @@ -590,9 +548,7 @@
>          }
>        }
>      }
> -
> -  info("MeshData \"boundary indicators\" transformed.");
> +    info("MeshData \"boundary indicators\" transformed.");
>    }
> -
>  }
> -
> +//-----------------------------------------------------------------------------
>
> === modified file 'dolfin/mesh/LocalMeshRefinement.h'
> --- dolfin/mesh/LocalMeshRefinement.h	2009-04-27 13:01:54 +0000
> +++ dolfin/mesh/LocalMeshRefinement.h	2010-02-07 21:53:42 +0000
> @@ -11,6 +11,7 @@
>  namespace dolfin
>  {
>
> +  // Forward declarations
>    class Cell;
>    class Edge;
>    class Mesh;
> @@ -23,19 +24,19 @@
>    public:
>
>      /// Refine simplicial mesh locally by edge bisection
> -    static void refineMeshByEdgeBisection(Mesh& mesh,
> +    static Mesh refineMeshByEdgeBisection(const Mesh& mesh,
>                                            MeshFunction<bool>& cell_marker,
>                                            bool refine_boundary = true);
>
>      /// Iteratively refine mesh locally by the longest edge bisection
> -    static void refineIterativelyByEdgeBisection(Mesh& mesh,
> -                                                 MeshFunction<bool>& cell_marker);
> +    static Mesh refineIterativelyByEdgeBisection(const Mesh& mesh,
> +                                              MeshFunction<bool>& cell_marker);
>
>      /// Recursively refine mesh locally by the longest edge bisection
>      /// Fast Rivara algorithm implementation with propagation MeshFunctions and
> -    ///   arrays for boundary indicators
> -    static void refineRecursivelyByEdgeBisection(Mesh& mesh,
> -                                                 MeshFunction<bool>& cell_marker);
> +    /// arrays for boundary indicators
> +    static Mesh refineRecursivelyByEdgeBisection(const Mesh& mesh,
> +                                              MeshFunction<bool>& cell_marker);
>
>    private:
>
> @@ -43,7 +44,7 @@
>      static void bisectEdgeOfSimplexCell(const Cell& cell, Edge& edge,
>                                          uint new_vertex,
>                                          MeshEditor& editor,
> -                                        uint& current_cell);
> +                                        uint current_cell);
>
>      /// Iteration of iterative algorithm
>      static bool iterationOfRefinement(Mesh& mesh,
> @@ -51,9 +52,9 @@
>                                        MeshFunction<uint>& bisected_edges);
>
>      /// Transform MeshData
> -    static void transformMeshData(Mesh& newmesh, Mesh& oldmesh,
> +    static void transformMeshData(Mesh& newmesh, const Mesh& oldmesh,
>                                    MeshFunction<uint>& cell_map,
> -		                  std::vector<int>& facet_map);
> +		                              std::vector<int>& facet_map);
>
>    };
>
>
> === modified file 'dolfin/mesh/RivaraRefinement.cpp'
> --- dolfin/mesh/RivaraRefinement.cpp	2009-10-08 19:11:44 +0000
> +++ dolfin/mesh/RivaraRefinement.cpp	2010-02-07 21:53:42 +0000
> @@ -2,27 +2,29 @@
>  // Licensed under the GNU LGPL Version 2.1.
>  //
>  // Modified by Bartosz Sawicki, 2009.
> +// Modified by Garth N. Wells, 2010.
> +//
> +// First added:  2008
> +// Last changed: 2010-02-07
>
>  #include <dolfin/log/dolfin_log.h>
>  #include <dolfin/mesh/Mesh.h>
>  #include <dolfin/mesh/MeshEditor.h>
> -#include <dolfin/mesh/MeshFunction.h>
>  #include <dolfin/mesh/Vertex.h>
>  #include <dolfin/mesh/Cell.h>
>  #include "RivaraRefinement.h"
>
> -
>  using namespace dolfin;
>
>  //-----------------------------------------------------------------------------
> -void RivaraRefinement::refine(Mesh& mesh,
> -			      MeshFunction<bool>& cell_marker,
> -			      MeshFunction<uint>& cell_map,
> -			      std::vector<int>& facet_map)
> +Mesh RivaraRefinement::refine(const Mesh& mesh,
> +			                        MeshFunction<bool>& cell_marker,
> +			                        MeshFunction<uint>& cell_map,
> +			                        std::vector<int>& facet_map)
>  {
>    info("Refining simplicial mesh by recursive Rivara bisection.");
>
> -  int dim = mesh.topology().dim();
> +  const uint dim = mesh.topology().dim();
>
>    // Create dynamic mesh and import data
>    DMesh dmesh;
> @@ -32,7 +34,7 @@
>    std::vector<bool> dmarked(mesh.num_cells());
>    for (CellIterator ci(mesh); !ci.end(); ++ci)
>    {
> -    if(cell_marker[ci->index()] == true)
> +    if(cell_marker[*ci] == true)
>        dmarked[ci->index()] = true;
>      else
>        dmarked[ci->index()] = false;
> @@ -58,11 +60,9 @@
>    std::vector<int> new2old_cell_arr;
>    std::vector<int> new2old_facet_arr;
>
> -  Mesh omesh;
> -
> -  dmesh.export_mesh(omesh, new2old_cell_arr, new2old_facet_arr);
> -
> -  mesh = omesh;
> +  // Refine mesh
> +  Mesh refined_mesh;
> +  dmesh.export_mesh(refined_mesh, new2old_cell_arr, new2old_facet_arr);
>
>    // Generate cell mesh function map
>    cell_map.init(mesh, dim);
> @@ -74,6 +74,8 @@
>    facet_map = new_facet_map;
>    for (uint i=0; i<new2old_facet_arr.size(); i++ )
>      facet_map[i] = new2old_facet_arr[i];
> +
> +  return refined_mesh;
>  }
>  //-----------------------------------------------------------------------------
>  RivaraRefinement::DVertex::DVertex() : id(0), cells(0), p(0.0, 0.0, 0.0)
> @@ -81,7 +83,8 @@
>    // Do nothing
>  }
>  //-----------------------------------------------------------------------------
> -RivaraRefinement::DCell::DCell() : id(0), parent_id(0), vertices(0), deleted(false), facets(0)
> +RivaraRefinement::DCell::DCell() : id(0), parent_id(0), vertices(0),
> +                                   deleted(false), facets(0)
>  {
>    // Do nothing
>  }
> @@ -96,15 +99,18 @@
>    // Delete allocated DVertices
>    for(std::list<DVertex* >::iterator it = vertices.begin();
>        it != vertices.end(); ++it)
> +  {
>      delete *it;
> -
> +  }
>    // Delete allocated DCells
>    for(std::list<DCell* >::iterator it = cells.begin();
>        it != cells.end(); ++it)
> +  {
>      delete *it;
> +  }
>  }
>  //-----------------------------------------------------------------------------
> -void RivaraRefinement::DMesh::import_mesh(Mesh& mesh)
> +void RivaraRefinement::DMesh::import_mesh(const Mesh& mesh)
>  {
>    cell_type = &(mesh.type());
>    dim = mesh.topology().dim();
> @@ -195,7 +201,6 @@
>      current_cell++;
>    }
>    editor.close();
> -
>  }
>  //-----------------------------------------------------------------------------
>  void RivaraRefinement::DMesh::number()
> @@ -339,7 +344,7 @@
>
>      if(c != dcell)
>      {
> -      int matches = 0;
> +      uint matches = 0;
>        for(uint i = 0; i < c->vertices.size(); i++)
>        {
>          if(c->vertices[i] == v1 || c->vertices[i] == v2)
> @@ -469,6 +474,3 @@
>    }
>  }
>  //-----------------------------------------------------------------------------
> -
> -
> -
>
> === modified file 'dolfin/mesh/RivaraRefinement.h'
> --- dolfin/mesh/RivaraRefinement.h	2009-04-27 13:01:54 +0000
> +++ dolfin/mesh/RivaraRefinement.h	2010-02-07 21:53:42 +0000
> @@ -2,6 +2,10 @@
>  // Licensed under the GNU LGPL Version 2.1.
>  //
>  // Modified by Bartosz Sawicki, 2009.
> +// Modified by Garth N. Wells, 2010.
> +//
> +// First added:  2008
> +// Last changed: 2010-02-07
>
>  #ifndef __RIVARAREFINEMENT_H
>  #define __RIVARAREFINEMENT_H
> @@ -15,6 +19,9 @@
>
>  namespace dolfin
>  {
> +  // Forward declarations
> +  //class CellType;
> +  //class Mesh
>
>    class RivaraRefinement
>    {
> @@ -22,10 +29,9 @@
>    public:
>
>      /// Refine simplicial mesh locally by recursive edge bisection
> -    static void refine(Mesh& mesh,
> -		       MeshFunction<bool>& cell_marker,
> -		       MeshFunction<uint>& cell_map,
> -		       std::vector<int>& facet_map);
> +    static Mesh refine(const Mesh& mesh, MeshFunction<bool>& cell_marker,
> +		                   MeshFunction<uint>& cell_map,
> +                       std::vector<int>& facet_map);
>
>    private:
>
> @@ -62,26 +68,22 @@
>        void add_vertex(DVertex* v);
>        void add_cell(DCell* c, std::vector<DVertex*> vs, int parent_id);
>        void remove_cell(DCell* c);
> -      void import_mesh(Mesh& mesh);
> +      void import_mesh(const Mesh& mesh);
>        void export_mesh(Mesh& mesh, std::vector<int>& new2old_cell, std::vector<int>& new2old_facet);
>        void number();
>        void bisect(DCell* dcell, DVertex* hangv, DVertex* hv0, DVertex* hv1);
>        void bisect_marked(std::vector<bool> marked_ids);
> -                        DCell* opposite(DCell* dcell, DVertex* v1, DVertex* v2);
> +      DCell* opposite(DCell* dcell, DVertex* v1, DVertex* v2);
>        void propagate_facets(DCell* dcell, DCell* c0, DCell* c1, uint ii, uint jj, DVertex* mv);
>
>        std::list<DVertex *> vertices;
>        std::list<DCell *> cells;
> -      CellType* cell_type;
> +      const CellType* cell_type;
>        uint dim;
>      };
>
> -
>    };
>
> -
> -
> -
>  }
>
>  #endif
>
> === modified file 'dolfin/mesh/UniformMeshRefinement.cpp'
> --- dolfin/mesh/UniformMeshRefinement.cpp	2009-05-02 20:55:48 +0000
> +++ dolfin/mesh/UniformMeshRefinement.cpp	2010-02-07 21:53:42 +0000
> @@ -1,8 +1,10 @@
>  // Copyright (C) 2006-2007 Anders Logg.
>  // Licensed under the GNU LGPL Version 2.1.
>  //
> +// Modified by Garth N. Wells, 2010
> +//
>  // First added:  2006-06-08
> -// Last changed: 2007-05-24
> +// Last changed: 2010-02-07
>
>  #include <dolfin/math/dolfin_math.h>
>  #include <dolfin/log/dolfin_log.h>
> @@ -19,13 +21,13 @@
>  using namespace dolfin;
>
>  //-----------------------------------------------------------------------------
> -void UniformMeshRefinement::refine(Mesh& mesh)
> +dolfin::Mesh UniformMeshRefinement::refine(const Mesh& mesh)
>  {
>    // Only know how to refine simplicial meshes
> -  refine_simplex(mesh);
> +  return refine_simplex(mesh);
>  }
>  //-----------------------------------------------------------------------------
> -void UniformMeshRefinement::refine_simplex(Mesh& mesh)
> +dolfin::Mesh UniformMeshRefinement::refine_simplex(const Mesh& mesh)
>  {
>    info(1, "Refining simplicial mesh uniformly.");
>
> @@ -42,7 +44,7 @@
>    Mesh refined_mesh;
>    MeshEditor editor;
>    editor.open(refined_mesh, cell_type.cell_type(),
> -	      mesh.topology().dim(), mesh.geometry().dim());
> +	            mesh.topology().dim(), mesh.geometry().dim());
>
>    // Get size of mesh
>    const uint num_vertices = mesh.size(0);
> @@ -67,8 +69,9 @@
>    for (CellIterator c(mesh); !c.end(); ++c)
>      cell_type.refine_cell(*c, editor, current_cell);
>
> -  // Overwrite old mesh with refined mesh
> +  // Close editor
>    editor.close();
> -  mesh = refined_mesh;
> +
> +  return refined_mesh;
>  }
>  //-----------------------------------------------------------------------------
>
> === modified file 'dolfin/mesh/UniformMeshRefinement.h'
> --- dolfin/mesh/UniformMeshRefinement.h	2009-04-27 07:46:35 +0000
> +++ dolfin/mesh/UniformMeshRefinement.h	2010-02-07 21:53:42 +0000
> @@ -1,8 +1,10 @@
>  // Copyright (C) 2006 Anders Logg.
>  // Licensed under the GNU LGPL Version 2.1.
>  //
> +// Modified by Garth N. Wells, 2010
> +//
>  // First added:  2006-06-07
> -// Last changed: 2006-06-16
> +// Last changed: 2010-02-07
>
>  #ifndef __UNIFORM_MESH_REFINEMENT_H
>  #define __UNIFORM_MESH_REFINEMENT_H
> @@ -21,10 +23,12 @@
>    public:
>
>      /// Refine mesh uniformly according to mesh type
> -    static void refine(Mesh& mesh);
> +    static Mesh refine(const Mesh& mesh);
> +
> +  private:
>
>      /// Refine simplicial mesh uniformly
> -    static void refine_simplex(Mesh& mesh);
> +    static Mesh refine_simplex(const Mesh& mesh);
>
>    };
>
>

Attachment: signature.asc
Description: Digital signature