v0.15.0
Loading...
Searching...
No Matches
Functions | Variables
damper_jacobian_test.cpp File Reference

Atom test testing calculation of element residual vectors and tangent matrices. More...

#include <BasicFiniteElements.hpp>

Go to the source code of this file.

Functions

int main (int argc, char *argv[])
 

Variables

static char help [] = "...\n\n"
 

Detailed Description

Atom test testing calculation of element residual vectors and tangent matrices.

Definition in file damper_jacobian_test.cpp.

Function Documentation

◆ main()

int main ( int  argc,
char *  argv[] 
)

Definition at line 17 of file damper_jacobian_test.cpp.

17 {
18
19 MoFEM::Core::Initialize(&argc, &argv, (char *)0, help);
20
21 try {
22
23 moab::Core mb_instance;
24 moab::Interface &moab = mb_instance;
25
26 {
27 PetscBool flg = PETSC_TRUE;
28 char mesh_file_name[255];
29 CHKERR PetscOptionsGetString(PETSC_NULLPTR, PETSC_NULLPTR, "-my_file",
30 mesh_file_name, 255, &flg);
31 ;
32 if (flg != PETSC_TRUE) {
33 SETERRQ(PETSC_COMM_SELF, 1, "*** ERROR -my_file (MESH FILE NEEDED)");
34 }
35 const char *option;
36 option = ""; //"PARALLEL=BCAST;";//;DEBUG_IO";
37 CHKERR moab.load_file(mesh_file_name, 0, option);
38 }
39
40 MoFEM::Core core(moab);
41 MoFEM::Interface &m_field = core;
42 BitRefLevel bit_level0;
43 bit_level0.set(0);
44 CHKERR m_field.getInterface<BitRefManager>()->setBitRefLevelByDim(
45 0, 3, bit_level0);
46 ;
47
48 // Define fields and finite elements
49 {
50
51 // Set approximation fields
52 {
53 // Seed all mesh entities to MoFEM database, those entities can be
54 // potentially used as finite elements or as entities which carry some
55 // approximation field.
56
57 bool check_if_spatial_field_exist =
58 m_field.check_field("SPATIAL_POSITION");
59 CHKERR m_field.add_field("SPATIAL_POSITION", H1,
60 AINSWORTH_LEGENDRE_BASE, 3, MB_TAG_SPARSE,
61 MF_ZERO);
62
63 // meshset consisting all entities in mesh
64 EntityHandle root_set = moab.get_root_set();
65 // add entities to field (root_mesh, i.e. on all mesh entities fields
66 // are approx.)
67 CHKERR m_field.add_ents_to_field_by_type(root_set, MBTET,
68 "SPATIAL_POSITION");
69
70 PetscBool flg;
71 PetscInt order;
72 CHKERR PetscOptionsGetInt(PETSC_NULLPTR, PETSC_NULLPTR, "-my_order", &order,
73 &flg);
74 ;
75 if (flg != PETSC_TRUE) {
76 order = 2;
77 }
78 if (order < 2) {
79 // SETERRQ()
80 }
81
82 CHKERR m_field.set_field_order(root_set, MBTET, "SPATIAL_POSITION",
83 order);
84 CHKERR m_field.set_field_order(root_set, MBTRI, "SPATIAL_POSITION",
85 order);
86 CHKERR m_field.set_field_order(root_set, MBEDGE, "SPATIAL_POSITION",
87 order);
88 CHKERR m_field.set_field_order(root_set, MBVERTEX, "SPATIAL_POSITION",
89 1);
90
91 CHKERR m_field.build_fields();
92
93 // Sett geometry approximation and initial spatial positions
94 // 10 node tets
95 if (!check_if_spatial_field_exist) {
96 Projection10NodeCoordsOnField ent_method_spatial(m_field,
97 "SPATIAL_POSITION");
98 CHKERR m_field.loop_dofs("SPATIAL_POSITION", ent_method_spatial);
99 }
100 }
101
102 // Set finite elements
103 {
104 CHKERR m_field.add_finite_element("DAMPER_FE", MF_ZERO);
105 CHKERR m_field.modify_finite_element_add_field_row("DAMPER_FE",
106 "SPATIAL_POSITION");
107 CHKERR m_field.modify_finite_element_add_field_col("DAMPER_FE",
108 "SPATIAL_POSITION");
109
110 CHKERR m_field.modify_finite_element_add_field_data("DAMPER_FE",
111 "SPATIAL_POSITION");
112
113 EntityHandle root_set = moab.get_root_set();
114 CHKERR m_field.add_ents_to_finite_element_by_type(root_set, MBTET,
115 "DAMPER_FE");
116
117 // build finite elemnts
118 CHKERR m_field.build_finite_elements();
119 // build adjacencies
120 CHKERR m_field.build_adjacencies(bit_level0);
121 }
122 }
123
124 // Create damper instance
125 KelvinVoigtDamper damper(m_field);
126 {
127
128 int id = 0;
129 KelvinVoigtDamper::BlockMaterialData &material_data =
130 damper.blockMaterialDataMap[id];
131 damper.constitutiveEquationMap.insert(
132 id,
133 new KelvinVoigtDamper::ConstitutiveEquation<adouble>(material_data));
134
135 // Set material parameters
136 CHKERR moab.get_entities_by_type(0, MBTET, material_data.tEts);
137 material_data.gBeta = 1;
138 material_data.vBeta = 0.3;
139
140 KelvinVoigtDamper::CommonData &common_data = damper.commonData;
141 common_data.spatialPositionName = "SPATIAL_POSITION";
142 common_data.spatialPositionNameDot = "DOT_SPATIAL_POSITION";
143
144 for (auto &&fe_ptr : {&damper.feRhs, &damper.feLhs}) {
145 // fe_ptr->getOpPtrVector().push_back(
146 // new OpCalculateVectorFieldValuesDot<3>(
147 // common_data.spatialPositionName,
148 // common_data.dataAtGaussTmpPtr));
149 fe_ptr->getOpPtrVector().push_back(
150 new OpCalculateVectorFieldGradient<3, 3>(
151 common_data.spatialPositionName,
152 common_data.gradDataAtGaussTmpPtr));
153 fe_ptr->getOpPtrVector().push_back(
154 new KelvinVoigtDamper::OpGetDataAtGaussPts(
155 common_data.spatialPositionName, common_data, false, true, false));
156 fe_ptr->getOpPtrVector().push_back(
157 new OpCalculateVectorFieldGradientDot<3, 3>(
158 common_data.spatialPositionName,
159 common_data.gradDataAtGaussTmpPtr));
160 fe_ptr->getOpPtrVector().push_back(
161 new KelvinVoigtDamper::OpGetDataAtGaussPts(
162 common_data.spatialPositionName, common_data, false, true, true));
163 }
164
165 // attach tags for each recorder
166 std::vector<int> tags;
167 tags.push_back(1);
168
169 KelvinVoigtDamper::ConstitutiveEquation<adouble> &ce =
170 damper.constitutiveEquationMap.at(id);
171
172 // Right hand side operators
173 damper.feRhs.getOpPtrVector().push_back(new KelvinVoigtDamper::OpJacobian(
174 "SPATIAL_POSITION", tags, ce, damper.commonData, true, false));
175 damper.feRhs.getOpPtrVector().push_back(
176 new KelvinVoigtDamper::OpRhsStress(damper.commonData));
177
178 // Left hand side operators
179 damper.feLhs.getOpPtrVector().push_back(new KelvinVoigtDamper::OpJacobian(
180 "SPATIAL_POSITION", tags, ce, damper.commonData, false, true));
181 damper.feLhs.getOpPtrVector().push_back(
182 new KelvinVoigtDamper::OpLhsdxdx(damper.commonData));
183 }
184
185 // Create dm instance
186 DM dm;
187 DMType dm_name = "DMDAMPER";
188 {
189 CHKERR DMRegister_MoFEM(dm_name);
190 CHKERR DMCreate(PETSC_COMM_WORLD, &dm);
191 CHKERR DMSetType(dm, dm_name);
192 CHKERR DMMoFEMCreateMoFEM(dm, &m_field, dm_name, bit_level0);
193 CHKERR DMSetFromOptions(dm);
194 // add elements to dm
195 CHKERR DMMoFEMAddElement(dm, "DAMPER_FE");
196 CHKERR DMSetUp(dm);
197 }
198
199 // Make calculations
200 Mat M;
201 Vec F, U_t;
202 {
203 CHKERR DMCreateGlobalVector_MoFEM(dm, &F);
204 CHKERR DMCreateGlobalVector_MoFEM(dm, &U_t);
205 CHKERR DMCreateMatrix_MoFEM(dm, &M);
206 CHKERR VecZeroEntries(F);
207 CHKERR VecGhostUpdateBegin(F, INSERT_VALUES, SCATTER_FORWARD);
208 CHKERR VecGhostUpdateEnd(F, INSERT_VALUES, SCATTER_FORWARD);
209 CHKERR VecZeroEntries(U_t);
210 CHKERR VecGhostUpdateBegin(U_t, INSERT_VALUES, SCATTER_FORWARD);
211 CHKERR VecGhostUpdateEnd(U_t, INSERT_VALUES, SCATTER_FORWARD);
212 CHKERR MatZeroEntries(M);
213 damper.feRhs.ts_F = F; // Set right hand side vector manually
214 damper.feRhs.ts_u_t = U_t;
215 damper.feRhs.ts_ctx = TSMethod::CTX_TSSETIFUNCTION;
216 CHKERR DMoFEMLoopFiniteElements(dm, "DAMPER_FE", &damper.feRhs);
217 damper.feLhs.ts_B = M; // Set matrix M
218 damper.feLhs.ts_a = 1.0; // Set time step parameter
219 damper.feLhs.ts_u_t = U_t;
220 CHKERR DMoFEMLoopFiniteElements(dm, "DAMPER_FE", &damper.feLhs);
221 CHKERR VecAssemblyBegin(F);
222 CHKERR VecAssemblyEnd(F);
223 CHKERR MatAssemblyBegin(M, MAT_FINAL_ASSEMBLY);
224 CHKERR MatAssemblyEnd(M, MAT_FINAL_ASSEMBLY);
225 }
226
227 // See results
228 {
229
230 PetscViewer viewer;
231 CHKERR PetscViewerASCIIOpen(PETSC_COMM_WORLD, "damper_jacobian_test.txt",
232 &viewer);
233 CHKERR PetscViewerDestroy(&viewer);
234 // MatView(M,PETSC_VIEWER_DRAW_WORLD);
235 // std::string wait;
236 // std::cin >> wait;
237 }
238
239 // Clean and destroy
240 {
241 CHKERR VecDestroy(&F);
242 CHKERR VecDestroy(&U_t);
243 CHKERR MatDestroy(&M);
244 CHKERR DMDestroy(&dm);
245 }
246 }
247 CATCH_ERRORS;
248
249 MoFEM::Core::Finalize();
250
251 return 0;
252}
help
static char help[]
Definition damper_jacobian_test.cpp:15
CATCH_ERRORS
#define CATCH_ERRORS
Catch errors.
Definition definitions.h:385
MF_ZERO
@ MF_ZERO
Definition definitions.h:111
AINSWORTH_LEGENDRE_BASE
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
Definition definitions.h:60
H1
@ H1
continuous field
Definition definitions.h:85
CHKERR
#define CHKERR
Inline error check.
Definition definitions.h:551
order
constexpr int order
Definition dg_projection.cpp:19
F
@ F
Definition free_surface.cpp:624
MoFEM::DMMoFEMAddElement
PetscErrorCode DMMoFEMAddElement(DM dm, std::string fe_name)
add element to dm
Definition DMMoFEM.cpp:488
MoFEM::DMMoFEMCreateMoFEM
PetscErrorCode DMMoFEMCreateMoFEM(DM dm, MoFEM::Interface *m_field_ptr, const char problem_name[], const MoFEM::BitRefLevel bit_level, const MoFEM::BitRefLevel bit_mask=MoFEM::BitRefLevel().set())
Must be called by user to set MoFEM data structures.
Definition DMMoFEM.cpp:114
MoFEM::DMCreateMatrix_MoFEM
PetscErrorCode DMCreateMatrix_MoFEM(DM dm, Mat *M)
Definition DMMoFEM.cpp:1187
MoFEM::DMRegister_MoFEM
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
Definition DMMoFEM.cpp:43
MoFEM::DMoFEMLoopFiniteElements
PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[], MoFEM::FEMethod *method, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr())
Executes FEMethod for finite elements in DM.
Definition DMMoFEM.cpp:576
MoFEM::DMCreateGlobalVector_MoFEM
PetscErrorCode DMCreateGlobalVector_MoFEM(DM dm, Vec *g)
DMShellSetCreateGlobalVector.
Definition DMMoFEM.cpp:1157
MoFEM::CoreInterface::add_finite_element
virtual MoFEMErrorCode add_finite_element(const std::string &fe_name, enum MoFEMTypes bh=MF_EXCL, int verb=DEFAULT_VERBOSITY)=0
add finite element
MoFEM::CoreInterface::build_finite_elements
virtual MoFEMErrorCode build_finite_elements(int verb=DEFAULT_VERBOSITY)=0
Build finite elements.
MoFEM::CoreInterface::modify_finite_element_add_field_col
virtual MoFEMErrorCode modify_finite_element_add_field_col(const std::string &fe_name, const std::string name_row)=0
set field col which finite element use
MoFEM::CoreInterface::add_ents_to_finite_element_by_type
virtual MoFEMErrorCode add_ents_to_finite_element_by_type(const EntityHandle entities, const EntityType type, const std::string name, const bool recursive=true)=0
add entities to finite element
MoFEM::CoreInterface::modify_finite_element_add_field_row
virtual MoFEMErrorCode modify_finite_element_add_field_row(const std::string &fe_name, const std::string name_row)=0
set field row which finite element use
MoFEM::CoreInterface::modify_finite_element_add_field_data
virtual MoFEMErrorCode modify_finite_element_add_field_data(const std::string &fe_name, const std::string name_field)=0
set finite element field data
MoFEM::CoreInterface::build_fields
virtual MoFEMErrorCode build_fields(int verb=DEFAULT_VERBOSITY)=0
MoFEM::CoreInterface::set_field_order
virtual MoFEMErrorCode set_field_order(const EntityHandle meshset, const EntityType type, const std::string &name, const ApproximationOrder order, int verb=DEFAULT_VERBOSITY)=0
Set order approximation of the entities in the field.
MoFEM::CoreInterface::add_ents_to_field_by_type
virtual MoFEMErrorCode add_ents_to_field_by_type(const Range &ents, const EntityType type, const std::string &name, int verb=DEFAULT_VERBOSITY)=0
Add entities to field meshset.
MoFEM::CoreInterface::check_field
virtual bool check_field(const std::string &name) const =0
check if field is in database
MoFEM::CoreInterface::loop_dofs
virtual MoFEMErrorCode loop_dofs(const Problem *problem_ptr, const std::string &field_name, RowColData rc, DofMethod &method, int lower_rank, int upper_rank, int verb=DEFAULT_VERBOSITY)=0
Make a loop over dofs.
mesh_file_name
char mesh_file_name[255]
Definition mesh_smoothing.cpp:22
EigenMatrix::Vec
const FTensor::Tensor2< T, Dim, Dim > Vec
Definition MatrixFunction.hpp:64
MoFEM::Types::BitRefLevel
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition Types.hpp:40
MoFEM::PetscOptionsGetInt
PetscErrorCode PetscOptionsGetInt(PetscOptions *, const char pre[], const char name[], PetscInt *ivalue, PetscBool *set)
Definition DeprecatedPetsc.hpp:142
MoFEM::PetscOptionsGetString
PetscErrorCode PetscOptionsGetString(PetscOptions *, const char pre[], const char name[], char str[], size_t size, PetscBool *set)
Definition DeprecatedPetsc.hpp:172
PlasticOps::M
FTensor::Index< 'M', 3 > M
Definition PlasticOps.hpp:117
KelvinVoigtDamper::BlockMaterialData::vBeta
double vBeta
Poisson ration spring alpha.
Definition KelvinVoigtDamper.hpp:31
KelvinVoigtDamper::BlockMaterialData::gBeta
double gBeta
Sheer modulus spring alpha.
Definition KelvinVoigtDamper.hpp:32
KelvinVoigtDamper::CommonData
Common data for nonlinear_elastic_elem model.
Definition KelvinVoigtDamper.hpp:169
KelvinVoigtDamper::CommonData::gradDataAtGaussTmpPtr
boost::shared_ptr< MatrixDouble > gradDataAtGaussTmpPtr
Definition KelvinVoigtDamper.hpp:179
KelvinVoigtDamper::OpRhsStress
Assemble internal force vector.
Definition KelvinVoigtDamper.hpp:564
KelvinVoigtDamper
Implementation of Kelvin Voigt Damper.
Definition KelvinVoigtDamper.hpp:20
MoFEM::BitRefManager
Managing BitRefLevels.
Definition BitRefManager.hpp:21
MoFEM::CoreInterface::build_adjacencies
virtual MoFEMErrorCode build_adjacencies(const Range &ents, int verb=DEFAULT_VERBOSITY)=0
build adjacencies
MoFEM::CoreInterface::add_field
virtual MoFEMErrorCode add_field(const std::string name, const FieldSpace space, const FieldApproximationBase base, const FieldCoefficientsNumber nb_of_coefficients, const TagType tag_type=MB_TAG_SPARSE, const enum MoFEMTypes bh=MF_EXCL, int verb=DEFAULT_VERBOSITY)=0
Add field.
MoFEM::CoreTmp< 0 >
Core (interface) class.
Definition Core.hpp:82
MoFEM::CoreTmp< 0 >::Initialize
static MoFEMErrorCode Initialize(int *argc, char ***args, const char file[], const char help[])
Initializes the MoFEM database PETSc, MOAB and MPI.
Definition Core.cpp:72
MoFEM::CoreTmp< 0 >::Finalize
static MoFEMErrorCode Finalize()
Checks for options to be called at the conclusion of the program.
Definition Core.cpp:118
MoFEM::DeprecatedCoreInterface
Deprecated interface functions.
Definition DeprecatedCoreInterface.hpp:16
MoFEM::OpCalculateVectorFieldGradientDot
Get field gradients time derivative at integration pts for scalar field rank 0, i....
Definition UserDataOperators.hpp:1759
MoFEM::OpCalculateVectorFieldGradient
Get field gradients at integration pts for scalar field rank 0, i.e. vector field.
Definition UserDataOperators.hpp:1745
MoFEM::Projection10NodeCoordsOnField
Projection of edge entities with one mid-node on hierarchical basis.
Definition Projection10NodeCoordsOnField.hpp:24
MoFEM::TSMethod::CTX_TSSETIFUNCTION
@ CTX_TSSETIFUNCTION
Setting up implicit function.
Definition LoopMethods.hpp:281
MoFEM::UnknownInterface::getInterface
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface reference to pointer of interface.
Definition UnknownInterface.hpp:93

Variable Documentation

◆ help

char help[] = "...\n\n"
static

Definition at line 15 of file damper_jacobian_test.cpp.

Generated by Doxygen 1.9.8 and hosted at University of Glasgow