Commit ae450cb3 authored by Katharina Heck's avatar Katharina Heck Committed by Timo Koch
Browse files

[feature] add spatialparams with base stuff for nonequilibrium

parent 11fe905e
// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
// vi: set et ts=4 sw=4 sts=4:
/*****************************************************************************
* See the file COPYING for full copying permissions. *
* *
* This program 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 2 of the License, or *
* (at your option) any later version. *
* *
* This program 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 this program. If not, see <http://www.gnu.org/licenses/>. *
*****************************************************************************/
/*!
* \file
* \ingroup SpatialParameters
* \brief base class for spatialparameters dealing with thermal and chemical nonequilibrium
*
*/
#ifndef DUMUX_FV_SPATIALPARAMS_NONEQUILIBRIUM_HH
#define DUMUX_FV_SPATIALPARAMS_NONEQUILIBRIUM_HH
#include <dumux/material/spatialparams/fv.hh>
// material laws for interfacial area
#include <dumux/material/fluidmatrixinteractions/2pia/efftoabslawia.hh>
#include <dumux/material/fluidmatrixinteractions/2pia/awnsurfacepolynomial2ndorder.hh>
#include <dumux/material/fluidmatrixinteractions/2pia/awnsurfacepcmaxfct.hh>
#include <dumux/material/fluidmatrixinteractions/2pia/awnsurfaceexpswpcto3.hh>
namespace Dumux {
/**
* \brief Definition of the spatial parameters for nonequilibrium
*/
template<class FVGridGeometry, class Scalar, class Implementation>
class FVNonEquilibriumSpatialParams
: public FVSpatialParams<FVGridGeometry,
Scalar,
Implementation>
{
using ParentType = FVSpatialParams<FVGridGeometry, Scalar, Implementation>;
using GridView = typename FVGridGeometry::GridView;
using FVElementGeometry = typename FVGridGeometry::LocalView;
using SubControlVolume = typename FVGridGeometry::SubControlVolume;
using Element = typename GridView::template Codim<0>::Entity;
using GlobalPosition = typename Element::Geometry::GlobalCoordinate;
public:
//! export the types used for interfacial area calculations
using AwnSurfaceParams = Scalar;
using AwsSurfaceParams = Scalar;
using AnsSurfaceParams = Scalar;
FVNonEquilibriumSpatialParams(std::shared_ptr<const FVGridGeometry> fvGridGeometry)
: ParentType(fvGridGeometry)
{ }
/*!\brief Return a reference to the container object for the
* parametrization of the surface between wetting and non-Wetting phase.
*
* The position is determined based on the coordinate of
* the vertex belonging to the considered sub control volume.
* \param element The finite element
* \param fvGeometry The finite volume geometry
* \param scvIdx The local index of the sub-control volume */
template<class ElementSolution>
AwnSurfaceParams aWettingNonWettingSurfaceParams(const Element &element,
const SubControlVolume &scv,
const ElementSolution &elemSol) const
{
DUNE_THROW(Dune::InvalidStateException,
"The spatial parameters do not provide a aWettingNonWettingSurfaceParams() method.");
}
/*!\brief Return a reference to the container object for the
* parametrization of the surface between non-Wetting and solid phase.
*
* The position is determined based on the coordinate of
* the vertex belonging to the considered sub control volume.
* \param element The finite element
* \param fvGeometry The finite volume geometry
* \param scvIdx The local index of the sub-control volume */
template<class ElementSolution>
AnsSurfaceParams aNonWettingSolidSurfaceParams(const Element &element,
const SubControlVolume &scv,
const ElementSolution &elemSol) const
{
DUNE_THROW(Dune::InvalidStateException,
"The spatial parameters do not provide a aNonWettingSolidSurfaceParams() method.");
}
/*!\brief Return a reference to the container object for the
* parametrization of the surface between non-Wetting and solid phase.
*
* The position is determined based on the coordinate of
* the vertex belonging to the considered sub control volume.
* \param element The finite element
* \param fvGeometry The finite volume geometry
* \param scvIdx The local index of the sub-control volume */
template<class ElementSolution>
AwsSurfaceParams aWettingSolidSurfaceParams(const Element &element,
const SubControlVolume &scv,
const ElementSolution &elemSol) const
{
DUNE_THROW(Dune::InvalidStateException,
"The spatial parameters do not provide a aWettingSolidSurfaceParams() method.");
}
/*!\brief Return the maximum capillary pressure for the given pc-Sw curve
*
* Of course physically there is no such thing as a maximum capillary pressure.
* The parametrization (VG/BC) goes to infinity and physically there is only one pressure
* for single phase conditions.
* Here, this is used for fitting the interfacial area surface: the capillary pressure,
* where the interfacial area is zero.
* Technically this value is obtained as the capillary pressure of saturation zero.
* This value of course only exists for the case of a regularized pc-Sw relation.
* \param element The finite element
* \param fvGeometry The finite volume geometry
* \param scvIdx The local index of the sub-control volume */
template<class ElementSolution>
const Scalar pcMax(const Element &element,
const SubControlVolume &scv,
const ElementSolution &elemSol) const
{ DUNE_THROW(Dune::InvalidStateException,
"The spatial parameters do not provide a aNonWettingSolidSurfaceParams() method.");
}
/*!\brief Return the characteristic length for the mass transfer.
*
* The position is determined based on the coordinate of
* the vertex belonging to the considered sub controle volume.
* \param element The finite element
* \param fvGeometry The finite volume geometry
* \param scvIdx The local index of the sub control volume */
template<class ElementSolution>
const Scalar characteristicLength(const Element & element,
const SubControlVolume &scv,
const ElementSolution &elemSol) const
{ return this->asImp_().characteristicLengthAtPos(scv.dofPosition()); }
/*!\brief Return the characteristic length for the mass transfer.
* \param globalPos The position in global coordinates.*/
const Scalar characteristicLengthAtPos(const GlobalPosition & globalPos) const
{
DUNE_THROW(Dune::InvalidStateException,
"The spatial parameters do not provide "
"a characteristicLengthAtPos() method.");
}
/*!\brief Return the pre factor the the energy transfer
*
* The position is determined based on the coordinate of
* the vertex belonging to the considered sub controle volume.
* \param element The finite element
* \param fvGeometry The finite volume geometry
* \param scvIdx The local index of the sub control volume */
template<class ElementSolution>
const Scalar factorEnergyTransfer(const Element &element,
const SubControlVolume &scv,
const ElementSolution &elemSol) const
{ return this->asImp_().factorEnergyTransferAtPos(scv.dofPosition()); }
/*!\brief Return the pre factor the the energy transfer
* \param globalPos The position in global coordinates.*/
const Scalar factorEnergyTransferAtPos(const GlobalPosition & globalPos) const
{
DUNE_THROW(Dune::InvalidStateException,
"The spatial parameters do not provide "
"a factorEnergyTransferAtPos() method.");
}
/*!\brief Return the pre factor the the mass transfer
*
* The position is determined based on the coordinate of
* the vertex belonging to the considered sub controle volume.
* \param element The finite element
* \param fvGeometry The finite volume geometry
* \param scvIdx The local index of the sub control volume */
template<class ElementSolution>
const Scalar factorMassTransfer(const Element &element,
const SubControlVolume &scv,
const ElementSolution &elemSol) const
{ return this->asImp_().factorMassTransferAtPos(scv.dofPosition()); }
/*!\brief Return the pre factor the the mass transfer
* \param globalPos The position in global coordinates.*/
const Scalar factorMassTransferAtPos(const GlobalPosition & globalPos) const
{
DUNE_THROW(Dune::InvalidStateException,
"The spatial parameters do not provide "
"a factorMassTransferAtPos() method.");
}
};
}
#endif // GUARDIAN
......@@ -68,10 +68,14 @@ class NonEquilibriumVolumeVariablesImplementation< Traits,
using ModelTraits = typename Traits::ModelTraits;
using Indices = typename ModelTraits::Indices;
using FS = typename Traits::FluidSystem;
static constexpr auto numEnergyEqFluid = ModelTraits::numEnergyEqFluid();
static constexpr auto numEnergyEqSolid = ModelTraits::numEnergyEqSolid();
static_assert((numEnergyEqFluid < 2), "This model is a specialization for a energy transfer of a fluid mixture and a solid");
static constexpr auto phase0Idx = FS::phase0Idx;
static constexpr auto phase1Idx = FS::phase1Idx;
static constexpr auto sPhaseIdx = FS::numPhases;
static_assert((numEnergyEqFluid < 3), "this model only has interfacial area relationships for maximum 2 fluids with one solid");
using DimLessNum = DimensionlessNumbers<Scalar>;
......@@ -97,7 +101,9 @@ public:
ParameterCache paramCache;
paramCache.updateAll(this->fluidState());
updateInterfacialArea(elemSol, this->fluidState(), paramCache, problem, element, scv);
updateDimLessNumbers(elemSol, this->fluidState(), paramCache, problem, element, scv);
if (numEnergyEqFluid == 2)
updateInterfacialArea(elemSol, this->fluidState(), paramCache, problem, element, scv);
}
/*!
......@@ -111,7 +117,7 @@ public:
* \param scv The sub-control volume
*/
template<class ElemSol, class Problem, class Element, class Scv>
void updateInterfacialArea(const ElemSol& elemSol,
void updateDimLessNumbers(const ElemSol& elemSol,
const FluidState& fluidState,
const ParameterCache& paramCache,
const Problem& problem,
......@@ -136,7 +142,7 @@ public:
const auto thermalConductivity = FluidSystem::thermalConductivity(fluidState, paramCache, phaseIdx);
const auto porosity = this->porosity();
reynoldsNumber_[phaseIdx] = DimLessNum::reynoldsNumber(darcyMagVelocity, characteristicLength_, kinematicViscosity);
reynoldsNumber_[phaseIdx] = DimLessNum::reynoldsNumber(darcyMagVelocity, characteristicLength_,kinematicViscosity);
prandtlNumber_[phaseIdx] = DimLessNum::prandtlNumber(dynamicViscosity, heatCapacity, thermalConductivity);
nusseltNumber_[phaseIdx] = DimLessNum::nusseltNumberForced(reynoldsNumber_[phaseIdx],
prandtlNumber_[phaseIdx],
......@@ -145,6 +151,80 @@ public:
}
}
/*!
* \brief Updates the volume specific interfacial area [m^2 / m^3] between the phases.
*
* \param elemSol A vector containing all primary variables connected to the element
* \param fluidState Container for all the secondary variables concerning the fluids
* \param paramCache The parameter cache corresponding to the fluid state
* \param problem The problem to be solved
* \param element An element which contains part of the control volume
* \param scv The sub-control volume
*/
template<class ElemSol, class Problem, class Element, class Scv>
void updateInterfacialArea(const ElemSol& elemSol,
const FluidState& fluidState,
const ParameterCache& paramCache,
const Problem& problem,
const Element& element,
const Scv& scv)
{
// obtain (standard) material parameters (needed for the residual saturations)
const auto& materialParams = problem.spatialParams().materialLawParams(element, scv, elemSol);
//obtain parameters for interfacial area constitutive relations
const auto& aWettingNonWettingSurfaceParams =problem.spatialParams().aWettingNonWettingSurfaceParams(element, scv, elemSol);
const auto& aNonWettingSolidSurfaceParams =problem.spatialParams().aNonWettingSolidSurfaceParams(element, scv, elemSol);
const Scalar pc = fluidState.pressure(phase1Idx) - fluidState.pressure(phase0Idx);
const Scalar Sw = fluidState.saturation(phase0Idx);
Scalar awn;
using AwnSurface = typename Problem::SpatialParams::AwnSurface;
awn = AwnSurface::interfacialArea(aWettingNonWettingSurfaceParams, materialParams, Sw, pc );
interfacialArea_[phase0Idx][phase1Idx] = awn;
interfacialArea_[phase1Idx][phase0Idx] = interfacialArea_[phase0Idx][phase1Idx];
interfacialArea_[phase0Idx][phase0Idx] = 0.;
using AnsSurface = typename Problem::SpatialParams::AnsSurface;
Scalar ans = AnsSurface::interfacialArea(aNonWettingSolidSurfaceParams, materialParams,Sw, pc);
// Switch for using a a_{wn} relations that has some "maximum capillary pressure" as parameter // That value is obtained by regularization of the pc(Sw) function.
#if USE_PCMAX
const Scalar pcMax = problem.spatialParams().pcMax(element, scv, elemSol);
//I know the solid surface from the pore network. But it is more consistent to use the fitvalue.
using AnsSurface = typename Problem::SpatialParams::AnsSurface;
solidSurface_ = AnsSurface::interfacialArea(aNonWettingSolidSurfaceParams, materialParams, /*Sw=*/0., pcMax);
const Scalar aws = solidSurface_ - ans;
interfacialArea_[phase0Idx][sPhaseIdx] = aws;
interfacialArea_[sPhaseIdx][phase0Idx] = interfacialArea_[phase0Idx][sPhaseIdx];
interfacialArea_[sPhaseIdx][sPhaseIdx] = 0.;
#else
using AwsSurface = typename Problem::SpatialParams::AwsSurface;
const auto& aWettingSolidSurfaceParams = problem.spatialParams().aWettingSolidSurfaceParams(element, scv, elemSol);
const auto aws = AwsSurface::interfacialArea(aWettingSolidSurfaceParams,materialParams, Sw, pc);
interfacialArea_[phase0Idx][sPhaseIdx] = aws ;
interfacialArea_[sPhaseIdx][phase0Idx] = interfacialArea_[phase0Idx][sPhaseIdx];
interfacialArea_[sPhaseIdx][sPhaseIdx] = 0.;
#endif
interfacialArea_[phase1Idx][sPhaseIdx] = ans;
interfacialArea_[sPhaseIdx][phase1Idx] = interfacialArea_[phase1Idx][sPhaseIdx];
interfacialArea_[phase1Idx][phase1Idx] = 0.;
}
/*!
* \brief The specific interfacial area between two fluid phases [m^2 / m^3]
* \note This is _only_ required by the kinetic mass/energy modules
*/
const Scalar interfacialArea(const unsigned int phaseIIdx, const unsigned int phaseJIdx) const
{
// there is no interfacial area between a phase and itself
assert(phaseIIdx not_eq phaseJIdx);
return interfacialArea_[phaseIIdx][phaseJIdx];
}
//! access function Reynolds Number
const Scalar reynoldsNumber(const unsigned int phaseIdx) const { return reynoldsNumber_[phaseIdx]; }
//! access function Prandtl Number
......@@ -168,6 +248,7 @@ private:
Scalar factorEnergyTransfer_;
Scalar factorMassTransfer_;
Scalar solidSurface_ ;
Scalar interfacialArea_[ModelTraits::numPhases()+numEnergyEqSolid][ModelTraits::numPhases()+numEnergyEqSolid];
};
////////////////////////////////////////////////////////////////////////////////////////////////////
......@@ -390,24 +471,6 @@ public:
Scalar awn;
// TODO can we delete this??? awn is overwritten anyway!
#define AwnRegul 0
// This regularizes the interfacial area between the fluid phases.
// This makes sure, that
// a) some saturation cannot be lost: Never leave two phase region.
// b) We cannot leave the fit region: no crazy (e.g. negative) values possible
// const Scalar Swr = aWettingNonWettingSurfaceParams.Swr() ;
// const Scalar Snr = aWettingNonWettingSurfaceParams.Snr() ;
// this just leads to a stalling newton error as soon as this kicks in.
// May be a spline or sth like this would help, but I do not which derivatives
// to specify.
#if AwnRegul
if(Sw < 5e-3 ) // or Sw > (1.-1e-5 )
awn = 0. ; // 10.; //
#endif
using AwnSurface = typename Problem::SpatialParams::AwnSurface;
awn = AwnSurface::interfacialArea(aWettingNonWettingSurfaceParams, materialParams, Sw, pc );
interfacialArea_[phase0Idx][phase1Idx] = awn;
......
......@@ -26,7 +26,7 @@
#define DUMUX_EVAPORATION_ATMOSPHERE_SPATIALPARAMS_HH
#include <dumux/porousmediumflow/properties.hh>
#include <dumux/material/spatialparams/fv.hh>
#include <dumux/material/spatialparams/fvnonequilibrium.hh>
#include <dumux/material/fluidmatrixinteractions/2p/linearmaterial.hh>
#include <dumux/material/fluidmatrixinteractions/2p/regularizedbrookscorey.hh>
#include <dumux/material/fluidmatrixinteractions/2p/efftoabslaw.hh>
......@@ -36,6 +36,12 @@
#include <dumux/material/fluidmatrixinteractions/mp/2poftadapter.hh>
#include <dumux/common/parameters.hh>
// material laws for interfacial area
#include <dumux/material/fluidmatrixinteractions/2pia/efftoabslawia.hh>
#include <dumux/material/fluidmatrixinteractions/2pia/awnsurfacepolynomial2ndorder.hh>
#include <dumux/material/fluidmatrixinteractions/2pia/awnsurfacepcmaxfct.hh>
#include <dumux/material/fluidmatrixinteractions/2pia/awnsurfaceexpswpcto3.hh>
namespace Dumux {
/**
......@@ -43,9 +49,9 @@ namespace Dumux {
*/
template<class TypeTag>
class EvaporationAtmosphereSpatialParams
: public FVSpatialParams<GetPropType<TypeTag, Properties::FVGridGeometry>,
GetPropType<TypeTag, Properties::Scalar>,
EvaporationAtmosphereSpatialParams<TypeTag>>
: public FVNonEquilibriumSpatialParams<GetPropType<TypeTag, Properties::FVGridGeometry>,
GetPropType<TypeTag, Properties::Scalar>,
EvaporationAtmosphereSpatialParams<TypeTag>>
{
using Scalar = GetPropType<TypeTag, Properties::Scalar>;
using FVGridGeometry = GetPropType<TypeTag, Properties::FVGridGeometry>;
......@@ -53,7 +59,7 @@ class EvaporationAtmosphereSpatialParams
using FVElementGeometry = typename FVGridGeometry::LocalView;
using SubControlVolume = typename FVElementGeometry::SubControlVolume;
using Element = typename GridView::template Codim<0>::Entity;
using ParentType = FVSpatialParams<FVGridGeometry, Scalar, EvaporationAtmosphereSpatialParams<TypeTag>>;
using ParentType = FVNonEquilibriumSpatialParams<FVGridGeometry, Scalar, EvaporationAtmosphereSpatialParams<TypeTag>>;
using GlobalPosition = Dune::FieldVector<Scalar, GridView::dimension>;
......@@ -69,13 +75,17 @@ public:
using PermeabilityType = Scalar;
//! export the material law type used
using MaterialLaw = TwoPAdapter<liquidPhaseIdx, EffToAbsLaw<RegularizedBrooksCorey<Scalar>>>;
//! export the types used for interfacial area calculations
using AwnSurface = GetPropType<TypeTag, Properties::AwnSurface>;
using AwsSurface = GetPropType<TypeTag, Properties::AwsSurface>;
using AnsSurface = GetPropType<TypeTag, Properties::AnsSurface>;
//! convenience aliases of the law parameters
using MaterialLawParams = typename MaterialLaw::Params;
//! export the types used for interfacial area calculations
using EffectiveIALawAws = AwnSurfacePolynomial2ndOrder<Scalar>;
using EffectiveIALawAwn = AwnSurfacePcMaxFct<Scalar>;
using EffectiveIALawAns = AwnSurfaceExpSwPcTo3<Scalar>;
using AwnSurface = EffToAbsLawIA<EffectiveIALawAwn, MaterialLawParams>;
using AwsSurface = EffToAbsLawIA<EffectiveIALawAws, MaterialLawParams>;
using AnsSurface = EffToAbsLawIA<EffectiveIALawAns, MaterialLawParams>;
using AwnSurfaceParams = typename AwnSurface::Params;
using AwsSurfaceParams = typename AwsSurface::Params;
using AnsSurfaceParams = typename AnsSurface::Params;
......@@ -284,21 +294,6 @@ public:
const ElementSolution &elemSol) const
{ return aWettingNonWettingSurfaceParams_.pcMax() ; }
/*!\brief Return the characteristic length for the mass transfer.
*
* The position is determined based on the coordinate of
* the vertex belonging to the considered sub controle volume.
* \param element The finite element
* \param fvGeometry The finite volume geometry
* \param scvIdx The local index of the sub control volume */
template<class ElementSolution>
const Scalar characteristicLength(const Element & element,
const SubControlVolume &scv,
const ElementSolution &elemSol) const
{ return characteristicLengthAtPos(scv.dofPosition()); }
/*!\brief Return the characteristic length for the mass transfer.
* \param globalPos The position in global coordinates.*/
const Scalar characteristicLengthAtPos(const GlobalPosition & globalPos) const
......@@ -310,19 +305,6 @@ public:
else DUNE_THROW(Dune::InvalidStateException, "You should not be here: x=" << globalPos[0] << " y= "<< globalPos[dimWorld-1]);
}
/*!\brief Return the pre factor the the energy transfer
*
* The position is determined based on the coordinate of
* the vertex belonging to the considered sub controle volume.
* \param element The finite element
* \param fvGeometry The finite volume geometry
* \param scvIdx The local index of the sub control volume */
template<class ElementSolution>
const Scalar factorEnergyTransfer(const Element &element,
const SubControlVolume &scv,
const ElementSolution &elemSol) const
{ return factorEnergyTransferAtPos(scv.dofPosition()); }
/*!\brief Return the pre factor the the energy transfer
* \param globalPos The position in global coordinates.*/
const Scalar factorEnergyTransferAtPos(const GlobalPosition & globalPos) const
......@@ -334,19 +316,6 @@ public:
else DUNE_THROW(Dune::InvalidStateException, "You should not be here: x=" << globalPos[0] << " y= "<< globalPos[dimWorld-1]);
}
/*!\brief Return the pre factor the the mass transfer
*
* The position is determined based on the coordinate of
* the vertex belonging to the considered sub controle volume.
* \param element The finite element
* \param fvGeometry The finite volume geometry
* \param scvIdx The local index of the sub control volume */
template<class ElementSolution>
const Scalar factorMassTransfer(const Element &element,
const SubControlVolume &scv,
const ElementSolution &elemSol) const
{ return factorMassTransferAtPos(scv.dofPosition()); }
/*!\brief Return the pre factor the the mass transfer
* \param globalPos The position in global coordinates.*/
const Scalar factorMassTransferAtPos(const GlobalPosition & globalPos) const
......
......@@ -26,6 +26,7 @@
#include <dune/common/parametertreeparser.hh>
#include <dumux/material/spatialparams/fvnonequilibrium.hh>
#include <dumux/material/fluidmatrixinteractions/2p/efftoabslaw.hh>
#include <dumux/material/fluidmatrixinteractions/2p/heatpipelaw.hh>
#include <dumux/material/fluidmatrixinteractions/2p/linearmaterial.hh>
......@@ -33,6 +34,12 @@
#include <dumux/porousmediumflow/properties.hh>
#include <dumux/material/spatialparams/fv.hh>
// material laws for interfacial area
#include <dumux/material/fluidmatrixinteractions/2pia/efftoabslawia.hh>
#include <dumux/material/fluidmatrixinteractions/2pia/awnsurfacepolynomial2ndorder.hh>
#include <dumux/material/fluidmatrixinteractions/2pia/awnsurfacepcmaxfct.hh>
#include <dumux/material/fluidmatrixinteractions/2pia/awnsurfaceexpswpcto3.hh>
namespace Dumux {
/**
......@@ -41,9 +48,9 @@ namespace Dumux {
*/
template<class TypeTag>
class CombustionSpatialParams
: public FVSpatialParams<GetPropType<TypeTag, Properties::FVGridGeometry>,
GetPropType<TypeTag, Properties::Scalar>,
CombustionSpatialParams<TypeTag>>
: public FVNonEquilibriumSpatialParams<GetPropType<TypeTag, Properties::FVGridGeometry>,
GetPropType<TypeTag, Properties::Scalar>,
CombustionSpatialParams<TypeTag>>
{
using Scalar = GetPropType<TypeTag, Properties::Scalar>;
using FVGridGeometry = GetPropType<TypeTag, Properties::FVGridGeometry>;
......@@ -51,7 +58,7 @@ class CombustionSpatialParams
using FVElementGeometry = typename FVGridGeometry::LocalView;
using SubControlVolume = typename FVElementGeometry::SubControlVolume;
using Element = typename GridView::template Codim<0>::Entity;
using ParentType = FVSpatialParams<FVGridGeometry, Scalar, CombustionSpatialParams<TypeTag>>;
using ParentType = FVNonEquilibriumSpatialParams<FVGridGeometry, Scalar, CombustionSpatialParams<TypeTag>>;
enum {dimWorld = GridView::dimensionworld};
using GlobalPosition = typename SubControlVolume::GlobalPosition;
......@@ -68,6 +75,18 @@ public:
using MaterialLaw = TwoPAdapter<wPhaseIdx, EffToAbsLaw<EffectiveLaw>>;
using MaterialLawParams = typename MaterialLaw::Params;
//! export the types used for interfacial area calculations
using EffectiveIALawAws = AwnSurfacePolynomial2ndOrder<Scalar>;
using EffectiveIALawAwn = AwnSurfacePcMaxFct<Scalar>;
using EffectiveIALawAns = AwnSurfaceExpSwPcTo3<Scalar>;
using AwnSurface = EffToAbsLawIA<EffectiveIALawAwn, MaterialLawParams>;
using AwsSurface = EffToAbsLawIA<EffectiveIALawAws, MaterialLawParams>;
using AnsSurface = EffToAbsLawIA<EffectiveIALawAns, MaterialLawParams>;
using AwnSurfaceParams = typename AwnSurface::Params;
using AwsSurfaceParams = typename AwsSurface::Params;
using AnsSurfaceParams = typename AnsSurface::Params;
CombustionSpatialParams(std::shared_ptr<const FVGridGeometry> fvGridGeometry) : ParentType(fvGridGeometry)
{
// this is the parameter value from file part
......@@ -166,59 +185,17 @@ public:
const MaterialLawParams& materialLawParamsAtPos(const GlobalPosition & globalPos) const
{ return materialParams_ ; }
/*!\brief Return the characteristic length for the mass transfer.
*
* The position is determined based on the coordinate of
* the vertex belonging to the considered sub controle volume.
* \param element The finite element
* \param fvGeometry The finite volume geometry
* \param scvIdx The local index of the sub control volume */
template<class ElementSolution>
const Scalar characteristicLength(const Element & element,
const SubControlVolume &scv,
const ElementSolution &elemSol) const
{ return characteristicLengthAtPos(scv.center()); }
/*!\brief Return the characteristic length for the mass transfer.
* \param globalPos The position in global coordinates.*/
const Scalar characteristicLengthAtPos(const GlobalPosition & globalPos) const
{ return characteristicLength_ ; }
/*!\brief Return the pre factor the the energy transfer
*
* The position is determined based on the coordinate of
* the vertex belonging to the considered sub controle volume.
* \param element The finite element
* \param fvGeometry The finite volume geometry
* \param scvIdx The local index of the sub control volume */
template<class ElementSolution>
const Scalar factorEnergyTransfer(const Element &element,
const SubControlVolume &scv,
const ElementSolution &elemSol) const
{ return factorEnergyTransferAtPos(scv.dofPosition()); }