v0.14.0
Loading...
Searching...
No Matches
NonlinearElasticElementInterface.hpp

Header file for NonlinearElasticElementInterface element implementation.

Header file for NonlinearElasticElementInterface element implementation

/** \file NonlinearElasticElementInterface.hpp
* \example NonlinearElasticElementInterface.hpp
\brief Header file for NonlinearElasticElementInterface element implementation
*/
#ifndef __NONLINEARELEMENTINTERFACE_HPP__
#define __NONLINEARELEMENTINTERFACE_HPP__
/** \brief Set of functions declaring elements and setting operators
* for generic element interface
*/
SmartPetscObj<DM> dM;
PetscBool isQuasiStatic;
PetscInt oRder;
BitRefLevel bIt;
boost::shared_ptr<NonlinearElasticElement> elasticElementPtr;
boost::shared_ptr<ElasticMaterials> elasticMaterialsPtr;
boost::shared_ptr<PostProcVolumeOnRefinedMesh> postProcMeshPtr;
string positionField;
string meshNodeField;
string postion_field,
string mesh_posi_field_name = "MESH_NODE_POSITIONS",
bool is_displacement_field = true,
PetscBool is_quasi_static = PETSC_TRUE)
: mField(m_field), positionField(postion_field),
meshNodeField(mesh_posi_field_name),
isDisplacementField(is_displacement_field),
oRder = 1;
}
MoFEMErrorCode getCommandLineParameters() {
isQuasiStatic = PETSC_FALSE;
oRder = 2;
CHKERR PetscOptionsGetBool(PETSC_NULL, "-is_quasi_static", &isQuasiStatic,
PETSC_NULL);
CHKERR PetscOptionsGetInt(PETSC_NULL, "-order", &oRder, PETSC_NULL);
};
MoFEMErrorCode addElementFields() {
auto simple = mField.getInterface<Simple>();
3);
CHKERR simple->addBoundaryField(positionField, H1,
CHKERR simple->setFieldOrder(positionField, oRder);
}
3);
CHKERR simple->setFieldOrder(meshNodeField, 2);
}
};
MoFEMErrorCode createElements() {
elasticElementPtr = boost::make_shared<NonlinearElasticElement>(mField, 2);
elasticMaterialsPtr = boost::make_shared<ElasticMaterials>(mField, "HOOKE");
CHKERR elasticMaterialsPtr->setBlocks(elasticElementPtr->setOfBlocks);
CHKERR addHOOpsVol(meshNodeField, elasticElementPtr->getLoopFeRhs(), true,
false, false, false);
CHKERR addHOOpsVol(meshNodeField, elasticElementPtr->getLoopFeLhs(), true,
false, false, false);
CHKERR addHOOpsVol(meshNodeField, elasticElementPtr->getLoopFeEnergy(),
true, false, false, false);
CHKERR elasticElementPtr->addElement("ELASTIC", positionField,
meshNodeField, false);
};
MoFEMErrorCode setOperators() {
auto &pipeline_rhs = elasticElementPtr->feRhs.getOpPtrVector();
auto &pipeline_lhs = elasticElementPtr->feLhs.getOpPtrVector();
pipeline_rhs.push_back(new OpSetBc(positionField, true, mBoundaryMarker));
pipeline_lhs.push_back(new OpSetBc(positionField, true, mBoundaryMarker));
pipeline_rhs.push_back(new OpUnSetBc(positionField));
pipeline_lhs.push_back(new OpUnSetBc(positionField));
}
BitRefLevel getBitRefLevel() { return bIt; };
MoFEMErrorCode addElementsToDM(SmartPetscObj<DM> dm) {
this->dM = dm;
CHKERR DMMoFEMAddElement(dM, "ELASTIC");
mField.getInterface<Simple>()->getOtherFiniteElements().push_back(
"ELASTIC");
};
MoFEMErrorCode setupSolverJacobianSNES() {
CHKERR DMMoFEMSNESSetJacobian(
dM, "ELASTIC", &elasticElementPtr->getLoopFeLhs(), NULL, NULL);
};
MoFEMErrorCode setupSolverFunctionSNES() {
CHKERR DMMoFEMSNESSetFunction(dM, "ELASTIC",
&elasticElementPtr->getLoopFeRhs(),
PETSC_NULL, PETSC_NULL);
};
MoFEMErrorCode setupSolverJacobianTS(const TSType type) {
auto &method = elasticElementPtr->getLoopFeLhs();
switch (type) {
case IM:
CHKERR DMMoFEMTSSetIJacobian(dM, "ELASTIC", &method, &method, &method);
break;
case IM2:
CHKERR DMMoFEMTSSetI2Jacobian(dM, "ELASTIC", &method, &method, &method);
break;
case EX:
CHKERR DMMoFEMTSSetRHSJacobian(dM, "ELASTIC", &method, &method, &method);
break;
default:
SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
"This TS is not yet implemented");
break;
}
};
MoFEMErrorCode setupSolverFunctionTS(const TSType type) {
auto &method = elasticElementPtr->getLoopFeRhs();
switch (type) {
case IM:
CHKERR DMMoFEMTSSetIFunction(dM, "ELASTIC", &method, &method, &method);
break;
case IM2:
CHKERR DMMoFEMTSSetI2Function(dM, "ELASTIC", &method, &method, &method);
break;
case EX:
CHKERR DMMoFEMTSSetRHSFunction(dM, "ELASTIC", &method, &method, &method);
break;
default:
SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not implemented");
break;
}
};
MoFEMErrorCode updateElementVariables() { return 0; };
MoFEMErrorCode postProcessElement(int step) {
if (elasticElementPtr->setOfBlocks.empty())
postProcMeshPtr = boost::make_shared<PostProcVolumeOnRefinedMesh>(mField);
if (mField.check_field("MESH_NODE_POSITIONS"))
CHKERR addHOOpsVol("MESH_NODE_POSITIONS", *postProcMeshPtr, true, false,
false, false);
CHKERR postProcMeshPtr->generateReferenceElementMesh();
CHKERR postProcMeshPtr->addFieldValuesPostProc(positionField);
CHKERR postProcMeshPtr->addFieldValuesPostProc(meshNodeField);
CHKERR postProcMeshPtr->addFieldValuesGradientPostProc(positionField);
for (auto &sit : elasticElementPtr->setOfBlocks) {
postProcMeshPtr->getOpPtrVector().push_back(new PostProcStress(
postProcMeshPtr->postProcMesh, postProcMeshPtr->mapGaussPts,
positionField, sit.second, postProcMeshPtr->commonData,
}
}
elasticElementPtr->getLoopFeEnergy().snes_ctx = SnesMethod::CTX_SNESNONE;
elasticElementPtr->getLoopFeEnergy().eNergy = 0;
// MOFEM_LOG("WORLD", Sev::inform) << "Loop energy\n";
CHKERR DMoFEMLoopFiniteElements(dM, "ELASTIC",
&elasticElementPtr->getLoopFeEnergy());
auto E = elasticElementPtr->getLoopFeEnergy().eNergy;
// Print elastic energy
MOFEM_LOG_C("WORLD", Sev::inform, "%d Time %3.2e Elastic energy %3.2e",
step, elasticElementPtr->getLoopFeRhs().ts_t, E);
CHKERR DMoFEMLoopFiniteElements(dM, "ELASTIC", postProcMeshPtr);
auto out_name = "out_vol_" + to_string(step) + ".h5m";
CHKERR postProcMeshPtr->writeFile(out_name);
};
};
#endif //__NONLINEARELEMENTINTERFACE_HPP__
#define MOFEM_LOG_C(channel, severity, format,...)
Definition: LogManager.hpp:311
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 .
Definition: definitions.h:60
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:447
@ 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 ...
Definition: definitions.h:346
@ MOFEM_NOT_IMPLEMENTED
Definition: definitions.h:32
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
#define CHKERR
Inline error check.
Definition: definitions.h:535
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:440
virtual bool check_field(const std::string &name) const =0
check if field is in database
PetscBool is_quasi_static
Definition: plastic.cpp:186
Set of functions declaring elements and setting operators for generic element interface.
Deprecated interface functions.
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.
Set of functions declaring elements and setting operators for generic element interface.
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