v0.9.0
Public Member Functions | Public Attributes | List of all members
BoneRemodeling::MonitorPostProc Struct Reference

#include <users_modules/bone_remodelling/src/Remodeling.hpp>

Inheritance diagram for BoneRemodeling::MonitorPostProc:
[legend]
Collaboration diagram for BoneRemodeling::MonitorPostProc:
[legend]

Public Member Functions

 MonitorPostProc (MoFEM::Interface &m_field, Remodeling::CommonData &common_data)
 
MoFEMErrorCode preProcess ()
 
MoFEMErrorCode operator() ()
 
MoFEMErrorCode postProcess ()
 

Public Attributes

MoFEM::InterfacemField
 
Remodeling::CommonDatacommonData
 
PostProcVolumeOnRefinedMesh postProc
 
PostProcVolumeOnRefinedMesh postProcElastic
 
PetscBool mass_postproc
 
PetscBool equilibrium_flg
 
double rate
 
bool iNit
 
int pRT
 

Detailed Description

Element used to post-process results at each time step

Definition at line 553 of file Remodeling.hpp.

Constructor & Destructor Documentation

◆ MonitorPostProc()

MonitorPostProc::MonitorPostProc ( MoFEM::Interface m_field,
Remodeling::CommonData common_data 
)
Examples
Remodeling.hpp.

Definition at line 1408 of file Remodeling.cpp.

1410  : mField(m_field), commonData(common_data), postProc(m_field),
1411  postProcElastic(m_field),
1412  // densityMaps(m_field),
1413  iNit(false) {}
bool iNit
Definition: Remodeling.hpp:563
MoFEM::Interface & mField
Definition: Remodeling.hpp:555
PostProcVolumeOnRefinedMesh postProc
Definition: Remodeling.hpp:557
PostProcVolumeOnRefinedMesh postProcElastic
Definition: Remodeling.hpp:558
Remodeling::CommonData & commonData
Definition: Remodeling.hpp:556

Member Function Documentation

◆ operator()()

MoFEMErrorCode MonitorPostProc::operator() ( )
Examples
Remodeling.hpp.

Definition at line 1420 of file Remodeling.cpp.

1420  {
1423 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508

◆ postProcess()

MoFEMErrorCode MonitorPostProc::postProcess ( )
Examples
Remodeling.hpp.

Definition at line 1425 of file Remodeling.cpp.

1425  {
1427 
1428  PetscBool mass_postproc = PETSC_FALSE;
1429  PetscBool equilibrium_flg = PETSC_FALSE;
1430  PetscBool analysis_mesh_flg = PETSC_FALSE;
1431  rate = 1;
1432 
1433  CHKERR PetscOptionsBegin(PETSC_COMM_WORLD, "", "Bone remodeling post-process",
1434  "none");
1435 
1436  CHKERR PetscOptionsBool(
1437  "-mass_postproc",
1438  "is used for testing, calculates mass and energy at each time step", "",
1439  mass_postproc, &mass_postproc, PETSC_NULL);
1440 
1441  CHKERR PetscOptionsBool("-analysis_mesh",
1442  "saves analysis mesh at each time step", "",
1443  analysis_mesh_flg, &analysis_mesh_flg, PETSC_NULL);
1444 
1445  CHKERR PetscOptionsScalar(
1446  "-equilibrium_stop_rate",
1447  "is used to stop calculations when equilibium state is achieved", "",
1448  rate, &rate, &equilibrium_flg);
1449  ierr = PetscOptionsEnd();
1450  CHKERRQ(ierr);
1451 
1452  if (!iNit) {
1453 
1454  CHKERR PetscOptionsBegin(PETSC_COMM_WORLD, "",
1455  "Bone remodeling post-process", "none");
1456 
1457  pRT = 1;
1458  CHKERR PetscOptionsInt("-my_output_prt",
1459  "frequncy how often results are dumped on hard disk",
1460  "", 1, &pRT, NULL);
1461  ierr = PetscOptionsEnd();
1462  CHKERRQ(ierr);
1463 
1465  CHKERR postProc.addFieldValuesPostProc("DISPLACEMENTS");
1466  CHKERR postProc.addFieldValuesPostProc("MESH_NODE_POSITIONS");
1467  if (!commonData.less_post_proc) {
1469  }
1470  // CHKERR postProc.addFieldValuesPostProc("RHO");
1471  // CHKERR postProc.addFieldValuesGradientPostProc("RHO");
1472 
1473  {
1474  boost::ptr_vector<MoFEM::ForcesAndSourcesCore::UserDataOperator>
1475  &list_of_operators = postProc.getOpPtrVector();
1476 
1477  list_of_operators.push_back(
1479  list_of_operators.push_back(new OpCalculateScalarFieldGradient<3>(
1480  "RHO", commonData.data.gradRhoPtr));
1481  // Get displacement gradient at integration points
1482  list_of_operators.push_back(new OpCalculateVectorFieldGradient<3, 3>(
1483  "DISPLACEMENTS", commonData.data.gradDispPtr));
1484  list_of_operators.push_back(new OpCalculateStress(commonData));
1485  list_of_operators.push_back(new OpPostProcStress(
1487  }
1488 
1491  CHKERR postProcElastic.addFieldValuesPostProc("MESH_NODE_POSITIONS");
1493  std::map<int, NonlinearElasticElement::BlockData>::iterator sit =
1494  commonData.elasticPtr->setOfBlocks.begin();
1495  for (; sit != commonData.elasticPtr->setOfBlocks.end(); sit++) {
1498  "DISPLACEMENTS", sit->second, postProcElastic.commonData));
1499  }
1500 
1501  iNit = true;
1502  }
1503 
1504  if (mass_postproc) {
1505  Vec mass_vec;
1506  Vec energ_vec;
1507  int mass_vec_ghost[] = {0};
1508  CHKERR VecCreateGhost(PETSC_COMM_WORLD, (!mField.get_comm_rank()) ? 1 : 0,
1509  1, 1, mass_vec_ghost, &mass_vec);
1510 
1511  CHKERR VecZeroEntries(mass_vec);
1512  CHKERR VecDuplicate(mass_vec, &energ_vec);
1513 
1514  Remodeling::Fe energyCalc(mField);
1515 
1516  energyCalc.getOpPtrVector().push_back(
1518  energyCalc.getOpPtrVector().push_back(new OpCalculateScalarFieldGradient<3>(
1519  "RHO", commonData.data.gradRhoPtr));
1520  energyCalc.getOpPtrVector().push_back(
1521  new OpCalculateVectorFieldGradient<3, 3>("DISPLACEMENTS",
1523  energyCalc.getOpPtrVector().push_back(new OpCalculateStress(commonData));
1524  energyCalc.getOpPtrVector().push_back(
1525  new OpMassAndEnergyCalculation("RHO", commonData, energ_vec, mass_vec));
1526  CHKERR DMoFEMLoopFiniteElements(commonData.dm, "FE_REMODELLING",
1527  &energyCalc);
1528 
1529  CHKERR VecAssemblyBegin(mass_vec);
1530  CHKERR VecAssemblyEnd(mass_vec);
1531  double mass_sum;
1532  CHKERR VecSum(mass_vec, &mass_sum);
1533 
1534  CHKERR VecAssemblyBegin(energ_vec);
1535  CHKERR VecAssemblyEnd(energ_vec);
1536  double energ_sum;
1537  CHKERR VecSum(energ_vec, &energ_sum);
1538 
1539  CHKERR PetscPrintf(PETSC_COMM_WORLD,
1540  "Time: %g Mass: %6.5g Elastic energy: %6.5g \n", ts_t,
1541  mass_sum, energ_sum);
1542  CHKERR VecDestroy(&mass_vec);
1543  CHKERR VecDestroy(&energ_vec);
1544 
1545  if (equilibrium_flg) {
1546  double equilibrium_rate =
1547  fabs(energ_sum - commonData.cUrrent_psi) / energ_sum;
1548  double equilibrium_mass_rate =
1549  fabs(mass_sum - commonData.cUrrent_mass) / mass_sum;
1550  if (equilibrium_rate < rate) {
1551  CHKERR PetscPrintf(
1552  PETSC_COMM_WORLD,
1553  "Energy equilibrium state is achieved! Difference = %0.6g %%. \n",
1554  equilibrium_rate * 100);
1555  }
1556  if (equilibrium_mass_rate < rate) {
1557  CHKERR PetscPrintf(
1558  PETSC_COMM_WORLD,
1559  "Mass equilibrium state is achieved! Difference = %0.6g %%. \n",
1560  equilibrium_mass_rate * 100);
1561  }
1562  commonData.cUrrent_psi = energ_sum;
1563  commonData.cUrrent_mass = mass_sum;
1564  }
1565  }
1566 
1567  int step;
1568 #if PETSC_VERSION_LE(3, 8, 0)
1569  CHKERR TSGetTimeStepNumber(ts, &step);
1570 #else
1571  CHKERR TSGetStepNumber(ts, &step);
1572 #endif
1573 
1574  if ((step) % pRT == 0) {
1575 
1576  ostringstream sss;
1577  sss << "out_" << step << ".h5m";
1578  CHKERR DMoFEMLoopFiniteElements(commonData.dm, "FE_REMODELLING", &postProc);
1579  CHKERR postProc.postProcMesh.write_file(sss.str().c_str(), "MOAB",
1580  "PARALLEL=WRITE_PART");
1581  if (analysis_mesh_flg) {
1582  ostringstream ttt;
1583  ttt << "analysis_mesh_" << step << ".h5m";
1584  CHKERR mField.get_moab().write_file(ttt.str().c_str(), "MOAB",
1585  "PARALLEL=WRITE_PART");
1586  }
1587 
1590  &postProcElastic);
1591  ostringstream ss;
1592  ss << "out_elastic_" << step << ".h5m";
1593  CHKERR postProcElastic.postProcMesh.write_file(ss.str().c_str(), "MOAB",
1594  "PARALLEL=WRITE_PART");
1595  }
1596  }
1598 }
PetscBool equilibrium_flg
Definition: Remodeling.hpp:561
MoFEMErrorCode addFieldValuesGradientPostProc(const std::string field_name, Vec v=PETSC_NULL)
Add operator to post-process L2 or H1 field gradient.
boost::shared_ptr< VectorDouble > rhoPtr
Ptr to density matrix container.
Definition: Remodeling.hpp:187
CommonData commonData
PetscBool less_post_proc
reduce file size
Definition: Remodeling.hpp:163
virtual moab::Interface & get_moab()=0
double cUrrent_psi
current free energy for evaluating equilibrium state
Definition: Remodeling.hpp:158
boost::shared_ptr< MatrixDouble > gradRhoPtr
Gradient of density field.
Definition: Remodeling.hpp:188
moab::Interface & postProcMesh
bool iNit
Definition: Remodeling.hpp:563
Get field gradients at integration pts for scalar filed rank 0, i.e. vector field.
double cUrrent_mass
current free energy for evaluating equilibrium state
Definition: Remodeling.hpp:160
double rate
Definition: Remodeling.hpp:562
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
PetscBool mass_postproc
Definition: Remodeling.hpp:560
PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[], MoFEM::FEMethod *method)
Executes FEMethod for finite elements in DM.
Definition: DMMMoFEM.cpp:496
boost::shared_ptr< MatrixDouble > gradDispPtr
Ptr to gradient of displacements matrix container.
Definition: Remodeling.hpp:185
virtual int get_comm_rank() const =0
MoFEMErrorCode addFieldValuesPostProc(const std::string field_name, Vec v=PETSC_NULL)
Add operator to post-process L2, H1, Hdiv, Hcurl field value.
MoFEM::Interface & mField
Definition: Remodeling.hpp:555
CHKERRQ(ierr)
int pRT
Definition: Remodeling.hpp:564
boost::ptr_vector< UserDataOperator > & getOpPtrVector()
Use to push back operator for row operator.
#define CHKERR
Inline error check.
Definition: definitions.h:596
PostProcVolumeOnRefinedMesh postProc
Definition: Remodeling.hpp:557
std::vector< EntityHandle > mapGaussPts
MoFEMErrorCode generateReferenceElementMesh()
Generate reference mesh on single element.
DM dm
Discretization manager.
Definition: Remodeling.hpp:127
Get field gradients at integration pts for scalar filed rank 0, i.e. vector field.
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407
PostProcVolumeOnRefinedMesh postProcElastic
Definition: Remodeling.hpp:558
Remodeling::CommonData & commonData
Definition: Remodeling.hpp:556
boost::shared_ptr< NonlinearElasticElement > elasticPtr
Definition: Remodeling.hpp:140
Get value at integration points for scalar field.

◆ preProcess()

MoFEMErrorCode MonitorPostProc::preProcess ( )
Examples
Remodeling.hpp.

Definition at line 1415 of file Remodeling.cpp.

1415  {
1418 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508

Member Data Documentation

◆ commonData

Remodeling::CommonData& BoneRemodeling::MonitorPostProc::commonData
Examples
Remodeling.hpp.

Definition at line 556 of file Remodeling.hpp.

◆ equilibrium_flg

PetscBool BoneRemodeling::MonitorPostProc::equilibrium_flg
Examples
Remodeling.hpp.

Definition at line 561 of file Remodeling.hpp.

◆ iNit

bool BoneRemodeling::MonitorPostProc::iNit
Examples
Remodeling.hpp.

Definition at line 563 of file Remodeling.hpp.

◆ mass_postproc

PetscBool BoneRemodeling::MonitorPostProc::mass_postproc
Examples
Remodeling.hpp.

Definition at line 560 of file Remodeling.hpp.

◆ mField

MoFEM::Interface& BoneRemodeling::MonitorPostProc::mField
Examples
Remodeling.hpp.

Definition at line 555 of file Remodeling.hpp.

◆ postProc

PostProcVolumeOnRefinedMesh BoneRemodeling::MonitorPostProc::postProc
Examples
Remodeling.hpp.

Definition at line 557 of file Remodeling.hpp.

◆ postProcElastic

PostProcVolumeOnRefinedMesh BoneRemodeling::MonitorPostProc::postProcElastic
Examples
Remodeling.hpp.

Definition at line 558 of file Remodeling.hpp.

◆ pRT

int BoneRemodeling::MonitorPostProc::pRT
Examples
Remodeling.hpp.

Definition at line 564 of file Remodeling.hpp.

◆ rate

double BoneRemodeling::MonitorPostProc::rate
Examples
Remodeling.hpp.

Definition at line 562 of file Remodeling.hpp.


The documentation for this struct was generated from the following files: