v0.14.0
HMHPMooneyRivlinWriggersEq63.cpp
Go to the documentation of this file.
1 /**
2  * @file HMHPMooneyRivlinWriggersEq63.cpp
3  * @brief Implement of MooneyRivlinWriggersEq63
4  *
5  * @copyright Copyright (c) 2024
6  *
7  */
8 
9 namespace EshelbianPlasticity {
10 
12 
13  static constexpr int numberOfActiveVariables = 9;
14  static constexpr int numberOfDependentVariables = 9;
15 
16  HMHPMooneyRivlinWriggersEq63(const double alpha, const double beta,
17  const double lambda)
20 
21  virtual OpJacobian *
22  returnOpJacobian(const int tag, const bool eval_rhs, const bool eval_lhs,
23  boost::shared_ptr<DataAtIntegrationPts> data_ptr,
24  boost::shared_ptr<PhysicalEquations> physics_ptr) {
25  return (new OpHMHH(tag, eval_rhs, eval_lhs, data_ptr, physics_ptr));
26  }
27 
30  CHKERR PetscOptionsBegin(PETSC_COMM_WORLD, "mooneyrivlin_", "", "none");
31 
32  alpha = 1;
33  CHKERR PetscOptionsScalar("-alpha", "Alpha", "", alpha, &alpha, PETSC_NULL);
34  beta = 1;
35  CHKERR PetscOptionsScalar("-beta", "Beta", "", beta, &beta, PETSC_NULL);
36 
37  lambda = 1;
38  CHKERR PetscOptionsScalar("-lambda", "Lambda", "", lambda, &lambda,
39  PETSC_NULL);
40 
41  epsilon = 0;
42  CHKERR PetscOptionsScalar("-epsilon", "Epsilon", "", epsilon, &epsilon,
43  PETSC_NULL);
44 
45  ierr = PetscOptionsEnd();
46  CHKERRG(ierr);
48  }
49 
50  double alpha;
51  double beta;
52  double lambda;
53  double epsilon;
54 
58 
63 
66 
71 
76 
78 
79  MoFEMErrorCode recordTape(const int tape, DTensor2Ptr *t_h_ptr) {
81 
91 
93 
95 
96  auto ih = get_h();
97  if (t_h_ptr)
98  ih(i, j) = (*t_h_ptr)(i, j);
99  else {
100  ih(i, j) = t_kd(i, j);
101  }
102 
103  auto r_P = get_P();
104 
105  enableMinMaxUsingAbs();
106  trace_on(tape);
107 
108  // Set active variables to ADOL-C
109  th(i, j) <<= get_h()(i, j);
110 
111  tF(i, I) = th(i, I);
114 
115  tCof(i, I) = detF * tInvF(I, i);
116 
117  A = tF(k, K) * tF(k, K);
118  B = tCof(k, K) * tCof(k, K);
119 
120  tBF(i, I) = 4 * alpha * (A * tF(i, I));
121  tBCof(i, I) = 4 * beta * (B * tCof(i, I));
122  tBj = (-12 * alpha - 24 * beta) / detF +
123  0.5 * (lambda / epsilon) *
124  (pow(detF, epsilon - 1) - pow(detF, -epsilon - 1));
125 
126  tP(i, I) = tBF(i, I);
127  tP(i, I) += (levi_civita(i, j, k) * tBCof(j, J)) *
128  (levi_civita(I, J, K) * tF(k, K));
129  tP(i, I) += tCof(i, I) * tBj;
130 
131  // Set dependent variables to ADOL-C
132  tP(i, j) >>= r_P(i, j);
133 
134  trace_off();
135 
137  }
138 };
139 
140 
141 
142 } // namespace EshelbianPlasticity
MoFEM::K
VectorDouble K
Definition: Projection10NodeCoordsOnField.cpp:125
EshelbianPlasticity::HMHPMooneyRivlinWriggersEq63::energy
adouble energy
Definition: HMHPMooneyRivlinWriggersEq63.cpp:72
EshelbianPlasticity::HMHPMooneyRivlinWriggersEq63::numberOfDependentVariables
static constexpr int numberOfDependentVariables
Definition: HMHPMooneyRivlinWriggersEq63.cpp:14
EshelbianPlasticity::HMHPMooneyRivlinWriggersEq63::tSigma
ATensor2 tSigma
Definition: HMHPMooneyRivlinWriggersEq63.cpp:65
EshelbianPlasticity::HMHPMooneyRivlinWriggersEq63::alpha
double alpha
Definition: HMHPMooneyRivlinWriggersEq63.cpp:50
MoFEM::Exceptions::MoFEMErrorCode
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
EshelbianPlasticity::HMHPMooneyRivlinWriggersEq63::tBCof
ATensor2 tBCof
Definition: HMHPMooneyRivlinWriggersEq63.cpp:69
EshelbianPlasticity::HMHPMooneyRivlinWriggersEq63::tP
ATensor2 tP
Definition: HMHPMooneyRivlinWriggersEq63.cpp:64
EshelbianPlasticity
Definition: CGGTonsorialBubbleBase.hpp:11
EshelbianPlasticity::HMHPMooneyRivlinWriggersEq63::tInvH
ATensor2 tInvH
Definition: HMHPMooneyRivlinWriggersEq63.cpp:61
EshelbianPlasticity::PhysicalEquations
Definition: EshelbianPlasticity.hpp:219
J
FTensor::Index< 'J', DIM1 > J
Definition: level_set.cpp:30
EshelbianPlasticity::HMHPMooneyRivlinWriggersEq63::recordTape
MoFEMErrorCode recordTape(const int tape, DTensor2Ptr *t_h_ptr)
Definition: HMHPMooneyRivlinWriggersEq63.cpp:79
FTensor::levi_civita
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
Definition: Levi_Civita.hpp:617
EshelbianPlasticity::HMHPMooneyRivlinWriggersEq63::beta
double beta
Definition: HMHPMooneyRivlinWriggersEq63.cpp:51
EshelbianPlasticity::HMHPMooneyRivlinWriggersEq63::tInvF
ATensor2 tInvF
Definition: HMHPMooneyRivlinWriggersEq63.cpp:62
FTensor::Tensor2< adouble, 3, 3 >
I
constexpr IntegrationType I
Definition: operators_tests.cpp:31
MoFEM::invertTensor3by3
MoFEMErrorCode invertTensor3by3(ublas::matrix< T, L, A > &jac_data, ublas::vector< T, A > &det_data, ublas::matrix< T, L, A > &inv_jac_data)
Calculate inverse of tensor rank 2 at integration points.
Definition: Templates.hpp:1588
FTensor::Number< 0 >
CHKERR
#define CHKERR
Inline error check.
Definition: definitions.h:548
EshelbianPlasticity::HMHPMooneyRivlinWriggersEq63::HMHPMooneyRivlinWriggersEq63
HMHPMooneyRivlinWriggersEq63(const double alpha, const double beta, const double lambda)
Definition: HMHPMooneyRivlinWriggersEq63.cpp:16
EshelbianPlasticity::HMHPMooneyRivlinWriggersEq63::th
ATensor2 th
Definition: HMHPMooneyRivlinWriggersEq63.cpp:55
EshelbianPlasticity::HMHPMooneyRivlinWriggersEq63::tH
ATensor2 tH
Definition: HMHPMooneyRivlinWriggersEq63.cpp:56
EshelbianPlasticity::HMHPMooneyRivlinWriggersEq63::tBF
ATensor2 tBF
Definition: HMHPMooneyRivlinWriggersEq63.cpp:68
EshelbianPlasticity::HMHPMooneyRivlinWriggersEq63::getOptions
MoFEMErrorCode getOptions()
Definition: HMHPMooneyRivlinWriggersEq63.cpp:28
EshelbianPlasticity::HMHPMooneyRivlinWriggersEq63::A
adouble A
Definition: HMHPMooneyRivlinWriggersEq63.cpp:74
EshelbianPlasticity::OpHMHH
Definition: EshelbianADOL-C.cpp:277
EshelbianPlasticity::HMHPMooneyRivlinWriggersEq63::detH
adouble detH
Definition: HMHPMooneyRivlinWriggersEq63.cpp:59
i
FTensor::Index< 'i', SPACE_DIM > i
Definition: hcurl_divergence_operator_2d.cpp:27
t_kd
constexpr auto t_kd
Definition: free_surface.cpp:137
FTensor::Index< 'i', 3 >
MoFEM::determinantTensor3by3
static auto determinantTensor3by3(T &t)
Calculate the determinant of a 3x3 matrix or a tensor of rank 2.
Definition: Templates.hpp:1540
adouble
EshelbianPlasticity::HMHPMooneyRivlinWriggersEq63::lambda
double lambda
Definition: HMHPMooneyRivlinWriggersEq63.cpp:52
EshelbianPlasticity::OpJacobian
Definition: EshelbianPlasticity.hpp:380
j
FTensor::Index< 'j', 3 > j
Definition: matrix_function.cpp:19
EshelbianPlasticity::HMHPMooneyRivlinWriggersEq63::returnOpJacobian
virtual OpJacobian * returnOpJacobian(const int tag, const bool eval_rhs, const bool eval_lhs, boost::shared_ptr< DataAtIntegrationPts > data_ptr, boost::shared_ptr< PhysicalEquations > physics_ptr)
Definition: HMHPMooneyRivlinWriggersEq63.cpp:22
MoFEM::Exceptions::ierr
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Exceptions.hpp:76
EshelbianPlasticity::HMHPMooneyRivlinWriggersEq63::detF
adouble detF
Definition: HMHPMooneyRivlinWriggersEq63.cpp:60
EshelbianPlasticity::HMHPMooneyRivlinWriggersEq63::tF
ATensor2 tF
Definition: HMHPMooneyRivlinWriggersEq63.cpp:57
EshelbianPlasticity::HMHPMooneyRivlinWriggersEq63::B
adouble B
Definition: HMHPMooneyRivlinWriggersEq63.cpp:75
FTensor::Kronecker_Delta_symmetric
Kronecker Delta class symmetric.
Definition: Kronecker_Delta.hpp:49
EshelbianPlasticity::HMHPMooneyRivlinWriggersEq63::tCof
ATensor2 tCof
Definition: HMHPMooneyRivlinWriggersEq63.cpp:67
EshelbianPlasticity::HMHPMooneyRivlinWriggersEq63::tBj
adouble tBj
Definition: HMHPMooneyRivlinWriggersEq63.cpp:70
EshelbianPlasticity::PhysicalEquations::get_h
DTensor2Ptr get_h()
Definition: EshelbianPlasticity.hpp:317
k
FTensor::Index< 'k', 3 > k
Definition: matrix_function.cpp:20
EshelbianPlasticity::HMHPMooneyRivlinWriggersEq63::phi
adouble phi
Definition: HMHPMooneyRivlinWriggersEq63.cpp:73
MoFEMFunctionReturn
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
EshelbianPlasticity::HMHPMooneyRivlinWriggersEq63
Definition: HMHPMooneyRivlinWriggersEq63.cpp:11
CHKERRG
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:496
EshelbianPlasticity::HMHPMooneyRivlinWriggersEq63::tPulledP
ATensor2 tPulledP
Definition: HMHPMooneyRivlinWriggersEq63.cpp:77
MoFEMFunctionBegin
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
EshelbianPlasticity::PhysicalEquations::get_P
DTensor2Ptr get_P()
Definition: EshelbianPlasticity.hpp:325
EshelbianPlasticity::HMHPMooneyRivlinWriggersEq63::epsilon
double epsilon
Definition: HMHPMooneyRivlinWriggersEq63.cpp:53
EshelbianPlasticity::HMHPMooneyRivlinWriggersEq63::numberOfActiveVariables
static constexpr int numberOfActiveVariables
Definition: HMHPMooneyRivlinWriggersEq63.cpp:13