v0.13.1
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
*/
/* This file is part of MoFEM.
* MoFEM is free software: you can redistribute it and/or modify it under
* the terms of the GNU Lesser General Public License as published by the
* Free Software Foundation, either version 3 of the License, or (at your
* option) any later version.
*
* MoFEM is distributed in the hope that it will be useful, but WITHOUT
* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
* FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
* License for more details.
*
* You should have received a copy of the GNU Lesser General Public
* License along with MoFEM. If not, see <http://www.gnu.org/licenses/>. */
#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;
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;
}
isQuasiStatic = PETSC_FALSE;
oRder = 2;
CHKERR PetscOptionsGetBool(PETSC_NULL, "-is_quasi_static", &isQuasiStatic,
PETSC_NULL);
CHKERR PetscOptionsGetInt(PETSC_NULL, "-order", &oRder, PETSC_NULL);
};
auto simple = mField.getInterface<Simple>();
3);
CHKERR simple->addBoundaryField(positionField, H1,
CHKERR simple->setFieldOrder(positionField, oRder);
}
3);
CHKERR simple->setFieldOrder(meshNodeField, 2);
}
};
elasticElementPtr = boost::make_shared<NonlinearElasticElement>(mField, 2);
elasticMaterialsPtr = boost::make_shared<ElasticMaterials>(mField, "HOOKE");
CHKERR elasticMaterialsPtr->setBlocks(elasticElementPtr->setOfBlocks);
false, false, false);
false, false, false);
true, false, false, false);
CHKERR elasticElementPtr->addElement("ELASTIC", positionField,
meshNodeField, false);
};
}
MoFEMErrorCode addElementsToDM(SmartPetscObj<DM> dm) {
this->dM = dm;
};
dM, "ELASTIC", &elasticElementPtr->getLoopFeLhs(), NULL, NULL);
};
&elasticElementPtr->getLoopFeRhs(),
PETSC_NULL, PETSC_NULL);
};
switch (type) {
case IM:
dM, "ELASTIC", &elasticElementPtr->getLoopFeLhs(), NULL, NULL);
break;
case IM2:
dM, "ELASTIC", &elasticElementPtr->getLoopFeLhs(), NULL, NULL);
break;
default:
SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not implemented");
break;
}
};
switch (type) {
case IM:
dM, "ELASTIC", &elasticElementPtr->getLoopFeLhs(), NULL, NULL);
break;
case IM2:
dM, "ELASTIC", &elasticElementPtr->getLoopFeLhs(), NULL, NULL);
break;
default:
SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not implemented");
break;
}
};
postProcMeshPtr = boost::make_shared<PostProcVolumeOnRefinedMesh>(mField);
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";
&elasticElementPtr->getLoopFeEnergy());
// Print elastic energy
MOFEM_LOG_C("WORLD", Sev::inform, "Elastic energy %6.4e\n",
elasticElementPtr->getLoopFeEnergy().eNergy);
auto out_name = "out_" + to_string(step) + ".h5m";
MOFEM_LOG_C("WORLD", Sev::inform, "out file %s\n", out_name.c_str());
CHKERR postProcMeshPtr->writeFile(out_name);
};
};
#endif //__NONLINEARELEMENTINTERFACE_HPP__
#define MOFEM_LOG_C(channel, severity, format,...)
Definition: LogManager.hpp:314
void simple(double P1[], double P2[], double P3[], double c[], const int N)
Definition: acoustic.cpp:69
constexpr bool is_quasi_static
Definition: contact.cpp:131
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:73
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:460
@ H1
continuous field
Definition: definitions.h:98
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
@ MOFEM_NOT_IMPLEMENTED
Definition: definitions.h:45
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
#define CHKERR
Inline error check.
Definition: definitions.h:548
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:453
PetscErrorCode DMMoFEMSNESSetFunction(DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
set SNES residual evaluation function
Definition: DMMMoFEM.cpp:677
PetscErrorCode DMMoFEMTSSetIFunction(DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
set TS implicit function evaluation function
Definition: DMMMoFEM.cpp:759
PetscErrorCode DMMoFEMAddElement(DM dm, const char fe_name[])
add element to dm
Definition: DMMMoFEM.cpp:462
PetscErrorCode DMMoFEMSNESSetJacobian(DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
set SNES Jacobian evaluation function
Definition: DMMMoFEM.cpp:718
PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[], MoFEM::FEMethod *method, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr())
Executes FEMethod for finite elements in DM.
Definition: DMMMoFEM.cpp:545
PetscErrorCode DMMoFEMTSSetIJacobian(DM dm, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
set TS Jacobian evaluation function
Definition: DMMMoFEM.cpp:812
PetscErrorCode DMMoFEMTSSetI2Jacobian(DM dm, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
set TS Jacobian evaluation function
Definition: DMMMoFEM.cpp:976
PetscErrorCode DMMoFEMTSSetI2Function(DM dm, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
set TS implicit function evaluation function
Definition: DMMMoFEM.cpp:934
virtual bool check_field(const std::string &name) const =0
check if field is in database
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:311
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:51
PetscErrorCode PetscOptionsGetInt(PetscOptions *, const char pre[], const char name[], PetscInt *ivalue, PetscBool *set)
PetscErrorCode PetscOptionsGetBool(PetscOptions *, const char pre[], const char name[], PetscBool *bval, PetscBool *set)
MoFEMErrorCode addHOOpsVol(const std::string field, E &e, bool h1, bool hcurl, bool hdiv, bool l2)
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.
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