v0.13.1
Public Member Functions | Public Attributes | List of all members
CohesiveElement::CohesiveInterfaceElement::PhysicalEquation Struct Reference

Constitutive (physical) equation for interface. More...

#include <users_modules/basic_finite_elements/cohesive_interface/src/CohesiveInterfaceElement.hpp>

Collaboration diagram for CohesiveElement::CohesiveInterfaceElement::PhysicalEquation:
[legend]

Public Member Functions

 PhysicalEquation (MoFEM::Interface &m_field)
 
virtual ~PhysicalEquation ()
 
MoFEMErrorCode iNitailise (const FEMethod *fe_method)
 Initialize history variable data. More...
 
double calcG (int gg, MatrixDouble gap_loc)
 Calculate gap opening. More...
 
MoFEMErrorCode getKappa (int nb_gauss_pts, const FEMethod *fe_method)
 Get pointer from the mesh to histoy variables \(\kappa\). More...
 
MoFEMErrorCode calcDglob (const double omega, MatrixDouble &R)
 Calculate stiffness material matrix. More...
 
MoFEMErrorCode calcOmega (const double kappa, double &omega)
 Calculate damage. More...
 
MoFEMErrorCode calcTangetDglob (const double omega, double g, const VectorDouble &gap_loc, MatrixDouble &R)
 Calculate tangent material stiffness. More...
 
virtual MoFEMErrorCode calculateTraction (VectorDouble &traction, int gg, CommonData &common_data, const FEMethod *fe_method)
 Calculate tractions. More...
 
virtual MoFEMErrorCode calculateTangentStiffeness (MatrixDouble &tangent_matrix, int gg, CommonData &common_data, const FEMethod *fe_method)
 Calculate tangent stiffness. More...
 
virtual MoFEMErrorCode updateHistory (CommonData &common_data, const FEMethod *fe_method)
 Update history variables when converged. More...
 

Public Attributes

MoFEM::InterfacemField
 
bool isInitialised
 
double h
 
double youngModulus
 
double beta
 
double ft
 
double Gf
 
Range pRisms
 
Tag thKappa
 
Tag thDamagedPrism
 
double E0
 
double g0
 
double kappa1
 
doublekappaPtr
 
int kappaSize
 
MatrixDouble Dglob
 
MatrixDouble Dloc
 

Detailed Description

Constitutive (physical) equation for interface.

This is linear degradation model. Material parameters are: strength \(f_t\), interface fracture energy \(G_f\), elastic material stiffness \(E\). Parameter \(\beta\) controls how interface opening is calculated.

Model parameter is interface penalty thickness \(h\).

Definition at line 51 of file CohesiveInterfaceElement.hpp.

Constructor & Destructor Documentation

◆ PhysicalEquation()

CohesiveElement::CohesiveInterfaceElement::PhysicalEquation::PhysicalEquation ( MoFEM::Interface m_field)

◆ ~PhysicalEquation()

virtual CohesiveElement::CohesiveInterfaceElement::PhysicalEquation::~PhysicalEquation ( )
virtual

Definition at line 58 of file CohesiveInterfaceElement.hpp.

58 {
59 }

Member Function Documentation

◆ calcDglob()

MoFEMErrorCode CohesiveElement::CohesiveInterfaceElement::PhysicalEquation::calcDglob ( const double  omega,
MatrixDouble &  R 
)

Calculate stiffness material matrix.

\[ \mathbf{D}_\textrm{loc} = (1-\Omega) \mathbf{I} E_0 \]

where \(E_0\) is initial interface penalty stiffness

\[ \mathbf{D}_\textrm{glob} = \mathbf{R}^\textrm{T} \mathbf{D}_\textrm{loc}\mathbf{R} \]

Definition at line 136 of file CohesiveInterfaceElement.hpp.

136 {
138 Dglob.resize(3,3);
139 Dloc.resize(3,3);
140 Dloc.clear();
141 double E = (1-omega)*E0;
142 Dloc(0,0) = E;
143 Dloc(1,1) = E;
144 Dloc(2,2) = E;
145 Dglob = prod( Dloc, R );
146 Dglob = prod( trans(R), Dglob );
148 }
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:447
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:440
constexpr double omega

◆ calcG()

double CohesiveElement::CohesiveInterfaceElement::PhysicalEquation::calcG ( int  gg,
MatrixDouble  gap_loc 
)

Calculate gap opening.

\[ g = \sqrt{ g_n^2 + \beta(g_{s1}^2 + g_{s2}^2)} \]

Definition at line 95 of file CohesiveInterfaceElement.hpp.

95 {
96 return sqrt(pow(gap_loc(gg,0),2)+beta*(pow(gap_loc(gg,1),2)+pow(gap_loc(gg,2),2)));
97 }

◆ calcOmega()

MoFEMErrorCode CohesiveElement::CohesiveInterfaceElement::PhysicalEquation::calcOmega ( const double  kappa,
double omega 
)

Calculate damage.

\[ \Omega = \frac{1}{2} \frac{(2 G_f E_0+f_t^2)\kappa}{(ft+E_0 \kappa)G_f} \]

Definition at line 157 of file CohesiveInterfaceElement.hpp.

◆ calcTangetDglob()

MoFEMErrorCode CohesiveElement::CohesiveInterfaceElement::PhysicalEquation::calcTangetDglob ( const double  omega,
double  g,
const VectorDouble &  gap_loc,
MatrixDouble &  R 
)

Calculate tangent material stiffness.

Definition at line 173 of file CohesiveInterfaceElement.hpp.

173 {
175 Dglob.resize(3,3);
176 Dloc.resize(3,3);
177 double domega = 0.5*(2*Gf*E0+ft*ft)/((ft+(g-ft/E0)*E0)*Gf) - 0.5*((g-ft/E0)*(2*Gf*E0+ft*ft)*E0)/(pow(ft+(g-ft/E0)*E0,2)*Gf);
178 Dloc.resize(3,3);
179 //r0
180 Dloc(0,0) = (1-omega)*E0 - domega*E0*gap_loc[0]*gap_loc[0]/g;
181 Dloc(0,1) = -domega*E0*gap_loc[0]*beta*gap_loc[1]/g;
182 Dloc(0,2) = -domega*E0*gap_loc[0]*beta*gap_loc[2]/g;
183 //r1
184 Dloc(1,0) = -domega*E0*gap_loc[1]*gap_loc[0]/g;
185 Dloc(1,1) = (1-omega)*E0 - domega*E0*gap_loc[1]*beta*gap_loc[1]/g;
186 Dloc(1,2) = -domega*E0*gap_loc[1]*beta*gap_loc[2]/g;
187 //r2
188 Dloc(2,0) = -domega*E0*gap_loc[2]*gap_loc[0]/g;
189 Dloc(2,1) = -domega*E0*gap_loc[2]*beta*gap_loc[1]/g;
190 Dloc(2,2) = (1-omega)*E0 - domega*E0*gap_loc[2]*beta*gap_loc[2]/g;
191 Dglob = prod(Dloc,R);
192 Dglob = prod(trans(R),Dglob);
194 }
constexpr double g

◆ calculateTangentStiffeness()

virtual MoFEMErrorCode CohesiveElement::CohesiveInterfaceElement::PhysicalEquation::calculateTangentStiffeness ( MatrixDouble &  tangent_matrix,
int  gg,
CommonData common_data,
const FEMethod fe_method 
)
virtual

Calculate tangent stiffness.

Definition at line 231 of file CohesiveInterfaceElement.hpp.

235 {
237
238 try {
239 if(!isInitialised) {
240 ierr = iNitailise(fe_method); CHKERRG(ierr);
241 isInitialised = true;
242 }
243 if(gg==0) {
244 ierr = getKappa(common_data.gapGlob.size1(),fe_method); CHKERRG(ierr);
245 }
246 double g = calcG(gg,common_data.gapLoc);
247 double kappa = fmax(g-g0,kappaPtr[gg]);
248 double omega = 0;
250 //std::cerr << gg << " " << omega << std::endl;
251 int iter;
252 ierr = SNESGetIterationNumber(fe_method->snes,&iter); CHKERRG(ierr);
253 if((kappa <= kappaPtr[gg])||(kappa>=kappa1)||(iter <= 1)) {
254 ierr = calcDglob(omega,common_data.R[gg]); CHKERRG(ierr);
255 } else {
256 ublas::matrix_row<MatrixDouble > g_loc(common_data.gapLoc,gg);
257 ierr = calcTangetDglob(omega,g,g_loc,common_data.R[gg]);
258 }
259 tangent_matrix.resize(3,3);
260 noalias(tangent_matrix) = Dglob;
261 //std::cerr << "t " << tangent_matrix << std::endl;
262 } catch (const std::exception& ex) {
263 std::ostringstream ss;
264 ss << "throw in method: " << ex.what() << std::endl;
265 SETERRQ(PETSC_COMM_SELF,1,ss.str().c_str());
266 }
268 }
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:483
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Exceptions.hpp:76
double calcG(int gg, MatrixDouble gap_loc)
Calculate gap opening.
MoFEMErrorCode iNitailise(const FEMethod *fe_method)
Initialize history variable data.
MoFEMErrorCode calcOmega(const double kappa, double &omega)
Calculate damage.
MoFEMErrorCode calcTangetDglob(const double omega, double g, const VectorDouble &gap_loc, MatrixDouble &R)
Calculate tangent material stiffness.
MoFEMErrorCode calcDglob(const double omega, MatrixDouble &R)
Calculate stiffness material matrix.
MoFEMErrorCode getKappa(int nb_gauss_pts, const FEMethod *fe_method)
Get pointer from the mesh to histoy variables .

◆ calculateTraction()

virtual MoFEMErrorCode CohesiveElement::CohesiveInterfaceElement::PhysicalEquation::calculateTraction ( VectorDouble &  traction,
int  gg,
CommonData common_data,
const FEMethod fe_method 
)
virtual

Calculate tractions.

\[ \mathbf{t} = \mathbf{D}_\textrm{glob}\mathbf{g} \]

Definition at line 203 of file CohesiveInterfaceElement.hpp.

207 {
209
210 if(!isInitialised) {
211 ierr = iNitailise(fe_method); CHKERRG(ierr);
212 isInitialised = true;
213 }
214 if(gg==0) {
215 ierr = getKappa(common_data.gapGlob.size1(),fe_method); CHKERRG(ierr);
216 }
217 double g = calcG(gg,common_data.gapLoc);
218 double kappa = fmax(g-g0,kappaPtr[gg]);
219 double omega = 0;
221 //std::cerr << gg << " " << omega << std::endl;
222 ierr = calcDglob(omega,common_data.R[gg]); CHKERRG(ierr);
223 traction.resize(3);
224 ublas::matrix_row<MatrixDouble > gap_glob(common_data.gapGlob,gg);
225 noalias(traction) = prod(Dglob,gap_glob);
227 }

◆ getKappa()

MoFEMErrorCode CohesiveElement::CohesiveInterfaceElement::PhysicalEquation::getKappa ( int  nb_gauss_pts,
const FEMethod fe_method 
)

Get pointer from the mesh to histoy variables \(\kappa\).

Definition at line 104 of file CohesiveInterfaceElement.hpp.

104 {
106 EntityHandle ent = fe_method->numeredEntFiniteElementPtr->getEnt();
107
108 rval = mField.get_moab().tag_get_by_ptr(thKappa,&ent,1,(const void **)&kappaPtr,&kappaSize);
109 if(rval != MB_SUCCESS || kappaSize != nb_gauss_pts) {
111 kappa.resize(nb_gauss_pts);
112 kappa.clear();
113 int tag_size[1];
114 tag_size[0] = nb_gauss_pts;
115 void const* tag_data[] = { &kappa[0] };
116 rval = mField.get_moab().tag_set_by_ptr(thKappa,&ent,1,tag_data,tag_size); CHKERRG(rval);
117 rval = mField.get_moab().tag_get_by_ptr(thKappa,&ent,1,(const void **)&kappaPtr,&kappaSize); CHKERRG(rval);
118 }
120 }
static MoFEMErrorCodeGeneric< moab::ErrorCode > rval
Definition: Exceptions.hpp:74
UBlasVector< double > VectorDouble
Definition: Types.hpp:68
virtual moab::Interface & get_moab()=0

◆ iNitailise()

MoFEMErrorCode CohesiveElement::CohesiveInterfaceElement::PhysicalEquation::iNitailise ( const FEMethod fe_method)

Initialize history variable data.

Create tag on the prism/interface to store damage history variable

Definition at line 72 of file CohesiveInterfaceElement.hpp.

72 {
74
75 double def_damaged = 0;
76 rval = mField.get_moab().tag_get_handle(
77 "DAMAGED_PRISM",1,MB_TYPE_INTEGER,thDamagedPrism,MB_TAG_CREAT|MB_TAG_SPARSE,&def_damaged
78 ); MOAB_THROW(rval);
79 const int def_len = 0;
80 rval = mField.get_moab().tag_get_handle("_KAPPA",def_len,MB_TYPE_DOUBLE,
81 thKappa,MB_TAG_CREAT|MB_TAG_SPARSE|MB_TAG_VARLEN,NULL); CHKERRG(rval);
83 g0 = ft/E0;
84 kappa1 = 2*Gf/ft;
86 }
#define MOAB_THROW(err)
Check error code of MoAB function and throw MoFEM exception.
Definition: definitions.h:541

◆ updateHistory()

virtual MoFEMErrorCode CohesiveElement::CohesiveInterfaceElement::PhysicalEquation::updateHistory ( CommonData common_data,
const FEMethod fe_method 
)
virtual

Update history variables when converged.

Definition at line 272 of file CohesiveInterfaceElement.hpp.

274 {
276
277
278 if(!isInitialised) {
279 ierr = iNitailise(fe_method); CHKERRG(ierr);
280 isInitialised = true;
281 }
282 ierr = getKappa(common_data.gapGlob.size1(),fe_method); CHKERRG(ierr);
283 bool all_gauss_pts_damaged = true;
284 for(unsigned int gg = 0;gg<common_data.gapGlob.size1();gg++) {
285 double omega = 0;
286 double g = calcG(gg,common_data.gapLoc);
287 double kappa = fmax(g-g0,kappaPtr[gg]);
288 kappaPtr[gg] = kappa;
290 //if(omega < 1.) {
291 all_gauss_pts_damaged = false;
292 //}
293 }
294 if(all_gauss_pts_damaged) {
295 EntityHandle ent = fe_method->numeredEntFiniteElementPtr->getEnt();
296 int set_prism_as_demaged = 1;
297 rval = mField.get_moab().tag_set_data(thDamagedPrism,&ent,1,&set_prism_as_demaged); CHKERRG(rval);
298 }
300 }

Member Data Documentation

◆ beta

double CohesiveElement::CohesiveInterfaceElement::PhysicalEquation::beta

Definition at line 61 of file CohesiveInterfaceElement.hpp.

◆ Dglob

MatrixDouble CohesiveElement::CohesiveInterfaceElement::PhysicalEquation::Dglob

Definition at line 122 of file CohesiveInterfaceElement.hpp.

◆ Dloc

MatrixDouble CohesiveElement::CohesiveInterfaceElement::PhysicalEquation::Dloc

Definition at line 122 of file CohesiveInterfaceElement.hpp.

◆ E0

double CohesiveElement::CohesiveInterfaceElement::PhysicalEquation::E0

Definition at line 65 of file CohesiveInterfaceElement.hpp.

◆ ft

double CohesiveElement::CohesiveInterfaceElement::PhysicalEquation::ft

Definition at line 61 of file CohesiveInterfaceElement.hpp.

◆ g0

double CohesiveElement::CohesiveInterfaceElement::PhysicalEquation::g0

Definition at line 65 of file CohesiveInterfaceElement.hpp.

◆ Gf

double CohesiveElement::CohesiveInterfaceElement::PhysicalEquation::Gf

Definition at line 61 of file CohesiveInterfaceElement.hpp.

◆ h

double CohesiveElement::CohesiveInterfaceElement::PhysicalEquation::h

Definition at line 61 of file CohesiveInterfaceElement.hpp.

◆ isInitialised

bool CohesiveElement::CohesiveInterfaceElement::PhysicalEquation::isInitialised

Definition at line 54 of file CohesiveInterfaceElement.hpp.

◆ kappa1

double CohesiveElement::CohesiveInterfaceElement::PhysicalEquation::kappa1

Definition at line 65 of file CohesiveInterfaceElement.hpp.

◆ kappaPtr

double* CohesiveElement::CohesiveInterfaceElement::PhysicalEquation::kappaPtr

Definition at line 99 of file CohesiveInterfaceElement.hpp.

◆ kappaSize

int CohesiveElement::CohesiveInterfaceElement::PhysicalEquation::kappaSize

Definition at line 100 of file CohesiveInterfaceElement.hpp.

◆ mField

MoFEM::Interface& CohesiveElement::CohesiveInterfaceElement::PhysicalEquation::mField

Definition at line 53 of file CohesiveInterfaceElement.hpp.

◆ pRisms

Range CohesiveElement::CohesiveInterfaceElement::PhysicalEquation::pRisms

Definition at line 62 of file CohesiveInterfaceElement.hpp.

◆ thDamagedPrism

Tag CohesiveElement::CohesiveInterfaceElement::PhysicalEquation::thDamagedPrism

Definition at line 63 of file CohesiveInterfaceElement.hpp.

◆ thKappa

Tag CohesiveElement::CohesiveInterfaceElement::PhysicalEquation::thKappa

Definition at line 63 of file CohesiveInterfaceElement.hpp.

◆ youngModulus

double CohesiveElement::CohesiveInterfaceElement::PhysicalEquation::youngModulus

Definition at line 61 of file CohesiveInterfaceElement.hpp.


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