Mercurial > fem-fenics-eugenio
changeset 86:66e4aa87c9a1
Rename the function accordingly to dolfin
* fem_init_mesh now becomes Mesh
* Makefile.in compile function with new name
| author | gedeone-octave <marcovass89@hotmail.it> |
|---|---|
| date | Mon, 05 Aug 2013 11:39:47 +0200 |
| parents | 8084ecfaa2b7 |
| children | d03627091414 |
| files | src/Makefile.in src/Mesh.cc src/fem_init_mesh.cc |
| diffstat | 3 files changed, 345 insertions(+), 346 deletions(-) [+] |
line wrap: on
line diff
--- a/src/Makefile.in Sat Aug 03 14:49:56 2013 +0200 +++ b/src/Makefile.in Mon Aug 05 11:39:47 2013 +0200 @@ -1,7 +1,7 @@ MKOCTFILE ?= mkoctfile OCTFILES=fem_init_env.oct\ - fem_init_mesh.oct \ + Mesh.oct \ fem_get_mesh.oct \ DirichletBC.oct \ Expression.oct \ @@ -19,11 +19,11 @@ fem_init_env.o: fem_init_env.cc mesh.h functionspace.h boundarycondition.h function.h coefficient.h $(MKOCTFILE) $(CPPFLAGS) -c fem_init_env.cc $(LDFLAGS) -o $@ -I. -fem_init_mesh.oct: fem_init_mesh.o - $(MKOCTFILE) $(CPPFLAGS) -s fem_init_mesh.o -o $@ $(LDFLAGS) $(LIBS) +Mesh.oct: Mesh.o + $(MKOCTFILE) $(CPPFLAGS) -s Mesh.o -o $@ $(LDFLAGS) $(LIBS) -fem_init_mesh.o: fem_init_mesh.cc mesh.h - $(MKOCTFILE) $(CPPFLAGS) -c fem_init_mesh.cc $(LDFLAGS) -o $@ -I. +Mesh.o: Mesh.cc mesh.h + $(MKOCTFILE) $(CPPFLAGS) -c Mesh.cc $(LDFLAGS) -o $@ -I. fem_get_mesh.oct: fem_get_mesh.o $(MKOCTFILE) $(CPPFLAGS) -s fem_get_mesh.o -o $@ $(LDFLAGS) $(LIBS)
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/src/Mesh.cc Mon Aug 05 11:39:47 2013 +0200 @@ -0,0 +1,340 @@ +/* + Copyright (C) 2013 Marco Vassallo + + This program is free software; you can redistribute it and/or modify it under + the terms of the GNU General Public License as published by the Free Software + Foundation; either version 2 of the License, or (at your option) any later + version. + + This program is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more + details. + + You should have received a copy of the GNU General Public License along with + this program; if not, see <http://www.gnu.org/licenses/>. +*/ + +#include "mesh.h" + +DEFUN_DLD (Mesh, args, ,"-*- texinfo -*-\n\ +@deftypefn {Function File} {[@var{mesh_out}]} = \ +Mesh (@var{mesh_in}) \n\ +The @var{mesh_in} should be either\n\ +@itemize @bullet \n\ +@item a string containing the name of the file where the mesh is stored\ +(compatible formats are the ones which are compatible with Fenics\n\ +@item a a PDE-tool like structure with matrix fields (p,e,t)\n\ +@end itemize\n\ +The output @var{mesh_out} is a representation of the\n\ +@var{mesh_in} which is compatible with fem-fenics\n\ +@seealso{fem_get_mesh}\n\ +@end deftypefn") +{ + int nargin = args.length (); + octave_value retval = 0; + + if (nargin < 1 || nargin > 1) + print_usage (); + else + { + if (!error_state) + { + if (args(0).is_string () == true) + { + std::string filename = args(0).string_value (); + //if the filename is not valid, dolfin takes care of it + + if (! mesh_type_loaded) + { + mesh::register_type (); + mesh_type_loaded = true; + mlock (); + } + retval = new mesh (filename); + } + + else + { + octave_scalar_map a = args(0).scalar_map_value (); + Array<double> p = a.contents ("p").matrix_value (); + Array<octave_idx_type> t = a.contents ("t").matrix_value (); + Array<octave_idx_type> e = a.contents ("e").matrix_value (); + + if (! error_state) + { + if (! mesh_type_loaded) + { + mesh::register_type (); + mesh_type_loaded = true; + mlock (); + } + + retval = new mesh (p, e, t); + } + else + error ("fem_init_mesh: the argument you provide is invalid"); + } + } + } + return retval; +} + +mesh::mesh (Array<double>& p, Array<octave_idx_type>& e, + Array<octave_idx_type>& t) +{ + + std::size_t D = p.rows (); + if (D < 2 || D > 3) + error ("mesh constructor: only 2D or 3D meshes are supported"); + else + { + dolfin::MeshEditor editor; + dolfin::Mesh msh; + editor.open (msh, D, D); + editor.init_vertices (p.cols ()); + editor.init_cells (t.cols ()); + + if (D == 2) + { + for (uint i = 0; i < p.cols (); ++i) + editor.add_vertex (i, p.xelem (0, i), p.xelem (1, i)); + + for (uint i = 0; i < t.cols (); ++i) + editor.add_cell (i, t.xelem (0, i) - 1, + t.xelem (1, i) - 1, t.xelem (2, i) - 1); + } + + if (D == 3) + { + for (uint i = 0; i < p.cols (); ++i) + editor.add_vertex (i, p.xelem (0, i), + p.xelem (1, i), p.xelem (2, i)); + + for (uint i = 0; i < t.cols (); ++i) + editor.add_cell (i, t.xelem (0, i) - 1, t.xelem (1, i) - 1, + t.xelem (2, i) - 1, t.xelem (3, i) - 1); + } + + editor.close (); + + // store information associated with e + msh.init (D - 1); + dolfin::MeshValueCollection<std::size_t> facet (D - 1); + std::size_t num_side_edges = e.cols (); + + if (D == 2) + { + for (uint i = 0; i < num_side_edges; ++i) + { + dolfin::Vertex v (msh, e.xelem (0, i) - 1); + for (dolfin::FacetIterator f (v); ! f.end (); ++f) + { + if ((*f).entities(0)[0] == e.xelem (0, i) - 1 + && (*f).entities(0)[1] == e.xelem (1, i) - 1 + || (*f).entities(0)[0] == e.xelem (1, i) - 1 + && (*f).entities(0)[1] == e.xelem (0, i) - 1) + { + facet.set_value ((*f).index (), e.xelem (4, i), msh); + break; + } + } + } + } + + + + if (D == 3) + { + for (uint i = 0; i < num_side_edges; ++i) + { + dolfin::Vertex v (msh, e.xelem (0, i) - 1); + for (dolfin::FacetIterator f (v); ! f.end (); ++f) + { + if ((*f).entities(0)[0] == e(0, i) - 1 + && (*f).entities(0)[1] == e.xelem (1, i) - 1 + && (*f).entities(0)[2] == e.xelem (2, i) - 1 + || (*f).entities(0)[0] == e.xelem (0, i) - 1 + && (*f).entities(0)[1] == e.xelem (2, i) - 1 + && (*f).entities(0)[2] == e.xelem (1, i) - 1 + || (*f).entities(0)[0] == e.xelem (1, i) - 1 + && (*f).entities(0)[1] == e.xelem (0, i) - 1 + && (*f).entities(0)[2] == e.xelem (2, i) - 1 + || (*f).entities(0)[0] == e.xelem (1, i) - 1 + && (*f).entities(0)[1] == e.xelem (2, i) - 1 + && (*f).entities(0)[2] == e.xelem (0, i) - 1 + || (*f).entities(0)[0] == e.xelem (2, i) - 1 + && (*f).entities(0)[1] == e.xelem (0, i) - 1 + && (*f).entities(0)[2] == e.xelem (1, i) - 1 + || (*f).entities(0)[0] == e.xelem (2, i) - 1 + && (*f).entities(0)[1] == e.xelem (1, i) - 1 + && (*f).entities(0)[2] == e.xelem (0, i) - 1) + { + facet.set_value ((*f).index (), e.xelem (9, i), msh); + break; + } + } + } + } + + *(msh.domains ().markers (D - 1)) = facet; + + // store information associated with t + dolfin::MeshValueCollection<std::size_t> cell (D); + std::size_t num_cells = t.cols (); + + if (D == 2) + { + for (uint i = 0; i < num_cells; ++i) + { + dolfin::Vertex v (msh, t.xelem (0, i) - 1); + for (dolfin::CellIterator f (v); ! f.end (); ++f) + { + if ((*f).entities(0)[0] == t.xelem (0, i) - 1 + && (*f).entities(0)[1] == t.xelem (1, i) - 1 + && (*f).entities(0)[2] == t.xelem (2, i) - 1 + || (*f).entities(0)[0] == t.xelem (0, i) - 1 + && (*f).entities(0)[1] == t.xelem (2, i) - 1 + && (*f).entities(0)[2] == t.xelem (1, i) - 1 + || (*f).entities(0)[0] == t.xelem (1, i) - 1 + && (*f).entities(0)[1] == t.xelem (0, i) - 1 + && (*f).entities(0)[2] == t.xelem (2, i) - 1 + || (*f).entities(0)[0] == t.xelem (1, i) - 1 + && (*f).entities(0)[1] == t.xelem (2, i) - 1 + && (*f).entities(0)[2] == t.xelem (0, i) - 1 + || (*f).entities(0)[0] == t.xelem (2, i) - 1 + && (*f).entities(0)[1] == t.xelem (0, i) - 1 + && (*f).entities(0)[2] == t.xelem (1, i) - 1 + || (*f).entities(0)[0] == t.xelem (2, i) - 1 + && (*f).entities(0)[1] == t.xelem (1, i) - 1 + && (*f).entities(0)[2] == t.xelem (0, i) - 1) + { + cell.set_value ((*f).index (), t.xelem (3, i), msh); + break; + } + } + } + } + + + + if (D == 3) + { + for (uint i = 0; i < num_cells; ++i) + { + dolfin::Vertex v (msh, t.xelem (0, i) - 1); + for (dolfin::CellIterator f (v); ! f.end (); ++f) + { + if ((*f).entities(0)[0] == t.xelem (0, i) - 1 + && (*f).entities(0)[1] == t.xelem (1, i) - 1 + && (*f).entities(0)[2] == t.xelem (2, i) - 1 + && (*f).entities(0)[3] == t.xelem (3, i) - 1 + || (*f).entities(0)[0] == t.xelem (0, i) - 1 + && (*f).entities(0)[1] == t.xelem (1, i) - 1 + && (*f).entities(0)[2] == t.xelem (3, i) - 1 + && (*f).entities(0)[3] == t.xelem (2, i) - 1 + || (*f).entities(0)[0] == t.xelem (0, i) - 1 + && (*f).entities(0)[1] == t.xelem (2, i) - 1 + && (*f).entities(0)[2] == t.xelem (1, i) - 1 + && (*f).entities(0)[3] == t.xelem (3, i) - 1 + || (*f).entities(0)[0] == t.xelem (0, i) - 1 + && (*f).entities(0)[1] == t.xelem (2, i) - 1 + && (*f).entities(0)[2] == t.xelem (3, i) - 1 + && (*f).entities(0)[3] == t.xelem (1, i) - 1 + || (*f).entities(0)[0] == t.xelem (0, i) - 1 + && (*f).entities(0)[1] == t.xelem (3, i) - 1 + && (*f).entities(0)[2] == t.xelem (1, i) - 1 + && (*f).entities(0)[3] == t.xelem (2, i) - 1 + || (*f).entities(0)[0] == t.xelem (0, i) - 1 + && (*f).entities(0)[1] == t.xelem (3, i) - 1 + && (*f).entities(0)[2] == t.xelem (2, i) - 1 + && (*f).entities(0)[3] == t.xelem (1, i) - 1 + + || (*f).entities(0)[0] == t.xelem (1, i) - 1 + && (*f).entities(0)[1] == t.xelem (0, i) - 1 + && (*f).entities(0)[2] == t.xelem (2, i) - 1 + && (*f).entities(0)[3] == t.xelem (3, i) - 1 + || (*f).entities(0)[0] == t.xelem (1, i) - 1 + && (*f).entities(0)[1] == t.xelem (0, i) - 1 + && (*f).entities(0)[2] == t.xelem (3, i) - 1 + && (*f).entities(0)[3] == t.xelem (2, i) - 1 + || (*f).entities(0)[0] == t.xelem (1, i) - 1 + && (*f).entities(0)[1] == t.xelem (2, i) - 1 + && (*f).entities(0)[2] == t.xelem (0, i) - 1 + && (*f).entities(0)[3] == t.xelem (3, i) - 1 + || (*f).entities(0)[0] == t.xelem (1, i) - 1 + && (*f).entities(0)[1] == t.xelem (2, i) - 1 + && (*f).entities(0)[2] == t.xelem (3, i) - 1 + && (*f).entities(0)[3] == t.xelem (0, i) - 1 + || (*f).entities(0)[0] == t.xelem (1, i) - 1 + && (*f).entities(0)[1] == t.xelem (3, i) - 1 + && (*f).entities(0)[2] == t.xelem (0, i) - 1 + && (*f).entities(0)[3] == t.xelem (2, i) - 1 + || (*f).entities(0)[0] == t.xelem (1, i) - 1 + && (*f).entities(0)[1] == t.xelem (3, i) - 1 + && (*f).entities(0)[2] == t.xelem (2, i) - 1 + && (*f).entities(0)[3] == t.xelem (0, i) - 1 + + || (*f).entities(0)[0] == t.xelem (2, i) - 1 + && (*f).entities(0)[1] == t.xelem (0, i) - 1 + && (*f).entities(0)[2] == t.xelem (1, i) - 1 + && (*f).entities(0)[3] == t.xelem (3, i) - 1 + || (*f).entities(0)[0] == t.xelem (2, i) - 1 + && (*f).entities(0)[1] == t.xelem (0, i) - 1 + && (*f).entities(0)[2] == t.xelem (3, i) - 1 + && (*f).entities(0)[3] == t.xelem (1, i) - 1 + || (*f).entities(0)[0] == t.xelem (2, i) - 1 + && (*f).entities(0)[1] == t.xelem (1, i) - 1 + && (*f).entities(0)[2] == t.xelem (0, i) - 1 + && (*f).entities(0)[3] == t.xelem (3, i) - 1 + || (*f).entities(0)[0] == t.xelem (2, i) - 1 + && (*f).entities(0)[1] == t.xelem (1, i) - 1 + && (*f).entities(0)[2] == t.xelem (3, i) - 1 + && (*f).entities(0)[3] == t.xelem (0, i) - 1 + || (*f).entities(0)[0] == t.xelem (2, i) - 1 + && (*f).entities(0)[1] == t.xelem (3, i) - 1 + && (*f).entities(0)[2] == t.xelem (0, i) - 1 + && (*f).entities(0)[3] == t.xelem (1, i) - 1 + || (*f).entities(0)[0] == t.xelem (2, i) - 1 + && (*f).entities(0)[1] == t.xelem (3, i) - 1 + && (*f).entities(0)[2] == t.xelem (1, i) - 1 + && (*f).entities(0)[3] == t.xelem (0, i) - 1 + + || (*f).entities(0)[0] == t.xelem (3, i) - 1 + && (*f).entities(0)[1] == t.xelem (0, i) - 1 + && (*f).entities(0)[2] == t.xelem (1, i) - 1 + && (*f).entities(0)[3] == t.xelem (2, i) - 1 + || (*f).entities(0)[0] == t.xelem (3, i) - 1 + && (*f).entities(0)[1] == t.xelem (0, i) - 1 + && (*f).entities(0)[2] == t.xelem (2, i) - 1 + && (*f).entities(0)[3] == t.xelem (1, i) - 1 + || (*f).entities(0)[0] == t.xelem (3, i) - 1 + && (*f).entities(0)[1] == t.xelem (1, i) - 1 + && (*f).entities(0)[2] == t.xelem (0, i) - 1 + && (*f).entities(0)[3] == t.xelem (2, i) - 1 + || (*f).entities(0)[0] == t.xelem (3, i) - 1 + && (*f).entities(0)[1] == t.xelem (1, i) - 1 + && (*f).entities(0)[2] == t.xelem (2, i) - 1 + && (*f).entities(0)[3] == t.xelem (0, i) - 1 + || (*f).entities(0)[0] == t.xelem (3, i) - 1 + && (*f).entities(0)[1] == t.xelem (2, i) - 1 + && (*f).entities(0)[2] == t.xelem (0, i) - 1 + && (*f).entities(0)[3] == t.xelem (1, i) - 1 + || (*f).entities(0)[0] == t.xelem (3, i) - 1 + && (*f).entities(0)[1] == t.xelem (2, i) - 1 + && (*f).entities(0)[2] == t.xelem (1, i) - 1 + && (*f).entities(0)[3] == t.xelem (0, i) - 1) + { + cell.set_value ((*f).index (), t.xelem (4, i), msh); + break; + } + } + } + } + + *(msh.domains ().markers (D)) = cell; + + boost::shared_ptr<const dolfin::Mesh> my_m (new dolfin::Mesh (msh)); + pmsh = my_m; + } +}
--- a/src/fem_init_mesh.cc Sat Aug 03 14:49:56 2013 +0200 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,341 +0,0 @@ -/* - Copyright (C) 2013 Marco Vassallo - - This program is free software; you can redistribute it and/or modify it under - the terms of the GNU General Public License as published by the Free Software - Foundation; either version 2 of the License, or (at your option) any later - version. - - This program is distributed in the hope that it will be useful, but WITHOUT - ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or - FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more - details. - - You should have received a copy of the GNU General Public License along with - this program; if not, see <http://www.gnu.org/licenses/>. -*/ - -#include "mesh.h" - -DEFUN_DLD (fem_init_mesh, args, ,"-*- texinfo -*-\n\ -@deftypefn {Function File} {[@var{mesh_out}]} = \ -fem_get_mesh (@var{mesh_in}) \n\ -The @var{mesh_in} should be either\n\ -@itemize @bullet \n\ -@item a string containing the name of the file where the mesh is stored\ -(compatible formats are the ones which are compatible with Fenics\n\ -@item a a PDE-tool like structure with matrix fields (p,e,t)\n\ -@end itemize\n\ -fem_init_mesh().\n\ -The output @var{mesh_out} is a representation of the\n\ -@var{mesh_in} which is compatible with fem-fenics\n\ -@seealso{fem_get_mesh}\n\ -@end deftypefn") -{ - int nargin = args.length (); - octave_value retval = 0; - - if (nargin < 1 || nargin > 1) - print_usage (); - else - { - if (!error_state) - { - if (args(0).is_string () == true) - { - std::string filename = args(0).string_value (); - //if the filename is not valid, dolfin takes care of it - - if (! mesh_type_loaded) - { - mesh::register_type (); - mesh_type_loaded = true; - mlock (); - } - retval = new mesh (filename); - } - - else - { - octave_scalar_map a = args(0).scalar_map_value (); - Array<double> p = a.contents ("p").matrix_value (); - Array<octave_idx_type> t = a.contents ("t").matrix_value (); - Array<octave_idx_type> e = a.contents ("e").matrix_value (); - - if (! error_state) - { - if (! mesh_type_loaded) - { - mesh::register_type (); - mesh_type_loaded = true; - mlock (); - } - - retval = new mesh (p, e, t); - } - else - error ("fem_init_mesh: the argument you provide is invalid"); - } - } - } - return retval; -} - -mesh::mesh (Array<double>& p, Array<octave_idx_type>& e, - Array<octave_idx_type>& t) -{ - - std::size_t D = p.rows (); - if (D < 2 || D > 3) - error ("mesh constructor: only 2D or 3D meshes are supported"); - else - { - dolfin::MeshEditor editor; - dolfin::Mesh msh; - editor.open (msh, D, D); - editor.init_vertices (p.cols ()); - editor.init_cells (t.cols ()); - - if (D == 2) - { - for (uint i = 0; i < p.cols (); ++i) - editor.add_vertex (i, p.xelem (0, i), p.xelem (1, i)); - - for (uint i = 0; i < t.cols (); ++i) - editor.add_cell (i, t.xelem (0, i) - 1, - t.xelem (1, i) - 1, t.xelem (2, i) - 1); - } - - if (D == 3) - { - for (uint i = 0; i < p.cols (); ++i) - editor.add_vertex (i, p.xelem (0, i), - p.xelem (1, i), p.xelem (2, i)); - - for (uint i = 0; i < t.cols (); ++i) - editor.add_cell (i, t.xelem (0, i) - 1, t.xelem (1, i) - 1, - t.xelem (2, i) - 1, t.xelem (3, i) - 1); - } - - editor.close (); - - // store information associated with e - msh.init (D - 1); - dolfin::MeshValueCollection<std::size_t> facet (D - 1); - std::size_t num_side_edges = e.cols (); - - if (D == 2) - { - for (uint i = 0; i < num_side_edges; ++i) - { - dolfin::Vertex v (msh, e.xelem (0, i) - 1); - for (dolfin::FacetIterator f (v); ! f.end (); ++f) - { - if ((*f).entities(0)[0] == e.xelem (0, i) - 1 - && (*f).entities(0)[1] == e.xelem (1, i) - 1 - || (*f).entities(0)[0] == e.xelem (1, i) - 1 - && (*f).entities(0)[1] == e.xelem (0, i) - 1) - { - facet.set_value ((*f).index (), e.xelem (4, i), msh); - break; - } - } - } - } - - - - if (D == 3) - { - for (uint i = 0; i < num_side_edges; ++i) - { - dolfin::Vertex v (msh, e.xelem (0, i) - 1); - for (dolfin::FacetIterator f (v); ! f.end (); ++f) - { - if ((*f).entities(0)[0] == e(0, i) - 1 - && (*f).entities(0)[1] == e.xelem (1, i) - 1 - && (*f).entities(0)[2] == e.xelem (2, i) - 1 - || (*f).entities(0)[0] == e.xelem (0, i) - 1 - && (*f).entities(0)[1] == e.xelem (2, i) - 1 - && (*f).entities(0)[2] == e.xelem (1, i) - 1 - || (*f).entities(0)[0] == e.xelem (1, i) - 1 - && (*f).entities(0)[1] == e.xelem (0, i) - 1 - && (*f).entities(0)[2] == e.xelem (2, i) - 1 - || (*f).entities(0)[0] == e.xelem (1, i) - 1 - && (*f).entities(0)[1] == e.xelem (2, i) - 1 - && (*f).entities(0)[2] == e.xelem (0, i) - 1 - || (*f).entities(0)[0] == e.xelem (2, i) - 1 - && (*f).entities(0)[1] == e.xelem (0, i) - 1 - && (*f).entities(0)[2] == e.xelem (1, i) - 1 - || (*f).entities(0)[0] == e.xelem (2, i) - 1 - && (*f).entities(0)[1] == e.xelem (1, i) - 1 - && (*f).entities(0)[2] == e.xelem (0, i) - 1) - { - facet.set_value ((*f).index (), e.xelem (9, i), msh); - break; - } - } - } - } - - *(msh.domains ().markers (D - 1)) = facet; - - // store information associated with t - dolfin::MeshValueCollection<std::size_t> cell (D); - std::size_t num_cells = t.cols (); - - if (D == 2) - { - for (uint i = 0; i < num_cells; ++i) - { - dolfin::Vertex v (msh, t.xelem (0, i) - 1); - for (dolfin::CellIterator f (v); ! f.end (); ++f) - { - if ((*f).entities(0)[0] == t.xelem (0, i) - 1 - && (*f).entities(0)[1] == t.xelem (1, i) - 1 - && (*f).entities(0)[2] == t.xelem (2, i) - 1 - || (*f).entities(0)[0] == t.xelem (0, i) - 1 - && (*f).entities(0)[1] == t.xelem (2, i) - 1 - && (*f).entities(0)[2] == t.xelem (1, i) - 1 - || (*f).entities(0)[0] == t.xelem (1, i) - 1 - && (*f).entities(0)[1] == t.xelem (0, i) - 1 - && (*f).entities(0)[2] == t.xelem (2, i) - 1 - || (*f).entities(0)[0] == t.xelem (1, i) - 1 - && (*f).entities(0)[1] == t.xelem (2, i) - 1 - && (*f).entities(0)[2] == t.xelem (0, i) - 1 - || (*f).entities(0)[0] == t.xelem (2, i) - 1 - && (*f).entities(0)[1] == t.xelem (0, i) - 1 - && (*f).entities(0)[2] == t.xelem (1, i) - 1 - || (*f).entities(0)[0] == t.xelem (2, i) - 1 - && (*f).entities(0)[1] == t.xelem (1, i) - 1 - && (*f).entities(0)[2] == t.xelem (0, i) - 1) - { - cell.set_value ((*f).index (), t.xelem (3, i), msh); - break; - } - } - } - } - - - - if (D == 3) - { - for (uint i = 0; i < num_cells; ++i) - { - dolfin::Vertex v (msh, t.xelem (0, i) - 1); - for (dolfin::CellIterator f (v); ! f.end (); ++f) - { - if ((*f).entities(0)[0] == t.xelem (0, i) - 1 - && (*f).entities(0)[1] == t.xelem (1, i) - 1 - && (*f).entities(0)[2] == t.xelem (2, i) - 1 - && (*f).entities(0)[3] == t.xelem (3, i) - 1 - || (*f).entities(0)[0] == t.xelem (0, i) - 1 - && (*f).entities(0)[1] == t.xelem (1, i) - 1 - && (*f).entities(0)[2] == t.xelem (3, i) - 1 - && (*f).entities(0)[3] == t.xelem (2, i) - 1 - || (*f).entities(0)[0] == t.xelem (0, i) - 1 - && (*f).entities(0)[1] == t.xelem (2, i) - 1 - && (*f).entities(0)[2] == t.xelem (1, i) - 1 - && (*f).entities(0)[3] == t.xelem (3, i) - 1 - || (*f).entities(0)[0] == t.xelem (0, i) - 1 - && (*f).entities(0)[1] == t.xelem (2, i) - 1 - && (*f).entities(0)[2] == t.xelem (3, i) - 1 - && (*f).entities(0)[3] == t.xelem (1, i) - 1 - || (*f).entities(0)[0] == t.xelem (0, i) - 1 - && (*f).entities(0)[1] == t.xelem (3, i) - 1 - && (*f).entities(0)[2] == t.xelem (1, i) - 1 - && (*f).entities(0)[3] == t.xelem (2, i) - 1 - || (*f).entities(0)[0] == t.xelem (0, i) - 1 - && (*f).entities(0)[1] == t.xelem (3, i) - 1 - && (*f).entities(0)[2] == t.xelem (2, i) - 1 - && (*f).entities(0)[3] == t.xelem (1, i) - 1 - - || (*f).entities(0)[0] == t.xelem (1, i) - 1 - && (*f).entities(0)[1] == t.xelem (0, i) - 1 - && (*f).entities(0)[2] == t.xelem (2, i) - 1 - && (*f).entities(0)[3] == t.xelem (3, i) - 1 - || (*f).entities(0)[0] == t.xelem (1, i) - 1 - && (*f).entities(0)[1] == t.xelem (0, i) - 1 - && (*f).entities(0)[2] == t.xelem (3, i) - 1 - && (*f).entities(0)[3] == t.xelem (2, i) - 1 - || (*f).entities(0)[0] == t.xelem (1, i) - 1 - && (*f).entities(0)[1] == t.xelem (2, i) - 1 - && (*f).entities(0)[2] == t.xelem (0, i) - 1 - && (*f).entities(0)[3] == t.xelem (3, i) - 1 - || (*f).entities(0)[0] == t.xelem (1, i) - 1 - && (*f).entities(0)[1] == t.xelem (2, i) - 1 - && (*f).entities(0)[2] == t.xelem (3, i) - 1 - && (*f).entities(0)[3] == t.xelem (0, i) - 1 - || (*f).entities(0)[0] == t.xelem (1, i) - 1 - && (*f).entities(0)[1] == t.xelem (3, i) - 1 - && (*f).entities(0)[2] == t.xelem (0, i) - 1 - && (*f).entities(0)[3] == t.xelem (2, i) - 1 - || (*f).entities(0)[0] == t.xelem (1, i) - 1 - && (*f).entities(0)[1] == t.xelem (3, i) - 1 - && (*f).entities(0)[2] == t.xelem (2, i) - 1 - && (*f).entities(0)[3] == t.xelem (0, i) - 1 - - || (*f).entities(0)[0] == t.xelem (2, i) - 1 - && (*f).entities(0)[1] == t.xelem (0, i) - 1 - && (*f).entities(0)[2] == t.xelem (1, i) - 1 - && (*f).entities(0)[3] == t.xelem (3, i) - 1 - || (*f).entities(0)[0] == t.xelem (2, i) - 1 - && (*f).entities(0)[1] == t.xelem (0, i) - 1 - && (*f).entities(0)[2] == t.xelem (3, i) - 1 - && (*f).entities(0)[3] == t.xelem (1, i) - 1 - || (*f).entities(0)[0] == t.xelem (2, i) - 1 - && (*f).entities(0)[1] == t.xelem (1, i) - 1 - && (*f).entities(0)[2] == t.xelem (0, i) - 1 - && (*f).entities(0)[3] == t.xelem (3, i) - 1 - || (*f).entities(0)[0] == t.xelem (2, i) - 1 - && (*f).entities(0)[1] == t.xelem (1, i) - 1 - && (*f).entities(0)[2] == t.xelem (3, i) - 1 - && (*f).entities(0)[3] == t.xelem (0, i) - 1 - || (*f).entities(0)[0] == t.xelem (2, i) - 1 - && (*f).entities(0)[1] == t.xelem (3, i) - 1 - && (*f).entities(0)[2] == t.xelem (0, i) - 1 - && (*f).entities(0)[3] == t.xelem (1, i) - 1 - || (*f).entities(0)[0] == t.xelem (2, i) - 1 - && (*f).entities(0)[1] == t.xelem (3, i) - 1 - && (*f).entities(0)[2] == t.xelem (1, i) - 1 - && (*f).entities(0)[3] == t.xelem (0, i) - 1 - - || (*f).entities(0)[0] == t.xelem (3, i) - 1 - && (*f).entities(0)[1] == t.xelem (0, i) - 1 - && (*f).entities(0)[2] == t.xelem (1, i) - 1 - && (*f).entities(0)[3] == t.xelem (2, i) - 1 - || (*f).entities(0)[0] == t.xelem (3, i) - 1 - && (*f).entities(0)[1] == t.xelem (0, i) - 1 - && (*f).entities(0)[2] == t.xelem (2, i) - 1 - && (*f).entities(0)[3] == t.xelem (1, i) - 1 - || (*f).entities(0)[0] == t.xelem (3, i) - 1 - && (*f).entities(0)[1] == t.xelem (1, i) - 1 - && (*f).entities(0)[2] == t.xelem (0, i) - 1 - && (*f).entities(0)[3] == t.xelem (2, i) - 1 - || (*f).entities(0)[0] == t.xelem (3, i) - 1 - && (*f).entities(0)[1] == t.xelem (1, i) - 1 - && (*f).entities(0)[2] == t.xelem (2, i) - 1 - && (*f).entities(0)[3] == t.xelem (0, i) - 1 - || (*f).entities(0)[0] == t.xelem (3, i) - 1 - && (*f).entities(0)[1] == t.xelem (2, i) - 1 - && (*f).entities(0)[2] == t.xelem (0, i) - 1 - && (*f).entities(0)[3] == t.xelem (1, i) - 1 - || (*f).entities(0)[0] == t.xelem (3, i) - 1 - && (*f).entities(0)[1] == t.xelem (2, i) - 1 - && (*f).entities(0)[2] == t.xelem (1, i) - 1 - && (*f).entities(0)[3] == t.xelem (0, i) - 1) - { - cell.set_value ((*f).index (), t.xelem (4, i), msh); - break; - } - } - } - } - - *(msh.domains ().markers (D)) = cell; - - boost::shared_ptr<const dolfin::Mesh> my_m (new dolfin::Mesh (msh)); - pmsh = my_m; - } -}