/home/dorival/mechsys/lib/sph/particle.h

Go to the documentation of this file.

/************************************************************************
 * MechSys - Open Library for Mechanical Systems                        *
 * Copyright (C) 2005 Dorival M. Pedroso, Raul Durand                   *
 * Copyright (C) 2009 Sergio Galindo                                    *
 *                                                                      *
 * This program is free software: you can redistribute it and/or modify *
 * it under the terms of the GNU General Public License as published by *
 * the Free Software Foundation, either version 3 of the License, or    *
 * any later version.                                                   *
 *                                                                      *
 * This program is distributed in the hope that it will be useful,      *
 * but WITHOUT ANY WARRANTY; without even the implied warranty of       *
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the         *
 * GNU General Public License for more details.                         *
 *                                                                      *
 * You should have received a copy of the GNU General Public License    *
 * along with this program. If not, see <http://www.gnu.org/licenses/>  *
 ************************************************************************/

#ifndef MECHSYS_SPH_PARTICLE_H
#define MECHSYS_SPH_PARTICLE_H

// Std lib
#include <iostream>

// MechSys
#include <mechsys/linalg/matvec.h>
#include <mechsys/dem/special_functions.h>

namespace SPH {

class Particle
{
public:

    // Constructor
    Particle(Vec3_t const & x0, Vec3_t const & v0, double density0, double h0, bool Fixed=false);


    // Data
    bool   IsFree;                  
    Vec3_t xo;                      
    Vec3_t x;                       
    Vec3_t xb;                      
    Vec3_t v;                       
    Vec3_t a;                       
    double Pressure;                
    double Density;                 
    double Densityb;                
    double Density0;                
    double dDensity;                
    double h;                       


    // Methods
    void Move (double dt);                                                  
    void StartAccel (Vec3_t acc = Vec3_t(0.0,0.0,0.0)) {a = acc;};          
    void Translate  (Vec3_t const & V) {x+=V; xb+=V;};                      
    void ResetDisplacements ();                                             
    double MaxDisplacement ();                                              

};

inline Particle::Particle(Vec3_t const & x0, Vec3_t const & v0, double density0, double h0,bool Fixed)
{
    x = x0;
    xb = x;
    v = v0;
    Density = density0;
    Densityb = Density;
    Density0 = Density;
    IsFree = !Fixed;
    h = h0;
}

inline void Particle::Move (double dt)
{
    if (IsFree)
    {
        // Evolve position and velocity
        Vec3_t xa;
        xa = 2*x - xb + a*dt*dt;
        v = 0.5*(xa - xb)/dt;
        xb = x;
        x = xa;

        // Evolve density
        double dens = Density;
        Density = Densityb + 2*dt*dDensity;
        //std::cout << Density << std::endl;
        Densityb = dens;
    }
}

inline void Particle::ResetDisplacements ()
{
    xo = x;
}

inline double Particle::MaxDisplacement ()
{
    return Distance(x,xo);
}

}; // namespace SPH

#endif // MECHSYS_SPH_PARTICLE_H