v0.13.1
Loading...
Searching...
No Matches
Public Member Functions | Public Attributes | List of all members
FractureMechanics::GriffithForceElement::SaveGriffithForceOnTags Struct Reference

#include <users_modules/fracture_mechanics/src/GriffithForceElement.hpp>

Collaboration diagram for FractureMechanics::GriffithForceElement::SaveGriffithForceOnTags:
[legend]

Public Member Functions

 SaveGriffithForceOnTags (MoFEM::Interface &m_field, Tag th, Range &tris, Tag *th_coords=nullptr)
 
MoFEMErrorCode operator() ()
 

Public Attributes

MoFEM::InterfacemField
 
Tag tH
 
RangetRis
 
Tag * thCoordsPtr
 
AuxFunctions< doubleauxFun
 

Detailed Description

Definition at line 530 of file GriffithForceElement.hpp.

Constructor & Destructor Documentation

◆ SaveGriffithForceOnTags()

FractureMechanics::GriffithForceElement::SaveGriffithForceOnTags::SaveGriffithForceOnTags ( MoFEM::Interface m_field,
Tag  th,
Range tris,
Tag *  th_coords = nullptr 
)
inline

Member Function Documentation

◆ operator()()

MoFEMErrorCode FractureMechanics::GriffithForceElement::SaveGriffithForceOnTags::operator() ( )
inline

Definition at line 544 of file GriffithForceElement.hpp.

544 {
546
547 Range verts;
548 CHKERR mField.get_moab().get_connectivity(tRis, verts, true);
549 {
550 std::vector<double> data(verts.size() * 3, 0);
551 CHKERR mField.get_moab().tag_set_data(tH, verts, &*data.begin());
552 }
553 MatrixDouble diff_n(3, 2);
554 CHKERR ShapeDiffMBTRI(&*diff_n.data().begin());
555
556 for (Range::iterator tit = tRis.begin(); tit != tRis.end(); ++tit) {
557 int num_nodes;
558 const EntityHandle *conn;
559 CHKERR mField.get_moab().get_connectivity(*tit, conn, num_nodes, true);
560 VectorDouble9 coords(9);
561
562 int nb_dofs = num_nodes * 3;
563 auxFun.referenceCoords.resize(9, false);
564 auxFun.currentCoords.resize(9, false);
565
566 if (thCoordsPtr) {
567 CHKERR mField.get_moab().tag_get_data(*thCoordsPtr, conn, num_nodes,
568 &*coords.begin());
569 } else {
570 CHKERR mField.get_moab().get_coords(conn, num_nodes,
571 &*coords.begin());
572 }
573 for (int dd = 0; dd != nb_dofs; dd++) {
574 auxFun.currentCoords[dd] = coords[dd];
575 }
576 if (thCoordsPtr)
577 CHKERR mField.get_moab().get_coords(conn, num_nodes,
578 &*coords.begin());
579 for (int dd = 0; dd != nb_dofs; dd++) {
580 auxFun.referenceCoords[dd] = coords[dd];
581 }
582
583 CHKERR auxFun.calculateGriffithForce(1, 0.5, diff_n);
584
585 double *tag_vals[3];
586 CHKERR mField.get_moab().tag_get_by_ptr(tH, conn, 3,
587 (const void **)tag_vals);
588 for (int nn = 0; nn != 3; ++nn) {
589 for (int dd = 0; dd != 3; ++dd) {
590 tag_vals[nn][dd] += auxFun.griffithForce[3 * nn + dd];
591 }
592 }
593 }
594
596 }
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
#define CHKERR
Inline error check.
Definition: definitions.h:535
PetscErrorCode ShapeDiffMBTRI(double *diffN)
calculate derivatives of shape functions
Definition: fem_tools.c:194
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
UBlasMatrix< double > MatrixDouble
Definition: Types.hpp:77
VectorBoundedArray< double, 9 > VectorDouble9
Definition: Types.hpp:96
MoFEMErrorCode calculateGriffithForce(const double gc, const double beta, const MatrixDouble &diff_n)
virtual moab::Interface & get_moab()=0

Member Data Documentation

◆ auxFun

AuxFunctions<double> FractureMechanics::GriffithForceElement::SaveGriffithForceOnTags::auxFun

Definition at line 542 of file GriffithForceElement.hpp.

◆ mField

MoFEM::Interface& FractureMechanics::GriffithForceElement::SaveGriffithForceOnTags::mField

Definition at line 532 of file GriffithForceElement.hpp.

◆ tH

Tag FractureMechanics::GriffithForceElement::SaveGriffithForceOnTags::tH

Definition at line 534 of file GriffithForceElement.hpp.

◆ thCoordsPtr

Tag* FractureMechanics::GriffithForceElement::SaveGriffithForceOnTags::thCoordsPtr

Definition at line 536 of file GriffithForceElement.hpp.

◆ tRis

Range& FractureMechanics::GriffithForceElement::SaveGriffithForceOnTags::tRis

Definition at line 535 of file GriffithForceElement.hpp.


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