#include <users_modules/basic_finite_elements/cohesive_interface/src/InterfaceGapArcLengthControl.hpp>
◆ ArcLengthIntElemFEMethod()
CohesiveElement::ArcLengthIntElemFEMethod::ArcLengthIntElemFEMethod |
( |
moab::Interface & |
moab, |
|
|
boost::shared_ptr< ArcLengthCtx > |
arcptr |
|
) |
| |
|
inline |
Definition at line 24 of file InterfaceGapArcLengthControl.hpp.
28 PetscInt ghosts[1] = { 0 };
30 if(pcomm->rank() == 0) {
37 for(Range::iterator pit = prisms.begin();pit!=prisms.end();pit++) {
54 double def_damaged = 0;
◆ ~ArcLengthIntElemFEMethod()
CohesiveElement::ArcLengthIntElemFEMethod::~ArcLengthIntElemFEMethod |
( |
| ) |
|
|
inline |
◆ calculate_db()
virtual MoFEMErrorCode CohesiveElement::ArcLengthIntElemFEMethod::calculate_db |
( |
| ) |
|
|
inlinevirtual |
Definition at line 145 of file InterfaceGapArcLengthControl.hpp.
150 NumeredDofEntityByLocalIdx::iterator dit,hi_dit;
151 dit = problemPtr->getNumeredRowDofsPtr()->get<PetscLocalIdx_mi_tag>().lower_bound(0);
152 hi_dit = problemPtr->getNumeredRowDofsPtr()->get<PetscLocalIdx_mi_tag>().upper_bound(
153 problemPtr->getNbLocalDofsRow()+problemPtr->getNbGhostDofsRow()
157 for(;dit!=hi_dit;dit++) {
158 if(dit->get()->getEntType() != MBVERTEX) {
159 array[dit->get()->getPetscLocalDofIdx()] = 0;
163 array[dit->get()->getPetscLocalDofIdx()] = +
arcPtr->alpha;
164 }
else if(
Nodes4.find(dit->get()->getEnt())!=
Nodes4.end()) {
165 array[dit->get()->getPetscLocalDofIdx()] = -
arcPtr->alpha;
166 }
else array[dit->get()->getPetscLocalDofIdx()] = 0;
◆ calculate_dx_and_dlambda()
MoFEMErrorCode CohesiveElement::ArcLengthIntElemFEMethod::calculate_dx_and_dlambda |
( |
Vec & |
x | ) |
|
|
inline |
◆ calculate_init_dlambda()
MoFEMErrorCode CohesiveElement::ArcLengthIntElemFEMethod::calculate_init_dlambda |
( |
double * |
dlambda | ) |
|
|
inline |
Definition at line 236 of file InterfaceGapArcLengthControl.hpp.
242 "\tInit dlambda = %6.4e s = %6.4e beta = %6.4e F_lambda2 = %6.4e\n",
247 std::ostringstream sss;
249 SETERRQ(PETSC_COMM_SELF,1,sss.str().c_str());
◆ calculate_lambda_int()
MoFEMErrorCode CohesiveElement::ArcLengthIntElemFEMethod::calculate_lambda_int |
( |
double & |
_lambda_int_ | ) |
|
|
inline |
Definition at line 103 of file InterfaceGapArcLengthControl.hpp.
106 NumeredDofEntityByLocalIdx::iterator dit,hi_dit;
107 dit = problemPtr->getNumeredRowDofsPtr()->get<PetscLocalIdx_mi_tag>().lower_bound(0);
108 hi_dit = problemPtr->getNumeredRowDofsPtr()->get<PetscLocalIdx_mi_tag>().upper_bound(problemPtr->getNbLocalDofsRow());
110 double *array_int_lambda;
116 array_int_lambda[0] = 0;
117 for(;dit!=hi_dit;dit++) {
118 if(dit->get()->getEntType() != MBVERTEX)
continue;
119 if(pcomm->rank() != dit->get()->getPart())
continue;
121 array_int_lambda[0] += array[dit->get()->getPetscLocalDofIdx()];
124 array_int_lambda[0] -= array[dit->get()->getPetscLocalDofIdx()];
◆ operator()()
MoFEMErrorCode CohesiveElement::ArcLengthIntElemFEMethod::operator() |
( |
| ) |
|
|
inline |
Definition at line 172 of file InterfaceGapArcLengthControl.hpp.
176 case CTX_SNESSETFUNCTION: {
180 PetscPrintf(PETSC_COMM_SELF,
"\tres_lambda = %6.4e lambda_int = %6.4e s = %6.4e\n",
184 case CTX_SNESSETJACOBIAN: {
◆ postProcess()
MoFEMErrorCode CohesiveElement::ArcLengthIntElemFEMethod::postProcess |
( |
| ) |
|
|
inline |
◆ preProcess()
MoFEMErrorCode CohesiveElement::ArcLengthIntElemFEMethod::preProcess |
( |
| ) |
|
|
inline |
◆ remove_damaged_prisms_nodes()
MoFEMErrorCode CohesiveElement::ArcLengthIntElemFEMethod::remove_damaged_prisms_nodes |
( |
| ) |
|
|
inline |
remove nodes of prims which are fully damaged
Definition at line 66 of file InterfaceGapArcLengthControl.hpp.
70 std::vector<int> is_prism_damaged(prisms.size());
72 Range::iterator pit = prisms.begin();
73 std::vector<int>::iterator vit = is_prism_damaged.begin();
74 for(;pit!=prisms.end();pit++,vit++) {
78 for(Range::iterator nit = nodes.begin();nit!=nodes.end();nit++) {
◆ set_dlambda_to_x()
MoFEMErrorCode CohesiveElement::ArcLengthIntElemFEMethod::set_dlambda_to_x |
( |
Vec & |
x, |
|
|
double |
dlambda |
|
) |
| |
|
inline |
Definition at line 255 of file InterfaceGapArcLengthControl.hpp.
258 if(
arcPtr->getPetscLocalDofIdx()!=-1) {
261 double lambda_old = array[
arcPtr->getPetscLocalDofIdx()];
262 if(!(dlambda == dlambda)) {
263 std::ostringstream sss;
265 SETERRQ(PETSC_COMM_SELF,1,sss.str().c_str());
267 array[
arcPtr->getPetscLocalDofIdx()] = lambda_old + dlambda;
268 PetscPrintf(PETSC_COMM_WORLD,
"\tlambda = %6.4e, %6.4e (%6.4e)\n",
269 lambda_old, array[
arcPtr->getPetscLocalDofIdx()], dlambda);
◆ arcPtr
boost::shared_ptr<ArcLengthCtx> CohesiveElement::ArcLengthIntElemFEMethod::arcPtr |
◆ Edges3
Range CohesiveElement::ArcLengthIntElemFEMethod::Edges3 |
◆ Edges4
Range CohesiveElement::ArcLengthIntElemFEMethod::Edges4 |
◆ Faces3
Range CohesiveElement::ArcLengthIntElemFEMethod::Faces3 |
◆ Faces4
Range CohesiveElement::ArcLengthIntElemFEMethod::Faces4 |
◆ GhostLambdaInt
Vec CohesiveElement::ArcLengthIntElemFEMethod::GhostLambdaInt |
◆ lambda_int
double CohesiveElement::ArcLengthIntElemFEMethod::lambda_int |
◆ mOab
moab::Interface& CohesiveElement::ArcLengthIntElemFEMethod::mOab |
◆ Nodes3
Range CohesiveElement::ArcLengthIntElemFEMethod::Nodes3 |
◆ Nodes4
Range CohesiveElement::ArcLengthIntElemFEMethod::Nodes4 |
◆ thDamagedPrism
Tag CohesiveElement::ArcLengthIntElemFEMethod::thDamagedPrism |
The documentation for this struct was generated from the following file: