diff --git a/src/functionObjects/solvers/electricPotential/electricPotential.C b/src/functionObjects/solvers/electricPotential/electricPotential.C index ff931d8583..d484a82e55 100644 --- a/src/functionObjects/solvers/electricPotential/electricPotential.C +++ b/src/functionObjects/solvers/electricPotential/electricPotential.C @@ -208,6 +208,7 @@ Foam::functionObjects::electricPotential::electricPotential ) ), fvOptions_(mesh_), + tol_(1), nCorr_(1), writeDerivedFields_(false), electricField_(false) @@ -258,6 +259,7 @@ bool Foam::functionObjects::electricPotential::read(const dictionary& dict) dict.readIfPresent("sigma", sigma_); dict.readIfPresent("epsilonr", epsilonr_); dict.readIfPresent("nCorr", nCorr_); + dict.readIfPresent("tolerance", tol_); dict.readIfPresent("writeDerivedFields", writeDerivedFields_); dict.readIfPresent("electricField", electricField_); @@ -346,6 +348,10 @@ bool Foam::functionObjects::electricPotential::execute() { Log << type() << " execute: " << name() << endl; + // Convergence monitor parameters + bool converged = false; + label iter = 0; + tmp tsigma = this->sigma(); const auto& sigma = tsigma(); @@ -362,7 +368,9 @@ bool Foam::functionObjects::electricPotential::execute() fvOptions_.constrain(eVEqn); - eVEqn.solve(); + ++iter; + converged = (eVEqn.solve().initialResidual() < tol_); + if (converged) break; } if (electricField_) @@ -371,6 +379,14 @@ bool Foam::functionObjects::electricPotential::execute() E == -fvc::grad(eV); } + if (converged) + { + Log << type() << ": " << name() << ": " + << eV.name() << " is converged." << nl + << tab << "initial-residual tolerance: " << tol_ << nl + << tab << "outer iteration: " << iter << nl; + } + Log << endl; return true; diff --git a/src/functionObjects/solvers/electricPotential/electricPotential.H b/src/functionObjects/solvers/electricPotential/electricPotential.H index af1b0e47fe..d41a2b37df 100644 --- a/src/functionObjects/solvers/electricPotential/electricPotential.H +++ b/src/functionObjects/solvers/electricPotential/electricPotential.H @@ -124,6 +124,7 @@ Usage electricField ; E ; fvOptions ; + tolerance ; // Inherited entries ... @@ -143,6 +144,7 @@ Usage electricField | Flag to calculate electric field | bool | no | false E | Name of electric field | word | no | electricPotential:E fvOptions | List of finite-volume options | dict | no | - + tolerance | Outer-loop initial-residual tolerance | scalar | no | 1 \endtable The inherited entries are elaborated in: @@ -218,6 +220,9 @@ class electricPotential //- Run-time selectable finite volume options fv::optionList fvOptions_; + //- Outer-loop initial-residual tolerance + scalar tol_; + //- Number of corrector iterations int nCorr_; diff --git a/src/functionObjects/solvers/energyTransport/energyTransport.C b/src/functionObjects/solvers/energyTransport/energyTransport.C index 2a0ac81eba..730a320f8e 100644 --- a/src/functionObjects/solvers/energyTransport/energyTransport.C +++ b/src/functionObjects/solvers/energyTransport/energyTransport.C @@ -26,11 +26,6 @@ License \*---------------------------------------------------------------------------*/ #include "energyTransport.H" -#include "surfaceFields.H" -#include "fvmDdt.H" -#include "fvmDiv.H" -#include "fvmLaplacian.H" -#include "fvmSup.H" #include "turbulentTransportModel.H" #include "turbulentFluidThermoModel.H" #include "addToRunTimeSelectionTable.H" @@ -188,23 +183,6 @@ Foam::functionObjects::energyTransport::energyTransport ) : fvMeshFunctionObject(name, runTime, dict), - fieldName_(dict.getOrDefault("field", "T")), - phiName_(dict.getOrDefault("phi", "phi")), - rhoName_(dict.getOrDefault("rho", "rho")), - nCorr_(0), - schemesField_("unknown-schemesField"), - fvOptions_(mesh_), - multiphaseThermo_(dict.subOrEmptyDict("phaseThermos")), - Cp_("Cp", dimEnergy/dimMass/dimTemperature, 0, dict), - kappa_ - ( - "kappa", - dimEnergy/dimTime/dimLength/dimTemperature, - 0, - dict - ), - rho_("rhoInf", dimDensity, 0, dict), - Prt_("Prt", dimless, 1, dict), rhoCp_ ( IOobject @@ -218,7 +196,25 @@ Foam::functionObjects::energyTransport::energyTransport ), mesh_, dimensionedScalar(dimEnergy/dimTemperature/dimVolume, Zero) - ) + ), + fvOptions_(mesh_), + multiphaseThermo_(dict.subOrEmptyDict("phaseThermos")), + Cp_("Cp", dimEnergy/dimMass/dimTemperature, 0, dict), + kappa_ + ( + "kappa", + dimEnergy/dimTime/dimLength/dimTemperature, + 0, + dict + ), + rho_("rhoInf", dimDensity, 0, dict), + Prt_("Prt", dimless, 1, dict), + fieldName_(dict.getOrDefault("field", "T")), + schemesField_("unknown-schemesField"), + phiName_(dict.getOrDefault("phi", "phi")), + rhoName_(dict.getOrDefault("rho", "rho")), + tol_(1), + nCorr_(0) { read(dict); @@ -305,17 +301,14 @@ Foam::functionObjects::energyTransport::energyTransport } -// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * // - -Foam::functionObjects::energyTransport::~energyTransport() -{} - - // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // bool Foam::functionObjects::energyTransport::read(const dictionary& dict) { - fvMeshFunctionObject::read(dict); + if (!fvMeshFunctionObject::read(dict)) + { + return false; + } dict.readIfPresent("phi", phiName_); dict.readIfPresent("rho", rhoName_); @@ -323,6 +316,7 @@ bool Foam::functionObjects::energyTransport::read(const dictionary& dict) schemesField_ = dict.getOrDefault("schemesField", fieldName_); dict.readIfPresent("nCorr", nCorr_); + dict.readIfPresent("tolerance", tol_); if (dict.found("fvOptions")) { @@ -355,12 +349,16 @@ bool Foam::functionObjects::energyTransport::execute() scalar relaxCoeff = 0; mesh_.relaxEquation(schemesField_, relaxCoeff); + // Convergence monitor parameters + bool converged = false; + int iter = 0; + if (phi.dimensions() == dimMass/dimTime) { rhoCp_ = rho()*Cp(); const surfaceScalarField rhoCpPhi(fvc::interpolate(Cp())*phi); - for (label i = 0; i <= nCorr_; i++) + for (int i = 0; i <= nCorr_; ++i) { fvScalarMatrix sEqn ( @@ -376,7 +374,9 @@ bool Foam::functionObjects::energyTransport::execute() fvOptions_.constrain(sEqn); - sEqn.solve(schemesField_); + ++iter; + converged = (sEqn.solve(schemesField_).initialResidual() < tol_); + if (converged) break; } } else if (phi.dimensions() == dimVolume/dimTime) @@ -393,7 +393,7 @@ bool Foam::functionObjects::energyTransport::execute() rhoCp ); - for (label i = 0; i <= nCorr_; i++) + for (int i = 0; i <= nCorr_; ++i) { fvScalarMatrix sEqn ( @@ -408,7 +408,9 @@ bool Foam::functionObjects::energyTransport::execute() fvOptions_.constrain(sEqn); - sEqn.solve(schemesField_); + ++iter; + converged = (sEqn.solve(schemesField_).initialResidual() < tol_); + if (converged) break; } } else @@ -419,6 +421,14 @@ bool Foam::functionObjects::energyTransport::execute() << dimVolume/dimTime << exit(FatalError); } + if (converged) + { + Log << type() << ": " << name() << ": " + << s.name() << " is converged." << nl + << tab << "initial-residual tolerance: " << tol_ << nl + << tab << "outer iteration: " << iter << nl; + } + Log << endl; return true; diff --git a/src/functionObjects/solvers/energyTransport/energyTransport.H b/src/functionObjects/solvers/energyTransport/energyTransport.H index aa3b9929f0..04512f81f3 100644 --- a/src/functionObjects/solvers/energyTransport/energyTransport.H +++ b/src/functionObjects/solvers/energyTransport/energyTransport.H @@ -30,151 +30,180 @@ Group grpSolversFunctionObjects Description - Evolves a simplified energy transport equation for incompressible flows. - It takes into account the inertia, conduction and convection terms plus - a source. + Computes the simplified energy transport equation in single-phase or + two-phase flow, considering incompressible cases: - - The field name must be temperature and its BC's specified in the time - directory. - - The turbulence model should be incompressible - - In order to use in a incompressible multi phase a list of thermal - properties are needed. See below + \f[ + \frac{\partial \rho \, C_p \, T}{\partial t} + + \nabla \cdot \left(\rho \, C_p \, \phi \, T \right) + - \nabla \cdot \left(\rho \, C_p \, \phi \right) \, T + - \nabla \cdot \left(\kappa_{eff} \, \nabla T \right) + = S_T + \f] + where: + \vartable + T | Scalar field + \rho | (Generic) Fluid density which is unity when not specified + C_p | Specific heat capacity at constant pressure + \phi | (Generic) Flux field + \kappa_{eff} | Effective thermal conductivity + S_T | Scalar field source term + \endvartable Usage - Example of function object specification to solve a energy transport - equation for a single phase flow plus a source term + Minimal example in \c system/controlDict.functions: \verbatim - functions + energyTransport1 { - energy - { - type energyTransport; - libs (energyTransportFunctionObjects); + // Mandatory entries + type energyTransport; + libs (solverFunctionObjects); - enabled true; - writeControl writeTime; - writeInterval 1; + // Optional entries + field ; + phi ; + rho ; + Cp ; + kappa ; + rhoInf ; + Prt ; + schemesField ; + tolerance ; + nCorr ; + fvOptions ; + phaseThermos ; - field T; - - // volumetric Flux - phi phi; - - // Thermal properties - Cp Cp [J/kg/K] 1e3; - kappa kappa [W/m/K] 0.0257; - rhoInf rho [kg/m^3] 1.2; - - write true; - - fvOptions - { - viscousDissipation - { - type viscousDissipation; - enabled true; - - viscousDissipationCoeffs - { - fields (T); - rhoInf $....rhoInf; - } - } - } - } + // Inherited entries + ... } \endverbatim - Example of function object specification to solve a energy transport - equation for a multiphase phase flow plus a source term - - equation: - \verbatim - functions - { - energy - { - type energyTransport; - libs (energyTransportFunctionObjects); - - enabled true; - writeControl writeTime; - writeInterval 1; - - field T; - - // rho field name - rho rho; - // mass flux for multiphase - phi rhoPhi; - - write true; - - // Thermal properties of the phases - phaseThermos - { - alpha.air - { - Cp 1e3; - kappa 0.0243; - } - alpha.mercury - { - Cp 140; - kappa 8.2; - } - alpha.oil - { - Cp 2e3; - kappa 0.2; - } - alpha.water - { - Cp 4e3; - kappa 0.6; - } - } - - - fvOptions - { - viscousDissipation - { - type viscousDissipation; - enabled true; - - viscousDissipationCoeffs - { - fields (T); - rho rho; //rho Field - } - } - } - } - } - \endverbatim - - Where the entries comprise: + where: \table - Property | Description | Required | Default value - type | Type name: energyTransport | yes | - field | Name of the scalar field | no | T - phi | Name of flux field | no | phi - rho | Name of density field | no | rho - nCorr | Number of correctors | no | 0 - schemesField | Name of field to specify schemes | no | field name - fvOptions | List of scalar sources | no | - Cp | Heat capacity for single phase | no | 0 - rhoInf | Density for single phase | no | 0 - kappa | Thermal conductivity for single phase | no | 0 - Prt | Turbulent Prandlt number | no | 1.0 - phaseThermos | Dictionary for multi-phase thermo |no | null - fvOptions | Opotional extra sources | no | null + Property | Description | Type | Reqd | Deflt + type | Type name: energyTransport | word | yes | - + libs | Library name: solverFunctionObjects | word | yes | - + field | Name of the passive-scalar field | word | no | s + phi | Name of flux field | word | no | phi + rho | Name of density field | word | no | rho + Cp | Specific heat capacity at constant pressure | scalar | no | 0 + kappa | Thermal conductivity | scalar | no | 0 + rhoInf | Fluid density | scalar | no | 0 + Prt | Turbulent Prandtl number | scalar | no | 1 + schemesField | Name of field to specify schemes | word | no | field + tolerance | Outer-loop initial-residual tolerance | scalar | no | 1 + nCorr | Number of outer-loop correctors | int | no | 0 + fvOptions | List of finite-volume options | dict | no | - + phaseThermos | Dictionary for multi-phase thermo | dict | no | null \endtable -See also - Foam::functionObjects::fvMeshFunctionObject + The inherited entries are elaborated in: + - \link fvMeshFunctionObject.H \endlink + - \link fvOption.H \endlink + + An example of function object specification to solve a energy transport + equation for a single phase flow plus a source term: + \verbatim + energyTransport1 + { + // Mandatory entries + type energyTransport; + libs (solverFunctionObjects); + + // Optional entries + field T; + phi phi; + Cp Cp [J/kg/K] 1e3; + kappa kappa [W/m/K] 0.0257; + rhoInf rho [kg/m^3] 1.2; + fvOptions + { + viscousDissipation + { + type viscousDissipation; + enabled true; + + viscousDissipationCoeffs + { + fields (T); + rhoInf $....rhoInf; + } + } + } + + // Inherited entries + enabled true; + writeControl writeTime; + writeInterval 1; + } + \endverbatim + + An example of function object specification to solve a energy transport + equation for a multiphase phase flow plus a source term: + \verbatim + energyTransport1 + { + // Mandatory entries + type energyTransport; + libs (solverFunctionObjects); + + // Optional entries + field T; + rho rho; + phi rhoPhi; + + // Thermal properties of the phases + phaseThermos + { + alpha.air + { + Cp 1e3; + kappa 0.0243; + } + alpha.mercury + { + Cp 140; + kappa 8.2; + } + alpha.oil + { + Cp 2e3; + kappa 0.2; + } + alpha.water + { + Cp 4e3; + kappa 0.6; + } + } + + fvOptions + { + viscousDissipation + { + type viscousDissipation; + enabled true; + + viscousDissipationCoeffs + { + fields (T); + rho rho; + } + } + } + + // Inherited entries + enabled true; + writeControl writeTime; + writeInterval 1; + } + \endverbatim + +Note + - The field name must be temperature and its boundary conditions + specified in the time directory. + - The turbulence model should be incompressible. SourceFiles energyTransport.C @@ -203,22 +232,10 @@ class energyTransport : public fvMeshFunctionObject { - // Private data + // Private Data - //- Name of the transport field. - word fieldName_; - - //- Name of flux field - word phiName_; - - //- Name of density field - word rhoName_; - - //- Number of corrector iterations (optional) - label nCorr_; - - //- Name of field whose schemes are used (optional) - word schemesField_; + //- Volumetric heat capacity field [J/m^3/K] + volScalarField rhoCp_; //- Run-time selectable finite volume options, e.g. sources, constraints fv::optionList fvOptions_; @@ -229,7 +246,7 @@ class energyTransport //- List of phase names wordList phaseNames_; - //- List of phase heat capacities + //- List of phase specific heat capacities at constant pressure PtrList Cps_; //- List of phase thermal diffusivity for temperature [J/m/s/K] @@ -238,7 +255,7 @@ class energyTransport //- Unallocated phase list UPtrList phases_; - //- Heat capacity for single phase flows + //- Specific heat capacity at constant pressure for single phase flows dimensionedScalar Cp_; //- Thermal diffusivity for temperature for single phase flows @@ -250,8 +267,23 @@ class energyTransport //- Turbulent Prandt number dimensionedScalar Prt_; - //- rhoCp - volScalarField rhoCp_; + //- Name of the transport field + word fieldName_; + + //- Name of field whose schemes are used + word schemesField_; + + //- Name of flux field + word phiName_; + + //- Name of density field + word rhoName_; + + //- Outer-loop initial-residual tolerance + scalar tol_; + + //- Number of corrector iterations + int nCorr_; // Private Member Functions @@ -262,21 +294,15 @@ class energyTransport //- Return the diffusivity field tmp kappaEff() const; - //- Return rho field + //- Return the density field, rho tmp rho() const; - //- Return Cp + //- Return the specific heat capacity at constant pressure field, Cp tmp Cp() const; - //- Return kappa + //- Return the thermal diffusivity field tmp kappa() const; - //- No copy construct - energyTransport(const energyTransport&) = delete; - - //- No copy assignment - void operator=(const energyTransport&) = delete; - public: @@ -296,7 +322,7 @@ public: //- Destructor - virtual ~energyTransport(); + virtual ~energyTransport() = default; // Member Functions diff --git a/src/functionObjects/solvers/scalarTransport/scalarTransport.C b/src/functionObjects/solvers/scalarTransport/scalarTransport.C index f418519024..99af527a7b 100644 --- a/src/functionObjects/solvers/scalarTransport/scalarTransport.C +++ b/src/functionObjects/solvers/scalarTransport/scalarTransport.C @@ -27,11 +27,6 @@ License \*---------------------------------------------------------------------------*/ #include "scalarTransport.H" -#include "surfaceFields.H" -#include "fvmDdt.H" -#include "fvmDiv.H" -#include "fvmLaplacian.H" -#include "fvmSup.H" #include "CMULES.H" #include "turbulentTransportModel.H" #include "turbulentFluidThermoModel.H" @@ -172,7 +167,9 @@ Foam::functionObjects::scalarTransport::scalarTransport ) : fvMeshFunctionObject(name, runTime, dict), + fvOptions_(mesh_), fieldName_(dict.getOrDefault("field", "s")), + schemesField_("unknown-schemesField"), phiName_(dict.getOrDefault("phi", "phi")), rhoName_(dict.getOrDefault("rho", "rho")), nutName_(dict.getOrDefault("nut", "none")), @@ -182,11 +179,12 @@ Foam::functionObjects::scalarTransport::scalarTransport dict.getOrDefault("phasePhiCompressed", "alphaPhiUn") ), D_(0), - constantD_(false), + alphaD_(1), + alphaDt_(1), + tol_(1), nCorr_(0), resetOnStartUp_(false), - schemesField_("unknown-schemesField"), - fvOptions_(mesh_), + constantD_(false), bounded01_(dict.getOrDefault("bounded01", true)) { read(dict); @@ -202,31 +200,30 @@ Foam::functionObjects::scalarTransport::scalarTransport } -// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * // - -Foam::functionObjects::scalarTransport::~scalarTransport() -{} - - // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // bool Foam::functionObjects::scalarTransport::read(const dictionary& dict) { - fvMeshFunctionObject::read(dict); + if (!fvMeshFunctionObject::read(dict)) + { + return false; + } dict.readIfPresent("phi", phiName_); dict.readIfPresent("rho", rhoName_); dict.readIfPresent("nut", nutName_); dict.readIfPresent("phase", phaseName_); - dict.readIfPresent("bounded01", bounded01_); + dict.readIfPresent("phasePhiCompressed", phasePhiCompressedName_); schemesField_ = dict.getOrDefault("schemesField", fieldName_); - constantD_ = dict.readIfPresent("D", D_); - alphaD_ = dict.getOrDefault("alphaD", 1); - alphaDt_ = dict.getOrDefault("alphaDt", 1); + dict.readIfPresent("alphaD", alphaD_); + dict.readIfPresent("alphaDt", alphaDt_); + dict.readIfPresent("tolerance", tol_); dict.readIfPresent("nCorr", nCorr_); dict.readIfPresent("resetOnStartUp", resetOnStartUp_); + constantD_ = dict.readIfPresent("D", D_); + dict.readIfPresent("bounded01", bounded01_); if (dict.found("fvOptions")) { @@ -256,6 +253,10 @@ bool Foam::functionObjects::scalarTransport::execute() scalar relaxCoeff = 0; mesh_.relaxEquation(schemesField_, relaxCoeff); + // Convergence monitor parameters + bool converged = false; + int iter = 0; + // Two phase scalar transport if (phaseName_ != "none") { @@ -272,7 +273,7 @@ bool Foam::functionObjects::scalarTransport::execute() // Solve tmp tTPhiUD; - for (label i = 0; i <= nCorr_; i++) + for (int i = 0; i <= nCorr_; ++i) { fvScalarMatrix sEqn ( @@ -285,9 +286,13 @@ bool Foam::functionObjects::scalarTransport::execute() sEqn.relax(relaxCoeff); fvOptions_.constrain(sEqn); - sEqn.solve(schemesField_); + + ++iter; + converged = (sEqn.solve(schemesField_).initialResidual() < tol_); tTPhiUD = sEqn.flux(); + + if (converged) break; } if (bounded01_) @@ -307,9 +312,8 @@ bool Foam::functionObjects::scalarTransport::execute() { const volScalarField& rho = lookupObject(rhoName_); - for (label i = 0; i <= nCorr_; i++) + for (int i = 0; i <= nCorr_; ++i) { - fvScalarMatrix sEqn ( fvm::ddt(rho, s) @@ -323,12 +327,14 @@ bool Foam::functionObjects::scalarTransport::execute() fvOptions_.constrain(sEqn); - sEqn.solve(schemesField_); + ++iter; + converged = (sEqn.solve(schemesField_).initialResidual() < tol_); + if (converged) break; } } else if (phi.dimensions() == dimVolume/dimTime) { - for (label i = 0; i <= nCorr_; i++) + for (int i = 0; i <= nCorr_; ++i) { fvScalarMatrix sEqn ( @@ -343,7 +349,9 @@ bool Foam::functionObjects::scalarTransport::execute() fvOptions_.constrain(sEqn); - sEqn.solve(schemesField_); + ++iter; + converged = (sEqn.solve(schemesField_).initialResidual() < tol_); + if (converged) break; } } else @@ -354,6 +362,14 @@ bool Foam::functionObjects::scalarTransport::execute() << dimVolume/dimTime << exit(FatalError); } + if (converged) + { + Log << type() << ": " << name() << ": " + << s.name() << " is converged." << nl + << tab << "initial-residual tolerance: " << tol_ << nl + << tab << "outer iteration: " << iter << nl; + } + Log << endl; return true; diff --git a/src/functionObjects/solvers/scalarTransport/scalarTransport.H b/src/functionObjects/solvers/scalarTransport/scalarTransport.H index 686d6cdba6..6a97eafb22 100644 --- a/src/functionObjects/solvers/scalarTransport/scalarTransport.H +++ b/src/functionObjects/solvers/scalarTransport/scalarTransport.H @@ -31,110 +31,134 @@ Group grpSolversFunctionObjects Description - Evolves a passive scalar transport equation. + Computes the transport equation for a passive scalar in single-phase or + two-phase flow, considering both incompressible and compressible cases: - - To specify the field name set the \c field entry - - To employ the same numerical schemes as another field set - the \c schemesField entry, - - The diffusivity can be set manually using the 'D' entry, retrieved - from the turbulence model or specified nut - - Alternatively if a turbulence model is available a turbulent diffusivity - may be constructed from the laminar and turbulent viscosities using the - optional diffusivity coefficients \c alphaD and \c alphaDt (which default - to 1): - \verbatim - D = alphaD*nu + alphaDt*nut - \endverbatim - - To specify a transport quantity within a phase enter phase. - - bounded01 bounds the transported scalar within 0 and 1. + \f[ + \frac{\partial \rho \, T}{\partial t} + + \nabla \cdot \left( \phi_\alpha \, T \right) + - \nabla \cdot (D_T \, \nabla T) + = \alpha \, S_T + \f] + + where: + \vartable + T | Passive scalar field + \rho | (Generic) Fluid density which is unity when not specified + \phi_\alpha | (Generic) Flux field + \alpha | Phase fraction which is unity for single-phase flows + D_T | Diffusivity representing the diffusive transport of T + S_T | Passive-scalar field source term + \endvartable Usage - Example of function object specification to solve a scalar transport - equation: + Minimal example in \c system/controlDict.functions: \verbatim - functions + scalarTransport1 { - scalar1 - { - type scalarTransport; - libs (solverFunctionObjects); + // Mandatory entries + type scalarTransport; + libs (solverFunctionObjects); - resetOnStartUp no; - region cabin; - field H2O; + // Optional entries + field ; + phi ; + rho ; + nut ; + phase ; + phasePhiCompressed ; + schemesField ; + bounded01 ; + D ; + alphaD ; + alphaDt ; + tolerance ; + nCorr ; + resetOnStartUp ; + fvOptions ; - - fvOptions - { - ... - } - } + // Inherited entries + ... } - \endverbatim - Example of function object specification to solve a residence time - in a two phase flow: - equation: - \verbatim - functions - { - sTransport - { - type scalarTransport; - libs (solverFunctionObjects); - - enabled true; - writeControl writeTime; - writeInterval 1; - - field s; - bounded01 false; - phase alpha.water; - - write true; - - fvOptions - { - unitySource - { - type scalarSemiImplicitSource; - enabled true; - - selectionMode all; - volumeMode specific; - - sources - { - s (1 0); - } - } - } - - resetOnStartUp false; - } - } - \endverbatim - - Where the entries comprise: + where: \table - Property | Description | Required | Default value - type | Type name: scalarTransport | yes | - field | Name of the scalar field | no | s - phi | Name of flux field | no | phi - rho | Name of density field | no | rho - phase | Name of the phase | no | none - nut | Name of the turbulence viscosity | no | none - D | Diffusion coefficient | no | auto generated - nCorr | Number of correctors | no | 0 - resetOnStartUp | Reset scalar to zero on start-up | no | no - schemesField | Name of field to specify schemes | no | field name - fvOptions | List of scalar sources | no | - bounded01 | Bounds scalar between 0-1 for multiphase | no | true - phasePhiCompressed | Compressed flux for VOF | no | alphaPhiUn + Property | Description | Type | Reqd | Deflt + type | Type name: scalarTransport | word | yes | - + libs | Library name: solverFunctionObjects | word | yes | - + field | Name of the passive-scalar field | word | no | s + phi | Name of flux field | word | no | phi + rho | Name of density field | word | no | rho + nut | Name of the turbulence viscosity | word | no | none + phase | Name of the phase | word | no | none + phasePhiCompressed | Name of compressed VOF flux | word | no | alphaPhiUn + schemesField | Name of field to specify schemes | word | no | field + bounded01 | Bounds scalar between 0-1 for multiphase | bool | no | true + D | Diffusion coefficient | scalar | no | - + alphaD | Laminar diffusivity coefficient | scalar | no | 1 + alphaDt | Turbulent diffusivity coefficient | scalar | no | 1 + tolerance | Outer-loop initial-residual tolerance | scalar | no | 1 + nCorr | Number of outer-loop correctors | int | no | 0 + resetOnStartUp | Flag to reset field to zero on start-up | bool | no | no + fvOptions | List of finite-volume options | dict | no | - \endtable -See also - Foam::functionObjects::fvMeshFunctionObject + The inherited entries are elaborated in: + - \link fvMeshFunctionObject.H \endlink + - \link fvOption.H \endlink + + An example of function object specification to solve a residence time + in a two-phase flow: + \verbatim + scalarTransport1 + { + // Mandatory entries + type scalarTransport; + libs (solverFunctionObjects); + + // Optional entries + field s; + bounded01 false; + phase alpha.water; + tolerance 1e-5; + resetOnStartUp false; + fvOptions + { + unitySource + { + type scalarSemiImplicitSource; + enabled true; + + selectionMode all; + volumeMode specific; + + sources + { + s (1 0); + } + } + } + + // Inherited entries + enabled true; + writeControl writeTime; + writeInterval 1; + } + \endverbatim + +Note + - To use the same numerical schemes as another field, + set the \c schemesField entry. + - The diffusivity can be set manually using the \c D entry, obtained + from the turbulence model or specified as `nut`. + - Alternatively, if a turbulence model is available, turbulent diffusivity + can be constructed from the laminar and turbulent viscosities using the + optional diffusivity coefficients \c alphaD and \c alphaDt + (which default to 1): + + \f[ + D = \alpha_D \, \nu + \alpha_{Dt} \, \nu_t + \f] SourceFiles scalarTransport.C @@ -163,49 +187,52 @@ class scalarTransport : public fvMeshFunctionObject { - // Private data + // Private Data + + //- Run-time selectable finite volume options, e.g. sources, constraints + fv::optionList fvOptions_; //- Name of the transport field. word fieldName_; - //- Name of flux field (optional) + //- Name of field whose schemes are used + word schemesField_; + + //- Name of flux field word phiName_; - //- Name of density field (optional) + //- Name of density field word rhoName_; - //- Name of turbulent viscosity field (optional) + //- Name of turbulent viscosity field word nutName_; - //- Name of phase field (optional) + //- Name of phase field word phaseName_; - //- Name of phase field compressed flux (optional) + //- Name of phase field compressed flux word phasePhiCompressedName_; - //- Diffusion coefficient (optional) + //- Diffusion coefficient scalar D_; - //- Flag to indicate whether a constant, uniform D_ is specified - bool constantD_; - - //- Laminar diffusion coefficient (optional) + //- Laminar diffusion coefficient scalar alphaD_; - //- Turbulent diffusion coefficient (optional) + //- Turbulent diffusion coefficient scalar alphaDt_; - //- Number of corrector iterations (optional) - label nCorr_; + //- Outer-loop initial-residual tolerance + scalar tol_; + + //- Number of corrector iterations + int nCorr_; //- Flag to reset the scalar to zero on start-up bool resetOnStartUp_; - //- Name of field whose schemes are used (optional) - word schemesField_; - - //- Run-time selectable finite volume options, e.g. sources, constraints - fv::optionList fvOptions_; + //- Flag to indicate whether a constant, uniform D_ is specified + bool constantD_; //- Bound scalar between 0-1 using MULES for multiphase case bool bounded01_; @@ -223,12 +250,6 @@ class scalarTransport const surfaceScalarField& phi ) const; - //- No copy construct - scalarTransport(const scalarTransport&) = delete; - - //- No copy assignment - void operator=(const scalarTransport&) = delete; - public: @@ -248,7 +269,7 @@ public: //- Destructor - virtual ~scalarTransport(); + virtual ~scalarTransport() = default; // Member Functions