v0.13.1
Loading...
Searching...
No Matches
Public Member Functions | Public Attributes | List of all members
ConeRigidBody Struct Reference

#include <users_modules/multifield_plasticity/src/RigidBodies.hpp>

Inheritance diagram for ConeRigidBody:
[legend]
Collaboration diagram for ConeRigidBody:
[legend]

Public Member Functions

 ConeRigidBody (VectorDouble c_coords, VectorDouble roller_disp, int id)
 
MoFEMErrorCode getRollerDataForTag ()
 
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 > 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)
 
Tensor1< double, 3 > getdGap (Tensor1< TPack3, 3 > &t_coords, Tensor1< TPack1, 3 > &t_normal)
 
template<typename T1 , typename T2 >
Tensor1< double, 3 > & getNormalImpl (Tensor1< T1, 3 > &t_coords, Tensor1< T2, 3 > &t_disp)
 
template<typename T1 >
double getGapImpl (Tensor1< T1, 3 > &t_coords)
 
template<typename T1 , typename T2 , typename T3 >
Tensor2< double, 3, 3 > getDiffNormalImpl (Tensor1< T1, 3 > &t_coords, Tensor1< T2, 3 > &t_disp, Tensor1< T3, 3 > &t_normal)
 
template<typename T1 , typename T2 >
Tensor1< double, 3 > getdGapImpl (Tensor1< T1, 3 > &t_coords, Tensor1< T2, 3 > &t_normal)
 
virtual MoFEMErrorCode getBodyOptions ()
 
- Public Member Functions inherited from RigidBodyData
 RigidBodyData (VectorDouble c_coords, VectorDouble roller_disp, int id)
 
 RigidBodyData ()=delete
 
virtual ~RigidBodyData ()
 
virtual Tensor1< double, 3 > getNormal (Tensor1< TPack3, 3 > &t_coords, Tensor1< TPack1, 3 > &t_disp)=0
 
virtual Tensor1< double, 3 > getNormal (Tensor1< TPack3, 3 > &t_coords, Tensor1< double, 3 > &t_disp)=0
 
virtual Tensor1< double, 3 > getNormal (Tensor1< double, 3 > &t_coords, Tensor1< double, 3 > &t_disp)=0
 
virtual Tensor2< double, 3, 3 > getDiffNormal (Tensor1< TPack3, 3 > &t_coords, Tensor1< TPack1, 3 > &t_disp, Tensor1< TPack1, 3 > &t_normal)=0
 
virtual double getGap (Tensor1< TPack3, 3 > &t_coords)=0
 
virtual Tensor1< double, 3 > getdGap (Tensor1< TPack3, 3 > &t_coords, Tensor1< TPack1, 3 > &t_normal)=0
 
virtual MoFEMErrorCode getBodyOptions ()=0
 
MoFEMErrorCode computeRotationMatrix ()
 
Tensor1< double, 3 > getBodyOffset ()
 
template<typename T1 , typename T2 >
Tensor1< double, 3 > & getPointCoords (Tensor1< T1, 3 > &t_coords, Tensor1< T2, 3 > &t_disp)
 
MoFEMErrorCode saveBasicDataOnTag (moab::Interface &moab_debug, EntityHandle &vertex)
 
virtual MoFEMErrorCode getRollerDataForTag ()=0
 

Public Attributes

shared_ptr< ArbitrarySurfaceDatacP
 
double sLope
 
double hEight
 
double rAdius
 
double smallRadius
 
double oFfset
 
double aNgle
 
- Public Attributes inherited from RigidBodyData
const int iD
 
Index< 'i', 3 > i
 
Index< 'j', 3 > j
 
Index< 'k', 3 > k
 
VectorDouble originCoords
 
VectorDouble rollerDisp
 
VectorDouble BodyDispScaled
 
VectorDouble BodyDirectionScaled
 
double gAp
 
Tensor1< double, 3 > tNormal
 
Tensor1< double, 3 > dGap
 
Tensor1< double, 3 > pointCoords
 
Tensor1< double, 3 > closestPoint
 
Tensor1< double, 3 > defaultOrientation
 
Tensor1< double, 3 > oRientation
 
Tensor2< double, 3, 3 > rotationMat
 
Tensor2< double, 3, 3 > diffNormal
 
boost::shared_ptr< TimeAccelerogrammethodOpForRollerPosition
 
boost::shared_ptr< TimeAccelerogrammethodOpForRollerDirection
 
int bodyType
 
array< double, 9 > dataForTags
 

Additional Inherited Members

- Public Types inherited from RigidBodyData
using TPack3 = PackPtr< double *, 3 >
 
using TPack1 = PackPtr< double *, 1 >
 

Detailed Description

Definition at line 990 of file RigidBodies.hpp.

Constructor & Destructor Documentation

◆ ConeRigidBody()

ConeRigidBody::ConeRigidBody ( VectorDouble  c_coords,
VectorDouble  roller_disp,
int  id 
)
inline

Definition at line 999 of file RigidBodies.hpp.

1000 : RigidBodyData(c_coords, roller_disp, id) {
1001 // , cP(R, r, 1)
1002 sLope = 1;
1003 hEight = 1;
1004 rAdius = 1;
1005 oFfset = 0;
1006 // compute rotation matrix here
1007 }
RigidBodyData()=delete

Member Function Documentation

◆ getBodyOptions()

virtual MoFEMErrorCode ConeRigidBody::getBodyOptions ( )
inlinevirtual

Implements RigidBodyData.

Definition at line 1209 of file RigidBodies.hpp.

1209 {
1211 PetscBool rflg;
1212 int nb_dirs = 3;
1213 double angle = 45;
1214 CHKERR PetscOptionsBegin(PETSC_COMM_WORLD, "", "", "");
1215 string param = "-radius" + to_string(iD);
1216 CHKERR PetscOptionsScalar(param.c_str(), "set cone radius", "", rAdius,
1217 &rAdius, PETSC_NULL);
1218 param = "-height" + to_string(iD);
1219 CHKERR PetscOptionsScalar(param.c_str(), "set cone height", "", hEight,
1220 &hEight, PETSC_NULL);
1221 param = "-angle" + to_string(iD);
1222 CHKERR PetscOptionsScalar(param.c_str(), "set cone half angle", "", angle,
1223 &angle, PETSC_NULL);
1224 param = "-direction" + to_string(iD);
1225 CHKERR PetscOptionsGetRealArray(NULL, NULL, param.c_str(), &oRientation(0),
1226 &nb_dirs, &rflg);
1227 ierr = PetscOptionsEnd();
1228 CHKERRG(ierr);
1229 aNgle = angle;
1230 sLope = std::tan(aNgle * M_PI / 180);
1231 oFfset = (rAdius / sLope) - hEight;
1233
1234 if (rflg)
1236
1237 cP = make_shared<ArbitrarySurfaceData>(1, sLope);
1238
1240 }
static PetscErrorCode ierr
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:483
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
#define CHKERR
Inline error check.
Definition: definitions.h:535
PetscErrorCode PetscOptionsGetRealArray(PetscOptions *, const char pre[], const char name[], PetscReal dval[], PetscInt *nmax, PetscBool *set)
shared_ptr< ArbitrarySurfaceData > cP
double smallRadius
Tensor1< double, 3 > oRientation
Definition: RigidBodies.hpp:44
MoFEMErrorCode computeRotationMatrix()
Definition: RigidBodies.hpp:91
const int iD
Definition: RigidBodies.hpp:30

◆ getdGap()

Tensor1< double, 3 > ConeRigidBody::getdGap ( Tensor1< TPack3, 3 > &  t_coords,
Tensor1< TPack1, 3 > &  t_normal 
)
inlinevirtual

Implements RigidBodyData.

Definition at line 1040 of file RigidBodies.hpp.

1041 {
1042 return getdGapImpl(t_coords, t_normal);
1043 }
Tensor1< double, 3 > getdGapImpl(Tensor1< T1, 3 > &t_coords, Tensor1< T2, 3 > &t_normal)

◆ getdGapImpl()

template<typename T1 , typename T2 >
Tensor1< double, 3 > ConeRigidBody::getdGapImpl ( Tensor1< T1, 3 > &  t_coords,
Tensor1< T2, 3 > &  t_normal 
)
inline

Definition at line 1204 of file RigidBodies.hpp.

1205 {
1206 return dGap;
1207 };
Tensor1< double, 3 > dGap
Definition: RigidBodies.hpp:40

◆ getDiffNormal()

Tensor2< double, 3, 3 > ConeRigidBody::getDiffNormal ( Tensor1< TPack3, 3 > &  t_coords,
Tensor1< TPack1, 3 > &  t_disp,
Tensor1< TPack1, 3 > &  t_normal 
)
inlinevirtual

Implements RigidBodyData.

Definition at line 1031 of file RigidBodies.hpp.

1033 {
1034 return getDiffNormalImpl(t_coords, t_disp, t_normal);
1035 }
Tensor2< double, 3, 3 > getDiffNormalImpl(Tensor1< T1, 3 > &t_coords, Tensor1< T2, 3 > &t_disp, Tensor1< T3, 3 > &t_normal)

◆ getDiffNormalImpl()

template<typename T1 , typename T2 , typename T3 >
Tensor2< double, 3, 3 > ConeRigidBody::getDiffNormalImpl ( Tensor1< T1, 3 > &  t_coords,
Tensor1< T2, 3 > &  t_disp,
Tensor1< T3, 3 > &  t_normal 
)
inline

Definition at line 1182 of file RigidBodies.hpp.

1184 {
1185 Tensor2<double, 3, 3> diff_normal;
1186 Tensor1<double, 3> norm_diff;
1187 const double clean_gap = gAp;
1188 const double eps = 1e-8;
1189 for (int ii = 0; ii != 3; ++ii) {
1190 Tensor1<double, 3> pert_disp{t_disp(0), t_disp(1), t_disp(2)};
1191 pert_disp(ii) += eps;
1192 auto pert_normal = getNormal(t_coords, pert_disp);
1193 norm_diff(i) = (pert_normal(i) - t_normal(i)) / eps;
1194 dGap(ii) = (gAp - clean_gap) / eps;
1195
1196 for (int jj = 0; jj != 3; ++jj) {
1197 diff_normal(jj, ii) = norm_diff(jj);
1198 }
1199 }
1200 return diff_normal;
1201 };
static const double eps
Tensor1< double, 3 > getNormal(Tensor1< TPack3, 3 > &t_coords, Tensor1< TPack1, 3 > &t_disp)
Index< 'i', 3 > i
Definition: RigidBodies.hpp:31

◆ getGap()

double ConeRigidBody::getGap ( Tensor1< TPack3, 3 > &  t_coords)
inlinevirtual

Implements RigidBodyData.

Definition at line 1037 of file RigidBodies.hpp.

1037 {
1038 return getGapImpl(t_coords);
1039 }
double getGapImpl(Tensor1< T1, 3 > &t_coords)

◆ getGapImpl()

template<typename T1 >
double ConeRigidBody::getGapImpl ( Tensor1< T1, 3 > &  t_coords)
inline

Definition at line 1177 of file RigidBodies.hpp.

1177 {
1178 return gAp;
1179 };

◆ getNormal() [1/3]

Tensor1< double, 3 > ConeRigidBody::getNormal ( Tensor1< double, 3 > &  t_coords,
Tensor1< double, 3 > &  t_disp 
)
inlinevirtual

Implements RigidBodyData.

Definition at line 1026 of file RigidBodies.hpp.

1027 {
1028 return getNormalImpl(t_coords, t_disp);
1029 }
Tensor1< double, 3 > & getNormalImpl(Tensor1< T1, 3 > &t_coords, Tensor1< T2, 3 > &t_disp)

◆ getNormal() [2/3]

Tensor1< double, 3 > ConeRigidBody::getNormal ( Tensor1< TPack3, 3 > &  t_coords,
Tensor1< double, 3 > &  t_disp 
)
inlinevirtual

Implements RigidBodyData.

Definition at line 1022 of file RigidBodies.hpp.

1023 {
1024 return getNormalImpl(t_coords, t_disp);
1025 }

◆ getNormal() [3/3]

Tensor1< double, 3 > ConeRigidBody::getNormal ( Tensor1< TPack3, 3 > &  t_coords,
Tensor1< TPack1, 3 > &  t_disp 
)
inlinevirtual

Implements RigidBodyData.

Definition at line 1018 of file RigidBodies.hpp.

1019 {
1020 return getNormalImpl(t_coords, t_disp);
1021 }

◆ getNormalImpl()

template<typename T1 , typename T2 >
Tensor1< double, 3 > & ConeRigidBody::getNormalImpl ( Tensor1< T1, 3 > &  t_coords,
Tensor1< T2, 3 > &  t_disp 
)
inline

Definition at line 1046 of file RigidBodies.hpp.

1047 {
1048
1049 auto t_d = getPointCoords(t_coords, t_disp);
1051 pointCoords(i) = rotationMat(i, j) * t_d(j);
1052 t_d(i) = pointCoords(i);
1053
1054 // Tensor1<double, 3> c_center(0., 0., 0.);
1055 // Tensor1<double, 3> q_minus_c(0., 0., 0.);
1056 // Tensor1<double, 3> K(0., 0., 0.);
1057 // auto Q = pointCoords;
1058 // auto P_Q = pointCoords;
1059 // auto Q_K = pointCoords;
1060
1061 const double half_height = hEight / 2.;
1062 // const double rad2 = rAdius * rAdius;
1063 // const double small_rad2 = smallRadius * smallRadius;
1064
1065 // tNormal(i) = 0.;
1066 // pointCoords(2) += (oFfset + half_height);
1067 // t_d(2) += oFfset + half_height;
1068 // if (t_d(2) > half_height) {
1069 // tNormal(2) = 1;
1070 // if (pow(t_d(1), 2) + pow(t_d(0), 2) < rad2) {
1071 // gAp = tNormal(i) * t_d(i) - half_height;
1072 // } else {
1073 // c_center(2) = half_height;
1074 // Q(2) = c_center(2); // Q
1075 // q_minus_c(i) = Q(i) - c_center(i);
1076 // q_minus_c.normalize();
1077 // K(i) = c_center(i) + rAdius * q_minus_c(i);
1078 // P_Q(i) = t_d(i) - Q(i);
1079 // Q_K(i) = Q(i) - K(i);
1080 // gAp = sqrt(P_Q.l2() * P_Q.l2() + Q_K.l2() * Q_K.l2());
1081 // }
1082 // } else if (t_d(2) < -half_height) {
1083 // tNormal(2) = -1;
1084 // if (pow(t_d(1), 2) + pow(t_d(0), 2) < small_rad2) {
1085 // gAp = tNormal(i) * t_d(i) - half_height;
1086 // } else {
1087 // c_center(2) = -half_height;
1088 // Q(2) = c_center(2); // Q
1089 // q_minus_c(i) = Q(i) - c_center(i);
1090 // q_minus_c.normalize();
1091 // K(i) = c_center(i) + smallRadius * q_minus_c(i);
1092 // P_Q(i) = t_d(i) - Q(i);
1093 // Q_K(i) = Q(i) - K(i);
1094 // gAp = sqrt(P_Q.l2() * P_Q.l2() + Q_K.l2() * Q_K.l2());
1095 // }
1096 // } else {
1097 // auto closest_pt_cone = cP->getClosestPointToCone(t_d);
1098 // closestPoint(i) = closest_pt_cone(i);
1099
1100 // tNormal(j) = t_d(j) - closestPoint(j);
1101 // auto test_normal = tNormal.l2();
1102 // if (cP->isPointInsideCone(t_d))
1103 // tNormal(i) *= -1;
1104 // tNormal.normalize();
1105
1106 // gAp = tNormal(i) * (t_d(i) - closestPoint(i));
1107 // }
1108
1109 auto clamp = [](auto x, auto min, auto max) {
1110 if (x < min)
1111 x = min;
1112 else if (x > max)
1113 x = max;
1114 return x;
1115 };
1116
1117 auto sd_capped_cone = [&](Tensor1<double, 3> p, double h, double r1,
1118 double r2) {
1119 typedef Tensor1<double, 2> vec2;
1121 vec2 cb;
1122 vec2 q(sqrt(p(i) * p(i)), p(2));
1123 vec2 k1(r2, h);
1124 vec2 k2(r2 - r1, 2.0 * h);
1125 vec2 ca(q(0) - std::min(q(0), (q(1) < 0.0) ? r1 : r2), abs(q(1)) - h);
1126 vec2 k1_q(k1(0) - q(0), k1(1) - q(1));
1127 cb(i) = q(i) - k1(i) +
1128 k2(i) * clamp((k1_q(i) * k2(i)) / (k2(i) * k2(i)), 0.0, 1.0);
1129 double s = (cb(0) < 0.0 && ca(1) < 0.0) ? -1.0 : 1.0;
1130 return s * sqrt(std::min((ca(i) * ca(i)), (cb(i) * cb(i))));
1131 };
1132
1133 const double eps = 1e-6;
1134 auto cone_distance =
1135 sd_capped_cone(pointCoords, half_height, rAdius, smallRadius);
1136
1137 auto get_d = [&](int a) {
1138 auto pt = pointCoords;
1139 pt(a) += eps;
1140 return (sd_capped_cone(pt, half_height, rAdius, smallRadius) -
1141 cone_distance) /
1142 eps;
1143 };
1144
1145 // for debugging
1146 // auto my_rand_mn = [](double M, double N) {
1147 // return M + (rand() / (RAND_MAX / (N - M)));
1148 // };
1149 // auto t_off = getBodyOffset();
1150 // std::ofstream afile("cone_pts.csv", std::ios::ate);
1151 // if (afile.is_open()) {
1152 // afile << "X,Y,Z\n";
1153 // Tensor1<double, 3> p(0, 0, 0);
1154 // for (int n = 0; n != 1e6; n++) {
1155 // double x = my_rand_mn(-3, 3);
1156 // double y = my_rand_mn(-3, 3);
1157 // double z = my_rand_mn(-3, 3);
1158 // Tensor1<double, 3> coords(x, y, z);
1159 // if (sd_capped_cone(coords, half_height, rAdius, smallRadius) < 0) {
1160
1161 // afile << x + t_off(0) << "," << y + t_off(1) << "," << z + t_off(2)
1162 // << "\n";
1163 // }
1164 // }
1165 // afile.close();
1166 // }
1167
1168 Tensor1<double, 3> my_normal(get_d(0), get_d(1), get_d(2));
1169 my_normal.normalize();
1170
1171 gAp = cone_distance;
1172 tNormal(i) = rotationMat(j, i) * my_normal(j);
1173
1174 return tNormal;
1175 };
static Index< 'p', 3 > p
constexpr double a
Definition: single.cpp:4
double h
Tensor1< double, 3 > pointCoords
Definition: RigidBodies.hpp:41
Tensor1< double, 3 > tNormal
Definition: RigidBodies.hpp:39
Tensor1< double, 3 > & getPointCoords(Tensor1< T1, 3 > &t_coords, Tensor1< T2, 3 > &t_disp)
Tensor2< double, 3, 3 > rotationMat
Definition: RigidBodies.hpp:45
Index< 'j', 3 > j
Definition: RigidBodies.hpp:32

◆ getRollerDataForTag()

MoFEMErrorCode ConeRigidBody::getRollerDataForTag ( )
inlinevirtual

Implements RigidBodyData.

Definition at line 1009 of file RigidBodies.hpp.

1009 {
1011 bodyType = CONE;
1012 dataForTags[0] = rAdius;
1013 dataForTags[2] = hEight;
1014 dataForTags[3] = aNgle;
1016 }
@ CONE
Definition: RigidBodies.hpp:21
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:447
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:440
array< double, 9 > dataForTags
Definition: RigidBodies.hpp:55

Member Data Documentation

◆ aNgle

double ConeRigidBody::aNgle

Definition at line 998 of file RigidBodies.hpp.

◆ cP

shared_ptr<ArbitrarySurfaceData> ConeRigidBody::cP

Definition at line 991 of file RigidBodies.hpp.

◆ hEight

double ConeRigidBody::hEight

Definition at line 993 of file RigidBodies.hpp.

◆ oFfset

double ConeRigidBody::oFfset

Definition at line 996 of file RigidBodies.hpp.

◆ rAdius

double ConeRigidBody::rAdius

Definition at line 994 of file RigidBodies.hpp.

◆ sLope

double ConeRigidBody::sLope

Definition at line 992 of file RigidBodies.hpp.

◆ smallRadius

double ConeRigidBody::smallRadius

Definition at line 995 of file RigidBodies.hpp.


The documentation for this struct was generated from the following file: