1070 {
1073
1074#ifndef NDEBUG
1077 "Topological data pointer is null");
1080 "Broken displacement data pointer is null");
1083 "Hybrid displacement pointer is null");
1086 "Adjoint hybrid displacement pointer is null");
1087#endif
1088
1089 const int nb_dofs = data.
getIndices().size();
1090 if (!nb_dofs)
1092
1093 const int nb_integration_pts = getGaussPts().size2();
1094 const int nb_base_functions = data.
getN().size2();
1095
1096#ifndef NDEBUG
1097 if (this->
nF.size() != nb_dofs)
1099 "Size of nF %ld != nb_dofs %d", this->
nF.size(), nb_dofs);
1100 if (data.
getDiffN().size1() != nb_integration_pts)
1102 "Differential of base functions should have the same number of "
1103 "integration points as the data");
1104 if (data.
getDiffN().size2() != nb_base_functions * 2)
1106 "Differential of base functions should have the same number of "
1107 "base functions as the data");
1108#endif
1109
1113
1114 auto &coords = getCoords();
1115
1116
1117 const double h = std::get<2>(Tools::getTricircumcenter3d(coords.data().data()));
1118
1119 auto t_w = getFTensor0IntegrationWeight();
1120 auto t_tangent1 = getFTensor1Tangent1AtGaussPts();
1121 auto t_tangent2 = getFTensor1Tangent2AtGaussPts();
1123 auto t_var_u_hybrid =
1126
1127 for (int gg = 0; gg != nb_integration_pts; ++gg) {
1131
1132 double area = std::sqrt(t_normal(
i) * t_normal(
i));
1134 t_da(
i) = t_normal(
i) / area;
1135 area /= 2.;
1137
1138 double tau_density = 0;
1140 auto t_u_broken =
1143 bd.getVarFlux(), nb_integration_pts);
1144 for (int ss = 0; ss != gg; ++ss) {
1145 ++t_u_broken;
1146 ++t_var_u_broken;
1147 }
1148
1149
1150
1151
1152 const double hybrid_hybrid = t_var_u_hybrid(
i) * t_u_hybrid(
i);
1153 const double broken_broken = t_var_u_broken(
i) * t_u_broken(
i);
1154 const double hybrid_broken = -t_var_u_hybrid(
i) * t_u_broken(
i);
1155 const double broken_hybrid = -t_var_u_broken(
i) * t_u_hybrid(
i);
1156 tau_density +=
1157 hybrid_hybrid + broken_broken + hybrid_broken + broken_hybrid;
1158 }
1159
1161 locJ += t_w * tau * area * tau_density;
1162
1163 auto t_nf = getFTensor1FromArray<SPACE_DIM, SPACE_DIM>(
nF);
1164 int rr = 0;
1165 for (; rr != nb_dofs /
SPACE_DIM; ++rr) {
1170
1172 t_w *
alphaTau * tau_density * (t_da(
i) * t_normal_dX(
i,
I)) /
h;
1173 ++t_diff_base;
1174 ++t_nf;
1175 }
1176 for (; rr != nb_base_functions; ++rr)
1177 ++t_diff_base;
1178
1179 ++t_w;
1180 ++t_tangent1;
1181 ++t_tangent2;
1182 ++t_u_hybrid;
1183 ++t_var_u_hybrid;
1184 }
1185
1187}
#define FTENSOR_INDEXES(DIM,...)
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#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', 3 > j
FTensor::Index< 'k', 3 > k
constexpr std::enable_if<(Dim0<=2 &&Dim1<=2), Tensor2_Expr< Levi_Civita< T >, T, Dim0, Dim1, i, j > >::type levi_civita(const Index< i, Dim0 > &, const Index< j, Dim1 > &)
levi_civita functions to make for easy adhoc use
auto getFTensor1FromMat(M &data, int rr=0, int cc=0)
Get tensor rank 1 (vector) form data matrix.
constexpr IntegrationType I
boost::shared_ptr< TopologicalData > topoData
MatrixDouble & getDiffN(const FieldApproximationBase base)
get derivatives of base functions
auto getFTensor1DiffN(const FieldApproximationBase base)
Get derivatives of base functions.
MatrixDouble & getN(const FieldApproximationBase base)
get base functions this return matrix (nb. of rows is equal to nb. of Gauss pts, nb....
const VectorInt & getIndices() const
Get global indices of degrees of freedom on entity.
VectorDouble nF
local right hand side vector
Public Member Functions inherited from EshelbianPlasticity::OpAssembleTopologicalObjectiveDerivativeImplBase< OpAssembleFace >