Integrate grad-grad operator.
1029 {
1031
1033
1040
1041 auto &row_hessian = row_data.getN(BaseDerivatives::SecondDerivative);
1042 auto &col_hessian = col_data.getN(BaseDerivatives::SecondDerivative);
1043
1044#ifndef NDEBUG
1047 "Wrong number of integration pts (%zu != %d)", row_hessian.size1(),
1049 }
1050 if (row_hessian.size2() !=
1053 "Wrong number of base functions (%zu != %d)",
1056 }
1059 "Wrong number of base functions (%zu < %d)", row_hessian.size2(),
1061 }
1064 "Wrong number of integration pts (%zu != %d)", col_hessian.size1(),
1066 }
1069 "Wrong number of base functions (%zu < %d)", col_hessian.size2(),
1071 }
1072#endif
1073
1074
1075 double vol = OpBase::getMeasure();
1076
1077
1078 auto t_w = OpBase::getFTensor0IntegrationWeight();
1079
1080 auto t_row_diff2 = getFTensor2SymmetricLowerFromPtr<SPACE_DIM>(
1081 &*row_hessian.data().begin());
1082
1083
1084
1085 using DL = DataLayoutTraits<DataLayout::GaussByCoeffs>;
1086 auto get_D_at_pts =
1087 MatrixSizeHelper<GetFTensor4DdgFromMatType<SPACE_DIM, SPACE_DIM, S, DL>,
1089 auto t_D_at_pts = get_D_at_pts();
1091
1092
1094
1095
1096 double a = t_w * vol;
1097
1098
1099 int rr = 0;
1101
1102 auto t_mat = getFTensor2FromPtr<FIELD_DIM, FIELD_DIM, FIELD_DIM>(
1104
1106 t_rowD(
k,
l) = t_D_at_pts(
i,
j,
k,
l) * (
a * t_row_diff2(
i,
j));
1107
1108
1109 auto t_col_diff2 =
1110 getFTensor2SymmetricLowerFromPtr<SPACE_DIM>(&col_hessian(gg, 0));
1111
1112
1114
1115
1116 t_mat(
I,
J) += (t_rowD(
i,
j) * t_col_diff2(
i,
j)) *
t_kd(
I,
J);
1117
1118
1119 ++t_col_diff2;
1120
1121
1122 ++t_mat;
1123 }
1124
1125
1126 ++t_row_diff2;
1127 }
1128
1130 ++t_row_diff2;
1131
1132
1133 ++t_w;
1134 ++t_D_at_pts;
1135 }
1136 }
1137
1139}
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
@ MOFEM_DATA_INCONSISTENCY
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
FTensor::Index< 'i', SPACE_DIM > i
FTensor::Index< 'J', DIM1 > J
FTensor::Index< 'l', 3 > l
FTensor::Index< 'j', 3 > j
FTensor::Index< 'k', 3 > k
DataLayoutTraits< DataLayout::GaussByCoeffs > DL
constexpr IntegrationType I
int nbRows
number of dofs on rows
int nbIntegrationPts
number of integration points
MatrixDouble locMat
local entity block matrix
int nbCols
number if dof on column
int nbRowBaseFunctions
number or row base functions