Commit bb32ab58 authored by GILLES Sebastien's avatar GILLES Sebastien

#1303 Make the non linear membrane test work correctly with the corrections given by Patrick.

parent 5df26ecf
......@@ -4166,7 +4166,7 @@
BE7E686A2065615100AA2FB3 /* Model.hpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.h; path = Model.hpp; sourceTree = "<group>"; };
BE7E686B2065615100AA2FB3 /* README */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = text; path = README; sourceTree = "<group>"; };
BE7E686C2065615100AA2FB3 /* Model.hxx */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.h; path = Model.hxx; sourceTree = "<group>"; };
BE7E686D2065615100AA2FB3 /* demo_input_nonlinear_membrane.lua */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = text; path = demo_input_nonlinear_membrane.lua; sourceTree = "<group>"; };
BE7E686D2065615100AA2FB3 /* demo.lua */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = text; path = demo.lua; sourceTree = "<group>"; };
BE7E686E2065615100AA2FB3 /* Model.cpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.cpp; path = Model.cpp; sourceTree = "<group>"; };
BE7E68702065615100AA2FB3 /* main.cpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.cpp; path = main.cpp; sourceTree = "<group>"; };
BE7E68752066728800AA2FB3 /* ExpectedResults.hpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.h; path = ExpectedResults.hpp; sourceTree = "<group>"; };
......@@ -8721,7 +8721,7 @@
children = (
BE7E68692065615100AA2FB3 /* CMakeLists.txt */,
BE7E686B2065615100AA2FB3 /* README */,
BE7E686D2065615100AA2FB3 /* demo_input_nonlinear_membrane.lua */,
BE7E686D2065615100AA2FB3 /* demo.lua */,
BE7E68702065615100AA2FB3 /* main.cpp */,
BE7E68682065615100AA2FB3 /* InputData.hpp */,
BE7E686E2065615100AA2FB3 /* Model.cpp */,
......@@ -11,82 +11,176 @@
#include "Test/Operators/VariationalInstances/NonlinearMembrane/ExpectedResults.hpp"
namespace MoReFEM
namespace MoReFEM::TestNS::NonLinearMembraneOperatorNS
{
namespace // anonymous
{
class FillMatrix
{
public:
FillMatrix(double factor,
expected_results_type<IsMatrixOrVector::matrix>& matrix);
void AddRow(std::vector<double>&& row) const;
private:
double factor_;
expected_results_type<IsMatrixOrVector::matrix>& matrix_;
};
void ApplyFactorToVector(double factor, expected_results_type<IsMatrixOrVector::vector>& vector);
// Matrix given as expected one is the local matrix computed in Matlab.
// In the test I want the global matrix (after assembling result); as it is a P1 matrix it is easy
// to apply the 'translation'.
void ConvertLocal2Global(expected_results_type<IsMatrixOrVector::matrix>& matrix);
// Same for vector.
void ConvertLocal2Global(expected_results_type<IsMatrixOrVector::vector>& vector);
} // namespace anonymous
namespace TestNS::NonLinearMembraneOperatorNS
expected_results_type<IsMatrixOrVector::matrix> GetExpectedMatrixP1P1()
{
// The values below have been computed independantly in Matlab by Dominique Chapelle.
expected_results_type<IsMatrixOrVector::matrix> ret;
FillMatrix fill_matrix(1.e9, ret);
fill_matrix.AddRow({ 3.1217, -2.8516, -0.2701, -0.5325, 0.4134, 0.1191, -0.3028, 0.2113, 0.0914 });
fill_matrix.AddRow({ -2.8516, 3.0204, -0.1688, 0.4339, -0.3198, -0.1141, 0.2183, -0.1612, -0.0571 });
fill_matrix.AddRow({ -0.2701, -0.1688, 0.4389, 0.0986, -0.0936, -0.0050, 0.0844, -0.0502, -0.0343 });
fill_matrix.AddRow({ -0.5325, 0.4339, 0.0986, 8.0441, -7.4416, -0.6025, 2.2840, -2.0490, -0.2350 });
fill_matrix.AddRow({ 0.4134, -0.3198, -0.0936, -7.4416, 8.0378, -0.5961, -2.0270, 2.0660, -0.0390 });
fill_matrix.AddRow({ 0.1191, -0.1141, -0.0050, -0.6025, -0.5961, 1.1986, -0.2571, -0.0170, 0.2740 });
fill_matrix.AddRow({ -0.3028, 0.2183, 0.0844, 2.2840, -2.0270, -0.2571, 4.1101, -3.7093, -0.4008 });
fill_matrix.AddRow({ 0.2113, -0.1612, -0.0502, -2.0490, 2.0660, -0.0170, -3.7093, 3.8471, -0.1378 });
fill_matrix.AddRow({ 0.0914, -0.0571, -0.0343, -0.2350, -0.0390, 0.2740, -0.4008, -0.1378, 0.5385 });
ConvertLocal2Global(ret);
expected_results_type<IsMatrixOrVector::matrix> GetExpectedMatrixP1P1()
return ret;
}
expected_results_type<IsMatrixOrVector::vector> GetExpectedVectorP1P1()
{
// The values below have been computed independantly in Matlab by Dominique Chapelle.
expected_results_type<IsMatrixOrVector::vector> ret
{
// The values below have been computed independantly in Matlab by Dominique Chapelle.
return expected_results_type<IsMatrixOrVector::matrix>
0.3342,
-0.2507,
-0.0835,
-3.7759,
4.2108,
-0.4349,
-1.6428,
1.7009,
-0.0582,
};
ApplyFactorToVector(1.e11, ret);
ConvertLocal2Global(ret);
return ret;
}
namespace // anonymous
{
FillMatrix::FillMatrix(double factor,
expected_results_type<IsMatrixOrVector::matrix>& matrix)
: factor_(factor),
matrix_(matrix)
{ }
void FillMatrix::AddRow(std::vector<double>&& row) const
{
for (auto& item : row)
item *= factor_;
matrix_.emplace_back(std::move(row));
}
void ApplyFactorToVector(double factor, expected_results_type<IsMatrixOrVector::vector>& vector)
{
for (auto& item : vector)
item *= factor;
}
void ConvertLocal2Global(expected_results_type<IsMatrixOrVector::matrix>& matrix)
{
assert(matrix.size() == 9ul);
std::vector<std::size_t> loc2glob_array { 0, 3, 6, 1, 4, 7, 2, 5, 8 };
expected_results_type<IsMatrixOrVector::matrix> reordered_matrix;
reordered_matrix.reserve(matrix.size());
std::vector<double> reordered_row;
reordered_row.reserve(matrix.size());
for (auto loc2glob : loc2glob_array)
{
const auto& original_row = matrix[loc2glob];
assert(original_row.size() == 9ul);
reordered_row.clear();
for (auto loc2glob_col : loc2glob_array)
{
3.1217173e+09, -5.3254715e+08, -3.0275440e+08, -2.8516195e+09, 4.1341198e+08,
2.1134163e+08, -2.7009781e+08, 1.1913516e+08, 9.1412770e+07
},
{
-5.3254715e+08, 8.0440892e+09, 2.2840439e+09, 4.3393421e+08, -7.4416177e+09,
-2.0490275e+09, 9.8612930e+07, -6.0247154e+08, -2.3501641e+08
},
{
-3.0275440e+08, 2.2840439e+09, 4.1100877e+09, 2.1831961e+08, -2.0269924e+09,
-3.7093376e+09, 8.4434793e+07, -2.5705148e+08, -4.0075010e+08
},
{
-2.8516195e+09, 4.3393421e+08, 2.1831961e+08, 3.0203942e+09, -3.1979658e+08,
-1.6118229e+08, -1.6877461e+08, -1.1413764e+08, -5.7137319e+07
},
{
4.1341198e+08, -7.4416177e+09, -2.0269924e+09, -3.1979658e+08, 8.0377623e+09,
2.0659864e+09, -9.3615406e+07, -5.9614462e+08, -3.8993961e+07
},
{
2.1134163e+08, -2.0490275e+09, -3.7093376e+09, -1.6118229e+08, 2.0659864e+09,
3.8470915e+09, -5.0159341e+07, -1.6958891e+07, -1.3775390e+08
},
{
-2.7009781e+08, 9.8612930e+07, 8.4434793e+07, -1.6877461e+08, -9.3615406e+07,
-5.0159341e+07, 4.3887242e+08, -4.9975241e+06, -3.4275451e+07
},
{
1.1913516e+08, -6.0247154e+08, -2.5705148e+08, -1.1413764e+08, -5.9614462e+08,
-1.6958891e+07, -4.9975241e+06, 1.1986162e+09, 2.7401037e+08
},
{
9.1412770e+07, -2.3501641e+08, -4.0075010e+08, -5.7137319e+07, -3.8993961e+07,
-1.3775390e+08, -3.4275451e+07, 2.7401037e+08, 5.3850401e+08
reordered_row.push_back(original_row[loc2glob_col]);
}
};
reordered_matrix.push_back(reordered_row);
}
matrix = reordered_matrix;
}
expected_results_type<IsMatrixOrVector::vector> GetExpectedVectorP1P1()
void ConvertLocal2Global(expected_results_type<IsMatrixOrVector::vector>& vector)
{
// The values below have been computed independantly in Matlab by Dominique Chapelle.
return expected_results_type<IsMatrixOrVector::vector>
assert(vector.size() == 9ul);
std::vector<std::size_t> loc2glob_array { 0, 3, 6, 1, 4, 7, 2, 5, 8 };
expected_results_type<IsMatrixOrVector::vector> reordered_vector;
reordered_vector.reserve(vector.size());
for (auto loc2glob : loc2glob_array)
{
3.3424261e+10,
-3.7758597e+11,
-1.6427512e+11,
-2.5069387e+10,
4.2107970e+11,
1.7009310e+11,
-8.3548746e+09,
-4.3493729e+10,
-5.8179833e+09
};
reordered_vector.push_back(vector[loc2glob]);
}
vector = reordered_vector;
}
} // namespace TestNS::NonLinearMembraneOperatorNS
} // namespace anonymous
} // namespace MoReFEM
} // namespace MoReFEM::TestNS::NonLinearMembraneOperatorNS
......@@ -86,7 +86,6 @@ namespace MoReFEM
domain_surface,
0.);
GVO::NonlinearMembrane stiffness_operator(felt_space_surface,
displacement_ptr,
displacement_ptr,
......@@ -109,13 +108,13 @@ namespace MoReFEM
__FILE__, __LINE__);
content[0] = 30.;
content[3] = -42.;
content[6] = -17.;
content[1] = 10.;
content[1] = -42.;
content[2] = -17.;
content[3] = 10.;
content[4] = 97.;
content[7] = 41.;
content[2] = 5.;
content[5] = -84.;
content[5] = 41.;
content[6] = 5.;
content[7] = -84.;
content[8] = -20.5;
}
......@@ -128,19 +127,17 @@ namespace MoReFEM
stiffness_operator.Assemble(std::make_tuple(std::ref(mat), std::ref(vec)),
previous_iteration);
matrix_tangent_stiffness.View(god_of_dof.GetMpi(), __FILE__, __LINE__);
/* BOOST_CHECK_NO_THROW */(CheckMatrix(god_of_dof,
matrix_tangent_stiffness,
GetExpectedMatrixP1P1(),
__FILE__, __LINE__,
1.e6)); // #1259 Obviously we need more precise reference values!
1.e6));
/* BOOST_CHECK_NO_THROW */(CheckVector(god_of_dof,
vector_stiffness_residual,
GetExpectedVectorP1P1(),
__FILE__, __LINE__,
1.e6)); // #1259 Obviously we need more precise reference values!
1.e7));
}
else if (do_assemble_into_matrix == assemble_into_matrix::yes
&& do_assemble_into_vector == assemble_into_vector::no)
......@@ -152,7 +149,7 @@ namespace MoReFEM
matrix_tangent_stiffness,
GetExpectedMatrixP1P1(),
__FILE__, __LINE__,
1.e6)); // #1259 Obviously we need more precise reference values!
1.e6));
}
else if (do_assemble_into_matrix == assemble_into_matrix::no
&& do_assemble_into_vector == assemble_into_vector::yes)
......@@ -164,7 +161,7 @@ namespace MoReFEM
vector_stiffness_residual,
GetExpectedVectorP1P1(),
__FILE__, __LINE__,
1.e6)); // #1259 Obviously we need more precise reference values!
1.e7));
}
......
......@@ -93,7 +93,7 @@ namespace // anonymous
return
environment.SubstituteValues("${MOREFEM_ROOT}/Sources/Test/Operators/"
"VariationalInstances/NonlinearMembrane/"
"demo_input_nonlinear_membrane.lua");
"demo.lua");
}
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment