dbpp_HLLFlux1D.h

Go to the documentation of this file.

#pragma once
 
// C++ include
#include <tuple>
// SfxBase19 include
#include "include/Sfx_CellFaceVariables.h"
 
namespace dbpp 
{
    class HLLFlux1D 
    {
    public:
        enum eFluxComp
        {
            FL1 = 0, 
            FR1 = 1, 
            FL2 = 2, 
            FR2 = 3, 
        };

        using ShockSpeed = std::pair<float64/*SL*/, float64/*SR*/>;
        using PhysicalFlux = std::tuple<float64/*FL1*/, float64/*FR1*/, float64/*FL2*/, float64/*FR2*/>;
        void setShockSpeed( const ShockSpeed& aShockSpeed) { m_shocSpeed = aShockSpeed; }
        void setphysicalFlux( const PhysicalFlux& aPhyF) { m_physFlux = aPhyF; }
        virtual std::pair<float64/**/, float64>
            calculFF( const Sfx::cellFaceVariables& aCellFaceVar) const
        {
            // compute the HLL parameters
//          HLLparameters w_prm;
//          w_prm.computeHLLParametres(aUL, aUR, w_sectionFlow->B());
 
            //  // call flux function instead of calling the function below
            //  // this is the physical flux at left side
            //  const double w_physFL2 = w_stVen1D->flux( Sfx::StateVector(UL1, UL2)); //FL2
            //  const double w_physFR2 = w_stVen1D->flux( Sfx::StateVector(UR1, UR2)); //FR2
            //}
            // St-Venant flux term (conservative form) and not consider hydrostatic pressure
//          const double FL2 = HydroTerms::EvaluationF2_I_1D(UL1, UL2, w_secti->B());
//          const double FR2 = HydroTerms::EvaluationF2_I_1D(UR1, UR2, B);
 
            // NOTE: 
            //   We solving the Riemann problem, isn't?
            //   Might as well call the Riemann solver
 
            double FF1{};
            double FF2{};
 
            //
            //  HLL Algorithm   
            //
 
            //  Calcul de Fi+1/2 (we are at cell face)
            if (m_shocSpeed.first > 0.) //SL
            {
                FF1 = std::get<FL1>(m_physFlux);  //FL1;
                FF2 = std::get<FL2>(m_physFlux);  //FL2;
            }
            else if (m_shocSpeed.second < 0.) //SR
            {
                FF1 = std::get<FR1>(m_physFlux);  //FR1;
                FF2 = std::get<FR2>(m_physFlux);  //FR2;
            }
            else
            {
                // HLL flux 
    /*          FF1[i] = (SR*FL1 - SL*FR1 + SL*SR*(UR1 - UL1)) / (SR - SL);
                FF2[i] = (SR*FL2 - SL*FR2 + SL*SR*(UR2 - UL2)) / (SR - SL);*/
                const auto [SL, SR] = m_shocSpeed;
                const auto [FL1, FR1, FL2, FR2] = m_physFlux;
 
                const auto& w_leftStateVar = aCellFaceVar.UL();
                const auto& w_rightStateVar = aCellFaceVar.UR();
 
                const auto UR1 = w_rightStateVar.A();
                const auto UL1 = w_leftStateVar.A();
 
                const auto UR2 = w_rightStateVar.Q();
                const auto UL2 = w_leftStateVar.Q();
 
                FF1 = (SR * FL1 - SL * FR1 + SL * SR * (UR1 - UL1)) / (SR - SL);
                FF2 = (SR * FL2 - SL * FR2 + SL * SR * (UR2 - UL2)) / (SR - SL);
            }
 
            // debugging purpose
            return std::pair{ FF1, FF2 }; // default for debugging purpose 
        }
    private:
        ShockSpeed m_shocSpeed;  
        PhysicalFlux m_physFlux; 
    };
}// End of namespace