7 #define BOOST_TEST_DYN_LINK
8 #define BOOST_TEST_MODULE "Test module for heat-conduction related proto operations"
10 #include <boost/assign.hpp>
11 #include <boost/test/unit_test.hpp>
72 const Uint nb_segments = 5 ;
73 const Uint nb_nodes = nb_segments + 1;
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();
92 std::vector<Uint> node_connectivity, starting_indices;
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);
104 comm_pattern.
insert(
"gid",gids->array(),
false);
108 lss.
create(comm_pattern, 1u, node_connectivity, starting_indices);
112 lss.
matrix()->print(
"utest-ufem-buildsparsity_heat_matrix_1D.plt");
121 const Uint nb_segments = 5 ;
122 const Uint nb_nodes = (nb_segments+1) * (nb_segments+1);
136 std::vector<Uint> node_connectivity, starting_indices;
143 comm_pattern.insert(
"gid",gids->array(),
false);
147 lss.
create(comm_pattern, 1u, node_connectivity, starting_indices);
151 lss.
matrix()->print(
"utest-ufem-buildsparsity_heat_matrix_2DQuads.plt");
160 const Uint nb_segments = 5 ;
161 const Uint nb_nodes = (nb_segments+1) * (nb_segments+1);
175 std::vector<Uint> node_connectivity, starting_indices;
182 comm_pattern.insert(
"gid",gids->array(),
false);
186 lss.
create(comm_pattern, 1u, node_connectivity, starting_indices);
190 lss.
matrix()->print(
"utest-ufem-buildsparsity_heat_matrix_2DTris.plt");
200 const Uint nb_segments = 10;
201 const Uint nb_nodes = (nb_segments+1) * (nb_segments+1) * (nb_segments+1);
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;
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.;
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;
235 blocks.create_mesh(mesh);
240 std::vector<Uint> node_connectivity, starting_indices;
247 comm_pattern.insert(
"gid",gids->array(),
false);
251 lss.
create(comm_pattern, 1u, node_connectivity, starting_indices);
254 lss.
matrix()->print(
"utest-ufem-buildsparsity_heat_matrix_3DHexaBlock.plt");
263 const Uint nb_segments = 4;
264 const Uint nb_nodes = (nb_segments+1) * (nb_segments+1) * (nb_segments+1);
277 blocks.create_mesh(mesh);
282 std::vector<Uint> node_connectivity, starting_indices;
289 comm_pattern.insert(
"gid",gids->array(),
false);
293 lss.
create(comm_pattern, 1u, node_connectivity, starting_indices);
296 lss.
matrix()->print(
"utest-ufem-buildsparsity_heat_matrix_3DHexaChannel.plt");
305 const Uint nb_segments = 5 ;
306 const Uint nb_nodes = nb_segments + 1;
319 boost::mpl::vector1<mesh::LagrangeP1::Line1D> allowed_elements;
333 lss_action->system_matrix +=
_A
337 << allocate_component<UFEM::BoundaryConditions>(
"BoundaryConditions")
338 << allocate_component<math::LSS::SolveLSS>(
"SolveLSS")
349 lss_action->options().set(
"regions", std::vector<URI>(1, mesh.
topology().
uri()));
354 lss.
matrix()->print(
"utest-ufem-buildsparsity_heat_matrix_1DHeat.plt");
359 BOOST_AUTO_TEST_SUITE_END()
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 >())
BOOST_AUTO_TEST_CASE(InitMPI)
static boost::proto::terminal< ZeroTag >::type _0
Placeholder for the zero matrix.
virtual mesh::Domain & create_domain(const std::string &name)
creates a domain in this model
Safe pointer to an object. This is the supported method for referring to components.
virtual physics::PhysModel & create_physics(const std::string &builder)
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)
UFEMBuildSparsityFixture()
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.
URI uri() const
Construct the full path.
Handle< LSS::Matrix > matrix()
Accessor to matrix.
static boost::proto::terminal< ElementQuadratureTag >::type element_quadrature
Use element_quadrature(expr1, expr2, ..., exprN) to evaluate a group of expressions.
const Field & coordinates() const
boost::shared_ptr< ElementsExpression< ExprT, ElementTypes > > elements_expression(ElementTypes, const ExprT &expr)
Manage a collection of UFEM solvers.
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.
boost::shared_ptr< List< Uint > > build_sparsity(const std::vector< Handle< Region > > ®ions, 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.
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
Basic Classes for Mesh applications used by COOLFluiD.
tuple model
Global confifuration.
Handle< Component > get_child(const std::string &name)
void init(int argc=0, char **args=0)
Top-level namespace for coolfluid.
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
unsigned int Uint
typedef for unsigned int
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
Handle< Component > handle()
Get a handle to the component.
boost::shared_ptr< NodesExpression< ExprT, boost::mpl::range_c< Uint, 1, 4 > > > nodes_expression(const ExprT &expr)
Dictionary & geometry_fields() const
static Comm & instance()
Return a reference to the current PE.
Base class for defining CF components.
void set(const std::string &pname, const boost::any &val)
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.
const Option & option(const std::string &pname) const
get a constant option from the list
Most basic kernel library.
Handle< common::Table< Real > > create_points(const Uint dimensions, const Uint nb_points)
Create the table that holds the points for the blocks.