8 #ifndef OPENFPM_NUMERICS_SRC_FINITEDIFFERENCE_EQ_UNIT_TEST_3D_HPP_
9 #define OPENFPM_NUMERICS_SRC_FINITEDIFFERENCE_EQ_UNIT_TEST_3D_HPP_
12 #include "Laplacian.hpp"
13 #include "FiniteDifference/eq.hpp"
14 #include "FiniteDifference/sum.hpp"
15 #include "FiniteDifference/mul.hpp"
16 #include "Grid/grid_dist_id.hpp"
17 #include "Decomposition/CartDecomposition.hpp"
18 #include "Vector/Vector.hpp"
19 #include "Solvers/umfpack_solver.hpp"
20 #include "data_type/aggregate.hpp"
22 BOOST_AUTO_TEST_SUITE( eq_test_suite_3d )
28 static const unsigned int dims = 3;
30 static const unsigned int nvar = 4;
33 static const bool boundary[];
54 static const unsigned int dims = 3;
56 static const unsigned int nvar = 4;
89 static float val() {
return 1.0;}
92 template<
typename solver_type,
typename l
id_nn_3d>
void lid_driven_cavity_3d()
94 #include "Equations/stoke_flow_eq_3d.hpp"
107 long int sz[] = {36,12,12};
109 szu[0] = (size_t)sz[0];
110 szu[1] = (size_t)sz[1];
111 szu[2] = (size_t)sz[2];
116 constexpr
int velocity = 0;
117 constexpr
int pressure = 1;
130 fd.impose(
ic_eq(),0.0, EQ_4, {0,0,0},{sz[0]-2,sz[1]-2,sz[2]-2},
true);
131 fd.impose(
Prs(), 0.0, EQ_4, {0,0,0},{0,0,0});
132 fd.impose(
vx_eq(),0.0, EQ_1, {1,0},{sz[0]-2,sz[1]-2,sz[2]-2});
133 fd.impose(
vy_eq(),0.0, EQ_2, {0,1},{sz[0]-2,sz[1]-2,sz[2]-2});
134 fd.impose(
vz_eq(),0.0, EQ_3, {0,0,1},{sz[0]-2,sz[1]-2,sz[2]-2});
138 fd.impose(
v_x(),0.0, EQ_1, {0,0,0}, {0,sz[1]-2,sz[2]-2});
139 fd.impose(
v_x(),0.0, EQ_1, {sz[0]-1,0,0},{sz[0]-1,sz[1]-2,sz[2]-2});
142 fd.impose(
avg_y_vx_f(),0.0, EQ_1, {0,-1,0}, {sz[0]-1,-1,sz[2]-2});
143 fd.impose(
avg_y_vx(),0.0, EQ_1, {0,sz[1]-1,0},{sz[0]-1,sz[1]-1,sz[2]-2});
146 fd.impose(
avg_z_vx_f(),0.0, EQ_1, {0,-1,-1}, {sz[0]-1,sz[1]-1,-1});
147 fd.impose(
avg_z_vx(),0.0, EQ_1, {0,-1,sz[2]-1},{sz[0]-1,sz[1]-1,sz[2]-1});
151 fd.impose(
avg_x_vy_f(),0.0, EQ_2, {-1,0,0}, {-1,sz[1]-1,sz[2]-2});
152 fd.impose(
avg_x_vy(),1.0, EQ_2, {sz[0]-1,0,0},{sz[0]-1,sz[1]-1,sz[2]-2});
155 fd.impose(
v_y(), 0.0, EQ_2, {0,0,0}, {sz[0]-2,0,sz[2]-2});
156 fd.impose(
v_y(), 0.0, EQ_2, {0,sz[1]-1,0},{sz[0]-2,sz[1]-1,sz[2]-2});
159 fd.impose(
avg_z_vy(),0.0, EQ_2, {-1,0,sz[2]-1}, {sz[0]-1,sz[1]-1,sz[2]-1});
160 fd.impose(
avg_z_vy_f(),0.0, EQ_2, {-1,0,-1}, {sz[0]-1,sz[1]-1,-1});
164 fd.impose(
avg_x_vz_f(),0.0, EQ_3, {-1,0,0}, {-1,sz[1]-2,sz[2]-1});
165 fd.impose(
avg_x_vz(),1.0, EQ_3, {sz[0]-1,0,0},{sz[0]-1,sz[1]-2,sz[2]-1});
168 fd.impose(
avg_y_vz(),0.0, EQ_3, {-1,sz[1]-1,0},{sz[0]-1,sz[1]-1,sz[2]-1});
169 fd.impose(
avg_y_vz_f(),0.0, EQ_3, {-1,-1,0}, {sz[0]-1,-1,sz[2]-1});
172 fd.impose(
v_z(),0.0, EQ_3, {0,0,0}, {sz[0]-2,sz[1]-2,0});
173 fd.impose(
v_z(),0.0, EQ_3, {0,0,sz[2]-1},{sz[0]-2,sz[1]-2,sz[2]-1});
178 fd.impose(
Prs(), 0.0, EQ_4, {-1,-1,-1},{-1,sz[1]-1,sz[2]-1});
179 fd.impose(
Prs(), 0.0, EQ_4, {sz[0]-1,-1,-1},{sz[0]-1,sz[1]-1,sz[2]-1});
182 fd.impose(
Prs(), 0.0, EQ_4, {0,sz[1]-1,-1}, {sz[0]-2,sz[1]-1,sz[2]-1});
183 fd.impose(
Prs(), 0.0, EQ_4, {0,-1 ,-1}, {sz[0]-2,-1, sz[2]-1});
186 fd.impose(
Prs(), 0.0, EQ_4, {0,0,sz[2]-1}, {sz[0]-2,sz[1]-2,sz[2]-1});
187 fd.impose(
Prs(), 0.0, EQ_4, {0,0,-1}, {sz[0]-2,sz[1]-2,-1});
190 fd.impose(
v_x(), 0.0, EQ_1, {-1,-1,-1},{-1,sz[1]-1,sz[2]-1});
191 fd.impose(
v_y(), 0.0, EQ_2, {-1,-1,-1},{sz[0]-1,-1,sz[2]-1});
192 fd.impose(
v_z(), 0.0, EQ_3, {-1,-1,-1},{sz[0]-1,sz[1]-1,-1});
195 auto x = solver.try_solve(fd.getA(),fd.getB());
198 fd.template copy<velocity,pressure>(x,{0,0},{sz[0]-1,sz[1]-1,sz[2]-1},g_dist);
200 std::string s = std::string(demangle(
typeid(solver_type).name()));
203 g_dist.write(s +
"lid_driven_cavity_3d_p" + std::to_string(v_cl.
getProcessingUnits()) +
"_grid");
207 std::string file1 = std::string(
"test/") + s +
"lid_driven_cavity_3d_p" + std::to_string(v_cl.
getProcessingUnits()) +
"_grid_" + std::to_string(v_cl.
getProcessUnitID()) +
"_test_osx.vtk";
214 std::string file1 = std::string(
"test/") + s +
"lid_driven_cavity_3d_p" + std::to_string(v_cl.
getProcessingUnits()) +
"_grid_" + std::to_string(v_cl.
getProcessUnitID()) +
"_test_GCC5.vtk";
217 #elif __GNUC__ == 6 || __GNUC__ == 7
219 std::string file1 = std::string(
"test/") + s +
"lid_driven_cavity_3d_p" + std::to_string(v_cl.
getProcessingUnits()) +
"_grid_" + std::to_string(v_cl.
getProcessUnitID()) +
"_test_GCC6.vtk";
224 std::string file1 = std::string(
"test/") + s +
"lid_driven_cavity_3d_p" + std::to_string(v_cl.
getProcessingUnits()) +
"_grid_" + std::to_string(v_cl.
getProcessUnitID()) +
"_test_GCC4.vtk";
231 std::cout <<
"File1: " << file1 << std::endl;
232 std::cout <<
"File2: " << file2 << std::endl;
237 bool test = compare(file1,file2);
238 BOOST_REQUIRE_EQUAL(test,
true);
245 BOOST_AUTO_TEST_CASE(lid_driven_cavity)
247 #if defined(HAVE_EIGEN) && defined(HAVE_SUITESPARSE)
255 BOOST_AUTO_TEST_SUITE_END()
grid_dist_id< 3, float, aggregate< float[3], float >, CartDecomposition< 3, float > > b_grid
type of base grid
[Definition of the system]
SparseMatrix< double, int > SparseMatrix_type
type of SparseMatrix for the linear solver
static const bool boundary[]
boundary at X and Y
SparseMatrix< double, int, PETSC_BASE > SparseMatrix_type
type of SparseMatrix for the linear solver
size_t getProcessUnitID()
Get the process unit id.
Sparse Matrix implementation stub object when OpenFPM is compiled with no linear algebra support...
Class that contain Padding information on each direction positive and Negative direction.
grid_dist_id< 3, float, aggregate< float[3], float >, CartDecomposition< 3, float > > b_grid
type of base grid
Implementation of VCluster class.
Specify the general caratteristic of system to solve.
This class decompose a space into sub-sub-domains and distribute them across processors.
Sparse Matrix implementation.
This is a distributed grid.
static const unsigned int nvar
number of fields in the system
void const_field
define that eta is a constant field
Vector< double > Vector_type
type of Vector for the linear solver
static const unsigned int dims
dimensionaly of the equation ( 3D problem ...)
float stype
type of space float, double, ...
It model an expression expr1 + ... exprn.
Vector< double, PETSC_BASE > Vector_type
type of Vector for the linear solver
static float val()
therutn the value of the constant
static const bool boundary[]
boundary at X and Y
size_t getProcessingUnits()
Get the total number of processors.
float stype
type of space float, double, ...