v0.15.0
Loading...
Searching...
No Matches
Public Types | Public Member Functions | Private Attributes | List of all members
OpFaceSideMaterialForce Struct Reference

#include "users_modules/eshelbian_plasticity/src/EshelbianOperators.hpp"

Inheritance diagram for OpFaceSideMaterialForce:
[legend]
Collaboration diagram for OpFaceSideMaterialForce:
[legend]

Public Types

using OP = VolumeElementForcesAndSourcesCoreOnSide::UserDataOperator
 

Public Member Functions

 OpFaceSideMaterialForce (boost::shared_ptr< DataAtIntegrationPts > data_ptr)
 
MoFEMErrorCode doWork (int side, EntityType type, EntData &data)
 Caluclate face material force and normal pressure at gauss points.
 

Private Attributes

boost::shared_ptr< DataAtIntegrationPtsdataAtPts
 data at integration pts
 

Detailed Description

Examples
mofem/users_modules/eshelbian_plasticity/src/impl/EshelbianPlasticity.cpp.

Definition at line 921 of file EshelbianOperators.hpp.

Member Typedef Documentation

◆ OP

using OpFaceSideMaterialForce::OP = VolumeElementForcesAndSourcesCoreOnSide::UserDataOperator

Definition at line 925 of file EshelbianOperators.hpp.

Constructor & Destructor Documentation

◆ OpFaceSideMaterialForce()

OpFaceSideMaterialForce::OpFaceSideMaterialForce ( boost::shared_ptr< DataAtIntegrationPts data_ptr)
inline

Definition at line 927 of file EshelbianOperators.hpp.

928 : OP(NOSPACE, OPSPACE), dataAtPts(data_ptr) {}
@ NOSPACE
Definition definitions.h:83
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
data at integration pts
VolumeElementForcesAndSourcesCoreOnSide::UserDataOperator OP

Member Function Documentation

◆ doWork()

MoFEMErrorCode OpFaceSideMaterialForce::doWork ( int  side,
EntityType  type,
EntData data 
)

Caluclate face material force and normal pressure at gauss points.

Parameters
side
type
data
Returns
MoFEMErrorCode

the full gradient \f$U=\nabla u\f$ on a surface from the symmetric part and the surface gradient.

Inputs:

  • t_strain : \(\varepsilon=\tfrac12(U+U^\top)\) (symmetric strain on S),
  • t_grad_u_gamma : \(u^\Gamma = U P\) (right-projected/surface gradient), with \(P=I-\mathbf N\otimes\mathbf N\),
  • t_normal : (possibly non‑unit) surface normal.

Procedure (pointwise on S): 1) Normalize the normal \(\mathbf n=\mathbf N/\|\mathbf N\|\). 2) Form the residual \(R=\varepsilon-\operatorname{sym}(u^\Gamma)\), where \(\operatorname{sym}(A)=\tfrac12(A+A^\top)\). 3) Recover the normal directional derivative (a vector) \(\mathbf v=\partial_{\mathbf n}u=2R\mathbf n-(\mathbf n^\top R\,\mathbf n)\,\mathbf n\). 4) Assemble the full gradient \(U = u^\Gamma + \mathbf v\otimes \mathbf n\).

Properties (sanity checks):

  • \(\tfrac12(U+U^\top)=\varepsilon\) (matches the given symmetric part),
  • \(U P = u^\Gamma\) (tangential/right-projected columns unchanged),
  • Only the normal column is updated via \(\mathbf v\otimes\mathbf n\).

Mapping to variables in this snippet:

  • \(\varepsilon \leftrightarrow\) t_strain,
  • \(u^\Gamma \leftrightarrow\) t_grad_u_gamma,
  • \(\mathbf N \leftrightarrow\) t_normal (normalized into t_N),
  • \(R \leftrightarrow\) t_R,
  • \(U \leftrightarrow\) t_grad_u.
Precondition
t_normal is nonzero; t_strain is symmetric.
Note
All indices use Einstein summation; computation is local to the surface point.
Examples
mofem/users_modules/eshelbian_plasticity/src/impl/EshelbianOperators.cpp.

Definition at line 3367 of file EshelbianOperators.cpp.

3368 {
3370
3383
3384 dataAtPts->faceMaterialForceAtPts.resize(3, getGaussPts().size2(), false);
3385 dataAtPts->normalPressureAtPts.resize(getGaussPts().size2(), false);
3386 if (getNinTheLoop() == 0) {
3387 dataAtPts->faceMaterialForceAtPts.clear();
3388 dataAtPts->normalPressureAtPts.clear();
3389 }
3390 auto loop_size = getLoopSize();
3391 if (loop_size == 1) {
3392 auto numebered_fe_ptr = getSidePtrFE()->numeredEntFiniteElementPtr;
3393 auto pstatus = numebered_fe_ptr->getPStatus();
3394 if (pstatus & (PSTATUS_SHARED | PSTATUS_MULTISHARED)) {
3395 loop_size = 2;
3396 }
3397 }
3398
3400
3401 auto t_normal = getFTensor1NormalsAtGaussPts();
3402 auto t_T = getFTensor1FromMat<SPACE_DIM>(
3403 dataAtPts->faceMaterialForceAtPts); //< face material force
3404 auto t_p =
3405 getFTensor0FromVec(dataAtPts->normalPressureAtPts); //< normal pressure
3406 auto t_P = getFTensor2FromMat<SPACE_DIM, SPACE_DIM>(dataAtPts->approxPAtPts);
3407 // auto t_grad_P = getFTensor3FromMat<SPACE_DIM, SPACE_DIM, SPACE_DIM>(
3408 // dataAtPts->gradPAtPts);
3409 auto t_u_gamma =
3410 getFTensor1FromMat<SPACE_DIM>(dataAtPts->hybridDispAtPts);
3411 auto t_grad_u_gamma =
3412 getFTensor2FromMat<SPACE_DIM, SPACE_DIM>(dataAtPts->gradHybridDispAtPts);
3413 auto t_strain =
3414 getFTensor2SymmetricFromMat<SPACE_DIM>(dataAtPts->logStretchTensorAtPts);
3415 auto t_omega = getFTensor1FromMat<3>(dataAtPts->rotAxisAtPts);
3416
3422
3423 auto next = [&]() {
3424 ++t_normal;
3425 ++t_P;
3426 // ++t_grad_P;
3427 ++t_omega;
3428 ++t_u_gamma;
3429 ++t_grad_u_gamma;
3430 ++t_strain;
3431 ++t_T;
3432 ++t_p;
3433 };
3434
3436 case GRIFFITH_FORCE:
3437 for (auto gg = 0; gg != getGaussPts().size2(); ++gg) {
3438 t_N(I) = t_normal(I);
3439 t_N.normalize();
3440
3441 t_A(i, j) = levi_civita(i, j, k) * t_omega(k);
3442 t_R(i, k) = t_kd(i, k) + t_A(i, k);
3443 t_grad_u(i, j) = t_R(i, j) + t_strain(i, j);
3444
3445 t_T(I) +=
3446 t_N(J) * (t_grad_u(i, I) * t_P(i, J)) / loop_size;
3447 // note that works only for Hooke material, for nonlinear material we need
3448 // strain energy expressed by stress
3449 t_T(I) -= t_N(I) * ((t_strain(i, K) * t_P(i, K)) / 2.) / loop_size;
3450
3451 t_p += t_N(I) * (t_N(J) * (t_grad_u_gamma(i, I) * t_P(i, J))) / loop_size;
3452
3453 next();
3454 }
3455 break;
3456 case GRIFFITH_SKELETON:
3457 for (auto gg = 0; gg != getGaussPts().size2(); ++gg) {
3458
3459 // Normalize the normal
3460 t_N(I) = t_normal(I);
3461 t_N.normalize();
3462
3463 // R = ε − sym(u^Γ)
3464 t_R(i, j) =
3465 t_strain(i, j) - 0.5 * (t_grad_u_gamma(i, j) + t_grad_u_gamma(j, i));
3466
3467 // U = u^Γ + [2 R N − (Nᵀ R N) N] ⊗ N
3468 t_grad_u(i, J) =
3469 t_grad_u_gamma(i, J) +
3470 (2 * t_R(i, K) * t_N(K) - (t_R(k, L) * t_N(k) * t_N(L)) * t_N(i)) *
3471 t_N(J);
3472
3473 t_T(I) += t_N(J) * (t_grad_u(i, I) * t_P(i, J)) / loop_size;
3474 // note that works only for Hooke material, for nonlinear material we need
3475 // strain energy expressed by stress
3476 t_T(I) -= t_N(I) * ((t_strain(i, K) * t_P(i, K)) / 2.) / loop_size;
3477
3478 // calculate nominal face pressure
3479 t_p += t_N(I) * (t_N(J) * (t_grad_u_gamma(i, I) * t_P(i, J))) / loop_size;
3480
3481 next();
3482 }
3483 break;
3484
3485 default:
3486 SETERRQ(PETSC_COMM_WORLD, MOFEM_NOT_IMPLEMENTED,
3487 "Grffith energy release "
3488 "selector not implemented");
3489 };
3490
3491#ifndef NDEBUG
3492 auto side_fe_ptr = getSidePtrFE();
3493 auto side_fe_mi_ptr = side_fe_ptr->numeredEntFiniteElementPtr;
3494 auto pstatus = side_fe_mi_ptr->getPStatus();
3495 if (pstatus) {
3496 auto owner = side_fe_mi_ptr->getOwnerProc();
3497 MOFEM_LOG("SELF", Sev::noisy)
3498 << "OpFaceSideMaterialForce: owner proc is not 0, owner proc: " << owner
3499 << " " << getPtrFE()->mField.get_comm_rank() << " n in the loop "
3500 << getNinTheLoop() << " loop size " << getLoopSize();
3501 }
3502#endif // NDEBUG
3503
3505}
#define FTENSOR_INDEX(DIM, I)
constexpr int SPACE_DIM
Kronecker Delta class symmetric.
Tensor1< T, Tensor_Dim > normalize()
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
@ MOFEM_NOT_IMPLEMENTED
Definition definitions.h:32
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
constexpr auto t_kd
#define MOFEM_LOG(channel, severity)
Log.
FTensor::Index< 'i', SPACE_DIM > i
const double n
refractive index of diffusive medium
FTensor::Index< 'J', DIM1 > J
Definition level_set.cpp:30
FTensor::Index< 'l', 3 > l
FTensor::Index< 'j', 3 > j
FTensor::Index< 'k', 3 > k
constexpr std::enable_if<(Dim0<=2 &&Dim1<=2), Tensor2_Expr< Levi_Civita< T >, T, Dim0, Dim1, i, j > >::type levi_civita(const Index< i, Dim0 > &, const Index< j, Dim1 > &)
levi_civita functions to make for easy adhoc use
static auto getFTensor0FromVec(ublas::vector< T, A > &data)
Get tensor rank 0 (scalar) form data vector.
constexpr IntegrationType I
FTensor::Index< 'm', 3 > m
const int N
Definition speed_test.cpp:3
static enum EnergyReleaseSelector energyReleaseSelector

Member Data Documentation

◆ dataAtPts

boost::shared_ptr<DataAtIntegrationPts> OpFaceSideMaterialForce::dataAtPts
private

data at integration pts

Definition at line 934 of file EshelbianOperators.hpp.


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