v0.15.0
Loading...
Searching...
No Matches
NonlinearElasticElementInterface.hpp
Go to the documentation of this file.
1/** \file NonlinearElasticElementInterface.hpp
2* \example NonlinearElasticElementInterface.hpp
3
4 \brief Header file for NonlinearElasticElementInterface element implementation
5*/
6
7
8
9#ifndef __NONLINEARELEMENTINTERFACE_HPP__
10#define __NONLINEARELEMENTINTERFACE_HPP__
11
12/** \brief Set of functions declaring elements and setting operators
13 * for generic element interface
14 */
16
18 SmartPetscObj<DM> dM;
19 PetscBool isQuasiStatic;
20
21 PetscInt oRder;
23 BitRefLevel bIt;
24 boost::shared_ptr<NonlinearElasticElement> elasticElementPtr;
25 boost::shared_ptr<ElasticMaterials> elasticMaterialsPtr;
26 boost::shared_ptr<PostProcVolumeOnRefinedMesh> postProcMeshPtr;
27
30
32 string postion_field,
33 string mesh_posi_field_name = "MESH_NODE_POSITIONS",
34 bool is_displacement_field = true,
35 PetscBool is_quasi_static = PETSC_TRUE)
36 : mField(m_field), positionField(postion_field),
37 meshNodeField(mesh_posi_field_name),
38 isDisplacementField(is_displacement_field),
40 oRder = 1;
41 }
42
44
45 MoFEMErrorCode getCommandLineParameters() {
47 isQuasiStatic = PETSC_FALSE;
48 oRder = 2;
49 CHKERR PetscOptionsGetBool(PETSC_NULLPTR, "-is_quasi_static", &isQuasiStatic,
50 PETSC_NULLPTR);
51 CHKERR PetscOptionsGetInt(PETSC_NULLPTR, "-order", &oRder, PETSC_NULLPTR);
52
54 };
55
56 MoFEMErrorCode addElementFields() {
58 auto simple = mField.getInterface<Simple>();
61 3);
62 CHKERR simple->addBoundaryField(positionField, H1,
64 CHKERR simple->setFieldOrder(positionField, oRder);
65 }
68 3);
69 CHKERR simple->setFieldOrder(meshNodeField, 2);
70 }
71
73 };
74
75 MoFEMErrorCode createElements() {
77
78 elasticElementPtr = boost::make_shared<NonlinearElasticElement>(mField, 2);
79 elasticMaterialsPtr = boost::make_shared<ElasticMaterials>(mField, "HOOKE");
80 CHKERR elasticMaterialsPtr->setBlocks(elasticElementPtr->setOfBlocks);
81
82 CHKERR addHOOpsVol(meshNodeField, elasticElementPtr->getLoopFeRhs(), true,
83 false, false, false);
84 CHKERR addHOOpsVol(meshNodeField, elasticElementPtr->getLoopFeLhs(), true,
85 false, false, false);
86 CHKERR addHOOpsVol(meshNodeField, elasticElementPtr->getLoopFeEnergy(),
87 true, false, false, false);
88 CHKERR elasticElementPtr->addElement("ELASTIC", positionField,
89 meshNodeField, false);
90
91
93 };
94
95 MoFEMErrorCode setOperators() {
97 auto &pipeline_rhs = elasticElementPtr->feRhs.getOpPtrVector();
98 auto &pipeline_lhs = elasticElementPtr->feLhs.getOpPtrVector();
99
100 pipeline_rhs.push_back(new OpSetBc(positionField, true, mBoundaryMarker));
101 pipeline_lhs.push_back(new OpSetBc(positionField, true, mBoundaryMarker));
102
105
106 pipeline_rhs.push_back(new OpUnSetBc(positionField));
107 pipeline_lhs.push_back(new OpUnSetBc(positionField));
109 }
110
111 BitRefLevel getBitRefLevel() { return bIt; };
112 MoFEMErrorCode addElementsToDM(SmartPetscObj<DM> dm) {
114 this->dM = dm;
115 CHKERR DMMoFEMAddElement(dM, "ELASTIC");
116 mField.getInterface<Simple>()->getOtherFiniteElements().push_back(
117 "ELASTIC");
118
120 };
121
122 MoFEMErrorCode setupSolverJacobianSNES() {
124
125 CHKERR DMMoFEMSNESSetJacobian(
126 dM, "ELASTIC", &elasticElementPtr->getLoopFeLhs(), NULL, NULL);
127
129 };
130 MoFEMErrorCode setupSolverFunctionSNES() {
132 CHKERR DMMoFEMSNESSetFunction(dM, "ELASTIC",
133 &elasticElementPtr->getLoopFeRhs(),
134 PETSC_NULLPTR, PETSC_NULLPTR);
136 };
137
138 MoFEMErrorCode setupSolverJacobianTS(const TSType type) {
140 auto &method = elasticElementPtr->getLoopFeLhs();
141 switch (type) {
142 case IM:
143 CHKERR DMMoFEMTSSetIJacobian(dM, "ELASTIC", &method, &method, &method);
144 break;
145 case IM2:
146 CHKERR DMMoFEMTSSetI2Jacobian(dM, "ELASTIC", &method, &method, &method);
147 break;
148 case EX:
149 CHKERR DMMoFEMTSSetRHSJacobian(dM, "ELASTIC", &method, &method, &method);
150 break;
151 default:
152 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
153 "This TS is not yet implemented");
154 break;
155 }
157 };
158
159 MoFEMErrorCode setupSolverFunctionTS(const TSType type) {
161 auto &method = elasticElementPtr->getLoopFeRhs();
162 switch (type) {
163 case IM:
164 CHKERR DMMoFEMTSSetIFunction(dM, "ELASTIC", &method, &method, &method);
165 break;
166 case IM2:
167 CHKERR DMMoFEMTSSetI2Function(dM, "ELASTIC", &method, &method, &method);
168 break;
169 case EX:
170 CHKERR DMMoFEMTSSetRHSFunction(dM, "ELASTIC", &method, &method, &method);
171 break;
172 default:
173 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not implemented");
174 break;
175 }
176
178 };
179
180 MoFEMErrorCode updateElementVariables() { return 0; };
181 MoFEMErrorCode postProcessElement(int step) {
183
184 if (elasticElementPtr->setOfBlocks.empty())
186
187 if (!postProcMeshPtr) {
188 postProcMeshPtr = boost::make_shared<PostProcVolumeOnRefinedMesh>(mField);
189
190 if (mField.check_field("MESH_NODE_POSITIONS"))
191 CHKERR addHOOpsVol("MESH_NODE_POSITIONS", *postProcMeshPtr, true, false,
192 false, false);
193 CHKERR postProcMeshPtr->generateReferenceElementMesh();
194 CHKERR postProcMeshPtr->addFieldValuesPostProc(positionField);
195 CHKERR postProcMeshPtr->addFieldValuesPostProc(meshNodeField);
196 CHKERR postProcMeshPtr->addFieldValuesGradientPostProc(positionField);
197
198 for (auto &sit : elasticElementPtr->setOfBlocks) {
199 postProcMeshPtr->getOpPtrVector().push_back(new PostProcStress(
200 postProcMeshPtr->postProcMesh, postProcMeshPtr->mapGaussPts,
201 positionField, sit.second, postProcMeshPtr->commonData,
203 }
204 }
205
206 elasticElementPtr->getLoopFeEnergy().snes_ctx = SnesMethod::CTX_SNESNONE;
207 elasticElementPtr->getLoopFeEnergy().eNergy = 0;
208 // MOFEM_LOG("WORLD", Sev::inform) << "Loop energy\n";
209 CHKERR DMoFEMLoopFiniteElements(dM, "ELASTIC",
210 &elasticElementPtr->getLoopFeEnergy());
211
212 auto E = elasticElementPtr->getLoopFeEnergy().eNergy;
213 // Print elastic energy
214 MOFEM_LOG_C("WORLD", Sev::inform, "%d Time %3.2e Elastic energy %3.2e",
215 step, elasticElementPtr->getLoopFeRhs().ts_t, E);
216
217 CHKERR DMoFEMLoopFiniteElements(dM, "ELASTIC", postProcMeshPtr);
218 auto out_name = "out_vol_" + to_string(step) + ".h5m";
219
220 CHKERR postProcMeshPtr->writeFile(out_name);
221
223 };
224};
225
226#endif //__NONLINEARELEMENTINTERFACE_HPP__
#define MOFEM_LOG_C(channel, severity, format,...)
void simple(double P1[], double P2[], double P3[], double c[], const int N)
Definition acoustic.cpp:69
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base nme:nme847.
Definition definitions.h:60
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
@ H1
continuous field
Definition definitions.h:85
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
@ MOFEM_NOT_IMPLEMENTED
Definition definitions.h:32
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
virtual bool check_field(const std::string &name) const =0
check if field is in database
PetscBool is_quasi_static
Definition plastic.cpp:143
Set of functions declaring elements and setting operators for generic element interface.
Deprecated interface functions.
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface reference to pointer of interface.
Set of functions declaring elements and setting operators for generic element interface.
NonlinearElasticElementInterface(MoFEM::Interface &m_field, string postion_field, string mesh_posi_field_name="MESH_NODE_POSITIONS", bool is_displacement_field=true, PetscBool is_quasi_static=PETSC_TRUE)
boost::shared_ptr< ElasticMaterials > elasticMaterialsPtr
MoFEMErrorCode setupSolverJacobianTS(const TSType type)
MoFEMErrorCode setupSolverFunctionTS(const TSType type)
boost::shared_ptr< NonlinearElasticElement > elasticElementPtr
MoFEMErrorCode addElementsToDM(SmartPetscObj< DM > dm)
boost::shared_ptr< PostProcVolumeOnRefinedMesh > postProcMeshPtr