v0.16.3
Loading...
Searching...
No Matches
Classes | Namespaces | Macros | Functions
EshelbianPlasticity.cpp File Reference
#include <IntegrationRules.hpp>
#include <MoFEM.hpp>
#include <CGGTonsorialBubbleBase.hpp>
#include <EshelbianPlasticity.hpp>
#include <EshelbianRestart.hpp>
#include <boost/math/constants/constants.hpp>
#include <cholesky.hpp>
#include <EshelbianAux.hpp>
#include <EshelbianContact.hpp>
#include <PlasticIncrementalOptimization.hpp>
#include <impl/PlasticIncrementalOptimizationInternal.hpp>
#include <EshelbianTopologicalDerivative.hpp>
#include <TSElasticPostStep.hpp>
#include <phg-quadrule/quad.h>
#include <queue>
#include "impl/CGGUserPolynomialBase.cpp"
#include "impl/EshelbianMonitor.cpp"
#include "impl/EshelbianTestingMonitor.cpp"
#include "impl/SetUpSchurImpl.cpp"
#include <impl/EshelbianFracture.cpp>

Go to the source code of this file.

Classes

struct  EshelbianPlasticity::SetIntegrationAtFrontVolume
 
struct  EshelbianPlasticity::SetIntegrationAtFrontVolume::Fe
 
struct  EshelbianPlasticity::SetIntegrationAtFrontFace
 
struct  EshelbianPlasticity::SetIntegrationAtFrontFace::Fe
 
struct  EshelbianPlasticity::OpApplyPlasticFlowIncrement
 
struct  EshelbianPlasticity::solve_elastic_setup
 

Namespaces

namespace  EshelbianPlasticity
 

Macros

#define SINGULARITY
 

Functions

static auto send_type (MoFEM::Interface &m_field, Range r, const EntityType type)
 
static auto get_entities_by_handle (MoFEM::Interface &m_field, const std::string block_name)
 
static auto get_range_from_block (MoFEM::Interface &m_field, const std::string block_name, int dim)
 
static auto get_range_from_block_map (MoFEM::Interface &m_field, const std::string block_name, int dim)
 
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)
 
auto EshelbianPlasticity::vol_rule (int o)
 
auto EshelbianPlasticity::face_rule (int o)
 
static MoFEMErrorCode EshelbianPlasticity::RelaxationResidualMonitor (TS ts, PetscInt, PetscReal, Vec, void *)
 
static MoFEMErrorCode EshelbianPlasticity::checkDynamicToleranceCompatibility (TS ts)
 

Macro Definition Documentation

◆ SINGULARITY

#define SINGULARITY

Definition at line 12 of file EshelbianPlasticity.cpp.

Function Documentation

◆ filter_owners()

static auto filter_owners ( MoFEM::Interface &  m_field,
Range  skin 
)
static
Examples
/home/lk58p/mofem_install/vanilla_dev_release/mofem-cephas/mofem/users_modules/eshelbian_plasticity/src/impl/EshelbianPlasticity.cpp.

Definition at line 184 of file EshelbianPlasticity.cpp.

184 {
185 Range owner_ents;
186 ParallelComm *pcomm =
187 ParallelComm::get_pcomm(&m_field.get_moab(), MYPCOMM_INDEX);
188 CHK_MOAB_THROW(pcomm->filter_pstatus(skin, PSTATUS_NOT_OWNED, PSTATUS_NOT, -1,
189 &owner_ents),
190 "filter_pstatus");
191 return owner_ents;
192};
#define MYPCOMM_INDEX
default communicator number PCOMM
#define CHK_MOAB_THROW(err, msg)
Check error code of MoAB function and throw MoFEM exception.
virtual moab::Interface & get_moab()=0

◆ filter_true_skin()

static auto filter_true_skin ( MoFEM::Interface &  m_field,
Range &&  skin 
)
static
Examples
/home/lk58p/mofem_install/vanilla_dev_release/mofem-cephas/mofem/users_modules/eshelbian_plasticity/src/impl/EshelbianPlasticity.cpp, mofem/atom_tests/tensor_divergence_operator.cpp, mofem/tutorials/adv-0_plasticity/plastic.cpp, mofem/tutorials/adv-5_poroelasticity/seepage.cpp, mofem/users_modules/adolc-plasticity/adolc_plasticity.cpp, plastic.cpp, thermo_elastic.cpp, and thermoplastic.cpp.

Definition at line 173 of file EshelbianPlasticity.cpp.

173 {
174 Range boundary_ents;
175 ParallelComm *pcomm =
176 ParallelComm::get_pcomm(&m_field.get_moab(), MYPCOMM_INDEX);
177 CHK_MOAB_THROW(pcomm->filter_pstatus(skin,
178 PSTATUS_SHARED | PSTATUS_MULTISHARED,
179 PSTATUS_NOT, -1, &boundary_ents),
180 "filter_pstatus");
181 return boundary_ents;
182};

◆ get_crack_front_edges()

static auto get_crack_front_edges ( MoFEM::Interface &  m_field,
Range  crack_faces 
)
static
Examples
/home/lk58p/mofem_install/vanilla_dev_release/mofem-cephas/mofem/users_modules/eshelbian_plasticity/src/impl/EshelbianPlasticity.cpp.

Definition at line 201 of file EshelbianPlasticity.cpp.

202 {
203 ParallelComm *pcomm =
204 ParallelComm::get_pcomm(&m_field.get_moab(), MYPCOMM_INDEX);
205 auto &moab = m_field.get_moab();
206 Range crack_skin_without_bdy;
207 if (pcomm->rank() == 0) {
208 Range crack_edges;
209 CHKERR moab.get_adjacencies(crack_faces, 1, true, crack_edges,
210 moab::Interface::UNION);
211 auto crack_skin = get_skin(m_field, crack_faces);
212 Range body_ents;
214 m_field.get_moab().get_entities_by_dimension(0, SPACE_DIM, body_ents),
215 "get_entities_by_dimension");
216 auto body_skin = get_skin(m_field, body_ents);
217 Range body_skin_edges;
218 CHK_MOAB_THROW(moab.get_adjacencies(body_skin, 1, true, body_skin_edges,
219 moab::Interface::UNION),
220 "get_adjacencies");
221 crack_skin_without_bdy = subtract(crack_skin, body_skin_edges);
222 auto front_edges_map = get_range_from_block_map(m_field, "FRONT", 1);
223 for (auto &m : front_edges_map) {
224 auto add_front = subtract(m.second, crack_edges);
225 auto i = intersect(m.second, crack_edges);
226 if (i.empty()) {
227 crack_skin_without_bdy.merge(add_front);
228 } else {
229 auto i_skin = get_skin(m_field, i);
230 Range adj_i_skin;
231 CHKERR moab.get_adjacencies(i_skin, 1, true, adj_i_skin,
232 moab::Interface::UNION);
233 adj_i_skin = subtract(intersect(adj_i_skin, m.second), crack_edges);
234 crack_skin_without_bdy.merge(adj_i_skin);
235 }
236 }
237 }
238 return send_type(m_field, crack_skin_without_bdy, MBEDGE);
239}
static auto send_type(MoFEM::Interface &m_field, Range r, const EntityType type)
static auto get_range_from_block_map(MoFEM::Interface &m_field, const std::string block_name, int dim)
static auto get_skin(MoFEM::Interface &m_field, Range body_ents)
constexpr int SPACE_DIM
#define CHKERR
Inline error check.
FTensor::Index< 'i', SPACE_DIM > i
FTensor::Index< 'm', 3 > m

◆ get_entities_by_handle()

static auto get_entities_by_handle ( MoFEM::Interface &  m_field,
const std::string  block_name 
)
static
Examples
/home/lk58p/mofem_install/vanilla_dev_release/mofem-cephas/mofem/users_modules/eshelbian_plasticity/src/impl/EshelbianPlasticity.cpp.

Definition at line 95 of file EshelbianPlasticity.cpp.

96 {
97 Range r;
98
99 auto mesh_mng = m_field.getInterface<MeshsetsManager>();
100 auto bcs = mesh_mng->getCubitMeshsetPtr(
101
102 std::regex((boost::format("%s(.*)") % block_name).str())
103
104 );
105
106 for (auto bc : bcs) {
107 auto meshset = bc->getMeshset();
108 CHK_MOAB_THROW(m_field.get_moab().get_entities_by_handle(meshset, r, true),
109 "get meshset ents");
110 }
111
112 return r;
113};
MoFEMErrorCode getCubitMeshsetPtr(const int ms_id, const CubitBCType cubit_bc_type, const CubitMeshSets **cubit_meshset_ptr) const
get cubit meshset
int r
Definition sdf.py:205
Interface for managing meshsets containing materials and boundary conditions.
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface reference to pointer of interface.

◆ get_range_from_block()

static auto get_range_from_block ( MoFEM::Interface &  m_field,
const std::string  block_name,
int  dim 
)
static

Definition at line 115 of file EshelbianPlasticity.cpp.

116 {
117 Range r;
118
119 auto mesh_mng = m_field.getInterface<MeshsetsManager>();
120 auto bcs = mesh_mng->getCubitMeshsetPtr(
121
122 std::regex((boost::format("%s(.*)") % block_name).str())
123
124 );
125
126 for (auto bc : bcs) {
127 Range faces;
128 CHK_MOAB_THROW(bc->getMeshsetIdEntitiesByDimension(m_field.get_moab(), dim,
129 faces, true),
130 "get meshset ents");
131 r.merge(faces);
132 }
133
134 return r;
135};

◆ get_range_from_block_map()

static auto get_range_from_block_map ( MoFEM::Interface &  m_field,
const std::string  block_name,
int  dim 
)
static
Examples
/home/lk58p/mofem_install/vanilla_dev_release/mofem-cephas/mofem/users_modules/eshelbian_plasticity/src/impl/EshelbianPlasticity.cpp.

Definition at line 137 of file EshelbianPlasticity.cpp.

138 {
139 std::map<std::string, Range> r;
140
141 auto mesh_mng = m_field.getInterface<MeshsetsManager>();
142 auto bcs = mesh_mng->getCubitMeshsetPtr(
143
144 std::regex((boost::format("%s(.*)") % block_name).str())
145
146 );
147
148 for (auto bc : bcs) {
149 Range faces;
150 CHK_MOAB_THROW(bc->getMeshsetIdEntitiesByDimension(m_field.get_moab(), dim,
151 faces, true),
152 "get meshset ents");
153 r[bc->getName()] = faces;
154 }
155
156 return r;
157}

◆ get_skin()

static auto get_skin ( MoFEM::Interface &  m_field,
Range  body_ents 
)
static

Definition at line 194 of file EshelbianPlasticity.cpp.

194 {
195 Skinner skin(&m_field.get_moab());
196 Range skin_ents;
197 CHK_MOAB_THROW(skin.find_skin(0, body_ents, false, skin_ents), "find_skin");
198 return skin_ents;
199};

◆ get_two_sides_of_crack_surface()

static auto get_two_sides_of_crack_surface ( MoFEM::Interface &  m_field,
Range  crack_faces 
)
static
Examples
/home/lk58p/mofem_install/vanilla_dev_release/mofem-cephas/mofem/users_modules/eshelbian_plasticity/src/impl/EshelbianPlasticity.cpp.

Definition at line 241 of file EshelbianPlasticity.cpp.

242 {
243
244 ParallelComm *pcomm =
245 ParallelComm::get_pcomm(&m_field.get_moab(), MYPCOMM_INDEX);
246
247 MOFEM_LOG("EP", Sev::noisy) << "get_two_sides_of_crack_surface";
248
249 if (!pcomm->rank()) {
250
251 auto impl = [&](auto &saids) {
253
254 auto &moab = m_field.get_moab();
255
256 auto get_adj = [&](auto e, auto dim) {
257 Range adj;
258 CHK_MOAB_THROW(m_field.get_moab().get_adjacencies(
259 e, dim, true, adj, moab::Interface::UNION),
260 "get adj");
261 return adj;
262 };
263
264 auto get_conn = [&](auto e) {
265 Range conn;
266 CHK_MOAB_THROW(m_field.get_moab().get_connectivity(e, conn, true),
267 "get connectivity");
268 return conn;
269 };
270
271 constexpr bool debug = false;
272 Range body_ents;
273 CHKERR m_field.get_moab().get_entities_by_dimension(0, SPACE_DIM,
274 body_ents);
275 auto body_skin = get_skin(m_field, body_ents);
276 auto body_skin_edges = get_adj(body_skin, 1);
277
278 auto crack_skin =
279 subtract(get_skin(m_field, crack_faces), body_skin_edges);
280 auto crack_skin_conn = get_conn(crack_skin);
281 auto crack_skin_conn_edges = get_adj(crack_skin_conn, 1);
282 auto crack_edges = get_adj(crack_faces, 1);
283 crack_edges = subtract(crack_edges, crack_skin);
284 auto all_tets = get_adj(crack_edges, 3);
285 crack_edges = subtract(crack_edges, crack_skin_conn_edges);
286 auto crack_conn = get_conn(crack_edges);
287 all_tets.merge(get_adj(crack_conn, 3));
288
289 if (debug) {
290 CHKERR save_range(m_field.get_moab(), "crack_faces.vtk", crack_faces);
291 CHKERR save_range(m_field.get_moab(), "all_crack_tets.vtk", all_tets);
292 CHKERR save_range(m_field.get_moab(), "crack_edges_all.vtk",
293 crack_edges);
294 }
295
296 if (crack_faces.size()) {
297 auto grow = [&](auto r) {
298 auto crack_faces_conn = get_conn(crack_faces);
299 Range v;
300 auto size_r = 0;
301 while (size_r != r.size() && r.size() > 0) {
302 size_r = r.size();
303 CHKERR moab.get_connectivity(r, v, true);
304 v = subtract(v, crack_faces_conn);
305 if (v.size()) {
306 CHKERR moab.get_adjacencies(v, SPACE_DIM, true, r,
307 moab::Interface::UNION);
308 r = intersect(r, all_tets);
309 }
310 if (r.empty()) {
311 break;
312 }
313 }
314 return r;
315 };
316
317 Range all_tets_ord = all_tets;
318 while (all_tets.size()) {
319 Range faces = get_adj(unite(saids.first, saids.second), 2);
320 faces = subtract(crack_faces, faces);
321 if (faces.size()) {
322 Range tets;
323 auto fit = faces.begin();
324 for (; fit != faces.end(); ++fit) {
325 tets = intersect(get_adj(Range(*fit, *fit), 3), all_tets);
326 if (tets.size() == 2) {
327 break;
328 }
329 }
330 if (tets.empty()) {
331 break;
332 } else {
333 saids.first.insert(tets[0]);
334 saids.first = grow(saids.first);
335 all_tets = subtract(all_tets, saids.first);
336 if (tets.size() == 2) {
337 saids.second.insert(tets[1]);
338 saids.second = grow(saids.second);
339 all_tets = subtract(all_tets, saids.second);
340 }
341 }
342 } else {
343 break;
344 }
345 }
346
347 saids.first = subtract(all_tets_ord, saids.second);
348 saids.second = subtract(all_tets_ord, saids.first);
349 }
350
352 };
353
354 std::pair<Range, Range> saids;
355 if (crack_faces.size())
356 CHK_THROW_MESSAGE(impl(saids), "get crack both sides");
357 return saids;
358 }
359
360 MOFEM_LOG("EP", Sev::noisy) << "get_two_sides_of_crack_surface <- done";
361
362 return std::pair<Range, Range>();
363}
#define CHK_THROW_MESSAGE(err, msg)
Check and throw MoFEM exception.
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define MOFEM_LOG(channel, severity)
Log.
const double v
phase velocity of light in medium (cm/ns)
static const bool debug
auto save_range

◆ save_range()

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

Definition at line 159 of file EshelbianPlasticity.cpp.

160 {}) {
162 auto out_meshset = get_temp_meshset_ptr(moab);
163 CHKERR moab.add_entities(*out_meshset, r);
164 if (r.size()) {
165 CHKERR moab.write_file(name.c_str(), "VTK", "", out_meshset->get_ptr(), 1,
166 tags.data(), tags.size());
167 } else {
168 MOFEM_LOG("SELF", Sev::warning) << "Empty range for " << name;
169 }
171};
auto get_temp_meshset_ptr(moab::Interface &moab)
Create smart pointer to temporary meshset.

◆ send_type()

static auto send_type ( MoFEM::Interface &  m_field,
Range  r,
const EntityType  type 
)
static
Examples
/home/lk58p/mofem_install/vanilla_dev_release/mofem-cephas/mofem/users_modules/eshelbian_plasticity/src/impl/EshelbianPlasticity.cpp.

Definition at line 49 of file EshelbianPlasticity.cpp.

50 {
51 ParallelComm *pcomm =
52 ParallelComm::get_pcomm(&m_field.get_moab(), MYPCOMM_INDEX);
53
54 auto dim = CN::Dimension(type);
55
56 std::vector<int> sendcounts(pcomm->size());
57 std::vector<int> displs(pcomm->size());
58 std::vector<int> sendbuf(r.size());
59 if (pcomm->rank() == 0) {
60 for (auto p = 1; p != pcomm->size(); p++) {
61 auto part_ents = m_field.getInterface<CommInterface>()
62 ->getPartEntities(m_field.get_moab(), p)
63 .subset_by_dimension(SPACE_DIM);
64 Range faces;
65 CHKERR m_field.get_moab().get_adjacencies(part_ents, dim, true, faces,
66 moab::Interface::UNION);
67 faces = intersect(faces, r);
68 sendcounts[p] = faces.size();
69 displs[p] = sendbuf.size();
70 for (auto f : faces) {
71 auto id = id_from_handle(f);
72 sendbuf.push_back(id);
73 }
74 }
75 }
76
77 int recv_data;
78 MPI_Scatter(sendcounts.data(), 1, MPI_INT, &recv_data, 1, MPI_INT, 0,
79 pcomm->comm());
80 std::vector<int> recvbuf(recv_data);
81 MPI_Scatterv(sendbuf.data(), sendcounts.data(), displs.data(), MPI_INT,
82 recvbuf.data(), recv_data, MPI_INT, 0, pcomm->comm());
83
84 if (pcomm->rank() > 0) {
85 Range r;
86 for (auto &f : recvbuf) {
87 r.insert(ent_form_type_and_id(type, f));
88 }
89 return r;
90 }
91
92 return r;
93}
std::string type
auto id_from_handle(const EntityHandle h)
auto ent_form_type_and_id(const EntityType type, const EntityID id)
get entity handle from type and id
Managing BitRefLevels.