MercuryDPM  Alpha
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
AxisymmetricIntersectionOfWalls Class Reference

A AxisymmetricIntersectionOfWalls is an axisymmetric wall, defined by rotating a twodimensional IntersectionOfWalls around a symmetry axis. More...

#include <AxisymmetricIntersectionOfWalls.h>

+ Inheritance diagram for AxisymmetricIntersectionOfWalls:

Public Member Functions

 AxisymmetricIntersectionOfWalls ()
 Default constructor. More...
 
 AxisymmetricIntersectionOfWalls (const AxisymmetricIntersectionOfWalls &p)
 Copy constructor. More...
 
 AxisymmetricIntersectionOfWalls (Vec3D position, Vec3D orientation, std::vector< normalAndPosition > walls, const ParticleSpecies *species)
 Constructor setting values. More...
 
 ~AxisymmetricIntersectionOfWalls ()
 Destructor. More...
 
AxisymmetricIntersectionOfWallsoperator= (const AxisymmetricIntersectionOfWalls &other)
 Copy assignment operator. More...
 
AxisymmetricIntersectionOfWallscopy () const final
 Wall copy method. It calls the copy constructor of this Wall, useful for polymorphism. More...
 
bool getDistanceAndNormal (const BaseParticle &P, Mdouble &distance, Vec3D &normal_return) const final
 Computes the distance from the wall for a given BaseParticle and returns true if there is a collision. If there is a collision, also return the normal vector. More...
 
void read (std::istream &is) final
 reads wall More...
 
void write (std::ostream &os) const final
 outputs wall More...
 
std::string getName () const final
 Returns the name of the object. More...
 
void convertLimits (Vec3D &min, Vec3D &max) const
 
void writeVTK (VTKContainer &vtk) const override
 
- Public Member Functions inherited from IntersectionOfWalls
 IntersectionOfWalls ()
 Default constructor. More...
 
 IntersectionOfWalls (const IntersectionOfWalls &other)
 Copy constructor. More...
 
 IntersectionOfWalls (std::vector< normalAndPosition > walls, const ParticleSpecies *species)
 Constructor setting values. More...
 
virtual ~IntersectionOfWalls ()
 Destructor. More...
 
IntersectionOfWallsoperator= (const IntersectionOfWalls &other)
 
IntersectionOfWallscopy () const override
 Wall copy method. It calls the copy constructor of this Wall, useful for polymorphism. More...
 
MERCURY_DEPRECATED void clear () override
 Removes all parts of the walls. More...
 
void setSpecies (const ParticleSpecies *species)
 sets species of subwalls as well More...
 
void setHandler (WallHandler *wallHandler) override
 A function which sets the WallHandler for this BaseWall. More...
 
void addObject (Vec3D normal, Vec3D point)
 Adds a wall to the set of infinite walls, given an outward normal vector s.t. normal*x=normal*point. More...
 
MERCURY_DEPRECATED void addObject (Vec3D normal, Mdouble position)
 Adds a wall to the set of finite walls, given an outward normal vector s. t. normal*x=position. More...
 
void createOpenPrism (std::vector< Vec3D > points, Vec3D prismAxis)
 Creates an open prism which is a polygon between the points, except the first and last point, and extends infinitely in the PrismAxis direction. More...
 
void createPrism (std::vector< Vec3D > points, Vec3D prismAxis)
 Creates an open prism which is a polygon between the points and extends infinitely in the PrismAxis direction. More...
 
void createOpenPrism (std::vector< Vec3D > points)
 Creates an open prism which is a polygon between the points, except the first and last point, and extends infinitely in the direction perpendicular to the first and second wall. More...
 
void createPrism (std::vector< Vec3D > points)
 Creates an open prism which is a polygon between the points and extends infinitely in the direction perpendicular to the first and second wall. More...
 
bool getDistanceAndNormal (const BaseParticle &p, Mdouble &distance, Vec3D &normal_return) const override
 Compute the distance from the wall for a given BaseParticle and return if there is a collision. If there is a collision, also return the normal vector. More...
 
bool getDistanceAndNormal (const Vec3D &postition, Mdouble wallInteractionRadius, Mdouble &distance, Vec3D &normal_return) const
 Compute the distance from the wall for a given BaseParticle and return if there is an interaction. If there is an interaction, also return the normal vector. More...
 
void move (const Vec3D &move) override
 Move the IntersectionOfWalls to a new position, which is a Vec3D from the old position. More...
 
void read (std::istream &is) override
 Reads an IntersectionOfWalls from an input stream, for example a restart file. More...
 
void write (std::ostream &os) const override
 Writes an IntersectionOfWalls to an output stream, for example a restart file. More...
 
std::string getName () const override
 Returns the name of the object, here the string "IntersectionOfWalls". More...
 
void writeVTK (VTKContainer &vtk) const override
 
- Public Member Functions inherited from BaseWall
 BaseWall ()
 Default constructor. It makes an empty BaseWall. More...
 
 BaseWall (const BaseWall &w)
 Copy constructor. More...
 
virtual ~BaseWall ()
 Default destructor. More...
 
WallHandlergetHandler () const
 A function which returns the WallHandler that handles this BaseWall. More...
 
void setIndSpecies (unsigned int indSpecies)
 Define the species of this wall using the index of the species in the SpeciesHandler in this DPMBase. More...
 
void setSpecies (const ParticleSpecies *species)
 Define the species of this wall. More...
 
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, const Vec3D normal, const Vec3D position) const
 
virtual std::vector
< BaseInteraction * > 
getInteractionWith (BaseParticle *p, Mdouble timeStamp, InteractionHandler *interactionHandler)
 
- Public Member Functions inherited from BaseInteractable
 BaseInteractable ()
 Default BaseInteractable constructor, it simply creates an empty BaseInteractable. More...
 
 BaseInteractable (const BaseInteractable &p)
 Copy constructor. It copies the BaseInteractable and all objects it contains. More...
 
virtual ~BaseInteractable ()
 Destructor, it simply destructs the BaseInteractable and all the objects it contains. More...
 
unsigned int getIndSpecies () const
 Returns the index of the Species of this BaseInteractable. 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...
 
const Vec3DgetPosition () const
 Returns the position of this BaseInteractable. More...
 
const Vec3DgetOrientation () const
 Returns the orientation of this BaseInteractable. More...
 
void setPosition (const Vec3D &position)
 Sets the position of this BaseInteractable. More...
 
void setOrientation (const Vec3D &orientation)
 Sets the orientation of this BaseInteractable. More...
 
void rotate (const Vec3D &rotate)
 Rotates this BaseInteractable. More...
 
const std::list
< BaseInteraction * > & 
getInteractions () const
 Returns a reference to the list of interactions in this BaseInteractable. 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< Vec3D(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...
 
- Public Member Functions inherited from BaseObject
 BaseObject ()
 Default constructor. More...
 
 BaseObject (const BaseObject &p)
 Copy constructor, copies all the objects BaseObject contains. More...
 
virtual ~BaseObject ()
 virtual destructor More...
 
virtual void moveInHandler (const unsigned int index)
 Except that it is virtual, it does the same thing as setIndex() does. More...
 
void setIndex (const unsigned int index)
 Allows one to assign an index to an object in the handler/container. More...
 
void setId (const unsigned int 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...
 

Additional Inherited Members

- Static Public Member Functions inherited from BaseWall
static void addToVTK (const std::vector< Vec3D > &points, VTKContainer &vtk)
 
- Public Attributes inherited from BaseWall
WallHandlerhandler_
 
- Protected Attributes inherited from IntersectionOfWalls
std::vector< InfiniteWallwallObjects_
 The wall "segments"/directions that together make up the finite wall. More...
 

Detailed Description

A AxisymmetricIntersectionOfWalls is an axisymmetric wall, defined by rotating a twodimensional IntersectionOfWalls around a symmetry axis.

It is defined by defining an axis position p and orientation o, around which the wall is axisymmetric, and an IntersectionOfWalls object w. To determine if a particle of touches an AxisymmetricIntersectionOfWalls, the particle's location is determined is a cylindrical coordinate system $(r,\theta,z)$, where r is the radial distance from the axis, z is the axial distance along the axis, and $\theta$ the angular position around the axis (which is ignored as the object is axisymmetric). A particle touches an AxisymmetricIntersectionOfWalls if it touches the IntersectionOfWalls w in the $(r,\theta,z)$ coordinate system.

The code below defines a cylindrical outer wall of radius r with an axis position p and orientation o and contact properties s:

ParticleSpecies s = setSpecies(speciesHandler.getObject(0));
Vec3D p = Vec3D(0.0, 0.0, 0.0);
Vec3D o = Vec3D(0.0, 0.0, 1.0);
Mdouble r = 1.0;
w.addObject(Vec3D(1.0, 0.0, 0.0), Vec3D(r, 0.0, 0.0));
wallHandler.copyAndAddObject(w);

For a demonstration on how to use this class, see Flow through a 3D hourglass/silo (shown in the image below).

Definition at line 60 of file AxisymmetricIntersectionOfWalls.h.

Constructor & Destructor Documentation

AxisymmetricIntersectionOfWalls::AxisymmetricIntersectionOfWalls ( )

Default constructor.

Definition at line 31 of file AxisymmetricIntersectionOfWalls.cc.

References DEBUG, and logger.

Referenced by copy().

32 {
33  logger(DEBUG, "AxisymmetricIntersectionOfWalls() finished");
34 }
Logger< MERCURY_LOGLEVEL > logger("MercuryKernel")
AxisymmetricIntersectionOfWalls::AxisymmetricIntersectionOfWalls ( const AxisymmetricIntersectionOfWalls other)

Copy constructor.

Parameters
[in]otherThe AxisymmetricIntersectionOfWalls that must be copied.

Definition at line 39 of file AxisymmetricIntersectionOfWalls.cc.

References DEBUG, and logger.

40  : IntersectionOfWalls(other)
41 {
42  logger(DEBUG, "AxisymmetricIntersectionOfWalls(const AxisymmetricIntersectionOfWalls &p) finished");
43 }
Logger< MERCURY_LOGLEVEL > logger("MercuryKernel")
IntersectionOfWalls()
Default constructor.
AxisymmetricIntersectionOfWalls::AxisymmetricIntersectionOfWalls ( Vec3D  position,
Vec3D  orientation,
std::vector< normalAndPosition walls,
const ParticleSpecies species 
)

Constructor setting values.

Definition at line 45 of file AxisymmetricIntersectionOfWalls.cc.

References BaseInteractable::setOrientation(), and BaseInteractable::setPosition().

46  : IntersectionOfWalls(walls, species)
47 {
48  setPosition(position);
49  setOrientation(orientation);
50 }
void setOrientation(const Vec3D &orientation)
Sets the orientation of this BaseInteractable.
IntersectionOfWalls()
Default constructor.
void setPosition(const Vec3D &position)
Sets the position of this BaseInteractable.
AxisymmetricIntersectionOfWalls::~AxisymmetricIntersectionOfWalls ( )

Destructor.

Definition at line 52 of file AxisymmetricIntersectionOfWalls.cc.

References DEBUG, and logger.

53 {
54  logger(DEBUG, "~AxisymmetricIntersectionOfWalls() finished.");
55 }
Logger< MERCURY_LOGLEVEL > logger("MercuryKernel")

Member Function Documentation

void AxisymmetricIntersectionOfWalls::convertLimits ( Vec3D min,
Vec3D max 
) const

converts XYZ limits into RZ limits, to properly limit the VTK plotting area.

Definition at line 146 of file AxisymmetricIntersectionOfWalls.cc.

References Vec3D::dot(), Vec3D::getLengthSquared(), BaseInteractable::getOrientation(), BaseInteractable::getPosition(), Vec3D::X, Vec3D::Y, Z, and Vec3D::Z.

Referenced by writeVTK().

146  {
147  //define the box in the rz plane
148  double x[2] = {min.X, max.X};
149  double y[2] = {min.Y, max.Y};
150  double z[2] = {min.Z, max.Z};
151  Vec3D o = getOrientation();
152  min.X = 0;
153  Mdouble maxXSquared = 0; //to be set
154  min.Y = 0;
155  max.Y = 0.001;
156  min.Z = Vec3D::dot(max,o); //to be set
157  max.Z = min.Z; //to be set
158  for (unsigned i=0; i<2; i++)
159  for (unsigned j=0; j<2; j++)
160  for (unsigned k=0; k<2; k++) {
161  Vec3D p = Vec3D(x[i],y[j],z[k])-getPosition();
162  double Z = Vec3D::dot(p,o);
163  double XSquared = Vec3D::getLengthSquared(p-Z*o);
164  if (XSquared>maxXSquared) maxXSquared=XSquared;
165  if (Z<min.Z) min.Z=Z;
166  if (Z>max.Z) max.Z=Z;
167  }
168  max.X = sqrt(maxXSquared);
169 }
Mdouble X
the vector components
Definition: Vector.h:52
double Mdouble
static Mdouble dot(const Vec3D &a, const Vec3D &b)
Calculates the dot product of two Vec3D: .
Definition: Vector.cc:167
Mdouble getLengthSquared() const
Calculates the squared length of this Vec3D: .
Definition: Vector.cc:300
const Vec3D & getPosition() const
Returns the position of this BaseInteractable.
const Vec3D & getOrientation() const
Returns the orientation of this BaseInteractable.
Mdouble Y
Definition: Vector.h:52
Implementation of a 3D vector (by Vitaliy).
Definition: Vector.h:45
Mdouble Z
Definition: Vector.h:52
AxisymmetricIntersectionOfWalls * AxisymmetricIntersectionOfWalls::copy ( ) const
finalvirtual

Wall copy method. It calls the copy constructor of this Wall, useful for polymorphism.

Returns
pointer to a IntersectionOfWalls object allocated using new.

Implements BaseWall.

Definition at line 75 of file AxisymmetricIntersectionOfWalls.cc.

References AxisymmetricIntersectionOfWalls().

Referenced by operator=().

76 {
77  return new AxisymmetricIntersectionOfWalls(*this);
78 }
bool AxisymmetricIntersectionOfWalls::getDistanceAndNormal ( const BaseParticle p,
Mdouble distance,
Vec3D normalReturn 
) const
finalvirtual

Computes the distance from the wall for a given BaseParticle and returns true if there is a collision. If there is a collision, also return the normal vector.

First, the particle is translated by the vector position_, then the distance normal and tangential to the orientation is computed. This normal and tangential direction is interpreted as the x and z coordinate. With the particle shifted into the XZ plane, the distance and normal is computed, as if the AxisymmetricIntersectionOfWalls would be a simple IntersectionOfWalls. Finally, the object and the normal is rotated back to the original position.

See also AxisymmetricIntersectionOfWalls for details.

Implements BaseWall.

Definition at line 90 of file AxisymmetricIntersectionOfWalls.cc.

References Vec3D::dot(), IntersectionOfWalls::getDistanceAndNormal(), BaseObject::getIndex(), Vec3D::getLength(), BaseInteractable::getOrientation(), BaseInteractable::getPosition(), BaseParticle::getWallInteractionRadius(), logger, WARN, Vec3D::X, and Vec3D::Z.

91 {
92  //transform to axisymmetric coordinates
93  //move the coordinate system to the axis origin, so pOrigin=(xhat,yhat,zhat)
94  Vec3D pOrigin = p.getPosition() -getPosition();
95  Mdouble normal = Vec3D::dot(pOrigin, getOrientation());
96  //tangential is the projection into the (xhat,yhat) plane
97  Vec3D tangentialUnitVector = pOrigin - normal * getOrientation();
98  Mdouble tangential = tangentialUnitVector.getLength();
99  if (tangential!=0.0)
100  tangentialUnitVector /= tangential;
101  else //in this case the tangential vector is irrelevant
102  logger(WARN, "Warning: Particle % is exactly on the symmetry axis of wall %", p.getIndex(), getIndex());
103  Vec3D transformedPosition = Vec3D(tangential, 0.0, normal); //now P=(r,phi,zhat) is cylindrical
104  Vec3D transformedNormal;
105  //determine wall distance, normal and contact in axissymmetric coordinates
106  //and transform from axisymmetric coordinates
107  if (!IntersectionOfWalls::getDistanceAndNormal(transformedPosition, p.getWallInteractionRadius(), distance, transformedNormal))
108  {
109  //if not in contact
110  return false;
111  }
112  else
113  {
114  //if in contact
115  normalReturn = transformedNormal.Z * getOrientation() + transformedNormal.X * tangentialUnitVector;
116  return true;
117  }
118 }
unsigned int getIndex() const
Returns the index of the object in the handler.
Definition: BaseObject.cc:108
Mdouble X
the vector components
Definition: Vector.h:52
Logger< MERCURY_LOGLEVEL > logger("MercuryKernel")
double Mdouble
bool getDistanceAndNormal(const BaseParticle &p, Mdouble &distance, Vec3D &normal_return) const override
Compute the distance from the wall for a given BaseParticle and return if there is a collision...
static Mdouble dot(const Vec3D &a, const Vec3D &b)
Calculates the dot product of two Vec3D: .
Definition: Vector.cc:167
const Vec3D & getPosition() const
Returns the position of this BaseInteractable.
const Vec3D & getOrientation() const
Returns the orientation of this BaseInteractable.
static Mdouble getLength(const Vec3D &a)
Calculates the length of a Vec3D: .
Definition: Vector.cc:414
Mdouble getWallInteractionRadius() const
Returns the interaction radius for interaction with walls. See also BaseParticle::getInteractionRadiu...
Implementation of a 3D vector (by Vitaliy).
Definition: Vector.h:45
Mdouble Z
Definition: Vector.h:52
std::string AxisymmetricIntersectionOfWalls::getName ( ) const
finalvirtual

Returns the name of the object.

Returns
The string "AxisymmetricIntersectionOfWalls".

Implements BaseObject.

Definition at line 141 of file AxisymmetricIntersectionOfWalls.cc.

142 {
143  return "AxisymmetricIntersectionOfWalls";
144 }
AxisymmetricIntersectionOfWalls & AxisymmetricIntersectionOfWalls::operator= ( const AxisymmetricIntersectionOfWalls other)

Copy assignment operator.

Parameters
[in]otherThe AxisymmetricIntersectionOfWalls that must be copied.

Definition at line 60 of file AxisymmetricIntersectionOfWalls.cc.

References copy().

61 {
62  if (this == &other)
63  {
64  return *this;
65  }
66  else
67  {
68  return *(other.copy());
69  }
70 }
AxisymmetricIntersectionOfWalls * copy() const final
Wall copy method. It calls the copy constructor of this Wall, useful for polymorphism.
void AxisymmetricIntersectionOfWalls::read ( std::istream &  is)
finalvirtual

reads wall

Parameters
[in]isThe input stream from which the AxisymmetricIntersectionOfWalls is read, usually a restart file.

Implements BaseInteractable.

Definition at line 124 of file AxisymmetricIntersectionOfWalls.cc.

References IntersectionOfWalls::read().

125 {
127 }
void read(std::istream &is) override
Reads an IntersectionOfWalls from an input stream, for example a restart file.
void AxisymmetricIntersectionOfWalls::write ( std::ostream &  os) const
finalvirtual

outputs wall

Parameters
[in]osThe output stream where the AxisymmetricIntersectionOfWalls must be written to, usually a restart file.

Implements BaseInteractable.

Definition at line 133 of file AxisymmetricIntersectionOfWalls.cc.

References IntersectionOfWalls::write().

134 {
136 }
void write(std::ostream &os) const override
Writes an IntersectionOfWalls to an output stream, for example a restart file.
void AxisymmetricIntersectionOfWalls::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)
Bug:
once the quaternions are implemented, we can orient these walls properly

Reimplemented from BaseWall.

Definition at line 171 of file AxisymmetricIntersectionOfWalls.cc.

References convertLimits(), mathsFunc::cos(), BaseHandler< T >::getDPMBase(), BaseWall::getHandler(), DPMBase::getMax(), DPMBase::getMin(), BaseInteractable::getPosition(), BaseWall::intersectVTK(), constants::pi, VTKContainer::points, mathsFunc::sin(), VTKContainer::triangleStrips, IntersectionOfWalls::wallObjects_, X, Vec3D::X, XY, Y, Vec3D::Y, and Vec3D::Z.

172 {
173  for (auto wall=wallObjects_.begin(); wall!=wallObjects_.end(); wall++)
174  {
175  //plot each of the intersecting walls
176  std::vector<Vec3D> myPoints;
177 
178  //first create a slice of non-rotated wall in the xz plane, 0<y<1
179  Vec3D min = getHandler()->getDPMBase()->getMin();
180  Vec3D max = getHandler()->getDPMBase()->getMax();
181  convertLimits(min, max);
182 
183  //create the basic slice for the first wall using the InfiniteWall routine
184  wall->createVTK (myPoints,min,max);
185 
186  //create intersections with the other walls, similar to the IntersectionOfWalls routine
187  for (auto other=wallObjects_.begin(); other!=wallObjects_.end(); other++)
188  {
189  if (other!=wall)
190  {
191  intersectVTK(myPoints, -other->getNormal(), other->getPosition());
192  }
193  }
194 
195  //only keep the y=0 values
196  std::vector<Vec3D> rzVec;
197  for (auto& p: myPoints)
198  {
199  if (p.Y==0)
200  {
201  rzVec.push_back(p);
202  }
203  }
204  if (rzVec.size()==0)
205  return;
206 
207  //create 60 points on the unit circle
208  unsigned nr = 60;
209  struct XY {
210  double X;
211  double Y;
212  };
213  std::vector<XY> xyVec;
214  for (unsigned ir=0; ir<nr; ir++) {
215  Mdouble angle = 2.0 * constants::pi * ir/nr;
216  xyVec.push_back({cos(angle),sin(angle)});
217  }
218 
219  //now create rings of points on the axisym. shape
221  unsigned nPoints = vtk.points.size();
222  Vec3D p;
223  for (auto rz : rzVec)
224  {
225  for (auto xy : xyVec)
226  {
227  p = getPosition();
228  p.X += rz.X * xy.X;
229  p.Y += rz.X * xy.Y;
230  p.Z += rz.Z;
231  vtk.points.push_back(p);
232  }
233  }
234 
235  //finally create the connectivity matri to plot shell-like triangle strips.
236  unsigned nz = rzVec.size();
237  unsigned nCells = vtk.triangleStrips.size();
238  vtk.triangleStrips.reserve(nCells+(nz-1));
239  for (unsigned iz=0; iz<nz-1; iz++) {
240  std::vector<double> cell;
241  cell.reserve(2*nr+2);
242  for (unsigned ir=0; ir<nr; ir++) {
243  cell.push_back(nPoints+ir+iz*nr);
244  cell.push_back(nPoints+ir+(iz+1)*nr);
245  }
246  cell.push_back(nPoints+iz*nr);
247  cell.push_back(nPoints+(iz+1)*nr);
248  vtk.triangleStrips.push_back(cell);
249  }
250  }
251 
252 }
Mdouble X
the vector components
Definition: Vector.h:52
void convertLimits(Vec3D &min, Vec3D &max) const
Vec3D getMin() const
Return the "bottom left" corner of the domain, a vector with xMin_, yMin_ and zMin_.
Definition: DPMBase.h:330
double Mdouble
void intersectVTK(std::vector< Vec3D > &points, const Vec3D normal, const Vec3D position) const
Definition: BaseWall.cc:163
const Vec3D & getPosition() const
Returns the position of this BaseInteractable.
std::vector< InfiniteWall > wallObjects_
The wall "segments"/directions that together make up the finite wall.
Mdouble cos(Mdouble x)
Definition: ExtendedMath.cc:60
WallHandler * getHandler() const
A function which returns the WallHandler that handles this BaseWall.
Definition: BaseWall.cc:85
Mdouble sin(Mdouble x)
Definition: ExtendedMath.cc:42
const Mdouble pi
Definition: ExtendedMath.h:42
std::vector< std::vector< double > > triangleStrips
Definition: BaseWall.h:36
Vec3D getMax() const
Return the "upper right" corner of the domain, a vector with xMin_, yMin_ and zMin_.
Definition: DPMBase.h:335
Mdouble Y
Definition: Vector.h:52
std::vector< Vec3D > points
Definition: BaseWall.h:35
Implementation of a 3D vector (by Vitaliy).
Definition: Vector.h:45
DPMBase * getDPMBase()
Gets the problem that is solved using this handler.
Definition: BaseHandler.h:543
Mdouble Z
Definition: Vector.h:52

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