Sfx_StVenant1DEquations.h

Go to the documentation of this file.

 
#pragma once
 
// C++ includes
#include <iostream>
#include <cstdlib>
// numeric
#include <valarray>
// SfxBase19 includes
#include "include/Sfx_SCL.h"
#include "include/Sfx_UniversalConstants.h"
// App includes
//#include "../SfxBase/Sfx_StVenant1DTerms.h"
//#include "../include/Swe_SectionGeometryHandler.h"
#include "../Nov_VS2019/dbpp_SectionFlow.h"
 
namespace Sfx { 
    class StateVectorField; 
}
namespace dbpp { class Uh; }
 
namespace Sfx
{
    //  'Model of'
    //  Conservative form of the St-Venant Equation in 
    //  One-dimension. We have two equations for the conservative 
    //  variables (mass and momentum). 

    class StVenant1D : public Sfx::SCLEquation
    {
    public:
        using valar64 = std::valarray<float64>;
       public:  
             StVenant1D( eFluxType aFtype = eFluxType::incomplete)
             : m_FluxType{ aFtype },
                 m_useFriction{false}
             {}

       SCLEquation* Clone() override { return new StVenant1D(*this); }

             float64 flux( const Sfx::StateVector& aState,
                 const eFluxType aFluxT = eFluxType::incomplete) override final;
             std::valarray<float64> flux( const Sfx::StateVectorField& aState,
                 const eFluxType aFluxT = eFluxType::incomplete);
             std::valarray<float64> flux(const dbpp::Uh& aGblNval,
                 const eFluxType aFluxT = eFluxType::incomplete);
             float64 flux(const Sfx::StateVector& aState, const dbpp::SectionFlow& aSect);
             std::string name() const override final { return std::string{ "StVenant1DNoFriction" }; }
             bool hasSourceTerm() const override final  { return true; }
             //bool useCompleteFlux() const override final { return false; };  convective and pressure (f=f_c+f_p)

             uint32 getNumberOfEquations() const override final { return 2; }
             uint32 getDimension() const override final { return 1; }
             eFluxType getFluxType() const noexcept { return m_FluxType;}
             bool useFriction() const noexcept { return m_useFriction;}
 
//protected:
        // not accessible from client side
             // noncopyable
             //StVenant1D(const StVenant1D& aOther) = default;
             //StVenant1D& operator= (const StVenant1D& aOther) = default;
 
       private:
             eFluxType m_FluxType; 
             bool m_useFriction;    
#if 0   
             // St-Venant flux without the pressure term
             inline float64 EvaluationFlux_I_1D( const Sfx::StateVector & aState)
             {
                 // u^2*h (without pressure term)
                 return (aState.Q()*aState.Q())/aState.A();
             }
             //inline float64 EvaluationFlux_I_1D( const SectFlow& aSect)
             //{
                // Sfx::StateVector w_sectVec = aSect.getStateVector();
 
                // // u^2*h (without pressure term)
                // return (w_sectVec.Q()*w_sectVec.Q())/w_sectVec.A();
             //}
 
             // St-Venant flux with the pressure term (total flux)
             inline float64 EvaluationFlux_C_1D( const Sfx::StateVector & aState, 
                 float64 B = 1.)  // section width
             {
                 // Flux including: u^2*h + gh^2/2 (with pressure term)
                 return EvaluationFlux_I_1D(aState) + 
                     Sfx::cGravity<double>*HydrostaticPressure(aState.A());
             }
             //inline float64 EvaluationFlux_C_1D( const SectFlow& aSect)
             //{
                // // Flux including: u^2*h + gh^2/2 (with pressure term)
                // return EvaluationFlux_I_1D(aSect) + 
                //   Sfx::cGravity<double>*HydrostaticPressure(aSect);
             //}
 
             // Hydrostatic pressure (gh^2/2.), (B = section width)
             inline float64 HydrostaticPressure( const float64 aArea, float64 B = 1.)
             {
                 return (aArea*aArea)/2.*B;
             }
             //inline float64 HydrostaticPressure( const SectFlow& aSect)
             //{
                // Sfx::StateVector w_sectVec = aSect.getStateVector();
                // const float64 w_sectWidth = SectionGeometryHandler::computeWidth(aSect);
                // return (w_sectVec.A()*w_sectVec.A())/2.*w_sectWidth;
             //}
#endif //0
             valar64 EvaluationFlux_I_1D( const valar64& U1, const valar64& U2, float64 B = 1.)
             {
                 return (U2 * U2) / U1;
             }

             valar64 CalculTermePressionHydrostatique1D( const valar64& A, float64 B = 1.)
             {
                 return ((A * A) / (2. * B));
             }

             valar64 EvaluationFlux_C_1D( const valar64& U1, const valar64& U2, float64 B = 1.)
             {
                 return EvaluationFlux_I_1D(U1, U2) + Sfx::cGravity<float64>*CalculTermePressionHydrostatique1D(U1, B);
             }
    };   
}  // End of namespace