|
| v0.14.0
|
#include <users_modules/bone_remodelling/src/Remodeling.hpp>
Element used to post-process results at each time step
Definition at line 553 of file Remodeling.hpp.
◆ MonitorPostProc()
◆ operator()()
◆ postProcess()
- Examples
- Remodeling.hpp.
Definition at line 1404 of file Remodeling.cpp.
1409 PetscBool analysis_mesh_flg = PETSC_FALSE;
1412 CHKERR PetscOptionsBegin(PETSC_COMM_WORLD,
"",
"Bone remodeling post-process",
1417 "is used for testing, calculates mass and energy at each time step",
"",
1420 CHKERR PetscOptionsBool(
"-analysis_mesh",
1421 "saves analysis mesh at each time step",
"",
1422 analysis_mesh_flg, &analysis_mesh_flg, PETSC_NULL);
1424 CHKERR PetscOptionsScalar(
1425 "-equilibrium_stop_rate",
1426 "is used to stop calculations when equilibium state is achieved",
"",
1428 ierr = PetscOptionsEnd();
1433 CHKERR PetscOptionsBegin(PETSC_COMM_WORLD,
"",
1434 "Bone remodeling post-process",
"none");
1437 CHKERR PetscOptionsInt(
"-my_output_prt",
1438 "frequncy how often results are dumped on hard disk",
1440 ierr = PetscOptionsEnd();
1455 list_of_operators.push_back(
1471 std::map<int, NonlinearElasticElement::BlockData>::iterator sit =
1485 int mass_vec_ghost[] = {0};
1487 1, 1, mass_vec_ghost, &mass_vec);
1489 CHKERR VecZeroEntries(mass_vec);
1490 CHKERR VecDuplicate(mass_vec, &energ_vec);
1492 Remodeling::Fe energyCalc(
mField);
1495 energyCalc.getOpPtrVector().push_back(
1499 energyCalc.getOpPtrVector().push_back(
1503 energyCalc.getOpPtrVector().push_back(
1504 new OpMassAndEnergyCalculation(
"RHO",
commonData, energ_vec, mass_vec));
1508 CHKERR VecAssemblyBegin(mass_vec);
1509 CHKERR VecAssemblyEnd(mass_vec);
1511 CHKERR VecSum(mass_vec, &mass_sum);
1513 CHKERR VecAssemblyBegin(energ_vec);
1514 CHKERR VecAssemblyEnd(energ_vec);
1516 CHKERR VecSum(energ_vec, &energ_sum);
1518 CHKERR PetscPrintf(PETSC_COMM_WORLD,
1519 "Time: %g Mass: %6.5g Elastic energy: %6.5g \n", ts_t,
1520 mass_sum, energ_sum);
1521 CHKERR VecDestroy(&mass_vec);
1522 CHKERR VecDestroy(&energ_vec);
1525 double equilibrium_rate =
1527 double equilibrium_mass_rate =
1529 if (equilibrium_rate <
rate) {
1532 "Energy equilibrium state is achieved! Difference = %0.6g %%. \n",
1533 equilibrium_rate * 100);
1535 if (equilibrium_mass_rate <
rate) {
1538 "Mass equilibrium state is achieved! Difference = %0.6g %%. \n",
1539 equilibrium_mass_rate * 100);
1547 #if PETSC_VERSION_LE(3, 8, 0)
1548 CHKERR TSGetTimeStepNumber(ts, &step);
1550 CHKERR TSGetStepNumber(ts, &step);
1553 if ((step) %
pRT == 0) {
1556 sss <<
"out_" << step <<
".h5m";
1559 "PARALLEL=WRITE_PART");
1560 if (analysis_mesh_flg) {
1562 ttt <<
"analysis_mesh_" << step <<
".h5m";
1564 "PARALLEL=WRITE_PART");
1571 ss <<
"out_elastic_" << step <<
".h5m";
1573 "PARALLEL=WRITE_PART");
◆ preProcess()
◆ commonData
◆ equilibrium_flg
PetscBool BoneRemodeling::MonitorPostProc::equilibrium_flg |
◆ iNit
bool BoneRemodeling::MonitorPostProc::iNit |
◆ mass_postproc
PetscBool BoneRemodeling::MonitorPostProc::mass_postproc |
◆ mField
◆ postProc
◆ postProcElastic
◆ pRT
int BoneRemodeling::MonitorPostProc::pRT |
◆ rate
double BoneRemodeling::MonitorPostProc::rate |
The documentation for this struct was generated from the following files:
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
moab::Interface & postProcMesh
MoFEMErrorCode addHOOpsVol(const std::string field, E &e, bool h1, bool hcurl, bool hdiv, bool l2)
virtual int get_comm_rank() const =0
double cUrrent_psi
current free energy for evaluating equilibrium state
double cUrrent_mass
current free energy for evaluating equilibrium state
boost::shared_ptr< VectorDouble > rhoPtr
Ptr to density matrix container.
boost::shared_ptr< MatrixDouble > gradDispPtr
Ptr to gradient of displacements matrix container.
Get field gradients at integration pts for scalar filed rank 0, i.e. vector field.
PetscBool less_post_proc
reduce file size
#define CHKERR
Inline error check.
virtual moab::Interface & get_moab()=0
PostProcVolumeOnRefinedMesh postProcElastic
boost::shared_ptr< NonlinearElasticElement > elasticPtr
PetscBool equilibrium_flg
boost::shared_ptr< MatrixDouble > gradRhoPtr
Gradient of density field.
Get value at integration points for scalar field.
MoFEMErrorCode generateReferenceElementMesh()
Generate reference mesh on single element.
MoFEMErrorCode addFieldValuesGradientPostProc(const std::string field_name, Vec v=PETSC_NULL)
Add operator to post-process L2 or H1 field gradient.
Get field gradients at integration pts for scalar filed rank 0, i.e. vector field.
PostProcVolumeOnRefinedMesh postProc
std::vector< EntityHandle > mapGaussPts
MoFEM::Interface & mField
DM dm
Discretization manager.
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
boost::ptr_deque< UserDataOperator > & getOpPtrVector()
Use to push back operator for row operator.
const FTensor::Tensor2< T, Dim, Dim > Vec
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[], MoFEM::FEMethod *method, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr())
Executes FEMethod for finite elements in DM.
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Remodeling::CommonData & commonData
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
MoFEMErrorCode addFieldValuesPostProc(const std::string field_name, Vec v=PETSC_NULL)
Add operator to post-process L2, H1, Hdiv, Hcurl field value.