Merge branch 'feature-casting-oxide' into 'develop'

ENH: interfaceOxideRate: new mass transfer model for icoReactingMultiphaseInterFoam solver

See merge request Development/openfoam!510

src/phaseSystemModels/multiphaseInter/phasesSystem/InterfaceCompositionModel/InterfaceCompositionModels.C

@@ -5,7 +5,7 @@
     \\  /    A nd           | www.openfoam.com
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
-    Copyright (C) 2017 OpenCFD Ltd.
+    Copyright (C) 2017-2021 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -51,6 +51,7 @@ License
 #include "kineticGasEvaporation.H"
 #include "Lee.H"
 #include "interfaceHeatResistance.H"
+#include "interfaceOxideRate.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -573,6 +574,78 @@ namespace Foam
             hPowerSolidThermoPhysics
         );
 
+
+        // interfaceOxideRate model definitions
+
+        // From pure phase (tabulated) to solid phase (const)
+        makeInterfacePureType
+        (
+            interfaceOxideRate,
+            heRhoThermo,
+            rhoThermo,
+            pureMixture,
+            tabulatedThermoPhysics,
+            heSolidThermo,
+            solidThermo,
+            pureMixture,
+            hConstSolidThermoPhysics
+        );
+
+        // From pure phase (rho const) to phase (rho const)
+        makeInterfacePureType
+        (
+            interfaceOxideRate,
+            heRhoThermo,
+            rhoThermo,
+            pureMixture,
+            constRhoHThermoPhysics,
+            heRhoThermo,
+            rhoThermo,
+            pureMixture,
+            constRhoHThermoPhysics
+        );
+
+        // From pure phase (ico) to solid phase (const)
+        makeInterfacePureType
+        (
+            interfaceOxideRate,
+            heRhoThermo,
+            rhoThermo,
+            pureMixture,
+            icoPoly8HThermoPhysics,
+            heSolidThermo,
+            solidThermo,
+            pureMixture,
+            hConstSolidThermoPhysics
+        );
+
+        // From pure phase (tabulated) to phase (rho const)
+        makeInterfacePureType
+        (
+            interfaceOxideRate,
+            heRhoThermo,
+            rhoThermo,
+            pureMixture,
+            tabulatedThermoPhysics,
+            heRhoThermo,
+            rhoThermo,
+            pureMixture,
+            constRhoHThermoPhysics
+        );
+
+        // From pure phase (ico) to phase (rho const)
+        makeInterfacePureType
+        (
+            interfaceOxideRate,
+            heRhoThermo,
+            rhoThermo,
+            pureMixture,
+            icoPoly8HThermoPhysics,
+            heRhoThermo,
+            rhoThermo,
+            pureMixture,
+            constRhoHThermoPhysics
+        );
 }
 
 

src/phaseSystemModels/multiphaseInter/phasesSystem/Make/files

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

src/phaseSystemModels/multiphaseInter/phasesSystem/derivedFvPatchFields/timeVaryingMassSorption/timeVaryingMassSorptionFvPatchScalarField.C

@@ -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
+    );
+}
+
+// ************************************************************************* //

src/phaseSystemModels/multiphaseInter/phasesSystem/derivedFvPatchFields/timeVaryingMassSorption/timeVaryingMassSorptionFvPatchScalarField.H

@@ -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
+
+// ************************************************************************* //

src/phaseSystemModels/multiphaseInter/phasesSystem/interfaceCompositionModel/interfaceCompositionModel.C

@@ -131,7 +131,7 @@ const Foam::word& Foam::interfaceCompositionModel::variable() const
 }
 
 
-bool Foam::interfaceCompositionModel::includeDivU()
+bool Foam::interfaceCompositionModel::includeDivU() const noexcept
 {
     return true;
 }

src/phaseSystemModels/multiphaseInter/phasesSystem/interfaceCompositionModel/interfaceCompositionModel.H

@@ -5,7 +5,7 @@
     \\  /    A nd           | www.openfoam.com
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
-    Copyright (C) 2017-2020 OpenCFD Ltd.
+    Copyright (C) 2017-2021 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -195,11 +195,11 @@ public:
         ) = 0;
 
         //- Reference value
-        virtual const dimensionedScalar& Tactivate() const = 0;
+        virtual const dimensionedScalar& Tactivate() const noexcept = 0;
 
         //- Add/subtract alpha*div(U) as a source term
         //- for alpha, substituting div(U) = mDot(1/rho1 - 1/rho2)
-        virtual bool includeDivU();
+        virtual bool includeDivU() const noexcept;
 
         //- Add volume change in pEq
         bool includeVolChange();

src/phaseSystemModels/multiphaseInter/phasesSystem/massTransferModels/Lee/Lee.C

@@ -5,7 +5,7 @@
     \\  /    A nd           | www.openfoam.com
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
-    Copyright (C) 2017-2020 OpenCFD Ltd.
+    Copyright (C) 2017-2021 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -170,20 +170,4 @@ Foam::meltingEvaporationModels::Lee<Thermo, OtherThermo>::KSu
 }
 
 
-template<class Thermo, class OtherThermo>
-const Foam::dimensionedScalar&
-Foam::meltingEvaporationModels::Lee<Thermo, OtherThermo>::Tactivate() const
-{
-    return Tactivate_;
-}
-
-
-template<class Thermo, class OtherThermo>
-bool
-Foam::meltingEvaporationModels::Lee<Thermo, OtherThermo>::includeDivU()
-{
-    return true;
-}
-
-
 // ************************************************************************* //

src/phaseSystemModels/multiphaseInter/phasesSystem/massTransferModels/Lee/Lee.H

@@ -5,7 +5,7 @@
     \\  /    A nd           | www.openfoam.com
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
-    Copyright (C) 2017-2020 OpenCFD Ltd.
+    Copyright (C) 2017-2021 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -165,11 +165,17 @@ public:
         );
 
         //- Return T transition between phases
-        virtual const dimensionedScalar& Tactivate() const;
+        virtual const dimensionedScalar& Tactivate() const noexcept
+        {
+            return Tactivate_;
+        }
 
         //- Add/subtract alpha*div(U) as a source term
         //- for alpha, substituting div(U) = mDot(1/rho1 - 1/rho2)
-        virtual bool includeDivU();
+        virtual bool includeDivU() const noexcept
+        {
+            return true;
+        }
 };
 
 

src/phaseSystemModels/multiphaseInter/phasesSystem/massTransferModels/interfaceHeatResistance/interfaceHeatResistance.C

@@ -5,8 +5,8 @@
     \\  /    A nd           | www.openfoam.com
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
-    Copyright (C) 2020 OpenCFD Ltd.
     Copyright (C) 2020 Henning Scheufler
+    Copyright (C) 2020-2021 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -316,22 +316,4 @@ Foam::meltingEvaporationModels::interfaceHeatResistance<Thermo, OtherThermo>
 }
 
 
-template<class Thermo, class OtherThermo>
-const Foam::dimensionedScalar&
-Foam::meltingEvaporationModels::interfaceHeatResistance<Thermo, OtherThermo>
-::Tactivate() const
-{
-    return Tactivate_;
-}
-
-
-template<class Thermo, class OtherThermo>
-bool
-Foam::meltingEvaporationModels::
-interfaceHeatResistance<Thermo, OtherThermo>::includeDivU()
-{
-    return true;
-}
-
-
 // ************************************************************************* //

src/phaseSystemModels/multiphaseInter/phasesSystem/massTransferModels/interfaceHeatResistance/interfaceHeatResistance.H

@@ -5,8 +5,8 @@
     \\  /    A nd           | www.openfoam.com
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
-    Copyright (C) 2020 OpenCFD Ltd.
     Copyright (C) 2020 Henning Scheufler
+    Copyright (C) 2020-2021 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -171,11 +171,17 @@ public:
         );
 
         //- Return Tactivate
-        virtual const dimensionedScalar& Tactivate() const;
+        virtual const dimensionedScalar& Tactivate() const noexcept
+        {
+            return Tactivate_;
+        }
 
         //- Add/subtract alpha*div(U) as a source term
         //- for alpha, substituting div(U) = mDot(1/rho1 - 1/rho2)
-        virtual bool includeDivU();
+        virtual bool includeDivU() const noexcept
+        {
+            return true;
+        }
 };
 
 

src/phaseSystemModels/multiphaseInter/phasesSystem/massTransferModels/interfaceOxideRate/interfaceOxideRate.C

@@ -0,0 +1,207 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  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 "interfaceOxideRate.H"
+#include "cutCellIso.H"
+#include "volPointInterpolation.H"
+#include "timeVaryingMassSorptionFvPatchScalarField.H"
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+template<class Thermo, class OtherThermo>
+Foam::meltingEvaporationModels::interfaceOxideRate<Thermo, OtherThermo>
+::interfaceOxideRate
+(
+    const dictionary& dict,
+    const phasePair& pair
+)
+:
+    InterfaceCompositionModel<Thermo, OtherThermo>(dict, pair),
+    C_
+    (
+        dimensionedScalar
+        (
+            dimDensity/dimTime,
+            dict.getCheck<scalar>("C", scalarMinMax::ge(0))
+        )
+    ),
+    Tliquidus_
+    (
+        dimensionedScalar
+        (
+            dimTemperature,
+            dict.getCheck<scalar>("Tliquidus", scalarMinMax::ge(0))
+        )
+    ),
+    Tsolidus_
+    (
+        dimensionedScalar
+        (
+            dimTemperature,
+            dict.getCheck<scalar>("Tsolidus", scalarMinMax::ge(0))
+        )
+    ),
+    oxideCrit_
+    (
+        dimensionedScalar
+        (
+            dimDensity,
+            dict.getCheck<scalar>("oxideCrit", scalarMinMax::ge(0))
+        )
+    ),
+    mDotOxide_
+    (
+        IOobject
+        (
+            "mDotOxide",
+            this->mesh_.time().timeName(),
+            this->mesh_,
+            IOobject::NO_READ,
+            IOobject::AUTO_WRITE
+        ),
+        this->mesh_,
+        dimensionedScalar(dimDensity/dimTime, Zero)
+    ),
+    isoAlpha_(dict.getOrDefault<scalar>("isoAlpha", 0.5))
+{}
+
+
+// * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
+
+template<class Thermo, class OtherThermo>
+Foam::tmp<Foam::volScalarField>
+Foam::meltingEvaporationModels::interfaceOxideRate<Thermo, OtherThermo>::Kexp
+(
+    const volScalarField& T
+)
+{
+    const volScalarField& from = this->pair().from();
+    const volScalarField& to = this->pair().to();
+
+    // (CSC:Eq. 2)
+    const fvMesh& mesh = this->mesh_;
+    scalarField ap
+    (
+        volPointInterpolation::New(mesh).interpolate(from)
+    );
+
+    cutCellIso cutCell(mesh, ap);
+
+    tmp<volScalarField> tSalpha = scalar(0)*from;
+    volScalarField& Salpha = tSalpha.ref();
+
+    forAll(Salpha, celli)
+    {
+        const label status = cutCell.calcSubCell(celli, isoAlpha_);
+        if (status == 0) // cell is cut
+        {
+            Salpha[celli] = scalar(1);
+        }
+    }
+
+    // (CSC:Eq. 5)
+    tmp<volScalarField> tSoxide =
+        max((oxideCrit_.value() - to)/oxideCrit_.value(), scalar(0));
+
+    // (CSC:Eq. 4)
+    tmp<volScalarField> tST =
+        Foam::exp
+        (
+            scalar(1)
+          - scalar(1)/max((T - Tsolidus_)/(Tliquidus_ - Tsolidus_),scalar(1e-6))
+        );
+
+    // (CSC:Eq. 6)
+    mDotOxide_ = C_*tSalpha*tSoxide*tST;
+
+    const volScalarField::Boundary& alphab = to.boundaryField();
+
+    forAll(alphab, patchi)
+    {
+        if (isA<timeVaryingMassSorptionFvPatchScalarField>(alphab[patchi]))
+        {
+            const auto& pp =
+                refCast<const timeVaryingMassSorptionFvPatchScalarField>
+                (
+                    alphab[patchi]
+                );
+            const labelUList& fc = mesh.boundary()[patchi].faceCells();
+            tmp<scalarField> tsb = pp.source();
+
+            auto tRhoto = tmp<volScalarField>::New
+            (
+                IOobject
+                (
+                    "tRhoto",
+                    mesh.time().timeName(),
+                    mesh
+                ),
+                mesh,
+                dimensionedScalar(dimDensity, Zero)
+            );
+            volScalarField& rhoto = tRhoto.ref();
+
+            rhoto = this->pair().to().rho();
+
+            forAll(fc, faceI)
+            {
+                const label cellI = fc[faceI];
+                const scalar rhoI = rhoto[cellI];
+                mDotOxide_[cellI] += rhoI*tsb()[faceI];
+            }
+        }
+    }
+
+    return tmp<volScalarField>::New(mDotOxide_);
+}
+
+
+template<class Thermo, class OtherThermo>
+Foam::tmp<Foam::volScalarField>
+Foam::meltingEvaporationModels::interfaceOxideRate<Thermo, OtherThermo>::KSp
+(
+    label variable,
+    const volScalarField& refValue
+)
+{
+    return nullptr;
+}
+
+
+template<class Thermo, class OtherThermo>
+Foam::tmp<Foam::volScalarField>
+Foam::meltingEvaporationModels::interfaceOxideRate<Thermo, OtherThermo>::KSu
+(
+    label variable,
+    const volScalarField& refValue
+)
+{
+    return nullptr;
+}
+
+
+// ************************************************************************* //

src/phaseSystemModels/multiphaseInter/phasesSystem/massTransferModels/interfaceOxideRate/interfaceOxideRate.H

@@ -0,0 +1,238 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  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::meltingEvaporationModels::interfaceOxideRate
+
+Description
+    The \c interfaceOxideRate is a simple model to calculate the formation
+    rate of oxide inclusions (\c mDotOxide), which is a local function of
+    the volume fraction of reducing agent (\c alpha), temperature (\c T)
+    and oxide-inclusion density (\c chi).
+
+    The oxide-inclusion formation rate is modelled as follows:
+
+    \f[
+        S_\chi = C_\chi S_\alpha S_T S_\rho
+    \f]
+
+    with
+
+    \f[
+        S_\alpha = \alpha (1 - \alpha)
+    \f]
+
+    \f[
+        S_T =
+        \exp{
+            1 - \frac{1}{max(T - T_{solidus}/(T_liquidus - T_solidus), 1e-6)}
+        }
+    \f]
+
+    \f[
+        S_\rho =
+            max \left(\frac{\chi_{crit} - \chi_{curr}}{\chi_{curr}}, 0\right)
+    \f]
+
+    where
+    \vartable
+      S_\chi   | Oxide-inclusion formation rate          [kg/(m^3 s)]
+      C_\chi   | Oxide-inclusion formation rate constant [kg/(m^3 s)]
+      S_\alpha | Formation factor due to volume fraction of reducing agent [-]
+      S_T      | Formation factor due to temperature     [-]
+      S_\rho   | Formation factor due to oxide-inclusion density [-]
+      \alpha   | Volume fraction of reducing agent       [-]
+      T        | Local temperature                       [K]
+      T_{solidus}  | Solidus temperature of reducing agent  [K]
+      T_{liquidus} | Liquidus temperature of reducing agent [K]
+      \chi_{crit}  | Critical oxide-inclusion density    [kg/m^3]
+      \chi_{curr}  | Current oxide-inclusion density     [kg/m^3]
+    \endvartable
+
+    References:
+    \verbatim
+        Oxide-inclusion model (tag:CSC):
+            Cao, L., Sun, F., Chen, T., Tang, Y., & Liao, D. (2018).
+            Quantitative prediction of oxide inclusion defects inside
+            the casting and on the walls during cast-filling processes.
+            International Journal of Heat and Mass Transfer, 119, 614-623.
+            DOI:10.1016/j.ijheatmasstransfer.2017.11.127
+    \endverbatim
+
+Usage
+    Minimal example by using \c constant/phaseProperties.massTransferModel:
+    \verbatim
+    massTransferModel
+    (
+        (liquid to oxide)
+        {
+            type                interfaceOxideRate;
+            C                   <scalar>;
+            Tliquidus           <scalar>;
+            Tsolidus            <scalar>;
+            oxideCrit           <scalar>;
+            isoAlpha            <scalar>;
+        }
+    );
+    \endverbatim
+
+    where the entries mean:
+    \table
+      Property  | Description                             | Type | Reqd | Deflt
+      type      | Type name: interfaceOxideRate           | word   | yes | -
+      C         | Oxide-inclusion formation rate constant | scalar | yes | -
+      Tliquidus | Liquidus temperature of reducing agent  | scalar | yes | -
+      Tsolidus  | Solidus temperature of reducing agent   | scalar | yes | -
+      oxideCrit | Critical oxide-inclusion density        | scalar | yes | -
+      isoAlpha  | Location of the source                  | scalar | no  | 0.5
+    \endtable
+
+Note
+  - \c oxideCrit should be determined experimentally (CSC:p. 616).
+  - \c C should be determined by practical production (CSC:p. 616).
+
+SourceFiles
+    interfaceOxideRate.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef meltingEvaporationModels_interfaceOxideRate_H
+#define meltingEvaporationModels_interfaceOxideRate_H
+
+#include "InterfaceCompositionModel.H"
+
+// * * * * * * * * * * * * * * * * *  * * * * * * * * * * * * * * * * * * * *//
+
+namespace Foam
+{
+
+// Forward Declarations
+class phasePair;
+
+namespace meltingEvaporationModels
+{
+
+/*---------------------------------------------------------------------------*\
+                    Class interfaceOxideRate Declaration
+\*---------------------------------------------------------------------------*/
+
+template<class Thermo, class OtherThermo>
+class interfaceOxideRate
+:
+    public InterfaceCompositionModel<Thermo, OtherThermo>
+{
+    // Private Data
+
+        //- Oxide-inclusion formation rate constant
+        const dimensionedScalar C_;
+
+        //- Liquidus temperature of reducing agent (e.g. a casting metal)
+        const dimensionedScalar Tliquidus_;
+
+        //- Solidus temperature of reducing agent (e.g. a casting metal)
+        const dimensionedScalar Tsolidus_;
+
+        //- Critical oxide-inclusion density
+        const dimensionedScalar oxideCrit_;
+
+        //- Oxide-inclusion formation rate
+        volScalarField mDotOxide_;
+
+        //- Interface Iso-value
+        scalar isoAlpha_;
+
+
+public:
+
+    //- Runtime type information
+    TypeName("interfaceOxideRate");
+
+
+    // Constructors
+
+        //- Construct from components
+        interfaceOxideRate
+        (
+            const dictionary& dict,
+            const phasePair& pair
+        );
+
+
+    //- Destructor
+    virtual ~interfaceOxideRate() = default;
+
+
+    // Member Functions
+
+        //- Explicit total mass transfer coefficient
+        virtual tmp<volScalarField> Kexp
+        (
+            const volScalarField& field
+        );
+
+        //- Implicit mass transfer coefficient
+        virtual tmp<volScalarField> KSp
+        (
+            label modelVariable,
+            const volScalarField& field
+        );
+
+        //- Explicit mass transfer coefficient
+        virtual tmp<volScalarField> KSu
+        (
+            label modelVariable,
+            const volScalarField& field
+        );
+
+        //- Return solidus temperature of reducing agent
+        virtual const dimensionedScalar& Tactivate() const noexcept
+        {
+            return Tsolidus_;
+        }
+
+        //- Add/subtract alpha*div(U) as a source term
+        //- for alpha, substituting div(U) = mDot(1/rho1 - 1/rho2)
+        virtual bool includeDivU() const noexcept
+        {
+            return true;
+        }
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace meltingEvaporationModels
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#ifdef NoRepository
+#   include "interfaceOxideRate.C"
+#endif
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+#endif
+
+// ************************************************************************* //

src/phaseSystemModels/multiphaseInter/phasesSystem/massTransferModels/kineticGasEvaporation/kineticGasEvaporation.C

@@ -5,7 +5,7 @@
     \\  /    A nd           | www.openfoam.com
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
-    Copyright (C) 2017-2020 OpenCFD Ltd.
+    Copyright (C) 2017-2021 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -319,21 +319,4 @@ Foam::meltingEvaporationModels::kineticGasEvaporation<Thermo, OtherThermo>::KSu
 }
 
 
-template<class Thermo, class OtherThermo>
-const Foam::dimensionedScalar&
-Foam::meltingEvaporationModels::kineticGasEvaporation<Thermo, OtherThermo>
-::Tactivate() const
-{
-    return Tactivate_;
-}
-
-
-template<class Thermo, class OtherThermo>
-bool
-Foam::meltingEvaporationModels::
-kineticGasEvaporation<Thermo, OtherThermo>::includeDivU()
-{
-    return true;
-}
-
 // ************************************************************************* //

src/phaseSystemModels/multiphaseInter/phasesSystem/massTransferModels/kineticGasEvaporation/kineticGasEvaporation.H

@@ -5,7 +5,7 @@
     \\  /    A nd           | www.openfoam.com
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
-    Copyright (C) 2017-2020 OpenCFD Ltd.
+    Copyright (C) 2017-2021 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -210,11 +210,17 @@ public:
         );
 
         //- Return Tactivate
-        virtual const dimensionedScalar& Tactivate() const;
+        virtual const dimensionedScalar& Tactivate() const noexcept
+        {
+            return Tactivate_;
+        }
 
         //- Add/subtract alpha*div(U) as a source term
         //- for alpha, substituting div(U) = mDot(1/rho1 - 1/rho2)
-        virtual bool includeDivU();
+        virtual bool includeDivU() const noexcept
+        {
+            return true;
+        }
 };
 
 

src/phaseSystemModels/multiphaseInter/phasesSystem/multiphaseSystem/multiphaseSystem.C

@@ -358,6 +358,8 @@ void Foam::multiphaseSystem::solveAlphas()
                  Su + fvm::Sp(Sp, alpha1)
             );
 
+            alpha1Eqn.boundaryManipulate(alpha1.boundaryFieldRef());
+
             alpha1Eqn.solve();
 
             phiAlpha += alpha1Eqn.flux();

tutorials/multiphase/icoReactingMultiphaseInterFoam/oxideFormation/0.orig/T

@@ -0,0 +1,43 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  v2106                                 |
+|   \\  /    A nd           | Website:  www.openfoam.com                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       volScalarField;
+    object      T;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+dimensions      [0 0 0 1 0 0 0];
+
+internalField   uniform 1013.15;
+
+boundaryField
+{
+    bottomInlet
+    {
+        type            fixedValue;
+        value           $internalField;
+    }
+
+    topOutlet
+    {
+        type            inletOutlet;
+        value           $internalField;
+        inletValue      $internalField;
+    }
+
+    "(left|right|front|back|bottomLeft|bottomRight)"
+    {
+        type            zeroGradient;
+    }
+}
+
+
+// ************************************************************************* //

tutorials/multiphase/icoReactingMultiphaseInterFoam/oxideFormation/0.orig/U

@@ -0,0 +1,43 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  v2106                                 |
+|   \\  /    A nd           | Website:  www.openfoam.com                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       volVectorField;
+    object      U;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+dimensions      [0 1 -1 0 0 0 0];
+
+internalField   uniform (0 0 0);
+
+boundaryField
+{
+    bottomInlet
+    {
+        type            fixedValue;
+        value           uniform (0 0 0.1);
+    }
+
+    topOutlet
+    {
+        type            pressureInletOutletVelocity;
+        inletValue      $internalField;
+        value           $internalField;
+    }
+
+    "(left|right|front|back|bottomLeft|bottomRight)"
+    {
+        type            noSlip;
+    }
+}
+
+
+// ************************************************************************* //

tutorials/multiphase/icoReactingMultiphaseInterFoam/oxideFormation/0.orig/alpha.gas

@@ -0,0 +1,43 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  v2106                                 |
+|   \\  /    A nd           | Website:  www.openfoam.com                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       volScalarField;
+    object      alpha.gas;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+dimensions      [0 0 0 0 0 0 0];
+
+internalField   uniform 1;
+
+boundaryField
+{
+    bottomInlet
+    {
+        type            fixedValue;
+        value           uniform 0;
+    }
+
+    topOutlet
+    {
+        type            inletOutlet;
+        inletValue      $internalField;
+        value           $internalField;
+    }
+
+    "(left|right|front|back|bottomLeft|bottomRight)"
+    {
+        type            zeroGradient;
+    }
+}
+
+
+// ************************************************************************* //

tutorials/multiphase/icoReactingMultiphaseInterFoam/oxideFormation/0.orig/alpha.liquid

@@ -0,0 +1,43 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  v2106                                 |
+|   \\  /    A nd           | Website:  www.openfoam.com                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       volScalarField;
+    object      alpha.liquid;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+dimensions      [0 0 0 0 0 0 0];
+
+internalField   uniform 0;
+
+boundaryField
+{
+    bottomInlet
+    {
+        type            fixedValue;
+        value           uniform 1;
+    }
+
+    topOutlet
+    {
+        type            inletOutlet;
+        inletValue      $internalField;
+        value           $internalField;
+    }
+
+    "(left|right|front|back|bottomLeft|bottomRight)"
+    {
+        type            zeroGradient;
+    }
+}
+
+
+// ************************************************************************* //

tutorials/multiphase/icoReactingMultiphaseInterFoam/oxideFormation/0.orig/alpha.oxide

@@ -0,0 +1,46 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  v2106                                 |
+|   \\  /    A nd           | Website:  www.openfoam.com                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       volScalarField;
+    object      alpha.liquid;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+dimensions      [0 0 0 0 0 0 0];
+
+internalField   uniform 0;
+
+boundaryField
+{
+    bottomInlet
+    {
+        type                fixedValue;
+        value               $internalField;
+    }
+
+    topOutlet
+    {
+        type                inletOutlet;
+        inletValue          $internalField;
+        value               $internalField;
+    }
+
+    "(left|right|front|back|bottomLeft|bottomRight)"
+    {
+        type            timeVaryingMassSorption;
+        kabs            1e1;
+        max             1e-3;
+        value           $internalField;
+    }
+}
+
+
+// ************************************************************************* //