DAkOmegaSSTFIML.H

Go to the documentation of this file.

/*---------------------------------------------------------------------------*\
 
    DAFoam  : Discrete Adjoint with OpenFOAM
    Version : v3
 
    Description:
    Child class for the kOmegaSST model with a betaFieldInversion field 
    multiplying to the omega transport equation. This betaFieldInversion term
    can then be trained to improve the kOmegaSST model. 
 
    This file is modified from OpenFOAM's source code
    src/TurbulenceModels/turbulenceModels/RAS/kOmega/kOmega.H
 
    OpenFOAM: The Open Source CFD Toolbox
 
    Copyright (C): 2011-2016 OpenFOAM Foundation
 
    OpenFOAM License:
 
        OpenFOAM 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
        (at your option) any later version.
    
        OpenFOAM 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 OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
 
\*---------------------------------------------------------------------------*/
 
#ifndef DAkOmegaSSTFIML_H
#define DAkOmegaSSTFIML_H
 
#include "DATurbulenceModel.H"
#include "addToRunTimeSelectionTable.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
namespace Foam
{
 
/*---------------------------------------------------------------------------*\
                       Class DAkOmegaSSTFIML Declaration
\*---------------------------------------------------------------------------*/
 
class DAkOmegaSSTFIML
    : public DATurbulenceModel
{
 
protected:
 
    dimensionedScalar alphaK1_;
    dimensionedScalar alphaK2_;
 
    dimensionedScalar alphaOmega1_;
    dimensionedScalar alphaOmega2_;
 
    dimensionedScalar gamma1_;
    dimensionedScalar gamma2_;
 
    dimensionedScalar beta1_;
    dimensionedScalar beta2_;
 
    dimensionedScalar betaStar_;
 
    dimensionedScalar a1_;
    dimensionedScalar b1_;
    dimensionedScalar c1_;
 
    Switch F3_;
 
    scalar* inputs_ = new scalar[mesh_.nCells() * 9];
    scalar* outputs_ = new scalar[mesh_.nCells()];
 
#ifdef CODI_AD_FORWARD
    double* inputsDouble_ = new double[mesh_.nCells() * 9];
    double* outputsDouble_ = new double[mesh_.nCells()];
#endif
 
    void calcBetaField();
 
 
    tmp<volScalarField> F1(const volScalarField& CDkOmega) const;
    tmp<volScalarField> F2() const;
    tmp<volScalarField> F3() const;
    tmp<volScalarField> F23() const;
 
    tmp<volScalarField> blend(
        const volScalarField& F1,
        const dimensionedScalar& psi1,
        const dimensionedScalar& psi2) const
    {
        return F1 * (psi1 - psi2) + psi2;
    }
 
    tmp<volScalarField::Internal> blend(
        const volScalarField::Internal& F1,
        const dimensionedScalar& psi1,
        const dimensionedScalar& psi2) const
    {
        return F1 * (psi1 - psi2) + psi2;
    }
 
    tmp<volScalarField> alphaK(const volScalarField& F1) const
    {
        return blend(F1, alphaK1_, alphaK2_);
    }
 
    tmp<volScalarField> alphaOmega(const volScalarField& F1) const
    {
        return blend(F1, alphaOmega1_, alphaOmega2_);
    }
 
    tmp<volScalarField::Internal> beta(
        const volScalarField::Internal& F1) const
    {
        return blend(F1, beta1_, beta2_);
    }
 
    tmp<volScalarField::Internal> gamma(
        const volScalarField::Internal& F1) const
    {
        return blend(F1, gamma1_, gamma2_);
    }
 
    //- Return the effective diffusivity for k
    tmp<volScalarField> DkEff(const volScalarField& F1) const
    {
        return tmp<volScalarField>(
            new volScalarField("DkEff", alphaK(F1) * nut_ + this->nu()));
    }
 
    //- Return the effective diffusivity for omega
    tmp<volScalarField> DomegaEff(const volScalarField& F1) const
    {
        return tmp<volScalarField>(
            new volScalarField(
                "DomegaEff",
                alphaOmega(F1) * nut_ + this->nu()));
    }
 
    //- Return k production rate
    tmp<volScalarField::Internal> Pk(
        const volScalarField::Internal& G) const;
 
    //- Return epsilon/k which for standard RAS is betaStar*omega
    tmp<volScalarField::Internal> epsilonByk(
        const volScalarField& F1,
        const volTensorField& gradU) const;
 
    //- Return G/nu
    tmp<volScalarField::Internal> GbyNu(
        const volScalarField::Internal& GbyNu0,
        const volScalarField::Internal& F2,
        const volScalarField::Internal& S2) const;
 
    tmp<fvScalarMatrix> kSource() const;
 
    tmp<fvScalarMatrix> omegaSource() const;
 
    tmp<fvScalarMatrix> Qsas(
        const volScalarField::Internal& S2,
        const volScalarField::Internal& gamma,
        const volScalarField::Internal& beta) const;
 
 
    volScalarField& omega_;
    volScalarField omegaRes_;
    volScalarField omegaResRef_;
    volScalarField omegaResPartDeriv_;
    volScalarField omegaRef_;
    volScalarField& k_;
    volScalarField kRes_;
    volScalarField kResRef_;
    volScalarField kResPartDeriv_;
    volScalarField kRef_;
 
    // Field inversion and machine learning fields
    volScalarField betaFieldInversion_;
    volScalarField betaFieldInversionML_;
    volScalarField QCriterion_;
    volScalarField& p_;
    volScalarField pGradAlongStream_;
    volScalarField turbulenceIntensity_;
    IOdictionary transportProperties_;
    volScalarField ReT_;
    volScalarField convectionTKE_;
    volScalarField tauRatio_;
    volScalarField pressureStress_;
    volScalarField curvature_;
    volScalarField UGradMisalignment_;
 
    //word eqnToModify_;
 
    //scalar betaFieldInversionWeight_;
 
    const volScalarField& y_;
 
    scalarList omegaNearWall_;
 
    label solveTurbState_ = 0;
 
    label printInterval_;
 
public:
    TypeName("kOmegaSSTFIML");
    // Constructors
 
    //- Construct from components
    DAkOmegaSSTFIML(
        const word modelType,
        const fvMesh& mesh,
        const DAOption& daOption);
 
    //- Destructor
    virtual ~DAkOmegaSSTFIML()
    {
        delete[] inputs_;
        delete[] outputs_;
#ifdef CODI_AD_FORWARD
        delete[] inputsDouble_;
        delete[] outputsDouble_;
#endif
    }
 
    // Member functions
 
    virtual void correctModelStates(wordList& modelStates) const;
 
    virtual void correctNut();
 
    virtual void correctBoundaryConditions();
 
    virtual void updateIntermediateVariables();
 
    virtual void correctStateResidualModelCon(List<List<word>>& stateCon) const;
 
    virtual void addModelResidualCon(HashTable<List<List<word>>>& allCon) const;
 
    virtual void calcResiduals(const dictionary& options);
 
    virtual void correct();
 
    void saveOmegaNearWall();
 
    void setOmegaNearWall();
 
    void correctOmegaBoundaryConditions();
 
#ifdef CODI_AD_REVERSE
 
    static void betaCompute(
        const double* x,
        size_t n,
        double* y,
        size_t m,
        codi::ExternalFunctionUserData* d)
    {
        DAUtility::pyCalcBetaInterface(x, n, y, m, DAUtility::pyCalcBeta);
    }
 
    static void betaJacVecProd(
        const double* x,
        double* x_b,
        size_t n,
        const double* y,
        const double* y_b,
        size_t m,
        codi::ExternalFunctionUserData* d)
    {
        DAUtility::pyCalcBetaJacVecProdInterface(x, x_b, n, y, y_b, m, DAUtility::pyCalcBetaJacVecProd);
    }
#endif
};
 
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
} // End namespace Foam
 
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
#endif
 
// ************************************************************************* //