ENH: timeVaryingMassSorption: new boundary condition for icoReacting solvers

Co-authored-by: Kutalmis Bercin <kutalmis.bercin@esi-group.com>
This commit is contained in:
sergio
2021-12-09 11:31:08 +00:00
committed by Andrew Heather
parent ea9d424e24
commit 437ea6917a
3 changed files with 533 additions and 0 deletions

View File

@ -25,4 +25,6 @@ $(surfaceTension)/surfaceTensionModel/surfaceTensionModel.C
$(surfaceTension)/constantSurfaceTensionCoefficient/constantSurfaceTensionCoefficient.C $(surfaceTension)/constantSurfaceTensionCoefficient/constantSurfaceTensionCoefficient.C
*/ */
derivedFvPatchFields/timeVaryingMassSorption/timeVaryingMassSorptionFvPatchScalarField.C
LIB = $(FOAM_LIBBIN)/libincompressibleMultiphaseSystems LIB = $(FOAM_LIBBIN)/libincompressibleMultiphaseSystems

View File

@ -0,0 +1,282 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | www.openfoam.com
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2021 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 "timeVaryingMassSorptionFvPatchScalarField.H"
#include "addToRunTimeSelectionTable.H"
#include "fvPatchFieldMapper.H"
#include "volFields.H"
#include "EulerDdtScheme.H"
#include "CrankNicolsonDdtScheme.H"
#include "backwardDdtScheme.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
const Foam::Enum
<
Foam::timeVaryingMassSorptionFvPatchScalarField::ddtSchemeType
>
Foam::timeVaryingMassSorptionFvPatchScalarField::ddtSchemeTypeNames_
({
{
ddtSchemeType::tsEuler,
fv::EulerDdtScheme<scalar>::typeName_()
},
{
ddtSchemeType::tsCrankNicolson,
fv::CrankNicolsonDdtScheme<scalar>::typeName_()
},
{
ddtSchemeType::tsBackward,
fv::backwardDdtScheme<scalar>::typeName_()
},
});
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::timeVaryingMassSorptionFvPatchScalarField::
timeVaryingMassSorptionFvPatchScalarField
(
const fvPatch& p,
const DimensionedField<scalar, volMesh>& iF
)
:
fixedValueFvPatchScalarField(p, iF),
kabs_(scalar(1)),
max_(scalar(1)),
kdes_(scalar(1))
{}
Foam::timeVaryingMassSorptionFvPatchScalarField::
timeVaryingMassSorptionFvPatchScalarField
(
const fvPatch& p,
const DimensionedField<scalar, volMesh>& iF,
const dictionary& dict
)
:
fixedValueFvPatchScalarField(p, iF, dict, false),
kabs_(dict.getCheck<scalar>("kabs", scalarMinMax::ge(0))),
max_(dict.getCheck<scalar>("max", scalarMinMax::ge(0))),
kdes_(dict.getCheckOrDefault<scalar>("kdes", 0, scalarMinMax::ge(0)))
{
if (dict.found("value"))
{
fvPatchScalarField::operator=
(
scalarField("value", dict, p.size())
);
}
else
{
fvPatchField<scalar>::operator=(Zero);
}
}
Foam::timeVaryingMassSorptionFvPatchScalarField::
timeVaryingMassSorptionFvPatchScalarField
(
const timeVaryingMassSorptionFvPatchScalarField& ptf,
const fvPatch& p,
const DimensionedField<scalar, volMesh>& iF,
const fvPatchFieldMapper& mapper
)
:
fixedValueFvPatchScalarField(ptf, p, iF, mapper),
kabs_(ptf.kabs_),
max_(ptf.max_),
kdes_(ptf.kdes_)
{}
Foam::timeVaryingMassSorptionFvPatchScalarField::
timeVaryingMassSorptionFvPatchScalarField
(
const timeVaryingMassSorptionFvPatchScalarField& ptf
)
:
fixedValueFvPatchScalarField(ptf),
kabs_(ptf.kabs_),
max_(ptf.max_),
kdes_(ptf.kdes_)
{}
Foam::timeVaryingMassSorptionFvPatchScalarField::
timeVaryingMassSorptionFvPatchScalarField
(
const timeVaryingMassSorptionFvPatchScalarField& ptf,
const DimensionedField<scalar, volMesh>& iF
)
:
fixedValueFvPatchScalarField(ptf, iF),
kabs_(ptf.kabs_),
max_(ptf.max_),
kdes_(ptf.kdes_)
{}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
void Foam::timeVaryingMassSorptionFvPatchScalarField::autoMap
(
const fvPatchFieldMapper& m
)
{
fixedValueFvPatchScalarField::autoMap(m);
}
void Foam::timeVaryingMassSorptionFvPatchScalarField::rmap
(
const fvPatchScalarField& ptf,
const labelList& addr
)
{
fixedValueFvPatchScalarField::rmap(ptf, addr);
}
Foam::tmp<Foam::scalarField>
Foam::timeVaryingMassSorptionFvPatchScalarField::source() const
{
auto tsource = tmp<scalarField>::New(patch().size(), Zero);
auto& source = tsource.ref();
const scalarField cp(*this);
const scalarField w(max(1 - cp/max_, scalar(0)));
source = -kabs_*w*max(patchInternalField() - cp, scalar(0));
source += kdes_*max(cp - patchInternalField(), scalar(0));
return tsource;
}
void Foam::timeVaryingMassSorptionFvPatchScalarField::updateCoeffs()
{
if (updated())
{
return;
}
const label patchi = patch().index();
const scalar dt = db().time().deltaTValue();
const auto& fld =
db().lookupObject<volScalarField>(this->internalField().name());
const volScalarField& fld0 = fld.oldTime();
// Lookup d/dt scheme from database
const word ddtSchemeName(fld.mesh().ddtScheme(fld.name()));
const ddtSchemeType& ddtScheme = ddtSchemeTypeNames_[ddtSchemeName];
const scalarField cp(*this);
const scalarField w(max(1 - cp/max_, scalar(0)));
tmp<scalarField> dfldp =
kabs_
*w
*max(patchInternalField() - cp, scalar(0))
*dt;
dfldp.ref() -=
kdes_
*max(cp - patchInternalField(), scalar(0))
*dt;
switch (ddtScheme)
{
case tsEuler:
case tsCrankNicolson:
{
operator==(fld0.boundaryField()[patchi] + dfldp);
break;
}
case tsBackward:
{
const scalar dt0 = db().time().deltaT0Value();
const scalar c = scalar(1) + dt/(dt + dt0);
const scalar c00 = dt*dt/(dt0*(dt + dt0));
const scalar c0 = c + c00;
operator==
(
(
c0*fld0.boundaryField()[patchi]
- c00*fld0.oldTime().boundaryField()[patchi]
+ dfldp
)/c
);
break;
}
default:
{
FatalErrorInFunction
<< ddtSchemeName << nl
<< " on patch " << this->patch().name()
<< " of field " << this->internalField().name()
<< " in file " << this->internalField().objectPath()
<< exit(FatalError);
}
}
fixedValueFvPatchScalarField::updateCoeffs();
}
void Foam::timeVaryingMassSorptionFvPatchScalarField::write(Ostream& os) const
{
fvPatchScalarField::write(os);
os.writeEntry("kabs", kabs_);
os.writeEntry("max", max_);
os.writeEntryIfDifferent<scalar>("kdes", scalar(0), kdes_);
writeEntry("value", os);
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
makePatchTypeField
(
fvPatchScalarField,
timeVaryingMassSorptionFvPatchScalarField
);
}
// ************************************************************************* //

View File

@ -0,0 +1,249 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | www.openfoam.com
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2021 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::timeVaryingMassSorptionFvPatchScalarField
Group
grpGenericBoundaryConditions
Description
This boundary condition provides a first order fixed-value condition
for a given scalar field to model time-dependent adsorption-desoprtion
processes to be used with the \c interfaceOxideRate mass model
\f[
\frac{d c}{d t} =
k_{abs} w (c_{int} - c_{p_{w}}) + k_{des} (c_{p_{w} - c_{int}})
\f]
\f[
w = max(1 - c_{p_{w}/max, 0);
\f]
where
\vartable
c_{int} | Concentration at cell
c_{p_{w}} | Concentration at wall
k_{abs} | Adsorption rate constant [1/s]
k_{des} | Desorption rate constant [1/s]
w | Weight function
max | Max concentration at wall
\endvartable
Usage
Example of the boundary condition specification:
\verbatim
<patchName>
{
// Mandatory entries
type timeVaryingMassSorption;
kbas <scalar>;
max <scalar>;
// Optional entries
kdes <scalar>;
// Inherited entries
...
}
\endverbatim
where the entries mean:
\table
Property | Description | Type | Reqd | Deflt
type | Type name: timeVaryingAdsorption | word | yes | -
kbas | Adsorption rate constant | scalar | yes | -
max | Maximum concentation at wall | scalar | yes | -
kdes | Desorption rate constant | scalar | no | 0
\endtable
The inherited entries are elaborated in:
- \link fixedValueFvPatchFields.H \endlink
SourceFiles
timeVaryingMassSorptionFvPatchScalarField.C
\*---------------------------------------------------------------------------*/
#ifndef timeVaryingMassSorptionFvPatchScalarField_H
#define timeVaryingMassSorptionFvPatchScalarField_H
#include "fixedValueFvPatchFields.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
/*---------------------------------------------------------------------------*\
Class timeVaryingMassSorptionFvPatchScalarField Declaration
\*---------------------------------------------------------------------------*/
class timeVaryingMassSorptionFvPatchScalarField
:
public fixedValueFvPatchScalarField
{
public:
// Public Enumeration
//- Enumeration defining the available ddt schemes
enum ddtSchemeType
{
tsEuler,
tsCrankNicolson,
tsBackward
};
private:
// Private Data
//- Time scheme type names
static const Enum<ddtSchemeType> ddtSchemeTypeNames_;
//- Adsorption rate constant
scalar kabs_;
//- Maximum level of adsorption of a given substance on patch
scalar max_;
//- Desorption rate constant
scalar kdes_;
public:
//- Runtime type information
TypeName("timeVaryingMassSorption");
// Constructors
//- Construct from patch and internal field
timeVaryingMassSorptionFvPatchScalarField
(
const fvPatch&,
const DimensionedField<scalar, volMesh>&
);
//- Construct from patch, internal field and dictionary
timeVaryingMassSorptionFvPatchScalarField
(
const fvPatch&,
const DimensionedField<scalar, volMesh>&,
const dictionary&
);
//- Construct by mapping given
//- timeVaryingMassSorptionFvPatchScalarField onto a new patch
timeVaryingMassSorptionFvPatchScalarField
(
const timeVaryingMassSorptionFvPatchScalarField&,
const fvPatch&,
const DimensionedField<scalar, volMesh>&,
const fvPatchFieldMapper&
);
//- Construct as copy
timeVaryingMassSorptionFvPatchScalarField
(
const timeVaryingMassSorptionFvPatchScalarField&
);
//- Construct and return a clone
virtual tmp<fvPatchScalarField> clone() const
{
return tmp<fvPatchScalarField>
(
new timeVaryingMassSorptionFvPatchScalarField(*this)
);
}
//- Construct as copy setting internal field reference
timeVaryingMassSorptionFvPatchScalarField
(
const timeVaryingMassSorptionFvPatchScalarField&,
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 timeVaryingMassSorptionFvPatchScalarField(*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&
);
// Help
//- Return source rate
tmp<scalarField> source() const;
// Evaluation
//- Update the coefficients associated with the patch field
virtual void updateCoeffs();
//- Write
virtual void write(Ostream&) const;
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //