From 473cf07ddefb33f81ea4e11ae93d99b64fbacaf8 Mon Sep 17 00:00:00 2001
From: Henry <Henry>
Date: Wed, 1 May 2013 15:32:36 +0100
Subject: [PATCH] MULES: Split files to separate the explicit and implicit
forms IMULES: MULES for implicit solution
---
.../interPhaseChangeFoam.C | 2 +-
.../multiphaseSystem/multiphaseSystem.C | 3 +
.../multiphaseMixture/multiphaseMixture.C | 5 +-
src/finiteVolume/Make/files | 1 +
.../fvMatrices/solvers/MULES/IMULES.C | 51 ++++
.../fvMatrices/solvers/MULES/IMULES.H | 87 +++++++
.../solvers/MULES/IMULESTemplates.C | 244 ++++++++++++++++++
.../fvMatrices/solvers/MULES/MULES.C | 32 ---
.../fvMatrices/solvers/MULES/MULES.H | 27 +-
.../fvMatrices/solvers/MULES/MULESTemplates.C | 216 ----------------
10 files changed, 395 insertions(+), 273 deletions(-)
create mode 100644 src/finiteVolume/fvMatrices/solvers/MULES/IMULES.C
create mode 100644 src/finiteVolume/fvMatrices/solvers/MULES/IMULES.H
create mode 100644 src/finiteVolume/fvMatrices/solvers/MULES/IMULESTemplates.C
diff --git a/applications/solvers/multiphase/interPhaseChangeFoam/interPhaseChangeFoam.C b/applications/solvers/multiphase/interPhaseChangeFoam/interPhaseChangeFoam.C
index bc6af055f0c..edfc34dfb83 100644
--- a/applications/solvers/multiphase/interPhaseChangeFoam/interPhaseChangeFoam.C
+++ b/applications/solvers/multiphase/interPhaseChangeFoam/interPhaseChangeFoam.C
@@ -41,7 +41,7 @@ Description
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
-#include "MULES.H"
+#include "IMULES.H"
#include "subCycle.H"
#include "interfaceProperties.H"
#include "phaseChangeTwoPhaseMixture.H"
diff --git a/applications/solvers/multiphase/multiphaseEulerFoam/multiphaseSystem/multiphaseSystem.C b/applications/solvers/multiphase/multiphaseEulerFoam/multiphaseSystem/multiphaseSystem.C
index acad49cd966..0cde87042d9 100644
--- a/applications/solvers/multiphase/multiphaseEulerFoam/multiphaseSystem/multiphaseSystem.C
+++ b/applications/solvers/multiphase/multiphaseEulerFoam/multiphaseSystem/multiphaseSystem.C
@@ -29,7 +29,10 @@ License
#include "Time.H"
#include "subCycle.H"
#include "MULES.H"
+#include "surfaceInterpolate.H"
+#include "fvcGrad.H"
#include "fvcSnGrad.H"
+#include "fvcDiv.H"
#include "fvcFlux.H"
#include "fvcAverage.H"
diff --git a/applications/solvers/multiphase/multiphaseInterFoam/multiphaseMixture/multiphaseMixture.C b/applications/solvers/multiphase/multiphaseInterFoam/multiphaseMixture/multiphaseMixture.C
index fb6283d1efe..b3c4ecb1115 100644
--- a/applications/solvers/multiphase/multiphaseInterFoam/multiphaseMixture/multiphaseMixture.C
+++ b/applications/solvers/multiphase/multiphaseInterFoam/multiphaseMixture/multiphaseMixture.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -28,7 +28,10 @@ License
#include "Time.H"
#include "subCycle.H"
#include "MULES.H"
+#include "surfaceInterpolate.H"
+#include "fvcGrad.H"
#include "fvcSnGrad.H"
+#include "fvcDiv.H"
#include "fvcFlux.H"
// * * * * * * * * * * * * * * * Static Member Data * * * * * * * * * * * * //
diff --git a/src/finiteVolume/Make/files b/src/finiteVolume/Make/files
index dd37136533e..ddf9564230f 100644
--- a/src/finiteVolume/Make/files
+++ b/src/finiteVolume/Make/files
@@ -216,6 +216,7 @@ fields/surfaceFields/surfaceFields.C
fvMatrices/fvMatrices.C
fvMatrices/fvScalarMatrix/fvScalarMatrix.C
fvMatrices/solvers/MULES/MULES.C
+fvMatrices/solvers/MULES/IMULES.C
fvMatrices/solvers/GAMGSymSolver/GAMGAgglomerations/faceAreaPairGAMGAgglomeration/faceAreaPairGAMGAgglomeration.C
interpolation = interpolation/interpolation
diff --git a/src/finiteVolume/fvMatrices/solvers/MULES/IMULES.C b/src/finiteVolume/fvMatrices/solvers/MULES/IMULES.C
new file mode 100644
index 00000000000..e05e46b7278
--- /dev/null
+++ b/src/finiteVolume/fvMatrices/solvers/MULES/IMULES.C
@@ -0,0 +1,51 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2011 OpenFOAM Foundation
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ 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/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "IMULES.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+void Foam::MULES::implicitSolve
+(
+ volScalarField& psi,
+ const surfaceScalarField& phi,
+ surfaceScalarField& phiPsi,
+ const scalar psiMax,
+ const scalar psiMin
+)
+{
+ implicitSolve
+ (
+ geometricOneField(),
+ psi,
+ phi,
+ phiPsi,
+ zeroField(), zeroField(),
+ psiMax, psiMin
+ );
+}
+
+
+// ************************************************************************* //
diff --git a/src/finiteVolume/fvMatrices/solvers/MULES/IMULES.H b/src/finiteVolume/fvMatrices/solvers/MULES/IMULES.H
new file mode 100644
index 00000000000..5306dab8433
--- /dev/null
+++ b/src/finiteVolume/fvMatrices/solvers/MULES/IMULES.H
@@ -0,0 +1,87 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2011 OpenFOAM Foundation
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ 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/>.
+
+Global
+ IMULES
+
+Description
+ IMULES: Multidimensional universal limiter with implicit solution.
+
+SourceFiles
+ IMULES.C
+ IMULESTemplates.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef IMULES_H
+#define IMULES_H
+
+#include "MULES.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace MULES
+{
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+template<class RhoType, class SpType, class SuType>
+void implicitSolve
+(
+ const RhoType& rho,
+ volScalarField& gamma,
+ const surfaceScalarField& phi,
+ surfaceScalarField& phiCorr,
+ const SpType& Sp,
+ const SuType& Su,
+ const scalar psiMax,
+ const scalar psiMin
+);
+
+void implicitSolve
+(
+ volScalarField& gamma,
+ const surfaceScalarField& phi,
+ surfaceScalarField& phiCorr,
+ const scalar psiMax,
+ const scalar psiMin
+);
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace MULES
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#ifdef NoRepository
+# include "IMULESTemplates.C"
+#endif
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/finiteVolume/fvMatrices/solvers/MULES/IMULESTemplates.C b/src/finiteVolume/fvMatrices/solvers/MULES/IMULESTemplates.C
new file mode 100644
index 00000000000..14749ff9b73
--- /dev/null
+++ b/src/finiteVolume/fvMatrices/solvers/MULES/IMULESTemplates.C
@@ -0,0 +1,244 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ 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/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "IMULES.H"
+#include "gaussConvectionScheme.H"
+#include "surfaceInterpolate.H"
+#include "fvmDdt.H"
+#include "fvmSup.H"
+#include "fvcDiv.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace MULES
+{
+ template<class RhoType>
+ inline tmp<surfaceScalarField> interpolate(const RhoType& rho)
+ {
+ notImplemented
+ (
+ "tmp<surfaceScalarField> interpolate(const RhoType& rho)"
+ );
+ return tmp<surfaceScalarField>(NULL);
+ }
+
+ template<>
+ inline tmp<surfaceScalarField> interpolate(const volScalarField& rho)
+ {
+ return fvc::interpolate(rho);
+ }
+}
+}
+
+
+template<class RhoType, class SpType, class SuType>
+void Foam::MULES::implicitSolve
+(
+ const RhoType& rho,
+ volScalarField& psi,
+ const surfaceScalarField& phi,
+ surfaceScalarField& phiPsi,
+ const SpType& Sp,
+ const SuType& Su,
+ const scalar psiMax,
+ const scalar psiMin
+)
+{
+ const fvMesh& mesh = psi.mesh();
+
+ const dictionary& MULEScontrols = mesh.solverDict(psi.name());
+
+ label maxIter
+ (
+ readLabel(MULEScontrols.lookup("maxIter"))
+ );
+
+ label nLimiterIter
+ (
+ readLabel(MULEScontrols.lookup("nLimiterIter"))
+ );
+
+ scalar maxUnboundedness
+ (
+ readScalar(MULEScontrols.lookup("maxUnboundedness"))
+ );
+
+ scalar CoCoeff
+ (
+ readScalar(MULEScontrols.lookup("CoCoeff"))
+ );
+
+ scalarField allCoLambda(mesh.nFaces());
+
+ {
+ slicedSurfaceScalarField CoLambda
+ (
+ IOobject
+ (
+ "CoLambda",
+ mesh.time().timeName(),
+ mesh,
+ IOobject::NO_READ,
+ IOobject::NO_WRITE,
+ false
+ ),
+ mesh,
+ dimless,
+ allCoLambda,
+ false // Use slices for the couples
+ );
+
+ if (phi.dimensions() == dimDensity*dimVelocity*dimArea)
+ {
+ tmp<surfaceScalarField> Cof =
+ mesh.time().deltaT()*mesh.surfaceInterpolation::deltaCoeffs()
+ *mag(phi/interpolate(rho))/mesh.magSf();
+
+ CoLambda == 1.0/max(CoCoeff*Cof, scalar(1));
+ }
+ else
+ {
+ tmp<surfaceScalarField> Cof =
+ mesh.time().deltaT()*mesh.surfaceInterpolation::deltaCoeffs()
+ *mag(phi)/mesh.magSf();
+
+ CoLambda == 1.0/max(CoCoeff*Cof, scalar(1));
+ }
+ }
+
+ scalarField allLambda(allCoLambda);
+ //scalarField allLambda(mesh.nFaces(), 1.0);
+
+ slicedSurfaceScalarField lambda
+ (
+ IOobject
+ (
+ "lambda",
+ mesh.time().timeName(),
+ mesh,
+ IOobject::NO_READ,
+ IOobject::NO_WRITE,
+ false
+ ),
+ mesh,
+ dimless,
+ allLambda,
+ false // Use slices for the couples
+ );
+
+ linear<scalar> CDs(mesh);
+ upwind<scalar> UDs(mesh, phi);
+ //fv::uncorrectedSnGrad<scalar> snGrads(mesh);
+
+ fvScalarMatrix psiConvectionDiffusion
+ (
+ fvm::ddt(rho, psi)
+ + fv::gaussConvectionScheme<scalar>(mesh, phi, UDs).fvmDiv(phi, psi)
+ //- fv::gaussLaplacianScheme<scalar, scalar>(mesh, CDs, snGrads)
+ //.fvmLaplacian(Dpsif, psi)
+ - fvm::Sp(Sp, psi)
+ - Su
+ );
+
+ surfaceScalarField phiBD(psiConvectionDiffusion.flux());
+
+ surfaceScalarField& phiCorr = phiPsi;
+ phiCorr -= phiBD;
+
+ for (label i=0; i<maxIter; i++)
+ {
+ if (i != 0 && i < 4)
+ {
+ allLambda = allCoLambda;
+ }
+
+ limiter
+ (
+ allLambda,
+ rho,
+ psi,
+ phiBD,
+ phiCorr,
+ Sp,
+ Su,
+ psiMax,
+ psiMin,
+ nLimiterIter
+ );
+
+ solve
+ (
+ psiConvectionDiffusion + fvc::div(lambda*phiCorr),
+ MULEScontrols
+ );
+
+ scalar maxPsiM1 = gMax(psi.internalField()) - 1.0;
+ scalar minPsi = gMin(psi.internalField());
+
+ scalar unboundedness = max(max(maxPsiM1, 0.0), -min(minPsi, 0.0));
+
+ if (unboundedness < maxUnboundedness)
+ {
+ break;
+ }
+ else
+ {
+ Info<< "MULES: max(" << psi.name() << " - 1) = " << maxPsiM1
+ << " min(" << psi.name() << ") = " << minPsi << endl;
+
+ phiBD = psiConvectionDiffusion.flux();
+
+ /*
+ word gammaScheme("div(phi,gamma)");
+ word gammarScheme("div(phirb,gamma)");
+
+ const surfaceScalarField& phir =
+ mesh.lookupObject<surfaceScalarField>("phir");
+
+ phiCorr =
+ fvc::flux
+ (
+ phi,
+ psi,
+ gammaScheme
+ )
+ + fvc::flux
+ (
+ -fvc::flux(-phir, scalar(1) - psi, gammarScheme),
+ psi,
+ gammarScheme
+ )
+ - phiBD;
+ */
+ }
+ }
+
+ phiPsi = psiConvectionDiffusion.flux() + lambda*phiCorr;
+}
+
+
+// ************************************************************************* //
diff --git a/src/finiteVolume/fvMatrices/solvers/MULES/MULES.C b/src/finiteVolume/fvMatrices/solvers/MULES/MULES.C
index 6f5b58c3e5d..63d37106d6f 100644
--- a/src/finiteVolume/fvMatrices/solvers/MULES/MULES.C
+++ b/src/finiteVolume/fvMatrices/solvers/MULES/MULES.C
@@ -24,17 +24,6 @@ License
\*---------------------------------------------------------------------------*/
#include "MULES.H"
-#include "upwind.H"
-#include "uncorrectedSnGrad.H"
-#include "gaussConvectionScheme.H"
-#include "gaussLaplacianScheme.H"
-#include "uncorrectedSnGrad.H"
-#include "surfaceInterpolate.H"
-#include "fvcSurfaceIntegrate.H"
-#include "slicedSurfaceFields.H"
-#include "syncTools.H"
-
-#include "fvm.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -59,27 +48,6 @@ void Foam::MULES::explicitSolve
}
-void Foam::MULES::implicitSolve
-(
- volScalarField& psi,
- const surfaceScalarField& phi,
- surfaceScalarField& phiPsi,
- const scalar psiMax,
- const scalar psiMin
-)
-{
- implicitSolve
- (
- geometricOneField(),
- psi,
- phi,
- phiPsi,
- zeroField(), zeroField(),
- psiMax, psiMin
- );
-}
-
-
void Foam::MULES::limitSum(UPtrList<scalarField>& phiPsiCorrs)
{
forAll(phiPsiCorrs[0], facei)
diff --git a/src/finiteVolume/fvMatrices/solvers/MULES/MULES.H b/src/finiteVolume/fvMatrices/solvers/MULES/MULES.H
index 283c308664e..2cdd99eef74 100644
--- a/src/finiteVolume/fvMatrices/solvers/MULES/MULES.H
+++ b/src/finiteVolume/fvMatrices/solvers/MULES/MULES.H
@@ -36,17 +36,20 @@ Description
SourceFiles
MULES.C
+ MULESTemplates.C
\*---------------------------------------------------------------------------*/
#ifndef MULES_H
#define MULES_H
-#include "volFields.H"
+#include "volFieldsFwd.H"
#include "surfaceFieldsFwd.H"
#include "primitiveFieldsFwd.H"
#include "geometricOneField.H"
#include "zero.H"
+#include "zeroField.H"
+#include "UPtrList.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -89,28 +92,6 @@ void explicitSolve
const scalar psiMin
);
-template<class RhoType, class SpType, class SuType>
-void implicitSolve
-(
- const RhoType& rho,
- volScalarField& gamma,
- const surfaceScalarField& phi,
- surfaceScalarField& phiCorr,
- const SpType& Sp,
- const SuType& Su,
- const scalar psiMax,
- const scalar psiMin
-);
-
-void implicitSolve
-(
- volScalarField& gamma,
- const surfaceScalarField& phi,
- surfaceScalarField& phiCorr,
- const scalar psiMax,
- const scalar psiMin
-);
-
template<class RhoType, class SpType, class SuType>
void limiter
(
diff --git a/src/finiteVolume/fvMatrices/solvers/MULES/MULESTemplates.C b/src/finiteVolume/fvMatrices/solvers/MULES/MULESTemplates.C
index d45cc25f298..a3cb9f17e9d 100644
--- a/src/finiteVolume/fvMatrices/solvers/MULES/MULESTemplates.C
+++ b/src/finiteVolume/fvMatrices/solvers/MULES/MULESTemplates.C
@@ -25,18 +25,11 @@ License
#include "MULES.H"
#include "upwind.H"
-#include "uncorrectedSnGrad.H"
-#include "gaussConvectionScheme.H"
-#include "gaussLaplacianScheme.H"
-#include "uncorrectedSnGrad.H"
-#include "surfaceInterpolate.H"
#include "fvcSurfaceIntegrate.H"
#include "slicedSurfaceFields.H"
#include "wedgeFvPatch.H"
#include "syncTools.H"
-#include "fvm.H"
-
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
template<class RhoType, class SpType, class SuType>
@@ -103,215 +96,6 @@ void Foam::MULES::explicitSolve
}
-namespace Foam
-{
-namespace MULES
-{
- template<class RhoType>
- inline tmp<surfaceScalarField> interpolate(const RhoType& rho)
- {
- notImplemented
- (
- "tmp<surfaceScalarField> interpolate(const RhoType& rho)"
- );
- return tmp<surfaceScalarField>(NULL);
- }
-
- template<>
- inline tmp<surfaceScalarField> interpolate(const volScalarField& rho)
- {
- return fvc::interpolate(rho);
- }
-}
-}
-
-
-template<class RhoType, class SpType, class SuType>
-void Foam::MULES::implicitSolve
-(
- const RhoType& rho,
- volScalarField& psi,
- const surfaceScalarField& phi,
- surfaceScalarField& phiPsi,
- const SpType& Sp,
- const SuType& Su,
- const scalar psiMax,
- const scalar psiMin
-)
-{
- const fvMesh& mesh = psi.mesh();
-
- const dictionary& MULEScontrols = mesh.solverDict(psi.name());
-
- label maxIter
- (
- readLabel(MULEScontrols.lookup("maxIter"))
- );
-
- label nLimiterIter
- (
- readLabel(MULEScontrols.lookup("nLimiterIter"))
- );
-
- scalar maxUnboundedness
- (
- readScalar(MULEScontrols.lookup("maxUnboundedness"))
- );
-
- scalar CoCoeff
- (
- readScalar(MULEScontrols.lookup("CoCoeff"))
- );
-
- scalarField allCoLambda(mesh.nFaces());
-
- {
- slicedSurfaceScalarField CoLambda
- (
- IOobject
- (
- "CoLambda",
- mesh.time().timeName(),
- mesh,
- IOobject::NO_READ,
- IOobject::NO_WRITE,
- false
- ),
- mesh,
- dimless,
- allCoLambda,
- false // Use slices for the couples
- );
-
- if (phi.dimensions() == dimDensity*dimVelocity*dimArea)
- {
- tmp<surfaceScalarField> Cof =
- mesh.time().deltaT()*mesh.surfaceInterpolation::deltaCoeffs()
- *mag(phi/interpolate(rho))/mesh.magSf();
-
- CoLambda == 1.0/max(CoCoeff*Cof, scalar(1));
- }
- else
- {
- tmp<surfaceScalarField> Cof =
- mesh.time().deltaT()*mesh.surfaceInterpolation::deltaCoeffs()
- *mag(phi)/mesh.magSf();
-
- CoLambda == 1.0/max(CoCoeff*Cof, scalar(1));
- }
- }
-
- scalarField allLambda(allCoLambda);
- //scalarField allLambda(mesh.nFaces(), 1.0);
-
- slicedSurfaceScalarField lambda
- (
- IOobject
- (
- "lambda",
- mesh.time().timeName(),
- mesh,
- IOobject::NO_READ,
- IOobject::NO_WRITE,
- false
- ),
- mesh,
- dimless,
- allLambda,
- false // Use slices for the couples
- );
-
- linear<scalar> CDs(mesh);
- upwind<scalar> UDs(mesh, phi);
- //fv::uncorrectedSnGrad<scalar> snGrads(mesh);
-
- fvScalarMatrix psiConvectionDiffusion
- (
- fvm::ddt(rho, psi)
- + fv::gaussConvectionScheme<scalar>(mesh, phi, UDs).fvmDiv(phi, psi)
- //- fv::gaussLaplacianScheme<scalar, scalar>(mesh, CDs, snGrads)
- //.fvmLaplacian(Dpsif, psi)
- - fvm::Sp(Sp, psi)
- - Su
- );
-
- surfaceScalarField phiBD(psiConvectionDiffusion.flux());
-
- surfaceScalarField& phiCorr = phiPsi;
- phiCorr -= phiBD;
-
- for (label i=0; i<maxIter; i++)
- {
- if (i != 0 && i < 4)
- {
- allLambda = allCoLambda;
- }
-
- limiter
- (
- allLambda,
- rho,
- psi,
- phiBD,
- phiCorr,
- Sp,
- Su,
- psiMax,
- psiMin,
- nLimiterIter
- );
-
- solve
- (
- psiConvectionDiffusion + fvc::div(lambda*phiCorr),
- MULEScontrols
- );
-
- scalar maxPsiM1 = gMax(psi.internalField()) - 1.0;
- scalar minPsi = gMin(psi.internalField());
-
- scalar unboundedness = max(max(maxPsiM1, 0.0), -min(minPsi, 0.0));
-
- if (unboundedness < maxUnboundedness)
- {
- break;
- }
- else
- {
- Info<< "MULES: max(" << psi.name() << " - 1) = " << maxPsiM1
- << " min(" << psi.name() << ") = " << minPsi << endl;
-
- phiBD = psiConvectionDiffusion.flux();
-
- /*
- word gammaScheme("div(phi,gamma)");
- word gammarScheme("div(phirb,gamma)");
-
- const surfaceScalarField& phir =
- mesh.lookupObject<surfaceScalarField>("phir");
-
- phiCorr =
- fvc::flux
- (
- phi,
- psi,
- gammaScheme
- )
- + fvc::flux
- (
- -fvc::flux(-phir, scalar(1) - psi, gammarScheme),
- psi,
- gammarScheme
- )
- - phiBD;
- */
- }
- }
-
- phiPsi = psiConvectionDiffusion.flux() + lambda*phiCorr;
-}
-
-
template<class RhoType, class SpType, class SuType>
void Foam::MULES::limiter
(
--
GitLab