Integrate grad-grad operator.
889 {
891
892 const size_t nb_row_dofs = row_data.getIndices().size();
893 const size_t nb_col_dofs = col_data.getIndices().size();
894
895 if (nb_row_dofs && nb_col_dofs) {
896
897 const bool diag = (row_data.getFieldEntities()[0]->getLocalUniqueId() ==
898 col_data.getFieldEntities()[0]->getLocalUniqueId());
899
904
905
906 double vol = OpBase::getMeasure();
907
908
909 auto t_w = OpBase::getFTensor0IntegrationWeight();
910
911
912 auto t_row_diff_base = row_data.getFTensor1DiffN<
SPACE_DIM>();
913
914
915
916 auto t_D = getFTensor4DdgFromMat<SPACE_DIM, SPACE_DIM, S>(*
matD);
917
918
919 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
920
921
923
924
925 double a = t_w * vol *
betaCoeff(t_coords(0), t_coords(1), t_coords(2));
926
927
928 int rr = 0;
930
931
932 auto t_m = OpBase::template getLocMat<SPACE_DIM>(
SPACE_DIM * rr);
933
935
936
937
938 t_rowD(
l,
j,
k) = t_D(
i,
j,
k,
l) * (
a * t_row_diff_base(
i));
939
940
941 auto t_col_diff_base = col_data.getFTensor1DiffN<
SPACE_DIM>(gg, 0);
942
943
945 for (int cc = 0; cc <= nb_cols; ++cc) {
946
947
948 t_m(
i,
j) += t_rowD(
i,
j,
k) * t_col_diff_base(
k);
949
950
951 ++t_col_diff_base;
952
953
954 ++t_m;
955 }
956
957
958 ++t_row_diff_base;
959 }
960
962 ++t_row_diff_base;
963
964
965 ++t_w;
966 ++t_D;
967 ++t_coords;
968 }
969
970
971 if (diag) {
973 auto t_m_rr = getFTensor2FromArray<SPACE_DIM, SPACE_DIM, SPACE_DIM>(
975 auto t_m_cc = getFTensor2FromArray<SPACE_DIM, SPACE_DIM>(
979 t_m_rr(
i,
j) = t_m_cc(
j,
i);
980 ++t_m_rr;
981 ++t_m_cc;
982 }
983 }
984 }
985
986 }
987
989}
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
FTensor::Index< 'l', 3 > l
FTensor::Index< 'j', 3 > j
FTensor::Index< 'k', 3 > k
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
FTensor::Index< 'i', SPACE_DIM > i
summit Index