From 3e2b64c08d1a88e38791c391ae73201c648eaade Mon Sep 17 00:00:00 2001 From: Henry Date: Wed, 29 Apr 2015 14:37:41 +0100 Subject: [PATCH] MULES: nLimiterIter and smoothLimiter are now user-input via the corresponding fvSolution sub-dict nLimiterIter: Number of iterations during limiter construction 3 (default) is sufficient for 3D simulations with a Courant number 0.5 or so For larger Courant numbers larger values may be needed but this is only relevant for IMULES and CMULES smoothLimiter: Coefficient to smooth the limiter to avoid "diamond" staggering patters seen in regions of low particle phase-fraction in fluidised-bed simulations. The default is 0 as it is not needed for all simulations. A value of 0.1 is appropriate for fluidised-bed simulations. The useful range is 0 -> 0.5. Values larger than 0.5 may cause excessive smearing of the solution. --- .../multiphaseMixtureThermo.C | 1 - .../interFoam/interMixingFoam/alphaEqns.H | 6 +- .../multiphaseSystem/multiphaseSystem.C | 1 - .../fvMatrices/solvers/MULES/CMULES.H | 6 +- .../solvers/MULES/CMULESTemplates.C | 59 ++++++++--------- .../solvers/MULES/IMULESTemplates.C | 8 +-- .../fvMatrices/solvers/MULES/MULES.H | 4 +- .../fvMatrices/solvers/MULES/MULESTemplates.C | 66 +++++++++++++++---- .../bubbleColumn/system/fvSolution | 2 +- .../damBreak4phase/system/fvSolution | 2 +- .../damBreak4phaseFine/system/fvSolution | 2 +- .../mixerVessel2D/system/fvSolution | 2 +- .../mixerVesselAMI2D/system/fvSolution | 2 +- .../laminar/damBreak4phase/system/fvSolution | 2 +- .../damBreak4phaseFine/system/fvSolution | 2 +- .../laminar/mixerVessel2D/system/fvSolution | 2 +- .../RAS/fluidisedBed/system/fvSolution | 4 +- .../laminar/fluidisedBed/system/fvSolution | 4 +- 18 files changed, 103 insertions(+), 72 deletions(-) diff --git a/applications/solvers/multiphase/compressibleMultiphaseInterFoam/multiphaseMixtureThermo/multiphaseMixtureThermo.C b/applications/solvers/multiphase/compressibleMultiphaseInterFoam/multiphaseMixtureThermo/multiphaseMixtureThermo.C index 8b0ea995e..59cd77aee 100644 --- a/applications/solvers/multiphase/compressibleMultiphaseInterFoam/multiphaseMixtureThermo/multiphaseMixtureThermo.C +++ b/applications/solvers/multiphase/compressibleMultiphaseInterFoam/multiphaseMixtureThermo/multiphaseMixtureThermo.C @@ -996,7 +996,6 @@ void Foam::multiphaseMixtureThermo::solveAlphas zeroField(), 1, 0, - 3, true ); diff --git a/applications/solvers/multiphase/interFoam/interMixingFoam/alphaEqns.H b/applications/solvers/multiphase/interFoam/interMixingFoam/alphaEqns.H index 6bdb9b66f..e0fe8ff6a 100644 --- a/applications/solvers/multiphase/interFoam/interMixingFoam/alphaEqns.H +++ b/applications/solvers/multiphase/interFoam/interMixingFoam/alphaEqns.H @@ -83,8 +83,7 @@ zeroField(), zeroField(), 1, - 0, - 3 + 0 ); // Create the complete flux for alpha2 @@ -125,8 +124,7 @@ zeroField(), zeroField(), 1, - 0, - 3 + 0 ); // Construct the limited fluxes diff --git a/applications/solvers/multiphase/multiphaseEulerFoam/multiphaseSystem/multiphaseSystem.C b/applications/solvers/multiphase/multiphaseEulerFoam/multiphaseSystem/multiphaseSystem.C index b514449b9..f600e096b 100644 --- a/applications/solvers/multiphase/multiphaseEulerFoam/multiphaseSystem/multiphaseSystem.C +++ b/applications/solvers/multiphase/multiphaseEulerFoam/multiphaseSystem/multiphaseSystem.C @@ -158,7 +158,6 @@ void Foam::multiphaseSystem::solveAlphas() zeroField(), 1, 0, - 3, true ); diff --git a/src/finiteVolume/fvMatrices/solvers/MULES/CMULES.H b/src/finiteVolume/fvMatrices/solvers/MULES/CMULES.H index 78a1d9f9d..d21200ddc 100644 --- a/src/finiteVolume/fvMatrices/solvers/MULES/CMULES.H +++ b/src/finiteVolume/fvMatrices/solvers/MULES/CMULES.H @@ -129,8 +129,7 @@ void limiterCorr const SpType& Sp, const SuType& Su, const scalar psiMax, - const scalar psiMin, - const label nLimiterIter + const scalar psiMin ); template @@ -144,8 +143,7 @@ void limitCorr const SpType& Sp, const SuType& Su, const scalar psiMax, - const scalar psiMin, - const label nLimiterIter + const scalar psiMin ); diff --git a/src/finiteVolume/fvMatrices/solvers/MULES/CMULESTemplates.C b/src/finiteVolume/fvMatrices/solvers/MULES/CMULESTemplates.C index 090c096d6..a3901be04 100644 --- a/src/finiteVolume/fvMatrices/solvers/MULES/CMULESTemplates.C +++ b/src/finiteVolume/fvMatrices/solvers/MULES/CMULESTemplates.C @@ -89,13 +89,6 @@ void Foam::MULES::correct const fvMesh& mesh = psi.mesh(); const scalar rDeltaT = 1.0/mesh.time().deltaTValue(); - const dictionary& MULEScontrols = mesh.solverDict(psi.name()); - - label nLimiterIter - ( - readLabel(MULEScontrols.lookup("nLimiterIter")) - ); - limitCorr ( rDeltaT, @@ -103,9 +96,10 @@ void Foam::MULES::correct psi, phi, phiCorr, - Sp, Su, - psiMax, psiMin, - nLimiterIter + Sp, + Su, + psiMax, + psiMin ); correct(rDeltaT, rho, psi, phi, phiCorr, Sp, Su); @@ -130,13 +124,6 @@ void Foam::MULES::LTScorrect const volScalarField& rDeltaT = mesh.objectRegistry::lookupObject("rSubDeltaT"); - const dictionary& MULEScontrols = mesh.solverDict(psi.name()); - - label nLimiterIter - ( - readLabel(MULEScontrols.lookup("nLimiterIter")) - ); - limitCorr ( rDeltaT, @@ -144,9 +131,10 @@ void Foam::MULES::LTScorrect psi, phi, phiCorr, - Sp, Su, - psiMax, psiMin, - nLimiterIter + Sp, + Su, + psiMax, + psiMin ); correct(rDeltaT, rho, psi, phi, phiCorr, Sp, Su); } @@ -164,8 +152,7 @@ void Foam::MULES::limiterCorr const SpType& Sp, const SuType& Su, const scalar psiMax, - const scalar psiMin, - const label nLimiterIter + const scalar psiMin ) { const scalarField& psiIf = psi; @@ -175,9 +162,19 @@ void Foam::MULES::limiterCorr const dictionary& MULEScontrols = mesh.solverDict(psi.name()); + label nLimiterIter + ( + readLabel(MULEScontrols.lookup("nLimiterIter")) + ); + + scalar smoothLimiter + ( + MULEScontrols.lookupOrDefault("smoothLimiter", 0) + ); + scalar extremaCoeff ( - MULEScontrols.lookupOrDefault("extremaCoeff", 0.0) + MULEScontrols.lookupOrDefault("extremaCoeff", 0) ); const labelUList& owner = mesh.owner(); @@ -294,9 +291,13 @@ void Foam::MULES::limiterCorr psiMaxn = min(psiMaxn + extremaCoeff*(psiMax - psiMin), psiMax); psiMinn = max(psiMinn - extremaCoeff*(psiMax - psiMin), psiMin); - // scalar smooth = 0.5; - // psiMaxn = min((1.0 - smooth)*psiIf + smooth*psiMaxn, psiMax); - // psiMinn = max((1.0 - smooth)*psiIf + smooth*psiMinn, psiMin); + if (smoothLimiter > SMALL) + { + psiMaxn = + min(smoothLimiter*psiIf + (1.0 - smoothLimiter)*psiMaxn, psiMax); + psiMinn = + max(smoothLimiter*psiIf + (1.0 - smoothLimiter)*psiMinn, psiMin); + } psiMaxn = V @@ -483,8 +484,7 @@ void Foam::MULES::limitCorr const SpType& Sp, const SuType& Su, const scalar psiMax, - const scalar psiMin, - const label nLimiterIter + const scalar psiMin ) { const fvMesh& mesh = psi.mesh(); @@ -519,8 +519,7 @@ void Foam::MULES::limitCorr Sp, Su, psiMax, - psiMin, - nLimiterIter + psiMin ); phiCorr *= lambda; diff --git a/src/finiteVolume/fvMatrices/solvers/MULES/IMULESTemplates.C b/src/finiteVolume/fvMatrices/solvers/MULES/IMULESTemplates.C index 542181e90..50e7eb921 100644 --- a/src/finiteVolume/fvMatrices/solvers/MULES/IMULESTemplates.C +++ b/src/finiteVolume/fvMatrices/solvers/MULES/IMULESTemplates.C @@ -77,11 +77,6 @@ void Foam::MULES::implicitSolve readLabel(MULEScontrols.lookup("maxIter")) ); - label nLimiterIter - ( - readLabel(MULEScontrols.lookup("nLimiterIter")) - ); - scalar maxUnboundedness ( readScalar(MULEScontrols.lookup("maxUnboundedness")) @@ -187,8 +182,7 @@ void Foam::MULES::implicitSolve Sp, Su, psiMax, - psiMin, - nLimiterIter + psiMin ); solve diff --git a/src/finiteVolume/fvMatrices/solvers/MULES/MULES.H b/src/finiteVolume/fvMatrices/solvers/MULES/MULES.H index f0525ddc9..65817dc1f 100644 --- a/src/finiteVolume/fvMatrices/solvers/MULES/MULES.H +++ b/src/finiteVolume/fvMatrices/solvers/MULES/MULES.H @@ -138,8 +138,7 @@ void limiter const SpType& Sp, const SuType& Su, const scalar psiMax, - const scalar psiMin, - const label nLimiterIter + const scalar psiMin ); template @@ -154,7 +153,6 @@ void limit const SuType& Su, const scalar psiMax, const scalar psiMin, - const label nLimiterIter, const bool returnCorr ); diff --git a/src/finiteVolume/fvMatrices/solvers/MULES/MULESTemplates.C b/src/finiteVolume/fvMatrices/solvers/MULES/MULESTemplates.C index a7a528281..84f9dba61 100644 --- a/src/finiteVolume/fvMatrices/solvers/MULES/MULESTemplates.C +++ b/src/finiteVolume/fvMatrices/solvers/MULES/MULESTemplates.C @@ -108,8 +108,23 @@ void Foam::MULES::explicitSolve { const fvMesh& mesh = psi.mesh(); const scalar rDeltaT = 1.0/mesh.time().deltaTValue(); + psi.correctBoundaryConditions(); - limit(rDeltaT, rho, psi, phi, phiPsi, Sp, Su, psiMax, psiMin, 3, false); + + limit + ( + rDeltaT, + rho, + psi, + phi, + phiPsi, + Sp, + Su, + psiMax, + psiMin, + false + ); + explicitSolve(rDeltaT, rho, psi, phiPsi, Sp, Su); } @@ -133,7 +148,21 @@ void Foam::MULES::explicitLTSSolve mesh.objectRegistry::lookupObject("rSubDeltaT"); psi.correctBoundaryConditions(); - limit(rDeltaT, rho, psi, phi, phiPsi, Sp, Su, psiMax, psiMin, 3, false); + + limit + ( + rDeltaT, + rho, + psi, + phi, + phiPsi, + Sp, + Su, + psiMax, + psiMin, + false + ); + explicitSolve(rDeltaT, rho, psi, phiPsi, Sp, Su); } @@ -150,17 +179,28 @@ void Foam::MULES::limiter const SpType& Sp, const SuType& Su, const scalar psiMax, - const scalar psiMin, - const label nLimiterIter + const scalar psiMin ) { const scalarField& psiIf = psi; const volScalarField::GeometricBoundaryField& psiBf = psi.boundaryField(); - const scalarField& psi0 = psi.oldTime(); - const fvMesh& mesh = psi.mesh(); + const dictionary& MULEScontrols = mesh.solverDict(psi.name()); + + label nLimiterIter + ( + MULEScontrols.lookupOrDefault