LiquidEvaporation.C

Go to the documentation of this file.

/*---------------------------------------------------------------------------*\
  =========                 |
  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   \\    /   O peration     | Website:  https://openfoam.org
    \\  /    A nd           | Copyright (C) 2011-2022 OpenFOAM Foundation
     \\/     M anipulation  |
-------------------------------------------------------------------------------
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 "LiquidEvaporation.H"
#include "specie.H"
#include "mathematicalConstants.H"
using namespace Foam::constant::mathematical;
 
// * * * * * * * * * * * * Protected Member Functions  * * * * * * * * * * * //
 
template<class CloudType>
Foam::tmp<Foam::scalarField> Foam::LiquidEvaporation<CloudType>::calcXc
(
    const label celli
) const
{
    scalarField Xc(this->owner().composition().carrier().Y().size());
 
    forAll(Xc, i)
    {
        Xc[i] =
            this->owner().composition().carrier().Y()[i][celli]
           /this->owner().composition().carrier().Wi(i);
    }
 
    return Xc/sum(Xc);
}
 
 
template<class CloudType>
Foam::scalar Foam::LiquidEvaporation<CloudType>::Sh
(
    const scalar Re,
    const scalar Sc
) const
{
    return 2.0 + 0.6*Foam::sqrt(Re)*cbrt(Sc);
}
 
 
// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
template<class CloudType>
Foam::LiquidEvaporation<CloudType>::LiquidEvaporation
(
    const dictionary& dict,
    CloudType& owner
)
:
    PhaseChangeModel<CloudType>(dict, owner, typeName),
    liquids_(owner.thermo().liquids()),
    condensation_
    (
        this->coeffDict().template lookupOrDefault<Switch>
        (
            "condensation",
            false
        )
    ),
    activeLiquids_(this->coeffDict().lookup("activeLiquids")),
    liqToCarrierMap_(activeLiquids_.size(), -1),
    liqToLiqMap_(activeLiquids_.size(), -1)
{
    if (activeLiquids_.size() == 0)
    {
        WarningInFunction
            << "Evaporation model selected, but no active liquids defined"
            << nl << endl;
    }
    else
    {
        Info<< "Participating liquid species:" << endl;
 
        // Determine mapping between liquid and carrier phase species
        forAll(activeLiquids_, i)
        {
            Info<< "    " << activeLiquids_[i] << endl;
            liqToCarrierMap_[i] =
                owner.composition().carrierId(activeLiquids_[i]);
        }
 
        // Determine mapping between model active liquids and global liquids
        const label idLiquid = owner.composition().idLiquid();
        forAll(activeLiquids_, i)
        {
            liqToLiqMap_[i] =
                owner.composition().localId(idLiquid, activeLiquids_[i]);
        }
    }
}
 
 
template<class CloudType>
Foam::LiquidEvaporation<CloudType>::LiquidEvaporation
(
    const LiquidEvaporation<CloudType>& pcm
)
:
    PhaseChangeModel<CloudType>(pcm),
    liquids_(pcm.owner().thermo().liquids()),
    activeLiquids_(pcm.activeLiquids_),
    liqToCarrierMap_(pcm.liqToCarrierMap_),
    liqToLiqMap_(pcm.liqToLiqMap_)
{}
 
 
// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
 
template<class CloudType>
Foam::LiquidEvaporation<CloudType>::~LiquidEvaporation()
{}
 
 
// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
template<class CloudType>
void Foam::LiquidEvaporation<CloudType>::calculate
(
    const typename CloudType::parcelType& p,
    const typename CloudType::parcelType::trackingData& td,
    const scalar dt,
    const scalar Re,
    const scalar Pr,
    const scalar d,
    const scalar nu,
    const scalar T,
    const scalar Ts,
    const scalar pc,
    const scalar Tc,
    const scalarField& X,
    scalarField& dMassPC
) const
{
    // immediately evaporate mass that has reached critical condition
    if ((liquids_.Tc(X) - T) < small)
    {
        if (debug)
        {
            WarningInFunction
                << "Parcel reached critical conditions: "
                << "evaporating all available mass" << endl;
        }
 
        forAll(activeLiquids_, i)
        {
            const label lid = liqToLiqMap_[i];
            dMassPC[lid] = great;
        }
 
        return;
    }
 
    // construct carrier phase species volume fractions for cell, celli
    const scalarField Xc(calcXc(p.cell()));
 
    // calculate mass transfer of each specie in liquid
    forAll(activeLiquids_, i)
    {
        const label gid = liqToCarrierMap_[i];
        const label lid = liqToLiqMap_[i];
 
        // vapour diffusivity [m^2/s]
        const scalar Dab = liquids_.properties()[lid].D(pc, Ts);
 
        // saturation pressure for species i [pa]
        // - carrier phase pressure assumed equal to the liquid vapour pressure
        //   close to the surface
        // NOTE: if pSat > pc then particle is superheated
        // calculated evaporation rate will be greater than that of a particle
        // at boiling point, but this is not a boiling model
        const scalar pSat = liquids_.properties()[lid].pv(pc, T);
 
        // Schmidt number
        const scalar Sc = nu/(Dab + rootVSmall);
 
        // Sherwood number
        const scalar Sh = this->Sh(Re, Sc);
 
        // mass transfer coefficient [m/s]
        const scalar kc = Sh*Dab/(d + rootVSmall);
 
        // vapour concentration at surface [kmol/m^3] at film temperature
        const scalar Cs = X[lid]*pSat/(RR*Ts);
 
        // vapour concentration in bulk gas [kmol/m^3] at film temperature
        const scalar Cinf = Xc[gid]*pc/(RR*Ts);
 
        // molar flux of vapour [kmol/m2/s]
        scalar Ni = kc*(Cs - Cinf);
 
        // limit if not permitting condensation
        if (!condensation_)
        {
            Ni = max(Ni, 0);
        }
 
        // mass transfer [kg]
        dMassPC[lid] += Ni*pi*sqr(d)*liquids_.properties()[lid].W()*dt;
    }
}
 
 
template<class CloudType>
Foam::scalar Foam::LiquidEvaporation<CloudType>::dh
(
    const label idc,
    const label idl,
    const scalar p,
    const scalar T
) const
{
    scalar dh = 0;
 
    typedef PhaseChangeModel<CloudType> parent;
    switch (parent::enthalpyTransfer_)
    {
        case (parent::etLatentHeat):
        {
            dh = liquids_.properties()[idl].hl(p, T);
            break;
        }
        case (parent::etEnthalpyDifference):
        {
            scalar hc = this->owner().composition().carrier().Ha(idc, p, T);
            scalar hp = liquids_.properties()[idl].Ha(p, T);
 
            dh = hc - hp;
            break;
        }
        default:
        {
            FatalErrorInFunction
                << "Unknown enthalpyTransfer type" << abort(FatalError);
        }
    }
 
    return dh;
}
 
 
template<class CloudType>
Foam::scalar Foam::LiquidEvaporation<CloudType>::Tvap
(
    const scalarField& X
) const
{
    return liquids_.Tpt(X);
}
 
 
template<class CloudType>
Foam::scalar Foam::LiquidEvaporation<CloudType>::TMax
(
    const scalar p,
    const scalarField& X
) const
{
    return liquids_.pvInvert(p, X);
}
 
 
// ************************************************************************* //