62 {
64
65 int nb_integration_pts = OP::getGaussPts().size2();
66
67 auto get_tag = [&]() {
69 for (
const auto &tag_range_pair : *(
physicalPtr->tagVsRangePtr)) {
70 if (tag_range_pair.second.find(DomainEleOp::getFEEntityHandle()) !=
71 tag_range_pair.second.end()) {
72 return tag_range_pair.first;
73 }
74 }
75 }
76#ifndef NDEBUG
79 "ADOL-C tag not found " +
81 }
82#endif
84 };
85
86 const int current_tag = get_tag();
87 auto *fe_ptr =
const_cast<FEMethod *
>(this->getFEMethod());
89
90 physicalPtr->matOpsDataPtr->getActiveDataPtr(
"F")->resize(DIM, DIM,
false);
91 physicalPtr->matOpsDataPtr->getDependentDataPtr(
"P")->resize(DIM, DIM,
false);
92 physicalPtr->matOpsDataPtr->getDependentDerivativesDataPtr(
"P_dF")->resize(
93 DIM * DIM, DIM * DIM, false);
94 auto mat_grad_ptr =
physicalPtr->matOpsDataPtr->getCommonDataPtr(
"grad");
95 auto mat_P_ptr =
physicalPtr->matOpsDataPtr->getCommonDataPtr(
"P");
96 auto mat_P_dF_ptr =
physicalPtr->matOpsDataPtr->getCommonDataPtr(
"P_dF");
97
98#ifndef NDEBUG
99 if (!mat_grad_ptr || !mat_P_ptr || !mat_P_dF_ptr) {
101 "Missing common data for ADOL-C evaluation");
102 }
103 if (mat_grad_ptr->size2() != DIM * DIM) {
105 "Inconsistent size of gradient matrix for ADOL-C evaluation");
106 }
107 if (mat_grad_ptr->size1() != nb_integration_pts) {
109 "Inconsistent size of gradient matrix data for ADOL-C evaluation "
110 "%zu != %d",
111 mat_grad_ptr->size1(), nb_integration_pts);
112 }
113#endif
114
116 auto get_grad_at_pts =
118 *mat_grad_ptr, nb_integration_pts);
119 auto get_P_at_pts =
121 *mat_P_ptr, nb_integration_pts);
122 auto get_P_dF_at_pts =
124 DL>::size(*mat_P_dF_ptr, nb_integration_pts);
125
129
130 auto t_grad_at_pts = get_grad_at_pts();
131 auto t_P_at_pts = get_P_at_pts();
132
133 auto next = [&]() {
134 ++t_grad_at_pts;
135 ++t_P_at_pts;
136 };
137
138 for (auto gg = 0; gg != nb_integration_pts; ++gg) {
139 auto t_F = getFTensor2FromPtr<DIM, DIM>(
140 physicalPtr->matOpsDataPtr->getActiveDataPtr(
"F")->data().data());
141 t_F(
i,
J) = t_grad_at_pts(
i,
J);
144 auto t_P = getFTensor2FromPtr<DIM, DIM>(
145 physicalPtr->matOpsDataPtr->getDependentDataPtr(
"P")->data().data());
146 t_P_at_pts(
i,
J) = t_P(
i,
J);
147 next();
148 }
149 }
150
156
157 auto t_grad_at_pts = get_grad_at_pts();
158 auto t_P_dF_at_pts = get_P_dF_at_pts();
159 auto next = [&]() {
160 ++t_grad_at_pts;
161 ++t_P_dF_at_pts;
162 };
163
164 for (auto gg = 0; gg != nb_integration_pts; ++gg) {
165 auto t_F = getFTensor2FromPtr<DIM, DIM>(
166 physicalPtr->matOpsDataPtr->getActiveDataPtr(
"F")->data().data());
167 t_F(
i,
J) = t_grad_at_pts(
i,
J);
170 auto t_P_dF = getFTensor4FromPtr<DIM, DIM, DIM, DIM>(
171 physicalPtr->matOpsDataPtr->getDependentDerivativesDataPtr(
"P_dF")
172 ->data()
173 .data());
174 t_P_dF_at_pts(
i,
J,
k,
L) = t_P_dF(
i,
J,
k,
L);
175 next();
176 }
177 }
178
182
183 auto t_grad_at_pts = get_grad_at_pts();
184 auto next = [&]() { ++t_grad_at_pts; };
185
186 for (int gg = 0; gg != nb_integration_pts; ++gg) {
187 auto t_F = getFTensor2FromPtr<DIM, DIM>(
188 physicalPtr->matOpsDataPtr->getActiveDataPtr(
"F")->data().data());
189 t_F(
i,
J) = t_grad_at_pts(
i,
J);
192 next();
193 }
194 }
195
197}
#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.
