91 const double s =
v.l2();
95 const double coef = 1. / (1. +
c);
96 Tensor2<double, 3, 3> t_v;
110 template <
typename T1,
typename T2>
124 std::array<double, 9> def_real;
125 std::fill(def_real.begin(), def_real.end(), 0);
129 CHKERR moab_debug.tag_get_handle(
"_ROLLER_ID", 1, MB_TYPE_INTEGER, tag0,
130 MB_TAG_CREAT | MB_TAG_SPARSE, &def_int);
131 CHKERR moab_debug.tag_set_data(tag0, &vertex, 1, &this->iD);
133 CHKERR moab_debug.tag_get_handle(
"_BODY_TYPE", 1, MB_TYPE_INTEGER, tag0,
134 MB_TAG_CREAT | MB_TAG_SPARSE, &def_int);
135 CHKERR moab_debug.tag_set_data(tag0, &vertex, 1, &this->bodyType);
137 CHKERR moab_debug.tag_get_handle(
"_ORIENTATION", 3, MB_TYPE_DOUBLE, tag0,
138 MB_TAG_CREAT | MB_TAG_SPARSE, &def_real);
141 CHKERR moab_debug.tag_get_handle(
"_ROLLER_DATA", 9, MB_TYPE_DOUBLE, tag0,
142 MB_TAG_CREAT | MB_TAG_SPARSE, &def_real);
143 CHKERR moab_debug.tag_set_data(tag0, &vertex, 1, this->dataForTags.data());
203 template <
typename T1,
typename T2>
223 template <
typename T1,
typename T2,
typename T3>
227 const Tensor2<double, 3, 3> t_diff_norm(0., 0., 0., 0., 0., 0., 0., 0., 0.);
231 template <
typename T1,
typename T2>
241 PetscInt nb_dirs = 3;
242 CHKERR PetscOptionsBegin(PETSC_COMM_WORLD,
"",
"",
"");
243 string param =
"-direction" + to_string(
iD);
246 param =
"-flip" + to_string(
iD);
247 CHKERR PetscOptionsBool(param.c_str(),
"flip the plane normal",
"",
fLip,
251 ierr = PetscOptionsEnd();
297 template <
typename T1,
typename T2>
307 template <
typename T1,
typename T2,
typename T3>
314 const double norm = t_d.l2();
315 const double norm3 = norm * norm * norm;
328 template <
typename T1,
typename T2>
341 CHKERR PetscOptionsBegin(PETSC_COMM_WORLD,
"",
"",
"");
343 string param =
"-radius" + to_string(
iD);
344 CHKERR PetscOptionsScalar(param.c_str(),
"set roller radius",
"",
rAdius,
346 ierr = PetscOptionsEnd();
395 template <
typename T1,
typename T2>
408 const double half_height =
hEight / 2.;
417 c_center(2) = half_height;
419 q_minus_c(
i) = Q(
i) - c_center(
i);
420 q_minus_c.normalize();
423 Q_K(
i) = Q(
i) -
K(
i);
424 gAp = sqrt(P_Q.l2() * P_Q.l2() + Q_K.l2() * Q_K.l2());
431 c_center(2) = -half_height;
433 q_minus_c(
i) = Q(
i) - c_center(
i);
434 q_minus_c.normalize();
437 Q_K(
i) = Q(
i) -
K(
i);
438 gAp = sqrt(P_Q.l2() * P_Q.l2() + Q_K.l2() * Q_K.l2());
502 template <
typename T1,
typename T2,
typename T3>
507 Tensor2<double, 3, 3> diff_normal;
509 const double clean_gap =
gAp;
510 const double eps = 1e-8;
511 for (
int ii = 0; ii != 3; ++ii) {
513 pert_disp(ii) +=
eps;
514 auto pert_normal =
getNormal(t_coords, pert_disp);
515 norm_diff(
i) = (pert_normal(
i) - t_normal(
i)) /
eps;
518 for (
int jj = 0; jj != 3; ++jj) {
519 diff_normal(jj, ii) = norm_diff(jj);
530 template <
typename T1,
typename T2>
543 PetscInt nb_dirs = 3;
544 CHKERR PetscOptionsBegin(PETSC_COMM_WORLD,
"",
"",
"");
545 string param =
"-direction" + to_string(
iD);
548 param =
"-radius" + to_string(
iD);
549 CHKERR PetscOptionsScalar(param.c_str(),
"set roller radius",
"",
rAdius,
551 param =
"-height" + to_string(
iD);
552 CHKERR PetscOptionsScalar(param.c_str(),
"set cylinder height",
"",
hEight,
556 ierr = PetscOptionsEnd();
607 template <
typename T1,
typename T2>
632 const double x = t_d(0);
633 const double y = t_d(1);
634 const double z = t_d(2);
636 const double prod = sqrt(pow(x, 2) + pow(y, 2));
637 const double prod2 = sqrt(pow(-RR + prod, 2) + pow(z, 2));
639 tNormal(0) = (x * (-RR + prod)) / (prod * prod2);
640 tNormal(1) = (y * (-RR + prod)) / (prod * prod2);
654 template <
typename T1,
typename T2,
typename T3>
658 Tensor2<double, 3, 3> diff_normal;
660 const double clean_gap =
gAp;
661 const double eps = 1e-8;
662 for (
int ii = 0; ii != 3; ++ii) {
664 pert_disp(ii) +=
eps;
665 auto pert_normal =
getNormal(t_coords, pert_disp);
666 norm_diff(
i) = (pert_normal(
i) - t_normal(
i)) /
eps;
669 for (
int jj = 0; jj != 3; ++jj) {
670 diff_normal(jj, ii) = norm_diff(jj);
676 template <
typename T1,
typename T2>
686 CHKERR PetscOptionsBegin(PETSC_COMM_WORLD,
"",
"",
"");
687 string param =
"-inner_radius" + to_string(
iD);
688 CHKERR PetscOptionsScalar(param.c_str(),
"set torus small radius",
"",
innerRadius,
690 param =
"-radius" + to_string(
iD);
691 CHKERR PetscOptionsScalar(param.c_str(),
"set torus large Radius",
"",
rAdius,
693 param =
"-direction" + to_string(
iD);
696 ierr = PetscOptionsEnd();
773 CHKERR PetscOptionsBegin(PETSC_COMM_WORLD,
"",
"",
"");
774 string param =
"-radius" + to_string(
iD);
775 CHKERR PetscOptionsScalar(param.c_str(),
"set roller radius",
"",
rAdius,
777 param =
"-fillet" + to_string(
iD);
778 CHKERR PetscOptionsScalar(param.c_str(),
"set torus small radius",
"",
780 param =
"-angle_a" + to_string(
iD);
781 CHKERR PetscOptionsScalar(param.c_str(),
"set roller angle of attack",
"",
783 param =
"-angle_b" + to_string(
iD);
784 CHKERR PetscOptionsScalar(param.c_str(),
"set roller back angle",
"",
787 param =
"-height" + to_string(
iD);
788 CHKERR PetscOptionsScalar(param.c_str(),
"set roller height (optional)",
"",
790 param =
"-direction" + to_string(
iD);
793 ierr = PetscOptionsEnd();
815 template <
typename T1,
typename T2>
822 const double inv_eps = 1e6;
823 array<double, 3> gaps;
824 const double half_height =
hEight / 2.;
826 double cone_apex_offset_a =
coneOffsetA + half_height;
827 double cone_apex_offset_b =
coneOffsetB - half_height;
831 double c1 = 1. / sqrt(1 + (slope * slope));
832 double c2 = c1 * slope;
833 auto prod = sqrt(
p(0) *
p(0) +
p(1) *
p(1));
834 auto prod2 = sqrt(1 + slope * slope);
838 return c1 * prod + c2 *
p(2);
847 auto torus_dist = [&]() {
848 const double prod = sqrt(pow(t_d(0), 2) + pow(t_d(1), 2));
849 const double prod2 = sqrt(pow(-
torusRadius + prod, 2) + pow(t_d(2), 2));
856 gaps[0] = torus_dist();
857 gaps[1] = t_d(2) < -cone_apex_offset_a ? inv_eps
858 : sd_infinite_cone(t_apx_1,
angleA);
859 gaps[2] = t_d(2) > -cone_apex_offset_b ? inv_eps
860 : sd_infinite_cone(t_apx_2, -
angleB);
862 int nb = distance(gaps.begin(), std::min_element(gaps.begin(), gaps.end()));
872 sd_infinite_cone(t_apx_1,
angleA);
877 sd_infinite_cone(t_apx_2, -
angleB);
923 template <
typename T1,
typename T2,
typename T3>
927 Tensor2<double, 3, 3> diff_normal;
929 const double clean_gap =
gAp;
930 const double eps = 1e-8;
931 for (
int ii = 0; ii != 3; ++ii) {
933 pert_disp(ii) +=
eps;
934 auto pert_normal =
getNormal(t_coords, pert_disp);
935 norm_diff(
i) = (pert_normal(
i) - t_normal(
i)) /
eps;
938 for (
int jj = 0; jj != 3; ++jj) {
939 diff_normal(jj, ii) = norm_diff(jj);
949 template <
typename T1,
typename T2>
957 shared_ptr<ArbitrarySurfaceData>
cP;
1011 template <
typename T1,
typename T2>
1026 const double half_height =
hEight / 2.;
1074 auto clamp = [](
auto x,
auto min,
auto max) {
1087 vec2 q(sqrt(
p(
i) *
p(
i)),
p(2));
1089 vec2 k2(r2 - r1, 2.0 *
h);
1090 vec2 ca(q(0) - std::min(q(0), (q(1) < 0.0) ? r1 : r2), abs(q(1)) -
h);
1091 vec2 k1_q(k1(0) - q(0), k1(1) - q(1));
1092 cb(
i) = q(
i) - k1(
i) +
1093 k2(
i) * clamp((k1_q(
i) * k2(
i)) / (k2(
i) * k2(
i)), 0.0, 1.0);
1094 double s = (cb(0) < 0.0 && ca(1) < 0.0) ? -1.0 : 1.0;
1095 return s * sqrt(std::min((ca(
i) * ca(
i)), (cb(
i) * cb(
i))));
1098 const double eps = 1e-6;
1099 auto cone_distance =
1102 auto get_d = [&](
int a) {
1134 my_normal.normalize();
1136 gAp = cone_distance;
1146 template <
typename T1,
typename T2,
typename T3>
1150 Tensor2<double, 3, 3> diff_normal;
1152 const double clean_gap =
gAp;
1153 const double eps = 1e-8;
1154 for (
int ii = 0; ii != 3; ++ii) {
1156 pert_disp(ii) +=
eps;
1157 auto pert_normal =
getNormal(t_coords, pert_disp);
1158 norm_diff(
i) = (pert_normal(
i) - t_normal(
i)) /
eps;
1161 for (
int jj = 0; jj != 3; ++jj) {
1162 diff_normal(jj, ii) = norm_diff(jj);
1168 template <
typename T1,
typename T2>
1179 CHKERR PetscOptionsBegin(PETSC_COMM_WORLD,
"",
"",
"");
1180 string param =
"-radius" + to_string(
iD);
1181 CHKERR PetscOptionsScalar(param.c_str(),
"set cone radius",
"",
rAdius,
1183 param =
"-height" + to_string(
iD);
1184 CHKERR PetscOptionsScalar(param.c_str(),
"set cone height",
"",
hEight,
1186 param =
"-angle" + to_string(
iD);
1187 CHKERR PetscOptionsScalar(param.c_str(),
"set cone half angle",
"", angle,
1188 &angle, PETSC_NULL);
1189 param =
"-direction" + to_string(
iD);
1192 ierr = PetscOptionsEnd();
1202 cP = make_shared<ArbitrarySurfaceData>(1,
sLope);
1255 template <
typename T1,
typename T2>
1271 template <
typename T1,
typename T2>
1277 template <
typename T1,
typename T2,
typename T3>
1281 Tensor2<double, 3, 3> diff_normal;
1283 const double clean_gap =
gAp;
1284 const double eps = 1e-8;
1285 for (
int ii = 0; ii != 3; ++ii) {
1287 pert_disp(ii) +=
eps;
1288 auto pert_normal =
getNormal(t_coords, pert_disp);
1289 norm_diff(
i) = (pert_normal(
i) - t_normal(
i)) /
eps;
1292 for (
int jj = 0; jj != 3; ++jj) {
1293 diff_normal(jj, ii) = norm_diff(jj);
1301 CHKERR PetscOptionsBegin(PETSC_COMM_WORLD,
"",
"",
"");
1303 string param =
"-body_file_name" + to_string(
iD);
1305 CHKERR PetscOptionsString(param.c_str(),
"file name",
"",
"mesh.vtk",
1310 "body surface should be vtk surface");
1315 make_shared<OrientedBoxTreeTool>(
bodyMesh,
"ROOTSETSURF",
true);
1317 char normal_tag_name[255] =
"Normals";
1320 "Normal tag name %s cannot be found. Please "
1321 "provide the correct name. or run ./calculate_normals");
1324 ierr = PetscOptionsEnd();
#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 ...
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
@ MOFEM_DATA_INCONSISTENCY
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
double v
phase velocity of light in medium (cm/ns)
Tensor2_Expr< Kronecker_Delta< T >, T, Dim0, Dim1, i, j > kronecker_delta(const Index< i, Dim0 > &, const Index< j, Dim1 > &)
Rank 2.
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
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
UBlasVector< double > VectorDouble
PetscErrorCode PetscOptionsGetRealArray(PetscOptions *, const char pre[], const char name[], PetscReal dval[], PetscInt *nmax, PetscBool *set)
DeprecatedCoreInterface Interface
double h
convective heat coefficient
shared_ptr< ArbitrarySurfaceData > cP
Tensor2< double, 3, 3 > getDiffNormalImpl(Tensor1< T1, 3 > &t_coords, Tensor1< T2, 3 > &t_disp, Tensor1< T3, 3 > &t_normal)
Tensor1< double, 3 > getdGapImpl(Tensor1< T1, 3 > &t_coords, Tensor1< T2, 3 > &t_normal)
Tensor1< double, 3 > getNormal(Tensor1< TPack3, 3 > &t_coords, Tensor1< double, 3 > &t_disp)
Tensor1< double, 3 > & getNormalImpl(Tensor1< T1, 3 > &t_coords, Tensor1< T2, 3 > &t_disp)
ConeRigidBody(VectorDouble c_coords, VectorDouble roller_disp, int id)
Tensor1< double, 3 > getNormal(Tensor1< TPack3, 3 > &t_coords, Tensor1< TPack1, 3 > &t_disp)
double getGapImpl(Tensor1< T1, 3 > &t_coords)
MoFEMErrorCode getRollerDataForTag()
Tensor1< double, 3 > getNormal(Tensor1< double, 3 > &t_coords, Tensor1< double, 3 > &t_disp)
Tensor2< double, 3, 3 > getDiffNormal(Tensor1< TPack3, 3 > &t_coords, Tensor1< TPack1, 3 > &t_disp, Tensor1< TPack1, 3 > &t_normal)
double getGap(Tensor1< TPack3, 3 > &t_coords)
virtual MoFEMErrorCode getBodyOptions()
Tensor1< double, 3 > getdGap(Tensor1< TPack3, 3 > &t_coords, Tensor1< TPack1, 3 > &t_normal)
Tensor1< double, 3 > getNormal(Tensor1< TPack3, 3 > &t_coords, Tensor1< TPack1, 3 > &t_disp)
Tensor1< double, 3 > getNormal(Tensor1< double, 3 > &t_coords, Tensor1< double, 3 > &t_disp)
Tensor1< double, 3 > getNormal(Tensor1< TPack3, 3 > &t_coords, Tensor1< double, 3 > &t_disp)
Tensor2< double, 3, 3 > getDiffNormal(Tensor1< TPack3, 3 > &t_coords, Tensor1< TPack1, 3 > &t_disp, Tensor1< TPack1, 3 > &t_normal)
double getGapImpl(Tensor1< T1, 3 > &t_coords)
Tensor1< double, 3 > & getNormalImpl(Tensor1< T1, 3 > &t_coords, Tensor1< T2, 3 > &t_disp)
Tensor1< double, 3 > getdGapImpl(Tensor1< T1, 3 > &t_coords, Tensor1< T2, 3 > &t_normal)
MoFEMErrorCode getRollerDataForTag()
Tensor1< double, 3 > getdGap(Tensor1< TPack3, 3 > &t_coords, Tensor1< TPack1, 3 > &t_normal)
MoFEMErrorCode getBodyOptions()
Tensor2< double, 3, 3 > getDiffNormalImpl(Tensor1< T1, 3 > &t_coords, Tensor1< T2, 3 > &t_disp, Tensor1< T3, 3 > &t_normal)
double getGap(Tensor1< TPack3, 3 > &t_coords)
CylinderRigidBody(VectorDouble c_coords, VectorDouble roller_disp, int id)
double getGap(Tensor1< TPack3, 3 > &t_coords)
PlaneRigidBody(VectorDouble c_coords, VectorDouble roller_disp, int id)
Tensor1< double, 3 > getNormal(Tensor1< TPack3, 3 > &t_coords, Tensor1< TPack1, 3 > &t_disp)
Tensor1< double, 3 > getNormal(Tensor1< TPack3, 3 > &t_coords, Tensor1< double, 3 > &t_disp)
Tensor1< double, 3 > getdGap(Tensor1< TPack3, 3 > &t_coords, Tensor1< TPack1, 3 > &t_normal)
MoFEMErrorCode getRollerDataForTag()
double getGapImpl(Tensor1< T1, 3 > &t_coords)
Tensor2< double, 3, 3 > getDiffNormal(Tensor1< TPack3, 3 > &t_coords, Tensor1< TPack1, 3 > &t_disp, Tensor1< TPack1, 3 > &t_normal)
Tensor2< double, 3, 3 > getDiffNormalImpl(Tensor1< T1, 3 > &t_coords, Tensor1< T2, 3 > &t_disp, Tensor1< T3, 3 > &t_normal)
Tensor1< double, 3 > getdGapImpl(Tensor1< T1, 3 > &t_coords, Tensor1< T2, 3 > &t_normal)
MoFEMErrorCode getBodyOptions()
Tensor1< double, 3 > & getNormalImpl(Tensor1< T1, 3 > &t_coords, Tensor1< T2, 3 > &t_disp)
Tensor1< double, 3 > getNormal(Tensor1< double, 3 > &t_coords, Tensor1< double, 3 > &t_disp)
PackPtr< double *, 1 > TPack1
virtual MoFEMErrorCode getRollerDataForTag()=0
Tensor1< double, 3 > pointCoords
Tensor1< double, 3 > tNormal
virtual Tensor2< double, 3, 3 > getDiffNormal(Tensor1< TPack3, 3 > &t_coords, Tensor1< TPack1, 3 > &t_disp, Tensor1< TPack1, 3 > &t_normal)=0
Tensor1< double, 3 > dGap
Tensor1< double, 3 > oRientation
string positionDataParamName
RigidBodyData(VectorDouble c_coords, VectorDouble roller_disp, int id)
VectorDouble BodyDispScaled
virtual MoFEMErrorCode getBodyOptions()=0
virtual Tensor1< double, 3 > getdGap(Tensor1< TPack3, 3 > &t_coords, Tensor1< TPack1, 3 > &t_normal)=0
virtual Tensor1< double, 3 > getNormal(Tensor1< TPack3, 3 > &t_coords, Tensor1< double, 3 > &t_disp)=0
Tensor1< double, 3 > & getPointCoords(Tensor1< T1, 3 > &t_coords, Tensor1< T2, 3 > &t_disp)
Tensor2< double, 3, 3 > rotationMat
array< double, 9 > dataForTags
MoFEMErrorCode saveBasicDataOnTag(moab::Interface &moab_debug, EntityHandle &vertex)
VectorDouble originCoords
PackPtr< double *, 3 > TPack3
virtual Tensor1< double, 3 > getNormal(Tensor1< TPack3, 3 > &t_coords, Tensor1< TPack1, 3 > &t_disp)=0
MoFEMErrorCode computeRotationMatrix()
virtual Tensor1< double, 3 > getNormal(Tensor1< double, 3 > &t_coords, Tensor1< double, 3 > &t_disp)=0
Tensor1< double, 3 > getBodyOffset()
virtual double getGap(Tensor1< TPack3, 3 > &t_coords)=0
Tensor1< double, 3 > closestPoint
Tensor1< double, 3 > defaultOrientation
Tensor2< double, 3, 3 > diffNormal
Tensor1< double, 3 > getdGap(Tensor1< TPack3, 3 > &t_coords, Tensor1< TPack1, 3 > &t_normal)
Tensor1< double, 3 > & getNormalImpl(Tensor1< T1, 3 > &t_coords, Tensor1< T2, 3 > &t_disp)
Tensor2< double, 3, 3 > getDiffNormalImpl(Tensor1< T1, 3 > &t_coords, Tensor1< T2, 3 > &t_disp, Tensor1< T3, 3 > &t_normal)
Tensor1< double, 3 > getdGapImpl(Tensor1< T1, 3 > &t_coords, Tensor1< T2, 3 > &t_normal)
MoFEMErrorCode getBodyOptions()
double getGap(Tensor1< TPack3, 3 > &t_coords)
Tensor1< double, 3 > getNormal(Tensor1< TPack3, 3 > &t_coords, Tensor1< TPack1, 3 > &t_disp)
MoFEMErrorCode getRollerDataForTag()
Tensor1< double, 3 > getNormal(Tensor1< TPack3, 3 > &t_coords, Tensor1< double, 3 > &t_disp)
Tensor1< double, 3 > getNormal(Tensor1< double, 3 > &t_coords, Tensor1< double, 3 > &t_disp)
RollerRigidBody(VectorDouble c_coords, VectorDouble roller_disp, int id)
double getGapImpl(Tensor1< T1, 3 > &t_coords)
Tensor2< double, 3, 3 > getDiffNormal(Tensor1< TPack3, 3 > &t_coords, Tensor1< TPack1, 3 > &t_disp, Tensor1< TPack1, 3 > &t_normal)
Tensor1< double, 3 > getdGap(Tensor1< TPack3, 3 > &t_coords, Tensor1< TPack1, 3 > &t_normal)
double getGap(Tensor1< TPack3, 3 > &t_coords)
Tensor1< double, 3 > getNormal(Tensor1< TPack3, 3 > &t_coords, Tensor1< double, 3 > &t_disp)
MoFEMErrorCode getBodyOptions()
Tensor1< double, 3 > getdGapImpl(Tensor1< T1, 3 > &t_coords, Tensor1< T2, 3 > &t_normal)
moab::Interface * bodyMesh
Tensor1< double, 3 > getNormal(Tensor1< TPack3, 3 > &t_coords, Tensor1< TPack1, 3 > &t_disp)
Tensor2< double, 3, 3 > getDiffNormal(Tensor1< TPack3, 3 > &t_coords, Tensor1< TPack1, 3 > &t_disp, Tensor1< TPack1, 3 > &t_normal)
shared_ptr< OrientedBoxTreeTool > orientedBox
Tensor1< double, 3 > getNormal(Tensor1< double, 3 > &t_coords, Tensor1< double, 3 > &t_disp)
Tensor2< double, 3, 3 > getDiffNormalImpl(Tensor1< T1, 3 > &t_coords, Tensor1< T2, 3 > &t_disp, Tensor1< T3, 3 > &t_normal)
MoFEMErrorCode getRollerDataForTag()
Tensor1< double, 3 > & getNormalImpl(Tensor1< T1, 3 > &t_coords, Tensor1< T2, 3 > &t_disp)
STLRigidBody(VectorDouble c_coords, VectorDouble roller_disp, int id)
double getGapImpl(Tensor1< T1, 3 > &t_coords)
double getGap(Tensor1< TPack3, 3 > &t_coords)
SphereRigidBody(VectorDouble c_coords, VectorDouble roller_disp, int id)
Tensor1< double, 3 > getNormal(Tensor1< TPack3, 3 > &t_coords, Tensor1< TPack1, 3 > &t_disp)
MoFEMErrorCode getRollerDataForTag()
Tensor2< double, 3, 3 > getDiffNormal(Tensor1< TPack3, 3 > &t_coords, Tensor1< TPack1, 3 > &t_disp, Tensor1< TPack1, 3 > &t_normal)
MoFEMErrorCode getBodyOptions()
double getGapImpl(Tensor1< T1, 3 > &t_coords)
Tensor2< double, 3, 3 > getDiffNormalImpl(Tensor1< T1, 3 > &t_coords, Tensor1< T2, 3 > &t_disp, Tensor1< T3, 3 > &t_normal)
Tensor1< double, 3 > getdGapImpl(Tensor1< T1, 3 > &t_coords, Tensor1< T2, 3 > &t_normal)
Tensor1< double, 3 > getNormal(Tensor1< double, 3 > &t_coords, Tensor1< double, 3 > &t_disp)
Tensor1< double, 3 > & getNormalImpl(Tensor1< T1, 3 > &t_coords, Tensor1< T2, 3 > &t_disp)
Tensor1< double, 3 > getdGap(Tensor1< TPack3, 3 > &t_coords, Tensor1< TPack1, 3 > &t_normal)
Tensor1< double, 3 > getNormal(Tensor1< TPack3, 3 > &t_coords, Tensor1< double, 3 > &t_disp)
double getGap(Tensor1< TPack3, 3 > &t_coords)
MoFEMErrorCode getRollerDataForTag()
Tensor1< double, 3 > getNormal(Tensor1< double, 3 > &t_coords, Tensor1< double, 3 > &t_disp)
Tensor2< double, 3, 3 > getDiffNormalImpl(Tensor1< T1, 3 > &t_coords, Tensor1< T2, 3 > &t_disp, Tensor1< T3, 3 > &t_normal)
Tensor1< double, 3 > getdGapImpl(Tensor1< T1, 3 > &t_coords, Tensor1< T2, 3 > &t_normal)
virtual MoFEMErrorCode getBodyOptions()
double getGapImpl(Tensor1< T1, 3 > &t_coords)
TorusRigidBody(VectorDouble c_coords, VectorDouble roller_disp, int id)
Tensor2< double, 3, 3 > getDiffNormal(Tensor1< TPack3, 3 > &t_coords, Tensor1< TPack1, 3 > &t_disp, Tensor1< TPack1, 3 > &t_normal)
Tensor1< double, 3 > getNormal(Tensor1< TPack3, 3 > &t_coords, Tensor1< double, 3 > &t_disp)
Tensor1< double, 3 > & getNormalImpl(Tensor1< T1, 3 > &t_coords, Tensor1< T2, 3 > &t_disp)
Tensor1< double, 3 > getNormal(Tensor1< TPack3, 3 > &t_coords, Tensor1< TPack1, 3 > &t_disp)
Tensor1< double, 3 > getdGap(Tensor1< TPack3, 3 > &t_coords, Tensor1< TPack1, 3 > &t_normal)