CFD Online Discussion Forums

CFD Online Discussion Forums (https://www.cfd-online.com/Forums/)

- OpenFOAM Programming & Development (https://www.cfd-online.com/Forums/openfoam-programming-development/)

- - Hybrid Scheme (https://www.cfd-online.com/Forums/openfoam-programming-development/129787-hybrid-scheme.html)

Dear all

I was trying to implement the hybrid scheme but I have faced some problems.
I will appreciate any suggestions,

for the implementation I have added three more parameter to the limiter function to be able to define the Peclect number. then I pick the Quick scheme as the template for new scheme.

now the code is reading as :

Code:

    scalar limiter

    (

        const scalar cdWeight,

        const scalar faceFlux,

        const typename LimiterFunc::phiType& phiP,

        const typename LimiterFunc::phiType& phiN,

        const typename LimiterFunc::gradPhiType& gradcP,

        const typename LimiterFunc::gradPhiType& gradcN,

        const vector& d,

        const vector& Sf,

        const scalar& magSf ,

       const scalar Diffu

    ) const

    {



    scalar Pe = (faceFlux/magSf)*(d & Sf)/(Diffu*magSf);

    scalar QLimiter;

    QLimiter = pos(1-0.5*mag(Pe));



        // Limit the limiter between upwind and downwind

        return  QLimiter;

    }

the code compiles well. but then when I want to run a case it gives me floating point error. if I put the Qlimiter 0 or 1 it runs without problem. but as I use a function of Pe it wont run.

best
Mahdi

Dear Mahdi,
sorry if stating the obvious: prevent magSf and or Diffu from being zero or add SMALL to the divider

Hope anyone else has a better solution
Best regards
Wouter

Dear Wouter

your suggestion will make it work. But then I'm just worried what will happen to boundaries ( possible place that Diffu and Magsf could be zero?!)

best
Mahdi

now it runs

Im reading the Diffu in limitedscheme class in this way

Quote:

const volScalarField& Diff =
phi.db().objectRegistry::lookupObject<volScalarFie ld>(Diffname_);
tmp<surfaceScalarField> tDiffu=fvc::interpolate(Diff);

const surfaceScalarField& Diffu =tDiffu();

best
mahdi