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

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

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

Public Member Functions

 STLRigidBody (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 >
Tensor1< double, 3 > getdGapImpl (Tensor1< T1, 3 > &t_coords, Tensor1< T2, 3 > &t_normal)
 
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)
 
MoFEMErrorCode getBodyOptions ()
 
- 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

moab::Core mbInstance
 
moab::Interface * bodyMesh
 
shared_ptr< OrientedBoxTreeTool > orientedBox
 
EntityHandle treeMeshset
 
Tag tagNormal
 
EntityHandle triangle
 
- 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 1208 of file RigidBodies.hpp.

Constructor & Destructor Documentation

◆ STLRigidBody()

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

Definition at line 1216 of file RigidBodies.hpp.

1217  : RigidBodyData(c_coords, roller_disp, id) {
1218  bodyMesh = &mbInstance;
1219  }
RigidBodyData()=delete
moab::Core mbInstance
moab::Interface * bodyMesh

Member Function Documentation

◆ getBodyOptions()

MoFEMErrorCode STLRigidBody::getBodyOptions ( )
virtual

Implements RigidBodyData.

Definition at line 1299 of file RigidBodies.hpp.

1299  {
1301  CHKERR PetscOptionsBegin(PETSC_COMM_WORLD, "", "", "");
1302  PetscBool rflg;
1303  string param = "-body_file_name" + to_string(iD);
1304  char mesh_file_name[255];
1305  CHKERR PetscOptionsString(param.c_str(), "file name", "", "mesh.vtk",
1306  mesh_file_name, 255, &rflg);
1307  if (string(mesh_file_name).substr(string(mesh_file_name).length() - 3) !=
1308  "vtk")
1309  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1310  "body surface should be vtk surface");
1311  bodyMesh->load_file(mesh_file_name, 0, "");
1312  Range tris;
1313  CHKERR bodyMesh->get_entities_by_dimension(0, 2, tris, false);
1314  orientedBox =
1315  make_shared<OrientedBoxTreeTool>(bodyMesh, "ROOTSETSURF", true);
1316  CHKERR orientedBox->build(tris, treeMeshset);
1317  char normal_tag_name[255] = "Normals";
1318  if (bodyMesh->tag_get_handle(normal_tag_name, tagNormal) != MB_SUCCESS) {
1319  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1320  "Normal tag name %s cannot be found. Please "
1321  "provide the correct name. or run ./calculate_normals");
1322  }
1323 
1324  ierr = PetscOptionsEnd();
1325  CHKERRG(ierr);
1327  }
#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
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:44
#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
char mesh_file_name[255]
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Exceptions.hpp:87
const int iD
Definition: RigidBodies.hpp:30
shared_ptr< OrientedBoxTreeTool > orientedBox
EntityHandle treeMeshset

◆ getdGap()

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

Implements RigidBodyData.

Definition at line 1250 of file RigidBodies.hpp.

1251  {
1252  return getdGapImpl(t_coords, t_normal);
1253  }
Tensor1< double, 3 > getdGapImpl(Tensor1< T1, 3 > &t_coords, Tensor1< T2, 3 > &t_normal)

◆ getdGapImpl()

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

Definition at line 1272 of file RigidBodies.hpp.

1273  {
1274  return dGap;
1275  };
Tensor1< double, 3 > dGap
Definition: RigidBodies.hpp:39

◆ getDiffNormal()

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

Implements RigidBodyData.

Definition at line 1240 of file RigidBodies.hpp.

1242  {
1243  return getDiffNormalImpl(t_coords, t_disp, t_normal);
1244  }
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> STLRigidBody::getDiffNormalImpl ( Tensor1< T1, 3 > &  t_coords,
Tensor1< T2, 3 > &  t_disp,
Tensor1< T3, 3 > &  t_normal 
)

Definition at line 1278 of file RigidBodies.hpp.

1280  {
1281  Tensor2<double, 3, 3> diff_normal;
1282  Tensor1<double, 3> norm_diff;
1283  const double clean_gap = gAp;
1284  const double eps = 1e-8;
1285  for (int ii = 0; ii != 3; ++ii) {
1286  Tensor1<double, 3> pert_disp{t_disp(0), t_disp(1), t_disp(2)};
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;
1290  dGap(ii) = (gAp - clean_gap) / eps;
1291 
1292  for (int jj = 0; jj != 3; ++jj) {
1293  diff_normal(jj, ii) = norm_diff(jj);
1294  }
1295  }
1296  return diff_normal;
1297  };
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 STLRigidBody::getGap ( Tensor1< TPack3, 3 > &  t_coords)
virtual

Implements RigidBodyData.

Definition at line 1246 of file RigidBodies.hpp.

1246  {
1247  return getGapImpl(t_coords);
1248  }
double getGapImpl(Tensor1< T1, 3 > &t_coords)

◆ getGapImpl()

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

Definition at line 1267 of file RigidBodies.hpp.

1267  {
1268  return gAp;
1269  };

◆ getNormal() [1/3]

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

Implements RigidBodyData.

Definition at line 1235 of file RigidBodies.hpp.

1236  {
1237  return getNormalImpl(t_coords, t_disp);
1238  }
Tensor1< double, 3 > & getNormalImpl(Tensor1< T1, 3 > &t_coords, Tensor1< T2, 3 > &t_disp)

◆ getNormal() [2/3]

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

Implements RigidBodyData.

Definition at line 1231 of file RigidBodies.hpp.

1232  {
1233  return getNormalImpl(t_coords, t_disp);
1234  }

◆ getNormal() [3/3]

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

Implements RigidBodyData.

Definition at line 1227 of file RigidBodies.hpp.

1228  {
1229  return getNormalImpl(t_coords, t_disp);
1230  }

◆ getNormalImpl()

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

Definition at line 1256 of file RigidBodies.hpp.

1257  {
1258  auto t_d = getPointCoords(t_coords, t_disp);
1259 
1260  CHKERR orientedBox->closest_to_location(&t_d(0), treeMeshset,
1261  &closestPoint(0), triangle);
1262  CHKERR bodyMesh->tag_get_data(tagNormal, &triangle, 1, &tNormal(0));
1263  gAp = tNormal(i) * (pointCoords(i) - closestPoint(i));
1264  return tNormal;
1265  };
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)
Tensor1< double, 3 > closestPoint
Definition: RigidBodies.hpp:41
EntityHandle triangle

◆ getRollerDataForTag()

MoFEMErrorCode STLRigidBody::getRollerDataForTag ( )
virtual

Implements RigidBodyData.

Definition at line 1221 of file RigidBodies.hpp.

1221  {
1223  bodyType = STL;
1225  }
@ STL
Definition: RigidBodies.hpp:24
#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

Member Data Documentation

◆ bodyMesh

moab::Interface* STLRigidBody::bodyMesh

Definition at line 1210 of file RigidBodies.hpp.

◆ mbInstance

moab::Core STLRigidBody::mbInstance

Definition at line 1209 of file RigidBodies.hpp.

◆ orientedBox

shared_ptr<OrientedBoxTreeTool> STLRigidBody::orientedBox

Definition at line 1211 of file RigidBodies.hpp.

◆ tagNormal

Tag STLRigidBody::tagNormal

Definition at line 1213 of file RigidBodies.hpp.

◆ treeMeshset

EntityHandle STLRigidBody::treeMeshset

Definition at line 1212 of file RigidBodies.hpp.

◆ triangle

EntityHandle STLRigidBody::triangle

Definition at line 1214 of file RigidBodies.hpp.


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