From f2cd716decdd3efdc763c62199248ed324130f63 Mon Sep 17 00:00:00 2001 From: Henry Weller Date: Sun, 30 Oct 2022 06:34:26 +0000 Subject: [PATCH] fvmLaplacian: Added laplacianCorrection functions and updated all thermal transport implementations Now that thermal transport is implemented as an energy implicit correction on an explicit temperature gradient formulation it is more efficient if the implicit correction contains only the implicit terms of the matrix and not the explicit non-orthogonal or anisotropic correction terms which are cancelled anyway when the evaluation of the matrix for the current state is subtracted. The new fvm::laplacianCorrection functions provide a convenient mechanism to efficiently evaluate only the implicit correction to the laplacian and is now used in all the thermophysicalTransportModels. --- .../StationaryPhaseModel.C | 16 +++- .../fluid/laminar/Fickian/Fickian.C | 2 +- .../fluid/laminar/Fourier/Fourier.C | 8 +- .../laminar/MaxwellStefan/MaxwellStefan.C | 2 +- .../eddyDiffusivity/eddyDiffusivity.C | 4 +- .../nonUnityLewisEddyDiffusivity.C | 2 +- .../solid/anisotropic/anisotropic.C | 12 ++- .../solid/isotropic/isotropic.C | 10 +-- .../finiteVolume/fvm/fvmLaplacian.C | 83 +++++++++++++++++-- .../finiteVolume/fvm/fvmLaplacian.H | 30 +++++++ 10 files changed, 137 insertions(+), 32 deletions(-) diff --git a/applications/solvers/multiphase/multiphaseEulerFoam/phaseSystems/phaseModel/StationaryPhaseModel/StationaryPhaseModel.C b/applications/solvers/multiphase/multiphaseEulerFoam/phaseSystems/phaseModel/StationaryPhaseModel/StationaryPhaseModel.C index b5c9b0607e..257f2451af 100644 --- a/applications/solvers/multiphase/multiphaseEulerFoam/phaseSystems/phaseModel/StationaryPhaseModel/StationaryPhaseModel.C +++ b/applications/solvers/multiphase/multiphaseEulerFoam/phaseSystems/phaseModel/StationaryPhaseModel/StationaryPhaseModel.C @@ -316,11 +316,21 @@ Foam::StationaryPhaseModel::divq(volScalarField& he) const { const volScalarField& alpha = *this; - return -fvm::laplacian + const surfaceScalarField alphaKappa ( - fvc::interpolate(alpha)*fvc::interpolate(this->thermo().alphahe()), - he + alpha.name() + '*' + this->thermo().kappa().name(), + fvc::interpolate(alpha)*fvc::interpolate(this->thermo().kappa()) ); + + // Return heat flux source as an implicit energy correction + // to the temperature gradient flux + return + -fvc::laplacian(alphaKappa, this->thermo().T()) + -fvm::laplacianCorrection + ( + alphaKappa/fvc::interpolate(this->thermo().Cpv()), + he + ); } diff --git a/src/ThermophysicalTransportModels/fluid/laminar/Fickian/Fickian.C b/src/ThermophysicalTransportModels/fluid/laminar/Fickian/Fickian.C index 2bfa8e5a4b..75fc64ca5e 100644 --- a/src/ThermophysicalTransportModels/fluid/laminar/Fickian/Fickian.C +++ b/src/ThermophysicalTransportModels/fluid/laminar/Fickian/Fickian.C @@ -309,7 +309,7 @@ tmp Fickian::divq const PtrList& Y = composition.Y(); tmpDivq.ref() -= - correction(fvm::laplacian(this->alpha()*this->alphaEff(), he)); + fvm::laplacianCorrection(this->alpha()*this->alphaEff(), he); surfaceScalarField sumJ ( diff --git a/src/ThermophysicalTransportModels/fluid/laminar/Fourier/Fourier.C b/src/ThermophysicalTransportModels/fluid/laminar/Fourier/Fourier.C index 62a627ba7d..9e2c6e0a20 100644 --- a/src/ThermophysicalTransportModels/fluid/laminar/Fourier/Fourier.C +++ b/src/ThermophysicalTransportModels/fluid/laminar/Fourier/Fourier.C @@ -130,8 +130,12 @@ Fourier::divq(volScalarField& he) const // Return heat flux source as an implicit energy correction // to the temperature gradient flux return - -correction(fvm::laplacian(this->alpha()*thermo.alphahe(), he)) - -fvc::laplacian(this->alpha()*thermo.kappa(), thermo.T()); + -fvc::laplacian(this->alpha()*thermo.kappa(), thermo.T()) + -fvm::laplacianCorrection + ( + this->alpha()*thermo.kappa()/thermo.Cpv(), + he + ); } diff --git a/src/ThermophysicalTransportModels/fluid/laminar/MaxwellStefan/MaxwellStefan.C b/src/ThermophysicalTransportModels/fluid/laminar/MaxwellStefan/MaxwellStefan.C index e15bca5f28..1edcfc90f4 100644 --- a/src/ThermophysicalTransportModels/fluid/laminar/MaxwellStefan/MaxwellStefan.C +++ b/src/ThermophysicalTransportModels/fluid/laminar/MaxwellStefan/MaxwellStefan.C @@ -310,7 +310,7 @@ tmp MaxwellStefan::divq const PtrList& Y = composition.Y(); tmpDivq.ref() -= - correction(fvm::laplacian(this->alpha()*this->alphaEff(), he)); + fvm::laplacianCorrection(this->alpha()*this->alphaEff(), he); surfaceScalarField sumJ ( diff --git a/src/ThermophysicalTransportModels/fluid/turbulence/eddyDiffusivity/eddyDiffusivity.C b/src/ThermophysicalTransportModels/fluid/turbulence/eddyDiffusivity/eddyDiffusivity.C index 57824b631d..fb3e160df0 100644 --- a/src/ThermophysicalTransportModels/fluid/turbulence/eddyDiffusivity/eddyDiffusivity.C +++ b/src/ThermophysicalTransportModels/fluid/turbulence/eddyDiffusivity/eddyDiffusivity.C @@ -163,8 +163,8 @@ eddyDiffusivity::divq // Return heat flux source as an implicit energy correction // to the temperature gradient flux return - -correction(fvm::laplacian(this->alpha()*this->alphaEff(), he)) - -fvc::laplacian(this->alpha()*this->kappaEff(), this->thermo().T()); + -fvc::laplacian(this->alpha()*this->kappaEff(), this->thermo().T()) + -fvm::laplacianCorrection(this->alpha()*this->alphaEff(), he); } diff --git a/src/ThermophysicalTransportModels/fluid/turbulence/nonUnityLewisEddyDiffusivity/nonUnityLewisEddyDiffusivity.C b/src/ThermophysicalTransportModels/fluid/turbulence/nonUnityLewisEddyDiffusivity/nonUnityLewisEddyDiffusivity.C index b3e9b7a802..491d33da70 100644 --- a/src/ThermophysicalTransportModels/fluid/turbulence/nonUnityLewisEddyDiffusivity/nonUnityLewisEddyDiffusivity.C +++ b/src/ThermophysicalTransportModels/fluid/turbulence/nonUnityLewisEddyDiffusivity/nonUnityLewisEddyDiffusivity.C @@ -162,7 +162,7 @@ nonUnityLewisEddyDiffusivity::divq const PtrList& Y = composition.Y(); tmpDivq.ref() -= - correction(fvm::laplacian(this->alpha()*this->alphaEff(), he)); + fvm::laplacianCorrection(this->alpha()*this->alphaEff(), he); surfaceScalarField hGradY ( diff --git a/src/ThermophysicalTransportModels/solid/anisotropic/anisotropic.C b/src/ThermophysicalTransportModels/solid/anisotropic/anisotropic.C index f8368c9bed..9a933d7626 100644 --- a/src/ThermophysicalTransportModels/solid/anisotropic/anisotropic.C +++ b/src/ThermophysicalTransportModels/solid/anisotropic/anisotropic.C @@ -401,19 +401,17 @@ Foam::solidThermophysicalTransportModels::anisotropic::divq { const solidThermo& thermo = this->thermo(); const volSymmTensorField Kappa(this->Kappa()); + const surfaceVectorField& Sf = thermo.mesh().Sf(); + const surfaceScalarField& magSf = thermo.mesh().magSf(); // Return heat flux source as an implicit energy correction // to the temperature gradient flux return -fvc::laplacian(Kappa, thermo.T()) - -correction + -fvm::laplacianCorrection ( - fvm::laplacian - ( - Kappa/thermo.Cv(), - e, - "laplacian(alphae,e)" - ) + (Sf & fvc::interpolate(Kappa/thermo.Cv()) & Sf)/sqr(magSf), + e ); } diff --git a/src/ThermophysicalTransportModels/solid/isotropic/isotropic.C b/src/ThermophysicalTransportModels/solid/isotropic/isotropic.C index be11430a5f..91f8633023 100644 --- a/src/ThermophysicalTransportModels/solid/isotropic/isotropic.C +++ b/src/ThermophysicalTransportModels/solid/isotropic/isotropic.C @@ -105,15 +105,7 @@ Foam::solidThermophysicalTransportModels::isotropic::divq // to the temperature gradient flux return -fvc::laplacian(kappa(), thermo.T()) - -correction - ( - fvm::laplacian - ( - kappa()/thermo.Cv(), - e, - "laplacian(alphae,e)" - ) - ); + -fvm::laplacianCorrection(kappa()/thermo.Cv(), e); } diff --git a/src/finiteVolume/finiteVolume/fvm/fvmLaplacian.C b/src/finiteVolume/finiteVolume/fvm/fvmLaplacian.C index 663b5f60c9..0c41f4835b 100644 --- a/src/finiteVolume/finiteVolume/fvm/fvmLaplacian.C +++ b/src/finiteVolume/finiteVolume/fvm/fvmLaplacian.C @@ -27,6 +27,7 @@ License #include "surfaceFields.H" #include "fvMatrix.H" #include "laplacianScheme.H" +#include "gaussLaplacianScheme.H" // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // @@ -229,9 +230,9 @@ laplacian const word& name ) { - tmp> Laplacian(fvm::laplacian(tgamma(), vf, name)); + tmp> tLaplacian(fvm::laplacian(tgamma(), vf, name)); tgamma.clear(); - return Laplacian; + return tLaplacian; } @@ -260,9 +261,9 @@ laplacian const GeometricField& vf ) { - tmp> Laplacian(fvm::laplacian(tgamma(), vf)); + tmp> tLaplacian(fvm::laplacian(tgamma(), vf)); tgamma.clear(); - return Laplacian; + return tLaplacian; } @@ -325,9 +326,79 @@ laplacian const GeometricField& vf ) { - tmp> tfvm(fvm::laplacian(tGamma(), vf)); + tmp> tLaplacian(fvm::laplacian(tGamma(), vf)); tGamma.clear(); - return tfvm; + return tLaplacian; +} + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +template +tmp> +laplacianCorrection +( + const GeometricField& gamma, + const GeometricField& vf +) +{ + return fvm::laplacianCorrection(fvc::interpolate(gamma), vf); +} + + +template +tmp> +laplacianCorrection +( + const tmp>& tgamma, + const GeometricField& vf +) +{ + tmp> tLaplacianCorrection + ( + fvm::laplacianCorrection(tgamma(), vf) + ); + tgamma.clear(); + return tLaplacianCorrection; +} + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +template +tmp> +laplacianCorrection +( + const GeometricField& gamma, + const GeometricField& vf +) +{ + return correction + ( + fv::gaussLaplacianScheme::fvmLaplacianUncorrected + ( + gamma*vf.mesh().magSf(), + vf.mesh().nonOrthDeltaCoeffs(), + vf + ) + ); +} + + +template +tmp> +laplacianCorrection +( + const tmp>& tGamma, + const GeometricField& vf +) +{ + tmp> tLaplacianCorrection + ( + fvm::laplacianCorrection(tGamma(), vf) + ); + tGamma.clear(); + return tLaplacianCorrection; } diff --git a/src/finiteVolume/finiteVolume/fvm/fvmLaplacian.H b/src/finiteVolume/finiteVolume/fvm/fvmLaplacian.H index e5652a7f4b..b831e65403 100644 --- a/src/finiteVolume/finiteVolume/fvm/fvmLaplacian.H +++ b/src/finiteVolume/finiteVolume/fvm/fvmLaplacian.H @@ -176,6 +176,36 @@ namespace fvm const tmp>&, const GeometricField& ); + + + template + tmp> laplacianCorrection + ( + const GeometricField&, + const GeometricField& + ); + + template + tmp> laplacianCorrection + ( + const tmp>&, + const GeometricField& + ); + + + template + tmp> laplacianCorrection + ( + const GeometricField&, + const GeometricField& + ); + + template + tmp> laplacianCorrection + ( + const tmp>&, + const GeometricField& + ); }