v0.15.0
Loading...
Searching...
No Matches
Monitor Struct Reference

[Push operators to pipeline] More...

#include "tutorials/vec-2/src/Monitor.hpp"

Inheritance diagram for Monitor:
[legend]
Collaboration diagram for Monitor:
[legend]

Public Member Functions

 Monitor (SmartPetscObj< DM > dm, MoFEM::Interface &m_field, boost::shared_ptr< PostProcEle > post_proc, boost::shared_ptr< PostProcFaceEle > post_proc_bdry, boost::shared_ptr< MatrixDouble > velocity_field_ptr, boost::shared_ptr< MatrixDouble > x2_field_ptr, boost::shared_ptr< MatrixDouble > geometry_field_ptr, std::array< double, 3 > pass_field_eval_coords, boost::shared_ptr< SetPtsData > pass_field_eval_data)
 
MoFEMErrorCode postProcess ()
 function is run at the end of loop
 
 Monitor (SmartPetscObj< DM > dm, boost::shared_ptr< PostProcEle > post_proc)
 
MoFEMErrorCode preProcess ()
 function is run at the beginning of loop
 
MoFEMErrorCode operator() ()
 function is run for every finite element
 
MoFEMErrorCode postProcess ()
 function is run at the end of loop
 
 Monitor (SmartPetscObj< DM > dm, boost::shared_ptr< PostProcEle > post_proc)
 
MoFEMErrorCode preProcess ()
 function is run at the beginning of loop
 
MoFEMErrorCode operator() ()
 function is run for every finite element
 
MoFEMErrorCode postProcess ()
 function is run at the end of loop
 
 Monitor (SmartPetscObj< DM > dm, std::pair< boost::shared_ptr< PostProcEleDomain >, boost::shared_ptr< PostProcEleBdy > > pair_post_proc, boost::shared_ptr< DomainEle > reaction_fe)
 
MoFEMErrorCode postProcess ()
 
 Monitor (SmartPetscObj< DM > dm, boost::shared_ptr< PostProcEle > post_proc)
 
MoFEMErrorCode postProcess ()
 function is run at the end of loop
 
 Monitor (SmartPetscObj< DM > dm, boost::shared_ptr< PostProcEle > post_proc)
 
MoFEMErrorCode postProcess ()
 function is run at the end of loop
 
 Monitor (SmartPetscObj< DM > dm, boost::shared_ptr< moab::Core > post_proc_mesh, boost::shared_ptr< PostProcEleDomainCont > post_proc, boost::shared_ptr< PostProcEleBdyCont > post_proc_edge, std::pair< boost::shared_ptr< BoundaryEle >, boost::shared_ptr< VectorDouble > > p)
 
MoFEMErrorCode postProcess ()
 function is run at the end of loop
 
 Monitor (SmartPetscObj< DM > dm, std::pair< boost::shared_ptr< PostProcEleDomain >, boost::shared_ptr< PostProcEleBdy > > pair_post_proc_fe, boost::shared_ptr< DomainEle > reaction_fe, std::vector< boost::shared_ptr< ScalingMethod > > smv)
 
MoFEMErrorCode postProcess ()
 function is run at the end of loop
 
- Public Member Functions inherited from MoFEM::FEMethod
MoFEMErrorCode query_interface (boost::typeindex::type_index type_index, UnknownInterface **iface) const
 
 FEMethod ()=default
 
auto getFEName () const
 get finite element name
 
auto getDataDofsPtr () const
 
auto getDataVectorDofsPtr () const
 
const FieldEntity_vector_viewgetDataFieldEnts () const
 
boost::shared_ptr< FieldEntity_vector_view > & getDataFieldEntsPtr () const
 
auto & getRowFieldEnts () const
 
auto & getRowFieldEntsPtr () const
 
auto & getColFieldEnts () const
 
auto & getColFieldEntsPtr () const
 
auto getRowDofsPtr () const
 
auto getColDofsPtr () const
 
auto getNumberOfNodes () const
 
EntityHandle getFEEntityHandle () const
 
MoFEMErrorCode getNodeData (const std::string field_name, VectorDouble &data, const bool reset_dofs=true)
 
- Public Member Functions inherited from MoFEM::BasicMethod
 BasicMethod ()
 
virtual ~BasicMethod ()=default
 
int getNinTheLoop () const
 get number of evaluated element in the loop
 
int getLoopSize () const
 get loop size
 
auto getLoHiFERank () const
 Get lo and hi processor rank of iterated entities.
 
auto getLoFERank () const
 Get upper rank in loop for iterating elements.
 
auto getHiFERank () const
 Get upper rank in loop for iterating elements.
 
unsigned int getFieldBitNumber (std::string field_name) const
 
MoFEMErrorCode copyBasicMethod (const BasicMethod &basic)
 Copy data from other base method to this base method.
 
boost::weak_ptr< CacheTuplegetCacheWeakPtr () const
 Get the cache weak ptr object.
 
- Public Member Functions inherited from MoFEM::KspMethod
MoFEMErrorCode query_interface (boost::typeindex::type_index type_index, UnknownInterface **iface) const
 
 KspMethod ()
 
virtual ~KspMethod ()=default
 
MoFEMErrorCode copyKsp (const KspMethod &ksp)
 copy data form another method
 
- Public Member Functions inherited from MoFEM::PetscData
 PetscData ()
 
virtual ~PetscData ()=default
 
MoFEMErrorCode copyPetscData (const PetscData &petsc_data)
 
- Public Member Functions inherited from MoFEM::UnknownInterface
template<class IFACE >
MoFEMErrorCode registerInterface (bool error_if_registration_failed=true)
 Register interface.
 
template<class IFACE >
MoFEMErrorCode getInterface (IFACE *&iface) const
 Get interface reference to pointer of interface.
 
template<class IFACE >
MoFEMErrorCode getInterface (IFACE **const iface) const
 Get interface pointer to pointer of interface.
 
template<class IFACE , typename boost::enable_if< boost::is_pointer< IFACE >, int >::type = 0>
IFACE getInterface () const
 Get interface pointer to pointer of interface.
 
template<class IFACE , typename boost::enable_if< boost::is_reference< IFACE >, int >::type = 0>
IFACE getInterface () const
 Get reference to interface.
 
template<class IFACE >
IFACE * getInterface () const
 Function returning pointer to interface.
 
virtual ~UnknownInterface ()=default
 
- Public Member Functions inherited from MoFEM::SnesMethod
MoFEMErrorCode query_interface (boost::typeindex::type_index type_index, UnknownInterface **iface) const
 
 SnesMethod ()
 
virtual ~SnesMethod ()=default
 
MoFEMErrorCode copySnes (const SnesMethod &snes)
 Copy snes data.
 
- Public Member Functions inherited from MoFEM::TSMethod
MoFEMErrorCode query_interface (boost::typeindex::type_index type_index, UnknownInterface **iface) const
 
 TSMethod ()
 
virtual ~TSMethod ()=default
 
MoFEMErrorCode copyTs (const TSMethod &ts)
 Copy TS solver data.
 
- Public Member Functions inherited from MoFEM::TaoMethod
MoFEMErrorCode query_interface (boost::typeindex::type_index type_index, UnknownInterface **iface) const
 
 TaoMethod ()
 
virtual ~TaoMethod ()=default
 
MoFEMErrorCode copyTao (const TaoMethod &tao)
 Copy TAO data.
 

Public Attributes

MoFEM::InterfacemField
 
- Public Attributes inherited from MoFEM::FEMethod
std::string feName
 Name of finite element.
 
boost::shared_ptr< const NumeredEntFiniteElementnumeredEntFiniteElementPtr
 
boost::function< bool(FEMethod *fe_method_ptr)> exeTestHook
 Tet if element to skip element.
 
- Public Attributes inherited from MoFEM::BasicMethod
int nInTheLoop
 number currently of processed method
 
int loopSize
 local number oe methods to process
 
std::pair< int, int > loHiFERank
 Llo and hi processor rank of iterated entities.
 
int rAnk
 processor rank
 
int sIze
 number of processors in communicator
 
const RefEntity_multiIndexrefinedEntitiesPtr
 container of mofem dof entities
 
const RefElement_multiIndexrefinedFiniteElementsPtr
 container of mofem finite element entities
 
const ProblemproblemPtr
 raw pointer to problem
 
const Field_multiIndexfieldsPtr
 raw pointer to fields container
 
const FieldEntity_multiIndexentitiesPtr
 raw pointer to container of field entities
 
const DofEntity_multiIndexdofsPtr
 raw pointer container of dofs
 
const FiniteElement_multiIndexfiniteElementsPtr
 raw pointer to container finite elements
 
const EntFiniteElement_multiIndexfiniteElementsEntitiesPtr
 
const FieldEntityEntFiniteElementAdjacencyMap_multiIndexadjacenciesPtr
 
boost::function< MoFEMErrorCode()> preProcessHook
 Hook function for pre-processing.
 
boost::function< MoFEMErrorCode()> postProcessHook
 Hook function for post-processing.
 
boost::function< MoFEMErrorCode()> operatorHook
 Hook function for operator.
 
boost::movelib::unique_ptr< boolvecAssembleSwitch
 
boost::movelib::unique_ptr< boolmatAssembleSwitch
 
boost::weak_ptr< CacheTuplecacheWeakPtr
 
- Public Attributes inherited from MoFEM::KspMethod
KSPContext ksp_ctx
 Context.
 
KSP ksp
 KSP solver.
 
Vec & ksp_f
 
Mat & ksp_A
 
Mat & ksp_B
 
- Public Attributes inherited from MoFEM::PetscData
Switches data_ctx
 
Vec f
 
Mat A
 
Mat B
 
Vec x
 
Vec dx
 
Vec x_t
 
Vec x_tt
 
- Public Attributes inherited from MoFEM::SnesMethod
SNESContext snes_ctx
 
SNES snes
 snes solver
 
Vec & snes_x
 state vector
 
Vec & snes_dx
 solution update
 
Vec & snes_f
 residual
 
Mat & snes_A
 jacobian matrix
 
Mat & snes_B
 preconditioner of jacobian matrix
 
- Public Attributes inherited from MoFEM::TSMethod
TS ts
 time solver
 
TSContext ts_ctx
 
PetscInt ts_step
 time step number
 
PetscReal ts_a
 shift for U_t (see PETSc Time Solver)
 
PetscReal ts_aa
 shift for U_tt shift for U_tt
 
PetscReal ts_t
 time
 
PetscReal ts_dt
 time step size
 
Vec & ts_u
 state vector
 
Vec & ts_u_t
 time derivative of state vector
 
Vec & ts_u_tt
 second time derivative of state vector
 
Vec & ts_F
 residual vector
 
Mat & ts_A
 
Mat & ts_B
 Preconditioner for ts_A.
 
- Public Attributes inherited from MoFEM::TaoMethod
TAOContext tao_ctx
 
Tao tao
 tao solver
 
Vec & tao_x
 
Vec & tao_f
 state vector
 
Mat & tao_A
 gradient vector
 
Mat & tao_B
 hessian matrix
 

Static Public Attributes

static constexpr int saveEveryNthStep = 1
 
- Static Public Attributes inherited from MoFEM::PetscData
static constexpr Switches CtxSetNone = PetscData::Switches(CTX_SET_NONE)
 
static constexpr Switches CtxSetF = PetscData::Switches(CTX_SET_F)
 
static constexpr Switches CtxSetA = PetscData::Switches(CTX_SET_A)
 
static constexpr Switches CtxSetB = PetscData::Switches(CTX_SET_B)
 
static constexpr Switches CtxSetX = PetscData::Switches(CTX_SET_X)
 
static constexpr Switches CtxSetDX = PetscData::Switches(CTX_SET_DX)
 
static constexpr Switches CtxSetX_T = PetscData::Switches(CTX_SET_X_T)
 
static constexpr Switches CtxSetX_TT = PetscData::Switches(CTX_SET_X_TT)
 
static constexpr Switches CtxSetTime = PetscData::Switches(CTX_SET_TIME)
 

Private Attributes

SmartPetscObj< DM > dM
 
boost::shared_ptr< PostProcElepostProc
 
boost::shared_ptr< PostProcFaceElepostProcBdy
 
boost::shared_ptr< MatrixDoublevelocityFieldPtr
 
boost::shared_ptr< MatrixDoublex2FieldPtr
 
boost::shared_ptr< MatrixDoublegeometryFieldPtr
 
std::array< double, 3 > fieldEvalCoords
 
boost::shared_ptr< SetPtsDatafieldEvalData
 
boost::shared_ptr< PostProcEleDomainpostProcDomain
 
boost::shared_ptr< PostProcEleBdypostProcSkin
 
boost::shared_ptr< DomainElereactionFe
 
PetscInt saveEvery = 1
 
boost::shared_ptr< moab::Core > postProcMesh
 
boost::shared_ptr< PostProcEleDomainContpostProc
 
boost::shared_ptr< PostProcEleBdyContpostProcEdge
 
boost::shared_ptr< BoundaryEleliftFE
 
boost::shared_ptr< VectorDoubleliftVec
 
boost::shared_ptr< PostProcEleDomaindomainPostProcFe
 
boost::shared_ptr< PostProcEleBdyskinPostProcFe
 
boost::shared_ptr< FEMethodreactionFE
 
boost::shared_ptr< BoundaryEleintegrateTraction
 
SmartPetscObj< Vec > fRes
 
VecOfTimeScalingMethods vecOfTimeScalingMethods
 

Additional Inherited Members

- Public Types inherited from MoFEM::KspMethod
enum  KSPContext { CTX_SETFUNCTION , CTX_OPERATORS , CTX_KSPNONE }
 pass information about context of KSP/DM for with finite element is computed More...
 
- Public Types inherited from MoFEM::PetscData
enum  DataContext {
  CTX_SET_NONE = 0 , CTX_SET_F = 1 << 0 , CTX_SET_A = 1 << 1 , CTX_SET_B = 1 << 2 ,
  CTX_SET_X = 1 << 3 , CTX_SET_DX = 1 << 4 , CTX_SET_X_T = 1 << 5 , CTX_SET_X_TT = 1 << 6 ,
  CTX_SET_TIME = 1 << 7
}
 
using Switches = std::bitset<8>
 
- Public Types inherited from MoFEM::SnesMethod
enum  SNESContext { CTX_SNESSETFUNCTION , CTX_SNESSETJACOBIAN , CTX_SNESNONE }
 
- Public Types inherited from MoFEM::TSMethod
enum  TSContext {
  CTX_TSSETRHSFUNCTION , CTX_TSSETRHSJACOBIAN , CTX_TSSETIFUNCTION , CTX_TSSETIJACOBIAN ,
  CTX_TSTSMONITORSET , CTX_TSNONE
}
 
- Public Types inherited from MoFEM::TaoMethod
enum  TAOContext { CTX_TAO_OBJECTIVE , CTX_TAO_GRADIENT , CTX_TAO_HESSIAN , CTX_TAO_NONE }
 
- Static Public Member Functions inherited from MoFEM::UnknownInterface
static MoFEMErrorCode getLibVersion (Version &version)
 Get library version.
 
static MoFEMErrorCode getFileVersion (moab::Interface &moab, Version &version)
 Get database major version.
 
static MoFEMErrorCode setFileVersion (moab::Interface &moab, Version version=Version(MoFEM_VERSION_MAJOR, MoFEM_VERSION_MINOR, MoFEM_VERSION_BUILD))
 Get database major version.
 
static MoFEMErrorCode getInterfaceVersion (Version &version)
 Get database major version.
 

Detailed Description

[Push operators to pipeline]

[Push operators to pip]

Monitor class for post-processing results.

Monitor solution.

Monitor solution

This functions is called by TS solver at the end of each step. It is used to output results to the hard drive.

This functions is called by TS solver at the end of each step. It is used to output results to the hard drive.

Monitor solution

This functions is called by TS kso at the end of each step. It is used

Examples
PlasticOpsMonitor.hpp, adolc_plasticity.cpp, dynamic_first_order_con_law.cpp, free_surface.cpp, heat_equation.cpp, photon_diffusion.cpp, plastic.cpp, reaction_diffusion.cpp, shallow_wave.cpp, and wave_equation.cpp.

Definition at line 774 of file dynamic_first_order_con_law.cpp.

Constructor & Destructor Documentation

◆ Monitor() [1/8]

Monitor::Monitor ( SmartPetscObj< DM > dm,
MoFEM::Interface & m_field,
boost::shared_ptr< PostProcEle > post_proc,
boost::shared_ptr< PostProcFaceEle > post_proc_bdry,
boost::shared_ptr< MatrixDouble > velocity_field_ptr,
boost::shared_ptr< MatrixDouble > x2_field_ptr,
boost::shared_ptr< MatrixDouble > geometry_field_ptr,
std::array< double, 3 > pass_field_eval_coords,
boost::shared_ptr< SetPtsData > pass_field_eval_data )
inline
Examples
adolc_plasticity.cpp, dynamic_first_order_con_law.cpp, free_surface.cpp, heat_equation.cpp, shallow_wave.cpp, and wave_equation.cpp.

Definition at line 776 of file dynamic_first_order_con_law.cpp.

784 : dM(dm), mField(m_field), postProc(post_proc),
785 postProcBdy(post_proc_bdry), velocityFieldPtr(velocity_field_ptr),
786 x2FieldPtr(x2_field_ptr), geometryFieldPtr(geometry_field_ptr),
787 fieldEvalCoords(pass_field_eval_coords),
788 fieldEvalData(pass_field_eval_data){};
boost::shared_ptr< PostProcFaceEle > postProcBdy
std::array< double, 3 > fieldEvalCoords
MoFEM::Interface & mField
boost::shared_ptr< MatrixDouble > geometryFieldPtr
SmartPetscObj< DM > dM
boost::shared_ptr< MatrixDouble > velocityFieldPtr
boost::shared_ptr< SetPtsData > fieldEvalData
boost::shared_ptr< MatrixDouble > x2FieldPtr
boost::shared_ptr< PostProcEle > postProc

◆ Monitor() [2/8]

Monitor::Monitor ( SmartPetscObj< DM > dm,
boost::shared_ptr< PostProcEle > post_proc )
inline

Definition at line 68 of file heat_equation.cpp.

69 : dM(dm), postProc(post_proc){};

◆ Monitor() [3/8]

Monitor::Monitor ( SmartPetscObj< DM > dm,
boost::shared_ptr< PostProcEle > post_proc )
inline

Definition at line 65 of file wave_equation.cpp.

66 : dM(dm), postProc(post_proc){};

◆ Monitor() [4/8]

Monitor::Monitor ( SmartPetscObj< DM > dm,
std::pair< boost::shared_ptr< PostProcEleDomain >, boost::shared_ptr< PostProcEleBdy > > pair_post_proc,
boost::shared_ptr< DomainEle > reaction_fe )
inline

Definition at line 7 of file Monitor.hpp.

12 : dM(dm), postProcDomain(pair_post_proc.first),
13 postProcSkin(pair_post_proc.second), reactionFe(reaction_fe) {
14 CHK_THROW_MESSAGE(PetscOptionsGetInt(PETSC_NULLPTR, "", "-save_every",
15 &saveEvery, PETSC_NULLPTR),
16 "Cannot get -save_every option");
17 }
#define CHK_THROW_MESSAGE(err, msg)
Check and throw MoFEM exception.
PetscErrorCode PetscOptionsGetInt(PetscOptions *, const char pre[], const char name[], PetscInt *ivalue, PetscBool *set)
PetscInt saveEvery
Definition Monitor.hpp:26
boost::shared_ptr< DomainEle > reactionFe
Definition Monitor.hpp:25
boost::shared_ptr< PostProcEleDomain > postProcDomain
Definition Monitor.hpp:23
boost::shared_ptr< PostProcEleBdy > postProcSkin
Definition Monitor.hpp:24

◆ Monitor() [5/8]

Monitor::Monitor ( SmartPetscObj< DM > dm,
boost::shared_ptr< PostProcEle > post_proc )
inline

Definition at line 309 of file nonlinear_dynamic_elastic.cpp.

310 : dM(dm), postProc(post_proc){};

◆ Monitor() [6/8]

Monitor::Monitor ( SmartPetscObj< DM > dm,
boost::shared_ptr< PostProcEle > post_proc )
inline

Definition at line 731 of file shallow_wave.cpp.

732 : dM(dm), postProc(post_proc){};

◆ Monitor() [7/8]

Monitor::Monitor ( SmartPetscObj< DM > dm,
boost::shared_ptr< moab::Core > post_proc_mesh,
boost::shared_ptr< PostProcEleDomainCont > post_proc,
boost::shared_ptr< PostProcEleBdyCont > post_proc_edge,
std::pair< boost::shared_ptr< BoundaryEle >, boost::shared_ptr< VectorDouble > > p )
inline

Definition at line 1989 of file free_surface.cpp.

1995 : dM(dm), postProcMesh(post_proc_mesh), postProc(post_proc),
1996 postProcEdge(post_proc_edge), liftFE(p.first), liftVec(p.second) {}
boost::shared_ptr< moab::Core > postProcMesh
boost::shared_ptr< VectorDouble > liftVec
boost::shared_ptr< BoundaryEle > liftFE
boost::shared_ptr< PostProcEleBdyCont > postProcEdge

◆ Monitor() [8/8]

Monitor::Monitor ( SmartPetscObj< DM > dm,
std::pair< boost::shared_ptr< PostProcEleDomain >, boost::shared_ptr< PostProcEleBdy > > pair_post_proc_fe,
boost::shared_ptr< DomainEle > reaction_fe,
std::vector< boost::shared_ptr< ScalingMethod > > smv )
inline

Definition at line 603 of file adolc_plasticity.cpp.

609 : dM(dm), reactionFE(reaction_fe), vecOfTimeScalingMethods(smv) {
611 domainPostProcFe = pair_post_proc_fe.first;
612 skinPostProcFe = pair_post_proc_fe.second;
613
614 MoFEM::Interface *m_field_ptr;
615 CHKERR DMoFEMGetInterfacePtr(dM, &m_field_ptr);
616#ifdef ADD_CONTACT
617 auto get_integrate_traction = [&]() {
618 auto integrate_traction = boost::make_shared<BoundaryEle>(*m_field_ptr);
619 auto common_data_ptr = boost::make_shared<ContactOps::CommonData>();
622 integrate_traction->getOpPtrVector(), {HDIV}, "GEOMETRY")),
623 "Apply transform");
624 // We have to integrate on curved face geometry, thus integration weight
625 // have to adjusted.
626 integrate_traction->getOpPtrVector().push_back(
628 integrate_traction->getRuleHook = [](int, int, int approx_order) {
629 return 2 * approx_order + 2 - 1;
630 };
631
635 integrate_traction->getOpPtrVector(), "SIGMA", 0)),
636 "push operators to calculate traction");
637
638 return integrate_traction;
639 };
640
641 integrateTraction = get_integrate_traction();
642#endif
643 }
constexpr int SPACE_DIM
#define CHKERR
Inline error check.
auto createDMVector(DM dm)
Get smart vector from DM.
Definition DMMoFEM.hpp:1102
PetscErrorCode DMoFEMGetInterfacePtr(DM dm, MoFEM::Interface **m_field_ptr)
Get pointer to MoFEM::Interface.
Definition DMMoFEM.cpp:410
MoFEMErrorCode opFactoryCalculateTraction(boost::ptr_deque< ForcesAndSourcesCore::UserDataOperator > &pip, std::string sigma, bool is_axisymmetric=false)
static constexpr int approx_order
Add operators pushing bases from local to physical configuration.
Deprecated interface functions.
boost::shared_ptr< PostProcEleDomain > domainPostProcFe
boost::shared_ptr< BoundaryEle > integrateTraction
boost::shared_ptr< FEMethod > reactionFE
VecOfTimeScalingMethods vecOfTimeScalingMethods
boost::shared_ptr< PostProcEleBdy > skinPostProcFe
SmartPetscObj< Vec > fRes

Member Function Documentation

◆ operator()() [1/2]

MoFEMErrorCode Monitor::operator() ( )
inlinevirtual

function is run for every finite element

It is used to calculate element local matrices and assembly. It can be used for post-processing.

Reimplemented from MoFEM::BasicMethod.

Examples
heat_equation.cpp, and wave_equation.cpp.

Definition at line 72 of file heat_equation.cpp.

72{ return 0; }

◆ operator()() [2/2]

MoFEMErrorCode Monitor::operator() ( )
inlinevirtual

function is run for every finite element

It is used to calculate element local matrices and assembly. It can be used for post-processing.

Reimplemented from MoFEM::BasicMethod.

Definition at line 69 of file wave_equation.cpp.

69{ return 0; }

◆ postProcess() [1/8]

MoFEMErrorCode Monitor::postProcess ( )
inlinevirtual

function is run at the end of loop

It is used to assembly matrices and vectors, calculating global variables, f.e. total internal energy, ect.

Iterating over dofs: Example1 iterating over dofs in row by name of the field for(IT_GET_FEROW_BY_NAME_DOFS_FOR_LOOP(this,"DISPLACEMENT",it)) { ... }

Reimplemented from MoFEM::BasicMethod.

Examples
adolc_plasticity.cpp, dynamic_first_order_con_law.cpp, free_surface.cpp, heat_equation.cpp, shallow_wave.cpp, and wave_equation.cpp.

Definition at line 789 of file dynamic_first_order_con_law.cpp.

789 {
791
792 auto *simple = mField.getInterface<Simple>();
793
795 ->evalFEAtThePoint<SPACE_DIM>(
796 fieldEvalCoords.data(), 1e-12, simple->getProblemName(),
797 simple->getDomainFEName(), fieldEvalData, mField.get_comm_rank(),
799
800 if (velocityFieldPtr->size1()) {
802 auto t_x2_field = getFTensor1FromMat<SPACE_DIM>(*x2FieldPtr);
804
805 double u_x = t_x2_field(0) - t_geom(0);
806 double u_y = t_x2_field(1) - t_geom(1);
807 double u_z = t_x2_field(2) - t_geom(2);
808
809 MOFEM_LOG("SYNC", Sev::inform)
810 << "Velocities x: " << t_vel(0) << " y: " << t_vel(1)
811 << " z: " << t_vel(2) << "\n";
812 MOFEM_LOG("SYNC", Sev::inform) << "Displacement x: " << u_x
813 << " y: " << u_y << " z: " << u_z << "\n";
814 }
815
817 std::regex((boost::format("%s(.*)") % "Data_Vertex").str()))) {
818 Range ents;
819 mField.get_moab().get_entities_by_dimension(m->getMeshset(), 0, ents,
820 true);
821 auto print_vets = [](boost::shared_ptr<FieldEntity> ent_ptr) {
823 if (!(ent_ptr->getPStatus() & PSTATUS_NOT_OWNED)) {
824 MOFEM_LOG("SYNC", Sev::inform)
825 << "Velocities: " << ent_ptr->getEntFieldData()[0] << " "
826 << ent_ptr->getEntFieldData()[1] << " "
827 << ent_ptr->getEntFieldData()[2] << "\n";
828 }
830 };
831 CHKERR mField.getInterface<FieldBlas>()->fieldLambdaOnEntities(
832 print_vets, "V", &ents);
833 }
835
836 PetscBool print_volume = PETSC_FALSE;
837 CHKERR PetscOptionsGetBool(PETSC_NULLPTR, "", "-print_volume", &print_volume,
838 PETSC_NULLPTR);
839
840 PetscBool print_skin = PETSC_FALSE;
841 CHKERR PetscOptionsGetBool(PETSC_NULLPTR, "", "-print_skin", &print_skin,
842 PETSC_NULLPTR);
843
844 int save_every_nth_step = 1;
845 CHKERR PetscOptionsGetInt(PETSC_NULLPTR, "", "-save_step",
846 &save_every_nth_step, PETSC_NULLPTR);
847 if (ts_step % save_every_nth_step == 0) {
848
849 if (print_volume) {
851 CHKERR postProc->writeFile(
852 "out_step_" + boost::lexical_cast<std::string>(ts_step) + ".h5m");
853 }
854
855 if (print_skin) {
857 CHKERR postProcBdy->writeFile(
858 "out_boundary_" + boost::lexical_cast<std::string>(ts_step) +
859 ".h5m");
860 }
861 }
863 }
#define MOFEM_LOG_SEVERITY_SYNC(comm, severity)
Synchronise "SYNC" on curtain severity level.
void simple(double P1[], double P2[], double P3[], double c[], const int N)
Definition acoustic.cpp:69
@ QUIET
@ MF_EXIST
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[], MoFEM::FEMethod *method, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr())
Executes FEMethod for finite elements in DM.
Definition DMMoFEM.cpp:576
#define MOFEM_LOG(channel, severity)
Log.
MoFEMErrorCode getCubitMeshsetPtr(const int ms_id, const CubitBCType cubit_bc_type, const CubitMeshSets **cubit_meshset_ptr) const
get cubit meshset
PetscErrorCode PetscOptionsGetBool(PetscOptions *, const char pre[], const char name[], PetscBool *bval, PetscBool *set)
FTensor::Tensor1< FTensor::PackPtr< T *, S >, Tensor_Dim > getFTensor1FromMat(ublas::matrix< T, L, A > &data)
Get tensor rank 1 (vector) form data matrix.
int save_every_nth_step
FTensor::Index< 'm', 3 > m
boost::weak_ptr< CacheTuple > getCacheWeakPtr() const
Get the cache weak ptr object.
virtual moab::Interface & get_moab()=0
virtual MPI_Comm & get_comm() const =0
virtual int get_comm_rank() const =0
Basic algebra on fields.
Definition FieldBlas.hpp:21
Field evaluator interface.
Interface for managing meshsets containing materials and boundary conditions.
Simple interface for fast problem set-up.
Definition Simple.hpp:27
PetscInt ts_step
time step number
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface reference to pointer of interface.

◆ postProcess() [2/8]

MoFEMErrorCode Monitor::postProcess ( )
inlinevirtual

function is run at the end of loop

It is used to assembly matrices and vectors, calculating global variables, f.e. total internal energy, ect.

Iterating over dofs: Example1 iterating over dofs in row by name of the field for(IT_GET_FEROW_BY_NAME_DOFS_FOR_LOOP(this,"DISPLACEMENT",it)) { ... }

Reimplemented from MoFEM::BasicMethod.

Definition at line 76 of file heat_equation.cpp.

76 {
78 if (ts_step % saveEveryNthStep == 0) {
80 CHKERR postProc->writeFile(
81 "out_level_" + boost::lexical_cast<std::string>(ts_step) + ".h5m");
82 }
84 }
static constexpr int saveEveryNthStep

◆ postProcess() [3/8]

MoFEMErrorCode Monitor::postProcess ( )
inlinevirtual

function is run at the end of loop

It is used to assembly matrices and vectors, calculating global variables, f.e. total internal energy, ect.

Iterating over dofs: Example1 iterating over dofs in row by name of the field for(IT_GET_FEROW_BY_NAME_DOFS_FOR_LOOP(this,"DISPLACEMENT",it)) { ... }

Reimplemented from MoFEM::BasicMethod.

Definition at line 73 of file wave_equation.cpp.

73 {
75 if (ts_step % saveEveryNthStep == 0) {
77 CHKERR postProc->writeFile(
78 "out_level_" + boost::lexical_cast<std::string>(ts_step) + ".h5m");
79 }
81 }

◆ postProcess() [4/8]

MoFEMErrorCode Monitor::postProcess ( )
virtual

postProcess

Reimplemented from MoFEM::BasicMethod.

◆ postProcess() [5/8]

MoFEMErrorCode Monitor::postProcess ( )
inlinevirtual

function is run at the end of loop

It is used to assembly matrices and vectors, calculating global variables, f.e. total internal energy, ect.

Iterating over dofs: Example1 iterating over dofs in row by name of the field for(IT_GET_FEROW_BY_NAME_DOFS_FOR_LOOP(this,"DISPLACEMENT",it)) { ... }

Reimplemented from MoFEM::BasicMethod.

Definition at line 311 of file nonlinear_dynamic_elastic.cpp.

311 {
313 constexpr int save_every_nth_step = 1;
314 if (ts_step % save_every_nth_step == 0) {
316 CHKERR postProc->writeFile(
317 "out_step_" + boost::lexical_cast<std::string>(ts_step) + ".h5m");
318 }
320 }

◆ postProcess() [6/8]

MoFEMErrorCode Monitor::postProcess ( )
inlinevirtual

function is run at the end of loop

It is used to assembly matrices and vectors, calculating global variables, f.e. total internal energy, ect.

Iterating over dofs: Example1 iterating over dofs in row by name of the field for(IT_GET_FEROW_BY_NAME_DOFS_FOR_LOOP(this,"DISPLACEMENT",it)) { ... }

Reimplemented from MoFEM::BasicMethod.

Definition at line 733 of file shallow_wave.cpp.

733 {
735 constexpr int save_every_nth_step = 50;
736 if (ts_step % save_every_nth_step == 0) {
738 CHKERR postProc->writeFile(
739 "out_step_" + boost::lexical_cast<std::string>(ts_step) + ".h5m");
740 MOFEM_LOG("SW", Sev::verbose)
741 << "writing vector in binary to vector.dat ...";
742 PetscViewer viewer;
743 PetscViewerBinaryOpen(PETSC_COMM_WORLD, "vector.dat", FILE_MODE_WRITE,
744 &viewer);
745 VecView(ts_u, viewer);
746 PetscViewerDestroy(&viewer);
747 }
749 }
Vec & ts_u
state vector

◆ postProcess() [7/8]

MoFEMErrorCode Monitor::postProcess ( )
inlinevirtual

function is run at the end of loop

It is used to assembly matrices and vectors, calculating global variables, f.e. total internal energy, ect.

Iterating over dofs: Example1 iterating over dofs in row by name of the field for(IT_GET_FEROW_BY_NAME_DOFS_FOR_LOOP(this,"DISPLACEMENT",it)) { ... }

Reimplemented from MoFEM::BasicMethod.

Definition at line 1997 of file free_surface.cpp.

1997 {
1999
2000 MOFEM_LOG("FS", Sev::verbose) << "Monitor";
2001 constexpr int save_every_nth_step = 1;
2002 if (ts_step % save_every_nth_step == 0) {
2003 MOFEM_LOG("FS", Sev::verbose) << "Mesh pre proc";
2004 MoFEM::Interface *m_field_ptr;
2005 CHKERR DMoFEMGetInterfacePtr(dM, &m_field_ptr);
2006 auto post_proc_begin =
2007 boost::make_shared<PostProcBrokenMeshInMoabBaseBegin>(*m_field_ptr,
2008 postProcMesh);
2009 auto post_proc_end = boost::make_shared<PostProcBrokenMeshInMoabBaseEnd>(
2010 *m_field_ptr, postProcMesh);
2011 CHKERR DMoFEMPreProcessFiniteElements(dM, post_proc_begin->getFEMethod());
2013 this->getCacheWeakPtr());
2015 this->getCacheWeakPtr());
2016 CHKERR DMoFEMPostProcessFiniteElements(dM, post_proc_end->getFEMethod());
2017 CHKERR post_proc_end->writeFile(
2018 "out_step_" + boost::lexical_cast<std::string>(ts_step) + ".h5m");
2019 MOFEM_LOG("FS", Sev::verbose) << "Mesh pre proc done";
2020 }
2021
2022 liftVec->resize(SPACE_DIM, false);
2023 liftVec->clear();
2025 MPI_Allreduce(MPI_IN_PLACE, &(*liftVec)[0], SPACE_DIM, MPI_DOUBLE, MPI_SUM,
2026 MPI_COMM_WORLD);
2027 MOFEM_LOG("FS", Sev::inform)
2028 << "Step " << ts_step << " time " << ts_t
2029 << " lift vec x: " << (*liftVec)[0] << " y: " << (*liftVec)[1];
2030
2032 }
PetscErrorCode DMoFEMPostProcessFiniteElements(DM dm, MoFEM::FEMethod *method)
execute finite element method for each element in dm (problem)
Definition DMMoFEM.cpp:546
PetscErrorCode DMoFEMPreProcessFiniteElements(DM dm, MoFEM::FEMethod *method)
execute finite element method for each element in dm (problem)
Definition DMMoFEM.cpp:536
PetscReal ts_t
time

◆ postProcess() [8/8]

MoFEMErrorCode Monitor::postProcess ( )
inlinevirtual

function is run at the end of loop

It is used to assembly matrices and vectors, calculating global variables, f.e. total internal energy, ect.

Iterating over dofs: Example1 iterating over dofs in row by name of the field for(IT_GET_FEROW_BY_NAME_DOFS_FOR_LOOP(this,"DISPLACEMENT",it)) { ... }

Reimplemented from MoFEM::BasicMethod.

Definition at line 644 of file adolc_plasticity.cpp.

644 {
646 int save_every_nth_step = 1;
647 CHKERR PetscOptionsGetInt(PETSC_NULLPTR, "", "-save_every_nth_step",
648 &save_every_nth_step, PETSC_NULLPTR);
649#ifdef ADD_CONTACT
650 auto print_traction = [&](const std::string msg) {
652 MoFEM::Interface *m_field_ptr;
653 CHKERR DMoFEMGetInterfacePtr(dM, &m_field_ptr);
654 if (!m_field_ptr->get_comm_rank()) {
655 const double *t_ptr;
656 CHKERR VecGetArrayRead(ContactOps::CommonData::totalTraction, &t_ptr);
657 MOFEM_LOG_C("CONTACT", Sev::inform, "%s time %3.4e %3.4e %3.4e %3.4e",
658 msg.c_str(), ts_t, t_ptr[0], t_ptr[1], t_ptr[2]);
660 &t_ptr);
661 }
663 };
664#endif
665
666 auto make_vtk = [&]() {
668 if (domainPostProcFe) {
671 CHKERR domainPostProcFe->writeFile(
672 "out_plastic_" + boost::lexical_cast<std::string>(ts_step) +
673 ".h5m");
674 }
675 if (skinPostProcFe) {
678 CHKERR skinPostProcFe->writeFile(
679 "out_skin_plastic_" + boost::lexical_cast<std::string>(ts_step) +
680 ".h5m");
681 }
683 };
684
685 if (!(ts_step % save_every_nth_step)) {
686 CHKERR make_vtk();
687 }
688 if (reactionFE) {
689 CHKERR VecZeroEntries(fRes);
690 reactionFE->f = fRes;
692 CHKERR VecAssemblyBegin(fRes);
693 CHKERR VecAssemblyEnd(fRes);
694 CHKERR VecGhostUpdateBegin(fRes, ADD_VALUES, SCATTER_REVERSE);
695 CHKERR VecGhostUpdateEnd(fRes, ADD_VALUES, SCATTER_REVERSE);
696
697 MoFEM::Interface *m_field_ptr;
698 CHKERR DMoFEMGetInterfacePtr(dM, &m_field_ptr);
700 *m_field_ptr, reactionFE, fRes)();
701
702 double nrm;
703 CHKERR VecNorm(fRes, NORM_2, &nrm);
704 MOFEM_LOG("PlasticPrb", Sev::verbose)
705 << "Residual norm " << nrm << " at time step " << ts_step;
706 }
707
708#ifdef ADD_CONTACT
709 auto calculate_traction = [&] {
716 };
717#endif
718
719 auto get_min_max_displacement = [&]() {
721 MoFEM::Interface *m_field_ptr;
722 CHKERR DMoFEMGetInterfacePtr(dM, &m_field_ptr);
723
724 std::array<double, 4> a_min = {DBL_MAX, DBL_MAX, DBL_MAX, 0};
725 std::array<double, 4> a_max = {-DBL_MAX, -DBL_MAX, -DBL_MAX, 0};
726
727 auto get_min_max = [&](boost::shared_ptr<FieldEntity> field_entity_ptr) {
729 int d = 0;
730 for (auto v : field_entity_ptr->getEntFieldData()) {
731 a_min[d] = std::min(a_min[d], v);
732 a_max[d] = std::max(a_max[d], v);
733 ++d;
734 }
736 };
737
738 a_min[SPACE_DIM] = 0;
739 a_max[SPACE_DIM] = 0;
740
741 Range verts;
742 CHKERR m_field_ptr->get_field_entities_by_type("U", MBVERTEX, verts);
743 CHKERR m_field_ptr->getInterface<FieldBlas>()->fieldLambdaOnEntities(
744 get_min_max, "U", &verts);
745
746 auto mpi_reduce = [&](auto &a, auto op) {
747 std::array<double, 3> a_mpi = {0, 0, 0};
748 MPI_Allreduce(a.data(), a_mpi.data(), 3, MPI_DOUBLE, op,
749 m_field_ptr->get_comm());
750 return a_mpi;
751 };
752
753 auto a_min_mpi = mpi_reduce(a_min, MPI_MIN);
754 auto a_max_mpi = mpi_reduce(a_max, MPI_MAX);
755
756 MOFEM_LOG("PlasticPrb", Sev::inform)
757 << "Min displacement " << a_min_mpi[0] << " " << a_min_mpi[1] << " "
758 << a_min_mpi[2];
759 MOFEM_LOG("PlasticPrb", Sev::inform)
760 << "Max displacement " << a_max_mpi[0] << " " << a_max_mpi[1] << " "
761 << a_max_mpi[2];
763 };
764 CHKERR get_min_max_displacement();
765#ifdef ADD_CONTACT
766 CHKERR calculate_traction();
767 CHKERR print_traction("Contact force");
768#endif
769 double scale = 1;
770 for (auto s : vecOfTimeScalingMethods) {
771 scale *= s->getScale(this->ts_t);
772 }
773
774 MOFEM_LOG("PlasticPrb", Sev::inform)
775 << "Time: " << this->ts_t << " scale: " << scale;
776
778 }
#define MOFEM_LOG_C(channel, severity, format,...)
double scale
constexpr double a
virtual MoFEMErrorCode get_field_entities_by_type(const std::string name, EntityType type, Range &ents) const =0
get entities in the field by type
const double v
phase velocity of light in medium (cm/ns)
static SmartPetscObj< Vec > totalTraction
Class (Function) to calculate residual side diagonal.
Definition Essential.hpp:49

◆ preProcess() [1/2]

MoFEMErrorCode Monitor::preProcess ( )
inlinevirtual

function is run at the beginning of loop

It is used to zeroing matrices and vectors, calculation of shape functions on reference element, preprocessing boundary conditions, etc.

Reimplemented from MoFEM::BasicMethod.

Examples
heat_equation.cpp, and wave_equation.cpp.

Definition at line 71 of file heat_equation.cpp.

71{ return 0; }

◆ preProcess() [2/2]

MoFEMErrorCode Monitor::preProcess ( )
inlinevirtual

function is run at the beginning of loop

It is used to zeroing matrices and vectors, calculation of shape functions on reference element, preprocessing boundary conditions, etc.

Reimplemented from MoFEM::BasicMethod.

Definition at line 68 of file wave_equation.cpp.

68{ return 0; }

Member Data Documentation

◆ dM

◆ domainPostProcFe

boost::shared_ptr<PostProcEleDomain> Monitor::domainPostProcFe
private
Examples
adolc_plasticity.cpp.

Definition at line 782 of file adolc_plasticity.cpp.

◆ fieldEvalCoords

std::array<double, 3> Monitor::fieldEvalCoords
private

◆ fieldEvalData

boost::shared_ptr<SetPtsData> Monitor::fieldEvalData
private

◆ fRes

SmartPetscObj<Vec> Monitor::fRes
private
Examples
adolc_plasticity.cpp.

Definition at line 788 of file adolc_plasticity.cpp.

◆ geometryFieldPtr

boost::shared_ptr<MatrixDouble> Monitor::geometryFieldPtr
private

◆ integrateTraction

boost::shared_ptr<BoundaryEle> Monitor::integrateTraction
private
Examples
adolc_plasticity.cpp.

Definition at line 786 of file adolc_plasticity.cpp.

◆ liftFE

boost::shared_ptr<BoundaryEle> Monitor::liftFE
private
Examples
free_surface.cpp.

Definition at line 2039 of file free_surface.cpp.

◆ liftVec

boost::shared_ptr<VectorDouble> Monitor::liftVec
private
Examples
free_surface.cpp.

Definition at line 2040 of file free_surface.cpp.

◆ mField

MoFEM::Interface& Monitor::mField

◆ postProc [1/2]

boost::shared_ptr< PostProcEle > Monitor::postProc
private

◆ postProc [2/2]

boost::shared_ptr<PostProcEleDomainCont> Monitor::postProc
private

Definition at line 2037 of file free_surface.cpp.

◆ postProcBdy

boost::shared_ptr<PostProcFaceEle> Monitor::postProcBdy
private

◆ postProcDomain

boost::shared_ptr<PostProcEleDomain> Monitor::postProcDomain
private

Definition at line 23 of file Monitor.hpp.

◆ postProcEdge

boost::shared_ptr<PostProcEleBdyCont> Monitor::postProcEdge
private
Examples
free_surface.cpp.

Definition at line 2038 of file free_surface.cpp.

◆ postProcMesh

boost::shared_ptr<moab::Core> Monitor::postProcMesh
private
Examples
free_surface.cpp.

Definition at line 2036 of file free_surface.cpp.

◆ postProcSkin

boost::shared_ptr<PostProcEleBdy> Monitor::postProcSkin
private

Definition at line 24 of file Monitor.hpp.

◆ reactionFE

boost::shared_ptr<FEMethod> Monitor::reactionFE
private
Examples
adolc_plasticity.cpp.

Definition at line 784 of file adolc_plasticity.cpp.

◆ reactionFe

boost::shared_ptr<DomainEle> Monitor::reactionFe
private

Definition at line 25 of file Monitor.hpp.

◆ saveEvery

PetscInt Monitor::saveEvery = 1
private

Definition at line 26 of file Monitor.hpp.

◆ saveEveryNthStep

static constexpr int Monitor::saveEveryNthStep = 1
staticconstexpr
Examples
heat_equation.cpp, and wave_equation.cpp.

Definition at line 74 of file heat_equation.cpp.

◆ skinPostProcFe

boost::shared_ptr<PostProcEleBdy> Monitor::skinPostProcFe
private
Examples
adolc_plasticity.cpp.

Definition at line 783 of file adolc_plasticity.cpp.

◆ vecOfTimeScalingMethods

VecOfTimeScalingMethods Monitor::vecOfTimeScalingMethods
private
Examples
adolc_plasticity.cpp.

Definition at line 789 of file adolc_plasticity.cpp.

◆ velocityFieldPtr

boost::shared_ptr<MatrixDouble> Monitor::velocityFieldPtr
private

◆ x2FieldPtr

boost::shared_ptr<MatrixDouble> Monitor::x2FieldPtr
private

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