Added more inputs for the regModel and variance obj (#565)

* Moved FIData util to preProcessing. * Added more reg inputs.
mdolab · Jan 18, 2024 · 458148b · 458148b
1 parent bcfadd7
commit 458148b
Show file tree

Hide file tree

Showing 29 changed files with 1,081 additions and 230 deletions.
diff --git a/Allclean b/Allclean
@@ -9,4 +9,5 @@ cd src/adjoint && ./Allclean && cd -
 cd src/pyDASolvers && ./Allclean && cd -
 cd src/utilities/coloring && ./Allclean && cd -
 cd src/utilities/postProcessing && ./Allclean && cd -
+cd src/utilities/preProcessing && ./Allclean && cd -
 cd tests && ./Allclean && cd -
diff --git a/Allmake b/Allmake
@@ -18,4 +18,5 @@ fi
 cd src/adjoint && ./Allmake $argm && cd -
 cd src/pyDASolvers && ./Allmake $argm && cd -
 cd src/utilities/coloring && ./Allmake $argm && cd -
-cd src/utilities/postProcessing && ./Allmake $argm && cd -
+cd src/utilities/postProcessing && ./Allmake $argm && cd -
+cd src/utilities/preProcessing && ./Allmake $argm && cd -
diff --git a/src/adjoint/DAModel/DAModel.C b/src/adjoint/DAModel/DAModel.C
@@ -314,6 +314,28 @@ void DAModel::getTurbProdOverDestruct(scalarList& PoD) const
 
 }
 
+void DAModel::getTurbConvOverProd(scalarList& CoP) const
+{
+    /*
+    Description: 
+        return the value of the convective over production term from the turbulence model
+    */
+
+#ifndef SolidDASolver
+    if (hasTurbulenceModel_)
+    {
+        DATurbulenceModel& daTurb = const_cast<DATurbulenceModel&>(
+            mesh_.thisDb().lookupObject<DATurbulenceModel>("DATurbulenceModel"));
+        daTurb.getTurbConvOverProd(CoP);
+    }
+
+    if (hasRadiationModel_)
+    {
+    }
+#endif
+
+}
+
 #ifdef CompressibleFlow
 const fluidThermo& DAModel::getThermo() const
 {

diff --git a/src/adjoint/DAModel/DAModel.H b/src/adjoint/DAModel/DAModel.H
@@ -111,6 +111,9 @@ public:
     /// return the value of the destruction term from the turbulence model
     void getTurbProdOverDestruct(scalarList& PoD) const;
 
+    /// return the value of the convective over production term from the turbulence model
+    void getTurbConvOverProd(scalarList& CoP) const;
+
 #ifndef SolidDASolver
     /// get a reference to DATurbulenceModel
     const DATurbulenceModel& getDATurbulenceModel() const;

diff --git a/src/adjoint/DAModel/DATurbulenceModel/DASpalartAllmaras.C b/src/adjoint/DAModel/DATurbulenceModel/DASpalartAllmaras.C
@@ -539,6 +539,27 @@ void DASpalartAllmaras::getTurbProdOverDestruct(scalarList& PoD) const
         PoD[cellI] = P[cellI] / (D[cellI] + 1e-16);
     }
 }
+
+void DASpalartAllmaras::getTurbConvOverProd(scalarList& CoP) const
+{
+    /*
+    Description:
+        Return the value of the convective over production term from the turbulence model 
+    */
+
+    const volScalarField chi(this->chi());
+    const volScalarField fv1(this->fv1(chi));
+
+    const volScalarField Stilda(this->Stilda(chi, fv1));
+
+    volScalarField P = Cb1_ * phase_ * rho_ * Stilda * nuTilda_;
+    volScalarField C = fvc::div(phaseRhoPhi_, nuTilda_);
+
+    forAll(P, cellI)
+    {
+        CoP[cellI] = C[cellI] / (P[cellI] + 1e-16);
+    }
+}
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 } // End namespace Foam

diff --git a/src/adjoint/DAModel/DATurbulenceModel/DASpalartAllmaras.H b/src/adjoint/DAModel/DATurbulenceModel/DASpalartAllmaras.H
@@ -149,6 +149,9 @@ public:
 
     /// return the value of the destruction term from the turbulence model
     virtual void getTurbProdOverDestruct(scalarList& PoD) const;
+
+    /// return the value of the convective over production term from the turbulence model
+    virtual void getTurbConvOverProd(scalarList& CoP) const;
 };
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

diff --git a/src/adjoint/DAModel/DATurbulenceModel/DASpalartAllmarasFv3.C b/src/adjoint/DAModel/DATurbulenceModel/DASpalartAllmarasFv3.C
@@ -767,6 +767,29 @@ void DASpalartAllmarasFv3::getTurbProdOverDestruct(scalarList& PoD) const
         PoD[cellI] = P[cellI] / (D[cellI] + 1e-16);
     }
 }
+
+void DASpalartAllmarasFv3::getTurbConvOverProd(scalarList& CoP) const
+{
+    /*
+    Description:
+        Return the value of the convective over production term from the turbulence model 
+    */
+
+    const volScalarField chi(this->chi());
+    const volScalarField fv1(this->fv1(chi));
+
+    const volScalarField Stilda(
+        this->fv3(chi, fv1) * ::sqrt(2.0) * mag(skew(fvc::grad(U_)))
+        + this->fv2(chi, fv1) * nuTilda_ / sqr(kappa_ * y_));
+
+    volScalarField P = Cb1_ * phase_ * rho_ * Stilda * nuTilda_;
+    volScalarField C = fvc::div(phaseRhoPhi_, nuTilda_);
+
+    forAll(P, cellI)
+    {
+        CoP[cellI] = C[cellI] / (P[cellI] + 1e-16);
+    }
+}
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 } // End namespace Foam

diff --git a/src/adjoint/DAModel/DATurbulenceModel/DASpalartAllmarasFv3.H b/src/adjoint/DAModel/DATurbulenceModel/DASpalartAllmarasFv3.H
@@ -171,6 +171,9 @@ public:
 
     /// return the value of the destruction term from the turbulence model
     virtual void getTurbProdOverDestruct(scalarList& PoD) const;
+
+    /// return the value of the convective over production term from the turbulence model
+    virtual void getTurbConvOverProd(scalarList& CoP) const;
 };
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

diff --git a/src/adjoint/DAModel/DATurbulenceModel/DATurbulenceModel.C b/src/adjoint/DAModel/DATurbulenceModel/DATurbulenceModel.C
@@ -504,6 +504,18 @@ void DATurbulenceModel::getTurbProdOverDestruct(scalarList& PoD) const
         << abort(FatalError);
 }
 
+void DATurbulenceModel::getTurbConvOverProd(scalarList& CoP) const
+{
+    /*
+    Description:
+        Return the value of the convective over production term from the turbulence model 
+    */
+
+    FatalErrorIn("DATurbulenceModel::getTurbConvOverProd")
+        << "Child class not implemented!"
+        << abort(FatalError);
+}
+
 void DATurbulenceModel::invTranProdNuTildaEqn(
     const volScalarField& mySource,
     volScalarField& pseudoNuTilda)

diff --git a/src/adjoint/DAModel/DATurbulenceModel/DATurbulenceModel.H b/src/adjoint/DAModel/DATurbulenceModel/DATurbulenceModel.H
@@ -214,9 +214,12 @@ public:
     /// return the value of the production term from the turbulence model
     virtual void getTurbProdTerm(scalarList& prodTerm) const;
 
-    /// return the value of the destruction term from the turbulence model
+    /// return the ratio of the production over destruction term from the turbulence model
     virtual void getTurbProdOverDestruct(scalarList& PoD) const;
 
+    /// return the value of the convective over production term from the turbulence model
+    virtual void getTurbConvOverProd(scalarList& CoP) const;
+
     /// dev terms
     tmp<volSymmTensorField> devRhoReff() const;
 

diff --git a/src/adjoint/DAModel/DATurbulenceModel/DAkEpsilon.C b/src/adjoint/DAModel/DATurbulenceModel/DAkEpsilon.C
@@ -704,6 +704,47 @@ void DAkEpsilon::getFvMatrixFields(
         lower = kEqn.lower();
     }
 }
+
+void DAkEpsilon::getTurbProdOverDestruct(scalarList& PoD) const
+{
+    /*
+    Description:
+        Return the value of the production over destruction term from the turbulence model 
+    */
+    tmp<volTensorField> tgradU = fvc::grad(U_);
+    volScalarField::Internal G(
+        "kEpsilon:G",
+        nut_.v() * (dev(twoSymm(tgradU().v())) && tgradU().v()));
+
+    volScalarField::Internal P = C1_ * phase_() * rho_() * G * epsilon_() / k_();
+    volScalarField::Internal D = C2_ * phase_() * rho_() * sqr(epsilon_()) / k_();
+
+    forAll(P, cellI)
+    {
+        PoD[cellI] = P[cellI] / (D[cellI] + 1e-16);
+    }
+}
+
+void DAkEpsilon::getTurbConvOverProd(scalarList& CoP) const
+{
+    /*
+    Description:
+        Return the value of the convective over production term from the turbulence model 
+    */
+
+    tmp<volTensorField> tgradU = fvc::grad(U_);
+    volScalarField::Internal G(
+        "kEpsilon:G",
+        nut_.v() * (dev(twoSymm(tgradU().v())) && tgradU().v()));
+
+    volScalarField::Internal P = C1_ * phase_() * rho_() * G * epsilon_() / k_();
+    volScalarField C = fvc::div(phaseRhoPhi_, epsilon_);
+
+    forAll(P, cellI)
+    {
+        CoP[cellI] = C[cellI] / (P[cellI] + 1e-16);
+    }
+}
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 } // End namespace Foam

diff --git a/src/adjoint/DAModel/DATurbulenceModel/DAkEpsilon.H b/src/adjoint/DAModel/DATurbulenceModel/DAkEpsilon.H
@@ -157,6 +157,12 @@ public:
         scalarField& diag,
         scalarField& upper,
         scalarField& lower);
+
+    /// return the value of the destruction term from the turbulence model
+    virtual void getTurbProdOverDestruct(scalarList& PoD) const;
+
+    /// return the value of the convective over production term from the turbulence model
+    virtual void getTurbConvOverProd(scalarList& CoP) const;
 };
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

diff --git a/src/adjoint/DAModel/DATurbulenceModel/DAkOmega.C b/src/adjoint/DAModel/DATurbulenceModel/DAkOmega.C
@@ -705,6 +705,43 @@ void DAkOmega::getFvMatrixFields(
         lower = kEqn.lower();
     }
 }
+
+void DAkOmega::getTurbProdOverDestruct(scalarList& PoD) const
+{
+    /*
+    Description:
+        Return the value of the production over destruction term from the turbulence model 
+    */
+    tmp<volTensorField> tgradU = fvc::grad(U_);
+    volScalarField G("kOmega:G", nut_ * (tgradU() && dev(twoSymm(tgradU()))));
+
+    volScalarField P = gamma_ * phase_ * rho_ * G * omega_ / k_;
+    volScalarField D = beta_ * phase_ * rho_ * sqr(omega_);
+
+    forAll(P, cellI)
+    {
+        PoD[cellI] = P[cellI] / (D[cellI] + 1e-16);
+    }
+}
+
+void DAkOmega::getTurbConvOverProd(scalarList& CoP) const
+{
+    /*
+    Description:
+        Return the value of the convective over production term from the turbulence model 
+    */
+
+    tmp<volTensorField> tgradU = fvc::grad(U_);
+    volScalarField G("kOmega:G", nut_ * (tgradU() && dev(twoSymm(tgradU()))));
+
+    volScalarField P = gamma_ * phase_ * rho_ * G * omega_ / k_;
+    volScalarField C = fvc::div(phaseRhoPhi_, omega_);
+
+    forAll(P, cellI)
+    {
+        CoP[cellI] = C[cellI] / (P[cellI] + 1e-16);
+    }
+}
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 } // End namespace Foam

diff --git a/src/adjoint/DAModel/DATurbulenceModel/DAkOmega.H b/src/adjoint/DAModel/DATurbulenceModel/DAkOmega.H
@@ -164,6 +164,12 @@ public:
         scalarField& diag,
         scalarField& upper,
         scalarField& lower);
+
+    /// return the value of the destruction term from the turbulence model
+    virtual void getTurbProdOverDestruct(scalarList& PoD) const;
+
+    /// return the value of the convective over production term from the turbulence model
+    virtual void getTurbConvOverProd(scalarList& CoP) const;
 };
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

diff --git a/src/adjoint/DAModel/DATurbulenceModel/DAkOmegaSST.C b/src/adjoint/DAModel/DATurbulenceModel/DAkOmegaSST.C
@@ -891,6 +891,65 @@ void DAkOmegaSST::getFvMatrixFields(
         lower = kEqn.lower();
     }
 }
+
+void DAkOmegaSST::getTurbProdOverDestruct(scalarList& PoD) const
+{
+    /*
+    Description:
+        Return the value of the production over destruction term from the turbulence model 
+    */
+    tmp<volTensorField> tgradU = fvc::grad(U_);
+    volScalarField S2(2 * magSqr(symm(tgradU())));
+    volScalarField::Internal GbyNu0((tgradU() && dev(twoSymm(tgradU()))));
+    volScalarField::Internal G("kOmegaSST:G", nut_ * GbyNu0);
+
+    volScalarField CDkOmega(
+        (scalar(2) * alphaOmega2_) * (fvc::grad(k_) & fvc::grad(omega_)) / omega_);
+
+    volScalarField F1(this->F1(CDkOmega));
+    volScalarField F23(this->F23());
+
+    volScalarField::Internal gamma(this->gamma(F1));
+    volScalarField::Internal beta(this->beta(F1));
+
+    volScalarField::Internal P = phase_() * rho_() * gamma * GbyNu(GbyNu0, F23(), S2());
+    volScalarField::Internal D = phase_() * rho_() * beta * sqr(omega_());
+
+    forAll(P, cellI)
+    {
+        PoD[cellI] = P[cellI] / (D[cellI] + 1e-16);
+    }
+}
+
+void DAkOmegaSST::getTurbConvOverProd(scalarList& CoP) const
+{
+    /*
+    Description:
+        Return the value of the convective over production term from the turbulence model 
+    */
+
+    tmp<volTensorField> tgradU = fvc::grad(U_);
+    volScalarField S2(2 * magSqr(symm(tgradU())));
+    volScalarField::Internal GbyNu0((tgradU() && dev(twoSymm(tgradU()))));
+    volScalarField::Internal G("kOmegaSST:G", nut_ * GbyNu0);
+
+    volScalarField CDkOmega(
+        (scalar(2) * alphaOmega2_) * (fvc::grad(k_) & fvc::grad(omega_)) / omega_);
+
+    volScalarField F1(this->F1(CDkOmega));
+    volScalarField F23(this->F23());
+
+    volScalarField::Internal gamma(this->gamma(F1));
+    volScalarField::Internal beta(this->beta(F1));
+
+    volScalarField::Internal P = phase_() * rho_() * gamma * GbyNu(GbyNu0, F23(), S2());
+    volScalarField C = fvc::div(phaseRhoPhi_, omega_);
+
+    forAll(P, cellI)
+    {
+        CoP[cellI] = C[cellI] / (P[cellI] + 1e-16);
+    }
+}
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 } // End namespace Foam

diff --git a/src/adjoint/DAModel/DATurbulenceModel/DAkOmegaSST.H b/src/adjoint/DAModel/DATurbulenceModel/DAkOmegaSST.H
@@ -238,6 +238,12 @@ public:
         scalarField& diag,
         scalarField& upper,
         scalarField& lower);
+
+    /// return the value of the destruction term from the turbulence model
+    virtual void getTurbProdOverDestruct(scalarList& PoD) const;
+
+    /// return the value of the convective over production term from the turbulence model
+    virtual void getTurbConvOverProd(scalarList& CoP) const;
 };
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //