Merge pull request #92 from gerlero/flow

Update flow models
gerlero · Jan 10, 2024 · 0e34871 · 0e34871
2 parents fd7feee + 9e22d9b
commit 0e34871
Show file tree

Hide file tree

Showing 26 changed files with 122 additions and 101 deletions.
diff --git a/libraries/unsaturatedFlowModels/BrooksAndCorey/BrooksAndCorey.C b/libraries/unsaturatedFlowModels/BrooksAndCorey/BrooksAndCorey.C
@@ -1,5 +1,6 @@
 #include "BrooksAndCorey.H"
 #include "constantFields.H"
+#include "dictionaries.H"
 #include "porousMedium.H"
 #include "fluidPhase.H"
 #include "phaseFractionField.H"
@@ -27,15 +28,17 @@ Foam::Pmt::unsaturatedFlowModels::BrooksAndCorey::BrooksAndCorey
 (
     const porousMedium& medium,
     const fluidPhase& phase,
-    const dictionary& coeffs
+    const phaseFractionField& frac,
+    const dictionary& transportProperties
 )
 :
     medium_{medium},
     phase_{phase},
-    pc0_{constantFields::read("pc0", medium.mesh(), dimPressure, coeffs)},
-    alpha_{constantFields::read("alpha", medium.mesh(), dimless, coeffs)},
-    n_{constantFields::read("n", medium.mesh(), dimless, coeffs)},
-    l_{constantFields::readOrDefault("l", alpha_, coeffs)}
+    frac_{frac},
+    pc0_{constantFields::read("pc0", medium.mesh(), dimPressure,dictionaries::subOrNullDictRef(transportProperties, "BrooksAndCoreyCoeffs"))},
+    alpha_{constantFields::read("alpha", medium.mesh(), dimless, dictionaries::subOrNullDictRef(transportProperties, "BrooksAndCoreyCoeffs"))},
+    n_{constantFields::read("n", medium.mesh(), dimless, dictionaries::subOrNullDictRef(transportProperties, "BrooksAndCoreyCoeffs"))},
+    l_{constantFields::readOrDefault("l", alpha_, dictionaries::subOrNullDictRef(transportProperties, "BrooksAndCoreyCoeffs"))}
 {
     Info<< nl
         << typeName << " model" << nl
@@ -49,17 +52,17 @@ Foam::Pmt::unsaturatedFlowModels::BrooksAndCorey::BrooksAndCorey
 }
 
 Foam::tmp<Foam::volScalarField>
-Foam::Pmt::unsaturatedFlowModels::BrooksAndCorey::C(const phaseFractionField& frac)
+Foam::Pmt::unsaturatedFlowModels::BrooksAndCorey::C()
 {
-    volScalarField p{-pc0_*pow(frac.eff(), -alpha_)};
+    volScalarField p{-pc0_*pow(frac_.eff(), -alpha_)};
 
-    return -neg(p + pc0_)*(frac.max() - frac.min())/(alpha_*p*pow(-p/pc0_, 1/alpha_));
+    return -neg(p + pc0_)*(frac_.max() - frac_.min())/(alpha_*p*pow(-p/pc0_, 1/alpha_));
 }
 
 Foam::tmp<Foam::volScalarField>
-Foam::Pmt::unsaturatedFlowModels::BrooksAndCorey::M(const phaseFractionField& frac)
+Foam::Pmt::unsaturatedFlowModels::BrooksAndCorey::M()
 {
-    volScalarField eff{frac.eff()};
+    volScalarField eff{frac_.eff()};
 
     return medium_.K()/phase_.mu()*pow(eff, n_ + l_ - 1)
          + neg0(1 - eff)*dimensionedScalar{dimArea/dimDynamicViscosity, One};

diff --git a/libraries/unsaturatedFlowModels/BrooksAndCorey/BrooksAndCorey.H b/libraries/unsaturatedFlowModels/BrooksAndCorey/BrooksAndCorey.H
@@ -32,16 +32,18 @@ public:
     (
         const porousMedium& medium,
         const fluidPhase& phase,
-        const dictionary& coeffs
+        const phaseFractionField& frac,
+        const dictionary& transportProperties
     );
 
-    tmp<volScalarField> C(const phaseFractionField& frac) override;
+    tmp<volScalarField> C() override;
 
-    tmp<volScalarField> M(const phaseFractionField& frac) override;
+    tmp<volScalarField> M() override;
 
 private:
     const porousMedium& medium_;
     const fluidPhase& phase_;
+    const phaseFractionField& frac_;
     const volScalarField pc0_;
     const volScalarField alpha_;
     const volScalarField n_;

diff --git a/libraries/unsaturatedFlowModels/LETd/LETd.C b/libraries/unsaturatedFlowModels/LETd/LETd.C
@@ -1,5 +1,6 @@
 #include "LETd.H"
 #include "constantFields.H"
+#include "dictionaries.H"
 #include "porousMedium.H"
 #include "phaseFractionField.H"
 
@@ -23,13 +24,15 @@ Foam::Pmt::unsaturatedFlowModels::LETd::LETd
 (
     const porousMedium& medium,
     const fluidPhase&,
-    const dictionary& coeffs
+    const phaseFractionField& frac,
+    const dictionary& transportProperties
 )
 :
-    Dwt_{constantFields::read("Dwt", medium.mesh(), dimViscosity, coeffs)},
-    L_{"L", dimless, coeffs},
-    E_{"E", dimless, coeffs},
-    T_{"T", dimless, coeffs}
+    frac_{frac},
+    Dwt_{constantFields::read("Dwt", medium.mesh(), dimViscosity, dictionaries::subOrNullDictRef(transportProperties, "LETCoeffs"))},
+    L_{"L", dimless, dictionaries::subOrNullDictRef(transportProperties, "LETCoeffs")},
+    E_{"E", dimless, dictionaries::subOrNullDictRef(transportProperties, "LETCoeffs")},
+    T_{"T", dimless, dictionaries::subOrNullDictRef(transportProperties, "LETCoeffs")}
 {
     Info<< nl
         << typeName << " model" << nl
@@ -43,9 +46,9 @@ Foam::Pmt::unsaturatedFlowModels::LETd::LETd
 }
 
 Foam::tmp<Foam::volScalarField>
-Foam::Pmt::unsaturatedFlowModels::LETd::D(const phaseFractionField& frac)
+Foam::Pmt::unsaturatedFlowModels::LETd::D()
 {
-    volScalarField Swp{frac.eff()};
+    volScalarField Swp{frac_.eff()};
 
     return Dwt_*pow(Swp, L_)/(pow(Swp, L_) + E_*pow(1 - Swp, T_));
 }
diff --git a/libraries/unsaturatedFlowModels/LETd/LETd.H b/libraries/unsaturatedFlowModels/LETd/LETd.H
@@ -28,12 +28,14 @@ public:
     (
         const porousMedium& medium,
         const fluidPhase& phase,
-        const dictionary& coeffs
+        const phaseFractionField& frac,
+        const dictionary& transportProperties
     );
 
-    tmp<volScalarField> D(const phaseFractionField& frac) override;
+    tmp<volScalarField> D() override;
 
 private:
+    const phaseFractionField& frac_;
     const volScalarField Dwt_;
     const dimensionedScalar L_;
     const dimensionedScalar E_;

diff --git a/libraries/unsaturatedFlowModels/LETxs/LETxs.C b/libraries/unsaturatedFlowModels/LETxs/LETxs.C
@@ -1,5 +1,6 @@
 #include "LETxs.H"
 #include "constantFields.H"
+#include "dictionaries.H"
 #include "porousMedium.H"
 #include "fluidPhase.H"
 #include "phaseFractionField.H"
@@ -24,18 +25,20 @@ Foam::Pmt::unsaturatedFlowModels::LETxs::LETxs
 (
     const porousMedium& medium,
     const fluidPhase& phase,
-    const dictionary& coeffs
+    const phaseFractionField& frac,
+    const dictionary& transportProperties
 )
 :
     medium_{medium},
     phase_{phase},
-    pc0_{constantFields::read("pc0", medium.mesh(), dimPressure, coeffs)},
-    Lw_{"Lw", dimless, coeffs},
-    Ew_{"Ew", dimless, coeffs},
-    Tw_{"Tw", dimless, coeffs},
-    Ls_{"Ls", dimless, coeffs},
-    Es_{"Es", dimless, coeffs},
-    Ts_{"Ts", dimless, coeffs}
+    frac_{frac},
+    pc0_{constantFields::read("pc0", medium.mesh(), dimPressure, dictionaries::subOrNullDictRef(transportProperties, "LETCoeffs"))},
+    Lw_{"Lw", dimless, dictionaries::subOrNullDictRef(transportProperties, "LETCoeffs")},
+    Ew_{"Ew", dimless, dictionaries::subOrNullDictRef(transportProperties, "LETCoeffs")},
+    Tw_{"Tw", dimless, dictionaries::subOrNullDictRef(transportProperties, "LETCoeffs")},
+    Ls_{"Ls", dimless, dictionaries::subOrNullDictRef(transportProperties, "LETCoeffs")},
+    Es_{"Es", dimless, dictionaries::subOrNullDictRef(transportProperties, "LETCoeffs")},
+    Ts_{"Ts", dimless, dictionaries::subOrNullDictRef(transportProperties, "LETCoeffs")}
 {
     Info<< nl
         << typeName << " model" << nl
@@ -52,13 +55,13 @@ Foam::Pmt::unsaturatedFlowModels::LETxs::LETxs
 }
 
 Foam::tmp<Foam::volScalarField>
-Foam::Pmt::unsaturatedFlowModels::LETxs::D(const phaseFractionField& frac)
+Foam::Pmt::unsaturatedFlowModels::LETxs::D()
 {
-    volScalarField Swp{frac.eff()};
-    auto Swir = frac.min()/frac.max();
+    volScalarField Swp{frac_.eff()};
+    auto Swir = frac_.min()/frac_.max();
 
     return
         medium_.K()/phase_.mu()*pc0_*Es_*pow(Swp,Lw_)*pow(Swp,Ts_)*pow(1 - Swp,Ls_)*(Ls_*Swp - Swp*Ts_ + Ts_)
-       /(frac.max()*Swp*(Swir - 1)*(Swp - 1)*pow(Es_*pow(Swp,Ts_) + pow(1 - Swp,Ls_),2)*(Ew_*pow(1 - Swp,Tw_) + pow(Swp,Lw_)));
+       /(frac_.max()*Swp*(Swir - 1)*(Swp - 1)*pow(Es_*pow(Swp,Ts_) + pow(1 - Swp,Ls_),2)*(Ew_*pow(1 - Swp,Tw_) + pow(Swp,Lw_)));
 
 }
diff --git a/libraries/unsaturatedFlowModels/LETxs/LETxs.H b/libraries/unsaturatedFlowModels/LETxs/LETxs.H
@@ -28,14 +28,16 @@ public:
     (
         const porousMedium& medium,
         const fluidPhase& phase,
-        const dictionary& coeffs
+        const phaseFractionField& frac,
+        const dictionary& transportProperties
     );
 
-    tmp<volScalarField> D(const phaseFractionField& frac) override;
+    tmp<volScalarField> D() override;
 
 private:
     const porousMedium& medium_;
     const fluidPhase& phase_;
+    const phaseFractionField& frac_;
     const volScalarField pc0_;
     const dimensionedScalar Lw_;
     const dimensionedScalar Ew_;

diff --git a/libraries/unsaturatedFlowModels/RichardsModel/RichardsModel.C b/libraries/unsaturatedFlowModels/RichardsModel/RichardsModel.C
@@ -14,7 +14,7 @@ defineRunTimeSelectionTable(RichardsModel, dictionary);
 }
 }
 
-Foam::tmp<Foam::volScalarField> Foam::Pmt::RichardsModel::D(const phaseFractionField& frac)
+Foam::tmp<Foam::volScalarField> Foam::Pmt::RichardsModel::D()
 {
-    return M(frac)/C(frac);
+    return M()/C();
 }
diff --git a/libraries/unsaturatedFlowModels/RichardsModel/RichardsModel.H b/libraries/unsaturatedFlowModels/RichardsModel/RichardsModel.H
@@ -32,23 +32,25 @@ public:
         (
             const porousMedium& medium,
             const fluidPhase& phase,
-            const dictionary& coeffs
+            const phaseFractionField& frac,
+            const dictionary& transportProperties
         ),
-        (medium, phase, coeffs)
+        (medium, phase, frac, transportProperties)
     );
 
     static autoPtr<RichardsModel> New
     (
         const porousMedium& medium,
         const fluidPhase& phase,
+        const phaseFractionField& frac,
         const dictionary& transportProperties
     );
 
-    virtual tmp<volScalarField> C(const phaseFractionField& frac) = 0;
+    virtual tmp<volScalarField> C() = 0;
 
-    virtual tmp<volScalarField> M(const phaseFractionField& frac) = 0;
+    virtual tmp<volScalarField> M() = 0;
 
-    tmp<volScalarField> D(const phaseFractionField& frac) override;
+    tmp<volScalarField> D() override;
 };
 
 }

diff --git a/libraries/unsaturatedFlowModels/RichardsModel/RichardsModelNew.C b/libraries/unsaturatedFlowModels/RichardsModel/RichardsModelNew.C
@@ -5,8 +5,9 @@ Foam::autoPtr<Foam::Pmt::RichardsModel> Foam::Pmt::RichardsModel::New
 (
     const porousMedium& medium,
     const fluidPhase& phase,
+    const phaseFractionField& frac,
     const dictionary& transportProperties
 )
 {
-    return Detail::unsaturatedFlowModelNew<RichardsModel>(medium, phase, transportProperties);
+    return Detail::unsaturatedFlowModelNew<RichardsModel>(medium, phase, frac, transportProperties);
 }
diff --git a/libraries/unsaturatedFlowModels/VanGenuchten/VanGenuchten.C b/libraries/unsaturatedFlowModels/VanGenuchten/VanGenuchten.C
@@ -1,5 +1,6 @@
 #include "VanGenuchten.H"
 #include "constantFields.H"
+#include "dictionaries.H"
 #include "porousMedium.H"
 #include "fluidPhase.H"
 #include "phaseFractionField.H"
@@ -27,26 +28,28 @@ Foam::Pmt::unsaturatedFlowModels::VanGenuchten::VanGenuchten
 (
     const porousMedium& medium,
     const fluidPhase& phase,
-    const dictionary& coeffs
+    const phaseFractionField& frac,
+    const dictionary& transportProperties
 )
 :
     medium_{medium},
     phase_{phase},
-    pc0_{constantFields::read("pc0", medium.mesh(), dimPressure, coeffs)},
+    frac_{frac},
+    pc0_{constantFields::read("pc0", medium.mesh(), dimPressure, dictionaries::subOrNullDictRef(transportProperties, "VanGenuchtenCoeffs"))},
     m_{
         [&]
         {
-            if (auto m = constantFields::readIfPresent("m", medium.mesh(), dimless, coeffs))
+            if (auto m = constantFields::readIfPresent("m", medium.mesh(), dimless, dictionaries::subOrNullDictRef(transportProperties, "VanGenuchtenCoeffs")))
             {
                 return m;
             }
 
-            auto n = constantFields::read("n", medium.mesh(), dimless, coeffs);
+            auto n = constantFields::read("n", medium.mesh(), dimless, dictionaries::subOrNullDictRef(transportProperties, "VanGenuchtenCoeffs"));
 
             return volScalarField::New("m", 1 - 1/n);
         }()
     },
-    l_{constantFields::readOrDefault("l", medium.mesh(), dimless, 0.5, coeffs)}
+    l_{constantFields::readOrDefault("l", medium.mesh(), dimless, 0.5, dictionaries::subOrNullDictRef(transportProperties, "VanGenuchtenCoeffs"))}
 {
     Info<< nl
         << typeName << " model" << nl
@@ -60,17 +63,17 @@ Foam::Pmt::unsaturatedFlowModels::VanGenuchten::VanGenuchten
 }
 
 Foam::tmp<Foam::volScalarField>
-Foam::Pmt::unsaturatedFlowModels::VanGenuchten::C(const phaseFractionField& frac)
+Foam::Pmt::unsaturatedFlowModels::VanGenuchten::C()
 {
-    volScalarField Se{frac.eff()};
+    volScalarField Se{frac_.eff()};
 
-    return 1/pc0_*m_/(1 - m_)*(frac.max() - frac.min())*pow(Se, 1/m_)*pow((1 - pow(Se, 1/m_)), m_);
+    return 1/pc0_*m_/(1 - m_)*(frac_.max() - frac_.min())*pow(Se, 1/m_)*pow((1 - pow(Se, 1/m_)), m_);
 }
 
 Foam::tmp<Foam::volScalarField>
-Foam::Pmt::unsaturatedFlowModels::VanGenuchten::M(const phaseFractionField& frac)
+Foam::Pmt::unsaturatedFlowModels::VanGenuchten::M()
 {
-    volScalarField Se{frac.eff()};
+    volScalarField Se{frac_.eff()};
 
     return
         pos(1 - Se)*medium_.K()/phase_.mu()*pow(Se, l_)*pow(1 - pow(1 - pow(Se, 1/m_), m_), 2)

diff --git a/libraries/unsaturatedFlowModels/VanGenuchten/VanGenuchten.H b/libraries/unsaturatedFlowModels/VanGenuchten/VanGenuchten.H
@@ -31,16 +31,18 @@ public:
     (
         const porousMedium& medium,
         const fluidPhase& phase,
-        const dictionary& coeffs
+        const phaseFractionField& frac,
+        const dictionary& transportProperties
     );
 
-    tmp<volScalarField> C(const phaseFractionField& frac) override;
+    tmp<volScalarField> C() override;
 
-    tmp<volScalarField> M(const phaseFractionField& frac) override;
+    tmp<volScalarField> M() override;
 
 private:
     const porousMedium& medium_;
     const fluidPhase& phase_;
+    const phaseFractionField& frac_;
     const volScalarField pc0_;
     const volScalarField m_;
     const volScalarField l_;

diff --git a/libraries/unsaturatedFlowModels/capillarityModel/capillarityModel.H b/libraries/unsaturatedFlowModels/capillarityModel/capillarityModel.H
@@ -30,14 +30,15 @@ public:
         (
             const porousMedium& medium,
             const fluidPhase& phase,
-            const dictionary& coeffs
+            const phaseFractionField& frac,
+            const dictionary& transportProperties
         ),
-        (medium, phase, coeffs)
+        (medium, phase, frac, transportProperties)
     );
 
     capillarityModel() = default;
 
-    virtual tmp<volScalarField> C(const phaseFractionField& frac) = 0;
+    virtual tmp<volScalarField> C() = 0;
 
     virtual ~capillarityModel() = default;
 

diff --git a/libraries/unsaturatedFlowModels/mixedRichardsModel/mixedRichardsModel.C b/libraries/unsaturatedFlowModels/mixedRichardsModel/mixedRichardsModel.C
@@ -28,13 +28,13 @@ Foam::Pmt::mixedRichardsModel::mixedRichardsModel
 }
 
 Foam::tmp<Foam::volScalarField>
-Foam::Pmt::mixedRichardsModel::C(const phaseFractionField& frac)
+Foam::Pmt::mixedRichardsModel::C()
 {
-    return capillarityModel_->C(frac);
+    return capillarityModel_->C();
 }
 
 Foam::tmp<Foam::volScalarField>
-Foam::Pmt::mixedRichardsModel::M(const phaseFractionField& frac)
+Foam::Pmt::mixedRichardsModel::M()
 {
-    return relativePermeabilityModel_->M(frac);
+    return relativePermeabilityModel_->M();
 }
diff --git a/libraries/unsaturatedFlowModels/mixedRichardsModel/mixedRichardsModel.H b/libraries/unsaturatedFlowModels/mixedRichardsModel/mixedRichardsModel.H
@@ -28,9 +28,9 @@ public:
         autoPtr<relativePermeabilityModel> relativePermeability
     );
 
-    tmp<volScalarField> C(const phaseFractionField& frac) override;
+    tmp<volScalarField> C() override;
 
-    tmp<volScalarField> M(const phaseFractionField& frac) override;
+    tmp<volScalarField> M() override;
 
 private:
     autoPtr<capillarityModel> capillarityModel_;