v0.13.0
PlasticOpsMonitor.hpp
Go to the documentation of this file.
1 /* This file is part of MoFEM.
2  * MoFEM is free software: you can redistribute it and/or modify it under
3  * the terms of the GNU Lesser General Public License as published by the
4  * Free Software Foundation, either version 3 of the License, or (at your
5  * option) any later version.
6  *
7  * MoFEM is distributed in the hope that it will be useful, but WITHOUT
8  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
9  * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
10  * License for more details.
11  *
12  * You should have received a copy of the GNU Lesser General Public
13  * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>. */
14 
15 /** \file PlasticOpsMonitor.hpp
16  * \example PlasticOpsMonitor.hpp
17  */
18 
19 namespace PlasticOps {
20 
22  const std::string field_name, moab::Interface &post_proc_mesh,
23  std::vector<EntityHandle> &map_gauss_pts,
24  boost::shared_ptr<CommonData> common_data_ptr)
25  : DomainEleOp(field_name, DomainEleOp::OPROW), postProcMesh(post_proc_mesh),
26  mapGaussPts(map_gauss_pts), commonDataPtr(common_data_ptr) {
27  // Opetor is only executed for vertices
28  std::fill(&doEntities[MBEDGE], &doEntities[MBMAXTYPE], false);
29 }
30 
31 //! [Postprocessing]
33  EntData &data) {
35 
36  std::array<double, 9> def;
37  std::fill(def.begin(), def.end(), 0);
38 
39  auto get_tag = [&](const std::string name, size_t size) {
40  Tag th;
41  CHKERR postProcMesh.tag_get_handle(name.c_str(), size, MB_TYPE_DOUBLE, th,
42  MB_TAG_CREAT | MB_TAG_SPARSE,
43  def.data());
44  return th;
45  };
46 
47  MatrixDouble3by3 mat(3, 3);
48 
49  auto set_matrix_3d = [&](auto &t) -> MatrixDouble3by3 & {
50  mat.clear();
51  for (size_t r = 0; r != SPACE_DIM; ++r)
52  for (size_t c = 0; c != SPACE_DIM; ++c)
53  mat(r, c) = t(r, c);
54  return mat;
55  };
56 
57  auto set_scalar = [&](auto t) -> MatrixDouble3by3 & {
58  mat.clear();
59  mat(0, 0) = t;
60  return mat;
61  };
62 
63  auto set_float_precision = [](const double x) {
64  if (std::abs(x) < std::numeric_limits<float>::epsilon())
65  return 0.;
66  else
67  return x;
68  };
69 
70  auto set_tag = [&](auto th, auto gg, MatrixDouble3by3 &mat) {
71  for (auto &v : mat.data())
72  v = set_float_precision(v);
73  return postProcMesh.tag_set_data(th, &mapGaussPts[gg], 1,
74  &*mat.data().begin());
75  };
76 
77  auto th_plastic_surface = get_tag("PLASTIC_SURFACE", 1);
78  auto th_hardening = get_tag("HARDENING", 1);
79  auto th_tau = get_tag("PLASTIC_MULTIPLIER", 1);
80  auto th_temperature = get_tag("TEMPERATURE", 1);
81  auto th_plastic_flow = get_tag("PLASTIC_FLOW", 9);
82  auto th_plastic_strain = get_tag("PLASTIC_STRAIN", 9);
83 
84  auto t_flow =
85  getFTensor2SymmetricFromMat<SPACE_DIM>(commonDataPtr->plasticFlow);
86  auto t_plastic_strain =
87  getFTensor2SymmetricFromMat<SPACE_DIM>(commonDataPtr->plasticStrain);
88  if (commonDataPtr->tempVal.size() != commonDataPtr->plasticSurface.size()) {
89  commonDataPtr->tempVal.resize(commonDataPtr->plasticSurface.size(), 0);
90  commonDataPtr->tempVal.clear();
91  }
92 
93  size_t gg = 0;
94  for (int gg = 0; gg != commonDataPtr->plasticSurface.size(); ++gg) {
95  const double temp = (commonDataPtr->tempVal)[gg];
96  const double tau = (commonDataPtr->plasticTau)[gg];
97  const double f = (commonDataPtr->plasticSurface)[gg];
98  const double h = hardening(tau, temp);
99  CHKERR set_tag(th_plastic_surface, gg, set_scalar(f));
100  CHKERR set_tag(th_hardening, gg, set_scalar(h));
101  CHKERR set_tag(th_tau, gg, set_scalar(tau));
102  CHKERR set_tag(th_temperature, gg, set_scalar(temp));
103  CHKERR set_tag(th_plastic_flow, gg, set_matrix_3d(t_flow));
104  CHKERR set_tag(th_plastic_strain, gg, set_matrix_3d(t_plastic_strain));
105  ++t_flow;
106  ++t_plastic_strain;
107  }
108 
110 }
111 
112 struct Monitor : public FEMethod {
113 
114  Monitor(SmartPetscObj<DM> &dm, boost::shared_ptr<PostProcEle> &post_proc_fe,
115  boost::shared_ptr<DomainEle> &reaction_fe,
116  std::tuple<SmartPetscObj<Vec>, SmartPetscObj<VecScatter>> ux_scatter,
117  std::tuple<SmartPetscObj<Vec>, SmartPetscObj<VecScatter>> uy_scatter,
118  std::tuple<SmartPetscObj<Vec>, SmartPetscObj<VecScatter>> uz_scatter)
119  : dM(dm), postProcFe(post_proc_fe), reactionFe(reaction_fe),
120  uXScatter(ux_scatter), uYScatter(uy_scatter), uZScatter(uz_scatter){};
121 
122  MoFEMErrorCode preProcess() { return 0; }
123  MoFEMErrorCode operator()() { return 0; }
124 
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  }
178 
179 private:
180  SmartPetscObj<DM> dM;
181  boost::shared_ptr<PostProcEle> postProcFe;
182  boost::shared_ptr<DomainEle> reactionFe;
183  std::tuple<SmartPetscObj<Vec>, SmartPetscObj<VecScatter>> uXScatter;
184  std::tuple<SmartPetscObj<Vec>, SmartPetscObj<VecScatter>> uYScatter;
185  std::tuple<SmartPetscObj<Vec>, SmartPetscObj<VecScatter>> uZScatter;
186 };
187 
188 }; // namespace PlasticOps
#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
void temp(int x, int y=10)
Definition: simple.cpp:4
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
double v
phase velocity of light in medium (cm/ns)
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
MatrixBoundedArray< double, 9 > MatrixDouble3by3
Definition: Types.hpp:116
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1965
auto hardening(long double tau)
const double r
rate factor
double h
convective heat coefficient
constexpr double t
plate stiffness
Definition: plate.cpp:76
std::array< bool, MBMAXTYPE > doEntities
If true operator is executed for entity.
Data on single entity (This is passed as argument to DataOperator::doWork)
SmartPetscObj< DM > dM
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)
std::tuple< SmartPetscObj< Vec >, SmartPetscObj< VecScatter > > uXScatter
std::tuple< SmartPetscObj< Vec >, SmartPetscObj< VecScatter > > uZScatter
MoFEMErrorCode preProcess()
MoFEMErrorCode postProcess()
boost::shared_ptr< PostProcEle > postProcFe
std::tuple< SmartPetscObj< Vec >, SmartPetscObj< VecScatter > > uYScatter
MoFEMErrorCode operator()()
boost::shared_ptr< DomainEle > reactionFe
boost::shared_ptr< CommonData > commonDataPtr
Definition: PlasticOps.hpp:289
std::vector< EntityHandle > & mapGaussPts
Definition: PlasticOps.hpp:288
OpPostProcPlastic(const std::string field_name, moab::Interface &post_proc_mesh, std::vector< EntityHandle > &map_gauss_pts, boost::shared_ptr< CommonData > common_data_ptr)
moab::Interface & postProcMesh
Definition: PlasticOps.hpp:287
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
[Postprocessing]