From 22e02bcced8b3b71d188279a7708948241acd06e Mon Sep 17 00:00:00 2001 From: Hrvoje Jasak Date: Fri, 20 Oct 2017 12:52:34 +0100 Subject: [PATCH] Convection-diffusion steady intertial ddt scheme --- .../steadyInertialDdtScheme.C | 123 +++++++++++++++++- .../steadyInertialDdtScheme.H | 30 ++++- 2 files changed, 144 insertions(+), 9 deletions(-) diff --git a/src/finiteVolume/finiteVolume/ddtSchemes/steadyInertialDdtScheme/steadyInertialDdtScheme.C b/src/finiteVolume/finiteVolume/ddtSchemes/steadyInertialDdtScheme/steadyInertialDdtScheme.C index ed798b4ce..0bfc58a41 100644 --- a/src/finiteVolume/finiteVolume/ddtSchemes/steadyInertialDdtScheme/steadyInertialDdtScheme.C +++ b/src/finiteVolume/finiteVolume/ddtSchemes/steadyInertialDdtScheme/steadyInertialDdtScheme.C @@ -26,6 +26,7 @@ License #include "steadyInertialDdtScheme.H" #include "surfaceInterpolate.H" #include "fvMatrices.H" +#include "autoPtr.H" // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // @@ -42,8 +43,52 @@ namespace fv template tmp steadyInertialDdtScheme::CorDeltaT() const { + const objectRegistry& registry = this->mesh(); + + autoPtr cofrDeltaTPtr; + + if (registry.foundObject(phiName_)) + { + cofrDeltaTPtr.reset(convectionCofrDeltaT().ptr()); + Info<< "Convection min ddt: " + << gMin(1/cofrDeltaTPtr().internalField()) + << endl; + } + + if + ( + registry.foundObject(nuName_) + || registry.foundObject(nuName_) + ) + { + if (cofrDeltaTPtr.valid()) + { + cofrDeltaTPtr() = + Foam::max(cofrDeltaTPtr(), diffusionCofrDeltaT()); + Info<< "Combined min ddt: " + << gMin(1/cofrDeltaTPtr().internalField()) + << endl; + } + else + { + cofrDeltaTPtr.reset(diffusionCofrDeltaT().ptr()); + Info<< "Diffusion min ddt: " + << gMin(1/cofrDeltaTPtr().internalField()) + << endl; + } + } + + if (cofrDeltaTPtr.empty()) + { + FatalErrorIn + ( + "steaddyInertialDdtScheme::CorDeltaT() const" + ) << "Cannot find phi or nu: " << phiName_ << " " << nuName_ + << abort(FatalError); + } + // Collect face delta t and pick the smallest for the cell - surfaceScalarField cofrDeltaT = CofrDeltaT(); + surfaceScalarField& cofrDeltaT = cofrDeltaTPtr(); tmp tcorDeltaT ( @@ -60,7 +105,6 @@ tmp steadyInertialDdtScheme::CorDeltaT() const zeroGradientFvPatchScalarField::typeName ) ); - volScalarField& corDeltaT = tcorDeltaT(); const unallocLabelList& owner = mesh().owner(); @@ -100,7 +144,8 @@ tmp steadyInertialDdtScheme::CorDeltaT() const template -tmp steadyInertialDdtScheme::CofrDeltaT() const +tmp +steadyInertialDdtScheme::convectionCofrDeltaT() const { const objectRegistry& registry = this->mesh(); @@ -143,8 +188,76 @@ tmp steadyInertialDdtScheme::CofrDeltaT() const } else { - FatalErrorIn("steaddyInertialDdtScheme::CofrDeltaT() const") - << "Incorrect dimensions of phi: " << phi.dimensions() + FatalErrorIn + ( + "steaddyInertialDdtScheme::convectionCofrDeltaT() const" + ) << "Incorrect dimensions of phi: " << phi.dimensions() + << abort(FatalError); + + return tmp(NULL); + } +} + + +template +tmp +steadyInertialDdtScheme::diffusionCofrDeltaT() const +{ + const objectRegistry& registry = this->mesh(); + + if (registry.foundObject(nuName_)) + { + const volScalarField& nu = + registry.lookupObject(nuName_); + + return diffusionCofrDeltaT(fvc::interpolate(nu)()); + } + else if (registry.foundObject(nuName_)) + { + const surfaceScalarField& nuf = + registry.lookupObject(nuName_); + + return diffusionCofrDeltaT(nuf); + } + else + { + FatalErrorIn + ( + "steaddyInertialDdtScheme::diffusionCofrDeltaT() const" + ) << "Cannot find nu" + << abort(FatalError); + + return tmp(NULL); + } +} + +template +tmp +steadyInertialDdtScheme::diffusionCofrDeltaT +( + const surfaceScalarField& nuf +) const +{ + const objectRegistry& registry = this->mesh(); + + if (nuf.dimensions() == dimensionSet(0, 2, -1, 0, 0)) + { + return nuf*sqr(mesh().surfaceInterpolation::deltaCoeffs())/maxCo_; + } + else if (nuf.dimensions() == dimensionSet(1, -1, -1, 0, 0)) + { + const volScalarField& rho = + registry.lookupObject(rhoName_); + + return nuf*sqr(mesh().surfaceInterpolation::deltaCoeffs())/ + (fvc::interpolate(rho)*maxCo_); + } + else + { + FatalErrorIn + ( + "steaddyInertialDdtScheme::diffusionCofrDeltaT() const" + ) << "Incorrect dimensions of nu: " << nuf.dimensions() << abort(FatalError); return tmp(NULL); diff --git a/src/finiteVolume/finiteVolume/ddtSchemes/steadyInertialDdtScheme/steadyInertialDdtScheme.H b/src/finiteVolume/finiteVolume/ddtSchemes/steadyInertialDdtScheme/steadyInertialDdtScheme.H index 3e3c47730..87599b48a 100644 --- a/src/finiteVolume/finiteVolume/ddtSchemes/steadyInertialDdtScheme/steadyInertialDdtScheme.H +++ b/src/finiteVolume/finiteVolume/ddtSchemes/steadyInertialDdtScheme/steadyInertialDdtScheme.H @@ -25,7 +25,8 @@ Class Foam::fv::steadyInertialDdtScheme Description - Stabilised local time-step first-order Euler implicit/explicit ddt. + Stabilised local time-step first-order Euler implicit/explicit ddt + for intertially relaxed steady-state simulations. The time-step is adjusted locally so that an advective equations remains diagonally dominant. @@ -35,6 +36,9 @@ Description See also CoEulerDdtScheme. +Author + Hrvoje Jasak, Wikki Ltd. All rights reserved. + SourceFiles steadyInertialDdtScheme.C @@ -66,8 +70,13 @@ class steadyInertialDdtScheme { // Private Data - //- Name of the flux field used to calculate the local time-step - word phiName_; + //- Name of the flux field used to calculate the local + // convective time-step + word phiName_; + + //- Name of the diffusivity field used to calculate the local + // diffusive time-step + word nuName_; //- Name of the density field used to obtain the volumetric flux // from the mass flux if required @@ -89,7 +98,19 @@ class steadyInertialDdtScheme tmp CorDeltaT() const; //- Return the reciprocal of the face-Courant-number limited time-step - tmp CofrDeltaT() const; + // for convective transport + tmp convectionCofrDeltaT() const; + + //- Return the reciprocal of the face-Courant-number limited time-step + // for diffusive transport + tmp diffusionCofrDeltaT() const; + + //- Return the reciprocal of the face-Courant-number limited time-step + // for diffusive transport ginev face diffusivity + tmp diffusionCofrDeltaT + ( + const surfaceScalarField& nuf + ) const; public: @@ -105,6 +126,7 @@ public: : ddtScheme(mesh, is), phiName_(is), + nuName_(is), rhoName_(is), maxCo_(readScalar(is)) {}