Detailed Description

Definition at line 51 of file TESReactionAdaptor.h.

#include <TESReactionAdaptor.h>

Inheritance diagram for ProcessLib::TES::TESFEMReactionAdaptorAdsorption:

Collaboration diagram for ProcessLib::TES::TESFEMReactionAdaptorAdsorption:

Public Member Functions
	TESFEMReactionAdaptorAdsorption (TESLocalAssemblerData const &data)

bool	checkBounds (std::vector< double > const &local_x, std::vector< double > const &local_x_prev_ts) override

ReactionRate	initReaction (const unsigned int_pt) override

void	preZerothTryAssemble () override

double	getReactionDampingFactor () const override

Public Member Functions inherited from ProcessLib::TES::TESFEMReactionAdaptor
virtual	~TESFEMReactionAdaptor ()=default

Private Member Functions
ReactionRate	initReaction_slowDownUndershootStrategy (const unsigned int_pt)

double	estimateAdsorptionEquilibrium (const double p_V0, const double C0) const

Private Attributes
double	_reaction_damping_factor = 1.0

std::vector< bool >	_bounds_violation

TESLocalAssemblerData const &	_d

Additional Inherited Members
Static Public Member Functions inherited from ProcessLib::TES::TESFEMReactionAdaptor
static std::unique_ptr< TESFEMReactionAdaptor >	newInstance (TESLocalAssemblerData const &data)

Constructor & Destructor Documentation

◆ TESFEMReactionAdaptorAdsorption()

ProcessLib::TES::TESFEMReactionAdaptorAdsorption::TESFEMReactionAdaptorAdsorption ( TESLocalAssemblerData const & data )

explicit

Definition at line 51 of file TESReactionAdaptor.cpp.

    : _bounds_violation(data.solid_density.size(), false), _d(data)
{
    assert(dynamic_cast<Adsorption::AdsorptionReaction const*>(
               data.ap.react_sys.get()) != nullptr &&
           "Reactive system has wrong type.");
    assert(!_bounds_violation.empty());
}

References _bounds_violation, ProcessLib::TES::TESLocalAssemblerData::ap, and ProcessLib::TES::AssemblyParams::react_sys.

Member Function Documentation

◆ checkBounds()

bool ProcessLib::TES::TESFEMReactionAdaptorAdsorption::checkBounds	(	std::vector< double > const &	local_x,
		std::vector< double > const &	local_x_prev_ts )

overridevirtual

Reimplemented from ProcessLib::TES::TESFEMReactionAdaptor.

Definition at line 194 of file TESReactionAdaptor.cpp.

{
    double alpha = 1.0;
 
    const double min_xmV = 1e-6;
    const std::size_t nnodes = local_x.size() / NODAL_DOF;
    const std::size_t xmV_offset = COMPONENT_ID_MASS_FRACTION * nnodes;
 
    for (std::size_t i = 0; i < nnodes; ++i)
    {
        auto const xnew = local_x[xmV_offset + i];
        auto const xold = local_x_prev_ts[xmV_offset + i];
 
        if (xnew < min_xmV)
        {
            const auto a = xold / (xold - xnew);
            alpha = std::min(alpha, a);
            _bounds_violation[i] = true;
        }
        else if (xnew > 1.0)
        {
            const auto a = xold / (xnew - xold);
            alpha = std::min(alpha, a);
            _bounds_violation[i] = true;
        }
        else
        {
            _bounds_violation[i] = false;
        }
    }
 
    assert(alpha > 0.0);
 
    if (alpha != 1.0)
    {
        if (alpha > 0.5)
        {
            alpha = 0.5;
        }
        if (alpha < 0.05)
        {
            alpha = 0.05;
        }
 
        if (_d.ap.number_of_try_of_iteration <= 3)
        {
            _reaction_damping_factor *= sqrt(alpha);
        }
        else
        {
            _reaction_damping_factor *= alpha;
        }
    }
 
    return alpha == 1.0;
}

References _bounds_violation, _d, _reaction_damping_factor, ProcessLib::TES::TESLocalAssemblerData::ap, ProcessLib::TES::COMPONENT_ID_MASS_FRACTION, ProcessLib::TES::NODAL_DOF, and ProcessLib::TES::AssemblyParams::number_of_try_of_iteration.

◆ estimateAdsorptionEquilibrium()

double ProcessLib::TES::TESFEMReactionAdaptorAdsorption::estimateAdsorptionEquilibrium	(	const double	p_V0,
		const double	C0 ) const

private

returns estimated equilibrium vapour pressure based on a local (i.e. no diffusion/advection) balance of adsorbate loading and vapour partial pressure

Definition at line 162 of file TESReactionAdaptor.cpp.

{
    auto f = [this, p_V0, C0](double pV) -> double
    {
        // pV0 := _p_V
        const double C_eq =
            _d.ap.react_sys->getEquilibriumLoading(pV, _d.T, _d.ap.M_react);
        return (pV - p_V0) * _d.ap.M_react /
                   MaterialLib::PhysicalConstant::IdealGasConstant / _d.T *
                   _d.ap.poro +
               (1.0 - _d.ap.poro) * (C_eq - C0) * _d.ap.rho_SR_dry;
    };
 
    // range where to search for roots of f
    const double C_eq0 =
        _d.ap.react_sys->getEquilibriumLoading(p_V0, _d.T, _d.ap.M_react);
    const double limit =
        (C_eq0 > C0)
            ? 1e-8
            : Adsorption::AdsorptionReaction::getEquilibriumVapourPressure(
                  _d.T);
 
    // search for roots
    auto rf = MathLib::Nonlinear::makeRegulaFalsi<MathLib::Nonlinear::Pegasus>(
        f, p_V0, limit);
    rf.step(3);
 
    // set vapour pressure
    return rf.getResult();
}

References _d, ProcessLib::TES::TESLocalAssemblerData::ap, Adsorption::AdsorptionReaction::getEquilibriumVapourPressure(), MaterialLib::PhysicalConstant::IdealGasConstant, ProcessLib::TES::AssemblyParams::M_react, MathLib::Nonlinear::makeRegulaFalsi(), ProcessLib::TES::AssemblyParams::poro, ProcessLib::TES::AssemblyParams::react_sys, ProcessLib::TES::AssemblyParams::rho_SR_dry, and ProcessLib::TES::TESLocalAssemblerData::T.

Referenced by initReaction_slowDownUndershootStrategy().

◆ getReactionDampingFactor()

double ProcessLib::TES::TESFEMReactionAdaptorAdsorption::getReactionDampingFactor ( ) const

inlineoverridevirtual

Reimplemented from ProcessLib::TES::TESFEMReactionAdaptor.

Definition at line 67 of file TESReactionAdaptor.h.

    {
        return _reaction_damping_factor;
    }

References _reaction_damping_factor.

◆ initReaction()

ReactionRate ProcessLib::TES::TESFEMReactionAdaptorAdsorption::initReaction ( const unsigned int_pt )

inlineoverridevirtual

Implements ProcessLib::TES::TESFEMReactionAdaptor.

Definition at line 59 of file TESReactionAdaptor.h.

    {
        return initReaction_slowDownUndershootStrategy(int_pt);
    }

References initReaction_slowDownUndershootStrategy().

◆ initReaction_slowDownUndershootStrategy()

ReactionRate ProcessLib::TES::TESFEMReactionAdaptorAdsorption::initReaction_slowDownUndershootStrategy ( const unsigned int_pt )

private

Definition at line 64 of file TESReactionAdaptor.cpp.

{
    assert(_d.ap.number_of_try_of_iteration <= 20);
 
    const double loading = Adsorption::AdsorptionReaction::getLoading(
        _d.solid_density_prev_ts[int_pt], _d.ap.rho_SR_dry);
 
    double react_rate_R =
        _d.ap.react_sys->getReactionRate(_d.p_V, _d.T, _d.ap.M_react, loading) *
        _d.ap.rho_SR_dry;
 
    // set reaction rate based on current damping factor
    react_rate_R = (_reaction_damping_factor > 1e-3)
                       ? _reaction_damping_factor * react_rate_R
                       : 0.0;
 
    if (_d.p_V <
            0.01 * Adsorption::AdsorptionReaction::getEquilibriumVapourPressure(
                       _d.T) &&
        react_rate_R > 0.0)
    {
        react_rate_R = 0.0;
    }
    else if (_d.p_V < 100.0 ||
             _d.p_V < 0.05 * Adsorption::AdsorptionReaction::
                                 getEquilibriumVapourPressure(_d.T))
    {
        // use equilibrium reaction for dry regime
 
        // in the case of zeroth try in zeroth iteration: _p_V and loading are
        // the values
        // at the end of the previous timestep
 
        const double pV_eq = estimateAdsorptionEquilibrium(_d.p_V, loading);
        // TODO [CL]: it would be more correct to subtract pV from the previous
        // timestep here
        const double delta_pV = pV_eq - _d.p_V;
        const double delta_rhoV =
            delta_pV * _d.ap.M_react /
            MaterialLib::PhysicalConstant::IdealGasConstant / _d.T * _d.ap.poro;
        const double delta_rhoSR = delta_rhoV / (_d.ap.poro - 1.0);
        double react_rate_R2 = delta_rhoSR / _d.ap.delta_t;
 
        if (_bounds_violation[int_pt])
        {
            react_rate_R2 *= 0.5;
        }
 
        // 0th try: make sure reaction is not slower than allowed by local
        // estimation
        // nth try: make sure reaction is not made faster by local estimation
        if ((_d.ap.number_of_try_of_iteration == 1 &&
             std::abs(react_rate_R2) > std::abs(react_rate_R)) ||
            (_d.ap.number_of_try_of_iteration > 1 &&
             std::abs(react_rate_R2) < std::abs(react_rate_R)))
        {
            react_rate_R = react_rate_R2;
        }
    }
 
    // smooth out readjustment of reaction rate
    if (_d.ap.iteration_in_current_timestep > 4)
    {
        if (_d.ap.iteration_in_current_timestep <= 9)
        {
            // update reaction rate for for five iterations
            const auto N = _d.ap.iteration_in_current_timestep - 4;
 
            // take average s.t. does not oscillate so much
            react_rate_R = 1.0 / (1.0 + N) *
                           (N * _d.reaction_rate[int_pt] + 1.0 * react_rate_R);
        }
        else
        {
            // afterwards no update anymore
            react_rate_R = _d.reaction_rate[int_pt];
        }
    }
 
    if (_d.ap.number_of_try_of_iteration > 1)
    {
        // assert that within tries reaction does not get faster
        // (e.g. due to switch equilibrium reaction <--> kinetic reaction)
 
        // factor of 0.9*N: in fact, even slow down reaction over tries
        const double r = std::pow(0.9, _d.ap.number_of_try_of_iteration - 1) *
                         _d.reaction_rate[int_pt];
        if (std::abs(react_rate_R) > std::abs(r))
        {
            react_rate_R = r;
        }
    }
 
    return {react_rate_R,
            _d.solid_density_prev_ts[int_pt] + react_rate_R * _d.ap.delta_t};
}

References _bounds_violation, _d, _reaction_damping_factor, ProcessLib::TES::TESLocalAssemblerData::ap, ProcessLib::TES::AssemblyParams::delta_t, estimateAdsorptionEquilibrium(), Adsorption::AdsorptionReaction::getEquilibriumVapourPressure(), Adsorption::AdsorptionReaction::getLoading(), MaterialLib::PhysicalConstant::IdealGasConstant, ProcessLib::TES::AssemblyParams::iteration_in_current_timestep, ProcessLib::TES::AssemblyParams::M_react, ProcessLib::TES::AssemblyParams::number_of_try_of_iteration, ProcessLib::TES::TESLocalAssemblerData::p_V, ProcessLib::TES::AssemblyParams::poro, ProcessLib::TES::AssemblyParams::react_sys, ProcessLib::TES::TESLocalAssemblerData::reaction_rate, ProcessLib::TES::AssemblyParams::rho_SR_dry, ProcessLib::TES::TESLocalAssemblerData::solid_density_prev_ts, and ProcessLib::TES::TESLocalAssemblerData::T.

Referenced by initReaction().

◆ preZerothTryAssemble()

void ProcessLib::TES::TESFEMReactionAdaptorAdsorption::preZerothTryAssemble ( )

overridevirtual

Reimplemented from ProcessLib::TES::TESFEMReactionAdaptor.

Definition at line 253 of file TESReactionAdaptor.cpp.

{
    if (_reaction_damping_factor < 1e-3)
    {
        _reaction_damping_factor = 1e-3;
    }
 
    _reaction_damping_factor = std::min(std::sqrt(_reaction_damping_factor),
                                        10.0 * _reaction_damping_factor);
}