cleaning after debug

classes/sph0/louis/src_cpp/FixedParticle.cpp

@@ -32,27 +32,15 @@ FixedParticle::update_vars()
         this->rho[2] = this->rho[0] + drho_dt * dt / 2.0; // prepare second step
         this->speed[1] = this->speed[0];
         this->coord[1] = this->coord[0];
-        double p = this->model.eqState->pressure(this->rho[1]);
-        double c = this->model.eqState->speed_of_sound(this->rho[1]);
-        this->p[1] = this->compute_pressure(this->rho[1]);
-        this->c[1] = this->compute_speedofsound(this->rho[1]); 
-        if(p != this->p[1] || c != this->c[1])
-        {
-            std::cout << "Pressure or speed of sound changed" << std::endl;
-            std::cout << "p: " << p << " c: " << c << std::endl;
-            std::cout << "p1: " << this->p[1] << " c1: " << this->c[1] << std::endl;
-        }
-
-
+        this->p[1] = this->model.eqState->pressure(this->rho[1]);
+        this->c[1] = this->model.eqState->speed_of_sound(this->rho[1]);
     }
     else // 2nd RK step
     {
         this->rho[2] = this->rho[2] + drho_dt * dt / 2.0;
         this->speed[2] = this->speed[1];
         this->coord[2] = this->coord[1];
-        // this->p[2] = this->model.eqState->pressure(this->rho[2]);
-        // this->c[2] = this->model.eqState->speed_of_sound(this->rho[2]);
-        this->p[2] = this->compute_pressure(this->rho[2]);
-        this->c[2] = this->compute_speedofsound(this->rho[2]);
+        this->p[2] = this->model.eqState->pressure(this->rho[2]);
+        this->c[2] = this->model.eqState->speed_of_sound(this->rho[2]);
     }
 }

classes/sph0/louis/src_cpp/MobileParticle.cpp

@@ -72,20 +72,16 @@ MobileParticle::update_vars()
         this->speed[2] = this->speed[0] + du_dt * dt / 2.0;
         this->coord[1] = this->coord[0] + this->speed[0] * dt;
         this->coord[2] = this->coord[0] + this->speed[0] * dt / 2.0;
-        this->p[1] = this->compute_pressure(this->rho[1]);
-        this->c[1] = this->compute_speedofsound(this->rho[1]);        
-        // this->p[1] = this->model.eqState->pressure(this->rho[1]);
-        // this->c[1] = this->model.eqState->speed_of_sound(this->rho[1]);
+        this->p[1] = this->model.eqState->pressure(this->rho[1]);
+        this->c[1] = this->model.eqState->speed_of_sound(this->rho[1]);
     }
     else // 2nd RK step
     {
         this->rho[2] = this->rho[2] + drho_dt * dt / 2.0;
         this->speed[2] = this->speed[2] + du_dt * dt / 2.0;
         this->coord[2] = this->coord[2] + this->speed[1] * dt / 2.0;
-        this->p[2] = this->compute_pressure(this->rho[2]);
-        this->c[2] = this->compute_speedofsound(this->rho[2]);
-        // this->p[2] = this->model.eqState->pressure(this->rho[2]);
-        // this->c[2] = this->model.eqState->speed_of_sound(this->rho[2]);
+        this->p[2] = this->model.eqState->pressure(this->rho[2]);
+        this->c[2] = this->model.eqState->speed_of_sound(this->rho[2]);
     }
 }
 

classes/sph0/louis/src_cpp/Model.cpp

@@ -200,12 +200,6 @@ Model::load_parameters(std::string const &param_path)
     file >> this->saveInt;
     file.close();
 
-    this->c_0 = c0;
-    this->rho_0 = rho0;
-    this->eqnState = eqnState;
-    this->state_gamma = gamma;
-    this->molMass = molMass;
-
     // create kernel
     switch (kernelKind)
     {

classes/sph0/louis/src_cpp/Model.h

@@ -28,23 +28,13 @@ class Model
     double alpha;         ///< weighting factor in the artificial viscosity formulation
     double beta;          ///< weighting factor in the artificial viscosity formulation
 
-    int eqnState;         ///< equation of state
-                          ///<   1 = ideal gas law
-                          ///<   2 = quasi-incompressible fluid
-    int state_gamma;      ///< power in eqn State 2.
-                          ///< often taken around 7
-    double molMass;       ///< Molar mass of the fluid for the prefect gas law
-
     int kernelCorrection; ///< correction of the kernel
                           ///<   0 = no correction
                           ///<   1 = correction enabled
     double maxTime;       ///< simulation time in seconds
     double saveInt;       ///< saving interval
     double h_0;           ///< initial smoothing length
 
-    double rho_0;         ///< density of the fluid at free surface
-    double c_0;           ///< speed of sound in normal conditions
-
     double timeStep;      ///< timestep (not constant)
     double currentTime;   ///< current time
     int RKstep;           ///< used to know in which RK iteration we are      [RB] (1 or 2)

classes/sph0/louis/src_cpp/Particle.cpp

@@ -21,70 +21,12 @@ Particle::load(std::ifstream &ufile, double h_0)
     this->speed[0] << u_x, u_y, u_z;
     this->rho[0] = rho;
     this->m = m;
-    this->h = h_0;
-    this->p[0] = this->compute_pressure(rho);
-    this->c[0] = this->compute_speedofsound(rho);   
-    // this->p[0] = this->model.eqState->pressure(rho);
-    // this->c[0] = this->model.eqState->speed_of_sound(rho);
+    this->h = h_0;  
+    this->p[0] = this->model.eqState->pressure(rho);
+    this->c[0] = this->model.eqState->speed_of_sound(rho);
     this->max_mu_ab = 0.0;
 }
 
-// calculates the pressure according to the equation of state chosen.
-// @param rho  : actual density
-
-double
-Particle::compute_pressure(double rho) const
-{
-    double pressure;
-    const double idealGasCst = 8.3144621;
-
-    switch (this->model.eqnState)
-    {
-    case LAW_IDEAL_GAS:
-    {
-        pressure = (rho / this->model.rho_0 - 1.0) *
-                   idealGasCst * 293.15 / this->model.molMass; // eq (3.24)
-        break;
-    }
-    case LAW_QINC_FLUID:
-    {
-        double B = this->model.c_0 * this->model.c_0 *
-                   this->model.rho_0 / this->model.state_gamma; // eq (3.27)
-        pressure = B * (pow(rho / this->model.rho_0, this->model.state_gamma) - 1.0);
-        break;
-    }
-    default:
-        throw std::runtime_error("Bad value for equation of state (1,2)");
-    }
-    return pressure;
-}
-
-/// Calculates the celerity according to the equation of state chosen.
-/// The equation used is @f[c = \sqrt{\frac{dp}{d\rho}}@f]
-/// @param rho  : actual density
-
-double
-Particle::compute_speedofsound(double rho) const
-{
-    double celerity;
-
-    switch (this->model.eqnState)
-    {
-    case LAW_IDEAL_GAS:
-        // 1 = considering the ideal gas law at 20 degrees C
-        celerity = this->model.c_0; // eq (3.36)
-        break;
-    case LAW_QINC_FLUID:
-        // 2 = considering a quasi-incompressible fluid
-        celerity = this->model.c_0 *
-                   pow(rho / this->model.rho_0, (this->model.state_gamma - 1) / 2); // eq (3.37)
-        break;
-    default:
-        throw std::runtime_error("Bad value for equation of state (1,2)");
-    }
-    return celerity;
-}
-
 /// Saves the state of a particle onto disk
 
 void

classes/sph0/louis/src_cpp/Particle.h

@@ -47,8 +47,6 @@ class Particle
     virtual void update_vars() = 0;
 
 protected:
-    double compute_pressure(double rho) const;
-    double compute_speedofsound(double rho) const;
     void gradW();
     void kernel_corr();
     void getNeighbours();

classes/sph0/louis/src_cpp/QtVTKHook.cpp

@@ -19,9 +19,7 @@ DisplayWindow::DisplayWindow(QWidget *parent) : QWidget(parent)
     resize(800, 600);
 
     setupGUI();
-
     addCube();
-
 }
 
 DisplayWindow::~DisplayWindow()

classes/sph0/louis/src_cpp/Sorter.cpp

@@ -95,12 +95,6 @@ Sorter::compute_hmax()
     for(auto &p : this->model.particles)
         if (p->h > hmax)
             hmax = p->h;
-
-    // Increase of h_max in order to have a security if h changes.
-    // This is done according to the equation of state used.
-    if (this->model.eqnState == LAW_IDEAL_GAS)
-        hmax = 1.1 * hmax;
-    else
-        hmax = 1.02 * hmax;
+    hmax *= this->model.eqState->h_factor();
     return hmax;
 }

classes/sph0/louis/src_cpp/Sorter.h

@@ -10,7 +10,6 @@
 class Sorter
 {
     Model &model;
-   //double h_max;    ///< maximum smoothing length
 
 public:
     double dx; ///< length of a side of a cube