diff --git a/src/TurbulenceModels/schemes/DEShybrid/DEShybrid.H b/src/TurbulenceModels/schemes/DEShybrid/DEShybrid.H index 015b6122d3..9a19dee379 100644 --- a/src/TurbulenceModels/schemes/DEShybrid/DEShybrid.H +++ b/src/TurbulenceModels/schemes/DEShybrid/DEShybrid.H @@ -6,7 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2015 OpenFOAM Foundation - Copyright (C) 2015-2021 OpenCFD Ltd. + Copyright (C) 2015-2022 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -226,32 +226,60 @@ class DEShybrid ) ); - const volScalarField factor + + const word factorName(IOobject::scopedName(typeName, "Factor")); + const fvMesh& mesh = this->mesh(); + + const IOobject factorIO ( - IOobject - ( - typeName + ":Factor", - this->mesh().time().timeName(), - this->mesh(), - IOobject::NO_READ, - IOobject::NO_WRITE - ), - max(sigmaMax_*tanh(pow(A, CH1_)), sigmaMin_) + factorName, + mesh.time().timeName(), + mesh, + IOobject::NO_READ, + IOobject::NO_WRITE ); if (blendedSchemeBaseName::debug) { - factor.write(); - } + auto* factorPtr = mesh.getObjectPtr(factorName); - return tmp - ( - new surfaceScalarField + if (!factorPtr) + { + factorPtr = + new volScalarField + ( + factorIO, + mesh, + dimensionedScalar(dimless, Zero) + ); + + const_cast(mesh).objectRegistry::store(factorPtr); + } + + auto& factor = *factorPtr; + + factor = max(sigmaMax_*tanh(pow(A, CH1_)), sigmaMin_); + + return tmp::New ( vf.name() + "BlendingFactor", fvc::interpolate(factor) - ) - ); + ); + } + else + { + const volScalarField factor + ( + factorIO, + max(sigmaMax_*tanh(pow(A, CH1_)), sigmaMin_) + ); + + return tmp::New + ( + vf.name() + "BlendingFactor", + fvc::interpolate(factor) + ); + } } diff --git a/src/TurbulenceModels/turbulenceModels/derivedFvPatchFields/wallFunctions/nutWallFunctions/nutUBlendedWallFunction/nutUBlendedWallFunctionFvPatchScalarField.C b/src/TurbulenceModels/turbulenceModels/derivedFvPatchFields/wallFunctions/nutWallFunctions/nutUBlendedWallFunction/nutUBlendedWallFunctionFvPatchScalarField.C index 38899885b4..47bdfbbd55 100644 --- a/src/TurbulenceModels/turbulenceModels/derivedFvPatchFields/wallFunctions/nutWallFunctions/nutUBlendedWallFunction/nutUBlendedWallFunctionFvPatchScalarField.C +++ b/src/TurbulenceModels/turbulenceModels/derivedFvPatchFields/wallFunctions/nutWallFunctions/nutUBlendedWallFunction/nutUBlendedWallFunctionFvPatchScalarField.C @@ -5,7 +5,7 @@ \\ / A nd | www.openfoam.com \\/ M anipulation | ------------------------------------------------------------------------------- - Copyright (C) 2016-2020 OpenCFD Ltd. + Copyright (C) 2016-2022 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -103,7 +103,8 @@ Foam::nutUBlendedWallFunctionFvPatchScalarField::calcUTau { const scalar yPlus = y[facei]*ut/nuw[facei]; const scalar uTauVis = magUp[facei]/yPlus; - const scalar uTauLog = kappa_*magUp[facei]/log(E_*yPlus); + const scalar uTauLog = + kappa_*magUp[facei]/log(max(E_*yPlus, 1 + 1e-4)); const scalar utNew = pow(pow(uTauVis, n_) + pow(uTauLog, n_), 1.0/n_); diff --git a/src/atmosphericModels/fvOptions/atmPlantCanopyUSource/atmPlantCanopyUSource.C b/src/atmosphericModels/fvOptions/atmPlantCanopyUSource/atmPlantCanopyUSource.C index f1044030bd..15f448ff57 100644 --- a/src/atmosphericModels/fvOptions/atmPlantCanopyUSource/atmPlantCanopyUSource.C +++ b/src/atmosphericModels/fvOptions/atmPlantCanopyUSource/atmPlantCanopyUSource.C @@ -6,7 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2020 ENERCON GmbH - Copyright (C) 2020-2021 OpenCFD Ltd. + Copyright (C) 2020-2022 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -28,6 +28,7 @@ License #include "atmPlantCanopyUSource.H" #include "addToRunTimeSelectionTable.H" +#include "fvmSup.H" // * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * // @@ -99,7 +100,7 @@ void Foam::fv::atmPlantCanopyUSource::addSup if (V_ > VSMALL) { // (SP:Eq. 42) - eqn -= (plantCd_*leafAreaDensity_*mag(U))*U; + eqn -= fvm::Sp(plantCd_*leafAreaDensity_*mag(U), U); } } @@ -115,7 +116,7 @@ void Foam::fv::atmPlantCanopyUSource::addSup if (V_ > VSMALL) { - eqn -= rho*(plantCd_*leafAreaDensity_*mag(U))*U; + eqn -= fvm::Sp(rho*plantCd_*leafAreaDensity_*mag(U), U); } } @@ -132,7 +133,7 @@ void Foam::fv::atmPlantCanopyUSource::addSup if (V_ > VSMALL) { - eqn -= alpha*rho*(plantCd_*leafAreaDensity_*mag(U))*U; + eqn -= fvm::Sp(alpha*rho*plantCd_*leafAreaDensity_*mag(U), U); } } diff --git a/src/finiteVolume/fields/fvPatchFields/derived/activePressureForceBaffleVelocity/activePressureForceBaffleVelocityFvPatchVectorField.C b/src/finiteVolume/fields/fvPatchFields/derived/activePressureForceBaffleVelocity/activePressureForceBaffleVelocityFvPatchVectorField.C index 2df16fccf5..e8148719a0 100644 --- a/src/finiteVolume/fields/fvPatchFields/derived/activePressureForceBaffleVelocity/activePressureForceBaffleVelocityFvPatchVectorField.C +++ b/src/finiteVolume/fields/fvPatchFields/derived/activePressureForceBaffleVelocity/activePressureForceBaffleVelocityFvPatchVectorField.C @@ -6,7 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2011-2016 OpenFOAM Foundation - Copyright (C) 2016-2020 OpenCFD Ltd. + Copyright (C) 2016-2022 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -260,11 +260,13 @@ void Foam::activePressureForceBaffleVelocityFvPatchVectorField::updateCoeffs() const labelList& nbrFaceCells = nbrPatch.patch().faceCells(); scalar valueDiff = 0; + scalar area = 0; // Add this side (p*area) forAll(cyclicFaceCells, facei) { valueDiff +=p[cyclicFaceCells[facei]]*mag(initCyclicSf_[facei]); + area += mag(initCyclicSf_[facei]); } // Remove other side @@ -275,7 +277,7 @@ void Foam::activePressureForceBaffleVelocityFvPatchVectorField::updateCoeffs() if (!fBased_) //pressure based then weighted by area { - valueDiff = valueDiff/gSum(patch().magSf()); + valueDiff = valueDiff/(area + VSMALL); } reduce(valueDiff, sumOp()); @@ -284,12 +286,13 @@ void Foam::activePressureForceBaffleVelocityFvPatchVectorField::updateCoeffs() { if (fBased_) { - Info<< "Force difference = " << valueDiff << endl; + Info<< "Force difference (threshold) = " << valueDiff + << "(" << minThresholdValue_ << ")" << endl; } else { - Info<< "Area-averaged pressure difference = " - << valueDiff << endl; + Info<< "Area-averaged pressure difference (threshold) = " + << valueDiff << "(" << minThresholdValue_ << ")" << endl; } } diff --git a/src/functionObjects/field/Make/options b/src/functionObjects/field/Make/options index a99cf3db84..cc5819e29c 100644 --- a/src/functionObjects/field/Make/options +++ b/src/functionObjects/field/Make/options @@ -14,6 +14,7 @@ EXE_INC = \ -I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \ -I$(LIB_SRC)/TurbulenceModels/incompressible/lnInclude \ -I$(LIB_SRC)/TurbulenceModels/compressible/lnInclude \ + -I$(LIB_SRC)/TurbulenceModels/schemes/lnInclude \ -I$(LIB_SRC)/transportModels/compressible/lnInclude \ -I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \ -I$(LIB_SRC)/thermophysicalModels/solidThermo/lnInclude \ @@ -39,6 +40,7 @@ LIB_LIBS = \ -lcompressibleTransportModels \ -lincompressibleTurbulenceModels \ -lcompressibleTurbulenceModels \ + -lturbulenceModelSchemes \ -lchemistryModel \ -lreactionThermophysicalModels \ -lpairPatchAgglomeration diff --git a/src/functionObjects/field/blendingFactor/blendingFactor.H b/src/functionObjects/field/blendingFactor/blendingFactor.H index ba5e09aa36..841994e144 100644 --- a/src/functionObjects/field/blendingFactor/blendingFactor.H +++ b/src/functionObjects/field/blendingFactor/blendingFactor.H @@ -6,7 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2013-2016 OpenFOAM Foundation - Copyright (C) 2016-2020 OpenCFD Ltd. + Copyright (C) 2016-2022 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -31,26 +31,14 @@ Group grpFieldFunctionObjects Description - Computes the blending coefficient employed by blended divergence schemes, - giving an indicator as to which of the schemes is active across the domain. + Computes blending factor as an indicator about which + of the schemes is active across the domain. - Blended schemes combine contributions from two schemes, i.e. \f$\phi_1\f$ - and \f$\phi_2\f$, using a weight field, i.e. \f$w\f$, such that the - effective scheme value, i.e. \f$\phi_{eff}\f$, is computed as follows: - - \f[ - \phi_{eff} = w \phi_1 + (1 - w) \phi_2 - \f] - - The weight field, i.e. \f$w\f$, is surface field and converted to a volume - field for easier post-processing by setting the cell value to one minus - the minimum of the face values. - - This conversion leads to blending indicator field whose values mean: + Blending factor values mean: \verbatim 0 = scheme 0 1 = scheme 1 - 0-1 = a blend between scheme 0 and scheme 1 + 0-1 = a blending factor between scheme 0 and scheme 1 \endverbatim Operands: @@ -66,30 +54,28 @@ Usage \verbatim blendingFactor1 { - // Mandatory entries (unmodifiable) + // Mandatory entries type blendingFactor; libs (fieldFunctionObjects); - - // Mandatory (inherited) entry (runtime modifiable) field ; - // Optional entries (runtime modifiable) + // Optional entries phi phi; tolerance 0.001; - // Optional (inherited) entries + // Inherited entries ... } \endverbatim where the entries mean: \table - Property | Description | Type | Req'd | Deflt - type | Type name: blendingFactor | word | yes | - - libs | Library name: fieldFunctionObjects | word | yes | - - field | Name of the operand field | word | yes | - - phi | Name of flux field | word | no | phi - tolerance | Tolerance for number of blended cells | scalar | no | 0.001 + Property | Description | Type | Reqd | Deflt + type | Type name: blendingFactor | word | yes | - + libs | Library name: fieldFunctionObjects | word | yes | - + field | Name of the operand field | word | yes | - + phi | Name of flux field | word | no | phi + tolerance | Tolerance for number of blended cells | scalar | no | 0.001 \endtable The inherited entries are elaborated in: @@ -99,13 +85,12 @@ Usage Usage by the \c postProcess utility is not available. -See also - - Foam::functionObject - - Foam::functionObjects::fieldExpression - - Foam::functionObjects::fvMeshFunctionObject - - Foam::functionObjects::writeFile - - Foam::CoBlended - - ExtendedCodeGuide::functionObjects::field::blendingFactor +Note + - The \c blendingFactor function object overwrites the \c DEShybrid:Factor + field internally when \c blendedSchemeBase debug flag is active. + However, users are allowed to write out the original \c DEShybrid:Factor + field by executing the \c writeObjects function object before + any \c blendingFactor function object execution. SourceFiles blendingFactor.C diff --git a/src/functionObjects/field/blendingFactor/blendingFactorTemplates.C b/src/functionObjects/field/blendingFactor/blendingFactorTemplates.C index 2eb39896c2..e238f4677e 100644 --- a/src/functionObjects/field/blendingFactor/blendingFactorTemplates.C +++ b/src/functionObjects/field/blendingFactor/blendingFactorTemplates.C @@ -6,7 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2013-2016 OpenFOAM Foundation - Copyright (C) 2016 OpenCFD Ltd. + Copyright (C) 2016-2022 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -30,6 +30,7 @@ License #include "boundedConvectionScheme.H" #include "blendedSchemeBase.H" #include "fvcCellReduce.H" +#include "DEShybrid.H" // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // @@ -72,12 +73,18 @@ void Foam::functionObjects::blendingFactor::calcBlendingFactor // Convert into vol field whose values represent the local face minima // Note: - // - factor applied to 1st scheme, and (1-factor) to 2nd scheme // - not using the store(...) mechanism due to need to correct BCs - volScalarField& indicator = - lookupObjectRef(resultName_); + auto& indicator = lookupObjectRef(resultName_); + + if (isA>(blendedScheme)) + { + indicator = fvc::cellReduce(factorf, minEqOp(), GREAT); + } + else + { + indicator = 1 - fvc::cellReduce(factorf, minEqOp(), GREAT); + } - indicator = 1 - fvc::cellReduce(factorf, minEqOp(), GREAT); indicator.correctBoundaryConditions(); } diff --git a/src/mesh/snappyHexMesh/snappyHexMeshDriver/snappyLayerDriver.C b/src/mesh/snappyHexMesh/snappyHexMeshDriver/snappyLayerDriver.C index 628a00400b..abb6ce6493 100644 --- a/src/mesh/snappyHexMesh/snappyHexMeshDriver/snappyLayerDriver.C +++ b/src/mesh/snappyHexMesh/snappyHexMeshDriver/snappyLayerDriver.C @@ -6,7 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2011-2015 OpenFOAM Foundation - Copyright (C) 2015-2021 OpenCFD Ltd. + Copyright (C) 2015-2022 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -5154,7 +5154,7 @@ void Foam::snappyLayerDriver::addLayers meshRefinement::updateList ( map.cellMap(), - 0, + label(0), cellNLayers ); forAll(cellNLayers, i) @@ -5237,7 +5237,7 @@ void Foam::snappyLayerDriver::addLayers meshRefinement::updateList ( newToOldPatchPoints, - 0, + label(0), basePatchNLayers ); meshRefinement::updateList @@ -5267,13 +5267,13 @@ void Foam::snappyLayerDriver::addLayers meshRefinement::updateList ( newToOldPatchPoints, - 0, + label(0), deltaNLayers ); meshRefinement::updateList ( newToOldPatchPoints, - 0, + label(0), nAddedLayers ); } diff --git a/src/phaseSystemModels/reactingEuler/multiphaseSystem/diameterModels/linearTsubDiameter/linearTsubDiameter.C b/src/phaseSystemModels/reactingEuler/multiphaseSystem/diameterModels/linearTsubDiameter/linearTsubDiameter.C index bfddaca6e3..fb91423137 100644 --- a/src/phaseSystemModels/reactingEuler/multiphaseSystem/diameterModels/linearTsubDiameter/linearTsubDiameter.C +++ b/src/phaseSystemModels/reactingEuler/multiphaseSystem/diameterModels/linearTsubDiameter/linearTsubDiameter.C @@ -6,7 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2018-2019 OpenFOAM Foundation - Copyright (C) 2020 OpenCFD Ltd. + Copyright (C) 2020-2022 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -126,7 +126,7 @@ void Foam::diameterModels::linearTsub::correct() const volScalarField Tsub ( - liquid.thermo().T() - satModel.Tsat(liquid.thermo().p()) + satModel.Tsat(liquid.thermo().p()) - liquid.thermo().T() ); d_ = max @@ -135,7 +135,7 @@ void Foam::diameterModels::linearTsub::correct() min ( d2_, - (d1_*(Tsub - Tsub2_) + d2_*(Tsub - Tsub1_))/(Tsub2_ - Tsub1_) + (d1_*(Tsub - Tsub2_) - d2_*(Tsub - Tsub1_))/(Tsub1_ - Tsub2_) ) ); } diff --git a/src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/TDACChemistryModel.C b/src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/TDACChemistryModel.C index b249e9ac54..8bfbca15fb 100644 --- a/src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/TDACChemistryModel.C +++ b/src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/TDACChemistryModel.C @@ -5,7 +5,7 @@ \\ / A nd | www.openfoam.com \\/ M anipulation | ------------------------------------------------------------------------------- - Copyright (C) 2016-2017 OpenFOAM Foundation + Copyright (C) 2016-2021 OpenFOAM Foundation Copyright (C) 2016-2021 OpenCFD Ltd. ------------------------------------------------------------------------------- License @@ -864,7 +864,7 @@ Foam::scalar Foam::TDACChemistryModel::solve << " " << nActiveSpecies/nAvg << endl; } - if (Pstream::parRun()) + if (reduced && Pstream::parRun()) { List active(composition.active()); Pstream::listCombineAllGather(active, orEqOp()); diff --git a/src/thermophysicalModels/thermophysicalPropertiesFvPatchFields/liquidProperties/humidityTemperatureCoupledMixed/humidityTemperatureCoupledMixedFvPatchScalarField.C b/src/thermophysicalModels/thermophysicalPropertiesFvPatchFields/liquidProperties/humidityTemperatureCoupledMixed/humidityTemperatureCoupledMixedFvPatchScalarField.C index 1baf7faece..61a150ce10 100644 --- a/src/thermophysicalModels/thermophysicalPropertiesFvPatchFields/liquidProperties/humidityTemperatureCoupledMixed/humidityTemperatureCoupledMixedFvPatchScalarField.C +++ b/src/thermophysicalModels/thermophysicalPropertiesFvPatchFields/liquidProperties/humidityTemperatureCoupledMixed/humidityTemperatureCoupledMixedFvPatchScalarField.C @@ -5,7 +5,7 @@ \\ / A nd | www.openfoam.com \\/ M anipulation | ------------------------------------------------------------------------------- - Copyright (C) 2015-2021 OpenCFD Ltd. + Copyright (C) 2015-2022 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -77,9 +77,10 @@ Foam::humidityTemperatureCoupledMixedFvPatchScalarField::htcCondensation const scalar Re ) const { - if (Tsat > 295 && Tsat < 373) + // (BLID:Eq. 10.52-10.53) + if (Tsat > 295.15 && Tsat < 373.15) { - return 51104 + 2044*Tsat; + return 51104 + 2044*(Tsat - 273.15); } else { diff --git a/src/thermophysicalModels/thermophysicalPropertiesFvPatchFields/liquidProperties/humidityTemperatureCoupledMixed/humidityTemperatureCoupledMixedFvPatchScalarField.H b/src/thermophysicalModels/thermophysicalPropertiesFvPatchFields/liquidProperties/humidityTemperatureCoupledMixed/humidityTemperatureCoupledMixedFvPatchScalarField.H index bc5d081cdb..64cfd0f070 100644 --- a/src/thermophysicalModels/thermophysicalPropertiesFvPatchFields/liquidProperties/humidityTemperatureCoupledMixed/humidityTemperatureCoupledMixedFvPatchScalarField.H +++ b/src/thermophysicalModels/thermophysicalPropertiesFvPatchFields/liquidProperties/humidityTemperatureCoupledMixed/humidityTemperatureCoupledMixedFvPatchScalarField.H @@ -5,7 +5,7 @@ \\ / A nd | www.openfoam.com \\/ M anipulation | ------------------------------------------------------------------------------- - Copyright (C) 2015-2021 OpenCFD Ltd. + Copyright (C) 2015-2022 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -52,13 +52,19 @@ Description It assumes a drop-wise type of condensation, whereby its heat transfer Nusselt number is calculated using: \f{eqnarray*}{ - 51104 + 2044 T & T > 295 & T < 373 \\ - 255510 & T > 373 & + 51104 + 2044 (T - 273.15) & T > 295 & T < 373 \\ + 255510 & T > 373 & \f} - Reference: - - T. Bergam, A.Lavine, F. Incropera and D. Dewitt. Heat and Mass Transfer. - 7th Edition. Chapter 10. + References: + \verbatim + Standard models (tag:BLID): + Bergman, T. L., Lavine, A. S., + Incropera, F. P., & Dewitt, D. P. (2011). + Fundamentals of heat and mass transfer. + John Wiley & Sons. 7th Edition. Chapter 10. + ISBN:9780470501979 + \endverbatim The mass transfer correlation used is: diff --git a/tutorials/incompressible/pimpleFoam/RAS/propeller/system/forces b/tutorials/incompressible/pimpleFoam/RAS/propeller/system/forces index 4f32e012ab..15cbfd0162 100644 --- a/tutorials/incompressible/pimpleFoam/RAS/propeller/system/forces +++ b/tutorials/incompressible/pimpleFoam/RAS/propeller/system/forces @@ -14,8 +14,6 @@ forces writeControl timeStep; timeInterval 1; - log yes; - patches ("propeller.*"); rho rhoInf; // Indicates incompressible log true; diff --git a/tutorials/lagrangian/coalChemistryFoam/simplifiedSiwek/constant/limestoneCloud1Properties b/tutorials/lagrangian/coalChemistryFoam/simplifiedSiwek/constant/limestoneCloud1Properties index 541fa8fe6d..6108d91d09 100644 --- a/tutorials/lagrangian/coalChemistryFoam/simplifiedSiwek/constant/limestoneCloud1Properties +++ b/tutorials/lagrangian/coalChemistryFoam/simplifiedSiwek/constant/limestoneCloud1Properties @@ -87,7 +87,7 @@ subModels { minValue 5e-06; maxValue 0.000565; - d 4.8e-05; + lambda 4.8e-05; n 0.5; } } diff --git a/tutorials/lagrangian/reactingParcelFoam/liquidFilmStepWithSprinkles/constant/reactingCloud1Properties b/tutorials/lagrangian/reactingParcelFoam/liquidFilmStepWithSprinkles/constant/reactingCloud1Properties index df0dea3cb8..7f178048ee 100644 --- a/tutorials/lagrangian/reactingParcelFoam/liquidFilmStepWithSprinkles/constant/reactingCloud1Properties +++ b/tutorials/lagrangian/reactingParcelFoam/liquidFilmStepWithSprinkles/constant/reactingCloud1Properties @@ -104,7 +104,7 @@ subModels { minValue 5e-04; maxValue 0.0012; - d 7.5e-05; + lambda 7.5e-05; n 0.5; } } diff --git a/tutorials/lagrangian/reactingParcelFoam/verticalChannelLTS/constant/thermo.incompressiblePoly b/tutorials/lagrangian/reactingParcelFoam/verticalChannelLTS/constant/thermo.incompressiblePoly index 9b364f90eb..ace33ab1cc 100644 --- a/tutorials/lagrangian/reactingParcelFoam/verticalChannelLTS/constant/thermo.incompressiblePoly +++ b/tutorials/lagrangian/reactingParcelFoam/verticalChannelLTS/constant/thermo.incompressiblePoly @@ -49,8 +49,8 @@ H2O } thermodynamics { - Hf -13423000; - Sf 10482; + Hf 0; + Sf 0; CpCoeffs<8> ( 1563.1 1.604 -0.0029334 3.2168e-06 -1.1571e-09 0 0 0 ); } transport diff --git a/tutorials/lagrangian/simpleReactingParcelFoam/verticalChannel/system/fvOptions b/tutorials/lagrangian/simpleReactingParcelFoam/verticalChannel/system/fvOptions new file mode 100644 index 0000000000..0145de64bb --- /dev/null +++ b/tutorials/lagrangian/simpleReactingParcelFoam/verticalChannel/system/fvOptions @@ -0,0 +1,26 @@ +/*--------------------------------*- C++ -*----------------------------------*\ +| ========= | | +| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox | +| \\ / O peration | Version: v2112 | +| \\ / A nd | Website: www.openfoam.com | +| \\/ M anipulation | | +\*---------------------------------------------------------------------------*/ +FoamFile +{ + version 2.0; + format ascii; + class dictionary; + object fvOptions; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +limitT +{ + type limitTemperature; + min 101; + max 600; + selectionMode all; +} + + +//************************************************************************** //