v0.9.1
Public Member Functions | Public Attributes | List of all members
EshelbianPlasticity::EshelbianCore Struct Reference

#include <users_modules/eshelbian_plasticty/src/EshelbianPlasticity.hpp>

Inheritance diagram for EshelbianPlasticity::EshelbianCore:
[legend]
Collaboration diagram for EshelbianPlasticity::EshelbianCore:
[legend]

Public Member Functions

MoFEMErrorCode query_interface (const MOFEMuuid &uuid, UnknownInterface **iface) const
 Getting interface of core database. More...
 
 EshelbianCore (MoFEM::Interface &m_field)
 
virtual ~EshelbianCore ()
 
MoFEMErrorCode getOptions ()
 
template<typename BC >
MoFEMErrorCode getBc (boost::shared_ptr< BC > &bc_vec_ptr, const std::string block_set_name, const int nb_attributes)
 
MoFEMErrorCode getSpatialDispBc ()
 
MoFEMErrorCode getSpatialRotationBc ()
 
MoFEMErrorCode getSpatialTractionBc ()
 
MoFEMErrorCode getTractionFreeBc (const EntityHandle meshset, boost::shared_ptr< TractionFreeBc > &bc_ptr, const std::string disp_block_set_name, const std::string rot_block_set_name)
 
MoFEMErrorCode getSpatialTractionFreeBc (const EntityHandle meshset=0)
 
MoFEMErrorCode addFields (const EntityHandle meshset=0)
 
MoFEMErrorCode addVolumeFiniteElement (const EntityHandle meshset=0)
 
MoFEMErrorCode addBoundaryFiniteElement (const EntityHandle meshset=0)
 
MoFEMErrorCode addDMs (const BitRefLevel bit=BitRefLevel().set(0))
 
MoFEMErrorCode addMaterial_HMHHStVenantKirchhoff (const int tape, const double lambda, const double mu, const double sigma_y)
 
MoFEMErrorCode addMaterial_HMHMooneyRivlin (const int tape, const double alpha, const double beta, const double lambda, const double sigma_y)
 
MoFEMErrorCode setBaseVolumeElementOps (const int tag, const bool do_rhs, const bool do_lhs, boost::shared_ptr< EpElement< VolumeElementForcesAndSourcesCore >> &fe)
 
MoFEMErrorCode setGenericVolumeElementOps (const int tag, const bool add_elastic, const bool add_material, boost::shared_ptr< EpElement< VolumeElementForcesAndSourcesCore >> &fe_rhs, boost::shared_ptr< EpElement< VolumeElementForcesAndSourcesCore >> &fe_lhs)
 
MoFEMErrorCode setGenericFaceElementOps (const bool add_elastic, const bool add_material, boost::shared_ptr< EpElement< FaceElementForcesAndSourcesCore >> &fe_rhs, boost::shared_ptr< EpElement< FaceElementForcesAndSourcesCore >> &fe_lhs)
 
MoFEMErrorCode setElasticElementOps (const int tag)
 
MoFEMErrorCode setElasticElementToTs (DM dm)
 
MoFEMErrorCode setUpTSElastic (TS ts, Mat m, Vec f)
 
MoFEMErrorCode solveElastic (TS ts, Vec x)
 
MoFEMErrorCode postProcessResults (const int tag, const std::string file)
 
- Public Member Functions inherited from MoFEM::UnknownInterface
virtual MoFEMErrorCode query_interface (const MOFEMuuid &uuid, UnknownInterface **iface) const =0
 
template<class IFACE >
MoFEMErrorCode registerInterface (const MOFEMuuid &uuid, bool error_if_registration_failed=true)
 Register interface. More...
 
template<class IFACE , bool VERIFY = false>
MoFEMErrorCode getInterface (const MOFEMuuid &uuid, IFACE *&iface) const
 Get interface by uuid and return reference to pointer of interface. More...
 
template<class IFACE >
MoFEMErrorCode getInterface (IFACE *&iface) const
 Get interface refernce to pointer of interface. More...
 
template<class IFACE >
MoFEMErrorCode getInterface (IFACE **const iface) const
 Get interface pointer to pointer of interface. More...
 
template<class IFACE , typename boost::enable_if< boost::is_pointer< IFACE >, int >::type = 0>
IFACE getInterface () const
 Get interface pointer to pointer of interface. More...
 
template<class IFACE , typename boost::enable_if< boost::is_reference< IFACE >, int >::type = 0>
IFACE getInterface () const
 Get reference to interface. More...
 
template<class IFACE >
IFACE * getInterface () const
 Function returning pointer to interface. More...
 
virtual ~UnknownInterface ()=default
 
virtual MoFEMErrorCode getLibVersion (Version &version) const
 Get library version. More...
 
virtual const MoFEMErrorCode getFileVersion (moab::Interface &moab, Version &version) const
 Get database major version. More...
 
virtual MoFEMErrorCode getInterfaceVersion (Version &version) const
 Get database major version. More...
 
template<>
MoFEMErrorCode getInterface (const MOFEMuuid &uuid, UnknownInterface *&iface) const
 

Public Attributes

MoFEM::InterfacemField
 
boost::shared_ptr< DataAtIntegrationPtsdataAtPts
 
boost::shared_ptr< PhysicalEquationsphysicalEquations
 
boost::shared_ptr< EpElement< VolumeElementForcesAndSourcesCore > > elasticFeRhs
 
boost::shared_ptr< EpElement< VolumeElementForcesAndSourcesCore > > elasticFeLhs
 
boost::shared_ptr< EpElement< FaceElementForcesAndSourcesCore > > elasticBcLhs
 
boost::shared_ptr< EpElement< FaceElementForcesAndSourcesCore > > elasticBcRhs
 
boost::shared_ptr< EpFEMethodschurAssembly
 
SmartPetscObj< DM > dM
 Coupled problem all fields. More...
 
SmartPetscObj< DM > dmElastic
 Elastic problem. More...
 
SmartPetscObj< DM > dmElasticSchurStreach
 Sub problem of dmElastic Schur. More...
 
SmartPetscObj< DM > dmElasticSchurBubble
 Sub problem of dmElastic Schur. More...
 
SmartPetscObj< DM > dmElasticSchurSpatialDisp
 Sub problem of dmElastic Schur. More...
 
SmartPetscObj< DM > dmElasticSchurOmega
 Sub problem of dmElastic Schur. More...
 
SmartPetscObj< DM > dmMaterial
 Material problem. More...
 
const std::string piolaStress
 
const std::string eshelbyStress
 
const std::string spatialDisp
 
const std::string materialDisp
 
const std::string streachTensor
 
const std::string rotAxis
 
const std::string materialGradient
 
const std::string tauField
 
const std::string lambdaField
 
const std::string bubbleField
 
const std::string elementVolumeName
 
const std::string naturalBcElement
 
const std::string essentialBcElement
 
int spaceOrder
 
int materialOrder
 
double alpha_u
 
double alpha_w
 
double preconditioner_eps
 
boost::shared_ptr< BcDispVecbcSpatialDispVecPtr
 
boost::shared_ptr< BcRotVecbcSpatialRotationVecPtr
 
boost::shared_ptr< TractionBcVecbcSpatialTraction
 
boost::shared_ptr< TractionFreeBcbcSpatialFreeTraction
 

Additional Inherited Members

- Protected Member Functions inherited from MoFEM::UnknownInterface
boost::typeindex::type_index getClassIdx (const MOFEMuuid &uid) const
 Get type name for interface Id. More...
 
MOFEMuuid getUId (const boost::typeindex::type_index &class_idx) const
 Get interface Id for class name. More...
 

Detailed Description

Examples
ep.cpp.

Definition at line 1178 of file EshelbianPlasticity.hpp.

Constructor & Destructor Documentation

◆ EshelbianCore()

EshelbianPlasticity::EshelbianCore::EshelbianCore ( MoFEM::Interface m_field)

Definition at line 45 of file EshelbianPlasticity.cpp.

46  : mField(m_field), piolaStress("P"), eshelbyStress("S"),
47  spatialDisp("w"), materialDisp("W"),
48  streachTensor("u"), rotAxis("omega"), materialGradient("G"),
49  tauField("TAU"), lambdaField("LAMBDA"), bubbleField("BUBBLE"),
50  elementVolumeName("EP"), naturalBcElement("NATURAL_BC"),
51  essentialBcElement("ESSENTIAL_BC") {
52 
53  ierr = getOptions();
54  CHKERRABORT(PETSC_COMM_WORLD, ierr);
55 }

◆ ~EshelbianCore()

EshelbianPlasticity::EshelbianCore::~EshelbianCore ( )
virtual

Definition at line 57 of file EshelbianPlasticity.cpp.

57  {
58 }

Member Function Documentation

◆ addBoundaryFiniteElement()

MoFEMErrorCode EshelbianPlasticity::EshelbianCore::addBoundaryFiniteElement ( const EntityHandle  meshset = 0)
Examples
ep.cpp.

Definition at line 261 of file EshelbianPlasticity.cpp.

261  {
263 
264  auto bc_elements_add_to_range = [&](const std::string disp_block_set_name,
265  Range &r) {
268  if (it->getName().compare(0, disp_block_set_name.length(),
269  disp_block_set_name) == 0) {
270  Range faces;
271  CHKERR it->getMeshsetIdEntitiesByDimension(mField.get_moab(), 2, faces,
272  true);
273  r.merge(faces);
274  }
275  }
277  };
278 
279  // set finite element fields
280  auto add_field_to_fe = [this](const std::string fe,
281  const std::string field_name) {
287  };
288 
289  Range natural_bc_elements;
290  CHKERR bc_elements_add_to_range("SPATIAL_DISP_BC", natural_bc_elements);
291  CHKERR bc_elements_add_to_range("SPATIAL_ROTATION_BC", natural_bc_elements);
292  Range essentail_bc_elements;
293  CHKERR bc_elements_add_to_range("SPATIAL_TRACTION_BC", essentail_bc_elements);
294 
296  CHKERR mField.add_ents_to_finite_element_by_type(natural_bc_elements, MBTRI,
298  CHKERR add_field_to_fe(naturalBcElement, piolaStress);
299  CHKERR add_field_to_fe(naturalBcElement, eshelbyStress);
301 
303  CHKERR mField.add_ents_to_finite_element_by_type(essentail_bc_elements, MBTRI,
305  CHKERR add_field_to_fe(essentialBcElement, piolaStress);
306  CHKERR add_field_to_fe(essentialBcElement, eshelbyStress);
308 
310 }
virtual MoFEMErrorCode modify_finite_element_add_field_col(const std::string &fe_name, const std::string &name_row)=0
set field col which finite element use
virtual MoFEMErrorCode add_ents_to_finite_element_by_type(const EntityHandle entities, const EntityType type, const std::string &name, const bool recursive=true)=0
add entities to finite element
virtual moab::Interface & get_moab()=0
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
virtual MoFEMErrorCode add_finite_element(const std::string &fe_name, enum MoFEMTypes bh=MF_EXCL, int verb=DEFAULT_VERBOSITY)=0
add finite element
#define _IT_CUBITMESHSETS_BY_SET_TYPE_FOR_LOOP_(MESHSET_MANAGER, CUBITBCTYPE, IT)
Iterator that loops over a specific Cubit MeshSet having a particular BC meshset in a moFEM field.
#define CHKERR
Inline error check.
Definition: definitions.h:601
virtual MoFEMErrorCode modify_finite_element_add_field_data(const std::string &fe_name, const std::string &name_filed)=0
set finite element field data
virtual MoFEMErrorCode modify_finite_element_add_field_row(const std::string &fe_name, const std::string &name_row)=0
set field row which finite element use
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412
virtual MoFEMErrorCode build_finite_elements(int verb=DEFAULT_VERBOSITY)=0
Build finite elementsBuild finite element data structures. Have to be run before problem and adjacenc...

◆ addDMs()

MoFEMErrorCode EshelbianPlasticity::EshelbianCore::addDMs ( const BitRefLevel  bit = BitRefLevel().set(0))
Examples
ep.cpp.

Definition at line 312 of file EshelbianPlasticity.cpp.

312  {
314 
315  // find adjacencies between finite elements and dofs
317 
318  // Create coupled problem
319  dM = createSmartDM(mField.get_comm(), "DMMOFEM");
320  CHKERR DMMoFEMCreateMoFEM(dM, &mField, "ESHELBY_PLASTICITY", bit,
321  BitRefLevel().set());
322  CHKERR DMMoFEMSetDestroyProblem(dM, PETSC_TRUE);
323  CHKERR DMMoFEMSetIsPartitioned(dM, PETSC_TRUE);
327  mField.getInterface<ProblemsManager>()->buildProblemFromFields = PETSC_TRUE;
328  CHKERR DMSetUp(dM);
329  mField.getInterface<ProblemsManager>()->buildProblemFromFields = PETSC_FALSE;
330 
331  auto remove_dofs_on_essential_spatial_stress_boundary =
332  [&](const std::string prb_name) {
334  for (int d : {0, 1, 2})
336  prb_name, piolaStress, (*bcSpatialFreeTraction)[d], d, d, NOISY,
337  true);
339  };
340  CHKERR remove_dofs_on_essential_spatial_stress_boundary("ESHELBY_PLASTICITY");
341 
342  // Create elastic sub-problem
343  dmElastic = createSmartDM(mField.get_comm(), "DMMOFEM");
344  CHKERR DMMoFEMCreateSubDM(dmElastic, dM, "ELASTIC_PROBLEM");
356  CHKERR DMSetUp(dmElastic);
357 
358  // Create elastic streach-problem
361  "ELASTIC_PROBLEM_STREACH_SCHUR");
372  CHKERR DMSetUp(dmElasticSchurStreach);
373 
374  // Create elastic bubble-problem
377  "ELASTIC_PROBLEM_BUBBLE_SCHUR");
387  CHKERR DMSetUp(dmElasticSchurBubble);
388 
389  // Create elastic omega-problem
392  "ELASTIC_PROBLEM_OMEGA_SCHUR");
401  CHKERR DMSetUp(dmElasticSchurOmega);
402 
403  // Create elastic tet_stress-problem
406  "ELASTIC_PROBLEM_SPATIAL_DISP_SCHUR");
410  elementVolumeName.c_str());
413  essentialBcElement.c_str());
417 
418  {
419  PetscSection section;
420  CHKERR mField.getInterface<ISManager>()->sectionCreate("ELASTIC_PROBLEM",
421  &section);
422  CHKERR DMSetDefaultSection(dmElastic, section);
423  CHKERR DMSetDefaultGlobalSection(dmElastic, section);
424  CHKERR PetscSectionDestroy(&section);
425  }
426 
428 }
SmartPetscObj< DM > dmElastic
Elastic problem.
Problem manager is used to build and partition problems.
virtual MoFEMErrorCode build_adjacencies(const Range &ents, int verb=DEFAULT_VERBOSITY)=0
build adjacencies
SmartPetscObj< DM > dmElasticSchurOmega
Sub problem of dmElastic Schur.
Section manager is used to create indexes and sectionsFIXME: ISManager is not properly testsed by ato...
Definition: ISManager.hpp:36
PetscErrorCode DMMoFEMSetSquareProblem(DM dm, PetscBool square_problem)
set squared problemIt if true is assumed that matrix has the same indexing on rows and columns....
Definition: DMMMoFEM.cpp:378
PetscErrorCode DMMoFEMAddElement(DM dm, const char fe_name[])
add element to dm
Definition: DMMMoFEM.cpp:411
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
MoFEMErrorCode removeDofsOnEntities(const std::string problem_name, const std::string field_name, const Range ents, const int lo_coeff=0, const int hi_coeff=MAX_DOFS_ON_ENTITY, int verb=VERBOSE, const bool debug=false)
Remove DOFs from problemRemove DOFs from problem which are on entities on the given range and given f...
PetscErrorCode DMMoFEMCreateSubDM(DM subdm, DM dm, const char problem_name[])
Must be called by user to set Sub DM MoFEM data structures.
Definition: DMMMoFEM.cpp:177
PetscErrorCode DMMoFEMAddSubFieldRow(DM dm, const char field_name[], EntityType lo_type=MBVERTEX, EntityType hi_type=MBMAXTYPE)
Definition: DMMMoFEM.cpp:198
MoFEMErrorCode getInterface(const MOFEMuuid &uuid, IFACE *&iface) const
Get interface by uuid and return reference to pointer of interface.
const Tensor1_Expr< const dTensor0< T, Dim, i >, typename promote< T, double >::V, Dim, i > d(const Tensor0< T * > &a, const Index< i, Dim > index, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition: dTensor0.hpp:27
SmartPetscObj< DM > dmElasticSchurStreach
Sub problem of dmElastic Schur.
SmartPetscObj< DM > dmElasticSchurSpatialDisp
Sub problem of dmElastic Schur.
boost::shared_ptr< TractionFreeBc > bcSpatialFreeTraction
PetscErrorCode DMMoFEMCreateMoFEM(DM dm, MoFEM::Interface *m_field_ptr, const char problem_name[], const MoFEM::BitRefLevel bit_level, const MoFEM::BitRefLevel bit_mask=MoFEM::BitRefLevel().set())
Must be called by user to set MoFEM data structures.
Definition: DMMMoFEM.cpp:105
#define CHKERR
Inline error check.
Definition: definitions.h:601
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:51
SmartPetscObj< DM > dM
Coupled problem all fields.
SmartPetscObj< DM > dmElasticSchurBubble
Sub problem of dmElastic Schur.
PetscErrorCode DMMoFEMSetDestroyProblem(DM dm, PetscBool destroy_problem)
Definition: DMMMoFEM.cpp:360
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412
virtual MPI_Comm & get_comm() const =0
auto createSmartDM
Creates smart DM object.
Definition: AuxPETSc.hpp:174
PetscErrorCode DMMoFEMSetIsPartitioned(DM dm, PetscBool is_partitioned)
Definition: DMMMoFEM.cpp:968

◆ addFields()

MoFEMErrorCode EshelbianPlasticity::EshelbianCore::addFields ( const EntityHandle  meshset = 0)
Examples
ep.cpp.

Definition at line 90 of file EshelbianPlasticity.cpp.

90  {
92 
93  Range tets;
94  CHKERR mField.get_moab().get_entities_by_type(meshset, MBTET, tets);
95  Range tets_skin_part;
96  Skinner skin(&mField.get_moab());
97  CHKERR skin.find_skin(0, tets, false, tets_skin_part);
98  ParallelComm *pcomm =
99  ParallelComm::get_pcomm(&mField.get_moab(), MYPCOMM_INDEX);
100  Range tets_skin;
101  CHKERR pcomm->filter_pstatus(tets_skin_part,
102  PSTATUS_SHARED | PSTATUS_MULTISHARED,
103  PSTATUS_NOT, -1, &tets_skin);
104 
105  auto subtract_faces_where_displacements_are_applied =
106  [&](const std::string disp_block_set_name) {
109  if (it->getName().compare(0, disp_block_set_name.length(),
110  disp_block_set_name) == 0) {
111  Range faces;
112  CHKERR it->getMeshsetIdEntitiesByDimension(mField.get_moab(), 2,
113  faces, true);
114  tets_skin = subtract(tets_skin, faces);
115  }
116  }
118  };
119  CHKERR subtract_faces_where_displacements_are_applied("SPATIAL_DISP_BC");
120  CHKERR subtract_faces_where_displacements_are_applied("SPATIAL_ROTATION_BC");
121  CHKERR subtract_faces_where_displacements_are_applied("SPATIAL_TRACTION_BC");
122 
123  Range faces;
124  CHKERR mField.get_moab().get_adjacencies(tets, 2, true, faces,
125  moab::Interface::UNION);
126  Range faces_not_on_the_skin = subtract(faces, tets_skin);
127 
128  auto add_hdiv_field = [&](const std::string field_name, const int order,
129  const int dim) {
132  MB_TAG_SPARSE, MF_ZERO);
133  CHKERR mField.add_ents_to_field_by_type(meshset, MBTET, field_name);
134  CHKERR mField.set_field_order(meshset, MBTET, field_name, order);
135  CHKERR mField.set_field_order(faces_not_on_the_skin, field_name, order);
136  CHKERR mField.set_field_order(tets_skin, field_name, 0);
138  };
139 
140  auto add_hdiv_rt_field = [&](const std::string field_name, const int order,
141  const int dim) {
144  MB_TAG_DENSE, MF_ZERO);
145  CHKERR mField.add_ents_to_field_by_type(meshset, MBTET, field_name);
146  CHKERR mField.set_field_order(meshset, MBTET, field_name, 0);
147  CHKERR mField.set_field_order(tets_skin, field_name, order);
149  };
150 
151  auto add_l2_field = [this, meshset](const std::string field_name,
152  const int order, const int dim) {
155  MB_TAG_DENSE, MF_ZERO);
156  CHKERR mField.add_ents_to_field_by_type(meshset, MBTET, field_name);
157  CHKERR mField.set_field_order(meshset, MBTET, field_name, order);
159  };
160 
161  auto add_h1_field = [this, meshset](const std::string field_name,
162  const int order, const int dim) {
165  MB_TAG_DENSE, MF_ZERO);
166  CHKERR mField.add_ents_to_field_by_type(meshset, MBTET, field_name);
167  CHKERR mField.set_field_order(meshset, MBVERTEX, field_name, 1);
168  CHKERR mField.set_field_order(meshset, MBEDGE, field_name, order);
169  CHKERR mField.set_field_order(meshset, MBTRI, field_name, order);
170  CHKERR mField.set_field_order(meshset, MBTET, field_name, order);
172  };
173 
174  auto add_bubble_field = [this, meshset](const std::string field_name,
175  const int order, const int dim) {
177  CHKERR mField.add_field(field_name, HDIV, USER_BASE, dim, MB_TAG_DENSE,
178  MF_ZERO);
179  // Modify field
180  auto field_ptr = mField.get_field_structure(field_name);
181  auto field_order_table =
182  const_cast<Field *>(field_ptr)->getFieldOrderTable();
183  auto get_cgg_bubble_order_zero = [](int p) { return 0; };
184  auto get_cgg_bubble_order_tet = [](int p) {
185  return NBVOLUMETET_CCG_BUBBLE(p);
186  };
187  field_order_table[MBVERTEX] = get_cgg_bubble_order_zero;
188  field_order_table[MBEDGE] = get_cgg_bubble_order_zero;
189  field_order_table[MBTRI] = get_cgg_bubble_order_zero;
190  field_order_table[MBTET] = get_cgg_bubble_order_tet;
191  CHKERR mField.add_ents_to_field_by_type(meshset, MBTET, field_name);
192  CHKERR mField.set_field_order(meshset, MBTRI, field_name, order);
193  CHKERR mField.set_field_order(meshset, MBTET, field_name, order);
195  };
196 
197  // spatial fields
198  CHKERR add_hdiv_field(piolaStress, spaceOrder, 3);
199  CHKERR add_bubble_field(bubbleField, spaceOrder, 1);
200  CHKERR add_l2_field(spatialDisp, spaceOrder - 1, 3);
201  CHKERR add_l2_field(rotAxis, spaceOrder - 1, 3);
202  CHKERR add_l2_field(streachTensor, spaceOrder, 6);
203 
204  // material fields
205  CHKERR add_hdiv_field(eshelbyStress, materialOrder, 3);
206  CHKERR add_l2_field(materialGradient, materialOrder - 1, 9);
207  // CHKERR add_l2_field(materialDisp, materialOrder - 1, 3);
208  // CHKERR add_l2_field(tauField, materialOrder - 1, 1);
209  // CHKERR add_l2_field(lambdaField, materialOrder - 1, 1);
210 
211  // Add history filedes
212  CHKERR add_l2_field(materialGradient + "0", materialOrder - 1, 9);
213 
215 
217 }
field with continuous normal traction
Definition: definitions.h:178
user implemented approximation base
Definition: definitions.h:159
virtual moab::Interface & get_moab()=0
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
virtual MoFEMErrorCode add_field(const std::string &name, const FieldSpace space, const FieldApproximationBase base, const FieldCoefficientsNumber nb_of_coefficients, const TagType tag_type=MB_TAG_SPARSE, const enum MoFEMTypes bh=MF_EXCL, int verb=DEFAULT_VERBOSITY)=0
Add field.
virtual MoFEMErrorCode set_field_order(const EntityHandle meshset, const EntityType type, const std::string &name, const ApproximationOrder order, int verb=DEFAULT_VERBOSITY)=0
Set order approximation of the entities in the field.
constexpr int order
virtual MoFEMErrorCode add_ents_to_field_by_type(const Range &ents, const EntityType type, const std::string &name, int verb=DEFAULT_VERBOSITY)=0
Add entities to field meshset.
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:151
FieldOrderTable & getFieldOrderTable()
Get the Field Order Table.
#define _IT_CUBITMESHSETS_BY_SET_TYPE_FOR_LOOP_(MESHSET_MANAGER, CUBITBCTYPE, IT)
Iterator that loops over a specific Cubit MeshSet having a particular BC meshset in a moFEM field.
#define CHKERR
Inline error check.
Definition: definitions.h:601
#define MYPCOMM_INDEX
default communicator number PCOMM
Definition: definitions.h:290
virtual MoFEMErrorCode build_fields(int verb=DEFAULT_VERBOSITY)=0
#define NBVOLUMETET_CCG_BUBBLE(P)
Bubble function for CGG H div space.
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412
continuous field
Definition: definitions.h:176
virtual const Field * get_field_structure(const std::string &name)=0
get field structure
field with C-1 continuity
Definition: definitions.h:179

◆ addMaterial_HMHHStVenantKirchhoff()

MoFEMErrorCode EshelbianPlasticity::EshelbianCore::addMaterial_HMHHStVenantKirchhoff ( const int  tape,
const double  lambda,
const double  mu,
const double  sigma_y 
)

Definition at line 644 of file EshelbianADOL-C.cpp.

646  {
648  physicalEquations = boost::shared_ptr<HMHStVenantKirchhoff>(
649  new HMHStVenantKirchhoff(lambda, mu, sigma_y));
650  CHKERR physicalEquations->recordTape(tape, nullptr);
652 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
boost::shared_ptr< PhysicalEquations > physicalEquations
#define CHKERR
Inline error check.
Definition: definitions.h:601
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412

◆ addMaterial_HMHMooneyRivlin()

MoFEMErrorCode EshelbianPlasticity::EshelbianCore::addMaterial_HMHMooneyRivlin ( const int  tape,
const double  alpha,
const double  beta,
const double  lambda,
const double  sigma_y 
)
Examples
ep.cpp.

Definition at line 654 of file EshelbianADOL-C.cpp.

656  {
658  physicalEquations = boost::shared_ptr<HMHPMooneyRivlinWriggersEq63>(
659  new HMHPMooneyRivlinWriggersEq63(alpha, beta, lambda, sigma_y));
660  CHKERR physicalEquations->recordTape(tape, nullptr);
662 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
boost::shared_ptr< PhysicalEquations > physicalEquations
#define CHKERR
Inline error check.
Definition: definitions.h:601
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412

◆ addVolumeFiniteElement()

MoFEMErrorCode EshelbianPlasticity::EshelbianCore::addVolumeFiniteElement ( const EntityHandle  meshset = 0)
Examples
ep.cpp.

Definition at line 220 of file EshelbianPlasticity.cpp.

220  {
222 
223  // set finite element fields
224  auto add_field_to_fe = [this](const std::string fe,
225  const std::string field_name) {
231  };
232 
237 
238  CHKERR add_field_to_fe(elementVolumeName, piolaStress);
239  CHKERR add_field_to_fe(elementVolumeName, bubbleField);
240  CHKERR add_field_to_fe(elementVolumeName, eshelbyStress);
241  CHKERR add_field_to_fe(elementVolumeName, streachTensor);
242  CHKERR add_field_to_fe(elementVolumeName, rotAxis);
243  CHKERR add_field_to_fe(elementVolumeName, spatialDisp);
244  CHKERR add_field_to_fe(elementVolumeName, streachTensor);
245  CHKERR add_field_to_fe(elementVolumeName, materialGradient);
246  // CHKERR add_field_to_fe(elementVolumeName, materialDisp);
247  // CHKERR add_field_to_fe(elementVolumeName, tauField);
248  // CHKERR add_field_to_fe(elementVolumeName, lambdaField);
249 
251  materialGradient + "0");
252  }
253 
254  // build finite elements data structures
256 
258 }
virtual MoFEMErrorCode modify_finite_element_add_field_col(const std::string &fe_name, const std::string &name_row)=0
set field col which finite element use
virtual MoFEMErrorCode add_ents_to_finite_element_by_type(const EntityHandle entities, const EntityType type, const std::string &name, const bool recursive=true)=0
add entities to finite element
virtual bool check_finite_element(const std::string &name) const =0
Check if finite element is in database.
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
virtual MoFEMErrorCode add_finite_element(const std::string &fe_name, enum MoFEMTypes bh=MF_EXCL, int verb=DEFAULT_VERBOSITY)=0
add finite element
#define CHKERR
Inline error check.
Definition: definitions.h:601
virtual MoFEMErrorCode modify_finite_element_add_field_data(const std::string &fe_name, const std::string &name_filed)=0
set finite element field data
virtual MoFEMErrorCode modify_finite_element_add_field_row(const std::string &fe_name, const std::string &name_row)=0
set field row which finite element use
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412
virtual MoFEMErrorCode build_finite_elements(int verb=DEFAULT_VERBOSITY)=0
Build finite elementsBuild finite element data structures. Have to be run before problem and adjacenc...

◆ getBc()

template<typename BC >
MoFEMErrorCode EshelbianPlasticity::EshelbianCore::getBc ( boost::shared_ptr< BC > &  bc_vec_ptr,
const std::string  block_set_name,
const int  nb_attributes 
)

Definition at line 1241 of file EshelbianPlasticity.hpp.

1243  {
1246  auto block_name = it->getName();
1247  if (block_name.compare(0, block_set_name.length(), block_set_name) == 0) {
1248  std::vector<double> block_attributes;
1249  CHKERR it->getAttributes(block_attributes);
1250  if (block_attributes.size() != nb_attributes) {
1251  SETERRQ3(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1252  "In block %s expected %d attributes, but given %d",
1253  it->getName().c_str(), nb_attributes,
1254  block_attributes.size());
1255  }
1256  Range faces;
1257  CHKERR it->getMeshsetIdEntitiesByDimension(mField.get_moab(), 2, faces,
1258  true);
1259  bc_vec_ptr->emplace_back(block_name, block_attributes, faces);
1260  }
1261  }
1263  }
virtual moab::Interface & get_moab()=0
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
#define _IT_CUBITMESHSETS_BY_SET_TYPE_FOR_LOOP_(MESHSET_MANAGER, CUBITBCTYPE, IT)
Iterator that loops over a specific Cubit MeshSet having a particular BC meshset in a moFEM field.
#define CHKERR
Inline error check.
Definition: definitions.h:601
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412

◆ getOptions()

MoFEMErrorCode EshelbianPlasticity::EshelbianCore::getOptions ( )

Definition at line 60 of file EshelbianPlasticity.cpp.

60  {
62  CHKERR PetscOptionsBegin(PETSC_COMM_WORLD, "", "Eshelbian plasticity",
63  "none");
64 
65  spaceOrder = 1;
66  CHKERR PetscOptionsInt("-space_order", "approximation oder for space", "",
67  spaceOrder, &spaceOrder, PETSC_NULL);
68  materialOrder = 1;
69  CHKERR PetscOptionsInt("-material_order", "approximation oder for material",
70  "", materialOrder, &materialOrder, PETSC_NULL);
71 
72  alpha_u = 0;
73  CHKERR PetscOptionsScalar("-viscosity_alpha_u", "viscosity", "", alpha_u,
74  &alpha_u, PETSC_NULL);
75 
76  alpha_w = 0;
77  CHKERR PetscOptionsScalar("-viscosity_alpha_w", "viscosity", "", alpha_w,
78  &alpha_w, PETSC_NULL);
79 
80  preconditioner_eps = 1e-3;
81  CHKERR PetscOptionsScalar("-preconditioner_eps", "preconditioner_eps", "",
83  PETSC_NULL);
84 
85  ierr = PetscOptionsEnd();
86  CHKERRG(ierr);
88 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:549
#define CHKERR
Inline error check.
Definition: definitions.h:601
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412

◆ getSpatialDispBc()

MoFEMErrorCode EshelbianPlasticity::EshelbianCore::getSpatialDispBc ( )
Examples
ep.cpp.

Definition at line 1265 of file EshelbianPlasticity.hpp.

1265  {
1266  bcSpatialDispVecPtr = boost::make_shared<BcDispVec>();
1267  return getBc(bcSpatialDispVecPtr, "SPATIAL_DISP_BC", 6);
1268  }
boost::shared_ptr< BcDispVec > bcSpatialDispVecPtr
MoFEMErrorCode getBc(boost::shared_ptr< BC > &bc_vec_ptr, const std::string block_set_name, const int nb_attributes)

◆ getSpatialRotationBc()

MoFEMErrorCode EshelbianPlasticity::EshelbianCore::getSpatialRotationBc ( )
Examples
ep.cpp.

Definition at line 1270 of file EshelbianPlasticity.hpp.

1270  {
1271  bcSpatialRotationVecPtr = boost::make_shared<BcRotVec>();
1272  return getBc(bcSpatialRotationVecPtr, "SPATIAL_ROTATION_BC", 4);
1273  }
boost::shared_ptr< BcRotVec > bcSpatialRotationVecPtr
MoFEMErrorCode getBc(boost::shared_ptr< BC > &bc_vec_ptr, const std::string block_set_name, const int nb_attributes)

◆ getSpatialTractionBc()

MoFEMErrorCode EshelbianPlasticity::EshelbianCore::getSpatialTractionBc ( )
Examples
ep.cpp.

Definition at line 1275 of file EshelbianPlasticity.hpp.

1275  {
1276  bcSpatialTraction = boost::make_shared<TractionBcVec>();
1277  return getBc(bcSpatialTraction, "SPATIAL_TRACTION_BC", 6);
1278  }
MoFEMErrorCode getBc(boost::shared_ptr< BC > &bc_vec_ptr, const std::string block_set_name, const int nb_attributes)
boost::shared_ptr< TractionBcVec > bcSpatialTraction

◆ getSpatialTractionFreeBc()

MoFEMErrorCode EshelbianPlasticity::EshelbianCore::getSpatialTractionFreeBc ( const EntityHandle  meshset = 0)
Examples
ep.cpp.

Definition at line 1285 of file EshelbianPlasticity.hpp.

1285  {
1287  boost::shared_ptr<TractionFreeBc>(new TractionFreeBc());
1288  return getTractionFreeBc(meshset, bcSpatialFreeTraction, "SPATIAL_DISP_BC",
1289  "SPATIAL_ROTATION_BC");
1290  }
std::vector< Range > TractionFreeBc
boost::shared_ptr< TractionFreeBc > bcSpatialFreeTraction
MoFEMErrorCode getTractionFreeBc(const EntityHandle meshset, boost::shared_ptr< TractionFreeBc > &bc_ptr, const std::string disp_block_set_name, const std::string rot_block_set_name)

◆ getTractionFreeBc()

MoFEMErrorCode EshelbianPlasticity::EshelbianCore::getTractionFreeBc ( const EntityHandle  meshset,
boost::shared_ptr< TractionFreeBc > &  bc_ptr,
const std::string  disp_block_set_name,
const std::string  rot_block_set_name 
)

Definition at line 461 of file EshelbianPlasticity.cpp.

464  {
466  Range tets;
467  CHKERR mField.get_moab().get_entities_by_type(meshset, MBTET, tets);
468  Range tets_skin_part;
469  Skinner skin(&mField.get_moab());
470  CHKERR skin.find_skin(0, tets, false, tets_skin_part);
471  ParallelComm *pcomm =
472  ParallelComm::get_pcomm(&mField.get_moab(), MYPCOMM_INDEX);
473  Range tets_skin;
474  CHKERR pcomm->filter_pstatus(tets_skin_part,
475  PSTATUS_SHARED | PSTATUS_MULTISHARED,
476  PSTATUS_NOT, -1, &tets_skin);
477 
478  bc_ptr->resize(3);
479  for (int dd = 0; dd != 3; ++dd)
480  (*bc_ptr)[dd] = tets_skin;
481 
483  if (it->getName().compare(0, disp_block_set_name.length(),
484  disp_block_set_name) == 0) {
485  std::vector<double> block_attributes;
486  CHKERR it->getAttributes(block_attributes);
487  if (block_attributes.size() != 6) {
488  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
489  "In block %s six attributes are required for given BC "
490  "blockset (3 values + "
491  "3 flags) != %d",
492  it->getName().c_str(), block_attributes.size());
493  }
494  Range faces;
495  CHKERR it->getMeshsetIdEntitiesByDimension(mField.get_moab(), 2, faces,
496  true);
497  if (block_attributes[3] != 0)
498  (*bc_ptr)[0] = subtract((*bc_ptr)[0], faces);
499  if (block_attributes[4] != 0)
500  (*bc_ptr)[1] = subtract((*bc_ptr)[1], faces);
501  if (block_attributes[5] != 0)
502  (*bc_ptr)[2] = subtract((*bc_ptr)[2], faces);
503  }
504  if (it->getName().compare(0, rot_block_set_name.length(),
505  rot_block_set_name) == 0) {
506  Range faces;
507  CHKERR it->getMeshsetIdEntitiesByDimension(mField.get_moab(), 2, faces,
508  true);
509  (*bc_ptr)[0] = subtract((*bc_ptr)[0], faces);
510  (*bc_ptr)[1] = subtract((*bc_ptr)[1], faces);
511  (*bc_ptr)[2] = subtract((*bc_ptr)[2], faces);
512  }
513  }
514 
515  // for (int dd = 0; dd != 3; ++dd) {
516  // EntityHandle meshset;
517  // CHKERR mField.get_moab().create_meshset(MESHSET_SET, meshset);
518  // CHKERR mField.get_moab().add_entities(meshset, (*bc_ptr)[dd]);
519  // std::string file_name = disp_block_set_name +
520  // "_traction_free_bc_" + boost::lexical_cast<std::string>(dd) + ".vtk";
521  // CHKERR mField.get_moab().write_file(file_name.c_str(), " VTK ", "",
522  // &meshset, 1);
523  // CHKERR mField.get_moab().delete_entities(&meshset, 1);
524  // }
525 
527 }
virtual moab::Interface & get_moab()=0
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
const Tensor2_symmetric_Expr< const ddTensor0< T, Dim, i, j >, typename promote< T, double >::V, Dim, i, j > dd(const Tensor0< T * > &a, const Index< i, Dim > index1, const Index< j, Dim > index2, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition: ddTensor0.hpp:33
#define _IT_CUBITMESHSETS_BY_SET_TYPE_FOR_LOOP_(MESHSET_MANAGER, CUBITBCTYPE, IT)
Iterator that loops over a specific Cubit MeshSet having a particular BC meshset in a moFEM field.
#define CHKERR
Inline error check.
Definition: definitions.h:601
#define MYPCOMM_INDEX
default communicator number PCOMM
Definition: definitions.h:290
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412

◆ postProcessResults()

MoFEMErrorCode EshelbianPlasticity::EshelbianCore::postProcessResults ( const int  tag,
const std::string  file 
)

Definition at line 1443 of file EshelbianPlasticity.cpp.

1444  {
1446 
1447  if (!dataAtPts) {
1448  dataAtPts =
1449  boost::shared_ptr<DataAtIntegrationPts>(new DataAtIntegrationPts());
1450  dataAtPts->approxPAtPts = boost::make_shared<MatrixDouble>();
1451  dataAtPts->approxSigmaAtPts = boost::make_shared<MatrixDouble>();
1452  dataAtPts->divSigmaAtPts = boost::make_shared<MatrixDouble>();
1453  dataAtPts->wAtPts = boost::make_shared<MatrixDouble>();
1454  dataAtPts->rotAxisAtPts = boost::make_shared<MatrixDouble>();
1455  dataAtPts->logStreachTensorAtPts = boost::make_shared<MatrixDouble>();
1456  dataAtPts->streachTensorAtPts = boost::make_shared<MatrixDouble>();
1457  dataAtPts->diffStreachTensorAtPts = boost::make_shared<MatrixDouble>();
1458  dataAtPts->rotMatAtPts = boost::make_shared<MatrixDouble>();
1459  dataAtPts->diffRotMatAtPts = boost::make_shared<MatrixDouble>();
1460  dataAtPts->hAtPts = boost::make_shared<MatrixDouble>();
1461  dataAtPts->GAtPts = boost::make_shared<MatrixDouble>();
1462  dataAtPts->G0AtPts = boost::make_shared<MatrixDouble>();
1463  }
1464 
1466  post_proc.getUserPolynomialBase() =
1467  boost::shared_ptr<BaseFunction>(new CGGUserPolynomialBase());
1468 
1469  CHKERR post_proc.generateReferenceElementMesh();
1470  post_proc.getOpPtrVector().push_back(new OpCalculateHVecTensorField<3, 3>(
1471  piolaStress, dataAtPts->approxPAtPts));
1472  post_proc.getOpPtrVector().push_back(new OpCalculateHTensorTensorField<3, 3>(
1473  bubbleField, dataAtPts->approxPAtPts, MBMAXTYPE));
1474  post_proc.getOpPtrVector().push_back(
1476  streachTensor, dataAtPts->logStreachTensorAtPts, MBTET));
1477  post_proc.getOpPtrVector().push_back(new OpCalculateVectorFieldValues<3>(
1478  rotAxis, dataAtPts->rotAxisAtPts, MBTET));
1479  post_proc.getOpPtrVector().push_back(new OpCalculateTensor2FieldValues<3, 3>(
1480  materialGradient, dataAtPts->GAtPts, MBTET));
1481  post_proc.getOpPtrVector().push_back(new OpCalculateVectorFieldValues<3>(
1482  spatialDisp, dataAtPts->wAtPts, MBTET));
1483 
1484  // evaluate derived quantities
1485  post_proc.getOpPtrVector().push_back(
1486  new OpCalculateRotationAndSpatialGradient(rotAxis, dataAtPts));
1487 
1488  // evaluate integration points
1489  post_proc.getOpPtrVector().push_back(physicalEquations->returnOpJacobian(
1490  spatialDisp, tag, true, false, dataAtPts, physicalEquations));
1491  post_proc.getOpPtrVector().push_back(new OpPostProcDataStructure(
1492  post_proc.postProcMesh, post_proc.mapGaussPts, spatialDisp, dataAtPts));
1493 
1494  CHKERR DMoFEMLoopFiniteElements(dM, elementVolumeName.c_str(), &post_proc);
1495  CHKERR post_proc.writeFile(file.c_str());
1497 }
Calculate tenor field using vectorial base, i.e. Hdiv/Hcurl.
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
Get values at integration pts for tensor filed rank 1, i.e. vector field.
PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[], MoFEM::FEMethod *method)
Executes FEMethod for finite elements in DM.
Definition: DMMMoFEM.cpp:493
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
Post processing.
boost::shared_ptr< PhysicalEquations > physicalEquations
Calculate tenor field using vectorial base, i.e. Hdiv/Hcurl.
#define CHKERR
Inline error check.
Definition: definitions.h:601
SmartPetscObj< DM > dM
Coupled problem all fields.
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412
Calculate symmetric tensor field values at integration pts.
Get values at integration pts for tensor filed rank 2, i.e. matrix field.

◆ query_interface()

MoFEMErrorCode EshelbianPlasticity::EshelbianCore::query_interface ( const MOFEMuuid uuid,
UnknownInterface **  iface 
) const

Getting interface of core database.

Parameters
uuidunique ID of interface
ifacereturned pointer to interface
Returns
error code

Definition at line 19 of file EshelbianPlasticity.cpp.

20  {
22  *iface = NULL;
24  *iface = const_cast<EshelbianCore *>(this);
26  }
27  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "unknown interface");
29 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:506
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:513
static const MOFEMuuid IDD_MOFEMEshelbianCrackInterface

◆ setBaseVolumeElementOps()

MoFEMErrorCode EshelbianPlasticity::EshelbianCore::setBaseVolumeElementOps ( const int  tag,
const bool  do_rhs,
const bool  do_lhs,
boost::shared_ptr< EpElement< VolumeElementForcesAndSourcesCore >> &  fe 
)

Definition at line 647 of file EshelbianPlasticity.cpp.

649  {
651  fe = boost::make_shared<EpElement<VolumeElementForcesAndSourcesCore>>(mField);
652 
653  fe->getUserPolynomialBase() =
654  boost::shared_ptr<BaseFunction>(new CGGUserPolynomialBase());
655 
656  // set integration rule
657  fe->getRuleHook = VolRule();
658 
659  if (!dataAtPts) {
660  dataAtPts =
661  boost::shared_ptr<DataAtIntegrationPts>(new DataAtIntegrationPts());
662  dataAtPts->approxPAtPts = boost::make_shared<MatrixDouble>();
663  dataAtPts->divPAtPts = boost::make_shared<MatrixDouble>();
664  dataAtPts->approxSigmaAtPts = boost::make_shared<MatrixDouble>();
665  dataAtPts->divSigmaAtPts = boost::make_shared<MatrixDouble>();
666  dataAtPts->wAtPts = boost::make_shared<MatrixDouble>();
667  dataAtPts->wDotAtPts = boost::make_shared<MatrixDouble>();
668  dataAtPts->rotAxisAtPts = boost::make_shared<MatrixDouble>();
669  dataAtPts->logStreachTensorAtPts = boost::make_shared<MatrixDouble>();
670  dataAtPts->streachTensorAtPts = boost::make_shared<MatrixDouble>();
671  dataAtPts->diffStreachTensorAtPts = boost::make_shared<MatrixDouble>();
672  dataAtPts->logStreachDotTensorAtPts = boost::make_shared<MatrixDouble>();
673  dataAtPts->rotMatAtPts = boost::make_shared<MatrixDouble>();
674  dataAtPts->diffRotMatAtPts = boost::make_shared<MatrixDouble>();
675  dataAtPts->hAtPts = boost::make_shared<MatrixDouble>();
676  dataAtPts->GAtPts = boost::make_shared<MatrixDouble>();
677  dataAtPts->G0AtPts = boost::make_shared<MatrixDouble>();
678  dataAtPts->physicsPtr = physicalEquations;
679  }
680 
681  // calculate fields values
682  fe->getOpPtrVector().push_back(new OpCalculateHVecTensorField<3, 3>(
683  piolaStress, dataAtPts->approxPAtPts));
684  fe->getOpPtrVector().push_back(new OpCalculateHTensorTensorField<3, 3>(
685  bubbleField, dataAtPts->approxPAtPts, MBMAXTYPE));
686  fe->getOpPtrVector().push_back(new OpCalculateHVecTensorDivergence<3, 3>(
687  piolaStress, dataAtPts->divPAtPts));
688  fe->getOpPtrVector().push_back(new OpCalculateHVecTensorField<3, 3>(
689  eshelbyStress, dataAtPts->approxSigmaAtPts));
690  fe->getOpPtrVector().push_back(new OpCalculateHVecTensorDivergence<3, 3>(
691  eshelbyStress, dataAtPts->divSigmaAtPts));
692  fe->getOpPtrVector().push_back(new OpCalculateTensor2SymmetricFieldValues<3>(
693  streachTensor, dataAtPts->logStreachTensorAtPts, MBTET));
694  fe->getOpPtrVector().push_back(new OpCalculateVectorFieldValues<3>(
695  rotAxis, dataAtPts->rotAxisAtPts, MBTET));
696  fe->getOpPtrVector().push_back(new OpCalculateTensor2FieldValues<3, 3>(
697  materialGradient, dataAtPts->GAtPts, MBTET));
698  fe->getOpPtrVector().push_back(new OpCalculateTensor2FieldValues<3, 3>(
699  materialGradient + "0", dataAtPts->G0AtPts, MBTET));
700  fe->getOpPtrVector().push_back(new OpCalculateVectorFieldValues<3>(
701  spatialDisp, dataAtPts->wAtPts, MBTET));
702 
703  // velocities
704  fe->getOpPtrVector().push_back(new OpCalculateVectorFieldValuesDot<3>(
705  spatialDisp, dataAtPts->wDotAtPts, MBTET));
706  fe->getOpPtrVector().push_back(
708  streachTensor, dataAtPts->logStreachDotTensorAtPts, MBTET));
709 
710  // calculate other derived quantities
711  fe->getOpPtrVector().push_back(
712  new OpCalculateRotationAndSpatialGradient(rotAxis, dataAtPts));
713 
714  // evaluate integration points
715  fe->getOpPtrVector().push_back(physicalEquations->returnOpJacobian(
716  spatialDisp, tag, do_rhs, do_lhs, dataAtPts, physicalEquations));
717 
719 }
Set integration rule.
Get time direvatives of values at integration pts for tensor filed rank 1, i.e. vector field.
Calculate tenor field using vectorial base, i.e. Hdiv/Hcurl.
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
Get values at integration pts for tensor filed rank 1, i.e. vector field.
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
boost::shared_ptr< PhysicalEquations > physicalEquations
Calculate divergence of tonsorial field using vectorial base.
Calculate tenor field using vectorial base, i.e. Hdiv/Hcurl.
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412
Calculate symmetric tensor field rates ant integratio pts.
Calculate symmetric tensor field values at integration pts.
Get values at integration pts for tensor filed rank 2, i.e. matrix field.

◆ setElasticElementOps()

MoFEMErrorCode EshelbianPlasticity::EshelbianCore::setElasticElementOps ( const int  tag)
Examples
ep.cpp.

Definition at line 831 of file EshelbianPlasticity.cpp.

831  {
834  elasticFeLhs);
837 }
MoFEMErrorCode setGenericFaceElementOps(const bool add_elastic, const bool add_material, boost::shared_ptr< EpElement< FaceElementForcesAndSourcesCore >> &fe_rhs, boost::shared_ptr< EpElement< FaceElementForcesAndSourcesCore >> &fe_lhs)
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
MoFEMErrorCode setGenericVolumeElementOps(const int tag, const bool add_elastic, const bool add_material, boost::shared_ptr< EpElement< VolumeElementForcesAndSourcesCore >> &fe_rhs, boost::shared_ptr< EpElement< VolumeElementForcesAndSourcesCore >> &fe_lhs)
boost::shared_ptr< EpElement< FaceElementForcesAndSourcesCore > > elasticBcRhs
#define CHKERR
Inline error check.
Definition: definitions.h:601
boost::shared_ptr< EpElement< VolumeElementForcesAndSourcesCore > > elasticFeLhs
boost::shared_ptr< EpElement< FaceElementForcesAndSourcesCore > > elasticBcLhs
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412
boost::shared_ptr< EpElement< VolumeElementForcesAndSourcesCore > > elasticFeRhs

◆ setElasticElementToTs()

MoFEMErrorCode EshelbianPlasticity::EshelbianCore::setElasticElementToTs ( DM  dm)
Examples
ep.cpp.

Definition at line 839 of file EshelbianPlasticity.cpp.

839  {
841  boost::shared_ptr<FEMethod> null;
842  boost::shared_ptr<EpElement<FeTractionBc>> spatial_traction_bc(
843  new EpElement<FeTractionBc>(mField, piolaStress, bcSpatialTraction));
844 
845  CHKERR DMMoFEMTSSetIFunction(dm, DM_NO_ELEMENT, null, spatial_traction_bc,
846  null);
850  spatial_traction_bc);
851 
852  schurAssembly = boost::make_shared<EpFEMethod>();
858 }
#define DM_NO_ELEMENT
Definition: DMMoFEM.hpp:22
boost::shared_ptr< EpFEMethod > schurAssembly
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
PetscErrorCode DMMoFEMTSSetIJacobian(DM dm, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
set TS Jacobian evaluation function
Definition: DMMMoFEM.cpp:758
boost::shared_ptr< EpElement< FaceElementForcesAndSourcesCore > > elasticBcRhs
#define CHKERR
Inline error check.
Definition: definitions.h:601
boost::shared_ptr< EpElement< VolumeElementForcesAndSourcesCore > > elasticFeLhs
boost::shared_ptr< EpElement< FaceElementForcesAndSourcesCore > > elasticBcLhs
boost::shared_ptr< TractionBcVec > bcSpatialTraction
PetscErrorCode DMMoFEMTSSetIFunction(DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
set TS implicit function evaluation function
Definition: DMMMoFEM.cpp:705
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412
boost::shared_ptr< EpElement< VolumeElementForcesAndSourcesCore > > elasticFeRhs

◆ setGenericFaceElementOps()

MoFEMErrorCode EshelbianPlasticity::EshelbianCore::setGenericFaceElementOps ( const bool  add_elastic,
const bool  add_material,
boost::shared_ptr< EpElement< FaceElementForcesAndSourcesCore >> &  fe_rhs,
boost::shared_ptr< EpElement< FaceElementForcesAndSourcesCore >> &  fe_lhs 
)

Definition at line 805 of file EshelbianPlasticity.cpp.

808  {
810 
811  fe_rhs =
812  boost::make_shared<EpElement<FaceElementForcesAndSourcesCore>>(mField);
813  fe_lhs =
814  boost::make_shared<EpElement<FaceElementForcesAndSourcesCore>>(mField);
815 
816  // set integration rule
817  fe_rhs->getRuleHook = FaceRule();
818  fe_lhs->getRuleHook = FaceRule();
819 
820  if (add_elastic) {
821  fe_rhs->getOpPtrVector().push_back(
822  new OpDispBc(piolaStress, dataAtPts, bcSpatialDispVecPtr));
823  fe_rhs->getOpPtrVector().push_back(
824  new OpRotationBc(piolaStress, dataAtPts, bcSpatialRotationVecPtr));
825 
826  }
827 
829 }
boost::shared_ptr< BcRotVec > bcSpatialRotationVecPtr
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
boost::shared_ptr< BcDispVec > bcSpatialDispVecPtr
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
Set integration rule to boundary elements.
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412

◆ setGenericVolumeElementOps()

MoFEMErrorCode EshelbianPlasticity::EshelbianCore::setGenericVolumeElementOps ( const int  tag,
const bool  add_elastic,
const bool  add_material,
boost::shared_ptr< EpElement< VolumeElementForcesAndSourcesCore >> &  fe_rhs,
boost::shared_ptr< EpElement< VolumeElementForcesAndSourcesCore >> &  fe_lhs 
)

Definition at line 721 of file EshelbianPlasticity.cpp.

724  {
726 
727  // Right hand side
728  CHKERR setBaseVolumeElementOps(tag, true, false, fe_rhs);
729 
730  // elastic
731  if (add_elastic) {
732  fe_rhs->getOpPtrVector().push_back(
733  new OpSpatialEquilibrium(spatialDisp, dataAtPts, alpha_w));
734  fe_rhs->getOpPtrVector().push_back(
735  new OpSpatialRotation(rotAxis, dataAtPts));
736  fe_rhs->getOpPtrVector().push_back(
737  new OpSpatialPhysical(streachTensor, dataAtPts, alpha_u));
738  fe_rhs->getOpPtrVector().push_back(
739  new OpSpatialConsistencyP(piolaStress, dataAtPts));
740  fe_rhs->getOpPtrVector().push_back(
741  new OpSpatialConsistencyBubble(bubbleField, dataAtPts));
742  fe_rhs->getOpPtrVector().push_back(
743  new OpSpatialConsistencyDivTerm(piolaStress, dataAtPts));
744  }
745 
746  // Left hand side
747  CHKERR setBaseVolumeElementOps(tag, true, true, fe_lhs);
748 
749  // elastic
750  if (add_elastic) {
751 
752  // Schur
753  fe_lhs->getOpPtrVector().push_back(
754  new OpSpatialSchurBegin(spatialDisp, dataAtPts));
755 
756  fe_lhs->getOpPtrVector().push_back(new OpSpatialPhysical_du_du(
758  fe_lhs->getOpPtrVector().push_back(new OpSpatialPhysical_du_dBubble(
760 
761  fe_lhs->getOpPtrVector().push_back(new OpSpatialEquilibrium_dw_dP(
763  fe_lhs->getOpPtrVector().push_back(new OpSpatialEquilibrium_dw_dw(
765  fe_lhs->getOpPtrVector().push_back(
766  new OpSpatialConsistency_dP_domega(piolaStress, rotAxis, dataAtPts));
767  fe_lhs->getOpPtrVector().push_back(new OpSpatialConsistency_dBubble_domega(
768  bubbleField, rotAxis, dataAtPts, false));
769 
770  fe_lhs->getOpPtrVector().push_back(new OpSpatialPhysical_du_dP(
772  fe_lhs->getOpPtrVector().push_back(new OpSpatialPhysical_du_domega(
773  streachTensor, rotAxis, dataAtPts, false));
774 
775  fe_lhs->getOpPtrVector().push_back(new OpSpatialRotation_domega_dP(
776  rotAxis, piolaStress, dataAtPts, false));
777  fe_lhs->getOpPtrVector().push_back(new OpSpatialRotation_domega_dBubble(
778  rotAxis, bubbleField, dataAtPts, false));
779  fe_lhs->getOpPtrVector().push_back(
780  new OpSpatialRotation_domega_domega(rotAxis, rotAxis, dataAtPts));
781  // fe_lhs->getOpPtrVector().push_back(
782  // new OpSpatialRotation_domega_du(rotAxis, streachTensor, dataAtPts));
783 
784  // Schur
785  dataAtPts->ooMatPtr = boost::make_shared<MatrixDouble>();
786  fe_lhs->getOpPtrVector().push_back(
787  new OpSpatialPreconditionMass(rotAxis, dataAtPts->ooMatPtr));
788  if (alpha_w < std::numeric_limits<double>::epsilon()) {
789  dataAtPts->wwMatPtr = boost::make_shared<MatrixDouble>();
790  fe_lhs->getOpPtrVector().push_back(
791  new OpSpatialPreconditionMass(spatialDisp, dataAtPts->wwMatPtr));
792  } else {
793  dataAtPts->wwMatPtr.reset();
794  }
795  fe_lhs->getOpPtrVector().push_back(
796  new OpSpatialSchurEnd(spatialDisp, dataAtPts, preconditioner_eps));
797 
798  if (add_material) {
799  }
800  }
801 
803 }
MoFEMErrorCode setBaseVolumeElementOps(const int tag, const bool do_rhs, const bool do_lhs, boost::shared_ptr< EpElement< VolumeElementForcesAndSourcesCore >> &fe)
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
#define CHKERR
Inline error check.
Definition: definitions.h:601
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412

◆ setUpTSElastic()

MoFEMErrorCode EshelbianPlasticity::EshelbianCore::setUpTSElastic ( TS  ts,
Mat  m,
Vec  f 
)
Examples
ep.cpp.

Definition at line 860 of file EshelbianPlasticity.cpp.

860  {
862  boost::shared_ptr<TsCtx> ts_ctx;
864  CHKERR TSSetIFunction(ts, f, PETSC_NULL, PETSC_NULL);
865  CHKERR TSSetIJacobian(ts, m, m, PETSC_NULL, PETSC_NULL);
866  CHKERR TSMonitorSet(ts, TsMonitorSet, ts_ctx.get(), PETSC_NULL);
867 
868  auto add_schur_streach_op = [this](auto &list, Mat S, AO aoS) {
870  for (auto &fe : list)
871  if (auto fe_cast = dynamic_cast<EpElementBase *>(fe.second))
872  fe_cast->addStreachSchurMatrix(S, aoS);
873  else
874  SETERRQ(mField.get_comm(), MOFEM_DATA_INCONSISTENCY, "No EpElement");
876  };
877 
878  auto add_schur_streach_pre = [this](auto &list, Mat S, AO aoS) {
880  for (auto &fe : list)
881  if (auto fe_cast = dynamic_cast<EpElementBase *>(fe.getSharedPtr().get()))
882  fe_cast->addStreachSchurMatrix(S, aoS);
883  else
884  SETERRQ(mField.get_comm(), MOFEM_DATA_INCONSISTENCY, "No EpElement");
886  };
887 
888  auto add_schur_bubble_op = [this](auto &list, Mat S, AO aoS) {
890  for (auto &fe : list)
891  if (auto fe_cast = dynamic_cast<EpElementBase *>(fe.second))
892  fe_cast->addBubbleSchurMatrix(S, aoS);
893  else
894  SETERRQ(mField.get_comm(), MOFEM_DATA_INCONSISTENCY, "No EpElement");
896  };
897 
898  auto add_schur_bubble_pre = [this](auto &list, Mat S, AO aoS) {
900  for (auto &fe : list)
901  if (auto fe_cast = dynamic_cast<EpElementBase *>(fe.getSharedPtr().get()))
902  fe_cast->addBubbleSchurMatrix(S, aoS);
903  else
904  SETERRQ(mField.get_comm(), MOFEM_DATA_INCONSISTENCY, "No EpElement");
906  };
907 
908  auto add_schur_spatial_disp_op = [this](auto &list, Mat S, AO aoS) {
910  for (auto &fe : list)
911  if (auto fe_cast = dynamic_cast<EpElementBase *>(fe.second))
912  fe_cast->addSpatialDispStressSchurMatrix(S, aoS);
913  else
914  SETERRQ(mField.get_comm(), MOFEM_DATA_INCONSISTENCY, "No EpElement");
916  };
917 
918  auto add_schur_spatial_disp_pre = [this](auto &list, Mat S, AO aoS) {
920  for (auto &fe : list)
921  if (auto fe_cast = dynamic_cast<EpElementBase *>(fe.getSharedPtr().get()))
922  fe_cast->addSpatialDispStressSchurMatrix(S, aoS);
923  else
924  SETERRQ(mField.get_comm(), MOFEM_DATA_INCONSISTENCY, "No EpElement");
926  };
927 
928  auto add_schur_omega_op = [this](auto &list, Mat S, AO aoS) {
930  for (auto &fe : list)
931  if (auto fe_cast = dynamic_cast<EpElementBase *>(fe.second))
932  fe_cast->addOmegaSchurMatrix(S, aoS);
933  else
934  SETERRQ(mField.get_comm(), MOFEM_DATA_INCONSISTENCY, "No EpElement");
936  };
937 
938  auto add_schur_omega_pre = [this](auto &list, Mat S, AO ao) {
940  for (auto &fe : list)
941  if (auto fe_cast = dynamic_cast<EpElementBase *>(fe.getSharedPtr().get()))
942  fe_cast->addOmegaSchurMatrix(S, ao);
943  else
944  SETERRQ(mField.get_comm(), MOFEM_DATA_INCONSISTENCY, "No EpElement");
946  };
947 
948  const MoFEM::Problem *schur_streach_prb_ptr;
949  CHKERR DMMoFEMGetProblemPtr(dmElasticSchurStreach, &schur_streach_prb_ptr);
950  if (auto sub_data = schur_streach_prb_ptr->subProblemData) {
951  Mat Suu;
952  CHKERR DMCreateMatrix(dmElasticSchurStreach, &Suu);
953  AO aoSuu;
954  CHKERR sub_data->getRowMap(&aoSuu);
955 
956  CHKERR add_schur_streach_op(ts_ctx->loops_to_do_IJacobian, Suu, aoSuu);
957  CHKERR add_schur_streach_pre(ts_ctx->preProcess_IJacobian, Suu, aoSuu);
958  CHKERR add_schur_streach_pre(ts_ctx->postProcess_IJacobian, Suu, aoSuu);
959 
960  CHKERR MatDestroy(&Suu);
961  CHKERR AODestroy(&aoSuu);
962 
963  const MoFEM::Problem *schur_bubble_prb_ptr;
964  CHKERR DMMoFEMGetProblemPtr(dmElasticSchurBubble, &schur_bubble_prb_ptr);
965  if (auto bubble_data = schur_bubble_prb_ptr->subProblemData) {
966  Mat SBubble;
967  CHKERR DMCreateMatrix(dmElasticSchurBubble, &SBubble);
968  AO aoSBubble;
969  CHKERR bubble_data->getRowMap(&aoSBubble);
970 
971  CHKERR add_schur_bubble_op(ts_ctx->loops_to_do_IJacobian, SBubble,
972  aoSBubble);
973  CHKERR add_schur_bubble_pre(ts_ctx->preProcess_IJacobian, SBubble,
974  aoSBubble);
975  CHKERR add_schur_bubble_pre(ts_ctx->postProcess_IJacobian, SBubble,
976  aoSBubble);
977 
978  CHKERR MatDestroy(&SBubble);
979  CHKERR AODestroy(&aoSBubble);
980 
981  const MoFEM::Problem *schur_omega_prb_ptr;
982  CHKERR DMMoFEMGetProblemPtr(dmElasticSchurOmega, &schur_omega_prb_ptr);
983  if (auto tet_stress_data = schur_omega_prb_ptr->subProblemData) {
984  Mat SOmega;
985  CHKERR DMCreateMatrix(dmElasticSchurOmega, &SOmega);
986  AO aoSOmega;
987  CHKERR tet_stress_data->getRowMap(&aoSOmega);
988 
989  CHKERR add_schur_omega_op(ts_ctx->loops_to_do_IJacobian, SOmega,
990  aoSOmega);
991  CHKERR add_schur_omega_pre(ts_ctx->preProcess_IJacobian, SOmega,
992  aoSOmega);
993  CHKERR add_schur_omega_pre(ts_ctx->postProcess_IJacobian, SOmega,
994  aoSOmega);
995 
996  const MoFEM::Problem *schur_spatial_disp_prb_ptr;
998  &schur_spatial_disp_prb_ptr);
999  if (auto spatial_disp_data =
1000  schur_spatial_disp_prb_ptr->subProblemData) {
1001 
1002  Mat Sw;
1003  CHKERR DMCreateMatrix(dmElasticSchurSpatialDisp, &Sw);
1004  AO aoSw;
1005  CHKERR spatial_disp_data->getRowMap(&aoSw);
1006 
1007  CHKERR add_schur_spatial_disp_op(ts_ctx->loops_to_do_IJacobian, Sw,
1008  aoSw);
1009  CHKERR add_schur_spatial_disp_pre(ts_ctx->preProcess_IJacobian, Sw,
1010  aoSw);
1011  CHKERR add_schur_spatial_disp_pre(ts_ctx->postProcess_IJacobian, Sw,
1012  aoSw);
1013 
1014  CHKERR MatDestroy(&Sw);
1015  CHKERR AODestroy(&aoSw);
1016  } else
1018  "Problem does not have sub-problem data");
1019 
1020  CHKERR MatDestroy(&SOmega);
1021  CHKERR AODestroy(&aoSOmega);
1022 
1023  } else
1025  "Problem does not have sub-problem data");
1026 
1027  } else
1029  "Problem does not have sub-problem data");
1030 
1031  } else
1033  "Problem does not have sub-problem data");
1034 
1035  struct Monitor : public FEMethod {
1036 
1037  using Ele = ForcesAndSourcesCore;
1038  using VolEle = VolumeElementForcesAndSourcesCore;
1040  using SetPtsData = FieldEvaluatorInterface::SetPtsData;
1041 
1042  EshelbianCore &eP;
1043  boost::shared_ptr<SetPtsData> dataFieldEval;
1044  boost::shared_ptr<VolEle> volPostProcEle;
1045  boost::shared_ptr<double> gEnergy;
1046 
1047  Monitor(EshelbianCore &ep)
1048  : eP(ep),
1049  dataFieldEval(ep.mField.getInterface<FieldEvaluatorInterface>()
1050  ->getData<VolEle>()),
1051  volPostProcEle(new VolEle(ep.mField)), gEnergy(new double) {
1052  ierr = ep.mField.getInterface<FieldEvaluatorInterface>()->buildTree3D(
1053  dataFieldEval, "EP");
1054  CHKERRABORT(PETSC_COMM_WORLD, ierr);
1055 
1056  auto no_rule = [](int, int, int) { return -1; };
1057 
1058  auto set_element_for_field_eval = [&]() {
1059  boost::shared_ptr<Ele> vol_ele(dataFieldEval->feMethodPtr.lock());
1060  vol_ele->getRuleHook = no_rule;
1061  vol_ele->getUserPolynomialBase() =
1062  boost::shared_ptr<BaseFunction>(new CGGUserPolynomialBase());
1063  vol_ele->getOpPtrVector().push_back(
1064  new OpCalculateHVecTensorField<3, 3>(eP.piolaStress,
1065  eP.dataAtPts->approxPAtPts));
1066  vol_ele->getOpPtrVector().push_back(
1068  eP.bubbleField, eP.dataAtPts->approxPAtPts, MBMAXTYPE));
1069  vol_ele->getOpPtrVector().push_back(
1071  eP.streachTensor, eP.dataAtPts->logStreachTensorAtPts, MBTET));
1072  vol_ele->getOpPtrVector().push_back(new OpCalculateVectorFieldValues<3>(
1073  eP.rotAxis, eP.dataAtPts->rotAxisAtPts, MBTET));
1074  vol_ele->getOpPtrVector().push_back(
1076  eP.materialGradient, eP.dataAtPts->GAtPts, MBTET));
1077  vol_ele->getOpPtrVector().push_back(new OpCalculateVectorFieldValues<3>(
1078  eP.spatialDisp, eP.dataAtPts->wAtPts, MBTET));
1079  vol_ele->getOpPtrVector().push_back(
1080  new OpCalculateRotationAndSpatialGradient(eP.rotAxis,
1081  eP.dataAtPts));
1082  };
1083 
1084 
1085  auto set_element_for_post_process = [&]() {
1086  volPostProcEle->getRuleHook = VolRule();
1087  volPostProcEle->getUserPolynomialBase() =
1088  boost::shared_ptr<BaseFunction>(new CGGUserPolynomialBase());
1089  volPostProcEle->getOpPtrVector().push_back(
1090  new OpCalculateHVecTensorField<3, 3>(eP.piolaStress,
1091  eP.dataAtPts->approxPAtPts));
1092  volPostProcEle->getOpPtrVector().push_back(
1094  eP.bubbleField, eP.dataAtPts->approxPAtPts, MBMAXTYPE));
1095  volPostProcEle->getOpPtrVector().push_back(
1097  eP.streachTensor, eP.dataAtPts->logStreachTensorAtPts, MBTET));
1098  volPostProcEle->getOpPtrVector().push_back(
1100  eP.rotAxis, eP.dataAtPts->rotAxisAtPts, MBTET));
1101  volPostProcEle->getOpPtrVector().push_back(
1103  eP.materialGradient, eP.dataAtPts->GAtPts, MBTET));
1104  volPostProcEle->getOpPtrVector().push_back(
1105  new OpCalculateVectorFieldValues<3>(eP.spatialDisp,
1106  eP.dataAtPts->wAtPts, MBTET));
1107  volPostProcEle->getOpPtrVector().push_back(
1108  new OpCalculateRotationAndSpatialGradient(eP.rotAxis,
1109  eP.dataAtPts));
1110  volPostProcEle->getOpPtrVector().push_back(
1111  new OpCalculateStrainEnergy(eP.spatialDisp, eP.dataAtPts, gEnergy));
1112  };
1113 
1114  set_element_for_field_eval();
1115  set_element_for_post_process();
1116  }
1117 
1118  MoFEMErrorCode preProcess() { return 0; }
1119 
1120  MoFEMErrorCode operator()() { return 0; }
1121 
1122  MoFEMErrorCode postProcess() {
1124 
1125  auto get_str_time = [](auto ts_t) {
1126  std::ostringstream ss;
1127  ss << boost::str(boost::format("%d") %
1128  static_cast<int>(std::ceil(ts_t * 1e6)));
1129  std::string s = ss.str();
1130  return s;
1131  };
1132 
1133  PetscViewer viewer;
1134  CHKERR PetscViewerBinaryOpen(
1135  PETSC_COMM_WORLD, ("restart_" + get_str_time(ts_t) + ".dat").c_str(),
1136  FILE_MODE_WRITE, &viewer);
1137  CHKERR VecView(ts_u, viewer);
1138  CHKERR PetscViewerDestroy(&viewer);
1139 
1140  CHKERR eP.postProcessResults(1, "out_sol_elastic_" + get_str_time(ts_t) +
1141  ".h5m");
1142 
1143  // Loop boundary elements with traction boundary conditions
1144  *gEnergy = 0;
1145  CHKERR eP.mField.loop_finite_elements(problemPtr->getName(), "EP",
1146  *volPostProcEle,nullptr);
1147 
1148  double body_energy;
1149  MPI_Allreduce(gEnergy.get(), &body_energy, 1, MPI_DOUBLE, MPI_SUM,
1150  eP.mField.get_comm());
1151  CHKERR PetscPrintf(eP.mField.get_comm(),
1152  "Step %d time %3.4g strain energy %3.6e\n", ts_step,
1153  ts_t, body_energy);
1154 
1155  auto post_proc_at_points = [&](std::array<double, 3> point,
1156  std::string str) {
1158 
1159  dataFieldEval->setEvalPoints(point.data(), point.size() / 3);
1160 
1161  struct OpPrint : public VolOp {
1162 
1163  EshelbianCore &eP;
1164  std::array<double, 3> point;
1165  std::string str;
1166 
1167  OpPrint(EshelbianCore &ep, std::array<double, 3> &point,
1168  std::string &str)
1169  : VolOp(ep.spatialDisp, VolOp::OPROW), eP(ep), point(point),
1170  str(str) {}
1171 
1172  MoFEMErrorCode doWork(int side, EntityType type,
1175  if (type == MBTET) {
1176  if (getGaussPts().size2()) {
1177 
1178  auto t_h = getFTensor2FromMat<3, 3>(*(eP.dataAtPts->hAtPts));
1179  auto t_approx_P =
1180  getFTensor2FromMat<3, 3>(*(eP.dataAtPts->approxPAtPts));
1181 
1185  const double jac = dEterminant(t_h);
1187  t_cauchy(i, j) = (1. / jac) * (t_approx_P(i, k) * t_h(j, k));
1188 
1189  auto add = [&]() {
1190  std::ostringstream s;
1191  s << str << " elem " << getFEEntityHandle() << " ";
1192  return s.str();
1193  };
1194 
1195  std::ostringstream print;
1196  print << add() << "comm rank " << eP.mField.get_comm_rank()
1197  << std::endl;
1198  print << add() << "point " << getVectorAdaptor(point.data(), 3)
1199  << std::endl;
1200  print << add() << "coords at gauss pts " << getCoordsAtGaussPts()
1201  << std::endl;
1202  print << add() << "w " << (*eP.dataAtPts->wAtPts) << std::endl;
1203  print << add() << "Piola " << (*eP.dataAtPts->approxPAtPts)
1204  << std::endl;
1205  print << add() << "Cauchy " << t_cauchy << std::endl;
1206  print << std::endl;
1207  CHKERR PetscSynchronizedPrintf(eP.mField.get_comm(), "%s",
1208  print.str().c_str());
1209 
1210  }
1211  }
1213  }
1214  };
1215 
1216  if (auto fe_ptr = dataFieldEval->feMethodPtr.lock()) {
1217 
1218  fe_ptr->getOpPtrVector().push_back(new OpPrint(eP, point, str));
1219  CHKERR eP.mField.getInterface<FieldEvaluatorInterface>()
1220  ->evalFEAtThePoint3D(point.data(), 1e-12, problemPtr->getName(),
1221  "EP", dataFieldEval,
1222  eP.mField.get_comm_rank(),
1223  eP.mField.get_comm_rank(), MF_EXIST, QUIET);
1224  fe_ptr->getOpPtrVector().pop_back();
1225  }
1226 
1228  };
1229 
1230  // Points for Cook beam
1231  std::array<double, 3> pointA = {48., 60., 5.};
1232  CHKERR post_proc_at_points(pointA, "Point A");
1233  PetscSynchronizedFlush(eP.mField.get_comm(), PETSC_STDOUT);
1234 
1235  std::array<double, 3> pointB = {48. / 2., 44. + (60. - 44.) / 2., 0.};
1236  CHKERR post_proc_at_points(pointB, "Point B");
1237  PetscSynchronizedFlush(eP.mField.get_comm(), PETSC_STDOUT);
1238 
1239  std::array<double, 3> pointC = {48. / 2., (44. - 0.) / 2., 0.};
1240  CHKERR post_proc_at_points(pointC, "Point C");
1241  PetscSynchronizedFlush(eP.mField.get_comm(), PETSC_STDOUT);
1242 
1244  }
1245  };
1246 
1247  boost::shared_ptr<FEMethod> monitor_ptr(new Monitor(*this));
1248  ts_ctx->get_loops_to_do_Monitor().push_back(
1250 
1251  CHKERR TSAppendOptionsPrefix(ts, "elastic_");
1252  CHKERR TSSetFromOptions(ts);
1254 }
SmartPetscObj< DM > dmElastic
Elastic problem.
Set integration rule.
structure for User Loop Methods on finite elementsIt can be used to calculate stiffness matrices,...
VolumeElementForcesAndSourcesCoreSwitch< 0 > VolumeElementForcesAndSourcesCore
Volume finite element default.
SmartPetscObj< DM > dmElasticSchurOmega
Sub problem of dmElastic Schur.
boost::shared_ptr< SubProblemData > subProblemData
Calculate tenor field using vectorial base, i.e. Hdiv/Hcurl.
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
IFACE getInterface() const
Get interface pointer to pointer of interface.
PetscErrorCode TsMonitorSet(TS ts, PetscInt step, PetscReal t, Vec u, void *ctx)
Set monitor for TS solver.
Definition: TsCtx.cpp:202
keeps basic data about problemThis is low level structure with information about problem,...
Get values at integration pts for tensor filed rank 1, i.e. vector field.
EshelbianCore(MoFEM::Interface &m_field)
Field evaluator interface.
SmartPetscObj< DM > dmElasticSchurStreach
Sub problem of dmElastic Schur.
static double dEterminant(T &t)
Calculate the determinant of a 3x3 matrix or a tensor of rank 2.
Definition: Templates.hpp:452
FTensor::Index< 'm', 2 > m
Definition: PlasticOps.hpp:67
SmartPetscObj< DM > dmElasticSchurSpatialDisp
Sub problem of dmElastic Schur.
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:66
auto getVectorAdaptor(T1 ptr, const size_t n)
Get Vector adaptor.
Definition: Templates.hpp:42
Calculate tenor field using vectorial base, i.e. Hdiv/Hcurl.
#define CHKERR
Inline error check.
Definition: definitions.h:601
FTensor::Index< 'k', 2 > k
Definition: ContactOps.hpp:28
ForcesAndSourcesCore::UserDataOperator UserDataOperator
PetscErrorCode DMMoFEMGetTsCtx(DM dm, MoFEM::TsCtx **ts_ctx)
get MoFEM::TsCtx data structure
Definition: DMMMoFEM.cpp:987
SmartPetscObj< DM > dmElasticSchurBubble
Sub problem of dmElastic Schur.
Data on single entity (This is passed as argument to DataOperator::doWork)
FTensor::Index< 'j', 2 > j
Definition: ContactOps.hpp:27
structure to get information form mofem into DataForcesAndSourcesCore
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412
PetscErrorCode DMMoFEMGetProblemPtr(DM dm, const MoFEM::Problem **problem_ptr)
Get pointer to problem data structure.
Definition: DMMMoFEM.cpp:348
virtual MPI_Comm & get_comm() const =0
Calculate symmetric tensor field values at integration pts.
Get values at integration pts for tensor filed rank 2, i.e. matrix field.
FTensor::Index< 'i', 2 > i
[Common data]
Definition: ContactOps.hpp:26

◆ solveElastic()

MoFEMErrorCode EshelbianPlasticity::EshelbianCore::solveElastic ( TS  ts,
Vec  x 
)
Examples
ep.cpp.

Definition at line 1256 of file EshelbianPlasticity.cpp.

1256  {
1258 
1259  CHKERR TSSetDM(ts, dmElastic);
1260 
1261  SNES snes;
1262  CHKERR TSGetSNES(ts, &snes);
1263 
1264  PetscViewerAndFormat *vf;
1265  CHKERR PetscViewerAndFormatCreate(PETSC_VIEWER_STDOUT_WORLD,
1266  PETSC_VIEWER_DEFAULT, &vf);
1267  CHKERR SNESMonitorSet(
1268  snes,
1269  (MoFEMErrorCode(*)(SNES, PetscInt, PetscReal, void *))SNESMonitorFields,
1270  vf, (MoFEMErrorCode(*)(void **))PetscViewerAndFormatDestroy);
1271 
1272  PetscSection section;
1273  CHKERR DMGetDefaultSection(dmElastic, &section);
1274  int num_fields;
1275  CHKERR PetscSectionGetNumFields(section, &num_fields);
1276  for (int ff = 0; ff != num_fields; ff++) {
1277  const char *field_name;
1278  CHKERR PetscSectionGetFieldName(section, ff, &field_name);
1279  CHKERR PetscPrintf(PETSC_COMM_WORLD, "Field %d name %s\n", ff, field_name);
1280  }
1281 
1282  CHKERR DMoFEMMeshToLocalVector(dmElastic, x, INSERT_VALUES, SCATTER_FORWARD);
1283 
1284  // PetscRandom rctx;
1285  // PetscRandomCreate(PETSC_COMM_WORLD, &rctx);
1286  // VecSetRandom(x, rctx);
1287  // VecScale(x, -1e-1);
1288  // PetscRandomDestroy(&rctx);
1289  // CHKERR VecView(x, PETSC_VIEWER_STDOUT_WORLD);
1290 
1291  CHKERR VecGhostUpdateBegin(x, INSERT_VALUES, SCATTER_FORWARD);
1292  CHKERR VecGhostUpdateEnd(x, INSERT_VALUES, SCATTER_FORWARD);
1293 
1294  // CHKERR VecView(x, PETSC_VIEWER_STDOUT_WORLD);
1295 
1296  // Adding field split solver
1297 
1298  KSP ksp;
1299  CHKERR SNESGetKSP(snes, &ksp);
1300  PC pc;
1301  CHKERR KSPGetPC(ksp, &pc);
1302  PetscBool is_uu_field_split;
1303  PetscObjectTypeCompare((PetscObject)pc, PCFIELDSPLIT, &is_uu_field_split);
1304  if (is_uu_field_split) {
1305 
1306  const MoFEM::Problem *schur_uu_ptr;
1308  if (auto uu_data = schur_uu_ptr->subProblemData) {
1309 
1310  const MoFEM::Problem *prb_ptr;
1312  map<std::string, IS> is_map;
1313  for (int ff = 0; ff != num_fields; ff++) {
1314  const char *field_name;
1315  CHKERR PetscSectionGetFieldName(section, ff, &field_name);
1316  CHKERR mField.getInterface<ISManager>()->isCreateProblemFieldAndRank(
1317  prb_ptr->getName(), ROW, field_name, 0, MAX_DOFS_ON_ENTITY,
1318  &is_map[field_name]);
1319  }
1320  // CHKERR mField.getInterface<ISManager>()
1321  // ->isCreateProblemFieldAndEntityType(
1322  // prb_ptr->getName(), ROW, piolaStress, MBTET, MBTET, 0,
1323  // MAX_DOFS_ON_ENTITY, &is_map["T_STRESS_ON_TETS"]);
1324  // CHKERR mField.getInterface<ISManager>()
1325  // ->isCreateProblemFieldAndEntityType(
1326  // prb_ptr->getName(), ROW, piolaStress, MBTRI, MBTRI, 0,
1327  // MAX_DOFS_ON_ENTITY, &is_map["T_STRESS_ON_TRIS"]);
1328 
1329  CHKERR uu_data->getRowIs(&is_map["E_IS_SUU"]);
1330 
1331  CHKERR PCFieldSplitSetIS(pc, NULL, is_map[streachTensor]);
1332  CHKERR PCFieldSplitSetIS(pc, NULL, is_map["E_IS_SUU"]);
1333 
1334  CHKERR PCFieldSplitSetSchurPre(pc, PC_FIELDSPLIT_SCHUR_PRE_USER,
1335  schurAssembly->Suu);
1336 
1337  CHKERR PCSetUp(pc);
1338  PetscInt n;
1339  KSP *uu_ksp;
1340  CHKERR PCFieldSplitGetSubKSP(pc, &n, &uu_ksp);
1341  PC bubble_pc;
1342  CHKERR KSPGetPC(uu_ksp[1], &bubble_pc);
1343  PetscBool is_bubble_field_split;
1344  PetscObjectTypeCompare((PetscObject)bubble_pc, PCFIELDSPLIT,
1345  &is_bubble_field_split);
1346  if (is_bubble_field_split) {
1347 
1348  const MoFEM::Problem *schur_bubble_ptr;
1349  CHKERR DMMoFEMGetProblemPtr(dmElasticSchurBubble, &schur_bubble_ptr);
1350  if (auto bubble_data = schur_bubble_ptr->subProblemData) {
1351 
1352  CHKERR bubble_data->getRowIs(&is_map["E_IS_BUBBLE"]);
1353 
1354  AO uu_ao;
1355  CHKERR uu_data->getRowMap(&uu_ao);
1356 
1357  CHKERR AOApplicationToPetscIS(uu_ao, is_map[bubbleField]);
1358  CHKERR PCFieldSplitSetIS(bubble_pc, NULL, is_map[bubbleField]);
1359  CHKERR PCFieldSplitSetIS(bubble_pc, NULL, is_map["E_IS_BUBBLE"]);
1360  CHKERR PCFieldSplitSetSchurPre(
1361  bubble_pc, PC_FIELDSPLIT_SCHUR_PRE_USER, schurAssembly->SBubble);
1362 
1363  CHKERR PCSetUp(bubble_pc);
1364  PetscInt bubble_n;
1365  KSP *bubble_ksp;
1366  CHKERR PCFieldSplitGetSubKSP(bubble_pc, &bubble_n, &bubble_ksp);
1367  PC omega_pc;
1368  CHKERR KSPGetPC(bubble_ksp[1], &omega_pc);
1369  PetscBool is_omega_field_split;
1370  PetscObjectTypeCompare((PetscObject)omega_pc, PCFIELDSPLIT,
1371  &is_omega_field_split);
1372 
1373  if (is_omega_field_split) {
1374 
1375  const MoFEM::Problem *schur_omega_ptr;
1376  CHKERR DMMoFEMGetProblemPtr(dmElasticSchurOmega, &schur_omega_ptr);
1377  if (auto omega_data = schur_omega_ptr->subProblemData) {
1378 
1379  AO bubble_ao;
1380  CHKERR bubble_data->getRowMap(&bubble_ao);
1381 
1382  CHKERR AOApplicationToPetscIS(uu_ao, is_map[rotAxis]);
1383  CHKERR AOApplicationToPetscIS(bubble_ao, is_map[rotAxis]);
1384  CHKERR omega_data->getRowIs(&is_map["E_IS_OMEGA"]);
1385 
1386  CHKERR PCFieldSplitSetIS(omega_pc, NULL, is_map[rotAxis]);
1387  CHKERR PCFieldSplitSetIS(omega_pc, NULL, is_map["E_IS_OMEGA"]);
1388  CHKERR PCFieldSplitSetSchurPre(omega_pc,
1389  PC_FIELDSPLIT_SCHUR_PRE_USER,
1390  schurAssembly->SOmega);
1391 
1392  CHKERR PCSetUp(omega_pc);
1393  PetscInt omega_n;
1394  KSP *omega_ksp;
1395  CHKERR PCFieldSplitGetSubKSP(omega_pc, &omega_n, &omega_ksp);
1396  PC w_pc;
1397  CHKERR KSPGetPC(omega_ksp[1], &w_pc);
1398  PetscBool is_w_field_split;
1399  PetscObjectTypeCompare((PetscObject)w_pc, PCFIELDSPLIT,
1400  &is_w_field_split);
1401  if (is_w_field_split) {
1402 
1403  const MoFEM::Problem *schur_w_ptr;
1405  &schur_w_ptr);
1406  if (auto w_data = schur_w_ptr->subProblemData) {
1407 
1408  AO omega_ao;
1409  CHKERR omega_data->getRowMap(&omega_ao);
1410 
1411  CHKERR AOApplicationToPetscIS(uu_ao, is_map[spatialDisp]);
1412  CHKERR AOApplicationToPetscIS(bubble_ao, is_map[spatialDisp]);
1413  CHKERR AOApplicationToPetscIS(omega_ao, is_map[spatialDisp]);
1414  CHKERR w_data->getRowIs(&is_map["E_IS_W"]);
1415 
1416  CHKERR PCFieldSplitSetIS(w_pc, NULL, is_map[spatialDisp]);
1417  CHKERR PCFieldSplitSetIS(w_pc, NULL, is_map["E_IS_W"]);
1418  CHKERR PCFieldSplitSetSchurPre(
1419  w_pc, PC_FIELDSPLIT_SCHUR_PRE_USER, schurAssembly->Sw);
1420 
1421  CHKERR AODestroy(&omega_ao);
1422  }
1423  }
1424 
1425  CHKERR PetscFree(omega_ksp);
1426  }
1427  }
1428  CHKERR PetscFree(bubble_ksp);
1429  CHKERR AODestroy(&uu_ao);
1430  }
1431  }
1432  CHKERR PetscFree(uu_ksp);
1433 
1434  for (auto &m : is_map)
1435  CHKERR ISDestroy(&m.second);
1436  }
1437  }
1438 
1439  CHKERR TSSolve(ts, x);
1441 }
SmartPetscObj< DM > dmElastic
Elastic problem.
boost::shared_ptr< EpFEMethod > schurAssembly
SmartPetscObj< DM > dmElasticSchurOmega
Sub problem of dmElastic Schur.
Section manager is used to create indexes and sectionsFIXME: ISManager is not properly testsed by ato...
Definition: ISManager.hpp:36
PetscErrorCode DMoFEMMeshToLocalVector(DM dm, Vec l, InsertMode mode, ScatterMode scatter_mode)
set local (or ghosted) vector values on mesh for partition only
Definition: DMMMoFEM.cpp:431
boost::shared_ptr< SubProblemData > subProblemData
FTensor::Index< 'n', 2 > n
Definition: PlasticOps.hpp:68
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
keeps basic data about problemThis is low level structure with information about problem,...
MoFEMErrorCode getInterface(const MOFEMuuid &uuid, IFACE *&iface) const
Get interface by uuid and return reference to pointer of interface.
SmartPetscObj< DM > dmElasticSchurStreach
Sub problem of dmElastic Schur.
FTensor::Index< 'm', 2 > m
Definition: PlasticOps.hpp:67
SmartPetscObj< DM > dmElasticSchurSpatialDisp
Sub problem of dmElastic Schur.
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:66
#define CHKERR
Inline error check.
Definition: definitions.h:601
SmartPetscObj< DM > dmElasticSchurBubble
Sub problem of dmElastic Schur.
#define MAX_DOFS_ON_ENTITY
Maximal number of DOFs on entity.
Definition: definitions.h:303
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412
PetscErrorCode DMMoFEMGetProblemPtr(DM dm, const MoFEM::Problem **problem_ptr)
Get pointer to problem data structure.
Definition: DMMMoFEM.cpp:348
std::string getName() const

Member Data Documentation

◆ alpha_u

double EshelbianPlasticity::EshelbianCore::alpha_u

Definition at line 1229 of file EshelbianPlasticity.hpp.

◆ alpha_w

double EshelbianPlasticity::EshelbianCore::alpha_w

Definition at line 1230 of file EshelbianPlasticity.hpp.

◆ bcSpatialDispVecPtr

boost::shared_ptr<BcDispVec> EshelbianPlasticity::EshelbianCore::bcSpatialDispVecPtr

Definition at line 1235 of file EshelbianPlasticity.hpp.

◆ bcSpatialFreeTraction

boost::shared_ptr<TractionFreeBc> EshelbianPlasticity::EshelbianCore::bcSpatialFreeTraction

Definition at line 1238 of file EshelbianPlasticity.hpp.

◆ bcSpatialRotationVecPtr

boost::shared_ptr<BcRotVec> EshelbianPlasticity::EshelbianCore::bcSpatialRotationVecPtr

Definition at line 1236 of file EshelbianPlasticity.hpp.

◆ bcSpatialTraction

boost::shared_ptr<TractionBcVec> EshelbianPlasticity::EshelbianCore::bcSpatialTraction

Definition at line 1237 of file EshelbianPlasticity.hpp.

◆ bubbleField

const std::string EshelbianPlasticity::EshelbianCore::bubbleField

Definition at line 1218 of file EshelbianPlasticity.hpp.

◆ dataAtPts

boost::shared_ptr<DataAtIntegrationPts> EshelbianPlasticity::EshelbianCore::dataAtPts

Definition at line 1191 of file EshelbianPlasticity.hpp.

◆ dM

SmartPetscObj<DM> EshelbianPlasticity::EshelbianCore::dM

Coupled problem all fields.

Definition at line 1200 of file EshelbianPlasticity.hpp.

◆ dmElastic

SmartPetscObj<DM> EshelbianPlasticity::EshelbianCore::dmElastic

Elastic problem.

Examples
ep.cpp.

Definition at line 1201 of file EshelbianPlasticity.hpp.

◆ dmElasticSchurBubble

SmartPetscObj<DM> EshelbianPlasticity::EshelbianCore::dmElasticSchurBubble

Sub problem of dmElastic Schur.

Definition at line 1203 of file EshelbianPlasticity.hpp.

◆ dmElasticSchurOmega

SmartPetscObj<DM> EshelbianPlasticity::EshelbianCore::dmElasticSchurOmega

Sub problem of dmElastic Schur.

Definition at line 1206 of file EshelbianPlasticity.hpp.

◆ dmElasticSchurSpatialDisp

SmartPetscObj<DM> EshelbianPlasticity::EshelbianCore::dmElasticSchurSpatialDisp

Sub problem of dmElastic Schur.

Definition at line 1205 of file EshelbianPlasticity.hpp.

◆ dmElasticSchurStreach

SmartPetscObj<DM> EshelbianPlasticity::EshelbianCore::dmElasticSchurStreach

Sub problem of dmElastic Schur.

Definition at line 1202 of file EshelbianPlasticity.hpp.

◆ dmMaterial

SmartPetscObj<DM> EshelbianPlasticity::EshelbianCore::dmMaterial

Material problem.

Definition at line 1207 of file EshelbianPlasticity.hpp.

◆ elasticBcLhs

boost::shared_ptr<EpElement<FaceElementForcesAndSourcesCore> > EshelbianPlasticity::EshelbianCore::elasticBcLhs

Definition at line 1196 of file EshelbianPlasticity.hpp.

◆ elasticBcRhs

boost::shared_ptr<EpElement<FaceElementForcesAndSourcesCore> > EshelbianPlasticity::EshelbianCore::elasticBcRhs

Definition at line 1197 of file EshelbianPlasticity.hpp.

◆ elasticFeLhs

boost::shared_ptr<EpElement<VolumeElementForcesAndSourcesCore> > EshelbianPlasticity::EshelbianCore::elasticFeLhs

Definition at line 1195 of file EshelbianPlasticity.hpp.

◆ elasticFeRhs

boost::shared_ptr<EpElement<VolumeElementForcesAndSourcesCore> > EshelbianPlasticity::EshelbianCore::elasticFeRhs

Definition at line 1194 of file EshelbianPlasticity.hpp.

◆ elementVolumeName

const std::string EshelbianPlasticity::EshelbianCore::elementVolumeName

Definition at line 1220 of file EshelbianPlasticity.hpp.

◆ eshelbyStress

const std::string EshelbianPlasticity::EshelbianCore::eshelbyStress

Definition at line 1210 of file EshelbianPlasticity.hpp.

◆ essentialBcElement

const std::string EshelbianPlasticity::EshelbianCore::essentialBcElement

Definition at line 1222 of file EshelbianPlasticity.hpp.

◆ lambdaField

const std::string EshelbianPlasticity::EshelbianCore::lambdaField

Definition at line 1217 of file EshelbianPlasticity.hpp.

◆ materialDisp

const std::string EshelbianPlasticity::EshelbianCore::materialDisp

Definition at line 1212 of file EshelbianPlasticity.hpp.

◆ materialGradient

const std::string EshelbianPlasticity::EshelbianCore::materialGradient

Definition at line 1215 of file EshelbianPlasticity.hpp.

◆ materialOrder

int EshelbianPlasticity::EshelbianCore::materialOrder

Definition at line 1228 of file EshelbianPlasticity.hpp.

◆ mField

MoFEM::Interface& EshelbianPlasticity::EshelbianCore::mField

Definition at line 1189 of file EshelbianPlasticity.hpp.

◆ naturalBcElement

const std::string EshelbianPlasticity::EshelbianCore::naturalBcElement

Definition at line 1221 of file EshelbianPlasticity.hpp.

◆ physicalEquations

boost::shared_ptr<PhysicalEquations> EshelbianPlasticity::EshelbianCore::physicalEquations

Definition at line 1192 of file EshelbianPlasticity.hpp.

◆ piolaStress

const std::string EshelbianPlasticity::EshelbianCore::piolaStress

Definition at line 1209 of file EshelbianPlasticity.hpp.

◆ preconditioner_eps

double EshelbianPlasticity::EshelbianCore::preconditioner_eps

Definition at line 1231 of file EshelbianPlasticity.hpp.

◆ rotAxis

const std::string EshelbianPlasticity::EshelbianCore::rotAxis

Definition at line 1214 of file EshelbianPlasticity.hpp.

◆ schurAssembly

boost::shared_ptr<EpFEMethod> EshelbianPlasticity::EshelbianCore::schurAssembly

Definition at line 1198 of file EshelbianPlasticity.hpp.

◆ spaceOrder

int EshelbianPlasticity::EshelbianCore::spaceOrder

Definition at line 1227 of file EshelbianPlasticity.hpp.

◆ spatialDisp

const std::string EshelbianPlasticity::EshelbianCore::spatialDisp

Definition at line 1211 of file EshelbianPlasticity.hpp.

◆ streachTensor

const std::string EshelbianPlasticity::EshelbianCore::streachTensor

Definition at line 1213 of file EshelbianPlasticity.hpp.

◆ tauField

const std::string EshelbianPlasticity::EshelbianCore::tauField

Definition at line 1216 of file EshelbianPlasticity.hpp.


The documentation for this struct was generated from the following files: