BaseInteractable.h

Go to the documentation of this file.

//Copyright (c) 2013-2020, The MercuryDPM Developers Team. All rights reserved.
//For the list of developers, see <http://www.MercuryDPM.org/Team>.
//
//Redistribution and use in source and binary forms, with or without
//modification, are permitted provided that the following conditions are met:
//  * Redistributions of source code must retain the above copyright
//    notice, this list of conditions and the following disclaimer.
//  * Redistributions in binary form must reproduce the above copyright
//    notice, this list of conditions and the following disclaimer in the
//    documentation and/or other materials provided with the distribution.
//  * Neither the name MercuryDPM nor the
//    names of its contributors may be used to endorse or promote products
//    derived from this software without specific prior written permission.
//
//THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND
//ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
//WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
//DISCLAIMED. IN NO EVENT SHALL THE MERCURYDPM DEVELOPERS TEAM BE LIABLE FOR ANY
//DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
//(INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
//LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND
//ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
//(INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
//SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.

#ifndef BASEINTERACTABLE_H
#define BASEINTERACTABLE_H

#include <functional>
#include <vector>
#include "BaseObject.h"
#include "Math/Vector.h"
#include "Math/Quaternion.h"

class BaseParticle;

class ParticleSpecies;

class BaseInteraction;

class InteractionHandler;

class BaseInteractable : public BaseObject
{
public:
    BaseInteractable();
    
    BaseInteractable(const BaseInteractable& p);
    
    ~BaseInteractable() override;
    
    void read(std::istream& is) override;
    
    void write(std::ostream& os) const override;
    
    unsigned int getIndSpecies() const
    { return indSpecies_; }
    
    virtual void setIndSpecies(unsigned int indSpecies)
    { indSpecies_ = indSpecies; }
    
    const ParticleSpecies* getSpecies() const
    {
        return species_;
    }
    
    void setSpecies(const ParticleSpecies* species);
    
    const Vec3D& getForce() const
    { return force_; }
    
    const Vec3D& getTorque() const
    { return torque_; }
    
    void setForce(const Vec3D& force)
    { force_ = force; }
    
    void setTorque(const Vec3D& torque)
    { torque_ = torque; }
    
    void addForce(const Vec3D& addForce);
    
    void addTorque(const Vec3D& addTorque);
    
    /*
     * This function sets the forces and torques acting on an interactable (particle/wall) to zero at teh beginning of computeOneTimeStep.
     *
     * For omp simulations, it further sets does the first step of summing up all forces (and torques) acting on an interactable.
     * It is step 1 of the add-by-reduction approach:
     *  1. Declare a vector forceOMP_, one element per thread, and set it to zero.
     *  2. When a new force needs to be added to the interactable, don't add it immediately.
     *     Add it to forceOMP_[threadNum] instead
     *  3. Add all forces in forceOMP_ together to force_.
     * This reduction approach allows step 2 to be run in parallel, while step 3 needs to be done in sequence.
     * It is efficient, since step 3 is much quicker than step 2.
     */
    void resetForceTorque(int numberOfOMPthreads);
    
    /*
     * This function does the last step of summing up all forces (and torques) acting on an interactable.
     * It is step 3 of the add-by-reduction approach:
     *  1. Declare a vector forceOMP_, one element per thread, and set it to zero.
     *  2. When a new force needs to be added to the interactable, don't add it immediately.
     *     Add it to forceOMP_[threadNum] instead
     *  3. Add all forces in forceOMP_ together to force_.
     * This reduction approach allows step 2 to be run in parallel, while step 3 needs to be done in sequence.
     * It is efficient, since step 3 is much quicker than step 2.
     */
    void sumForceTorqueOMP();
    
    const Vec3D& getPosition() const
    { return position_; }
    
    const Quaternion& getOrientation() const
    { return orientation_; }
    
    void setPosition(const Vec3D& position)
    { position_ = position; }
    
    void setOrientationViaNormal(Vec3D normal);
    
    void setOrientationViaEuler(Vec3D eulerAngle);
    
    void setOrientation(const Quaternion& orientation)
    { orientation_ = orientation; }
    
    virtual void move(const Vec3D& move);
    
    virtual void rotate(const Vec3D& angularVelocityDt);
    
    const std::vector<BaseInteraction*>& getInteractions() const
    { return interactions_; }
    
    void addInteraction(BaseInteraction* I);
    
    bool removeInteraction(BaseInteraction* I);
    
    void copyInteractionsForPeriodicParticles(const BaseInteractable& p);
    
    void setVelocity(const Vec3D& velocity);
    
    void setAngularVelocity(const Vec3D& angularVelocity);
    
    void addVelocity(const Vec3D& velocity)
    { velocity_ += velocity; }
    
    void addAngularVelocity(const Vec3D& angularVelocity);
    
    virtual const Vec3D& getVelocity() const;
    
    virtual const Vec3D& getAngularVelocity() const;
    
    void setPrescribedPosition(const std::function<Vec3D(double)>& prescribedPosition);
    
    void applyPrescribedPosition(double time);
    
    void setPrescribedVelocity(const std::function<Vec3D(double)>& prescribedVelocity);
    
    void applyPrescribedVelocity(double time);
    
    void setPrescribedOrientation(const std::function<Quaternion(double)>& prescribedOrientation);
    
    void applyPrescribedOrientation(double time);
    
    void setPrescribedAngularVelocity(const std::function<Vec3D(double)>& prescribedAngularVelocity);
    
    void applyPrescribedAngularVelocity(double time);
    
    virtual BaseInteraction* getInteractionWith(BaseParticle* P, unsigned timeStamp, InteractionHandler* interactionHandler)=0;
    
    virtual const Vec3D getVelocityAtContact(const Vec3D& contact) const;
    
    void integrateBeforeForceComputation(double time, double timeStep);
    
    void integrateAfterForceComputation(double time, double timeStep);
    
    virtual bool isFixed() const =0;

    virtual Mdouble getInvMass() const
    { return 0.0; }

    virtual Mdouble getCurvature(const Vec3D& labFixedCoordinates) const
    { return 0.0; }

private:
    std::function<Vec3D(double)> prescribedPosition_;
    
    std::function<Vec3D(double)> prescribedVelocity_;
    
    std::function<Quaternion(double)> prescribedOrientation_;
    
    std::function<Vec3D(double)> prescribedAngularVelocity_;
    
    Vec3D position_;
    
    Quaternion orientation_;
    
    Vec3D angularVelocity_;
    
    Vec3D force_;
    
    Vec3D torque_;
    
    /*
     * This variable is needed in omp simulations to sum up all forces (and torques) acting on an interactable.
     * It is used in an add-by-reduction approach:
     *  1. Declare a vector forceOMP_, one element per thread, and set it to zero.
     *  2. When a new force needs to be added to the interactable, don't add it immediately.
     *     Add it to forceOMP_[threadNum] instead
     *  3. Add all forces in forceOMP_ together to force_.
     * This reduction approach allows step 2 to be run in parallel, while step 3 needs to be done in sequence.
     * It is efficient, since step 3 is much quicker than step 2.
     */
    std::vector<Vec3D> forceOMP_;
    
    std::vector<Vec3D> torqueOMP_;

    const ParticleSpecies* species_;
    
    unsigned int indSpecies_;
    
    Vec3D velocity_;
    
    std::vector<BaseInteraction*> interactions_;
};

#endif