75 {
77
78 int nb_integration_pts = OP::getGaussPts().size2();
79
80 auto get_tag = [&]() {
82 for (
const auto &tag_range_pair : *(
physicalPtr->tagVsRangePtr)) {
83 if (tag_range_pair.second.find(DomainEleOp::getFEEntityHandle()) !=
84 tag_range_pair.second.end()) {
85 return tag_range_pair.first;
86 }
87 }
88 }
89#ifndef NDEBUG
92 "ADOL-C tag not found " +
94 }
95#endif
97 };
98
99 const int current_tag = get_tag();
100 auto *fe_ptr =
const_cast<FEMethod *
>(this->getFEMethod());
102
103 physicalPtr->matOpsDataPtr->getActiveDataPtr(
"F")->resize(DIM, DIM,
false);
104 physicalPtr->matOpsDataPtr->getDependentDataPtr(
"P")->resize(DIM, DIM,
false);
105 physicalPtr->matOpsDataPtr->getDependentDerivativesDataPtr(
"P_dF")->resize(
106 DIM * DIM, DIM * DIM, false);
107 auto mat_grad_ptr =
physicalPtr->matOpsDataPtr->getCommonDataPtr(
"grad");
108 auto mat_P_ptr =
physicalPtr->matOpsDataPtr->getCommonDataPtr(
"P");
109 auto mat_P_dF_ptr =
physicalPtr->matOpsDataPtr->getCommonDataPtr(
"P_dF");
110
111#ifndef NDEBUG
112 if (!mat_grad_ptr || !mat_P_ptr || !mat_P_dF_ptr) {
114 "Missing common data for ADOL-C evaluation");
115 }
116 if (mat_grad_ptr->size2() != DIM * DIM) {
118 "Inconsistent size of gradient matrix for ADOL-C evaluation");
119 }
120 if (mat_grad_ptr->size1() != nb_integration_pts) {
122 "Inconsistent size of gradient matrix data for ADOL-C evaluation "
123 "%zu != %d",
124 mat_grad_ptr->size1(), nb_integration_pts);
125 }
126#endif
127
129 auto get_grad_at_pts =
131 *mat_grad_ptr, nb_integration_pts);
132 auto get_P_at_pts =
134 *mat_P_ptr, nb_integration_pts);
135 auto get_P_dF_at_pts =
137 DL>::size(*mat_P_dF_ptr, nb_integration_pts);
138
142
143 auto t_grad_at_pts = get_grad_at_pts();
144 auto t_P_at_pts = get_P_at_pts();
145
146 auto next = [&]() {
147 ++t_grad_at_pts;
148 ++t_P_at_pts;
149 };
150
152
153 for (auto gg = 0; gg != nb_integration_pts; ++gg) {
154 auto t_F = getFTensor2FromPtr<DIM, DIM>(
155 physicalPtr->matOpsDataPtr->getActiveDataPtr(
"F")->data().data());
156 t_F(
i,
J) = t_grad_at_pts(
i,
J);
159 auto t_P = getFTensor2FromPtr<DIM, DIM>(
160 physicalPtr->matOpsDataPtr->getDependentDataPtr(
"P")->data().data());
161 t_P_at_pts(
i,
J) = t_P(
i,
J);
162 next();
163 }
164 }
165
171
172 auto t_grad_at_pts = get_grad_at_pts();
173 auto t_P_dF_at_pts = get_P_dF_at_pts();
175
176 auto next = [&]() {
177 ++t_grad_at_pts;
178 ++t_P_dF_at_pts;
179 };
180
181 for (auto gg = 0; gg != nb_integration_pts; ++gg) {
182 auto t_F = getFTensor2FromPtr<DIM, DIM>(
183 physicalPtr->matOpsDataPtr->getActiveDataPtr(
"F")->data().data());
184 t_F(
i,
J) = t_grad_at_pts(
i,
J);
187 auto t_P_dF = getFTensor4FromPtr<DIM, DIM, DIM, DIM>(
188 physicalPtr->matOpsDataPtr->getDependentDerivativesDataPtr(
"P_dF")
189 ->data()
190 .data());
191 t_P_dF_at_pts(
i,
J,
k,
L) = t_P_dF(
i,
J,
k,
L);
192 next();
193 }
194 }
195
199
200 auto t_grad_at_pts = get_grad_at_pts();
201 auto next = [&]() { ++t_grad_at_pts; };
202
203 for (int gg = 0; gg != nb_integration_pts; ++gg) {
204 auto t_F = getFTensor2FromPtr<DIM, DIM>(
205 physicalPtr->matOpsDataPtr->getActiveDataPtr(
"F")->data().data());
206 t_F(
i,
J) = t_grad_at_pts(
i,
J);
209 next();
210 }
211 }
212
214}
#define FTENSOR_INDEX(DIM, I)
#define CHK_THROW_MESSAGE(err, msg)
Check and throw MoFEM exception.
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
@ MOFEM_OPERATION_UNSUCCESSFUL
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
FTensor::Index< 'i', SPACE_DIM > i
FTensor::Index< 'J', DIM1 > J
FTensor::Index< 'k', 3 > k
DataLayoutTraits< DataLayout::GaussByCoeffs > DL
decltype(GetFTensor4FromMatImpl< Tensor_Dim0, Tensor_Dim1, Tensor_Dim2, Tensor_Dim3, S, DL, M >::get(std::declval< M & >(), 0, 0)) GetFTensor4FromMatType
decltype(GetFTensor2FromMatImpl< Tensor_Dim0, Tensor_Dim1, S, DL, M >::get(std::declval< M & >(), 0, 0)) GetFTensor2FromMatType
Structure for user loop methods on finite elements.
EntityHandle getFEEntityHandle() const
Get the entity handle of the current finite element.
