v0.14.0
PlasticOpsMonitor.hpp
/** \file PlasticOpsMonitor.hpp
* \example PlasticOpsMonitor.hpp
*/
namespace PlasticOps {
template <int DIM> struct Monitor : public FEMethod {
SmartPetscObj<DM> dm,
std::pair<boost::shared_ptr<PostProcEle>,
boost::shared_ptr<SkinPostProcEle>>
pair_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,
std::array<double, DIM> pass_field_eval_coords,
boost::shared_ptr<SetPtsData> pass_field_eval_data,
boost::shared_ptr<std::map<std::string, boost::shared_ptr<VectorDouble>>>
scalar_field_ptrs,
boost::shared_ptr<std::map<std::string, boost::shared_ptr<MatrixDouble>>>
vec_field_ptrs,
boost::shared_ptr<std::map<std::string, boost::shared_ptr<MatrixDouble>>>
sym_tensor_field_ptrs,
boost::shared_ptr<std::map<std::string, boost::shared_ptr<MatrixDouble>>>
tensor_field_ptrs)
: dM(dm), reactionFe(reaction_fe), uXScatter(ux_scatter),
uYScatter(uy_scatter), uZScatter(uz_scatter),
fieldEvalCoords(pass_field_eval_coords),
fieldEvalData(pass_field_eval_data), scalarFieldPtrs(scalar_field_ptrs),
vecFieldPtrs(vec_field_ptrs), symTensorFieldPtrs(sym_tensor_field_ptrs),
tensorFieldPtrs(tensor_field_ptrs) {
postProcFe = pair_post_proc_fe.first;
skinPostProcFe = pair_post_proc_fe.second;
};
MoFEMErrorCode preProcess() { return 0; }
MoFEMErrorCode operator()() { return 0; }
private:
SmartPetscObj<DM> dM;
boost::shared_ptr<PostProcEle> postProcFe;
boost::shared_ptr<SkinPostProcEle> skinPostProcFe;
boost::shared_ptr<DomainEle> reactionFe;
std::tuple<SmartPetscObj<Vec>, SmartPetscObj<VecScatter>> uXScatter;
std::tuple<SmartPetscObj<Vec>, SmartPetscObj<VecScatter>> uYScatter;
std::tuple<SmartPetscObj<Vec>, SmartPetscObj<VecScatter>> uZScatter;
std::array<double, DIM> fieldEvalCoords;
boost::shared_ptr<SetPtsData> fieldEvalData;
boost::shared_ptr<std::map<std::string, boost::shared_ptr<VectorDouble>>>
boost::shared_ptr<std::map<std::string, boost::shared_ptr<MatrixDouble>>>
boost::shared_ptr<std::map<std::string, boost::shared_ptr<MatrixDouble>>>
boost::shared_ptr<std::map<std::string, boost::shared_ptr<MatrixDouble>>>
protected:
};
MoFEM::Interface *m_field_ptr;
CHKERR DMoFEMGetInterfacePtr(dM, &m_field_ptr);
auto *simple = m_field_ptr->getInterface<Simple>();
CHKERR m_field_ptr->getInterface<FieldEvaluatorInterface>()
->evalFEAtThePoint<DIM>(
fieldEvalCoords.data(), 1e-12, simple->getProblemName(),
simple->getDomainFEName(), fieldEvalData,
m_field_ptr->get_comm_rank(), m_field_ptr->get_comm_rank(),
getCacheWeakPtr().lock(), MF_EXIST, QUIET);
auto process_scalar_field =
[](const std::string label,
const boost::shared_ptr<VectorDouble> scalarFieldPtr) {
if (scalarFieldPtr->size()) {
auto t_scalar_holder = getFTensor0FromVec(*scalarFieldPtr);
MOFEM_LOG("SYNC", Sev::inform)
<< "For " << label << " field, " << t_scalar_holder;
}
};
auto process_vector_field =
[](const std::string label,
const boost::shared_ptr<MatrixDouble> vecFieldPtr) {
if (vecFieldPtr->size1()) {
auto t_vec_holder = getFTensor1FromMat<DIM>(*vecFieldPtr);
std::string vec_holder_string = "";
for (int i = 0; i < DIM; i++) {
vec_holder_string +=
" " + boost::lexical_cast<std::string>(t_vec_holder(i));
}
MOFEM_LOG("SYNC", Sev::inform)
<< "For " << label << " field," << vec_holder_string;
}
};
auto process_sym_tensor_field =
[](const std::string label,
const boost::shared_ptr<MatrixDouble> symTensorFieldPtr) {
if (symTensorFieldPtr->size1()) {
auto t_sym_tensor_holder =
getFTensor2SymmetricFromMat<DIM>(*symTensorFieldPtr);
std::string sym_tensor_holder_string = "";
for (int i = 0; i < DIM; i++) {
for (int j = i; j < DIM; j++) {
sym_tensor_holder_string +=
", entry " + boost::lexical_cast<std::string>(i) +
boost::lexical_cast<std::string>(j) + " = ";
sym_tensor_holder_string +=
boost::lexical_cast<std::string>(t_sym_tensor_holder(i, j));
}
}
MOFEM_LOG("SYNC", Sev::inform)
<< "For " << label << " field" << sym_tensor_holder_string;
}
};
auto process_tensor_field =
[](const std::string label,
const boost::shared_ptr<MatrixDouble> tensorFieldPtr) {
if (tensorFieldPtr->size1()) {
auto t_tensor_holder =
getFTensor2FromMat<DIM, DIM>(*tensorFieldPtr);
std::string tensor_holder_string = "";
for (int i = 0; i < DIM; i++) {
for (int j = 0; j < DIM; j++) {
tensor_holder_string +=
", entry " + boost::lexical_cast<std::string>(i) +
boost::lexical_cast<std::string>(j) + " = ";
tensor_holder_string +=
boost::lexical_cast<std::string>(t_tensor_holder(i, j));
}
}
MOFEM_LOG("SYNC", Sev::inform)
<< "For " << label << " field" << tensor_holder_string;
}
};
auto process_fields = [](auto &fieldPtr, auto process_field) {
for (const auto &pair : fieldPtr) {
const std::string &label = pair.first;
const boost::shared_ptr<MatrixDouble> &ptr = pair.second;
if (ptr) {
process_field(label, ptr);
}
}
};
if (scalarFieldPtrs) {
for (const auto &pair : *scalarFieldPtrs) {
const std::string &label = pair.first;
const boost::shared_ptr<VectorDouble> &ptr = pair.second;
if (ptr) {
process_scalar_field(label, ptr);
}
}
}
if (vecFieldPtrs) {
for (const auto &pair : *vecFieldPtrs) {
const std::string &label = pair.first;
const boost::shared_ptr<MatrixDouble> &ptr = pair.second;
if (ptr) {
process_vector_field(label, ptr);
}
}
}
if (symTensorFieldPtrs) {
for (const auto &pair : *symTensorFieldPtrs) {
const std::string &label = pair.first;
const boost::shared_ptr<MatrixDouble> &ptr = pair.second;
if (ptr) {
process_sym_tensor_field(label, ptr);
}
}
}
if (tensorFieldPtrs) {
for (const auto &pair : *tensorFieldPtrs) {
const std::string &label = pair.first;
const boost::shared_ptr<MatrixDouble> &ptr = pair.second;
if (ptr) {
process_tensor_field(label, ptr);
}
}
}
}
MOFEM_LOG_SEVERITY_SYNC(m_field_ptr->get_comm(), Sev::inform);
auto make_vtk = [&]() {
if (postProcFe) {
CHKERR DMoFEMLoopFiniteElements(dM, "dFE", postProcFe, getCacheWeakPtr());
CHKERR postProcFe->writeFile(
"out_plastic_" + boost::lexical_cast<std::string>(ts_step) + ".h5m");
}
if (skinPostProcFe) {
CHKERR DMoFEMLoopFiniteElements(dM, "bFE", skinPostProcFe,
getCacheWeakPtr());
CHKERR skinPostProcFe->writeFile(
"out_skin_plastic_" + boost::lexical_cast<std::string>(ts_step) +
".h5m");
}
};
auto calculate_reaction = [&]() {
auto r = createDMVector(dM);
reactionFe->f = r;
CHKERR VecZeroEntries(r);
CHKERR DMoFEMLoopFiniteElements(dM, "dFE", reactionFe, getCacheWeakPtr());
#ifndef NDEBUG
auto post_proc_residual = [&](auto dm, auto f_res, auto out_name) {
auto post_proc_fe =
boost::make_shared<PostProcBrokenMeshInMoab<DomainEle>>(*m_field_ptr);
using OpPPMap = OpPostProcMapInMoab<DIM, DIM>;
auto u_vec = boost::make_shared<MatrixDouble>();
post_proc_fe->getOpPtrVector().push_back(
new OpCalculateVectorFieldValues<DIM>("U", u_vec, f_res));
post_proc_fe->getOpPtrVector().push_back(
new OpPPMap(
post_proc_fe->getPostProcMesh(), post_proc_fe->getMapGaussPts(),
{},
{{"RES", u_vec}},
{}, {})
);
CHKERR DMoFEMLoopFiniteElements(dM, "dFE", post_proc_fe);
post_proc_fe->writeFile("res.h5m");
};
CHKERR post_proc_residual(dM, r, "reaction");
#endif // NDEBUG
};
auto print_max_min = [&](auto &tuple, const std::string msg) {
CHKERR VecScatterBegin(std::get<1>(tuple), ts_u, std::get<0>(tuple),
INSERT_VALUES, SCATTER_FORWARD);
CHKERR VecScatterEnd(std::get<1>(tuple), ts_u, std::get<0>(tuple),
INSERT_VALUES, SCATTER_FORWARD);
double max, min;
CHKERR VecMax(std::get<0>(tuple), PETSC_NULL, &max);
CHKERR VecMin(std::get<0>(tuple), PETSC_NULL, &min);
MOFEM_LOG_C("PLASTICITY", Sev::inform, "%s time %3.4e min %3.4e max %3.4e",
msg.c_str(), ts_t, min, max);
};
int se = 1;
CHKERR PetscOptionsGetInt(PETSC_NULL, "", "-save_every", &se, PETSC_NULL);
if (!(ts_step % se)) {
CHKERR make_vtk();
}
if (reactionFe)
CHKERR calculate_reaction();
CHKERR print_max_min(uXScatter, "Ux");
CHKERR print_max_min(uYScatter, "Uy");
if (DIM == 3)
CHKERR print_max_min(uZScatter, "Uz");
}
}; // namespace PlasticOps
MoFEM::UnknownInterface::getInterface
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface reference to pointer of interface.
Definition: UnknownInterface.hpp:93
PlasticOps::Monitor::fieldEvalData
boost::shared_ptr< SetPtsData > fieldEvalData
Definition: PlasticOpsMonitor.hpp:52
MOFEM_LOG_SEVERITY_SYNC
#define MOFEM_LOG_SEVERITY_SYNC(comm, severity)
Synchronise "SYNC" on curtain severity level.
Definition: LogManager.hpp:352
PlasticOps::Monitor::operator()
MoFEMErrorCode operator()()
Definition: PlasticOpsMonitor.hpp:41
PlasticOps::Monitor::preProcess
MoFEMErrorCode preProcess()
Definition: PlasticOpsMonitor.hpp:40
PlasticOps::Monitor::uYScatter
std::tuple< SmartPetscObj< Vec >, SmartPetscObj< VecScatter > > uYScatter
Definition: PlasticOpsMonitor.hpp:49
PlasticOps::Monitor::scalarFieldPtrs
boost::shared_ptr< std::map< std::string, boost::shared_ptr< VectorDouble > > > scalarFieldPtrs
Definition: PlasticOpsMonitor.hpp:54
MoFEM::CoreInterface::get_comm
virtual MPI_Comm & get_comm() const =0
do_eval_field
PetscBool do_eval_field
Evaluate field.
Definition: plastic.cpp:119
MoFEM::Exceptions::MoFEMErrorCode
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
PlasticOps::Monitor::symTensorFieldPtrs
boost::shared_ptr< std::map< std::string, boost::shared_ptr< MatrixDouble > > > symTensorFieldPtrs
Definition: PlasticOpsMonitor.hpp:58
MoFEM::CoreInterface::get_comm_rank
virtual int get_comm_rank() const =0
PlasticOps::Monitor::uXScatter
std::tuple< SmartPetscObj< Vec >, SmartPetscObj< VecScatter > > uXScatter
Definition: PlasticOpsMonitor.hpp:48
PlasticOps::Monitor::postProcFe
boost::shared_ptr< PostProcEle > postProcFe
Definition: PlasticOpsMonitor.hpp:45
sdf.r
int r
Definition: sdf.py:8
MoFEM::DeprecatedCoreInterface
Deprecated interface functions.
Definition: DeprecatedCoreInterface.hpp:16
PlasticOps::Monitor::skinPostProcFe
boost::shared_ptr< SkinPostProcEle > skinPostProcFe
Definition: PlasticOpsMonitor.hpp:46
FEMethod
CHKERR
#define CHKERR
Inline error check.
Definition: definitions.h:548
MoFEM::createDMVector
auto createDMVector(DM dm)
Get smart vector from DM.
Definition: DMMoFEM.hpp:1099
PlasticOps::Monitor::tensorFieldPtrs
boost::shared_ptr< std::map< std::string, boost::shared_ptr< MatrixDouble > > > tensorFieldPtrs
Definition: PlasticOpsMonitor.hpp:60
MOFEM_LOG_C
#define MOFEM_LOG_C(channel, severity, format,...)
Definition: LogManager.hpp:311
simple
void simple(double P1[], double P2[], double P3[], double c[], const int N)
Definition: acoustic.cpp:69
PlasticOps::Monitor::fieldEvalCoords
std::array< double, DIM > fieldEvalCoords
Definition: PlasticOpsMonitor.hpp:51
OpPPMap
OpPostProcMapInMoab< SPACE_DIM, SPACE_DIM > OpPPMap
Definition: photon_diffusion.cpp:29
MoFEM::getFTensor0FromVec
static auto getFTensor0FromVec(ublas::vector< T, A > &data)
Get tensor rank 0 (scalar) form data vector.
Definition: Templates.hpp:135
PlasticOps::Monitor::vecFieldPtrs
boost::shared_ptr< std::map< std::string, boost::shared_ptr< MatrixDouble > > > vecFieldPtrs
Definition: PlasticOpsMonitor.hpp:56
i
FTensor::Index< 'i', SPACE_DIM > i
Definition: hcurl_divergence_operator_2d.cpp:27
PlasticOps::Monitor::postProcess
MoFEMErrorCode postProcess()
Definition: PlasticOpsMonitor.hpp:66
PlasticOps
Definition: plastic.cpp:182
PlasticOps::Monitor::uZScatter
std::tuple< SmartPetscObj< Vec >, SmartPetscObj< VecScatter > > uZScatter
Definition: PlasticOpsMonitor.hpp:50
MOFEM_LOG
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:308
MoFEM::DMoFEMGetInterfacePtr
PetscErrorCode DMoFEMGetInterfacePtr(DM dm, MoFEM::Interface **m_field_ptr)
Get pointer to MoFEM::Interface.
Definition: DMMoFEM.cpp:414
j
FTensor::Index< 'j', 3 > j
Definition: matrix_function.cpp:19
Monitor
[Push operators to pipeline]
Definition: dynamic_first_order_con_law.cpp:774
PlasticOps::Monitor::Monitor
Monitor(SmartPetscObj< DM > dm, std::pair< boost::shared_ptr< PostProcEle >, boost::shared_ptr< SkinPostProcEle >> pair_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, std::array< double, DIM > pass_field_eval_coords, boost::shared_ptr< SetPtsData > pass_field_eval_data, boost::shared_ptr< std::map< std::string, boost::shared_ptr< VectorDouble >>> scalar_field_ptrs, boost::shared_ptr< std::map< std::string, boost::shared_ptr< MatrixDouble >>> vec_field_ptrs, boost::shared_ptr< std::map< std::string, boost::shared_ptr< MatrixDouble >>> sym_tensor_field_ptrs, boost::shared_ptr< std::map< std::string, boost::shared_ptr< MatrixDouble >>> tensor_field_ptrs)
Definition: PlasticOpsMonitor.hpp:11
QUIET
@ QUIET
Definition: definitions.h:221
PlasticOps::Monitor::dM
SmartPetscObj< DM > dM
Definition: PlasticOpsMonitor.hpp:44
MoFEM::DMoFEMLoopFiniteElements
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:586
MF_EXIST
@ MF_EXIST
Definition: definitions.h:113
MoFEM::PetscOptionsGetInt
PetscErrorCode PetscOptionsGetInt(PetscOptions *, const char pre[], const char name[], PetscInt *ivalue, PetscBool *set)
Definition: DeprecatedPetsc.hpp:142
MoFEMFunctionReturn
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
PlasticOps::Monitor::reactionFe
boost::shared_ptr< DomainEle > reactionFe
Definition: PlasticOpsMonitor.hpp:47
MoFEMFunctionBegin
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
MoFEM::OpPostProcMapInMoab
Post post-proc data at points from hash maps.
Definition: PostProcBrokenMeshInMoabBase.hpp:709