v0.14.0
HookeInternalStressElement.hpp
Go to the documentation of this file.
1 /**
2  * \file HookeInternalStressElement.hpp
3  * \example HookeInternalStressElement.hpp
4  *
5  * \brief Operators and data structures for calculating internal stress
6  *
7  */
8 
9 /*
10  * This file is part of MoFEM.
11  * MoFEM is free software: you can redistribute it and/or modify it under
12  * the terms of the GNU Lesser General Public License as published by the
13  * Free Software Foundation, either version 3 of the License, or (at your
14  * option) any later version.
15  *
16  * MoFEM is distributed in the hope that it will be useful, but WITHOUT
17  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
19  * License for more details.
20  *
21  * You should have received a copy of the GNU Lesser General Public
22  * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>. */
23 
24 #ifndef __HOOKE_INTERNAL_STRESS_ELEMENT_HPP__
25 #define __HOOKE_INTERNAL_STRESS_ELEMENT_HPP__
26 
27 #include <BasicFiniteElements.hpp>
28 
30 
31  struct DataAtIntegrationPts : public HookeElement::DataAtIntegrationPts {
32 
33  boost::shared_ptr<MatrixDouble> internalStressMat;
34  boost::shared_ptr<MatrixDouble> actualStressMat;
35  boost::shared_ptr<MatrixDouble> deviatoricStressMat;
36  boost::shared_ptr<MatrixDouble> hydrostaticStressMat;
37 
38  boost::shared_ptr<MatrixDouble> spatPosMat;
39  boost::shared_ptr<MatrixDouble> meshNodePosMat;
40 
42 
43  internalStressMat = boost::shared_ptr<MatrixDouble>(new MatrixDouble());
44  actualStressMat = boost::shared_ptr<MatrixDouble>(new MatrixDouble());
45  deviatoricStressMat = boost::shared_ptr<MatrixDouble>(new MatrixDouble());
47  boost::shared_ptr<MatrixDouble>(new MatrixDouble());
48 
49  spatPosMat = boost::shared_ptr<MatrixDouble>(new MatrixDouble());
50  meshNodePosMat = boost::shared_ptr<MatrixDouble>(new MatrixDouble());
51  }
52  };
53 
55 
56  OpGetInternalStress(const std::string row_field,
57  const std::string col_field,
58  boost::shared_ptr<DataAtIntegrationPts> data_at_pts,
59  moab::Interface &input_mesh, char *stress_tag_name)
60  : VolUserDataOperator(row_field, col_field, OPROW, true),
61  dataAtPts(data_at_pts), inputMesh(input_mesh),
62  stressTagName(stress_tag_name) {
63  std::fill(&doEntities[MBEDGE], &doEntities[MBMAXTYPE], false);
64  };
65 
66  MoFEMErrorCode doWork(int row_side, EntityType row_type, EntData &row_data);
67 
68  protected:
69  boost::shared_ptr<DataAtIntegrationPts> dataAtPts;
72  };
73 
74  struct OpInternalStrain_dx : public OpAssemble {
75 
76  OpInternalStrain_dx(const std::string row_field,
77  boost::shared_ptr<DataAtIntegrationPts> data_at_pts)
78  : OpAssemble(row_field, row_field, data_at_pts, OPROW, true),
79  dataAtPts(data_at_pts){};
80 
81  MoFEMErrorCode iNtegrate(EntData &row_data);
82 
83  protected:
84  boost::shared_ptr<DataAtIntegrationPts> dataAtPts;
85  };
86 
87  template <int S = 0>
89 
90  typedef boost::function<
91 
93 
95 
96  )
97 
98  >
100 
102  const std::string row_field, const std::string col_field,
103  boost::shared_ptr<DataAtIntegrationPts> data_at_pts,
104  StrainFunction strain_fun)
105  : VolUserDataOperator(row_field, col_field, OPROW, true),
106  dataAtPts(data_at_pts), strainFun(strain_fun) {
107  std::fill(&doEntities[MBEDGE], &doEntities[MBMAXTYPE], false);
108  }
109 
110  MoFEMErrorCode doWork(int row_side, EntityType row_type, EntData &row_data);
111 
112  protected:
114  boost::shared_ptr<DataAtIntegrationPts> dataAtPts;
115  };
116 
118  boost::shared_ptr<DataAtIntegrationPts> dataAtPts;
119  map<int, BlockData>
120  &blockSetsPtr; // FIXME: (works only with the first block)
122  bool isALE;
124  double scaleFactor;
125  bool saveMean;
126 
127  OpSaveStress(const string row_field, const string col_field,
128  boost::shared_ptr<DataAtIntegrationPts> data_at_pts,
129  map<int, BlockData> &block_sets_ptr,
130  moab::Interface &output_mesh, double scale_factor,
131  bool save_mean = false, bool is_ale = false,
132  bool is_field_disp = true)
133  : VolUserDataOperator(row_field, col_field, OPROW, true),
134  dataAtPts(data_at_pts), blockSetsPtr(block_sets_ptr),
135  outputMesh(output_mesh), scaleFactor(scale_factor), isALE(is_ale),
136  isFieldDisp(is_field_disp), saveMean(save_mean){};
137 
138  MoFEMErrorCode doWork(int row_side, EntityType row_type, EntData &row_data);
139  };
140 
141  template <class ELEMENT>
143  boost::shared_ptr<DataAtIntegrationPts> dataAtPts;
144  map<int, BlockData>
145  &blockSetsPtr; // FIXME: (works only with the first block)
147  std::vector<EntityHandle> &mapGaussPts;
148  bool isALE;
150 
151  OpPostProcHookeElement(const string row_field,
152  boost::shared_ptr<DataAtIntegrationPts> data_at_pts,
153  map<int, BlockData> &block_sets_ptr,
154  moab::Interface &post_proc_mesh,
155  std::vector<EntityHandle> &map_gauss_pts,
156  bool is_ale = false, bool is_field_disp = true)
157  : ELEMENT::UserDataOperator(row_field, UserDataOperator::OPROW),
158  dataAtPts(data_at_pts), blockSetsPtr(block_sets_ptr),
159  postProcMesh(post_proc_mesh), mapGaussPts(map_gauss_pts),
160  isALE(is_ale), isFieldDisp(is_field_disp) {}
161 
162  MoFEMErrorCode doWork(int side, EntityType type,
164  };
165 };
166 
168  int row_side, EntityType row_type, EntData &row_data) {
170 
171  const int nb_integration_pts = row_data.getN().size1();
172 
173  const EntityHandle ent = getFEEntityHandle();
174 
175  const int val_num = 9 * nb_integration_pts;
176  std::vector<double> def_vals(val_num, 0.0);
177  Tag th_internal_stress;
178  CHKERR inputMesh.tag_get_handle(
179  stressTagName, val_num, MB_TYPE_DOUBLE, th_internal_stress,
180  MB_TAG_CREAT | MB_TAG_SPARSE, &*def_vals.begin());
181 
182  dataAtPts->internalStressMat->resize(9, nb_integration_pts, false);
183 
184  CHKERR inputMesh.tag_get_data(
185  th_internal_stress, &ent, 1,
186  &*(dataAtPts->internalStressMat->data().begin()));
187 
189 }
190 
198 
199  auto get_tensor1 = [](VectorDouble &v, const int r) {
201  &v(r + 0), &v(r + 1), &v(r + 2));
202  };
203 
204  const int nb_integration_pts = getGaussPts().size2();
205 
206  auto t_internal_stress =
207  getFTensor2FromMat<3, 3>(*dataAtPts->internalStressMat);
208 
209  // get element volume
210  double vol = getVolume();
211  auto t_w = getFTensor0IntegrationWeight();
212 
213  nF.resize(nbRows, false);
214  nF.clear();
215 
216  // get derivatives of base functions on rows
217  auto t_row_diff_base = row_data.getFTensor1DiffN<3>();
218  const int row_nb_base_fun = row_data.getN().size2();
219 
220  for (size_t gg = 0; gg != nb_integration_pts; ++gg) {
221 
222  // calculate scalar weight times element volume
223  double a = t_w * vol;
224  auto t_nf = get_tensor1(nF, 0);
225 
226  int rr = 0;
227  for (; rr != nbRows / 3; ++rr) {
228  t_nf(i) += a * t_row_diff_base(j) * t_internal_stress(i, j);
229  ++t_row_diff_base;
230  ++t_nf;
231  }
232 
233  for (; rr != row_nb_base_fun; ++rr)
234  ++t_row_diff_base;
235 
236  ++t_w;
237  ++t_internal_stress;
238  }
239 
241 }
242 
243 template <int S>
246  int row_side, EntityType row_type, EntData &row_data) {
252 
253  auto tensor_to_tensor = [](const auto &t1, auto &t2) {
254  t2(0, 0) = t1(0, 0);
255  t2(1, 1) = t1(1, 1);
256  t2(2, 2) = t1(2, 2);
257  t2(0, 1) = t2(1, 0) = t1(1, 0);
258  t2(0, 2) = t2(2, 0) = t1(2, 0);
259  t2(1, 2) = t2(2, 1) = t1(2, 1);
260  };
261 
262  const int nb_integration_pts = getGaussPts().size2();
263  auto t_coords = getFTensor1CoordsAtGaussPts();
264 
265  // elastic stiffness tensor (4th rank tensor with minor and major
266  // symmetry)
268  MAT_TO_DDG(dataAtPts->stiffnessMat));
269 
270  dataAtPts->internalStressMat->resize(9, nb_integration_pts, false);
271  auto t_internal_stress =
272  getFTensor2FromMat<3, 3>(*dataAtPts->internalStressMat);
273 
274  FTensor::Tensor2_symmetric<double, 3> t_internal_stress_symm;
275 
276  for (size_t gg = 0; gg != nb_integration_pts; ++gg) {
277 
278  auto t_fun_strain = strainFun(t_coords);
279  t_internal_stress_symm(i, j) = -t_D(i, j, k, l) * t_fun_strain(k, l);
280 
281  tensor_to_tensor(t_internal_stress_symm, t_internal_stress);
282 
283  ++t_coords;
284  ++t_D;
285  ++t_internal_stress;
286  }
287 
289 }
290 
292  int row_side, EntityType row_type, EntData &row_data) {
294 
295  if (row_type != MBVERTEX) {
297  }
298 
299  const EntityHandle ent = getFEEntityHandle();
300 
301  auto tensor_to_tensor = [](const auto &t1, auto &t2) {
302  t2(0, 0) = t1(0, 0);
303  t2(1, 1) = t1(1, 1);
304  t2(2, 2) = t1(2, 2);
305  t2(0, 1) = t2(1, 0) = t1(1, 0);
306  t2(0, 2) = t2(2, 0) = t1(2, 0);
307  t2(1, 2) = t2(2, 1) = t1(2, 1);
308  };
309 
310  auto tensor_to_vector = [](const auto &t, auto &v) {
311  v(0) = t(0, 0);
312  v(1) = t(1, 1);
313  v(2) = t(2, 2);
314  v(3) = t(0, 1);
315  v(4) = t(0, 2);
316  v(5) = t(1, 2);
317  };
318 
319  auto get_tag_handle = [&](auto name, auto size) {
320  Tag th;
321  std::vector<double> def_vals(size, 0.0);
322  CHKERR outputMesh.tag_get_handle(name, size, MB_TYPE_DOUBLE, th,
323  MB_TAG_CREAT | MB_TAG_SPARSE,
324  def_vals.data());
325  return th;
326  };
327 
328  const int nb_integration_pts = row_data.getN().size1();
329 
330  auto th_internal_stress =
331  get_tag_handle("INTERNAL_STRESS", 9 * nb_integration_pts);
332  auto th_actual_stress =
333  get_tag_handle("ACTUAL_STRESS", 9 * nb_integration_pts);
334  auto th_deviatoric_stress =
335  get_tag_handle("DEVIATORIC_STRESS", 9 * nb_integration_pts);
336  auto th_hydrostatic_stress =
337  get_tag_handle("HYDROSTATIC_STRESS", 9 * nb_integration_pts);
338 
339  auto th_internal_stress_mean = get_tag_handle("MED_INTERNAL_STRESS", 9);
340  auto th_actual_stress_mean = get_tag_handle("MED_ACTUAL_STRESS", 9);
341 
346 
347  auto t_h = getFTensor2FromMat<3, 3>(*dataAtPts->hMat);
348  auto t_H = getFTensor2FromMat<3, 3>(*dataAtPts->HMat);
349 
350  dataAtPts->stiffnessMat->resize(36, 1, false);
352  MAT_TO_DDG(dataAtPts->stiffnessMat));
353  for (auto &m : (blockSetsPtr)) {
354  const double young = m.second.E;
355  const double poisson = m.second.PoissonRatio;
356 
357  const double coefficient = young / ((1 + poisson) * (1 - 2 * poisson));
358 
359  t_D(i, j, k, l) = 0.;
360  t_D(0, 0, 0, 0) = t_D(1, 1, 1, 1) = t_D(2, 2, 2, 2) = 1 - poisson;
361  t_D(0, 1, 0, 1) = t_D(0, 2, 0, 2) = t_D(1, 2, 1, 2) =
362  0.5 * (1 - 2 * poisson);
363  t_D(0, 0, 1, 1) = t_D(1, 1, 0, 0) = t_D(0, 0, 2, 2) = t_D(2, 2, 0, 0) =
364  t_D(1, 1, 2, 2) = t_D(2, 2, 1, 1) = poisson;
365  t_D(i, j, k, l) *= coefficient;
366 
367  break; // FIXME: calculates only first block
368  }
369 
370  double detH = 0.;
373 
375  FTensor::Tensor2_symmetric<double, 3> t_small_strain_symm;
376 
378 
379  auto t_internal_stress =
380  getFTensor2FromMat<3, 3>(*dataAtPts->internalStressMat);
381 
382  dataAtPts->actualStressMat->resize(9, nb_integration_pts, false);
383  auto t_actual_stress = getFTensor2FromMat<3, 3>(*dataAtPts->actualStressMat);
384 
385  dataAtPts->deviatoricStressMat->resize(9, nb_integration_pts, false);
386  auto t_deviatoric_stress =
387  getFTensor2FromMat<3, 3>(*dataAtPts->deviatoricStressMat);
388  dataAtPts->hydrostaticStressMat->resize(9, nb_integration_pts, false);
389  auto t_hydrostatic_stress =
390  getFTensor2FromMat<3, 3>(*dataAtPts->hydrostaticStressMat);
391 
392  FTensor::Tensor2<double, 3, 3> t_internal_stress_mean;
393  FTensor::Tensor2<double, 3, 3> t_actual_stress_mean;
394 
395  t_internal_stress_mean(i, j) = 0.;
396  t_actual_stress_mean(i, j) = 0.;
397 
398  auto t_w = getFTensor0IntegrationWeight();
399 
400  for (int gg = 0; gg != nb_integration_pts; ++gg) {
401 
402  if (!isALE) {
403  t_small_strain_symm(i, j) = (t_h(i, j) || t_h(j, i)) / 2.;
404  } else {
405  CHKERR determinantTensor3by3(t_H, detH);
406  CHKERR invertTensor3by3(t_H, detH, t_invH);
407  t_F(i, j) = t_h(i, k) * t_invH(k, j);
408  t_small_strain_symm(i, j) = (t_F(i, j) || t_F(j, i)) / 2.;
409  ++t_H;
410  }
411 
412  if (isALE || !isFieldDisp) {
413  for (auto ii : {0, 1, 2}) {
414  t_small_strain_symm(ii, ii) -= 1.;
415  }
416  }
417 
418  // symmetric tensors need improvement
419  t_stress_symm(i, j) = t_D(i, j, k, l) * t_small_strain_symm(k, l);
420  tensor_to_tensor(t_stress_symm, t_stress);
421  t_actual_stress(i, j) = t_stress(i, j) + t_internal_stress(i, j);
422 
423  t_actual_stress(i, j) *= scaleFactor;
424  t_internal_stress(i, j) *= scaleFactor;
425 
426  double hydrostatic_pressure =
427  (t_actual_stress(0, 0) + t_actual_stress(1, 1) +
428  t_actual_stress(2, 2)) /
429  3.;
430 
431  t_hydrostatic_stress(i, j) = 0.;
432  for (auto ii : {0, 1, 2}) {
433  t_hydrostatic_stress(ii, ii) += hydrostatic_pressure;
434  }
435 
436  t_deviatoric_stress(i, j) = t_actual_stress(i, j);
437  for (auto ii : {0, 1, 2}) {
438  t_deviatoric_stress(ii, ii) -= hydrostatic_pressure;
439  }
440 
441  t_actual_stress_mean(i, j) += t_w * t_actual_stress(i, j);
442  t_internal_stress_mean(i, j) += t_w * t_internal_stress(i, j);
443 
444  ++t_w;
445  ++t_h;
446  ++t_actual_stress;
447  ++t_internal_stress;
448  }
449 
450  if (saveMean) {
451  VectorDouble vec_actual_stress_mean;
452  vec_actual_stress_mean.resize(9, false);
453  vec_actual_stress_mean.clear();
454 
455  VectorDouble vec_internal_stress_mean;
456  vec_internal_stress_mean.resize(9, false);
457  vec_internal_stress_mean.clear();
458 
459  tensor_to_vector(t_actual_stress_mean, vec_actual_stress_mean);
460  tensor_to_vector(t_internal_stress_mean, vec_internal_stress_mean);
461 
462  CHKERR outputMesh.tag_set_data(th_actual_stress_mean, &ent, 1,
463  &*(vec_actual_stress_mean.data().begin()));
464  CHKERR outputMesh.tag_set_data(th_internal_stress_mean, &ent, 1,
465  &*(vec_internal_stress_mean.data().begin()));
466  } else {
467  CHKERR outputMesh.tag_set_data(
468  th_internal_stress, &ent, 1,
469  &*(dataAtPts->internalStressMat->data().begin()));
470  CHKERR outputMesh.tag_set_data(
471  th_actual_stress, &ent, 1,
472  &*(dataAtPts->actualStressMat->data().begin()));
473 
474  CHKERR outputMesh.tag_set_data(
475  th_hydrostatic_stress, &ent, 1,
476  &*(dataAtPts->hydrostaticStressMat->data().begin()));
477  CHKERR outputMesh.tag_set_data(
478  th_deviatoric_stress, &ent, 1,
479  &*(dataAtPts->deviatoricStressMat->data().begin()));
480  }
481 
483 }
484 
485 template <class ELEMENT>
488  int side, EntityType type, EntitiesFieldData::EntData &data) {
490 
491  if (type != MBVERTEX) {
493  }
494 
495  auto tensor_to_tensor = [](const auto &t1, auto &t2) {
496  t2(0, 0) = t1(0, 0);
497  t2(1, 1) = t1(1, 1);
498  t2(2, 2) = t1(2, 2);
499  t2(0, 1) = t2(1, 0) = t1(1, 0);
500  t2(0, 2) = t2(2, 0) = t1(2, 0);
501  t2(1, 2) = t2(2, 1) = t1(2, 1);
502  };
503 
504  auto get_tag_handle = [&](auto name, auto size) {
505  Tag th;
506  std::vector<double> def_vals(size, 0.0);
507  CHKERR postProcMesh.tag_get_handle(name, size, MB_TYPE_DOUBLE, th,
508  MB_TAG_CREAT | MB_TAG_SPARSE,
509  def_vals.data());
510  return th;
511  };
512 
513  auto th_stress = get_tag_handle("STRESS", 9);
514  auto th_strain = get_tag_handle("STRAIN", 9);
515  auto th_psi = get_tag_handle("ENERGY", 1);
516  auto th_disp = get_tag_handle("DISPLACEMENT", 3);
517 
518  const int nb_integration_pts = mapGaussPts.size();
519 
524 
525  auto t_h = getFTensor2FromMat<3, 3>(*dataAtPts->hMat);
526  auto t_H = getFTensor2FromMat<3, 3>(*dataAtPts->HMat);
527 
528  dataAtPts->stiffnessMat->resize(36, 1, false);
530  MAT_TO_DDG(dataAtPts->stiffnessMat));
531  for (auto &m : (blockSetsPtr)) {
532  const double young = m.second.E;
533  const double poisson = m.second.PoissonRatio;
534 
535  const double coefficient = young / ((1 + poisson) * (1 - 2 * poisson));
536 
537  t_D(i, j, k, l) = 0.;
538  t_D(0, 0, 0, 0) = t_D(1, 1, 1, 1) = t_D(2, 2, 2, 2) = 1 - poisson;
539  t_D(0, 1, 0, 1) = t_D(0, 2, 0, 2) = t_D(1, 2, 1, 2) =
540  0.5 * (1 - 2 * poisson);
541  t_D(0, 0, 1, 1) = t_D(1, 1, 0, 0) = t_D(0, 0, 2, 2) = t_D(2, 2, 0, 0) =
542  t_D(1, 1, 2, 2) = t_D(2, 2, 1, 1) = poisson;
543  t_D(i, j, k, l) *= coefficient;
544 
545  break; // FIXME: calculates only first block
546  }
547 
548  double detH = 0.;
552  FTensor::Tensor2<double, 3, 3> t_small_strain;
553 
555  FTensor::Tensor2_symmetric<double, 3> t_small_strain_symm;
556 
557  auto t_spat_pos = getFTensor1FromMat<3>(*dataAtPts->spatPosMat);
558  auto t_mesh_node_pos = getFTensor1FromMat<3>(*dataAtPts->meshNodePosMat);
560 
561  for (int gg = 0; gg != nb_integration_pts; ++gg) {
562 
563  if (isFieldDisp) {
564  t_h(0, 0) += 1;
565  t_h(1, 1) += 1;
566  t_h(2, 2) += 1;
567  }
568 
569  if (!isALE) {
570  t_small_strain_symm(i, j) = (t_h(i, j) || t_h(j, i)) / 2.;
571  } else {
572  CHKERR determinantTensor3by3(t_H, detH);
573  CHKERR invertTensor3by3(t_H, detH, t_invH);
574  t_F(i, j) = t_h(i, k) * t_invH(k, j);
575  t_small_strain_symm(i, j) = (t_F(i, j) || t_F(j, i)) / 2.;
576  ++t_H;
577  }
578 
579  t_small_strain_symm(0, 0) -= 1;
580  t_small_strain_symm(1, 1) -= 1;
581  t_small_strain_symm(2, 2) -= 1;
582 
583  // symmetric tensors need improvement
584  t_stress_symm(i, j) = t_D(i, j, k, l) * t_small_strain_symm(k, l);
585  tensor_to_tensor(t_stress_symm, t_stress);
586  tensor_to_tensor(t_small_strain_symm, t_small_strain);
587 
588  t_disp(i) = t_spat_pos(i) - t_mesh_node_pos(i);
589 
590  const double psi = 0.5 * t_stress_symm(i, j) * t_small_strain_symm(i, j);
591 
592  CHKERR postProcMesh.tag_set_data(th_psi, &mapGaussPts[gg], 1, &psi);
593  CHKERR postProcMesh.tag_set_data(th_stress, &mapGaussPts[gg], 1,
594  &t_stress(0, 0));
595  CHKERR postProcMesh.tag_set_data(th_strain, &mapGaussPts[gg], 1,
596  &t_small_strain(0, 0));
597  CHKERR postProcMesh.tag_set_data(th_disp, &mapGaussPts[gg], 1, &t_disp(0));
598 
599  ++t_h;
600  ++t_spat_pos;
601  ++t_mesh_node_pos;
602  }
603 
605 }
606 
607 #endif
HookeInternalStressElement::OpGetInternalStress
Definition: HookeInternalStressElement.hpp:54
MoFEMFunctionReturnHot
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:460
HookeInternalStressElement::OpGetInternalStress::stressTagName
char * stressTagName
Definition: HookeInternalStressElement.hpp:71
MoFEM::EntitiesFieldData::EntData
Data on single entity (This is passed as argument to DataOperator::doWork)
Definition: EntitiesFieldData.hpp:128
HookeInternalStressElement::OpSaveStress::saveMean
bool saveMean
Definition: HookeInternalStressElement.hpp:125
FTensor::Tensor1
Definition: Tensor1_value.hpp:8
EntityHandle
HookeInternalStressElement::OpSaveStress
Definition: HookeInternalStressElement.hpp:117
HookeInternalStressElement::OpGetAnalyticalInternalStress::StrainFunction
boost::function< FTensor::Tensor2_symmetric< double, 3 > FTensor::Tensor1< FTensor::PackPtr< double *, 3 >, 3 > &t_coords) > StrainFunction
Definition: HookeInternalStressElement.hpp:99
ELEMENT
HookeInternalStressElement::OpSaveStress::isALE
bool isALE
Definition: HookeInternalStressElement.hpp:122
MoFEM::Exceptions::MoFEMErrorCode
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
MoFEM::Types::MatrixDouble
UBlasMatrix< double > MatrixDouble
Definition: Types.hpp:77
HookeInternalStressElement::DataAtIntegrationPts::internalStressMat
boost::shared_ptr< MatrixDouble > internalStressMat
Definition: HookeInternalStressElement.hpp:33
HookeInternalStressElement::OpSaveStress::doWork
MoFEMErrorCode doWork(int row_side, EntityType row_type, EntData &row_data)
Operator for linear form, usually to calculate values on right hand side.
Definition: HookeInternalStressElement.hpp:291
MoFEM::th
Tag th
Definition: Projection10NodeCoordsOnField.cpp:122
HookeInternalStressElement::OpSaveStress::dataAtPts
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
Definition: HookeInternalStressElement.hpp:118
HookeInternalStressElement::OpSaveStress::scaleFactor
double scaleFactor
Definition: HookeInternalStressElement.hpp:124
HookeInternalStressElement::OpGetAnalyticalInternalStress::dataAtPts
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
Definition: HookeInternalStressElement.hpp:114
BasicFiniteElements.hpp
FTensor::Tensor2_symmetric
Definition: Tensor2_symmetric_value.hpp:13
HookeInternalStressElement::OpSaveStress::OpSaveStress
OpSaveStress(const string row_field, const string col_field, boost::shared_ptr< DataAtIntegrationPts > data_at_pts, map< int, BlockData > &block_sets_ptr, moab::Interface &output_mesh, double scale_factor, bool save_mean=false, bool is_ale=false, bool is_field_disp=true)
Definition: HookeInternalStressElement.hpp:127
HookeInternalStressElement::OpGetInternalStress::OpGetInternalStress
OpGetInternalStress(const std::string row_field, const std::string col_field, boost::shared_ptr< DataAtIntegrationPts > data_at_pts, moab::Interface &input_mesh, char *stress_tag_name)
Definition: HookeInternalStressElement.hpp:56
sdf.r
int r
Definition: sdf.py:8
FTensor::Tensor2< double, 3, 3 >
HookeInternalStressElement::OpPostProcHookeElement::blockSetsPtr
map< int, BlockData > & blockSetsPtr
Definition: HookeInternalStressElement.hpp:145
MoFEM::invertTensor3by3
MoFEMErrorCode invertTensor3by3(ublas::matrix< T, L, A > &jac_data, ublas::vector< T, A > &det_data, ublas::matrix< T, L, A > &inv_jac_data)
Calculate inverse of tensor rank 2 at integration points.
Definition: Templates.hpp:1588
MoFEM::Interface
DeprecatedCoreInterface Interface
Definition: Interface.hpp:2010
HookeInternalStressElement::OpPostProcHookeElement::OpPostProcHookeElement
OpPostProcHookeElement(const string row_field, boost::shared_ptr< DataAtIntegrationPts > data_at_pts, map< int, BlockData > &block_sets_ptr, moab::Interface &post_proc_mesh, std::vector< EntityHandle > &map_gauss_pts, bool is_ale=false, bool is_field_disp=true)
Definition: HookeInternalStressElement.hpp:151
CHKERR
#define CHKERR
Inline error check.
Definition: definitions.h:548
a
constexpr double a
Definition: approx_sphere.cpp:30
MoFEM::ForcesAndSourcesCore::UserDataOperator
Definition: ForcesAndSourcesCore.hpp:549
HookeInternalStressElement::DataAtIntegrationPts::deviatoricStressMat
boost::shared_ptr< MatrixDouble > deviatoricStressMat
Definition: HookeInternalStressElement.hpp:35
convert.type
type
Definition: convert.py:64
HookeInternalStressElement::OpGetInternalStress::dataAtPts
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
Definition: HookeInternalStressElement.hpp:69
HookeInternalStressElement::OpSaveStress::isFieldDisp
bool isFieldDisp
Definition: HookeInternalStressElement.hpp:123
HookeElement
HookeInternalStressElement::OpInternalStrain_dx::iNtegrate
MoFEMErrorCode iNtegrate(EntData &row_data)
Definition: HookeInternalStressElement.hpp:192
HookeInternalStressElement::OpInternalStrain_dx::OpInternalStrain_dx
OpInternalStrain_dx(const std::string row_field, boost::shared_ptr< DataAtIntegrationPts > data_at_pts)
Definition: HookeInternalStressElement.hpp:76
HookeInternalStressElement::OpGetInternalStress::inputMesh
moab::Interface & inputMesh
Definition: HookeInternalStressElement.hpp:70
FTensor::PackPtr
Definition: FTensor.hpp:54
HookeInternalStressElement::OpSaveStress::blockSetsPtr
map< int, BlockData > & blockSetsPtr
Definition: HookeInternalStressElement.hpp:120
HookeInternalStressElement::OpPostProcHookeElement::mapGaussPts
std::vector< EntityHandle > & mapGaussPts
Definition: HookeInternalStressElement.hpp:147
t
constexpr double t
plate stiffness
Definition: plate.cpp:58
i
FTensor::Index< 'i', SPACE_DIM > i
Definition: hcurl_divergence_operator_2d.cpp:27
HookeInternalStressElement::OpInternalStrain_dx::dataAtPts
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
Definition: HookeInternalStressElement.hpp:84
HookeInternalStressElement::OpPostProcHookeElement
Definition: HookeInternalStressElement.hpp:142
HookeInternalStressElement::OpPostProcHookeElement::isFieldDisp
bool isFieldDisp
Definition: HookeInternalStressElement.hpp:149
MoFEM::VolumeElementForcesAndSourcesCore::UserDataOperator
Definition: VolumeElementForcesAndSourcesCore.hpp:108
HookeInternalStressElement::OpPostProcHookeElement::postProcMesh
moab::Interface & postProcMesh
Definition: HookeInternalStressElement.hpp:146
FTensor::Index< 'i', 3 >
HookeInternalStressElement::OpGetAnalyticalInternalStress::OpGetAnalyticalInternalStress
OpGetAnalyticalInternalStress(const std::string row_field, const std::string col_field, boost::shared_ptr< DataAtIntegrationPts > data_at_pts, StrainFunction strain_fun)
Definition: HookeInternalStressElement.hpp:101
MoFEM::determinantTensor3by3
static auto determinantTensor3by3(T &t)
Calculate the determinant of a 3x3 matrix or a tensor of rank 2.
Definition: Templates.hpp:1540
v
const double v
phase velocity of light in medium (cm/ns)
Definition: initial_diffusion.cpp:40
HookeInternalStressElement::OpSaveStress::outputMesh
moab::Interface & outputMesh
Definition: HookeInternalStressElement.hpp:121
MoFEM::DataOperator::doEntities
std::array< bool, MBMAXTYPE > doEntities
If true operator is executed for entity.
Definition: DataOperators.hpp:85
HookeInternalStressElement::DataAtIntegrationPts::spatPosMat
boost::shared_ptr< MatrixDouble > spatPosMat
Definition: HookeInternalStressElement.hpp:38
MoFEM::EntitiesFieldData::EntData::getFTensor1DiffN
FTensor::Tensor1< FTensor::PackPtr< double *, Tensor_Dim >, Tensor_Dim > getFTensor1DiffN(const FieldApproximationBase base)
Get derivatives of base functions.
Definition: EntitiesFieldData.cpp:526
MoFEM::ForcesAndSourcesCore::UserDataOperator::getFEEntityHandle
EntityHandle getFEEntityHandle() const
Return finite element entity handle.
Definition: ForcesAndSourcesCore.hpp:1004
HookeInternalStressElement::DataAtIntegrationPts
Definition: HookeInternalStressElement.hpp:31
MoFEM::EntitiesFieldData::EntData::getN
MatrixDouble & getN(const FieldApproximationBase base)
get base functions this return matrix (nb. of rows is equal to nb. of Gauss pts, nb....
Definition: EntitiesFieldData.hpp:1318
HookeInternalStressElement::OpGetAnalyticalInternalStress
Definition: HookeInternalStressElement.hpp:88
HookeInternalStressElement::OpPostProcHookeElement::isALE
bool isALE
Definition: HookeInternalStressElement.hpp:148
MAT_TO_DDG
#define MAT_TO_DDG(SM)
Definition: HookeElement.hpp:142
UserDataOperator
ForcesAndSourcesCore::UserDataOperator UserDataOperator
Definition: HookeElement.hpp:75
j
FTensor::Index< 'j', 3 > j
Definition: matrix_function.cpp:19
FTensor::Ddg
Definition: Ddg_value.hpp:7
HookeInternalStressElement::DataAtIntegrationPts::hydrostaticStressMat
boost::shared_ptr< MatrixDouble > hydrostaticStressMat
Definition: HookeInternalStressElement.hpp:36
HookeInternalStressElement::OpPostProcHookeElement::doWork
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Definition: HookeInternalStressElement.hpp:487
MoFEM::Types::VectorDouble
UBlasVector< double > VectorDouble
Definition: Types.hpp:68
HookeInternalStressElement::OpInternalStrain_dx
Definition: HookeInternalStressElement.hpp:74
HookeInternalStressElement::OpGetInternalStress::doWork
MoFEMErrorCode doWork(int row_side, EntityType row_type, EntData &row_data)
Operator for linear form, usually to calculate values on right hand side.
Definition: HookeInternalStressElement.hpp:167
m
FTensor::Index< 'm', 3 > m
Definition: shallow_wave.cpp:80
OpAssemble
Definition: HookeElement.hpp:340
HookeInternalStressElement::DataAtIntegrationPts::actualStressMat
boost::shared_ptr< MatrixDouble > actualStressMat
Definition: HookeInternalStressElement.hpp:34
k
FTensor::Index< 'k', 3 > k
Definition: matrix_function.cpp:20
HookeInternalStressElement
Definition: HookeInternalStressElement.hpp:29
MoFEMFunctionReturn
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
HookeInternalStressElement::DataAtIntegrationPts::DataAtIntegrationPts
DataAtIntegrationPts()
Definition: HookeInternalStressElement.hpp:41
HookeInternalStressElement::DataAtIntegrationPts::meshNodePosMat
boost::shared_ptr< MatrixDouble > meshNodePosMat
Definition: HookeInternalStressElement.hpp:39
HookeInternalStressElement::OpGetAnalyticalInternalStress::strainFun
StrainFunction strainFun
Definition: HookeInternalStressElement.hpp:113
HookeInternalStressElement::OpGetAnalyticalInternalStress::doWork
MoFEMErrorCode doWork(int row_side, EntityType row_type, EntData &row_data)
Operator for linear form, usually to calculate values on right hand side.
Definition: HookeInternalStressElement.hpp:245
MoFEMFunctionBegin
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
l
FTensor::Index< 'l', 3 > l
Definition: matrix_function.cpp:21
MoFEM::ForcesAndSourcesCore::UserDataOperator::OPROW
@ OPROW
operator doWork function is executed on FE rows
Definition: ForcesAndSourcesCore.hpp:567
HookeInternalStressElement::OpPostProcHookeElement::dataAtPts
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
Definition: HookeInternalStressElement.hpp:143