v0.13.0
Public Member Functions | Private Attributes | List of all members
PlasticOps::Monitor Struct Reference

#include <users_modules/tutorials/adv-0/src/PlasticOpsMonitor.hpp>

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

Public Member Functions

 Monitor (SmartPetscObj< DM > &dm, boost::shared_ptr< PostProcEle > &post_proc_fe, boost::shared_ptr< DomainEle > &reaction_fe, std::tuple< SmartPetscObj< Vec >, SmartPetscObj< VecScatter >> ux_scatter, std::tuple< SmartPetscObj< Vec >, SmartPetscObj< VecScatter >> uy_scatter, std::tuple< SmartPetscObj< Vec >, SmartPetscObj< VecScatter >> uz_scatter)
 
MoFEMErrorCode preProcess ()
 
MoFEMErrorCode operator() ()
 
MoFEMErrorCode postProcess ()
 

Private Attributes

SmartPetscObj< DM > dM
 
boost::shared_ptr< PostProcElepostProcFe
 
boost::shared_ptr< DomainElereactionFe
 
std::tuple< SmartPetscObj< Vec >, SmartPetscObj< VecScatter > > uXScatter
 
std::tuple< SmartPetscObj< Vec >, SmartPetscObj< VecScatter > > uYScatter
 
std::tuple< SmartPetscObj< Vec >, SmartPetscObj< VecScatter > > uZScatter
 

Detailed Description

Definition at line 112 of file PlasticOpsMonitor.hpp.

Constructor & Destructor Documentation

◆ Monitor()

PlasticOps::Monitor::Monitor ( SmartPetscObj< DM > &  dm,
boost::shared_ptr< PostProcEle > &  post_proc_fe,
boost::shared_ptr< DomainEle > &  reaction_fe,
std::tuple< SmartPetscObj< Vec >, SmartPetscObj< VecScatter >>  ux_scatter,
std::tuple< SmartPetscObj< Vec >, SmartPetscObj< VecScatter >>  uy_scatter,
std::tuple< SmartPetscObj< Vec >, SmartPetscObj< VecScatter >>  uz_scatter 
)
Examples
PlasticOpsMonitor.hpp.

Definition at line 114 of file PlasticOpsMonitor.hpp.

119  : dM(dm), postProcFe(post_proc_fe), reactionFe(reaction_fe),
120  uXScatter(ux_scatter), uYScatter(uy_scatter), uZScatter(uz_scatter){};
SmartPetscObj< DM > dM
std::tuple< SmartPetscObj< Vec >, SmartPetscObj< VecScatter > > uXScatter
std::tuple< SmartPetscObj< Vec >, SmartPetscObj< VecScatter > > uZScatter
boost::shared_ptr< PostProcEle > postProcFe
std::tuple< SmartPetscObj< Vec >, SmartPetscObj< VecScatter > > uYScatter
boost::shared_ptr< DomainEle > reactionFe

Member Function Documentation

◆ operator()()

MoFEMErrorCode PlasticOps::Monitor::operator() ( )
Examples
PlasticOpsMonitor.hpp.

Definition at line 123 of file PlasticOpsMonitor.hpp.

123 { return 0; }

◆ postProcess()

MoFEMErrorCode PlasticOps::Monitor::postProcess ( )
Examples
PlasticOpsMonitor.hpp.

Definition at line 125 of file PlasticOpsMonitor.hpp.

125  {
127 
128  auto make_vtk = [&]() {
131  CHKERR postProcFe->writeFile(
132  "out_plastic_" + boost::lexical_cast<std::string>(ts_step) + ".h5m");
134  };
135 
136  auto calculate_reaction = [&]() {
138  auto r = smartCreateDMVector(dM);
139  reactionFe->f = r;
140  CHKERR VecZeroEntries(r);
142  CHKERR VecGhostUpdateBegin(r, ADD_VALUES, SCATTER_REVERSE);
143  CHKERR VecGhostUpdateEnd(r, ADD_VALUES, SCATTER_REVERSE);
144  CHKERR VecAssemblyBegin(r);
145  CHKERR VecAssemblyEnd(r);
146 
147  double sum;
148  CHKERR VecSum(r, &sum);
149  MOFEM_LOG_C("EXAMPLE", Sev::inform, "reaction time %3.4e %3.4e", ts_t,
150  sum);
151 
153  };
154 
155  auto print_max_min = [&](auto &tuple, const std::string msg) {
157  CHKERR VecScatterBegin(std::get<1>(tuple), ts_u, std::get<0>(tuple),
158  INSERT_VALUES, SCATTER_FORWARD);
159  CHKERR VecScatterEnd(std::get<1>(tuple), ts_u, std::get<0>(tuple),
160  INSERT_VALUES, SCATTER_FORWARD);
161  double max, min;
162  CHKERR VecMax(std::get<0>(tuple), PETSC_NULL, &max);
163  CHKERR VecMin(std::get<0>(tuple), PETSC_NULL, &min);
164  MOFEM_LOG_C("EXAMPLE", Sev::inform, "%s time %3.4e min %3.4e max %3.4e",
165  msg.c_str(), ts_t, min, max);
167  };
168 
169  CHKERR make_vtk();
170  CHKERR calculate_reaction();
171  CHKERR print_max_min(uXScatter, "Ux");
172  CHKERR print_max_min(uYScatter, "Uy");
173  if (SPACE_DIM == 3)
174  CHKERR print_max_min(uZScatter, "Uz");
175 
177  }
#define MOFEM_LOG_C(channel, severity, format,...)
Definition: LogManager.hpp:314
constexpr int SPACE_DIM
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
#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
auto smartCreateDMVector
Get smart vector from DM.
Definition: DMMoFEM.hpp:956
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:544
const double r
rate factor

◆ preProcess()

MoFEMErrorCode PlasticOps::Monitor::preProcess ( )
Examples
PlasticOpsMonitor.hpp.

Definition at line 122 of file PlasticOpsMonitor.hpp.

122 { return 0; }

Member Data Documentation

◆ dM

SmartPetscObj<DM> PlasticOps::Monitor::dM
private
Examples
PlasticOpsMonitor.hpp.

Definition at line 180 of file PlasticOpsMonitor.hpp.

◆ postProcFe

boost::shared_ptr<PostProcEle> PlasticOps::Monitor::postProcFe
private
Examples
PlasticOpsMonitor.hpp.

Definition at line 181 of file PlasticOpsMonitor.hpp.

◆ reactionFe

boost::shared_ptr<DomainEle> PlasticOps::Monitor::reactionFe
private
Examples
PlasticOpsMonitor.hpp.

Definition at line 182 of file PlasticOpsMonitor.hpp.

◆ uXScatter

std::tuple<SmartPetscObj<Vec>, SmartPetscObj<VecScatter> > PlasticOps::Monitor::uXScatter
private
Examples
PlasticOpsMonitor.hpp.

Definition at line 183 of file PlasticOpsMonitor.hpp.

◆ uYScatter

std::tuple<SmartPetscObj<Vec>, SmartPetscObj<VecScatter> > PlasticOps::Monitor::uYScatter
private
Examples
PlasticOpsMonitor.hpp.

Definition at line 184 of file PlasticOpsMonitor.hpp.

◆ uZScatter

std::tuple<SmartPetscObj<Vec>, SmartPetscObj<VecScatter> > PlasticOps::Monitor::uZScatter
private
Examples
PlasticOpsMonitor.hpp.

Definition at line 185 of file PlasticOpsMonitor.hpp.


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