v0.14.0
Classes | Namespaces | Macros | Functions
EshelbianPlasticity.cpp File Reference

Implementation of automatic differentiation. More...

#include <MoFEM.hpp>
#include <CGGTonsorialBubbleBase.hpp>
#include <EshelbianPlasticity.hpp>
#include <boost/math/constants/constants.hpp>
#include <cholesky.hpp>
#include <EshelbianAux.hpp>
#include <EshelbianContact.hpp>
#include <phg-quadrule/quad.h>
#include "impl/EshelbianMonitor.cpp"
#include "impl/SetUpSchurImpl.cpp"
#include "impl/TSElasticPostStep.cpp"
#include <impl/EshelbianContact.cpp>

Go to the source code of this file.

Classes

struct  EshelbianPlasticity::VolUserDataOperatorStabAssembly
 
struct  EshelbianPlasticity::FaceUserDataOperatorStabAssembly
 
struct  EshelbianPlasticity::DynamicRelaxationTimeScale
 
struct  EshelbianPlasticity::SetIntegrationAtFrontVolume
 
struct  EshelbianPlasticity::SetIntegrationAtFrontVolume::Fe
 
struct  EshelbianPlasticity::SetIntegrationAtFrontFace
 
struct  EshelbianPlasticity::SetIntegrationAtFrontFace::Fe
 
struct  EshelbianPlasticity::VolRule
 Set integration rule on element. More...
 
struct  EshelbianPlasticity::FaceRule
 
struct  EshelbianPlasticity::solve_elastic_setup
 

Namespaces

 EshelbianPlasticity
 

Macros

#define SINGULARITY
 

Functions

static auto send_type (MoFEM::Interface &m_field, Range r, const EntityType type)
 
static auto get_range_from_block (MoFEM::Interface &m_field, const std::string block_name, int dim)
 
static auto get_block_meshset (MoFEM::Interface &m_field, const int ms_id, const unsigned int cubit_bc_type)
 
static auto save_range (moab::Interface &moab, const std::string name, const Range r, std::vector< Tag > tags={})
 
static auto filter_true_skin (MoFEM::Interface &m_field, Range &&skin)
 
static auto filter_owners (MoFEM::Interface &m_field, Range skin)
 
static auto get_skin (MoFEM::Interface &m_field, Range body_ents)
 
static auto get_crack_front_edges (MoFEM::Interface &m_field, Range crack_faces)
 
static auto get_two_sides_of_crack_surface (MoFEM::Interface &m_field, Range crack_faces)
 

Detailed Description

Implementation of automatic differentiation.

Definition in file EshelbianPlasticity.cpp.

Macro Definition Documentation

◆ SINGULARITY

#define SINGULARITY

Definition at line 8 of file EshelbianPlasticity.cpp.

Function Documentation

◆ filter_owners()

static auto filter_owners ( MoFEM::Interface m_field,
Range  skin 
)
static
Examples
EshelbianPlasticity.cpp.

Definition at line 153 of file EshelbianPlasticity.cpp.

153  {
154  Range owner_ents;
155  ParallelComm *pcomm =
156  ParallelComm::get_pcomm(&m_field.get_moab(), MYPCOMM_INDEX);
157  CHK_MOAB_THROW(pcomm->filter_pstatus(skin, PSTATUS_NOT_OWNED, PSTATUS_NOT, -1,
158  &owner_ents),
159  "filter_pstatus");
160  return owner_ents;
161 };

◆ filter_true_skin()

static auto filter_true_skin ( MoFEM::Interface m_field,
Range &&  skin 
)
static
Examples
adolc_plasticity.cpp, EshelbianPlasticity.cpp, plastic.cpp, seepage.cpp, tensor_divergence_operator.cpp, and thermo_elastic.cpp.

Definition at line 142 of file EshelbianPlasticity.cpp.

142  {
143  Range boundary_ents;
144  ParallelComm *pcomm =
145  ParallelComm::get_pcomm(&m_field.get_moab(), MYPCOMM_INDEX);
146  CHK_MOAB_THROW(pcomm->filter_pstatus(skin,
147  PSTATUS_SHARED | PSTATUS_MULTISHARED,
148  PSTATUS_NOT, -1, &boundary_ents),
149  "filter_pstatus");
150  return boundary_ents;
151 };

◆ get_block_meshset()

static auto get_block_meshset ( MoFEM::Interface m_field,
const int  ms_id,
const unsigned int  cubit_bc_type 
)
static
Examples
EshelbianPlasticity.cpp.

Definition at line 120 of file EshelbianPlasticity.cpp.

121  {
122  auto mesh_mng = m_field.getInterface<MeshsetsManager>();
123  EntityHandle meshset;
124  CHKERR mesh_mng->getMeshset(ms_id, cubit_bc_type, meshset);
125  return meshset;
126 };

◆ get_crack_front_edges()

static auto get_crack_front_edges ( MoFEM::Interface m_field,
Range  crack_faces 
)
static
Examples
EshelbianPlasticity.cpp.

Definition at line 170 of file EshelbianPlasticity.cpp.

171  {
172  ParallelComm *pcomm =
173  ParallelComm::get_pcomm(&m_field.get_moab(), MYPCOMM_INDEX);
174  auto &moab = m_field.get_moab();
175  Range crack_skin;
176  if (pcomm->rank() == 0) {
177  crack_skin = get_skin(m_field, crack_faces);
178  Range body_ents;
180  m_field.get_moab().get_entities_by_dimension(0, SPACE_DIM, body_ents),
181  "get_entities_by_dimension");
182  auto body_skin = get_skin(m_field, body_ents);
183  Range body_skin_edges;
184  CHK_MOAB_THROW(moab.get_adjacencies(body_skin, 1, true, body_skin_edges,
185  moab::Interface::UNION),
186  "get_adjacencies");
187  crack_skin = subtract(crack_skin, body_skin_edges);
188  }
189  return send_type(m_field, crack_skin, MBEDGE);
190 }

◆ get_range_from_block()

static auto get_range_from_block ( MoFEM::Interface m_field,
const std::string  block_name,
int  dim 
)
static
Examples
EshelbianPlasticity.cpp.

Definition at line 98 of file EshelbianPlasticity.cpp.

99  {
100  Range r;
101 
102  auto mesh_mng = m_field.getInterface<MeshsetsManager>();
103  auto bcs = mesh_mng->getCubitMeshsetPtr(
104 
105  std::regex((boost::format("%s(.*)") % block_name).str())
106 
107  );
108 
109  for (auto bc : bcs) {
110  Range faces;
111  CHK_MOAB_THROW(bc->getMeshsetIdEntitiesByDimension(m_field.get_moab(), dim,
112  faces, true),
113  "get meshset ents");
114  r.merge(faces);
115  }
116 
117  return r;
118 };

◆ get_skin()

static auto get_skin ( MoFEM::Interface m_field,
Range  body_ents 
)
static
Examples
adolc_plasticity.cpp, EshelbianPlasticity.cpp, plastic.cpp, plate.cpp, seepage.cpp, tensor_divergence_operator.cpp, test_cache_on_entities.cpp, and thermo_elastic.cpp.

Definition at line 163 of file EshelbianPlasticity.cpp.

163  {
164  Skinner skin(&m_field.get_moab());
165  Range skin_ents;
166  CHK_MOAB_THROW(skin.find_skin(0, body_ents, false, skin_ents), "find_skin");
167  return skin_ents;
168 };

◆ get_two_sides_of_crack_surface()

static auto get_two_sides_of_crack_surface ( MoFEM::Interface m_field,
Range  crack_faces 
)
static
Examples
EshelbianPlasticity.cpp.

Definition at line 192 of file EshelbianPlasticity.cpp.

193  {
194 
195  ParallelComm *pcomm =
196  ParallelComm::get_pcomm(&m_field.get_moab(), MYPCOMM_INDEX);
197 
198  MOFEM_LOG("EP", Sev::inform) << "get_two_sides_of_crack_surface";
199 
200  if (!pcomm->rank()) {
201 
202  auto impl = [&](auto &saids) {
204 
205  auto &moab = m_field.get_moab();
206 
207  auto get_adj = [&](auto e, auto dim) {
208  Range adj;
209  CHK_MOAB_THROW(m_field.get_moab().get_adjacencies(
210  e, dim, true, adj, moab::Interface::UNION),
211  "get adj");
212  return adj;
213  };
214 
215  auto get_conn = [&](auto e) {
216  Range conn;
217  CHK_MOAB_THROW(m_field.get_moab().get_connectivity(e, conn, true),
218  "get connectivity");
219  return conn;
220  };
221 
222  constexpr bool debug = false;
223  Range body_ents;
224  CHKERR m_field.get_moab().get_entities_by_dimension(0, SPACE_DIM,
225  body_ents);
226  auto body_skin = get_skin(m_field, body_ents);
227  auto body_skin_edges = get_adj(body_skin, 1);
228 
229  auto crack_skin =
230  subtract(get_skin(m_field, crack_faces), body_skin_edges);
231  auto crack_skin_conn = get_conn(crack_skin);
232  auto crack_skin_conn_edges = get_adj(crack_skin_conn, 1);
233  auto crack_edges = get_adj(crack_faces, 1);
234  crack_edges = subtract(crack_edges, crack_skin);
235  auto all_tets = get_adj(crack_edges, 3);
236  crack_edges = subtract(crack_edges, crack_skin_conn_edges);
237  auto crack_conn = get_conn(crack_edges);
238  all_tets.merge(get_adj(crack_conn, 3));
239 
240  if (debug) {
241  CHKERR save_range(m_field.get_moab(), "crack_faces.vtk", crack_faces);
242  CHKERR save_range(m_field.get_moab(), "all_crack_tets.vtk", all_tets);
243  CHKERR save_range(m_field.get_moab(), "crack_edges_all.vtk",
244  crack_edges);
245  }
246 
247  if (crack_faces.size()) {
248  auto grow = [&](auto r) {
249  auto crack_faces_conn = get_conn(crack_faces);
250  Range v;
251  auto size_r = 0;
252  while (size_r != r.size() && r.size() > 0) {
253  size_r = r.size();
254  CHKERR moab.get_connectivity(r, v, true);
255  v = subtract(v, crack_faces_conn);
256  if (v.size()) {
257  CHKERR moab.get_adjacencies(v, SPACE_DIM, true, r,
258  moab::Interface::UNION);
259  r = intersect(r, all_tets);
260  }
261  if (r.empty()) {
262  break;
263  }
264  }
265  return r;
266  };
267 
268  Range all_tets_ord = all_tets;
269  while (all_tets.size()) {
270  Range faces = get_adj(unite(saids.first, saids.second), 2);
271  faces = subtract(crack_faces, faces);
272  if (faces.size()) {
273  Range tets;
274  auto fit = faces.begin();
275  for (;fit != faces.end(); ++fit) {
276  tets = intersect(get_adj(Range(*fit, *fit), 3), all_tets);
277  if (tets.size() == 2) {
278  break;
279  }
280  }
281  if (tets.empty()) {
282  break;
283  } else {
284  saids.first.insert(tets[0]);
285  saids.first = grow(saids.first);
286  all_tets = subtract(all_tets, saids.first);
287  if (tets.size() == 2) {
288  saids.second.insert(tets[1]);
289  saids.second = grow(saids.second);
290  all_tets = subtract(all_tets, saids.second);
291  }
292  }
293  } else {
294  break;
295  }
296  }
297 
298  saids.first = subtract(all_tets_ord, saids.second);
299  saids.second = subtract(all_tets_ord, saids.first);
300 
301  }
302 
304  };
305 
306  std::pair<Range, Range> saids;
307  if (crack_faces.size())
308  CHK_THROW_MESSAGE(impl(saids), "get crack both sides");
309  return saids;
310  }
311 
312  MOFEM_LOG("EP", Sev::inform) << "get_two_sides_of_crack_surface <- done";
313 
314  return std::pair<Range, Range>();
315 }

◆ save_range()

static auto save_range ( moab::Interface &  moab,
const std::string  name,
const Range  r,
std::vector< Tag >  tags = {} 
)
static
Examples
EshelbianPlasticity.cpp.

Definition at line 128 of file EshelbianPlasticity.cpp.

129  {}) {
131  auto out_meshset = get_temp_meshset_ptr(moab);
132  CHKERR moab.add_entities(*out_meshset, r);
133  if (r.size()) {
134  CHKERR moab.write_file(name.c_str(), "VTK", "", out_meshset->get_ptr(), 1,
135  tags.data(), tags.size());
136  } else {
137  MOFEM_LOG("SELF", Sev::warning) << "Empty range for " << name;
138  }
140 };

◆ send_type()

static auto send_type ( MoFEM::Interface m_field,
Range  r,
const EntityType  type 
)
static
Examples
EshelbianPlasticity.cpp.

Definition at line 52 of file EshelbianPlasticity.cpp.

53  {
54  ParallelComm *pcomm =
55  ParallelComm::get_pcomm(&m_field.get_moab(), MYPCOMM_INDEX);
56 
57  auto dim = CN::Dimension(type);
58 
59  std::vector<int> sendcounts(pcomm->size());
60  std::vector<int> displs(pcomm->size());
61  std::vector<int> sendbuf(r.size());
62  if (pcomm->rank() == 0) {
63  for (auto p = 1; p != pcomm->size(); p++) {
64  auto part_ents = m_field.getInterface<CommInterface>()
65  ->getPartEntities(m_field.get_moab(), p)
66  .subset_by_dimension(SPACE_DIM);
67  Range faces;
68  CHKERR m_field.get_moab().get_adjacencies(part_ents, dim, true, faces,
69  moab::Interface::UNION);
70  faces = intersect(faces, r);
71  sendcounts[p] = faces.size();
72  displs[p] = sendbuf.size();
73  for (auto f : faces) {
74  auto id = id_from_handle(f);
75  sendbuf.push_back(id);
76  }
77  }
78  }
79 
80  int recv_data;
81  MPI_Scatter(sendcounts.data(), 1, MPI_INT, &recv_data, 1, MPI_INT, 0,
82  pcomm->comm());
83  std::vector<int> recvbuf(recv_data);
84  MPI_Scatterv(sendbuf.data(), sendcounts.data(), displs.data(), MPI_INT,
85  recvbuf.data(), recv_data, MPI_INT, 0, pcomm->comm());
86 
87  if (pcomm->rank() > 0) {
88  Range r;
89  for (auto &f : recvbuf) {
90  r.insert(ent_form_type_and_id(type, f));
91  }
92  return r;
93  }
94 
95  return r;
96 }
MoFEM::UnknownInterface::getInterface
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface reference to pointer of interface.
Definition: UnknownInterface.hpp:93
CHK_MOAB_THROW
#define CHK_MOAB_THROW(err, msg)
Check error code of MoAB function and throw MoFEM exception.
Definition: definitions.h:589
MoFEM::id_from_handle
auto id_from_handle(const EntityHandle h)
Definition: Templates.hpp:1890
MYPCOMM_INDEX
#define MYPCOMM_INDEX
default communicator number PCOMM
Definition: definitions.h:228
save_range
static auto save_range(moab::Interface &moab, const std::string name, const Range r, std::vector< Tag > tags={})
Definition: EshelbianPlasticity.cpp:128
EntityHandle
CHK_THROW_MESSAGE
#define CHK_THROW_MESSAGE(err, msg)
Check and throw MoFEM exception.
Definition: definitions.h:609
sdf.r
int r
Definition: sdf.py:8
CHKERR
#define CHKERR
Inline error check.
Definition: definitions.h:548
send_type
static auto send_type(MoFEM::Interface &m_field, Range r, const EntityType type)
Definition: EshelbianPlasticity.cpp:52
MoFEM::CoreInterface::get_moab
virtual moab::Interface & get_moab()=0
SPACE_DIM
constexpr int SPACE_DIM
Definition: child_and_parent.cpp:16
MoFEM::get_temp_meshset_ptr
auto get_temp_meshset_ptr(moab::Interface &moab)
Create smart pointer to temporary meshset.
Definition: Templates.hpp:1886
convert.type
type
Definition: convert.py:64
get_skin
static auto get_skin(MoFEM::Interface &m_field, Range body_ents)
Definition: EshelbianPlasticity.cpp:163
debug
static const bool debug
Definition: dm_create_subdm.cpp:12
v
const double v
phase velocity of light in medium (cm/ns)
Definition: initial_diffusion.cpp:40
Range
MOFEM_LOG
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:308
MoFEM::CommInterface
Managing BitRefLevels.
Definition: CommInterface.hpp:21
HenckyOps::f
auto f
Definition: HenckyOps.hpp:15
MoFEM::MeshsetsManager
Interface for managing meshsets containing materials and boundary conditions.
Definition: MeshsetsManager.hpp:104
MoFEM::MeshsetsManager::getCubitMeshsetPtr
MoFEMErrorCode getCubitMeshsetPtr(const int ms_id, const CubitBCType cubit_bc_type, const CubitMeshSets **cubit_meshset_ptr) const
get cubit meshset
Definition: MeshsetsManager.cpp:589
MoFEM::ent_form_type_and_id
auto ent_form_type_and_id(const EntityType type, const EntityID id)
get entity handle from type and id
Definition: Templates.hpp:1906
MoFEMFunctionReturn
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
MoFEMFunctionBegin
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359