A infinite wall fills the half-space {point: (position_-point)*normal_<=0}. More...

#include <LevelSetWall.h>

Inheritance diagram for LevelSetWall:

Public Types
enum class	Shape { Sphere , Cube , Diamond , FourSided , Cylinder }

Public Member Functions
	LevelSetWall (Shape s, double radius, ParticleSpecies *sp=nullptr)

	~LevelSetWall () override
	Default destructor. More...

LevelSetWall *	copy () const override
	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 override

void	read (std::istream &is) override
	Reads LevelSetWall from a restart file. More...

std::string	getName () const override
	Returns the name of the object, in this case the string "LevelSetWall". More...

void	writeVTK (VTKContainer &vtk) const override

void	writeToFile (int n, double radiusContact) const

Public Member Functions inherited from BaseWall
	BaseWall ()
	Default constructor. More...

	BaseWall (const BaseWall &w)
	Copy constructor. More...

	~BaseWall () override
	Default destructor. More...

void	write (std::ostream &os) const override
	Function that writes a BaseWall to an output stream, usually a restart file. 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...

WallHandler *	getHandler () 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 BaseInteraction *	getInteractionWithSuperQuad (SuperQuadricParticle p, unsigned timeStamp, InteractionHandler interactionHandler)

void	getVTK (std::vector< Vec3D > &points, std::vector< std::vector< double >> &triangleStrips)

const Vec3D	getAxis () const

BaseInteraction *	getInteractionWith (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)

BaseWall *	getRenderedWall (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 ParticleSpecies *	getSpecies () const
	Returns a pointer to the species of this BaseInteractable. More...

void	setSpecies (const ParticleSpecies *species)
	Sets the species of this BaseInteractable. More...

const Vec3D &	getForce () const
	Returns the force on this BaseInteractable. More...

const Vec3D &	getTorque () 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 Vec3D &	getPosition () const
	Returns the position of this BaseInteractable. More...

const Quaternion &	getOrientation () 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...

virtual void	rotate (const Vec3D &angularVelocityDt)
	Rotates this BaseInteractable. 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 Vec3D &	getVelocity () const
	Returns the velocity of this interactable. More...

virtual const Vec3D &	getAngularVelocity () 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 Member Functions
bool	getDistanceAndNormalLabCoordinates (Vec3D position, Mdouble interactionRadius, Mdouble &distance, Vec3D &normal) const

void	setShapeSphere ()

void	setShapeCube ()

void	setShapeDiamond ()

void	setShapeCylinder ()

void	setShapeFourSided ()

void	createVTKSphere ()

void	createVTKCube ()

void	createVTKDiamond ()

void	createVTK ()

Private Attributes
double	levelSet_ [2 N+1][2 N+1][2 *N+1]

double	radius_ = 1

VTKContainer	vtkLabFrame_

Static Private Attributes
static const int	N = 10

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

A infinite wall fills the half-space {point: (position_-point)*normal_<=0}.

This is a class defining walls. It defines the interaction of regular walls and periodic walls with particles as defined in Particle Modifications:

Thus, the surface of the wall is a plane through position position_ with normal_ the outward unit normal vector of the wall (pointing away from the particles, into the wall). Please note that this wall is infinite and straight.

A particle touches an infinite wall if (position_-point)*normal_<=radius.

Member Enumeration Documentation

◆ Shape

enum LevelSetWall::Shape

strong

Enumerator
Sphere
Cube
Diamond
FourSided
Cylinder

     {
         Sphere,
         Cube,
         Diamond,
         FourSided,
         Cylinder
     };

Constructor & Destructor Documentation

◆ LevelSetWall()

LevelSetWall::LevelSetWall	(	Shape	s,
		double	radius,
		ParticleSpecies *	sp = `nullptr`
	)

                                                                       : radius_(radius)
 {
     setSpecies(sp);
     if (s == Shape::Sphere)
     {
         setShapeSphere();
         createVTKSphere();
     }
     else if (s == Shape::Cube)
     {
         setShapeCube();
         createVTKCube();
     }
     else if (s == Shape::Diamond)
     {
         setShapeDiamond();
         createVTK();
     }
     else if (s == Shape::Cylinder)
     {
         setShapeCylinder();
         createVTK();
     }
     else if (s == Shape::FourSided)
     {
         setShapeFourSided();
         createVTK();
     }
     else
     {
         logger(ERROR, "Shape unknown");
     }
 }

References createVTK(), createVTKCube(), createVTKSphere(), Cube, Cylinder, Diamond, ERROR, FourSided, logger, s, setShapeCube(), setShapeCylinder(), setShapeDiamond(), setShapeFourSided(), setShapeSphere(), BaseWall::setSpecies(), and Sphere.

Referenced by copy().

◆ ~LevelSetWall()

LevelSetWall::~LevelSetWall ( )

override

Default destructor.

 {
     logger(DEBUG, "LevelSetWall::~LevelSetWall finished");
 }

References DEBUG, and logger.

Member Function Documentation

◆ copy()

LevelSetWall * LevelSetWall::copy ( ) const

overridevirtual

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

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

Implements BaseWall.

 {
     return new LevelSetWall(*this);
 }

References LevelSetWall().

◆ createVTK()

void LevelSetWall::createVTK ( )

private

 {
     Mdouble distance;
     Vec3D normal;
     InfiniteWall wall;
  
     Mdouble interactionRadius = radius_/(2*N);
     for (int i=-N; i<N; ++i) {
         for (int j=-N; j<N; ++j) {
             for (int k=-N; k<N; ++k) {
                 Vec3D min = interactionRadius*Vec3D(2*i,2*j,2*k);
                 Vec3D max = min + interactionRadius*Vec3D(2,2,2);
                 Vec3D position = min + interactionRadius*Vec3D(1,1,1);
                 if (getDistanceAndNormalLabCoordinates(position, 2*interactionRadius, distance, normal)) {
                     wall.set(normal,position+distance*normal);
                     std::vector<Vec3D> points;
                     wall.createVTK(points,min,max);
                     if (!points.empty()) {
                         wall.addToVTK(points, vtkLabFrame_);
                     }
                 }
             }
         }
     }
 }

References BaseWall::addToVTK(), InfiniteWall::createVTK(), getDistanceAndNormalLabCoordinates(), i, j, k, max, min, N, WallFunction::normal(), radius_, InfiniteWall::set(), and vtkLabFrame_.

Referenced by LevelSetWall().

◆ createVTKCube()

void LevelSetWall::createVTKCube ( )

private

 {
     Mdouble r = radius_ * (N - 1) / N;
     vtkLabFrame_.points.resize(8);
     vtkLabFrame_.points[0] = r * Vec3D(-1, -1, -1);
     vtkLabFrame_.points[1] = r * Vec3D(-1, -1, 1);
     vtkLabFrame_.points[2] = r * Vec3D(-1, 1, -1);
     vtkLabFrame_.points[3] = r * Vec3D(-1, 1, 1);
     vtkLabFrame_.points[4] = r * Vec3D(1, -1, -1);
     vtkLabFrame_.points[5] = r * Vec3D(1, -1, 1);
     vtkLabFrame_.points[6] = r * Vec3D(1, 1, -1);
     vtkLabFrame_.points[7] = r * Vec3D(1, 1, 1);
     vtkLabFrame_.triangleStrips.resize(6);
     vtkLabFrame_.triangleStrips[0] = {0, 1, 2, 3};
     vtkLabFrame_.triangleStrips[1] = {4, 5, 6, 7};
     vtkLabFrame_.triangleStrips[2] = {0, 1, 4, 5};
     vtkLabFrame_.triangleStrips[3] = {2, 3, 6, 7};
     vtkLabFrame_.triangleStrips[4] = {0, 2, 4, 6};
     vtkLabFrame_.triangleStrips[5] = {1, 3, 5, 7};
 }

References N, VTKContainer::points, UniformPSDSelfTest::r, radius_, VTKContainer::triangleStrips, and vtkLabFrame_.

Referenced by LevelSetWall().

◆ createVTKDiamond()

void LevelSetWall::createVTKDiamond ( )

private

 {
     vtkLabFrame_.points.resize(6);
     vtkLabFrame_.points[0] = radius_ * Vec3D(-1, 0, 0);
     vtkLabFrame_.points[1] = radius_ * Vec3D(1, 0, 0);
     vtkLabFrame_.points[2] = radius_ * Vec3D(0, -1, 0);
     vtkLabFrame_.points[3] = radius_ * Vec3D(0, 1, 0);
     vtkLabFrame_.points[4] = radius_ * Vec3D(0, 0, -1);
     vtkLabFrame_.points[5] = radius_ * Vec3D(0, 0, 1);
     vtkLabFrame_.triangleStrips.resize(8);
     vtkLabFrame_.triangleStrips[0] = {0, 3, 4};
     vtkLabFrame_.triangleStrips[1] = {0, 3, 5};
     vtkLabFrame_.triangleStrips[2] = {0, 2, 4};
     vtkLabFrame_.triangleStrips[3] = {0, 2, 5};
     vtkLabFrame_.triangleStrips[4] = {1, 3, 4};
     vtkLabFrame_.triangleStrips[5] = {1, 3, 5};
     vtkLabFrame_.triangleStrips[6] = {1, 2, 4};
     vtkLabFrame_.triangleStrips[7] = {1, 2, 5};
 }

References VTKContainer::points, radius_, VTKContainer::triangleStrips, and vtkLabFrame_.

◆ createVTKSphere()

void LevelSetWall::createVTKSphere ( )

private

 {
     vtkLabFrame_.points.clear();
     vtkLabFrame_.triangleStrips.clear();
     
     const unsigned nr = 40;
     const unsigned nz = 40;
     
     std::array<Mdouble, nr> s, c;
     for (unsigned i = 0; i < nr; ++i)
     {
         s[i] = sin(2.0 * constants::pi * i / nr);
         c[i] = cos(2.0 * constants::pi * i / nr);
     }
     std::array<Mdouble, nz> h, r;
     for (unsigned j = 0; j < nz; ++j)
     {
         r[j] = radius_ * sin(2.0 * constants::pi * j / nz);
         h[j] = radius_ * cos(2.0 * constants::pi * j / nz);
     }
     
     for (unsigned j = 0; j < nz; ++j)
     {
         for (unsigned i = 0; i < nr; ++i)
         {
             vtkLabFrame_.points.emplace_back(r[j] * s[i], r[j] * c[i], h[j]);
         }
     }
     
     vtkLabFrame_.triangleStrips.reserve(nz - 1);
     for (unsigned j = 0; j < nz - 1; ++j)
     {
         std::vector<double> cell;
         cell.reserve(2 * nr + 2);
         for (unsigned i = 0; i < nr; ++i)
         {
             cell.push_back(i + j * nr);
             cell.push_back(i + (j + 1) * nr);
         }
         cell.push_back(j * nr);
         cell.push_back((j + 1) * nr);
         vtkLabFrame_.triangleStrips.push_back(cell);
     }
 }

References calibrate::c, cos(), i, j, Mesh_Parameters::nz, constants::pi, VTKContainer::points, UniformPSDSelfTest::r, radius_, s, sin(), VTKContainer::triangleStrips, and vtkLabFrame_.

Referenced by LevelSetWall().

◆ getDistanceAndNormal()

bool LevelSetWall::getDistanceAndNormal	(	const BaseParticle &	p,
		Mdouble &	distance,
		Vec3D &	normal_return
	)		const

overridevirtual

     \brief Returns the distance of the wall to the particle.
    &zwj;/

Mdouble getDistance(Vec3D position) const;

/*!
   \brief 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.

 \param[in] otherPosition The position to which the distance must be computed to.
 \return The distance of the wall to the particle.
&zwj;/

Mdouble LevelSetWall::getDistance(Vec3D otherPosition) const { return getOrientation().getDistance(otherPosition,getPosition()); }

/*!

Parameters

[in]	p	BaseParticle for which the distance to the wall must be computed.
[out]	distance	Distance between the particle and the wall.
[out]	normal_return	The normal of this wall, will only be set if there is a collision.

Returns: A boolean value for whether or not there is a collision.

First the distance is checked. If there is no collision, this function will return false and give the distance. If there is a collision, the function will return true and give the distance and the normal vector of this wall. Since this function should be called before calculating any Particle-Wall interactions, it can also be used to set the normal vector in case of curved walls.

Todo:: do this for all walls

Implements BaseWall.

 {
     //transform coordinates into position-orientation frame
     Vec3D position = p.getPosition() - getPosition();
     getOrientation().rotateBack(position);
     if (getDistanceAndNormalLabCoordinates(position, p.getWallInteractionRadius(this), distance, normal_return))
     {
         getOrientation().rotate(normal_return);
         return true;
     }
     return false;
     //return getDistanceAndNormal(p.getPosition(), distance, normal_return, p.getRadius(), p.getInteractionRadius());
 }

References getDistanceAndNormalLabCoordinates(), BaseInteractable::getOrientation(), BaseInteractable::getPosition(), p, Quaternion::rotate(), and Quaternion::rotateBack().

◆ getDistanceAndNormalLabCoordinates()

bool LevelSetWall::getDistanceAndNormalLabCoordinates	(	Vec3D	position,
		Mdouble	interactionRadius,
		Mdouble &	distance,
		Vec3D &	normal
	)		const

private

 {
     //scale values [-radius,radius] to [0,2*N]
     Vec3D positionScaled = position / radius_ * static_cast<double>(N) + Vec3D(1, 1, 1) * static_cast<double>(N);
     // index of level set value used for interpolation (the discrete level-set value left of the interpolation point)
     int i = std::max(std::min((int) positionScaled.X, 2 * N - 1), 0);
     int j = std::max(std::min((int) positionScaled.Y, 2 * N - 1), 0);
     int k = std::max(std::min((int) positionScaled.Z, 2 * N - 1), 0);
     // scaled difference between position and discrete level-set value, in [0,1) for interpolation
     double x = positionScaled.X - i;
     double y = positionScaled.Y - j;
     double z = positionScaled.Z - k;
     // trilinear interpolation
     distance = levelSet_[i][j][k] * (1 - x) * (1 - y) * (1 - z)
                + levelSet_[i + 1][j][k] * x * (1 - y) * (1 - z)
                + levelSet_[i][j + 1][k] * (1 - x) * y * (1 - z)
                + levelSet_[i + 1][j + 1][k] * x * y * (1 - z)
                + levelSet_[i][j][k + 1] * (1 - x) * (1 - y) * z
                + levelSet_[i + 1][j][k + 1] * x * (1 - y) * z
                + levelSet_[i][j + 1][k + 1] * (1 - x) * y * z
                + levelSet_[i + 1][j + 1][k + 1] * x * y * z;
     if (distance < interactionRadius)
     {
         normal.X = -levelSet_[i][j][k] * (1 - y) * (1 - z)
                    + levelSet_[i + 1][j][k] * (1 - y) * (1 - z)
                    - levelSet_[i][j + 1][k] * y * (1 - z)
                    + levelSet_[i + 1][j + 1][k] * y * (1 - z)
                    - levelSet_[i][j][k + 1] * (1 - y) * z
                    + levelSet_[i + 1][j][k + 1] * (1 - y) * z
                    - levelSet_[i][j + 1][k + 1] * y * z
                    + levelSet_[i + 1][j + 1][k + 1] * y * z;
         normal.Y = -levelSet_[i][j][k] * (1 - x) * (1 - z)
                    - levelSet_[i + 1][j][k] * x * (1 - z)
                    + levelSet_[i][j + 1][k] * (1 - x) * (1 - z)
                    + levelSet_[i + 1][j + 1][k] * x * (1 - z)
                    - levelSet_[i][j][k + 1] * (1 - x) * z
                    - levelSet_[i + 1][j][k + 1] * x * z
                    + levelSet_[i][j + 1][k + 1] * (1 - x) * z
                    + levelSet_[i + 1][j + 1][k + 1] * x * z;
         normal.Z = -levelSet_[i][j][k] * (1 - x) * (1 - y)
                    - levelSet_[i + 1][j][k] * x * (1 - y)
                    - levelSet_[i][j + 1][k] * (1 - x) * y
                    - levelSet_[i + 1][j + 1][k] * x * y
                    + levelSet_[i][j][k + 1] * (1 - x) * (1 - y)
                    + levelSet_[i + 1][j][k + 1] * x * (1 - y)
                    + levelSet_[i][j + 1][k + 1] * (1 - x) * y
                    + levelSet_[i + 1][j + 1][k + 1] * x * y;
         Mdouble length = normal.getLength();
         if (length != 0)
             normal /= -length;
         else
             normal = {1,0,0};
         return true;
     }
     return false;
 }

References i, j, k, levelSet_, max, min, N, WallFunction::normal(), radius_, plotDoE::x, Vec3D::X, y, Vec3D::Y, and Vec3D::Z.

Referenced by createVTK(), getDistanceAndNormal(), and writeToFile().

◆ getName()

std::string LevelSetWall::getName ( ) const

overridevirtual

Returns the name of the object, in this case the string "LevelSetWall".

Returns: The string "LevelSetWall", which is the name of this class.

Implements BaseObject.

 {
     return "LevelSetWall";
 }

◆ read()

void LevelSetWall::read ( std::istream & is )

overridevirtual

Reads LevelSetWall from a restart file.

Parameters

[in] is The input stream from which the LevelSetWall is read. Only needed for backward compatibility.

Todo:

Reimplemented from BaseWall.

 {
     BaseWall::read(is);
     std::string dummy;
     is >> dummy;
 }

References BaseWall::read(), and oomph::Global_string_for_annotation::string().

◆ setShapeCube()

void LevelSetWall::setShapeCube ( )

private

 {
     logger(INFO, "Creating a cube");
     radius_ *= N / (N - 1);
     Mdouble h = radius_ / N;
     for (int k = -N; k <= N; ++k)
     {
         for (int j = -N; j <= N; ++j)
         {
             for (int i = -N; i <= N; ++i)
             {
                 std::array<int, 3> sort;
                 sort[0] = abs(i);
                 sort[1] = abs(j);
                 sort[2] = abs(k);
                 std::sort(sort.begin(), sort.end());
                 if (sort[1] < N)
                 {
                     // face
                     levelSet_[i + N][j + N][k + N] = h * (sort[2] - (N - 1));
                 }
                 else if (sort[0] < N)
                 {
                     // edge
                     levelSet_[i + N][j + N][k + N] = h * 2;
                 }
                 else
                 {
                     //vertex
                     levelSet_[i + N][j + N][k + N] = h * 3;
                 }
                 //edge and vertex distances longer than real, otherwise trilinear extrapolation is problematic.
             }
         }
     }
 }

References abs(), i, INFO, j, k, levelSet_, logger, N, and radius_.

Referenced by LevelSetWall().

◆ setShapeCylinder()

void LevelSetWall::setShapeCylinder ( )

private

 {
     logger(INFO, "Creating a cylinder");
     for (int k = -N; k <= N; ++k)
     {
         for (int j = -N; j <= N; ++j)
         {
             for (int i = -N; i <= N; ++i)
             {
                 levelSet_[i + N][j + N][k + N] = radius_*std::fmax(sqrt(i * i + j * j) / N - 1,fabs(k)/N-1);
             }
         }
     }
 }

References boost::multiprecision::fabs(), Eigen::bfloat16_impl::fmax(), i, INFO, j, k, levelSet_, logger, N, radius_, and sqrt().

Referenced by LevelSetWall().

◆ setShapeDiamond()

void LevelSetWall::setShapeDiamond ( )

private

 {
     logger(INFO, "Creating a diamond square");
     for (int k = -N; k <= N; ++k)
     {
         for (int j = -N; j <= N; ++j)
         {
             for (int i = -N; i <= N; ++i)
             {
                 levelSet_[i + N][j + N][k + N] = radius_ / sqrt(3) / N * ((double) (abs(i) + abs(j) + abs(k)) - N);
             }
         }
     }
 }

References abs(), i, INFO, j, k, levelSet_, logger, N, radius_, and sqrt().

Referenced by LevelSetWall().

◆ setShapeFourSided()

void LevelSetWall::setShapeFourSided ( )

private

 {
     //open filestream
     std::string fileName = "fourSided.txt";
     std::ifstream is(fileName.c_str(), std::ios::in);
     if (is.fail())
     {
         helpers::writeToFile("fourSided.m",
                 "%% parameters\n"
                 "n=10; % 2n+1 is the number of cells in the level set per dimension\n"
                 "w=1; % [-w,w] is the width of the shape in z-direction\n"
                 "% define a few points on the xy-cross section of the shape that will be interpolated\n"
                 "% we define a quarter of the cross section only; thus the first point has differ \n"
                 "% from the last point by a quarter rotation around the origin\n"
                 "x0=[1 1.8 1.8 1.4 1];\n"
                 "y0=[-1 0 1 1.4 1];\n"
                 "\n"
                 "%% create curve, and plot\n"
                 "% define a parameter t0, such that the parametrisation is x0(t0), y0(t0)\n"
                 "t0=1:length(x0);\n"
                 "% define a more fine-grained parameter t1, and interpolate, obtaining a \n"
                 "% finer parametrisation x1(t1), y1(t1)\n"
                 "t1= linspace(1,length(x0));\n"
                 "x1=interpn(t0,x0,t1,'spline');\n"
                 "y1=interpn(t0,y0,t1,'spline');\n"
                 "% complete the shape by repeating it four times, each rotated by 90 deg\n"
                 "% the points (x,y) are on the shape \n"
                 "x=[x1 -y1 -x1 y1];\n"
                 "y=[y1 x1 -y1 -x1];\n"
                 "% plot the resulting shape\n"
                 "figure(1)\n"
                 "plot(x,y,'-',x0,y0,'o')\n"
                 "axis equal\n"
                 "\n"
                 "%% create level set, and plot\n"
                 "% compute the size of the domain needed to contain the shape ...\n"
                 "radius = max([max(abs(x)),max(abs(y)),abs(w)]);\n"
                 "% .. then divide by the number of cells to get the cell size\n"
                 "radius = n/floor(n/radius);\n"
                 "% create xLS,yLS,zLS values of the level set function dist(xLS,yLS,zLS)\n"
                 "LS = linspace(-radius,radius,2*n+1);\n"
                 "[xLS,yLS,zLS]=meshgrid(LS);\n"
                 "% define the level set function by computing the distance of (xLS,yLS,zLS)\n"
                 "% to the nearest point on the shape (x,y,+-w)\n"
                 "dist = zeros(size(xLS));\n"
                 "for i=1:length(xLS(:))\n"
                 "   dist(i) = sqrt(min((x-xLS(i)).^2+(y-yLS(i)).^2));\n"
                 "   % make it a signed distance\n"
                 "   if inpolygon(xLS(i),yLS(i),x,y)\n"
                 "       dist(i) = -dist(i);\n"
                 "   end\n"
                 "   % compute the signed distance as the max of signed distances in xy and z\n"
                 "   dist(i) = max(dist(i),abs(zLS(i))-w);\n"
                 "end\n"
                 "% rescale the distance\n"
                 "dist = dist/radius;\n"
                 "\n"
                 "% plot the i-th cross-section in xy (i=n+1 is in the center)\n"
                 "figure(2)\n"
                 "i=n+4\n"
                 "mesh(xLS(:,:,i),yLS(:,:,i),dist(:,:,i),'FaceColor','none')\n"
                 "hold on\n"
                 "contourf(xLS(:,:,i),yLS(:,:,i),dist(:,:,i),[0 0])\n"
                 "plot(x0,y0,'ro')\n"
                 "hold off\n"
                 "view(0,90)\n"
                 "axis equal\n"
                 "title(['z=' num2str(unique(zLS(:,:,i)))])\n"
                 "xlabel('x'), ylabel('y'), zlabel('dist')\n"
                 "\n"
                 "% plot the i-th cross-section in yz (i=n+1 is in the center)\n"
                 "figure(3)\n"
                 "mesh(squeeze(yLS(:,i,:)),squeeze(zLS(:,i,:)),squeeze(dist(:,i,:)),'FaceColor','none')\n"
                 "hold on\n"
                 "contourf(squeeze(yLS(:,i,:)),squeeze(zLS(:,i,:)),squeeze(dist(:,i,:)),[0 0])\n"
                 "hold off\n"
                 "view(90,0)\n"
                 "axis equal\n"
                 "title(['x=' num2str(unique(xLS(:,i,:)))])\n"
                 "xlabel('y'), ylabel('z'), zlabel('dist')\n"
                 "\n"
                 "% plot the i-th cross-section in xz (i=n+1 is in the center)\n"
                 "figure(4)\n"
                 "mesh(squeeze(xLS(i,:,:)),squeeze(zLS(i,:,:)),squeeze(dist(i,:,:)),'FaceColor','none')\n"
                 "hold on\n"
                 "contourf(squeeze(xLS(i,:,:)),squeeze(zLS(i,:,:)),squeeze(dist(i,:,:)),[0 0])\n"
                 "hold off\n"
                 "view(90,0)\n"
                 "axis equal\n"
                 "title(['y=' num2str(unique(yLS(i,:,:)))])\n"
                 "xlabel('x'), ylabel('z'), zlabel('dist')\n"
                 "\n"
                 "%% write n and dist values to file (read into MercuryDPM)\n"
                 "fid = fopen('fourSided.txt','w');\n"
                 "fprintf(fid,'%d\\n',n);\n"
                 "fprintf(fid,'%f ',dist);\n"
                 "fclose(fid);");
         logger(ERROR, "readFromFile: file % could not be opened; make sure you run fourSided.m first.", fileName);
     }
     
     int n;
     is >> n;
     if (n != N)
     {
         logger(ERROR, "readFromFile: level set size does not match %", n);
     }
     
     Mdouble unitDistance;
     for (int k = -N; k <= N; ++k)
     {
         for (int j = -N; j <= N; ++j)
         {
             for (int i = -N; i <= N; ++i)
             {
                 is >> unitDistance;
                 levelSet_[i + N][j + N][k + N] = radius_ * unitDistance;
             }
         }
     }
 }

References ERROR, UniformPSDSelfTest::fileName, i, j, k, levelSet_, logger, n, N, radius_, oomph::Global_string_for_annotation::string(), and helpers::writeToFile().

Referenced by LevelSetWall().

◆ setShapeSphere()

void LevelSetWall::setShapeSphere ( )

private

 {
     logger(INFO, "Creating a sphere");
     for (int k = -N; k <= N; ++k)
     {
         for (int j = -N; j <= N; ++j)
         {
             for (int i = -N; i <= N; ++i)
             {
                 levelSet_[i + N][j + N][k + N] = radius_ * sqrt(i * i + j * j + k * k) / N - radius_;
             }
         }
     }
 }

References i, INFO, j, k, levelSet_, logger, N, radius_, and sqrt().

Referenced by LevelSetWall().

◆ writeToFile()

void LevelSetWall::writeToFile	(	int	n,
		double	radiusContact
	)		const

Parameters

n	how many points should be used to interpolate
radiusContact	how much should be added to the radius

 {
     std::stringstream s;
     s << "Values of level set function\n";
     Vec3D normal;
     Mdouble distance;
     double meshSize = (radius_ + extraRadius) / n;
     for (int k = -n; k <= n; ++k)
     {
         for (int j = -n; j <= n; ++j)
         {
             for (int i = -n; i <= n; ++i)
             {
                 Vec3D position = meshSize * Vec3D(i, j, k);
                 getDistanceAndNormalLabCoordinates(position, radius_ + extraRadius, distance, normal);
                 s << position << ' '
                   << distance << ' '
                   << normal << '\n';
             }
         }
     }
     helpers::writeToFile("LevelSet.txt", s.str());
     helpers::writeToFile("LevelSet.m", "close all\n"
                                        "data=importdata('LevelSet.txt');\n"
                                        "N = nthroot(size(data.data,1),3);\n"
                                        "x = reshape(data.data(:,1),N,N,N);\n"
                                        "y = reshape(data.data(:,2),N,N,N);\n"
                                        "z = reshape(data.data(:,3),N,N,N);\n"
                                        "d = reshape(data.data(:,4),N,N,N);\n"
                                        "nx= reshape(data.data(:,5),N,N,N);\n"
                                        "ny= reshape(data.data(:,6),N,N,N);\n"
                                        "nz= reshape(data.data(:,7),N,N,N);\n"
                                        "mid = (N-1)/2;\n"
                                        "surf(x(:,:,mid),y(:,:,mid),d(:,:,mid),'FaceColor','none')\n"
                                        "hold on\n"
                                        "quiver(x(:,:,mid),y(:,:,mid),nx(:,:,mid),ny(:,:,mid))\n"
                                        "legend('distance','normal');\n"
                                        "if (min(min(d(:,:,mid)))<0 && max(max(d(:,:,mid)))>0)\n"
                                        "    contourf(x(:,:,mid),y(:,:,mid),d(:,:,mid),[0 0])\n"
                                        "    legend('distance','normal','distance>=0'); set(legend,'Location','best');\n"
                                        "end\n"
                                        "hold off\n"
                                        "xlabel('x'); ylabel('y'); title('cross-section at z=0');\n"
                                        "view(0,0); axis equal;");
     logger(INFO, "Run LevelSet.m to view level set");
 }

References getDistanceAndNormalLabCoordinates(), i, INFO, j, k, logger, n, WallFunction::normal(), radius_, s, and helpers::writeToFile().

◆ writeVTK()

void LevelSetWall::writeVTK ( VTKContainer & vtk ) const

overridevirtual

Adds the vtk wall representation to the VTK container

Reimplemented from BaseWall.

 {
     vtk.triangleStrips = vtkLabFrame_.triangleStrips;
     vtk.points = vtkLabFrame_.points;
     for (Vec3D& p : vtk.points)
     {
         getOrientation().rotate(p);
         p += getPosition();
     }
 }

References BaseInteractable::getOrientation(), BaseInteractable::getPosition(), p, VTKContainer::points, Quaternion::rotate(), VTKContainer::triangleStrips, and vtkLabFrame_.

Member Data Documentation

◆ levelSet_

double LevelSetWall::levelSet_[2 *N+1][2 *N+1][2 *N+1]

private

Referenced by getDistanceAndNormalLabCoordinates(), setShapeCube(), setShapeCylinder(), setShapeDiamond(), setShapeFourSided(), and setShapeSphere().

◆ N

const int LevelSetWall::N = 10

staticprivate

Referenced by createVTK(), createVTKCube(), getDistanceAndNormalLabCoordinates(), setShapeCube(), setShapeCylinder(), setShapeDiamond(), setShapeFourSided(), and setShapeSphere().

◆ radius_

double LevelSetWall::radius_ = 1

private

Referenced by createVTK(), createVTKCube(), createVTKDiamond(), createVTKSphere(), getDistanceAndNormalLabCoordinates(), setShapeCube(), setShapeCylinder(), setShapeDiamond(), setShapeFourSided(), setShapeSphere(), and writeToFile().

◆ vtkLabFrame_

VTKContainer LevelSetWall::vtkLabFrame_

private

Referenced by createVTK(), createVTKCube(), createVTKDiamond(), createVTKSphere(), and writeVTK().

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

Public Types

Public Member Functions

Private Member Functions

Private Attributes

Static Private Attributes

Additional Inherited Members

Detailed Description

Member Enumeration Documentation

◆ Shape

Constructor & Destructor Documentation

◆ LevelSetWall()

◆ ~LevelSetWall()

Member Function Documentation

◆ copy()

◆ createVTK()

◆ createVTKCube()

◆ createVTKDiamond()

◆ createVTKSphere()

◆ getDistanceAndNormal()

◆ getDistanceAndNormalLabCoordinates()

◆ getName()

◆ read()

◆ setShapeCube()

◆ setShapeCylinder()

◆ setShapeDiamond()

◆ setShapeFourSided()

◆ setShapeSphere()

◆ writeToFile()

◆ writeVTK()

Member Data Documentation

◆ levelSet_

◆ N

◆ radius_

◆ vtkLabFrame_