v0.13.0
Public Member Functions | Public Attributes | List of all members
RollerRigidBody Struct Reference

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

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

Public Member Functions

 RollerRigidBody (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)
 
MoFEMErrorCode getBodyOptions ()
 
template<typename T1 , typename T2 >
Tensor1< double, 3 > & getNormalImpl (Tensor1< T1, 3 > &t_coords, Tensor1< T2, 3 > &t_disp)
 
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 >
double getGapImpl (Tensor1< T1, 3 > &t_coords)
 
template<typename T1 , typename T2 >
Tensor1< double, 3 > getdGapImpl (Tensor1< T1, 3 > &t_coords, Tensor1< T2, 3 > &t_normal)
 
- Public Member Functions inherited from RigidBodyData
 RigidBodyData (VectorDouble c_coords, VectorDouble roller_disp, int id)
 
 RigidBodyData ()=delete
 
virtual ~RigidBodyData ()
 
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)
 

Public Attributes

double angleA
 
double angleB
 
double rAdius
 
double fIllet
 
double torusRadius
 
double hEight
 
double coneRadiusA
 
double coneTopRadiusA
 
double coneOffsetA
 
double coneRadiusB
 
double coneTopRadiusB
 
double coneOffsetB
 
- 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
 
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
 
string positionDataParamName
 
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 708 of file RigidBodies.hpp.

Constructor & Destructor Documentation

◆ RollerRigidBody()

RollerRigidBody::RollerRigidBody ( VectorDouble  c_coords,
VectorDouble  roller_disp,
int  id 
)

Definition at line 725 of file RigidBodies.hpp.

726  : RigidBodyData(c_coords, roller_disp, id) {}
RigidBodyData()=delete

Member Function Documentation

◆ getBodyOptions()

MoFEMErrorCode RollerRigidBody::getBodyOptions ( )
virtual

Implements RigidBodyData.

Definition at line 768 of file RigidBodies.hpp.

768  {
770  PetscBool rflg;
771  int nb_dirs = 3;
772 
773  CHKERR PetscOptionsBegin(PETSC_COMM_WORLD, "", "", "");
774  string param = "-radius" + to_string(iD);
775  CHKERR PetscOptionsScalar(param.c_str(), "set roller radius", "", rAdius,
776  &rAdius, PETSC_NULL);
777  param = "-fillet" + to_string(iD);
778  CHKERR PetscOptionsScalar(param.c_str(), "set torus small radius", "",
779  fIllet, &fIllet, PETSC_NULL);
780  param = "-angle_a" + to_string(iD);
781  CHKERR PetscOptionsScalar(param.c_str(), "set roller angle of attack", "",
782  angleA, &angleA, PETSC_NULL);
783  param = "-angle_b" + to_string(iD);
784  CHKERR PetscOptionsScalar(param.c_str(), "set roller back angle", "",
785  angleB, &angleB, PETSC_NULL);
786  hEight = 3 * fIllet;
787  param = "-height" + to_string(iD);
788  CHKERR PetscOptionsScalar(param.c_str(), "set roller height (optional)", "",
789  hEight, &hEight, PETSC_NULL);
790  param = "-direction" + to_string(iD);
791  CHKERR PetscOptionsGetRealArray(NULL, NULL, param.c_str(), &oRientation(0),
792  &nb_dirs, &rflg);
793  ierr = PetscOptionsEnd();
794  CHKERRG(ierr);
795  if (rflg)
797 
798  // cP = make_shared<ArbitrarySurfaceData>(rAdius, fIllet);
799 
801  angleA = std::tan(angleA * M_PI / 180);
802  angleB = std::tan(angleB * M_PI / 180);
803 
804  coneRadiusA = torusRadius + fIllet / sqrt(1 + angleA * angleA);
806  coneOffsetA = -fIllet * angleA / sqrt(1 + angleA * angleA) - hEight / 2.;
807 
808  coneTopRadiusB = torusRadius + fIllet / sqrt(1 + angleB * angleB);
810  coneOffsetB = fIllet * angleB / sqrt(1 + angleB * angleB) + hEight / 2.;
811 
813  }
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:496
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
#define CHKERR
Inline error check.
Definition: definitions.h:548
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Exceptions.hpp:87
PetscErrorCode PetscOptionsGetRealArray(PetscOptions *, const char pre[], const char name[], PetscReal dval[], PetscInt *nmax, PetscBool *set)
Tensor1< double, 3 > oRientation
Definition: RigidBodies.hpp:43
MoFEMErrorCode computeRotationMatrix()
Definition: RigidBodies.hpp:85
const int iD
Definition: RigidBodies.hpp:30

◆ getdGap()

Tensor1<double, 3> RollerRigidBody::getdGap ( Tensor1< TPack3, 3 > &  t_coords,
Tensor1< TPack1, 3 > &  t_normal 
)
virtual

Implements RigidBodyData.

Definition at line 763 of file RigidBodies.hpp.

764  {
765  return getdGapImpl(t_coords, t_normal);
766  }
Tensor1< double, 3 > getdGapImpl(Tensor1< T1, 3 > &t_coords, Tensor1< T2, 3 > &t_normal)

◆ getdGapImpl()

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

Definition at line 950 of file RigidBodies.hpp.

951  {
952  return dGap;
953  };
Tensor1< double, 3 > dGap
Definition: RigidBodies.hpp:39

◆ getDiffNormal()

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

Implements RigidBodyData.

Definition at line 753 of file RigidBodies.hpp.

755  {
756  return getDiffNormalImpl(t_coords, t_disp, t_normal);
757  }
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> RollerRigidBody::getDiffNormalImpl ( Tensor1< T1, 3 > &  t_coords,
Tensor1< T2, 3 > &  t_disp,
Tensor1< T3, 3 > &  t_normal 
)

Definition at line 924 of file RigidBodies.hpp.

926  {
927  Tensor2<double, 3, 3> diff_normal;
928  Tensor1<double, 3> norm_diff;
929  const double clean_gap = gAp;
930  const double eps = 1e-8;
931  for (int ii = 0; ii != 3; ++ii) {
932  Tensor1<double, 3> pert_disp{t_disp(0), t_disp(1), t_disp(2)};
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;
936  dGap(ii) = (gAp - clean_gap) / eps;
937 
938  for (int jj = 0; jj != 3; ++jj) {
939  diff_normal(jj, ii) = norm_diff(jj);
940  }
941  }
942  return diff_normal;
943  };
static const double eps
Index< 'i', 3 > i
Definition: RigidBodies.hpp:31
Tensor1< double, 3 > getNormal(Tensor1< TPack3, 3 > &t_coords, Tensor1< TPack1, 3 > &t_disp)

◆ getGap()

double RollerRigidBody::getGap ( Tensor1< TPack3, 3 > &  t_coords)
virtual

Implements RigidBodyData.

Definition at line 759 of file RigidBodies.hpp.

759  {
760  return getGapImpl(t_coords);
761  }
double getGapImpl(Tensor1< T1, 3 > &t_coords)

◆ getGapImpl()

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

Definition at line 945 of file RigidBodies.hpp.

945  {
946  return gAp;
947  };

◆ getNormal() [1/3]

Tensor1<double, 3> RollerRigidBody::getNormal ( Tensor1< double, 3 > &  t_coords,
Tensor1< double, 3 > &  t_disp 
)
virtual

Implements RigidBodyData.

Definition at line 748 of file RigidBodies.hpp.

749  {
750  return getNormalImpl(t_coords, t_disp);
751  }
Tensor1< double, 3 > & getNormalImpl(Tensor1< T1, 3 > &t_coords, Tensor1< T2, 3 > &t_disp)

◆ getNormal() [2/3]

Tensor1<double, 3> RollerRigidBody::getNormal ( Tensor1< TPack3, 3 > &  t_coords,
Tensor1< double, 3 > &  t_disp 
)
virtual

Implements RigidBodyData.

Definition at line 744 of file RigidBodies.hpp.

745  {
746  return getNormalImpl(t_coords, t_disp);
747  }

◆ getNormal() [3/3]

Tensor1<double, 3> RollerRigidBody::getNormal ( Tensor1< TPack3, 3 > &  t_coords,
Tensor1< TPack1, 3 > &  t_disp 
)
virtual

Implements RigidBodyData.

Definition at line 740 of file RigidBodies.hpp.

741  {
742  return getNormalImpl(t_coords, t_disp);
743  }

◆ getNormalImpl()

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

Definition at line 816 of file RigidBodies.hpp.

817  {
818  auto t_d = getPointCoords(t_coords, t_disp);
819  pointCoords(i) = rotationMat(i, j) * t_d(j);
820  t_d(i) = pointCoords(i);
821 
822  const double inv_eps = 1e6;
823  array<double, 3> gaps;
824  const double half_height = hEight / 2.;
825 
826  double cone_apex_offset_a = coneOffsetA + half_height;
827  double cone_apex_offset_b = coneOffsetB - half_height;
828 
829  // https://computergraphics.stackexchange.com/questions/7485/glsl-shapes-signed-distance-field-implementation-explanation
830  auto sd_infinite_cone = [&](Tensor1<double, 3> p, double slope) {
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);
835  tNormal(0) = p(0) / (prod2 * prod);
836  tNormal(1) = p(1) / (prod2 * prod);
837  tNormal(2) = slope / prod2;
838  return c1 * prod + c2 * p(2);
839  };
840 
841  const double full_h_a = -cone_apex_offset_a + coneRadiusA / angleA;
842  const double full_h_b = cone_apex_offset_b + coneTopRadiusB / angleB;
843 
844  Tensor1<double, 3> t_apx_1(t_d(0), t_d(1), t_d(2) - full_h_a);
845  Tensor1<double, 3> t_apx_2(t_d(0), t_d(1), t_d(2) + full_h_b);
846 
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));
850  tNormal(0) = (t_d(0) * (-torusRadius + prod)) / (prod * prod2);
851  tNormal(1) = (t_d(1) * (-torusRadius + prod)) / (prod * prod2);
852  tNormal(2) = t_d(2) / prod2;
853  return -fIllet + prod2;
854  };
855 
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);
861 
862  int nb = distance(gaps.begin(), std::min_element(gaps.begin(), gaps.end()));
863  gAp = gaps[nb];
864 
865  switch (nb) {
866  case 0: {
867  torus_dist();
868  break;
869  }
870  case 1: {
871  // compute normal for cone A
872  sd_infinite_cone(t_apx_1, angleA);
873  break;
874  }
875  case 2: {
876  // compute normal for cone b
877  sd_infinite_cone(t_apx_2, -angleB);
878  } break;
879  }
880 
881  // // for debugging
882  // auto my_rand_mn = [](double M, double N) {
883  // return M + (rand() / (RAND_MAX / (N - M)));
884  // };
885  // auto t_off = getBodyOffset();
886  // std::ofstream afile("cone_pts.csv", std::ios::ate);
887  // if (afile.is_open()) {
888  // afile << "x,y,z\n";
889  // Tensor1<double, 3> p(0, 0, 0);
890  // for (int n = 0; n != 1e6; n++) {
891  // double x = my_rand_mn(-1.5, 1.5);
892  // double y = my_rand_mn(-1.5, 1.5);
893  // double z = my_rand_mn(-1, 1);
894  // Tensor1<double, 3> coords(x, y, z);
895 
896  // t_d(i) = coords(i);
897  // t_dc1(i) = coords(i);
898  // t_dc2(i) = coords(i);
899  // t_dc1(2) -= full_h_a;
900  // t_dc2(2) += full_h_b;
901 
902  // gaps[0] = torus_dist();
903  // gaps[1] = t_d(2) < -cone_apex_offset_a ? inv_eps :
904  // sd_infinite_cone(t_dc1, angleA); gaps[2] = t_d(2) >
905  // -cone_apex_offset_b ? inv_eps : sd_infinite_cone(t_dc2, -angleB);
906 
907  // if (gaps[0] < 0 || gaps[1] < 0 || gaps[2] < 0) {
908 
909  // afile << x + t_off(0) << "," << y + t_off(1) << "," << z + t_off(2)
910  // << "\n";
911  // }
912  // }
913  // afile.close();
914  // }
915 
916  tNormal.normalize();
917  auto normal_copy = tNormal;
918  tNormal(i) = rotationMat(j, i) * normal_copy(j);
919 
920  return tNormal;
921  };
static Index< 'p', 3 > p
Tensor1< double, 3 > pointCoords
Definition: RigidBodies.hpp:40
Tensor1< double, 3 > tNormal
Definition: RigidBodies.hpp:38
Tensor1< double, 3 > & getPointCoords(Tensor1< T1, 3 > &t_coords, Tensor1< T2, 3 > &t_disp)
Tensor2< double, 3, 3 > rotationMat
Definition: RigidBodies.hpp:44
Index< 'j', 3 > j
Definition: RigidBodies.hpp:32

◆ getRollerDataForTag()

MoFEMErrorCode RollerRigidBody::getRollerDataForTag ( )
virtual

Implements RigidBodyData.

Definition at line 728 of file RigidBodies.hpp.

728  {
730  bodyType = ROLLER;
731  dataForTags[0] = rAdius;
732  dataForTags[1] = fIllet;
733  dataForTags[2] = hEight;
734  dataForTags[4] = angleA;
735  dataForTags[5] = angleB;
736 
738  }
@ ROLLER
Definition: RigidBodies.hpp:23
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:460
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:453
array< double, 9 > dataForTags
Definition: RigidBodies.hpp:50

Member Data Documentation

◆ angleA

double RollerRigidBody::angleA

Definition at line 710 of file RigidBodies.hpp.

◆ angleB

double RollerRigidBody::angleB

Definition at line 711 of file RigidBodies.hpp.

◆ coneOffsetA

double RollerRigidBody::coneOffsetA

Definition at line 719 of file RigidBodies.hpp.

◆ coneOffsetB

double RollerRigidBody::coneOffsetB

Definition at line 723 of file RigidBodies.hpp.

◆ coneRadiusA

double RollerRigidBody::coneRadiusA

Definition at line 717 of file RigidBodies.hpp.

◆ coneRadiusB

double RollerRigidBody::coneRadiusB

Definition at line 721 of file RigidBodies.hpp.

◆ coneTopRadiusA

double RollerRigidBody::coneTopRadiusA

Definition at line 718 of file RigidBodies.hpp.

◆ coneTopRadiusB

double RollerRigidBody::coneTopRadiusB

Definition at line 722 of file RigidBodies.hpp.

◆ fIllet

double RollerRigidBody::fIllet

Definition at line 713 of file RigidBodies.hpp.

◆ hEight

double RollerRigidBody::hEight

Definition at line 715 of file RigidBodies.hpp.

◆ rAdius

double RollerRigidBody::rAdius

Definition at line 712 of file RigidBodies.hpp.

◆ torusRadius

double RollerRigidBody::torusRadius

Definition at line 714 of file RigidBodies.hpp.


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