Screw Class Reference

This function defines an Archimedes' screw in the z-direction from a (constant) starting point, a (constant) length L, a (constant) radius r, a (constant) number or revelations N and a (constant) rotation speed (rev/s) More...

#include <Screw.h>

+ Inheritance diagram for Screw:

Public Member Functions

 Screw ()
 Default constructor: make a screw with default parameters. More...
 
 Screw (const Screw &other)
 Copy constructor, copies another Screw. More...
 
 Screw (Vec3D start, Mdouble l, Mdouble r, Mdouble n, Mdouble omega, Mdouble thickness, ScrewType screwType=ScrewType::doubleHelix)
 Constructor in which all parameters of the screw are set. More...
 
 ~Screw () override
 Default destructor. More...
 
Screwcopy () const final
 Copy this screw and return a pointer to the copy. More...
 
bool getDistanceAndNormal (const BaseParticle &P, Mdouble &distance, Vec3D &normal_return) const final
 Compute the distance from the Screw for a given BaseParticle and return if there is a collision. If there is a collision, also return the normal vector of the interaction point. More...
 
bool getDistanceAndNormalLabCoordinates (Vec3D position, Mdouble wallInteractionRadius, Mdouble &distance, Vec3D &normal_return) const
 
void move_time (Mdouble dt)
 Rotate the Screw for a period dt, so that the offset_ changes with omega_*dt. More...
 
void rotate (const Vec3D &angularVelocityDt) override
 
void read (std::istream &is) override
 Reads a Screw from an input stream, for example a restart file. More...
 
void oldRead (std::istream &is)
 Reads a Screw in the old style from an input stream, for example a restart file old style. More...
 
void write (std::ostream &os) const override
 Writes this Screw to an output stream, for example a restart file. More...
 
std::string getName () const final
 Returns the name of the object, here the string "Screw". More...
 
void writeVTK (VTKContainer &vtk) const override
 
void writeVTK (std::string filename) const
 
- Public Member Functions inherited from BaseWall
 BaseWall ()
 Default constructor. More...
 
 BaseWall (const BaseWall &w)
 Copy constructor. More...
 
 ~BaseWall () override
 Default destructor. More...
 
virtual bool getDistanceNormalOverlap (const BaseParticle &P, Mdouble &distance, Vec3D &normal_return, Mdouble &overlap) const
 
virtual bool getDistanceNormalOverlapSuperquadric (const SuperQuadricParticle &p, Mdouble &distance, Vec3D &normal_return, Mdouble &overlap) const
 
virtual Vec3D getFurthestPointSuperQuadric (const Vec3D &normalBodyFixed, const Vec3D &axes, Mdouble eps1, Mdouble eps2) const
 
virtual void setHandler (WallHandler *handler)
 A function which sets the WallHandler for this BaseWall. More...
 
WallHandlergetHandler () const
 A function which returns the WallHandler that handles this BaseWall. More...
 
void setIndSpecies (unsigned int indSpecies) override
 Define the species of this wall using the index of the species in the SpeciesHandler in this DPMBase. More...
 
void setSpecies (const ParticleSpecies *species)
 Defines the species of the current wall. More...
 
bool isFixed () const override
 
void setForceControl (Vec3D forceGoal, Vec3D gainFactor, Vec3D baseVelocity={0, 0, 0})
 Slowly adjusts the force on a wall towards a specified goal, by adjusting (prescribing) the velocity of the wall. More...
 
virtual bool isLocal (Vec3D &min, Vec3D &max) const
 
bool getLinePlaneIntersect (Vec3D &intersect, const Vec3D &p0, const Vec3D &p1, const Vec3D &n, const Vec3D &p)
 
bool isInsideWallVTK (const Vec3D &point, const Vec3D &normal, const Vec3D &position) const
 
void projectOntoWallVTK (Vec3D &point0, const Vec3D &point1, const Vec3D &normal, const Vec3D &position) const
 
void intersectVTK (std::vector< Vec3D > &points, Vec3D normal, Vec3D position) const
 
virtual BaseInteractiongetInteractionWithSuperQuad (SuperQuadricParticle *p, unsigned timeStamp, InteractionHandler *interactionHandler)
 
void getVTK (std::vector< Vec3D > &points, std::vector< std::vector< double >> &triangleStrips)
 
const Vec3D getAxis () const
 
BaseInteractiongetInteractionWith (BaseParticle *p, unsigned timeStamp, InteractionHandler *interactionHandler) override
 Returns the interaction between this wall and a given particle, nullptr if there is no interaction. More...
 
virtual void actionsOnRestart ()
 No implementation but can be overidden in its derived classes. More...
 
virtual void actionsAfterParticleGhostUpdate ()
 No implementation but can be overidden in its derived classes. More...
 
virtual void handleParticleAddition (unsigned int id, BaseParticle *p)
 Handles the addition of particles to the particleHandler. More...
 
virtual void handleParticleRemoval (unsigned int id)
 Handles the addition of particles to the particleHandler. More...
 
virtual void checkInteractions (InteractionHandler *interactionHandler, unsigned int timeStamp)
 Check if all interactions are valid. More...
 
bool getVTKVisibility () const
 
void setVTKVisibility (bool vtkVisibility)
 
void addRenderedWall (BaseWall *w)
 
BaseWallgetRenderedWall (size_t i) const
 
std::vector< BaseWall * > getRenderedWalls () const
 
void removeRenderedWalls ()
 
void renderWall (VTKContainer &vtk)
 
void addParticlesAtWall (unsigned numElements=50)
 
void setVelocityControl (Vec3D forceGoal, Vec3D gainFactor, Vec3D baseVelocity)
 
virtual void writeWallDetailsVTK (VTKData &data) const
 
virtual void computeWear ()
 
- Public Member Functions inherited from BaseInteractable
 BaseInteractable ()
 Default BaseInteractable constructor. More...
 
 BaseInteractable (const BaseInteractable &p)
 Copy constructor. More...
 
 ~BaseInteractable () override
 Destructor, it simply destructs the BaseInteractable and all the objects it contains. More...
 
unsigned int getIndSpecies () const
 Returns the index of the species associated with the interactable object. More...
 
const ParticleSpeciesgetSpecies () const
 Returns a pointer to the species of this BaseInteractable. More...
 
void setSpecies (const ParticleSpecies *species)
 Sets the species of this BaseInteractable. More...
 
const Vec3DgetForce () const
 Returns the force on this BaseInteractable. More...
 
const Vec3DgetTorque () const
 Returns the torque on this BaseInteractable. More...
 
void setForce (const Vec3D &force)
 Sets the force on this BaseInteractable. More...
 
void setTorque (const Vec3D &torque)
 Sets the torque on this BaseInteractable. More...
 
void addForce (const Vec3D &addForce)
 Adds an amount to the force on this BaseInteractable. More...
 
void addTorque (const Vec3D &addTorque)
 Adds an amount to the torque on this BaseInteractable. More...
 
virtual void resetForceTorque (int numberOfOMPthreads)
 
void sumForceTorqueOMP ()
 
const Vec3DgetPosition () const
 Returns the position of this BaseInteractable. More...
 
const QuaterniongetOrientation () const
 Returns the orientation of this BaseInteractable. More...
 
virtual void setPosition (const Vec3D &position)
 Sets the position of this BaseInteractable. More...
 
void setOrientationViaNormal (Vec3D normal)
 Sets the orientation of this BaseInteractable by defining the vector that results from the rotation of the (1,0,0) vector. More...
 
void setOrientationViaEuler (Vec3D eulerAngle)
 Sets the orientation of this BaseInteractable by defining the euler angles. More...
 
virtual void setOrientation (const Quaternion &orientation)
 Sets the orientation of this BaseInteractable. More...
 
virtual void move (const Vec3D &move)
 Moves this BaseInteractable by adding an amount to the position. More...
 
const std::vector< BaseInteraction * > & getInteractions () const
 Returns a list of interactions which belong to this interactable. More...
 
void addInteraction (BaseInteraction *I)
 Adds an interaction to this BaseInteractable. More...
 
bool removeInteraction (BaseInteraction *I)
 Removes an interaction from this BaseInteractable. More...
 
void copyInteractionsForPeriodicParticles (const BaseInteractable &p)
 Copies interactions to this BaseInteractable whenever a periodic copy made. More...
 
void setVelocity (const Vec3D &velocity)
 set the velocity of the BaseInteractable. More...
 
void setAngularVelocity (const Vec3D &angularVelocity)
 set the angular velocity of the BaseInteractble. More...
 
void addVelocity (const Vec3D &velocity)
 adds an increment to the velocity. More...
 
void addAngularVelocity (const Vec3D &angularVelocity)
 add an increment to the angular velocity. More...
 
virtual const Vec3DgetVelocity () const
 Returns the velocity of this interactable. More...
 
virtual const Vec3DgetAngularVelocity () const
 Returns the angular velocity of this interactable. More...
 
void setPrescribedPosition (const std::function< Vec3D(double)> &prescribedPosition)
 Allows the position of an infinite mass interactable to be prescribed. More...
 
void applyPrescribedPosition (double time)
 Computes the position from the user defined prescribed position function. More...
 
void setPrescribedVelocity (const std::function< Vec3D(double)> &prescribedVelocity)
 Allows the velocity of an infinite mass interactable to be prescribed. More...
 
void applyPrescribedVelocity (double time)
 Computes the velocity from the user defined prescribed velocity function. More...
 
void setPrescribedOrientation (const std::function< Quaternion(double)> &prescribedOrientation)
 Allows the orientation of the infinite mass interactbale to be prescribed. More...
 
void applyPrescribedOrientation (double time)
 Computes the orientation from the user defined prescribed orientation function. More...
 
void setPrescribedAngularVelocity (const std::function< Vec3D(double)> &prescribedAngularVelocity)
 Allows the angular velocity of the infinite mass interactable to be prescribed. More...
 
void applyPrescribedAngularVelocity (double time)
 Computes the angular velocity from the user defined prescribed angular velocity. More...
 
virtual const Vec3D getVelocityAtContact (const Vec3D &contact) const
 Returns the velocity at the contact point, use by many force laws. More...
 
void integrateBeforeForceComputation (double time, double timeStep)
 This is part of integrate routine for objects with infinite mass. More...
 
void integrateAfterForceComputation (double time, double timeStep)
 This is part of the integration routine for objects with infinite mass. More...
 
virtual Mdouble getInvMass () const
 
virtual Mdouble getCurvature (const Vec3D &labFixedCoordinates) const
 
virtual bool isFaceContact (const Vec3D &normal) const
 
- Public Member Functions inherited from BaseObject
 BaseObject ()=default
 Default constructor. More...
 
 BaseObject (const BaseObject &p)=default
 Copy constructor, copies all the objects BaseObject contains. More...
 
virtual ~BaseObject ()=default
 virtual destructor More...
 
virtual void moveInHandler (unsigned int index)
 Except that it is virtual, it does the same thing as setIndex() does. More...
 
void setIndex (unsigned int index)
 Allows one to assign an index to an object in the handler/container. More...
 
void setId (unsigned long id)
 Assigns a unique identifier to each object in the handler (container) which remains constant even after the object is deleted from the container/handler. More...
 
unsigned int getIndex () const
 Returns the index of the object in the handler. More...
 
unsigned int getId () const
 Returns the unique identifier of any particular object. More...
 
void setGroupId (unsigned groupId)
 
unsigned getGroupId () const
 

Private Attributes

Vec3D start_
 The centre of the lower end of the screw. More...
 
Mdouble l_
 The length of the Screw. More...
 
Mdouble maxR_
 The outer radius of the Screw. More...
 
Mdouble n_
 The number of revelations. More...
 
Mdouble omega_
 Rotation speed in rev/s. More...
 
Mdouble offset_
 The angle that describes how much the Screw has turned, going from 0 to 1 for a rotation. More...
 
Mdouble thickness_
 The thickness of the Screw. More...
 
ScrewType screwType_
 Single or double helix screw. More...
 

Additional Inherited Members

- Static Public Member Functions inherited from BaseWall
static void addToVTK (const std::vector< Vec3D > &points, VTKContainer &vtk)
 Takes the points provided and adds a triangle strip connecting these points to the vtk container. More...
 

Detailed Description

This function defines an Archimedes' screw in the z-direction from a (constant) starting point, a (constant) length L, a (constant) radius r, a (constant) number or revelations N and a (constant) rotation speed (rev/s)

q is a new coordinate going from 0 to 1 and t is the time, x=xs+r*cos(2*pi*(offset+N*q)), y=ys+r*sin(2*pi*(offset+N*q)), z=zs+q*L

Todo:
IFCD: Can these details about class Screw be made more clear? I don't understand them.

Constructor & Destructor Documentation

◆ Screw() [1/3]

Screw::Screw ( )

Default constructor: make a screw with default parameters.

Make a Screw which is centered in the origin, has a length of 1, one revelation, a radius of 1, turns with 1 revelation per second, is infinitely thin and starts at its normal initial point.

42 {
43  start_.setZero();
44  l_ = 1.0;
45  maxR_ = 1.0;
46  n_ = 1.0;
47  omega_ = 1.0;
48  offset_ = 0.0;
49  thickness_ = 0.0;
51  setOrientationViaNormal({0, 0, 1});//default screw is in z-direction
52  logger(DEBUG, "Screw() constructor finished.");
53 }
Logger< MERCURYDPM_LOGLEVEL > logger("MercuryKernel")
Definition of different loggers with certain modules. A user can define its own custom logger here.
@ DEBUG
void setOrientationViaNormal(Vec3D normal)
Sets the orientation of this BaseInteractable by defining the vector that results from the rotation o...
Definition: BaseInteractable.cc:199
Mdouble thickness_
The thickness of the Screw.
Definition: Screw.h:140
Vec3D start_
The centre of the lower end of the screw.
Definition: Screw.h:116
Mdouble l_
The length of the Screw.
Definition: Screw.h:120
Mdouble maxR_
The outer radius of the Screw.
Definition: Screw.h:124
Mdouble n_
The number of revelations.
Definition: Screw.h:128
ScrewType screwType_
Single or double helix screw.
Definition: Screw.h:144
Mdouble omega_
Rotation speed in rev/s.
Definition: Screw.h:132
Mdouble offset_
The angle that describes how much the Screw has turned, going from 0 to 1 for a rotation.
Definition: Screw.h:136
void setZero()
Sets all elements to zero.
Definition: Vector.cc:43

References DEBUG, doubleHelix, l_, logger, maxR_, n_, offset_, omega_, screwType_, BaseInteractable::setOrientationViaNormal(), Vec3D::setZero(), start_, and thickness_.

Referenced by copy().

◆ Screw() [2/3]

Screw::Screw ( const Screw other)

Copy constructor, copies another Screw.

Parameters
[in]otherThe Screw that has to be copied.
59  : BaseWall(other)
60 {
61  start_ = other.start_;
62  l_ = other.l_;
63  maxR_ = other.maxR_;
64  n_ = other.n_;
65  omega_ = other.omega_;
66  thickness_ = other.thickness_;
67  offset_ = other.offset_;
68  screwType_ = other.screwType_;
69  logger(DEBUG, "Screw(const Screw&) copy constructor finished.");
70 }
BaseWall()
Default constructor.
Definition: BaseWall.cc:36

References DEBUG, l_, logger, maxR_, n_, offset_, omega_, screwType_, start_, and thickness_.

◆ Screw() [3/3]

Screw::Screw ( Vec3D  start,
Mdouble  l,
Mdouble  r,
Mdouble  n,
Mdouble  omega,
Mdouble  thickness,
ScrewType  screwType = ScrewType::doubleHelix 
)

Constructor in which all parameters of the screw are set.

Parameters
[in]startA Vec3D which denotes the centre of the lower end of the Screw.
[in]lThe length of the Screw, must be positive.
[in]rThe radius of the Screw, must be positive.
[in]nThe number of revelations of the Screw, must be positive.
[in]omegaThe rotation speed of the Screw in rev/s.
[in]thicknessThe thickness of the Screw, must be non-negative.

Make a Screw by assigning all input parameters to the data-members of this class, and setting the offset_ to 0.

83 {
84  start_ = start;
85  l_ = l;
86  maxR_ = r;
87  n_ = n;
88  omega_ = omega;
89  thickness_ = thickness;
90  offset_ = 0.0;
91  screwType_ = screwType;
92  logger(DEBUG, "Screw(Vec3D, Mdouble, Mdouble, Mdouble, Mdouble, Mdouble) constructor finished.");
93 }
const unsigned n
Definition: CG3DPackingUnitTest.cpp:32

References DEBUG, l_, logger, maxR_, n, n_, offset_, omega_, screwType_, start_, and thickness_.

◆ ~Screw()

Screw::~Screw ( )
override

Default destructor.

96 {
97  logger(DEBUG, "~Screw() finished, destroyed the Screw.");
98 }

References DEBUG, and logger.

Member Function Documentation

◆ copy()

Screw * Screw::copy ( ) const
finalvirtual

Copy this screw and return a pointer to the copy.

Returns
A pointer to a copy of this Screw.

Implements BaseWall.

104 {
105  return new Screw(*this);
106 }
Screw()
Default constructor: make a screw with default parameters.
Definition: Screw.cc:41

References Screw().

◆ getDistanceAndNormal()

bool Screw::getDistanceAndNormal ( const BaseParticle P,
Mdouble distance,
Vec3D normal_return 
) const
finalvirtual

Compute the distance from the Screw for a given BaseParticle and return if there is a collision. If there is a collision, also return the normal vector of the interaction point.

Todo:
do this for all walls

Implements BaseWall.

109 {
110  //transform coordinates into position-orientation frame
111  Vec3D position = P.getPosition() - getPosition();
112  getOrientation().rotateBack(position);
115  if (getDistanceAndNormalLabCoordinates(position, P.getRadius() + s->getInteractionDistance(), distance,
116  normal_return))
117  {
118  getOrientation().rotate(normal_return);
119  return true;
120  }
121  else
122  {
123  return false;
124  }
125 }
DPMBase * getDPMBase()
Gets the problem that is solved using this handler.
Definition: BaseHandler.h:725
const Quaternion & getOrientation() const
Returns the orientation of this BaseInteractable.
Definition: BaseInteractable.h:230
const ParticleSpecies * getSpecies() const
Returns a pointer to the species of this BaseInteractable.
Definition: BaseInteractable.h:108
const Vec3D & getPosition() const
Returns the position of this BaseInteractable.
Definition: BaseInteractable.h:218
BaseSpecies is the class from which all other species are derived.
Definition: BaseSpecies.h:50
Mdouble getInteractionDistance() const
returns the largest separation distance at which adhesive short-range forces can occur.
Definition: BaseSpecies.h:146
WallHandler * getHandler() const
A function which returns the WallHandler that handles this BaseWall.
Definition: BaseWall.cc:134
SpeciesHandler speciesHandler
A handler to that stores the species type i.e. LinearViscoelasticSpecies, etc.
Definition: DPMBase.h:1427
void rotate(Vec3D &position) const
Definition: Quaternion.cc:563
void rotateBack(Vec3D &position) const
Definition: Quaternion.cc:610
bool getDistanceAndNormalLabCoordinates(Vec3D position, Mdouble wallInteractionRadius, Mdouble &distance, Vec3D &normal_return) const
Definition: Screw.cc:138
std::enable_if<!std::is_pointer< typename U::MixedSpeciesType >::value, typename U::MixedSpeciesType * >::type getMixedObject(const U *S, const U *T)
Definition: SpeciesHandler.h:74
Definition: Vector.h:51
double P
Uniform pressure.
Definition: TwenteMeshGluing.cpp:73

References getDistanceAndNormalLabCoordinates(), BaseHandler< T >::getDPMBase(), BaseWall::getHandler(), BaseSpecies::getInteractionDistance(), SpeciesHandler::getMixedObject(), BaseInteractable::getOrientation(), BaseInteractable::getPosition(), BaseInteractable::getSpecies(), Global_Physical_Variables::P, Quaternion::rotate(), Quaternion::rotateBack(), and DPMBase::speciesHandler.

◆ getDistanceAndNormalLabCoordinates()

bool Screw::getDistanceAndNormalLabCoordinates ( Vec3D  position,
Mdouble  wallInteractionRadius,
Mdouble distance,
Vec3D normal_return 
) const
Parameters
[in]pBaseParticle we want to calculate the distance and whether it collided of.
[out]distanceThe distance of the BaseParticle to this wall.
[out]normal_returnIf there was a collision, the normal vector to this wall will be placed here.
Returns
A boolean which says whether or not there was a collision.

This function computes whether or not there is a collision between a given BaseParticle and this Screw. If there is a collision, this function also computes the distance between the BaseParticle and Screw and the normal of the IntersectionOfWalls at the intersection point.

Todo:
Make this function readable and explain the steps in the details.
140 {
141  // compute the square of the radial distance (in yz-plane) between particle and screw.
142  const Mdouble RSquared = square(position.Y - start_.Y) + square(position.Z - start_.Z);
143  const Mdouble X = position.X - start_.X;
144 
145  //first do a simple check if particle is within the cylindrical hull of the screw
146  if (RSquared > square(maxR_ + wallInteractionRadius + thickness_)
147  || X > l_ + wallInteractionRadius + thickness_
148  || X < -wallInteractionRadius - thickness_)
149  {
150  return false;
151  }
152 
153  //else compute radial distance, and angle in yz-plane
154  const Mdouble R = sqrt(RSquared);
155  const Mdouble A = atan2(position.Z - start_.Z, position.Y - start_.Y);
156 
157  //after subtracting the start position and transforming the particle position from (XYZ) into (XRA)
158  //coordinates, we compute the distance to the wall at, located at (r,a,z)=(q*L,r,2*pi*(offset+N*q+k/2)), 0<q<1.
159 
160  //To find the contact point we have to minimize (with respect to r and q)
161  //distance^2=(x-x0-r*cos(2*pi*(offset+N*q)))^2+(y-y0-r*sin(2*pi*(offset+N*q)))^2+(z-z0-q*L)^2
162  //Using polar coordinates (i.e. x-x0=R*cos(A), y-y0=R*sin(A) and Z=z-z0)
163  //distance^2=R^2+r^2-2*R*r*cos(A-2*pi*(offset+N*q))+(Z-q*L)^2
164 
165  //Assumption: d(distance)/dr=0 at minDistance (should there be also a q-derivative?)
166  //Differentiate with respect to r and solve for zero:
167  //0=2*r-2*R*cos(A-2*pi*(offset+N*q)
168  //r=R*cos(A-2*pi*(offset+N*q))
169 
170  //Substitue back
171  //distance^2=R^2+R^2*cos^2(A-2*pi*(offset+N*q))-2*R^2*cos^2(A-2*pi*(offset+N*q))+(Z-q*L)^2
172  //distance^2=R^2*sin^2(A-2*pi*(offset+N*q))+(Z-q*L)^2
173 
174  //So we have to minimize:
175  //distance^2=R^2*sin^2(A-2*pi*(offset+N*q))^2 + (Z-q*L)^2 = f(q)
176  //f'(q)=(-2*pi*N)*R^2*sin(2*A-4*pi*(offset+N*q)) + 2*L*(Z-q*L) (D[Sin[x]^2,x]=Sin[2x])
177  //f''(q)=(4*pi^2*N^2)*R^2*cos(2*A-4*pi*(offset+N*q)) - 2*L*L
178  //For this we use the Euler algoritm
179 
180  Mdouble q; //Current guess
181  Mdouble dd; //Derivative at current guess
182  Mdouble ddd; //Second derivative at current guess
183 
184  //Set initial q to the closest position on the screw with the same angle as A
185  //The initial guess will be in the minimum of the sin closest to q0
186  //Minima of the sin are at
187  //A-2*pi*(offset+N*q)=k*pi (k=integer)
188  //q=A/(2*pi*N)-k/(2*N)-offset/N (k=integer)
189  {
190  const Mdouble q0 = X / l_; //assume closest point q0 is at same z-location as the particle
191  const Mdouble k = round(A / constants::pi - 2.0 * (offset_ + n_ * q0)); // k: |A-a(q0)-k*pi|=min
192  q = A / (2.0 * constants::pi * n_) - k / (2.0 * n_) - offset_ / n_; // q: a(q)=A
193  }
194 
195  //Now apply Newton's method
196  unsigned count = 0;
197  const unsigned maxCount = 30;
198  do
199  {
200  Mdouble arg = 2.0 * A - 4.0 * constants::pi * (n_ * q + offset_);
201  dd = -2.0 * constants::pi * n_ * RSquared * sin(arg) - 2.0 * l_ * (X - q * l_);
202  ddd = 8.0 * constants::sqr_pi * n_ * n_ * RSquared * cos(arg) + 2.0 * l_ * l_;
203  q -= dd / ddd;
204  if(++count>maxCount) {
205  logger(WARN,"Screw % contact detection did not converge (err=%); continuing with approximate contact point",getId(),fabs(dd/ddd));
206  break;
207  }
208  } while (fabs(dd / ddd) > 5e-14);
209 
210  //Calculate r
211  Mdouble r = R * cos(2.0 * constants::pi * (offset_ + n_ * q) - A);
212 
213  //If the particle touches the single screw center
214  if (screwType_==ScrewType::singleHelix && r>0) {
215  distance = R-thickness_;
216  if (distance>wallInteractionRadius) {
217  return false;
218  } else {
219  normal_return = Vec3D(0,-position.Y/R,-position.Z/R);
220  return true;
221  }
222  }
223 
224  //Check if the location is actually on the screw:
225  //First posibility is that the radius is too large:
226  Mdouble distanceSquared = 0;
227  if (fabs(r) > maxR_) //Left boundary of the coil
228  {
229  r = mathsFunc::sign(r) * maxR_;
230  distanceSquared = mathsFunc::square(R-maxR_);
231  count = 0;
232  //This case reduces to the coil problem
233  do
234  {
235  dd = -4.0 * R * r * constants::pi * n_ * sin(A - 2.0 * constants::pi * (n_ * q + offset_)) -
236  2.0 * l_ * (X - q * l_);
237  ddd = 8.0 * R * r * constants::sqr_pi * n_ * n_ * cos(A - 2.0 * constants::pi * (n_ * q + offset_)) +
238  2.0 * l_ * l_;
239  q -= dd / ddd;
240  if(++count>maxCount) {
241  logger(WARN,"Screw % contact detection at left boundary did not converge (err=%); continuing with approximate contact point",getId(),fabs(dd/ddd));
242  break;
243  }
244  } while (fabs(dd / ddd) > 1e-13);
245  }
246  //Second possibility is that it occured before the start of after the end
247  if (q < 0)
248  {
249  q = 0;
250  r = R * cos(A - 2.0 * constants::pi * (offset_ + q * n_));
251  if (fabs(r) > maxR_)
252  {
253  r = mathsFunc::sign(r) * maxR_;
254  }
255  }
256  else if (q > 1)
257  {
258  q = 1;
259  r = R * cos(A - 2.0 * constants::pi * (offset_ + q * n_));
260  if (fabs(r) > maxR_)
261  {
262  r = mathsFunc::sign(r) * maxR_;
263  }
264  }
265 
266  //compute squared distance between particle and contact point
267  distanceSquared += square(X - q * l_) + square(R*sin(A - 2 * constants::pi * (offset_ + n_ * q)));
268 
269  //If distance is too large there is no contact
270  if (distanceSquared >= square(wallInteractionRadius + thickness_))
271  {
272  return false;
273  }
274 
275  //Else
276  Vec3D ContactPoint;
277  distance = sqrt(distanceSquared) - thickness_;
278  ContactPoint.Y = start_.Y + r * cos(2.0 * constants::pi * (offset_ + n_ * q));
279  ContactPoint.Z = start_.Z + r * sin(2.0 * constants::pi * (offset_ + n_ * q));
280  ContactPoint.X = start_.X + q * l_;
281  normal_return = (ContactPoint - position);
282  normal_return.normalise();
283  return true;
284 }
@ WARN
@ R
Definition: StatisticsVector.h:42
@ A
Definition: StatisticsVector.h:42
@ X
Definition: StatisticsVector.h:42
unsigned int getId() const
Returns the unique identifier of any particular object.
Definition: BaseObject.h:125
Mdouble Y
Definition: Vector.h:66
Mdouble Z
Definition: Vector.h:66
Mdouble X
the vector components
Definition: Vector.h:66
void normalise()
Makes this Vec3D unit length.
Definition: Vector.cc:123
const Mdouble sqr_pi
Definition: ExtendedMath.h:47
const Mdouble pi
Definition: ExtendedMath.h:45
Mdouble round(Mdouble value, unsigned int precision)
rounds a floating point number with a given precision
Definition: MathHelpers.cc:28
Mdouble cos(Mdouble x)
Definition: ExtendedMath.cc:64
T square(const T val)
squares a number
Definition: ExtendedMath.h:106
int sign(T val)
This is a sign function, it returns -1 for negative numbers, 1 for positive numbers and 0 for 0.
Definition: ExtendedMath.h:97
Mdouble sin(Mdouble x)
Definition: ExtendedMath.cc:44

References A, mathsFunc::cos(), BaseObject::getId(), l_, logger, maxR_, n_, Vec3D::normalise(), offset_, constants::pi, R, helpers::round(), screwType_, mathsFunc::sign(), mathsFunc::sin(), singleHelix, constants::sqr_pi, mathsFunc::square(), start_, thickness_, WARN, Vec3D::X, X, Vec3D::Y, and Vec3D::Z.

Referenced by getDistanceAndNormal().

◆ getName()

std::string Screw::getName ( ) const
finalvirtual

Returns the name of the object, here the string "Screw".

Returns
The string "Screw".

Implements BaseObject.

361 {
362  return "Screw";
363 }

Referenced by writeVTK().

◆ move_time()

void Screw::move_time ( Mdouble  dt)

Rotate the Screw for a period dt, so that the offset_ changes with omega_*dt.

Parameters
[in]dtThe time for which the Screw has to be turned.
290 {
291  //offset_ += omega_ * dt;
292 }

◆ oldRead()

void Screw::oldRead ( std::istream &  is)

Reads a Screw in the old style from an input stream, for example a restart file old style.

Parameters
[in,out]isInput stream from which the Screw must be read.

Read the Screw in old style, please note that the thickness is not read in this function, so it has either to be set manually or it is 0.0 from the default constructor.

331 {
332  std::string dummy;
333  is >> dummy >> start_
334  >> dummy >> l_
335  >> dummy >> maxR_
336  >> dummy >> n_
337  >> dummy >> omega_
338  >> dummy >> offset_;
339 }

References l_, maxR_, n_, offset_, omega_, and start_.

◆ read()

void Screw::read ( std::istream &  is)
overridevirtual

Reads a Screw from an input stream, for example a restart file.

Parameters
[in,out]isInput stream from which the Screw must be read.

Reimplemented from BaseWall.

309 {
310  BaseWall::read(is);
311  std::string dummy;
312  unsigned screwType;
313  is >> dummy >> start_
314  >> dummy >> l_
315  >> dummy >> maxR_
316  >> dummy >> n_
317  >> dummy >> omega_
318  >> dummy >> thickness_
319  >> dummy >> offset_
320  >> dummy >> screwType;
321  screwType_ = static_cast<ScrewType>(screwType);
322 }
ScrewType
Definition: Screw.h:35
void read(std::istream &is) override
Function that reads a BaseWall from an input stream, usually a restart file.
Definition: BaseWall.cc:78

References l_, maxR_, n_, offset_, omega_, BaseWall::read(), screwType_, start_, and thickness_.

◆ rotate()

void Screw::rotate ( const Vec3D angularVelocityDt)
overridevirtual
Todo:
the move and rotate functions should only pass the time step, as teh velocity can be accessed directly by the object
Parameters
angularVelocityDt

Reimplemented from BaseInteractable.

299 {
300  BaseInteractable::rotate(angularVelocityDt);
302 }
virtual void rotate(const Vec3D &angularVelocityDt)
Rotates this BaseInteractable.
Definition: BaseInteractable.cc:230
Mdouble getTimeStep() const
Returns the simulation time step.
Definition: DPMBase.cc:1250

References BaseHandler< T >::getDPMBase(), BaseWall::getHandler(), DPMBase::getTimeStep(), offset_, omega_, and BaseInteractable::rotate().

◆ write()

void Screw::write ( std::ostream &  os) const
overridevirtual

Writes this Screw to an output stream, for example a restart file.

Parameters
[in,out]osOutput stream to which the Screw must be written.

Reimplemented from BaseWall.

345 {
346  BaseWall::write(os);
347  os << " start " << start_
348  << " Length " << l_
349  << " Radius " << maxR_
350  << " Revolutions " << n_
351  << " Omega " << omega_
352  << " Thickness " << thickness_
353  << " Offset " << offset_
354  << " ScrewType " << static_cast<unsigned>(screwType_);
355 }
void write(std::ostream &os) const override
Function that writes a BaseWall to an output stream, usually a restart file.
Definition: BaseWall.cc:100

References l_, maxR_, n_, offset_, omega_, screwType_, start_, thickness_, and BaseWall::write().

◆ writeVTK() [1/2]

void Screw::writeVTK ( std::string  filename) const
420 {
421  VTKContainer vtk;
422  writeVTK(vtk);
423 
424  std::stringstream file;
425  file << "# vtk DataFile Version 2.0\n"
426  << getName() << "\n"
427  "ASCII\n"
428  "DATASET UNSTRUCTURED_GRID\n"
429  "POINTS " << vtk.points.size() << " double\n";
430  for (const auto& vertex : vtk.points)
431  file << vertex << '\n';
432  file << "\nCELLS " << vtk.triangleStrips.size() << ' ' << 4 * vtk.triangleStrips.size() << "\n";
433  for (const auto& face : vtk.triangleStrips)
434  file << "3 " << face[0] << ' ' << face[1] << ' ' << face[2] << '\n';
435  file << "\nCELL_TYPES " << vtk.triangleStrips.size() << "\n";
436  for (const auto& face UNUSED : vtk.triangleStrips)
437  file << "5\n";
438  helpers::writeToFile(filename, file.str());
439 }
#define UNUSED
Definition: GeneralDefine.h:39
void writeVTK(VTKContainer &vtk) const override
Definition: Screw.cc:365
std::string getName() const final
Returns the name of the object, here the string "Screw".
Definition: Screw.cc:360
bool writeToFile(std::string filename, std::string filecontent)
Writes a string to a file.
Definition: FileIOHelpers.cc:58
Definition: BaseWall.h:38
std::vector< std::vector< double > > triangleStrips
Definition: BaseWall.h:40
std::vector< Vec3D > points
Definition: BaseWall.h:39

References getName(), VTKContainer::points, VTKContainer::triangleStrips, UNUSED, helpers::writeToFile(), and writeVTK().

◆ writeVTK() [2/2]

void Screw::writeVTK ( VTKContainer vtk) const
overridevirtual

adds extra information to the points and triangleStrips vectors needed to plot the wall in vtk format

Parameters
pointsCoordinates of the vertices of the triangulated surfaces (in the VTK file this is called POINTS)
triangleStripsIndices of three vertices forming one triangulated surface (in the VTK file this is called CELL)

Reimplemented from BaseWall.

366 {
367  //number of points in radial direction (for one side of the screw)
368  unsigned nr = 5;
369  //number of points in axial direction
370  unsigned nz = static_cast<unsigned int>(99 * fabs(n_));
371 
372  unsigned long nPoints = vtk.points.size();
373  Vec3D contactPoint;
374  // either one or two helices
375  for (Mdouble offset = offset_; offset<=offset_+(screwType_==ScrewType::doubleHelix?0.5:0); offset+=0.5)
376  {
377  for (unsigned iz = 0; iz < nz; iz++)
378  {
379  for (unsigned ir = 0; ir < nr; ir++)
380  {
381  double q = (double) iz / nz;
382  double r = (double) ir / (nr - 1) * maxR_;
383  contactPoint.Y = start_.Y - r * cos(2.0 * constants::pi * (offset + n_ * q));
384  contactPoint.Z = start_.Z - r * sin(2.0 * constants::pi * (offset + n_ * q));
385  contactPoint.X = start_.X + q * l_ - thickness_;
386  getOrientation().rotate(contactPoint);
387  contactPoint += getPosition();
388  vtk.points.push_back(contactPoint);
389  }
390  for (unsigned ir = 0; ir < nr; ir++)
391  {
392  double q = (double) iz / nz;
393  double r = (double) (nr - 1 - ir) / (nr - 1) * maxR_;
394  contactPoint.Y = start_.Y - r * cos(2.0 * constants::pi * (offset + n_ * q));
395  contactPoint.Z = start_.Z - r * sin(2.0 * constants::pi * (offset + n_ * q));
396  contactPoint.X = start_.X + q * l_ + thickness_;
397  getOrientation().rotate(contactPoint);
398  contactPoint += getPosition();
399  vtk.points.push_back(contactPoint);
400  }
401  }
402  }
403 
404  for (unsigned iz = 0; iz < (screwType_==ScrewType::doubleHelix?2:1)*nz-1; iz++)
405  {
406  //skip step that would connect the two screw parts
407  if (iz==nz-1) ++iz;
408  std::vector<double> cell;
409  cell.reserve(2 * nr);
410  for (unsigned ir = 0; ir < 2*nr; ir++)
411  {
412  cell.push_back(nPoints + ir + 2*iz * nr);
413  cell.push_back(nPoints + ir + 2*(iz + 1) * nr);
414  }
415  vtk.triangleStrips.push_back(cell);
416  }
417 }

References mathsFunc::cos(), doubleHelix, BaseInteractable::getOrientation(), BaseInteractable::getPosition(), l_, maxR_, n_, offset_, constants::pi, VTKContainer::points, Quaternion::rotate(), screwType_, mathsFunc::sin(), start_, thickness_, VTKContainer::triangleStrips, Vec3D::X, Vec3D::Y, and Vec3D::Z.

Referenced by TriangulatedScrewSelfTest::setupInitialConditions(), and writeVTK().

Member Data Documentation

◆ l_

Mdouble Screw::l_
private

The length of the Screw.

Referenced by getDistanceAndNormalLabCoordinates(), oldRead(), read(), Screw(), write(), and writeVTK().

◆ maxR_

Mdouble Screw::maxR_
private

The outer radius of the Screw.

Referenced by getDistanceAndNormalLabCoordinates(), oldRead(), read(), Screw(), write(), and writeVTK().

◆ n_

Mdouble Screw::n_
private

The number of revelations.

Referenced by getDistanceAndNormalLabCoordinates(), oldRead(), read(), Screw(), write(), and writeVTK().

◆ offset_

Mdouble Screw::offset_
private

The angle that describes how much the Screw has turned, going from 0 to 1 for a rotation.

Referenced by getDistanceAndNormalLabCoordinates(), oldRead(), read(), rotate(), Screw(), write(), and writeVTK().

◆ omega_

Mdouble Screw::omega_
private

Rotation speed in rev/s.

Referenced by oldRead(), read(), rotate(), Screw(), and write().

◆ screwType_

ScrewType Screw::screwType_
private

Single or double helix screw.

Referenced by getDistanceAndNormalLabCoordinates(), read(), Screw(), write(), and writeVTK().

◆ start_

Vec3D Screw::start_
private

The centre of the lower end of the screw.

Referenced by getDistanceAndNormalLabCoordinates(), oldRead(), read(), Screw(), write(), and writeVTK().

◆ thickness_

Mdouble Screw::thickness_
private

The thickness of the Screw.

Referenced by getDistanceAndNormalLabCoordinates(), read(), Screw(), write(), and writeVTK().


The documentation for this class was generated from the following files: