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

Inheritance diagram for CylinderRigidBody:

Collaboration diagram for CylinderRigidBody:

Public Member Functions
	CylinderRigidBody (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 , 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)

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
double	rAdius

double	hEight

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< TimeAccelerogram >	methodOpForRollerPosition

boost::shared_ptr< TimeAccelerogram >	methodOpForRollerDirection

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 384 of file RigidBodies.hpp.

Constructor & Destructor Documentation

◆ CylinderRigidBody()

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

inline

Definition at line 387 of file RigidBodies.hpp.

388 : RigidBodyData(c_coords, roller_disp, id), rAdius(1), hEight(INFINITY) {}

CylinderRigidBody::hEight

double hEight

Definition: RigidBodies.hpp:386

CylinderRigidBody::rAdius

double rAdius

Definition: RigidBodies.hpp:385

RigidBodyData::RigidBodyData

RigidBodyData()=delete

Member Function Documentation

◆ getBodyOptions()

MoFEMErrorCode CylinderRigidBody::getBodyOptions ( )

inlinevirtual

Implements RigidBodyData.

Definition at line 572 of file RigidBodies.hpp.

                                  {
    MoFEMFunctionBegin;
    PetscBool rflg;
    PetscInt nb_dirs = 3;
    CHKERR PetscOptionsBegin(PETSC_COMM_WORLD, "", "", "");
    string param = "-direction" + to_string(iD);
    CHKERR PetscOptionsGetRealArray(NULL, NULL, param.c_str(), &oRientation(0),
                                    &nb_dirs, &rflg);
    param = "-radius" + to_string(iD);
    CHKERR PetscOptionsScalar(param.c_str(), "set roller radius", "", rAdius,
                              &rAdius, PETSC_NULL);
    param = "-height" + to_string(iD);
    CHKERR PetscOptionsScalar(param.c_str(), "set cylinder height", "", hEight,
                              &hEight, PETSC_NULL);
    if (rflg)
      CHKERR computeRotationMatrix();
    ierr = PetscOptionsEnd();
    CHKERRG(ierr);
    MoFEMFunctionReturn(0);
  }

◆ getdGap()

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

inlinevirtual

Implements RigidBodyData.

Definition at line 421 of file RigidBodies.hpp.

                                                           {
    return getdGapImpl(t_coords, t_normal);
  }

◆ getdGapImpl()

template<typename T1 , typename T2 >

Tensor1< double, 3 > CylinderRigidBody::getdGapImpl	(	Tensor1< T1, 3 > &	t_coords,
		Tensor1< T2, 3 > &	t_normal
	)

inline

Definition at line 563 of file RigidBodies.hpp.

                                                                  {
    // dGap(i) =
    //     pointCoords(j) * diffNormal(i, j) +
    //     tNormal(j) * (kronecker_delta(i, j) - 2. * rAdius * diffNormal(i,
    //     j));
    return dGap;
  };

◆ getDiffNormal()

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

inlinevirtual

Implements RigidBodyData.

Definition at line 411 of file RigidBodies.hpp.

                                                                    {
    return getDiffNormalImpl(t_coords, t_disp, t_normal);
  }

◆ getDiffNormalImpl()

template<typename T1 , typename T2 , typename T3 >

Tensor2< double, 3, 3 > CylinderRigidBody::getDiffNormalImpl	(	Tensor1< T1, 3 > &	t_coords,
		Tensor1< T2, 3 > &	t_disp,
		Tensor1< T3, 3 > &	t_normal
	)

inline

Definition at line 535 of file RigidBodies.hpp.

                                                                           {
 
    Tensor2<double, 3, 3> diff_normal;
    Tensor1<double, 3> norm_diff;
    const double clean_gap = gAp;
    const double eps = 1e-8;
    for (int ii = 0; ii != 3; ++ii) {
      Tensor1<double, 3> pert_disp{t_disp(0), t_disp(1), t_disp(2)};
      pert_disp(ii) += eps;
      auto pert_normal = getNormal(t_coords, pert_disp);
      norm_diff(i) = (pert_normal(i) - t_normal(i)) / eps;
      dGap(ii) = (gAp - clean_gap) / eps;
 
      for (int jj = 0; jj != 3; ++jj) {
        diff_normal(jj, ii) = norm_diff(jj);
      }
    }
    return diff_normal;
  };

◆ getGap()

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

inlinevirtual

Implements RigidBodyData.

Definition at line 417 of file RigidBodies.hpp.

                                                     {
    return getGapImpl(t_coords);
  }

◆ getGapImpl()

template<typename T1 >

double CylinderRigidBody::getGapImpl ( Tensor1< T1, 3 > & t_coords )

inline

Definition at line 557 of file RigidBodies.hpp.

                                                                            {
    // gAp = tNormal(i) * (pointCoords(i) - tNormal(i) * rAdius);
    return gAp;
  };

◆ getNormal() [1/3]

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

inlinevirtual

Implements RigidBodyData.

Definition at line 406 of file RigidBodies.hpp.

                                                           {
    return getNormalImpl(t_coords, t_disp);
  }

◆ getNormal() [2/3]

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

inlinevirtual

Implements RigidBodyData.

Definition at line 402 of file RigidBodies.hpp.

                                                           {
    return getNormalImpl(t_coords, t_disp);
  }

◆ getNormal() [3/3]

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

inlinevirtual

Implements RigidBodyData.

Definition at line 398 of file RigidBodies.hpp.

                                                           {
    return getNormalImpl(t_coords, t_disp);
  }

◆ getNormalImpl()

template<typename T1 , typename T2 >

Tensor1< double, 3 > & CylinderRigidBody::getNormalImpl	(	Tensor1< T1, 3 > &	t_coords,
		Tensor1< T2, 3 > &	t_disp
	)

inline

Definition at line 427 of file RigidBodies.hpp.

                                                                   {
    auto t_d = getPointCoords(t_coords, t_disp);
    CHKERR computeRotationMatrix();
    pointCoords(i) = rotationMat(i, j) * t_d(j);
    Tensor1<double, 3> c_center(0., 0., 0.);
    Tensor1<double, 3> q_minus_c(0., 0., 0.);
    Tensor1<double, 3> K(0., 0., 0.);
    auto Q = pointCoords;
    auto P_Q = pointCoords;
    auto Q_K = pointCoords;
 
    // https://www.geometrictools.com/Documentation/DistanceToCircle3.pdf
    const double half_height = hEight / 2.;
    const double rad2 = rAdius * rAdius;
    tNormal(i) = 0.;
 
    if (pointCoords(2) > half_height) {
      tNormal(2) = 1;
      if (pow(pointCoords(1), 2) + pow(pointCoords(0), 2) < rad2) {
        gAp = tNormal(i) * pointCoords(i) - half_height;
      } else {
        c_center(2) = half_height;
        Q(2) = c_center(2); // Q
        q_minus_c(i) = Q(i) - c_center(i);
        q_minus_c.normalize();
        K(i) = c_center(i) + rAdius * q_minus_c(i);
        P_Q(i) = pointCoords(i) - Q(i);
        Q_K(i) = Q(i) - K(i);
        gAp = sqrt(P_Q.l2() * P_Q.l2() + Q_K.l2() * Q_K.l2());
      }
    } else if (pointCoords(2) < -half_height) {
      tNormal(2) = -1;
      if (pow(pointCoords(1), 2) + pow(pointCoords(0), 2) < rad2) {
        gAp = tNormal(i) * pointCoords(i) - half_height;
      } else {
        c_center(2) = -half_height;
        Q(2) = c_center(2); // Q
        q_minus_c(i) = Q(i) - c_center(i);
        q_minus_c.normalize();
        K(i) = c_center(i) + rAdius * q_minus_c(i);
        P_Q(i) = pointCoords(i) - Q(i);
        Q_K(i) = Q(i) - K(i);
        gAp = sqrt(P_Q.l2() * P_Q.l2() + Q_K.l2() * Q_K.l2());
      }
    } else {
      pointCoords(2) = 0;
      tNormal(i) = pointCoords(i);
      tNormal.normalize();
      gAp = tNormal(i) * (pointCoords(i) - tNormal(i) * rAdius);
    }
 
    // auto length = [&](auto v) { return sqrt(v(0) * v(0) + v(1) * v(1)); };
 
    // auto sd_cylinder = [&](Tensor1<double, 3> p) {
    //   typedef Tensor1<double, 2> vec2;
    //   Index<'i', 2> i;
    //   vec2 d(length(p) - rAdius, abs(p(2)) - half_height);
    //   vec2 mx_d(max(d(0), 0.0), max(d(1), 0.0));
    //   return min(max(d(0), d(1)), 0.0) + length(mx_d);
    // };
 
    // const double eps = 1e-6;
    // auto cyl_distance = sd_cylinder(pointCoords);
    // gAp = cyl_distance;
 
    // auto get_d = [&](int a) {
    //   auto pt = pointCoords;
    //   pt(a) += eps;
    //   return (sd_cylinder(pt) - cyl_distance) / eps;
    // };
 
    // for (int dd = 0; dd != 3; ++dd)
    //   tNormal(dd) = get_d(dd);
 
    // for debugging
    // auto my_rand_mn = [](double M, double N) {
    //   return M + (rand() / (RAND_MAX / (N - M)));
    // };
    // auto t_off = getBodyOffset();
    // std::ofstream afile("cylinder_pts.csv", std::ios::ate);
    // if (afile.is_open()) {
    //   afile << "X,Y,Z\n";
    //   Tensor1<double, 3> p(0, 0, 0);
    //   for (int n = 0; n != 1e6; n++) {
    //     double x = my_rand_mn(-300, 300);
    //     double y = my_rand_mn(-300, 300);
    //     double z = my_rand_mn(-300, 300);
    //     Tensor1<double, 3> coords(x, y, z);
    //     if (sd_cylinder(coords) < 0) {
 
    //       afile << x + t_off(0) << "," << y + t_off(1) << "," << z + t_off(2)
    //       << "\n";
    //     }
    //   }
    //   afile.close();
    // }
    // cout << "kill process " << endl;
 
    tNormal.normalize();
 
    auto normal_copy = tNormal;
    tNormal(i) = rotationMat(j, i) * normal_copy(j);
 
    return tNormal;
  };

◆ getRollerDataForTag()

MoFEMErrorCode CylinderRigidBody::getRollerDataForTag ( )

inlinevirtual

Implements RigidBodyData.

Definition at line 390 of file RigidBodies.hpp.

                                       {
    MoFEMFunctionBeginHot;
    bodyType = CYLINDER;
    dataForTags[0] = rAdius;
    dataForTags[2] = hEight;
    MoFEMFunctionReturnHot(0);
  }

Member Data Documentation

◆ hEight

double CylinderRigidBody::hEight

Definition at line 386 of file RigidBodies.hpp.

◆ rAdius

double CylinderRigidBody::rAdius

Definition at line 385 of file RigidBodies.hpp.

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

users_modules/multifield_plasticity/src/RigidBodies.hpp

Public Member Functions

Public Attributes

Additional Inherited Members

Detailed Description

Constructor & Destructor Documentation

◆ CylinderRigidBody()

Member Function Documentation

◆ getBodyOptions()

◆ getdGap()

◆ getdGapImpl()

◆ getDiffNormal()

◆ getDiffNormalImpl()

◆ getGap()

◆ getGapImpl()

◆ getNormal() [1/3]

◆ getNormal() [2/3]

◆ getNormal() [3/3]

◆ getNormalImpl()

◆ getRollerDataForTag()

Member Data Documentation

◆ hEight

◆ rAdius