1074 {
1077
1078#ifndef NDEBUG
1081 "Topological data pointer is null");
1084 "Broken displacement data pointer is null");
1087 "Hybrid displacement pointer is null");
1090 "Adjoint hybrid displacement pointer is null");
1091#endif
1092
1093 const int nb_dofs = data.
getIndices().size();
1094 if (!nb_dofs)
1096
1097 const int nb_integration_pts = getGaussPts().size2();
1098 const int nb_base_functions = data.
getN().size2();
1099
1100#ifndef NDEBUG
1101 if (this->
nF.size() != nb_dofs)
1103 "Size of nF %ld != nb_dofs %d", this->
nF.size(), nb_dofs);
1104 if (data.
getDiffN().size1() != nb_integration_pts)
1106 "Differential of base functions should have the same number of "
1107 "integration points as the data");
1108 if (data.
getDiffN().size2() != nb_base_functions * 2)
1110 "Differential of base functions should have the same number of "
1111 "base functions as the data");
1112#endif
1113
1117
1118 auto &coords = getCoords();
1119
1120
1121 const double h = std::get<2>(Tools::getTricircumcenter3d(coords.data().data()));
1122
1123 auto t_w = getFTensor0IntegrationWeight();
1124 auto t_tangent1 = getFTensor1Tangent1AtGaussPts();
1125 auto t_tangent2 = getFTensor1Tangent2AtGaussPts();
1127 auto t_var_u_hybrid =
1130
1131 for (int gg = 0; gg != nb_integration_pts; ++gg) {
1135
1136 double area = std::sqrt(t_normal(
i) * t_normal(
i));
1138 t_da(
i) = t_normal(
i) / area;
1139 area /= 2.;
1141
1142 double tau_density = 0;
1144 auto t_u_broken =
1147 bd.getVarFlux(), nb_integration_pts);
1148 for (int ss = 0; ss != gg; ++ss) {
1149 ++t_u_broken;
1150 ++t_var_u_broken;
1151 }
1152
1153
1154
1155
1156 const double hybrid_hybrid = t_var_u_hybrid(
i) * t_u_hybrid(
i);
1157 const double broken_broken = t_var_u_broken(
i) * t_u_broken(
i);
1158 const double hybrid_broken = -t_var_u_hybrid(
i) * t_u_broken(
i);
1159 const double broken_hybrid = -t_var_u_broken(
i) * t_u_hybrid(
i);
1160 tau_density +=
1161 hybrid_hybrid + broken_broken + hybrid_broken + broken_hybrid;
1162 }
1163
1165 locJ += t_w * tau * area * tau_density;
1166
1167 auto t_nf = getFTensor1FromArray<SPACE_DIM, SPACE_DIM>(
nF);
1168 int rr = 0;
1169 for (; rr != nb_dofs /
SPACE_DIM; ++rr) {
1174
1176 t_w *
alphaTau * tau_density * (t_da(
i) * t_normal_dX(
i,
I)) /
h;
1177 ++t_diff_base;
1178 ++t_nf;
1179 }
1180 for (; rr != nb_base_functions; ++rr)
1181 ++t_diff_base;
1182
1183 ++t_w;
1184 ++t_tangent1;
1185 ++t_tangent2;
1186 ++t_u_hybrid;
1187 ++t_var_u_hybrid;
1188 }
1189
1191}
#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 >