COOLFluiD  Release kernel
COOLFluiD is a Collaborative Simulation Environment (CSE) focused on complex MultiPhysics simulations.
utest-mesh-construction.cpp
Go to the documentation of this file.
1 // Copyright (C) 2010-2013 von Karman Institute for Fluid Dynamics, Belgium
2 //
3 // This software is distributed under the terms of the
4 // GNU Lesser General Public License version 3 (LGPLv3).
5 // See doc/lgpl.txt and doc/gpl.txt for the license text.
6 
7 #define BOOST_TEST_DYN_LINK
8 #define BOOST_TEST_MODULE "Tests mesh construction"
9 
10 #include <iostream>
11 #include <boost/test/unit_test.hpp>
12 #include <boost/foreach.hpp>
13 
14 #include "common/Log.hpp"
15 #include "common/Core.hpp"
17 #include "common/Group.hpp"
18 #include "common/List.hpp"
19 
20 #include "mesh/Connectivity.hpp"
21 #include "mesh/Mesh.hpp"
22 #include "mesh/Region.hpp"
23 #include "mesh/Elements.hpp"
24 #include "mesh/Dictionary.hpp"
25 #include "mesh/Space.hpp"
26 #include "mesh/Field.hpp"
27 #include "mesh/ElementData.hpp"
28 #include "mesh/ElementType.hpp"
30 
31 #include "mesh/MeshWriter.hpp"
32 
33 using namespace std;
34 using namespace boost;
35 using namespace cf3;
36 using namespace cf3::mesh;
37 using namespace cf3::common;
38 
40  public boost::iterator_facade<ElementIterator, // iterator
41  Entity, // Value
42  boost::random_access_traversal_tag, // search direction
43  Entity // return type of dereference
44  >
45 {
46  typedef boost::iterator_facade<ElementIterator,
47  Entity,
48  boost::random_access_traversal_tag,
49  Entity
50  > BaseT;
51 public:
52 
53  typedef BaseT::difference_type difference_type;
54 
58  explicit ElementIterator(Entities& entities,
59  const Uint element_idx=0)
60  : element(entities,element_idx) {}
61 
62 private:
63  friend class boost::iterator_core_access;
64 
65  bool equal(ElementIterator const& other) const { return (element == other.element); }
66 
67  void increment()
68  {
69  ++element.idx;
70  }
71 
72  void decrement()
73  {
74  --element.idx;
75  }
76 
77  void advance(const difference_type n) { element.idx += n; }
78 
79  template <typename T2>
80  difference_type distance_to(ElementIterator const& other) const
81  {
82  return other.element.idx - element.idx;
83  }
84 
85 public:
86 
88  const Entity& dereference() const { return element; }
89 
90 private:
91 
92  Entity element;
93 };
94 
96 {
97  ElementIterator begin(entities);
98  return boost::make_iterator_range(begin,begin+static_cast<ElementIterator::difference_type>(entities.size()));
99 }
100 
101 
102 
104 {
105 
106  SpaceElement(Space& component, const Uint index=0) :
107  comp( &component ),
108  idx(index)
109  {
110  cf3_assert(idx<comp->support().size());
111  }
112 
115 
117  const ShapeFunction& element_type() const { return comp->shape_function(); }
118 
120  Dictionary& dict() const { return comp->dict(); }
121 
122  Entity support() const { return Entity(comp->support(),idx); }
123 
124  Connectivity::ConstRow field_indices() const { return comp->connectivity()[idx]; }
125  Uint glb_idx() const { return comp->support().glb_idx()[idx]; }
126  Uint rank() const { return comp->support().rank()[idx]; }
127  bool is_ghost() const { return comp->support().is_ghost(idx); }
128  RealMatrix get_coordinates() const { return comp->get_coordinates(idx); }
129  void put_coordinates(RealMatrix& coordinates) const { return comp->put_coordinates(coordinates,idx); }
130  void allocate_coordinates(RealMatrix& coordinates) const { return comp->allocate_coordinates(coordinates); }
131 };
132 
133 
135 
137 {
140  {
141  // uncomment if you want to use arguments to the test executable
142  //int* argc = &boost::unit_test::framework::master_test_suite().argc;
143  //char*** argv = &boost::unit_test::framework::master_test_suite().argv;
144  }
145 
148  {
149  }
150 
152 
154 
156  RealVector create_coord(const Real& x, const Real& y) {
157  RealVector coordVec(2);
158  coordVec[XX]=x;
159  coordVec[YY]=y;
160  return coordVec;
161  }
162 
164  std::vector<Uint> create_quad(const Uint& A, const Uint& B, const Uint& C, const Uint& D) {
165  Uint quad[] = {A,B,C,D};
166  std::vector<Uint> quadVec;
167  quadVec.assign(quad,quad+4);
168  return quadVec;
169  }
170 
172  std::vector<Uint> create_triag(const Uint& A, const Uint& B, const Uint& C) {
173  Uint triag[] = {A,B,C};
174  std::vector<Uint> triagVec;
175  triagVec.assign(triag,triag+3);
176  return triagVec;
177  }
178 
179  std::vector<Uint> create_quad_p2(const Uint A, const Uint B, const Uint C, const Uint D,
180  const Uint E, const Uint F, const Uint G, const Uint H, const Uint I)
181  {
182  Uint quad[] = {A,B,C,D,E,F,G,H,I};
183  std::vector<Uint> quadVec;
184  quadVec.assign(quad,quad+9);
185  return quadVec;
186  }
187 
189  std::vector<Uint> create_triag_p2(const Uint A, const Uint B, const Uint C,
190  const Uint D, const Uint E, const Uint F)
191  {
192  Uint triag[] = {A,B,C,D,E,F};
193  std::vector<Uint> triagVec;
194  triagVec.assign(triag,triag+6);
195  return triagVec;
196  }
197 
199 
200 };
201 
203 
204 BOOST_FIXTURE_TEST_SUITE( MeshConstruction_TestSuite, MeshConstruction_Fixture )
205 
206 
208 BOOST_AUTO_TEST_CASE( P1_2D_MeshConstruction )
209 {
210  AssertionManager::instance().AssertionThrows = true;
211 
212  const Uint dim=2;
213 
214  // Create root and mesh component
215  Component& root = Core::instance().root();
216 
217  Mesh& mesh = *root.create_component<Mesh>( "mesh" ) ;
218 
219  // create regions
220  Region& superRegion = mesh.topology().create_region("superRegion");
221  Dictionary& nodes = mesh.geometry_fields();
222  mesh.initialize_nodes(0,dim);
223  BOOST_CHECK_EQUAL(nodes.coordinates().row_size() , dim);
224 
225  Elements& quadRegion = superRegion.create_elements("cf3.mesh.LagrangeP1.Quad2D",nodes);
226  Elements& triagRegion = superRegion.create_elements("cf3.mesh.LagrangeP1.Triag2D",nodes);
227 
228  Table<Uint>::Buffer qTableBuffer = quadRegion.geometry_space().connectivity().create_buffer();
229  Table<Uint>::Buffer tTableBuffer = triagRegion.geometry_space().connectivity().create_buffer();
230  Table<Real>::Buffer coordinatesBuffer = nodes.coordinates().create_buffer();
231 
232  // Mesh of quads and triangles with node and element numbering:
233  //
234  // 5----4----6
235  // | |\ 3 |
236  // | | \ |
237  // | 1 | \ |
238  // | | 2 \|
239  // 3----2----7
240  // | |\ 1 |
241  // | | \ |
242  // | 0 | \ |
243  // | | 0 \|
244  // 0----1----8
245 
246  // fill coordinates in the buffer
247  coordinatesBuffer.add_row(create_coord( 0.0 , 0.0 )); // 0
248  coordinatesBuffer.add_row(create_coord( 1.0 , 0.0 )); // 1
249  coordinatesBuffer.add_row(create_coord( 1.0 , 1.0 )); // 2
250  coordinatesBuffer.add_row(create_coord( 0.0 , 1.0 )); // 3
251  coordinatesBuffer.add_row(create_coord( 1.0 , 2.0 )); // 4
252  coordinatesBuffer.add_row(create_coord( 0.0 , 2.0 )); // 5
253  coordinatesBuffer.add_row(create_coord( 2.0 , 2.0 )); // 6
254  coordinatesBuffer.add_row(create_coord( 2.0 , 1.0 )); // 7
255  coordinatesBuffer.add_row(create_coord( 2.0 , 0.0 )); // 8
256 
257  // fill connectivity in the buffer
258  qTableBuffer.add_row(create_quad( 0 , 1 , 2 , 3 ));
259  qTableBuffer.add_row(create_quad( 3 , 2 , 4 , 5 ));
260 
261  tTableBuffer.add_row(create_triag( 1 , 8 , 2 ));
262  tTableBuffer.add_row(create_triag( 8 , 7 , 2 ));
263  tTableBuffer.add_row(create_triag( 2 , 7 , 4 ));
264  tTableBuffer.add_row(create_triag( 7 , 6 , 4 ));
265 
266  // flush buffers into the table.
267  // This causes the table and array to be resized and filled.
268  coordinatesBuffer.flush();
269  qTableBuffer.flush();
270  tTableBuffer.flush();
271 
272  nodes.resize(nodes.coordinates().size());
273  quadRegion.resize(quadRegion.size());
274  triagRegion.resize(triagRegion.size());
275 
276  mesh.raise_mesh_loaded();
277  BOOST_CHECK(true);
278 
279  // check if coordinates match
280  Uint elem=1;
281  Uint node=2;
282 
283  Table<Uint>::ConstRow nodesRef = triagRegion.geometry_space().connectivity()[elem];
284  Table<Real>::Row coordRef = triagRegion.geometry_fields().coordinates()[nodesRef[node]];
285  BOOST_CHECK_EQUAL(coordRef[0],1.0);
286  BOOST_CHECK_EQUAL(coordRef[1],1.0);
287 
288  // calculate all volumes of a region
289  BOOST_FOREACH( Elements& region, find_components_recursively<Elements>(superRegion))
290  {
291  const ElementType& elementType = region.element_type();
292  const Uint nbRows = region.geometry_space().connectivity().size();
293  std::vector<Real> volumes(nbRows);
294  const Table<Real>& region_coordinates = region.geometry_fields().coordinates();
295  const Table<Uint>& region_connTable = region.geometry_space().connectivity();
296 
297  // the loop
298  RealMatrix elementCoordinates(elementType.nb_nodes(), elementType.dimension());
299  for (Uint iElem=0; iElem<nbRows; ++iElem)
300  {
301  fill(elementCoordinates, region_coordinates, region_connTable[iElem]);
302 
303  volumes[iElem]=elementType.volume(elementCoordinates);
304 
305  // check
306  if(elementType.shape() == GeoShape::QUAD)
307  BOOST_CHECK_EQUAL(volumes[iElem],1.0);
308  if(elementType.shape() == GeoShape::TRIAG)
309  BOOST_CHECK_EQUAL(volumes[iElem],0.5);
310  }
311  }
312 
313 // BOOST_FOREACH(Table<Real>::Row node , elem_coord)
314 // {
315 // CFinfo << "node = ";
316 // for (Uint j=0; j<node.size(); j++) {
317 // CFinfo << node[j] << " ";
318 // }
319 // CFinfo << "\n" << CFflush;
320 // }
321 
322  boost::shared_ptr< MeshWriter > meshwriter = build_component_abstract_type<MeshWriter>("cf3.mesh.gmsh.Writer","meshwriter");
323  meshwriter->write_from_to(mesh,"p1-mesh.msh");
324 
325 
326  Entity entity(quadRegion,0);
327 
328 
329  SpaceElement space_elem(quadRegion.space(nodes),0);
330  BOOST_CHECK_EQUAL( space_elem.field_indices()[0] , 0 );
331 
332  RealMatrix geom_coords = space_elem.support().get_coordinates();
333  RealMatrix space_coords = space_elem.get_coordinates();
334 
335  BOOST_CHECK( geom_coords == space_coords );
336 
337 
338  for(SpaceElement elem(quadRegion.space(nodes)); elem.idx<quadRegion.size(); ++elem.idx)
339  {
340  CFinfo << "indices = ";
341  boost_foreach(Uint field_idx, elem.field_indices())
342  CFinfo << field_idx << " ";
343  CFinfo << CFendl;
344  }
345  ElementIterator it(quadRegion);
346  boost_foreach( const Entity& entity, make_range(quadRegion) )
347  {
348  std::cout << entity << std::endl;
349  }
350 
351 }
352 
354 
355 BOOST_AUTO_TEST_CASE( P2_2D_MeshConstruction )
356 {
357  AssertionManager::instance().AssertionThrows = true;
358 
359  const Uint dim=2;
360 
361  // Create root and mesh component
362  boost::shared_ptr<Component> root = boost::static_pointer_cast<Component>(allocate_component<Group>("root"));
363 
364  Mesh& mesh = *root->create_component<Mesh>( "mesh" );
365 
366  // create regions
367  Region& superRegion = mesh.topology().create_region("superRegion");
368  Dictionary& nodes = mesh.geometry_fields();
369  mesh.initialize_nodes(0,dim);
370  BOOST_CHECK_EQUAL(nodes.coordinates().row_size() , dim);
371  Elements& quadRegion = superRegion.create_elements("cf3.mesh.LagrangeP2.Quad2D",nodes);
372  Elements& triagRegion = superRegion.create_elements("cf3.mesh.LagrangeP2.Triag2D",nodes);
373 
374  Table<Uint>::Buffer qTableBuffer = quadRegion.geometry_space().connectivity().create_buffer();
375  Table<Uint>::Buffer tTableBuffer = triagRegion.geometry_space().connectivity().create_buffer();
376  Table<Real>::Buffer coordinatesBuffer = nodes.coordinates().create_buffer();
377 
378  // Mesh of quads and triangles with node numbering and element numbering in brackets:
379  //
380  // Y ^
381  // |
382  // 2 - 5--15---4--24--6
383  // | | |\ |
384  // | | | \ (3)|
385  // 1.5- 16 17 14 22 23
386  // | | | \ |
387  // | | (1) |(2) \ |
388  // | | | \|
389  // 1 - 3--11---2--21--7
390  // | | |\ |
391  // | | | \ (1)|
392  // 0.5- 12 13 10 19 20
393  // | | | \ |
394  // | | (0) | (0)\ |
395  // | | | \|
396  // 0 - 0---9---1--18--8
397  // |
398  // o--|---|---|--|--|------> X
399  // 0 0.5 1 1.5 2
400 
401  // fill coordinates in the buffer of the P1 mesh
402  coordinatesBuffer.add_row(create_coord( 0.0 , 0.0 )); // 0
403  coordinatesBuffer.add_row(create_coord( 1.0 , 0.0 )); // 1
404  coordinatesBuffer.add_row(create_coord( 1.0 , 1.0 )); // 2
405  coordinatesBuffer.add_row(create_coord( 0.0 , 1.0 )); // 3
406  coordinatesBuffer.add_row(create_coord( 1.0 , 2.0 )); // 4
407  coordinatesBuffer.add_row(create_coord( 0.0 , 2.0 )); // 5
408  coordinatesBuffer.add_row(create_coord( 2.0 , 2.0 )); // 6
409  coordinatesBuffer.add_row(create_coord( 2.0 , 1.0 )); // 7
410  coordinatesBuffer.add_row(create_coord( 2.0 , 0.0 )); // 8
411 
412  // enrich
413  coordinatesBuffer.add_row(create_coord( 0.5 , 0.0 )); // 9
414  coordinatesBuffer.add_row(create_coord( 1.0 , 0.5 )); // 10
415  coordinatesBuffer.add_row(create_coord( 0.5 , 1.0 )); // 11
416  coordinatesBuffer.add_row(create_coord( 0.0 , 0.5 )); // 12
417  coordinatesBuffer.add_row(create_coord( 0.5 , 0.5 )); // 13
418  coordinatesBuffer.add_row(create_coord( 1.0 , 1.5 )); // 14
419  coordinatesBuffer.add_row(create_coord( 0.5 , 2.0 )); // 15
420  coordinatesBuffer.add_row(create_coord( 0.0 , 1.5 )); // 16
421  coordinatesBuffer.add_row(create_coord( 0.5 , 1.5 )); // 17
422  coordinatesBuffer.add_row(create_coord( 1.5 , 0.0 )); // 18
423  coordinatesBuffer.add_row(create_coord( 1.5 , 0.5 )); // 19
424  coordinatesBuffer.add_row(create_coord( 2.0 , 0.5 )); // 20
425  coordinatesBuffer.add_row(create_coord( 1.5 , 1.0 )); // 21
426  coordinatesBuffer.add_row(create_coord( 1.5 , 1.5 )); // 22
427  coordinatesBuffer.add_row(create_coord( 2.0 , 1.5 )); // 23
428  coordinatesBuffer.add_row(create_coord( 1.5 , 2.0 )); // 24
429 
430 
431  // enrich
432 
433  // fill connectivity in the buffer
434  qTableBuffer.add_row(create_quad_p2( 0 , 1 , 2 , 3 , 9 , 10, 11, 12, 13 ));
435  qTableBuffer.add_row(create_quad_p2( 3 , 2 , 4 , 5 , 11, 14, 15, 16, 17 ));
436 
437  tTableBuffer.add_row(create_triag_p2( 1 , 8 , 2 , 18 , 19 , 10 ));
438  tTableBuffer.add_row(create_triag_p2( 8 , 7 , 2 , 20 , 21 , 19 ));
439  tTableBuffer.add_row(create_triag_p2( 2 , 7 , 4 , 21 , 22 , 14 ));
440  tTableBuffer.add_row(create_triag_p2( 7 , 6 , 4 , 23 , 24 , 22 ));
441 
442  // flush buffers into the table.
443  // This causes the table and array to be resized and filled.
444  coordinatesBuffer.flush();
445  qTableBuffer.flush();
446  tTableBuffer.flush();
447 
448  nodes.resize(nodes.coordinates().size());
449  quadRegion.resize(quadRegion.size());
450  triagRegion.resize(triagRegion.size());
451 
452  mesh.raise_mesh_loaded();
453  BOOST_CHECK(true);
454 
455  // check if coordinates match
456  Uint elem=1;
457  Uint node=2;
458 
459  Table<Uint>::ConstRow nodesRef = triagRegion.geometry_space().connectivity()[elem];
460  Table<Real>::Row coordRef = triagRegion.geometry_fields().coordinates()[nodesRef[node]];
461  BOOST_CHECK_EQUAL(coordRef[0],1.0);
462  BOOST_CHECK_EQUAL(coordRef[1],1.0);
463 
464 // // calculate all volumes of a region
465 // BOOST_FOREACH( Elements& region, find_components_recursively<Elements>(superRegion))
466 // {
467 // const ElementType& elementType = region.element_type();
468 // const Uint nbRows = region.geometry_space().connectivity().size();
469 // std::vector<Real> volumes(nbRows);
470 // const Table<Real>& region_coordinates = region.coordinates();
471 // const Table<Uint>& region_connTable = region.geometry_space().connectivity();
472 // // the loop
473 // ElementType::NodesT elementCoordinates(elementType.nb_nodes(), elementType.dimension());
474 // for (Uint iElem=0; iElem<nbRows; ++iElem)
475 // {
476 // elementCoordinates.fill(region_coordinates, region_connTable[iElem]);
477 //
478 // volumes[iElem]=elementType.compute_volume(elementCoordinates);
479 //
480 // // check
481 // if(elementType.shape() == GeoShape::QUAD)
482 // BOOST_CHECK_EQUAL(volumes[iElem],1.0);
483 // if(elementType.shape() == GeoShape::TRIAG)
484 // BOOST_CHECK_EQUAL(volumes[iElem],0.5);
485 // }
486 // }
487 //
488 // // BOOST_FOREACH(Table<Real>::Row node , elem_coord)
489 // // {
490 // // CFinfo << "node = ";
491 // // for (Uint j=0; j<node.size(); j++) {
492 // // CFinfo << node[j] << " ";
493 // // }
494 // // CFinfo << "\n" << CFflush;
495 // // }
496 
497 
498  boost::shared_ptr< MeshWriter > meshwriter = build_component_abstract_type<MeshWriter>("cf3.mesh.gmsh.Writer","meshwriter");
499  meshwriter->write_from_to(mesh,"p2-mesh.msh");
500 
501 }
502 
504 
505 BOOST_AUTO_TEST_SUITE_END()
506 
507 
#define F(x, y, z)
#define CFinfo
these are always defined
Definition: Log.hpp:104
const ShapeFunction & shape_function() const
Definition: Space.cpp:102
RealMatrix get_coordinates() const
Definition: Entities.cpp:213
external boost library namespace
Space & geometry_space() const
Definition: Entities.hpp:94
void initialize_nodes(const Uint nb_nodes, const Uint dimension)
will among others set the coordinate dimension for the nodes
Definition: Mesh.cpp:127
difference_type distance_to(ElementIterator const &other) const
Connectivity::ConstRow field_indices() const
std::vector< Uint > create_triag_p2(const Uint A, const Uint B, const Uint C, const Uint D, const Uint E, const Uint F)
create a Uint vector with 3 node ID's
#define cf3_assert(a)
Definition: Assertions.hpp:93
ElementType & element_type() const
return the elementType
Definition: Entities.cpp:116
Region & create_region(const std::string &name)
Definition: Region.cpp:50
STL namespace.
#define Mesh_API
Definition: LibMesh.hpp:24
Buffer create_buffer(const size_t buffersize=16384)
Definition: Table.hpp:102
bool is_ghost(const Uint idx) const
Definition: Entities.cpp:197
#define boost_foreach
lowercase version of BOOST_FOREACH
Definition: Foreach.hpp:16
const Field & coordinates() const
Definition: Dictionary.cpp:481
common::List< Uint > & rank()
Definition: Entities.hpp:71
void advance(const difference_type n)
tuple root
Definition: coolfluid.py:24
Entities & support() const
Access the geometric support.
Definition: Space.hpp:102
#define CFendl
Definition: Log.hpp:109
boost::proto::terminal< SFOp< NodesOp > >::type const nodes
~MeshConstruction_Fixture()
common tear-down for each test case
Dictionary & geometry_fields() const
Const access to the coordinates.
Definition: Entities.hpp:63
Definition: Defs.hpp:17
Eigen::Matrix< Real, Eigen::Dynamic, Eigen::Dynamic > RealMatrix
Dynamic sized matrix of Real scalars.
Definition: MatrixTypes.hpp:22
Uint size() const
Definition: Table.hpp:127
SpaceElement(Space &component, const Uint index=0)
ElementIterator(Entities &entities, const Uint element_idx=0)
Uint size() const
return the number of elements
Definition: Entities.cpp:161
RealMatrix get_coordinates() const
Entity support() const
const ShapeFunction & element_type() const
return the elementType
common::List< Uint > & glb_idx()
Mutable access to the list of nodes.
Definition: Entities.hpp:66
Dictionary & dict() const
Const access to the coordinates.
Basic Classes for Mesh applications used by COOLFluiD.
virtual Real volume(const RealMatrix &nodes) const =0
const Entity & dereference() const
dereferencing
Eigen::Matrix< Real, Eigen::Dynamic, 1 > RealVector
Dynamic sized column vector.
Definition: MatrixTypes.hpp:25
BaseT::difference_type difference_type
GeoShape::Type shape() const
Definition: ElementType.hpp:64
bool equal(ElementIterator const &other) const
boost::iterator_facade< ElementIterator, Entity, boost::random_access_traversal_tag, Entity > BaseT
void allocate_coordinates(RealMatrix &coordinates) const
Allocate element coordinates.
Definition: Space.cpp:172
Top-level namespace for coolfluid.
Definition: Action.cpp:18
std::vector< Uint > create_quad(const Uint &A, const Uint &B, const Uint &C, const Uint &D)
create a Uint vector with 4 node ID's
void resize(const Uint nb_elem)
Definition: Entities.cpp:246
BOOST_AUTO_TEST_CASE(P1_2D_MeshConstruction)
void fill(NodeValuesT &to_fill, const common::Table< Real > &data_array, const RowT &element_row, const Uint start=0)
Fill STL-vector like per-node data storage.
Definition: ElementData.hpp:28
Component holding a 2 dimensional array of a templated type.
Definition: Table.hpp:45
RealMatrix get_coordinates(const Uint elem_idx) const
Lookup element coordinates.
Definition: Space.cpp:179
RealVector create_coord(const Real &x, const Real &y)
possibly common functions used on the tests below
void raise_mesh_loaded()
Definition: Mesh.cpp:479
Uint row_size(Uint i=0) const
Definition: Table.hpp:133
boost::proto::terminal< SFOp< CoordinatesOp > >::type const coordinates
boost::iterator_range< ElementIterator > make_range(Entities &entities)
Low storage struct to uniquely identify one element.
Definition: Entities.hpp:141
unsigned int Uint
typedef for unsigned int
Definition: CF.hpp:90
MeshConstruction_Fixture()
common setup for each test case
void allocate_coordinates(RealMatrix &coordinates) const
void resize(const Uint size)
Resize the contained fields.
Definition: Dictionary.cpp:106
Uint nb_nodes() const
Definition: ElementType.hpp:73
Region & topology() const
Definition: Mesh.hpp:51
std::vector< Uint > create_triag(const Uint &A, const Uint &B, const Uint &C)
create a Uint vector with 3 node ID's
Definition: Defs.hpp:17
Dictionary & dict() const
Dictionary this space belongs to.
Definition: Space.hpp:92
Dictionary & geometry_fields() const
Definition: Mesh.cpp:339
Base class for defining CF components.
Definition: Component.hpp:82
Space component class.
Definition: Space.hpp:59
Handle< Component > create_component(const std::string &name, const std::string &builder)
Build a (sub)component of this component using the extended type_name of the component.
Definition: Component.cpp:568
Most basic kernel library.
Definition: Action.cpp:19
std::vector< Uint > create_quad_p2(const Uint A, const Uint B, const Uint C, const Uint D, const Uint E, const Uint F, const Uint G, const Uint H, const Uint I)
void put_coordinates(RealMatrix &coordinates, const Uint elem_idx) const
Lookup element coordinates.
Definition: Space.cpp:153
Uint dimension() const
Definition: ElementType.hpp:82
Connectivity & connectivity()
connectivity table to dictionary entries
Definition: Space.hpp:110
void put_coordinates(RealMatrix &coordinates) const
Send comments to:
COOLFluiD Web Admin