diff --git a/src/thermoTools/Make/files b/src/thermoTools/Make/files index b66d55ed1360452d2d996dab7fdfdfd0eda6de9d..b1a3a6cef0934ad1fd6944ae7157dde775aa8965 100644 --- a/src/thermoTools/Make/files +++ b/src/thermoTools/Make/files @@ -20,5 +20,6 @@ $(BCs)/lumpedMassWallTemperature/lumpedMassWallTemperatureFvPatchScalarField.C $(BCs)/wallFunctions/alphatWallFunctions/alphatWallFunction/alphatWallFunctionFvPatchScalarField.C $(BCs)/wallFunctions/alphatWallFunctions/alphatJayatillekeWallFunction/alphatJayatillekeWallFunctionFvPatchScalarField.C +$(BCs)/wallFunctions/sorptionWallFunction/sorptionWallFunctionFvPatchScalarField.C LIB = $(FOAM_LIBBIN)/libthermoTools diff --git a/src/thermoTools/derivedFvPatchFields/wallFunctions/sorptionWallFunction/sorptionWallFunctionFvPatchScalarField.C b/src/thermoTools/derivedFvPatchFields/wallFunctions/sorptionWallFunction/sorptionWallFunctionFvPatchScalarField.C new file mode 100644 index 0000000000000000000000000000000000000000..d289f63444821e90b7d16cd2444dd40acff6a341 --- /dev/null +++ b/src/thermoTools/derivedFvPatchFields/wallFunctions/sorptionWallFunction/sorptionWallFunctionFvPatchScalarField.C @@ -0,0 +1,429 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | www.openfoam.com + \\/ M anipulation | +------------------------------------------------------------------------------- + Copyright (C) 2022 OpenCFD Ltd. +------------------------------------------------------------------------------- +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 "sorptionWallFunctionFvPatchScalarField.H" +#include "addToRunTimeSelectionTable.H" +#include "fvPatchFieldMapper.H" +#include "volFields.H" +#include "wallDist.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +// * * * * * * * * * * * * * * * Local Functions * * * * * * * * * * * * * * // + +//- Estimate the y* at the intersection of the two sublayers +static scalar calcYStarLam +( + const scalar kappa, + const scalar E, + const scalar Sc, + const scalar Sct, + const scalar Pc +) +{ + scalar ypl = 11; + + for (int iter = 0; iter < 10; ++iter) + { + // (F:Eq. 5.5) + ypl = (log(max(E*ypl, scalar(1)))/kappa + Pc)*Sct/Sc; + } + + return ypl; +} + + +// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // + +tmp<scalarField> sorptionWallFunctionFvPatchScalarField::yPlus() const +{ + // Calculate fields of interest + const label patchi = patch().index(); + + const auto& k = db().lookupObject<volScalarField>(kName_); + tmp<scalarField> tkwc = k.boundaryField()[patchi].patchInternalField(); + const scalarField& kwc = tkwc.cref(); + + const auto& nu = db().lookupObject<volScalarField>(nuName_); + tmp<scalarField> tnuwc = nu.boundaryField()[patchi].patchInternalField(); + const scalarField& nuwc = tnuwc.cref(); + + const volScalarField& y = wallDist::New(internalField().mesh()).y(); + tmp<scalarField> tywc = y.boundaryField()[patchi].patchInternalField(); + const scalarField& ywc = tywc.cref(); + + + // Calculate the empirical constant given by (Jayatilleke, 1966) (FDC:Eq. 6) + const scalar Pc = + 9.24*(pow(Sc_/Sct_, 0.75) - 1)*(1 + 0.28*exp(-0.007*Sc_/Sct_)); + const scalar Cmu25 = pow025(wallCoeffs_.Cmu()); + const scalar kappa = wallCoeffs_.kappa(); + const scalar E = wallCoeffs_.E(); + + auto tyPlus = tmp<scalarField>::New(patch().size(), Zero); + auto& yPlus = tyPlus.ref(); + + forAll(yPlus, facei) + { + // (FDC:Eq. 3) + const scalar yStar = Cmu25*sqrt(kwc[facei])*ywc[facei]/nuwc[facei]; + + // (FDC:Eq. 4) + const scalar yPlusVis = Sc_*yStar; + + // (FDC:Eq. 5) + const scalar yPlusLog = Sct_*(log(max(E*yStar, 1 + 1e-4))/kappa + Pc); + + switch (blender_) + { + case blenderType::EXPONENTIAL: + { + // (FDC:Eq. 2) + const scalar Gamma = + 0.01*pow4(Sc_*yStar)/(1 + 5*pow3(Sc_)*yStar); + const scalar invGamma = scalar(1)/max(Gamma, ROOTVSMALL); + + // (FDC:Eq. 1) + yPlus[facei] = yPlusVis*exp(-Gamma) + yPlusLog*exp(-invGamma); + break; + } + + case blenderType::STEPWISE: + { + static const scalar yStarLam = + calcYStarLam(kappa, E, Sc_, Sct_, Pc); + + // (F:Eq. 5.3) + if (yStar < yStarLam) + { + yPlus[facei] = yPlusVis; + } + else + { + yPlus[facei] = yPlusLog; + } + break; + } + + case blenderType::BINOMIAL: + { + yPlus[facei] = + pow + ( + pow(yPlusVis, n_) + pow(yPlusLog, n_), + scalar(1)/n_ + ); + break; + } + + case blenderType::MAX: + { + yPlus[facei] = max(yPlusVis, yPlusLog); + break; + } + + case blenderType::TANH: + { + const scalar phiTanh = tanh(pow4(0.1*yStar)); + const scalar b1 = yPlusVis + yPlusLog; + const scalar b2 = + pow(pow(yPlusVis, 1.2) + pow(yPlusLog, 1.2), 1.0/1.2); + + yPlus[facei] = phiTanh*b1 + (1 - phiTanh)*b2; + break; + } + } + } + + return tyPlus; +} + + +tmp<scalarField> sorptionWallFunctionFvPatchScalarField::flux() const +{ + // Calculate fields of interest + const label patchi = patch().index(); + + const auto& k = db().lookupObject<volScalarField>(kName_); + tmp<scalarField> tkwc = k.boundaryField()[patchi].patchInternalField(); + + const volScalarField& y = wallDist::New(internalField().mesh()).y(); + tmp<scalarField> tywc = y.boundaryField()[patchi].patchInternalField(); + + + // Calculate mass-transfer coefficient field (FDC:Eqs. 8-9) + tmp<scalarField> ta = + laminar_ + ? D_/tywc + : pow025(wallCoeffs_.Cmu())*sqrt(tkwc)/yPlus(); + + + // Calculate wall-surface concentration + const auto& Csurf = *this; + + // Calculate wall-adjacent concentration + const scalar t = db().time().timeOutputValue(); + tmp<scalarField> tkAbs = kAbsPtr_->value(t); + tmp<scalarField> tCstar = Csurf/tkAbs; + + // Calculate near-wall cell concentration + const word& fldName = internalField().name(); + const auto& C = db().lookupObject<volScalarField>(fldName); + tmp<scalarField> tCwc = C.boundaryField()[patchi].patchInternalField(); + + // Return adsorption or absorption/permeation flux + // normalised by the mass-transfer coefficient (FDC:Fig. 1) + return (tCstar - tCwc)/ta; +} + + +void sorptionWallFunctionFvPatchScalarField::writeLocalEntries +( + Ostream& os +) const +{ + wallFunctionBlenders::writeEntries(os); + os.writeEntryIfDifferent<bool>("laminar", false, laminar_); + os.writeEntry("Sc", Sc_); + os.writeEntry("Sct", Sct_); + os.writeEntryIfDifferent<scalar>("D", -1, D_); + wallCoeffs_.writeEntries(os); + os.writeEntryIfDifferent<word>("k", "k", kName_); + os.writeEntryIfDifferent<word>("nu", "nu", nuName_); + + if (kAbsPtr_) + { + kAbsPtr_->writeData(os); + } +} + + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +sorptionWallFunctionFvPatchScalarField::sorptionWallFunctionFvPatchScalarField +( + const fvPatch& p, + const DimensionedField<scalar, volMesh>& iF +) +: + fixedGradientFvPatchScalarField(p, iF), + wallFunctionBlenders(), + laminar_(false), + kAbsPtr_(nullptr), + Sc_(1), + Sct_(1), + D_(-1), + kName_("k"), + nuName_("nu"), + wallCoeffs_() +{} + + +sorptionWallFunctionFvPatchScalarField::sorptionWallFunctionFvPatchScalarField +( + const sorptionWallFunctionFvPatchScalarField& ptf, + const fvPatch& p, + const DimensionedField<scalar, volMesh>& iF, + const fvPatchFieldMapper& mapper +) +: + fixedGradientFvPatchScalarField(ptf, p, iF, mapper), + wallFunctionBlenders(ptf), + laminar_(ptf.laminar_), + kAbsPtr_(ptf.kAbsPtr_.clone(patch().patch())), + Sc_(ptf.Sc_), + Sct_(ptf.Sct_), + D_(ptf.D_), + kName_(ptf.kName_), + nuName_(ptf.nuName_), + wallCoeffs_(ptf.wallCoeffs_) +{} + + +sorptionWallFunctionFvPatchScalarField::sorptionWallFunctionFvPatchScalarField +( + const fvPatch& p, + const DimensionedField<scalar, volMesh>& iF, + const dictionary& dict +) +: + fixedGradientFvPatchScalarField(p, iF), + wallFunctionBlenders(dict, blenderType::STEPWISE, scalar(2)), + laminar_(dict.getOrDefault<bool>("laminar", false)), + kAbsPtr_(PatchFunction1<scalar>::New(p.patch(), "kAbs", dict)), + Sc_(dict.getCheck<scalar>("Sc", scalarMinMax::ge(0))), + Sct_(dict.getCheck<scalar>("Sct", scalarMinMax::ge(0))), + D_(dict.getOrDefault<scalar>("D", -1)), + kName_(dict.getOrDefault<word>("k", "k")), + nuName_(dict.getOrDefault<word>("nu", "nu")), + wallCoeffs_(dict) +{ + if (laminar_) + { + if (D_ < 0) + { + FatalIOErrorInFunction(dict) + << "Molecular diffusion coefficient cannot be non-positive. " + << "D = " << D_ + << exit(FatalIOError); + } + } + + if (!kAbsPtr_) + { + FatalIOErrorInFunction(dict) + << "Adsorption or absorption coefficient is not set." + << exit(FatalIOError); + } + + if (dict.found("value") && dict.found("gradient")) + { + fvPatchField<scalar>::operator = + ( + Field<scalar>("value", dict, p.size()) + ); + gradient() = Field<scalar>("gradient", dict, p.size()); + } + else + { + fvPatchField<scalar>::operator=(patchInternalField()); + gradient() = 0.0; + } +} + + +sorptionWallFunctionFvPatchScalarField::sorptionWallFunctionFvPatchScalarField +( + const sorptionWallFunctionFvPatchScalarField& swfpsf +) +: + fixedGradientFvPatchScalarField(swfpsf), + wallFunctionBlenders(swfpsf), + laminar_(swfpsf.laminar_), + kAbsPtr_(swfpsf.kAbsPtr_.clone(patch().patch())), + Sc_(swfpsf.Sc_), + Sct_(swfpsf.Sct_), + D_(swfpsf.D_), + kName_(swfpsf.kName_), + nuName_(swfpsf.nuName_), + wallCoeffs_(swfpsf.wallCoeffs_) +{} + + +sorptionWallFunctionFvPatchScalarField::sorptionWallFunctionFvPatchScalarField +( + const sorptionWallFunctionFvPatchScalarField& swfpsf, + const DimensionedField<scalar, volMesh>& iF +) +: + fixedGradientFvPatchScalarField(swfpsf, iF), + wallFunctionBlenders(swfpsf), + laminar_(swfpsf.laminar_), + kAbsPtr_(swfpsf.kAbsPtr_.clone(patch().patch())), + Sc_(swfpsf.Sc_), + Sct_(swfpsf.Sct_), + D_(swfpsf.D_), + kName_(swfpsf.kName_), + nuName_(swfpsf.nuName_), + wallCoeffs_(swfpsf.wallCoeffs_) +{} + + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + +void sorptionWallFunctionFvPatchScalarField::autoMap +( + const fvPatchFieldMapper& mapper +) +{ + fixedGradientFvPatchScalarField::autoMap(mapper); + + if (kAbsPtr_) + { + kAbsPtr_->autoMap(mapper); + } +} + + +void sorptionWallFunctionFvPatchScalarField::rmap +( + const fvPatchScalarField& ptf, + const labelList& addr +) +{ + fixedGradientFvPatchScalarField::rmap(ptf, addr); + + const auto& swfptf = + refCast<const sorptionWallFunctionFvPatchScalarField>(ptf); + + if (kAbsPtr_) + { + kAbsPtr_->rmap(swfptf.kAbsPtr_(), addr); + } +} + + +void sorptionWallFunctionFvPatchScalarField::updateCoeffs() +{ + if (updated()) + { + return; + } + + gradient() = flux()/patch().deltaCoeffs(); + + fixedGradientFvPatchScalarField::updateCoeffs(); +} + + +void sorptionWallFunctionFvPatchScalarField::write(Ostream& os) const +{ + fixedGradientFvPatchScalarField::write(os); + + writeLocalEntries(os); + + writeEntry("value", os); +} + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +makePatchTypeField +( + fvPatchScalarField, + sorptionWallFunctionFvPatchScalarField +); + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace Foam + +// ************************************************************************* // diff --git a/src/thermoTools/derivedFvPatchFields/wallFunctions/sorptionWallFunction/sorptionWallFunctionFvPatchScalarField.H b/src/thermoTools/derivedFvPatchFields/wallFunctions/sorptionWallFunction/sorptionWallFunctionFvPatchScalarField.H new file mode 100644 index 0000000000000000000000000000000000000000..bcc7c2522ea1f7c72101b5b0546d7ec93b93fc64 --- /dev/null +++ b/src/thermoTools/derivedFvPatchFields/wallFunctions/sorptionWallFunction/sorptionWallFunctionFvPatchScalarField.H @@ -0,0 +1,323 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | www.openfoam.com + \\/ M anipulation | +------------------------------------------------------------------------------- + Copyright (C) 2022 OpenCFD Ltd. +------------------------------------------------------------------------------- +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/>. + +Class + Foam::sorptionWallFunctionFvPatchScalarField + +Group + grpWallFunctions + +Description + The \c sorptionWallFunction is a wall boundary condition to + specify scalar/concentration gradient for turbulent and laminar flows. + + The governing equation of the boundary condition is: + + \f[ + \nabla C = \frac{C^* - C_p}{\Delta_y} = \frac{F}{a \Delta_y} + \f] + + with + + \f[ + C^* = \frac{C_{surf}}{K} + \f] + + and with the mass-transfer coefficient is calculated for turbulent flows + + \f[ + a = \frac{C_\mu^{0.25} k_p^{0.5}}{y^+_{blended}} + \f] + + or for laminar-flow and molecular-diffusion-only states + + \f[ + a = \frac{D_m}{y_1} + \f] + + where + \vartable + \nabla C | Gradient of concentration + C^* | Wall-adjacent concentration + C_p | Near-wall cell concentration + \Delta_y | First-cell centre wall distance + F | Flux of concentration + a | Mass-transfer coefficient + C_{surf} | Wall-surface concentration + K | Adsorption or absorption/permeation coefficient + C_\mu | Empirical model coefficient + k_p | Turbulent kinetic energy in near-wall cell + y^+_{blended} | Non-dimensional blended near-wall cell height + D_m | Molecular-diffusion coefficient + y_1 | First-cell centre wall distance + \endvartable + + Required fields: + \verbatim + x | Arbitrary scalar field, e.g. species, passive scalars etc. + \endverbatim + + Reference: + \verbatim + Standard model for exponential blending (tag:FDC): + Foat, T., Drodge, J., Charleson, A., Whatmore, B., + Pownall, S., Glover, P., ... & Marr, A. (2022). + Predicting vapour transport from semi-volatile organic + compounds concealed within permeable packaging. + International Journal of Heat and Mass Transfer, 183, 122012. + DOI:10.1016/j.ijheatmasstransfer.2021.122012 + + Standard model for stepwise blending (tag:F): + Foat, T. (2021). + Modelling vapour transport in indoor environments for + improved detection of explosives using dogs. + Doctoral dissertation. University of Southampton. + URI:http://eprints.soton.ac.uk/id/eprint/456709 + \endverbatim + +Usage + Example of the boundary condition specification: + \verbatim + <patchName> + { + // Mandatory entries + type sorptionWallFunction; + Sc <scalar>; + Sct <scalar>; + kAbs <PatchFunction1<scalar>>; + + // Optional entries + laminar <bool>; + D <scalar>; + kName <word>; + nuName <word>; + + // Inherited entries + Cmu <scalar>; + kappa <scalar>; + E <scalar>; + blending <word>; + ... + } + \endverbatim + + where the entries mean: + \table + Property | Description | Type | Reqd | Deflt + type | Type name: sorptionWallFunction | word | yes | - + Sc | Schmidt number | scalar | yes | - + Sct | Turbulent Schmidt number | scalar | yes | - + kAbs | Adsorption or absorption/permeation coefficient <!-- + --> | PatchFunction1\<scalar\> | yes | - + laminar | Flag to calculate mass-transfer coefficient under the <!-- + --> laminar-flow or molecular-diffusion-only states <!-- + --> | bool | no | false + kName | Name of operand turbulent kinetic energy field | word | no | k + nuName | Name of operand kinematic viscosity field | word | no | nu + \endtable + + The inherited entries are elaborated in: + - \link fixedGradientFvPatchField.H \endlink + - \link wallFunctionCoefficients.H \endlink + - \link wallFunctionBlenders.H \endlink + - \link PatchFunction1.H \endlink + +SourceFiles + sorptionWallFunctionFvPatchScalarField.C + +\*---------------------------------------------------------------------------*/ + +#ifndef Foam_sorptionWallFunctionFvPatchScalarFields_H +#define Foam_sorptionWallFunctionFvPatchScalarFields_H + +#include "fvPatchFields.H" +#include "fixedGradientFvPatchFields.H" +#include "wallFunctionCoefficients.H" +#include "wallFunctionBlenders.H" +#include "PatchFunction1.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +/*---------------------------------------------------------------------------*\ + Class sorptionWallFunctionFvPatchScalarField Declaration +\*---------------------------------------------------------------------------*/ + +class sorptionWallFunctionFvPatchScalarField +: + public fixedGradientFvPatchScalarField, + private wallFunctionBlenders +{ + //- Flag to calculate mass-transfer coefficient + //- under the laminar-flow or molecular-diffusion-only states + bool laminar_; + + //- Adsorption or absorption/permeation coefficient + autoPtr<PatchFunction1<scalar>> kAbsPtr_; + + //- Schmidt number + scalar Sc_; + + //- Turbulent Schmidt number + scalar Sct_; + + //- Molecular diffusion coefficient + scalar D_; + + //- Name of operand turbulent kinetic energy field + word kName_; + + //- Name of operand kinematic viscosity field + word nuName_; + + //- Standard wall-function coefficients + wallFunctionCoefficients wallCoeffs_; + + + // Private Member Functions + + //- Return the non-dimensional near-wall cell height field + //- blended between the viscous and inertial sublayers + tmp<scalarField> yPlus() const; + + //- Return the flux normalised by the mass-transfer coefficient + tmp<scalarField> flux() const; + + //- Write local wall-function variables + void writeLocalEntries(Ostream&) const; + + +public: + + //- Runtime type information + TypeName("sorptionWallFunction"); + + + // Constructors + + //- Construct from patch and internal field + sorptionWallFunctionFvPatchScalarField + ( + const fvPatch&, + const DimensionedField<scalar, volMesh>& + ); + + //- Construct from patch, internal field and dictionary + sorptionWallFunctionFvPatchScalarField + ( + const fvPatch&, + const DimensionedField<scalar, volMesh>&, + const dictionary& + ); + + //- Construct by mapping given + //- sorptionWallFunctionFvPatchScalarField onto + //- a new patch + sorptionWallFunctionFvPatchScalarField + ( + const sorptionWallFunctionFvPatchScalarField&, + const fvPatch&, + const DimensionedField<scalar, volMesh>&, + const fvPatchFieldMapper& + ); + + //- Construct as copy + sorptionWallFunctionFvPatchScalarField + ( + const sorptionWallFunctionFvPatchScalarField& + ); + + //- Construct and return a clone + virtual tmp<fvPatchScalarField> clone() const + { + return tmp<fvPatchScalarField> + ( + new sorptionWallFunctionFvPatchScalarField(*this) + ); + } + + //- Construct as copy setting internal field reference + sorptionWallFunctionFvPatchScalarField + ( + const sorptionWallFunctionFvPatchScalarField&, + const DimensionedField<scalar, volMesh>& + ); + + //- Construct and return a clone setting internal field reference + virtual tmp<fvPatchScalarField> clone + ( + const DimensionedField<scalar, volMesh>& iF + ) const + { + return tmp<fvPatchScalarField> + ( + new sorptionWallFunctionFvPatchScalarField + ( + *this, + iF + ) + ); + } + + + // Member Functions + + // Mapping + + //- Map (and resize as needed) from self given a mapping object + virtual void autoMap(const fvPatchFieldMapper&); + + //- Reverse map the given fvPatchField onto this fvPatchField + virtual void rmap + ( + const fvPatchScalarField&, + const labelList& + ); + + + // Evaluation + + //- Update the coefficients associated with the patch field + virtual void updateCoeffs(); + + + // I-O + + //- Write + virtual void write(Ostream&) const; +}; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* //