COOLFluiD  Release kernel
COOLFluiD is a Collaborative Simulation Environment (CSE) focused on complex MultiPhysics simulations.
utest-ufem-buildsparsity.cpp
Go to the documentation of this file.
1 // Copyright (C) 2010-2011 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 "Test module for heat-conduction related proto operations"
9 
10 #include <boost/assign.hpp>
11 #include <boost/test/unit_test.hpp>
12 
13 #include "common/Core.hpp"
14 #include "common/Environment.hpp"
15 
17 
18 #include "math/LSS/System.hpp"
19 
20 #include "mesh/Domain.hpp"
22 
23 #include "solver/Model.hpp"
24 
27 
29 #include "mesh/MeshGenerator.hpp"
30 
31 #include "UFEM/LSSAction.hpp"
32 #include "UFEM/Solver.hpp"
33 #include "UFEM/SparsityBuilder.hpp"
34 #include "UFEM/Tags.hpp"
35 #include "math/LSS/SolveLSS.hpp"
36 
37 using namespace cf3;
38 using namespace cf3::solver;
39 using namespace cf3::solver::actions;
40 using namespace cf3::solver::actions::Proto;
41 using namespace cf3::common;
42 using namespace cf3::math;
43 using namespace cf3::math::Consts;
44 using namespace cf3::mesh;
45 
46 using namespace boost::assign;
47 
49 {
51  root( Core::instance().root() )
52  {
53  }
54 
56 };
57 
58 BOOST_FIXTURE_TEST_SUITE( UFEMBuildSparsitySuite, UFEMBuildSparsityFixture )
59 
61 {
62  common::PE::Comm::instance().init(boost::unit_test::framework::master_test_suite().argc, boost::unit_test::framework::master_test_suite().argv);
63  BOOST_CHECK_EQUAL(common::PE::Comm::instance().size(), 1);
64 }
65 
66 BOOST_AUTO_TEST_CASE( Sparsity1D )
67 {
68  Core::instance().environment().options().set("log_level", 4u);
69 
70  // Parameters
71  Real length = 5.;
72  const Uint nb_segments = 5 ;
73  const Uint nb_nodes = nb_segments + 1;
74 
75  // Setup a model
76  Model& model = *root.create_component<Model>("Model");
77  Domain& domain = model.create_domain("Domain");
78 
79  LSS::System& lss = *model.create_component<LSS::System>("LSS");
80  lss.options().option("matrix_builder").change_value(std::string("cf3.math.LSS.TrilinosFEVbrMatrix"));
81 
82  // Setup mesh
83  // Mesh& mesh = *domain.create_component<Mesh>("Mesh");
84  // Tools::MeshGeneration::create_line(mesh, length, nb_segments);
85  boost::shared_ptr<MeshGenerator> create_line = build_component_abstract_type<MeshGenerator>("cf3.mesh.SimpleMeshGenerator","create_line");
86  create_line->options().set("mesh",domain.uri()/"Mesh");
87  create_line->options().set("lengths",std::vector<Real>(DIM_1D, length));
88  create_line->options().set("nb_cells",std::vector<Uint>(DIM_1D, nb_segments));
89  Mesh& mesh = create_line->generate();
90 
91  // Setup sparsity
92  std::vector<Uint> node_connectivity, starting_indices;
93  Handle< List<Uint> > gids = domain.create_component< List<Uint> >("GIDs");
94  Handle< List<Uint> > ranks = domain.create_component< List<Uint> >("Ranks");
95  Handle< List<int> > used_node_map = domain.create_component< List<int> >("used_node_map");
96  UFEM::build_sparsity(std::vector< Handle<Region> >(1, mesh.topology().handle<Region>()), mesh.geometry_fields(), node_connectivity, starting_indices, *gids, *ranks, *used_node_map);
97 
98  // Check result
99  BOOST_CHECK_EQUAL(starting_indices[0], 0u);
100  for(Uint i = 2; i != nb_nodes; ++i)
101  BOOST_CHECK_EQUAL(starting_indices[i] - starting_indices[i-1], 3);
102 
103  PE::CommPattern& comm_pattern = *domain.create_component<PE::CommPattern>("CommPattern");
104  comm_pattern.insert("gid",gids->array(),false);
105  comm_pattern.setup(Handle<PE::CommWrapper>(comm_pattern.get_child("gid")),ranks->array());
106 
107  // Create the LSS
108  lss.create(comm_pattern, 1u, node_connectivity, starting_indices);
109 
110 
111  // Write the matrix
112  lss.matrix()->print("utest-ufem-buildsparsity_heat_matrix_1D.plt");
113 }
114 
115 BOOST_AUTO_TEST_CASE( Sparsity2DQuads )
116 {
117  Core::instance().environment().options().set("log_level", 4u);
118 
119  // Parameters
120  Real length = 5.;
121  const Uint nb_segments = 5 ;
122  const Uint nb_nodes = (nb_segments+1) * (nb_segments+1);
123 
124  // Setup a model
125  Model& model = *root.create_component<Model>("Model");
126  Domain& domain = model.create_domain("Domain");
127 
128  LSS::System& lss = *model.create_component<LSS::System>("LSS");
129  lss.options().option("matrix_builder").change_value(std::string("cf3.math.LSS.TrilinosFEVbrMatrix"));
130 
131  // Setup mesh
132  Mesh& mesh = *domain.create_component<Mesh>("Mesh");
133  Tools::MeshGeneration::create_rectangle(mesh, length, length, nb_segments, nb_segments);
134 
135  // Setup sparsity
136  std::vector<Uint> node_connectivity, starting_indices;
137  Handle< List<Uint> > gids = domain.create_component< List<Uint> >("GIDs");
138  Handle< List<Uint> > ranks = domain.create_component< List<Uint> >("Ranks");
139  Handle< List<int> > used_node_map = domain.create_component< List<int> >("used_node_map");
140  UFEM::build_sparsity(std::vector< Handle<Region> >(1, mesh.topology().handle<Region>()), mesh.geometry_fields(), node_connectivity, starting_indices, *gids, *ranks, *used_node_map);
141 
142  PE::CommPattern& comm_pattern = *domain.create_component<PE::CommPattern>("CommPattern");
143  comm_pattern.insert("gid",gids->array(),false);
144  comm_pattern.setup(Handle<PE::CommWrapper>(comm_pattern.get_child("gid")),ranks->array());
145 
146  // Create the LSS
147  lss.create(comm_pattern, 1u, node_connectivity, starting_indices);
148 
149 
150  // Write the matrix
151  lss.matrix()->print("utest-ufem-buildsparsity_heat_matrix_2DQuads.plt");
152 }
153 
154 BOOST_AUTO_TEST_CASE( Sparsity2DTris )
155 {
156  Core::instance().environment().options().set("log_level", 4u);
157 
158  // Parameters
159  Real length = 5.;
160  const Uint nb_segments = 5 ;
161  const Uint nb_nodes = (nb_segments+1) * (nb_segments+1);
162 
163  // Setup a model
164  Model& model = *root.create_component<Model>("Model");
165  Domain& domain = model.create_domain("Domain");
166 
167  LSS::System& lss = *model.create_component<LSS::System>("LSS");
168  lss.options().option("matrix_builder").change_value(std::string("cf3.math.LSS.TrilinosFEVbrMatrix"));
169 
170  // Setup mesh
171  Mesh& mesh = *domain.create_component<Mesh>("Mesh");
172  Tools::MeshGeneration::create_rectangle_tris(mesh, length, length, nb_segments, nb_segments);
173 
174  // Setup sparsity
175  std::vector<Uint> node_connectivity, starting_indices;
176  Handle< List<Uint> > gids = domain.create_component< List<Uint> >("GIDs");
177  Handle< List<Uint> > ranks = domain.create_component< List<Uint> >("Ranks");
178  Handle< List<int> > used_node_map = domain.create_component< List<int> >("used_node_map");
179  UFEM::build_sparsity(std::vector< Handle<Region> >(1, mesh.topology().handle<Region>()), mesh.geometry_fields(), node_connectivity, starting_indices, *gids, *ranks, *used_node_map);
180 
181  PE::CommPattern& comm_pattern = *domain.create_component<PE::CommPattern>("CommPattern");
182  comm_pattern.insert("gid",gids->array(),false);
183  comm_pattern.setup(Handle<PE::CommWrapper>(comm_pattern.get_child("gid")),ranks->array());
184 
185  // Create the LSS
186  lss.create(comm_pattern, 1u, node_connectivity, starting_indices);
187 
188 
189  // Write the matrix
190  lss.matrix()->print("utest-ufem-buildsparsity_heat_matrix_2DTris.plt");
191 }
192 
193 // Single block, meshed with the blockmesher
194 BOOST_AUTO_TEST_CASE( Sparsity3DHexaBlock )
195 {
196  Core::instance().environment().options().set("log_level", 4u);
197 
198  // Parameters
199  Real length = 5.;
200  const Uint nb_segments = 10;
201  const Uint nb_nodes = (nb_segments+1) * (nb_segments+1) * (nb_segments+1);
202 
203  // Setup a model
204  Model& model = *root.create_component<Model>("Model");
205  Domain& domain = model.create_domain("Domain");
206 
207  LSS::System& lss = *model.create_component<LSS::System>("LSS");
208  lss.options().option("matrix_builder").change_value(std::string("cf3.math.LSS.TrilinosFEVbrMatrix"));
209 
210  // Setup mesh
211  Mesh& mesh = *domain.create_component<Mesh>("Mesh");
212 
213  BlockMesh::BlockArrays& blocks = *domain.create_component<BlockMesh::BlockArrays>("blocks");
214 
215  *blocks.create_points(3, 8) << 0. << 0. << 0.
216  << length << 0. << 0.
217  << 0. << length << 0.
218  << length << length << 0.
219  << 0. << 0. << length
220  << length << 0. << length
221  << 0. << length << length
222  << length << length << length;
223 
224  *blocks.create_blocks(1) << 0 << 1 << 3 << 2 << 4 << 5 << 7 << 6;
225  *blocks.create_block_subdivisions() << nb_segments << nb_segments << nb_segments;
226  *blocks.create_block_gradings() << 1. << 1. << 1. << 1. << 1. << 1. << 1. << 1. << 1. << 1. << 1. << 1.;
227 
228  *blocks.create_patch("bottomWall", 1) << 0 << 2 << 3 << 1;
229  *blocks.create_patch("topWall", 1) << 4 << 5 << 7 << 6;
230  *blocks.create_patch("side1", 1) << 1 << 3 << 7 << 5;
231  *blocks.create_patch("side2", 1) << 0 << 1 << 5 << 4;
232  *blocks.create_patch("side3", 1) << 0 << 4 << 6 << 2;
233  *blocks.create_patch("side4", 1) << 2 << 6 << 7 << 3;
234 
235  blocks.create_mesh(mesh);
236 
237  BOOST_CHECK_EQUAL(nb_nodes, mesh.geometry_fields().coordinates().size());
238 
239  // Setup sparsity
240  std::vector<Uint> node_connectivity, starting_indices;
241  Handle< List<Uint> > gids = domain.create_component< List<Uint> >("GIDs");
242  Handle< List<Uint> > ranks = domain.create_component< List<Uint> >("Ranks");
243  Handle< List<int> > used_node_map = domain.create_component< List<int> >("used_node_map");
244  UFEM::build_sparsity(std::vector< Handle<Region> >(1, mesh.topology().handle<Region>()), mesh.geometry_fields(), node_connectivity, starting_indices, *gids, *ranks, *used_node_map);
245 
246  PE::CommPattern& comm_pattern = *domain.create_component<PE::CommPattern>("CommPattern");
247  comm_pattern.insert("gid",gids->array(),false);
248  comm_pattern.setup(Handle<PE::CommWrapper>(comm_pattern.get_child("gid")),ranks->array());
249 
250  // Create the LSS
251  lss.create(comm_pattern, 1u, node_connectivity, starting_indices);
252 
253  // Write the matrix
254  lss.matrix()->print("utest-ufem-buildsparsity_heat_matrix_3DHexaBlock.plt");
255 }
256 
257 BOOST_AUTO_TEST_CASE( Sparsity3DHexaChannel )
258 {
259  Core::instance().environment().options().set("log_level", 4u);
260 
261  // Parameters
262  Real length = 5.;
263  const Uint nb_segments = 4;
264  const Uint nb_nodes = (nb_segments+1) * (nb_segments+1) * (nb_segments+1);
265 
266  // Setup a model
267  Model& model = *root.create_component<Model>("Model");
268  Domain& domain = model.create_domain("Domain");
269 
270  LSS::System& lss = *model.create_component<LSS::System>("LSS");
271  lss.options().option("matrix_builder").change_value(std::string("cf3.math.LSS.TrilinosFEVbrMatrix"));
272 
273  // Setup mesh
274  Mesh& mesh = *domain.create_component<Mesh>("Mesh");
275  BlockMesh::BlockArrays& blocks = *domain.create_component<BlockMesh::BlockArrays>("blocks");
276  Tools::MeshGeneration::create_channel_3d(blocks, length, length/8., length, nb_segments, nb_segments/2, nb_segments, 1.);
277  blocks.create_mesh(mesh);
278 
279  BOOST_CHECK_EQUAL(nb_nodes, mesh.geometry_fields().coordinates().size());
280 
281  // Setup sparsity
282  std::vector<Uint> node_connectivity, starting_indices;
283  Handle< List<Uint> > gids = domain.create_component< List<Uint> >("GIDs");
284  Handle< List<Uint> > ranks = domain.create_component< List<Uint> >("Ranks");
285  Handle< List<int> > used_node_map = domain.create_component< List<int> >("used_node_map");
286  UFEM::build_sparsity(std::vector< Handle<Region> >(1, mesh.topology().handle<Region>()), mesh.geometry_fields(), node_connectivity, starting_indices, *gids, *ranks, *used_node_map);
287 
288  PE::CommPattern& comm_pattern = *domain.create_component<PE::CommPattern>("CommPattern");
289  comm_pattern.insert("gid",gids->array(),false);
290  comm_pattern.setup(Handle<PE::CommWrapper>(comm_pattern.get_child("gid")),ranks->array());
291 
292  // Create the LSS
293  lss.create(comm_pattern, 1u, node_connectivity, starting_indices);
294 
295  // Write the matrix
296  lss.matrix()->print("utest-ufem-buildsparsity_heat_matrix_3DHexaChannel.plt");
297 }
298 
299 BOOST_AUTO_TEST_CASE( Heat1DComponent )
300 {
301  Core::instance().environment().options().set("log_level", 4u);
302 
303  // Parameters
304  Real length = 5.;
305  const Uint nb_segments = 5 ;
306  const Uint nb_nodes = nb_segments + 1;
307 
308  // Setup a model
309  Model& model = *root.create_component<Model>("Model");
310  Domain& domain = model.create_domain("Domain");
311  UFEM::Solver& solver = *model.create_component<UFEM::Solver>("Solver");
312 
313  Handle<UFEM::LSSAction> lss_action(solver.add_direct_solver("cf3.UFEM.LSSAction"));
314 
315  // Proto placeholders
316  FieldVariable<0, ScalarField> temperature("Temperature", UFEM::Tags::solution());
317 
318  // Allowed elements (reducing this list improves compile times)
319  boost::mpl::vector1<mesh::LagrangeP1::Line1D> allowed_elements;
320 
321  // add the top-level actions (assembly, BC and solve)
322  *lss_action
324  (
325  "Assembly",
327  (
328  allowed_elements,
329  group
330  (
331  _A = _0,
332  element_quadrature( _A(temperature) += transpose(nabla(temperature)) * nabla(temperature) ),
333  lss_action->system_matrix += _A
334  )
335  )
336  )
337  << allocate_component<UFEM::BoundaryConditions>("BoundaryConditions")
338  << allocate_component<math::LSS::SolveLSS>("SolveLSS")
339  << create_proto_action("Increment", nodes_expression(temperature += lss_action->solution(temperature)))
340  << create_proto_action("Output", nodes_expression(_cout << "T(" << coordinates(0,0) << ") = " << temperature << "\n"));
341 
342  // Setup physics
343  model.create_physics("cf3.physics.DynamicModel");
344 
345  // Setup mesh
346  Mesh& mesh = *domain.create_component<Mesh>("Mesh");
347  Tools::MeshGeneration::create_line(mesh, length, nb_segments);
348 
349  lss_action->options().set("regions", std::vector<URI>(1, mesh.topology().uri()));
350 
351  LSS::System& lss = lss_action->create_lss();
352 
353  // Write the matrix
354  lss.matrix()->print("utest-ufem-buildsparsity_heat_matrix_1DHeat.plt");
355 }
356 
358 
359 BOOST_AUTO_TEST_SUITE_END()
360 
361 
void create(cf3::common::PE::CommPattern &cp, Uint neq, std::vector< Uint > &node_connectivity, std::vector< Uint > &starting_indices, const std::vector< Uint > &periodic_links_nodes=std::vector< Uint >(), const std::vector< bool > &periodic_links_active=std::vector< bool >())
Definition: System.cpp:82
BOOST_AUTO_TEST_CASE(InitMPI)
static boost::proto::terminal< ZeroTag >::type _0
Placeholder for the zero matrix.
Definition: Terminals.hpp:209
void create_line(Mesh &mesh, const Real x_len, const Uint x_segments)
Create a 1D line mesh.
virtual mesh::Domain & create_domain(const std::string &name)
creates a domain in this model
Definition: Model.cpp:201
Safe pointer to an object. This is the supported method for referring to components.
Definition: Handle.hpp:39
virtual physics::PhysModel & create_physics(const std::string &builder)
Definition: Model.cpp:182
static boost::proto::terminal< ElementSystemMatrix< boost::mpl::int_< 0 > > >::type const _A
Some predefined element matrices (more can be user-defined, but you have to change the number in the ...
Handle< common::Action > add_direct_solver(const std::string &builder_name)
Definition: Solver.cpp:128
This header collects all the headers needed for the linear system solver, also including configure-ti...
Basic Classes for Mathematical applications used by COOLFluiD.
Parallel Communication Pattern. This class provides functionality to collect communication. For efficiency it works such a way that you submit your request via the constructor or the add/remove/move magic triangle and then call setup to modify the commpattern. The data needed to be kept synchronous can be registered via the insert function. The word node here means any kind of "point of storage", in this context it is not directly related with the computational mesh.
void setup(const Handle< CommWrapper > &gid, std::vector< Uint > &rank)
Basic Classes for Solver applications used by CF.
Definition: Action.cpp:29
URI uri() const
Construct the full path.
Definition: Component.cpp:248
Handle< LSS::Matrix > matrix()
Accessor to matrix.
Definition: System.hpp:145
static boost::proto::terminal< ElementQuadratureTag >::type element_quadrature
Use element_quadrature(expr1, expr2, ..., exprN) to evaluate a group of expressions.
const Field & coordinates() const
Definition: Dictionary.cpp:481
void create_rectangle_tris(Mesh &mesh, const Real x_len, const Real y_len, const Uint x_segments, const Uint y_segments)
Create a triangulated version of a 2D rectangular grid.
tuple root
Definition: coolfluid.py:24
boost::shared_ptr< ElementsExpression< ExprT, ElementTypes > > elements_expression(ElementTypes, const ExprT &expr)
Definition: Expression.hpp:292
void create_rectangle(Mesh &mesh, const Real x_len, const Real y_len, const Uint x_segments, const Uint y_segments)
Create a rectangular, 2D, quad-only mesh. No buffer for creation.
Manage a collection of UFEM solvers.
Definition: Solver.hpp:31
void insert(const std::string &name, T *&data, const int size, const unsigned int stride=1, const bool needs_update=true)
static boost::proto::terminal< std::ostream & >::type _cout
Wrap std::cout.
Definition: Terminals.hpp:219
boost::shared_ptr< List< Uint > > build_sparsity(const std::vector< Handle< Region > > &regions, const Dictionary &dictionary, std::vector< Uint > &node_connectivity, std::vector< Uint > &start_indices, List< Uint > &gids, List< Uint > &ranks, List< int > &used_node_map)
Static functions for mathematical constants.
Definition: Consts.hpp:25
Uint size() const
Definition: Table.hpp:127
boost::proto::terminal< TransposeFunction >::type const transpose
boost::proto::terminal< SFOp< NablaOp > >::type const nabla
virtual void change_value(const boost::any &value)
change the value of this option
Definition: Option.cpp:89
Basic Classes for Mesh applications used by COOLFluiD.
tuple model
Global confifuration.
Handle< Component > get_child(const std::string &name)
Definition: Component.cpp:441
void init(int argc=0, char **args=0)
Definition: Comm.cpp:80
Top-level namespace for coolfluid.
Definition: Action.cpp:18
void create_channel_3d(BlockArrays &blocks, const Real length, const Real half_height, const Real width, const Uint x_segs, const Uint y_segs_half, const Uint z_segs, const Real ratio)
boost::proto::terminal< SFOp< CoordinatesOp > >::type const coordinates
static boost::proto::terminal< ExpressionGroupTag >::type group
Use group(expr1, expr2, ..., exprN) to evaluate a group of expressions.
common::Environment & environment() const
Definition: Core.cpp:168
unsigned int Uint
typedef for unsigned int
Definition: CF.hpp:90
boost::shared_ptr< ProtoAction > create_proto_action(const std::string &name, const boost::shared_ptr< Expression > &expression)
Create a new ProtoAction, immediatly setting the expression.
Region & topology() const
Definition: Mesh.hpp:51
static Core & instance()
Definition: Core.cpp:37
Handle< Component > handle()
Get a handle to the component.
Definition: Component.hpp:179
boost::shared_ptr< NodesExpression< ExprT, boost::mpl::range_c< Uint, 1, 4 > > > nodes_expression(const ExprT &expr)
Definition: Expression.hpp:308
OptionList & options()
Definition: Component.cpp:856
Dictionary & geometry_fields() const
Definition: Mesh.cpp:339
static Comm & instance()
Return a reference to the current PE.
Definition: Comm.cpp:44
Base class for defining CF components.
Definition: Component.hpp:82
void set(const std::string &pname, const boost::any &val)
Definition: OptionList.cpp:132
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
const Option & option(const std::string &pname) const
get a constant option from the list
Definition: OptionList.cpp:52
Most basic kernel library.
Definition: Action.cpp:19
Handle< common::Table< Real > > create_points(const Uint dimensions, const Uint nb_points)
Create the table that holds the points for the blocks.
Definition: BlockData.cpp:1287
Send comments to:
COOLFluiD Web Admin