v0.8.23
Public Member Functions | List of all members
EshelbianPlasticity::OpHMHH Struct Reference
Inheritance diagram for EshelbianPlasticity::OpHMHH:
[legend]
Collaboration diagram for EshelbianPlasticity::OpHMHH:
[legend]

Public Member Functions

 OpHMHH (const std::string &field_name, const int tag, const bool eval_rhs, const bool eval_lhs, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, boost::shared_ptr< PhysicalEquations > &physics_ptr)
 
MoFEMErrorCode evaluateRhs (EntData &data)
 
MoFEMErrorCode evaluateLhs (EntData &data)
 
- Public Member Functions inherited from EshelbianPlasticity::OpJacobian
 OpJacobian (const std::string &field_name, const int tag, const bool eval_rhs, const bool eval_lhs, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, boost::shared_ptr< PhysicalEquations > &physics_ptr)
 
MoFEMErrorCode doWork (int side, EntityType type, EntData &data)
 

Additional Inherited Members

- Public Attributes inherited from EshelbianPlasticity::OpJacobian
const int tAg
 adol-c tape More...
 
const bool evalRhs
 
const bool evalLhs
 
boost::shared_ptr< DataAtIntegrationPtsdataAtPts
 data at integration pts More...
 
boost::shared_ptr< PhysicalEquationsphysicsPtr
 material physical equations More...
 

Detailed Description

Definition at line 18 of file EshelbianADOL-C.cpp.

Constructor & Destructor Documentation

◆ OpHMHH()

EshelbianPlasticity::OpHMHH::OpHMHH ( const std::string &  field_name,
const int  tag,
const bool  eval_rhs,
const bool  eval_lhs,
boost::shared_ptr< DataAtIntegrationPts > &  data_ptr,
boost::shared_ptr< PhysicalEquations > &  physics_ptr 
)

Definition at line 19 of file EshelbianADOL-C.cpp.

22  : OpJacobian(field_name, tag, eval_rhs, eval_lhs, data_ptr, physics_ptr) {
23  }
OpJacobian(const std::string &field_name, const int tag, const bool eval_rhs, const bool eval_lhs, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, boost::shared_ptr< PhysicalEquations > &physics_ptr)

Member Function Documentation

◆ evaluateLhs()

MoFEMErrorCode EshelbianPlasticity::OpHMHH::evaluateLhs ( EntData data)
virtual

Implements EshelbianPlasticity::OpJacobian.

Definition at line 98 of file EshelbianADOL-C.cpp.

98  {
103  const int number_of_active_variables = physicsPtr->activeVariables.size();
104  const int number_of_dependent_variables =
105  physicsPtr->dependentVariables.size();
106  std::vector<double *> jac_ptr(number_of_dependent_variables);
107  for (unsigned int n = 0; n != number_of_dependent_variables; ++n) {
108  jac_ptr[n] =
109  &(physicsPtr
110  ->dependentVariablesDirevatives[n * number_of_active_variables]);
111  }
112 
113  const auto nb_integration_pts = ent_data.getN().size1();
114 
115  auto create_data_mat =
116  [nb_integration_pts](boost::shared_ptr<MatrixDouble> &m_ptr) {
117  if (!m_ptr)
118  m_ptr = boost::make_shared<MatrixDouble>();
119  m_ptr->resize(9, nb_integration_pts, false);
120  };
121 
122  auto create_data_ten =
123  [nb_integration_pts](boost::shared_ptr<MatrixDouble> &m_ptr) {
124  if (!m_ptr)
125  m_ptr = boost::make_shared<MatrixDouble>();
126  m_ptr->resize(27, nb_integration_pts, false);
127  };
128 
129  create_data_ten(dataAtPts->P_dh0);
130  create_data_ten(dataAtPts->P_dh1);
131  create_data_ten(dataAtPts->P_dh2);
132  create_data_ten(dataAtPts->P_dH0);
133  create_data_ten(dataAtPts->P_dH1);
134  create_data_ten(dataAtPts->P_dH2);
135  create_data_ten(dataAtPts->Sigma_dh0);
136  create_data_ten(dataAtPts->Sigma_dh1);
137  create_data_ten(dataAtPts->Sigma_dh2);
138  create_data_ten(dataAtPts->Sigma_dH0);
139  create_data_ten(dataAtPts->Sigma_dH1);
140  create_data_ten(dataAtPts->Sigma_dH2);
141  create_data_mat(dataAtPts->phi_dh);
142  create_data_mat(dataAtPts->phi_dH);
143  create_data_ten(dataAtPts->Flow_dh0);
144  create_data_ten(dataAtPts->Flow_dh1);
145  create_data_ten(dataAtPts->Flow_dh2);
146  create_data_ten(dataAtPts->Flow_dH0);
147  create_data_ten(dataAtPts->Flow_dH1);
148  create_data_ten(dataAtPts->Flow_dH2);
149 
150  auto ix_grad = getFTensor2FromMat<3, 3>(*dataAtPts->xGradAtPts);
151  auto iu = getFTensor2SymmetricFromMat<3>(*dataAtPts->streachTensorAtPts);
152  auto iG = getFTensor2FromMat<3, 3>(*(dataAtPts->GAtPts));
153 
154  auto r_P_dh0 = getFTensor3FromMat(*(dataAtPts->P_dh0));
155  auto r_P_dh1 = getFTensor3FromMat(*(dataAtPts->P_dh1));
156  auto r_P_dh2 = getFTensor3FromMat(*(dataAtPts->P_dh2));
157  auto r_P_dH0 = getFTensor3FromMat(*(dataAtPts->P_dH0));
158  auto r_P_dH1 = getFTensor3FromMat(*(dataAtPts->P_dH1));
159  auto r_P_dH2 = getFTensor3FromMat(*(dataAtPts->P_dH2));
160  auto r_Sigma_dh0 = getFTensor3FromMat(*(dataAtPts->Sigma_dh0));
161  auto r_Sigma_dh1 = getFTensor3FromMat(*(dataAtPts->Sigma_dh1));
162  auto r_Sigma_dh2 = getFTensor3FromMat(*(dataAtPts->Sigma_dh2));
163  auto r_Sigma_dH0 = getFTensor3FromMat(*(dataAtPts->Sigma_dH0));
164  auto r_Sigma_dH1 = getFTensor3FromMat(*(dataAtPts->Sigma_dH1));
165  auto r_Sigma_dH2 = getFTensor3FromMat(*(dataAtPts->Sigma_dH2));
166  auto r_phi_dh = getFTensor2FromMat<3, 3>(*(dataAtPts->phi_dh));
167  auto r_phi_dH = getFTensor2FromMat<3, 3>(*(dataAtPts->phi_dH));
168  auto r_Flow_dh0 = getFTensor3FromMat(*(dataAtPts->Flow_dh0));
169  auto r_Flow_dh1 = getFTensor3FromMat(*(dataAtPts->Flow_dh1));
170  auto r_Flow_dh2 = getFTensor3FromMat(*(dataAtPts->Flow_dh2));
171  auto r_Flow_dH0 = getFTensor3FromMat(*(dataAtPts->Flow_dH0));
172  auto r_Flow_dH1 = getFTensor3FromMat(*(dataAtPts->Flow_dH1));
173  auto r_Flow_dH2 = getFTensor3FromMat(*(dataAtPts->Flow_dH2));
174 
175  for (unsigned int gg = 0; gg != nb_integration_pts; ++gg) {
176 
177  auto t_h = physicsPtr->get_h();
178  t_h(i, j) = iu(i, k) * ix_grad(k, j) + ix_grad(i, j);
179 
180  physicsPtr->get_H()(i, j) = iG(i, j);
181  for (int ii = 0; ii != 3; ++ii)
182  physicsPtr->get_H()(ii, ii) += 1;
183 
184  // play recorder for jacobians
185  int r = ::jacobian(tAg, number_of_dependent_variables,
186  number_of_active_variables,
187  &physicsPtr->activeVariables[0], &jac_ptr[0]);
188  if (r < 0) {
189  SETERRQ(PETSC_COMM_SELF, MOFEM_OPERATION_UNSUCCESSFUL,
190  "ADOL-C function evaluation with error");
191  }
192 
193  r_P_dh0(i, j, k) = physicsPtr->get_P_dh0()(i, j, k);
194  r_P_dh1(i, j, k) = physicsPtr->get_P_dh1()(i, j, k);
195  r_P_dh2(i, j, k) = physicsPtr->get_P_dh2()(i, j, k);
196  r_P_dH0(i, j, k) = physicsPtr->get_P_dH0()(i, j, k);
197  r_P_dH1(i, j, k) = physicsPtr->get_P_dH1()(i, j, k);
198  r_P_dH2(i, j, k) = physicsPtr->get_P_dH2()(i, j, k);
199  r_Sigma_dh0(i, j, k) = physicsPtr->get_Sigma_dh0()(i, j, k);
200  r_Sigma_dh1(i, j, k) = physicsPtr->get_Sigma_dh1()(i, j, k);
201  r_Sigma_dh2(i, j, k) = physicsPtr->get_Sigma_dh2()(i, j, k);
202  r_Sigma_dH0(i, j, k) = physicsPtr->get_Sigma_dH0()(i, j, k);
203  r_Sigma_dH1(i, j, k) = physicsPtr->get_Sigma_dH1()(i, j, k);
204  r_Sigma_dH2(i, j, k) = physicsPtr->get_Sigma_dH2()(i, j, k);
205  r_phi_dh(i, j) = physicsPtr->get_Phi_dh()(i, j);
206  r_phi_dH(i, j) = physicsPtr->get_Phi_dH()(i, j);
207  r_Flow_dh0(i, j, k) = physicsPtr->get_Flow_dh0()(i, j, k);
208  r_Flow_dh1(i, j, k) = physicsPtr->get_Flow_dh1()(i, j, k);
209  r_Flow_dh2(i, j, k) = physicsPtr->get_Flow_dh2()(i, j, k);
210  r_Flow_dH0(i, j, k) = physicsPtr->get_Flow_dH0()(i, j, k);
211  r_Flow_dH1(i, j, k) = physicsPtr->get_Flow_dH1()(i, j, k);
212  r_Flow_dH2(i, j, k) = physicsPtr->get_Flow_dH2()(i, j, k);
213 
214  ++iu;
215  ++ix_grad;
216  ++iG;
217  ++r_P_dh0;
218  ++r_P_dh1;
219  ++r_P_dh2;
220  ++r_P_dH0;
221  ++r_P_dH1;
222  ++r_P_dH2;
223  ++r_Sigma_dh0;
224  ++r_Sigma_dh1;
225  ++r_Sigma_dh2;
226  ++r_Sigma_dH0;
227  ++r_Sigma_dH1;
228  ++r_Sigma_dH2;
229  ++r_phi_dh;
230  ++r_phi_dH;
231  ++r_Flow_dh0;
232  ++r_Flow_dh1;
233  ++r_Flow_dh2;
234  ++r_Flow_dH0;
235  ++r_Flow_dH1;
236  ++r_Flow_dH2;
237  }
239 }
boost::shared_ptr< PhysicalEquations > physicsPtr
material physical equations
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:476
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
data at integration pts
FTensor::Tensor3< FTensor::PackPtr< double *, 1 >, 3, 3, 3 > getFTensor3FromMat(MatrixDouble &m)
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:406

◆ evaluateRhs()

MoFEMErrorCode EshelbianPlasticity::OpHMHH::evaluateRhs ( EntData data)
virtual

Implements EshelbianPlasticity::OpJacobian.

Definition at line 29 of file EshelbianADOL-C.cpp.

29  {
34 
35  const auto nb_integration_pts = ent_data.getN().size1();
36  auto ix_grad = getFTensor2FromMat<3, 3>(*dataAtPts->xGradAtPts);
37  auto iu = getFTensor2SymmetricFromMat<3>(*dataAtPts->streachTensorAtPts);
38  // FIXME: that should work with material streach
39  auto iG = getFTensor2FromMat<3, 3>(*dataAtPts->GAtPts);
40 
41  auto create_data_vec =
42  [nb_integration_pts](boost::shared_ptr<VectorDouble> &v_ptr) {
43  if (!v_ptr)
44  v_ptr = boost::make_shared<VectorDouble>();
45  v_ptr->resize(nb_integration_pts, false);
46  };
47  auto create_data_mat =
48  [nb_integration_pts](boost::shared_ptr<MatrixDouble> &m_ptr) {
49  if (!m_ptr)
50  m_ptr = boost::make_shared<MatrixDouble>();
51  m_ptr->resize(9, nb_integration_pts, false);
52  };
53 
54  create_data_mat(dataAtPts->PAtPts);
55  create_data_mat(dataAtPts->SigmaAtPts);
56  create_data_vec(dataAtPts->phiAtPts);
57  create_data_mat(dataAtPts->flowAtPts);
58 
59  auto r_P = getFTensor2FromMat<3, 3>(*(dataAtPts->PAtPts));
60  auto r_Sigma = getFTensor2FromMat<3, 3>(*(dataAtPts->SigmaAtPts));
61  auto r_phi = getFTensor0FromVec(*(dataAtPts->phiAtPts));
62  auto r_flow = getFTensor2FromMat<3, 3>(*(dataAtPts->flowAtPts));
63 
64  for (unsigned int gg = 0; gg != nb_integration_pts; ++gg) {
65 
66  auto t_h = physicsPtr->get_h();
67  t_h(i, j) = iu(i, k) * ix_grad(k, j) + ix_grad(i, j);
68 
69  physicsPtr->get_H()(i, j) = iG(i, j);
70  for (int ii = 0; ii != 3; ++ii)
71  physicsPtr->get_H()(ii, ii) += 1;
72 
73  int r = ::function(tAg, physicsPtr->dependentVariables.size(),
74  physicsPtr->activeVariables.size(),
75  &physicsPtr->activeVariables[0],
76  &physicsPtr->dependentVariables[0]);
77  if (r < 0) { // function is locally analytic
78  SETERRQ1(PETSC_COMM_SELF, MOFEM_OPERATION_UNSUCCESSFUL,
79  "ADOL-C function evaluation with error r = %d", r);
80  }
81 
82  r_P(i, j) = physicsPtr->get_P()(i, j);
83  r_Sigma(i, j) = physicsPtr->get_Sigma()(i, j);
84  r_phi = physicsPtr->get_Phi();
85  r_flow(i, j) = physicsPtr->get_Flow()(i, j);
86 
87  ++iu;
88  ++ix_grad;
89  ++iG;
90  ++r_P;
91  ++r_Sigma;
92  ++r_phi;
93  ++r_flow;
94  }
96 }
boost::shared_ptr< PhysicalEquations > physicsPtr
material physical equations
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:476
static FTensor::Tensor0< FTensor::PackPtr< double *, 1 > > getFTensor0FromVec(ublas::vector< T, A > &data)
Get tensor rank 0 (scalar) form data vectorExample how to use it.
Definition: Templates.hpp:143
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
data at integration pts
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:406

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