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

[volumevariables][model][test] add solidstate and solidsystem in all models

The volumevariables derive now from the energyvolumevariables which set the solid propererties needed for nonisothermal models in the new solidstate
Each model now also includes a solidstate (inert or compositional) and each problem a solidsystem.
For nonreactive models the default is set to constant component where properties can be specified in the input file.
The fluidsystems now do not have solid properties anymore. Each solid component now has all thermal properties which might depend on temperature.
parent 0ad80436
......@@ -65,7 +65,7 @@ class FouriersLawNonEquilibriumImplementation<TypeTag, DiscretizationMethod::box
enum { dimWorld = GridView::dimensionworld} ;
enum { numPhases = GET_PROP_TYPE(TypeTag, ModelTraits)::numPhases()} ;
enum { numEnergyEqFluid = GET_PROP_VALUE(TypeTag, NumEnergyEqFluid) };
enum {sPhaseIdx = FluidSystem::sPhaseIdx};
enum {sPhaseIdx = FluidSystem::numPhases};
using DimWorldMatrix = Dune::FieldMatrix<Scalar, dimWorld, dimWorld>;
......
......@@ -64,7 +64,7 @@ public:
/*!
* \brief The mass density \f$\mathrm{[kg/m^3]}\f$ of CaO.
*/
static Scalar density()
static Scalar solidDensity(Scalar temperature)
{
return 3370;
}
......@@ -72,10 +72,18 @@ public:
/*!
* \brief The specific heat capacity \f$\mathrm{[J/kg K]}\f$ of CaO.
*/
static Scalar heatCapacity()
static Scalar solidHeatCapacity(Scalar temperature)
{
return 934; //Nagel et al. (2014) : 934 J/kgK
}
/*!
* \brief The thermal conductivity \f$\mathrm{[W/(m K)]}\f$ of the porous material.
*/
static Scalar solidThermalConductivity(Scalar temperature)
{
return 0.4; //Nagel et al. (2014)
}
};
} // end namespace Components
......
......@@ -65,7 +65,7 @@ public:
/*!
* \brief The mass density \f$\mathrm{[kg/m^3]}\f$ of CaO2H2.
*/
static Scalar density()
static Scalar solidDensity(Scalar temperature)
{
return 2200.0; //at 293 K ; Shao et al. (2013)
}
......@@ -73,10 +73,18 @@ public:
/*!
* \brief The specific heat capacity \f$\mathrm{[J/kgK]}\f$ of CaO2H2.
*/
static Scalar heatCapacity()
static Scalar solidHeatCapacity(Scalar temperature)
{
return 1530; //Nagel et al. (2014) : 1530 J/kgK
}
/*!
* \brief The thermal conductivity \f$\mathrm{[W/(m K)]}\f$ of the porous material.
*/
static Scalar solidThermalConductivity(Scalar temperature)
{
return 0.4; //Nagel et al. (2014)
}
};
} // end namespace Components
......
......@@ -76,7 +76,7 @@ public:
/*!
* \brief The mass density \f$\mathrm{[kg/m^3]}\f$ of NaCl.
*/
static Scalar density()
static Scalar solidDensity(Scalar temperature)
{
return 2165.0;
}
......@@ -84,7 +84,7 @@ public:
/*!
* \brief The specific heat capacity \f$\mathrm{[J/molK]}\f$ of NaCl.
*/
static Scalar heatCapacity()
static Scalar solidHeatCapacity(Scalar temperature)
{
return 50.50;
}
......
......@@ -54,36 +54,10 @@ public:
* \param pressure pressure of component in \f$\mathrm{[Pa]}\f$
*/
template<class C = Component>
static Scalar solidDensity(Scalar temperature, Scalar pressure)
static Scalar solidDensity(Scalar temperature)
{
static_assert(AlwaysFalse<C>::value, "Mandatory function not implemented: solidDensity(t,p)");
DUNE_THROW(Dune::NotImplemented, "solidDensity(t,p)");
}
/*!
* \brief Specific enthalpy \f$\mathrm{[J/kg]}\f$ of the pure component in solid.
*
* \param temperature temperature of component in \f$\mathrm{[K]}\f$
* \param pressure pressure of component in \f$\mathrm{[Pa]}\f$
*/
template<class C = Component>
static const Scalar solidEnthalpy(Scalar temperature, Scalar pressure)
{
static_assert(AlwaysFalse<C>::value, "Mandatory function not implemented: solidEnthalpy(t,p)");
DUNE_THROW(Dune::NotImplemented, "solidEnthalpy(t,p)");
}
/*!
* \brief Specific internal energy \f$\mathrm{[J/kg]}\f$ of the pure component in solid.
*
* \param temperature temperature of component in \f$\mathrm{[K]}\f$
* \param pressure pressure of component in \f$\mathrm{[Pa]}\f$
*/
template<class C = Component>
static const Scalar solidInternalEnergy(Scalar temperature, Scalar pressure)
{
static_assert(AlwaysFalse<C>::value, "Mandatory function not implemented: solidInternalEnergy(t,p)");
DUNE_THROW(Dune::NotImplemented, "solidInternalEnergy(t,p)");
static_assert(AlwaysFalse<C>::value, "Mandatory function not implemented: solidDensity(t)");
DUNE_THROW(Dune::NotImplemented, "solidDensity(t)");
}
/*!
......@@ -92,10 +66,10 @@ public:
* \param pressure pressure of component in \f$\mathrm{[Pa]}\f$
*/
template<class C = Component>
static Scalar solidThermalConductivity(Scalar temperature, Scalar pressure)
static Scalar solidThermalConductivity(Scalar temperature)
{
static_assert(AlwaysFalse<C>::value, "Mandatory function not implemented: solidThermalConductivity(t,p)");
DUNE_THROW(Dune::NotImplemented, "solidThermalConductivity(t,p)");
static_assert(AlwaysFalse<C>::value, "Mandatory function not implemented: solidThermalConductivity(t)");
DUNE_THROW(Dune::NotImplemented, "solidThermalConductivity(t)");
}
/*!
......@@ -104,10 +78,10 @@ public:
* \param pressure pressure of component in \f$\mathrm{[Pa]}\f$
*/
template<class C = Component>
static Scalar solidHeatCapacity(Scalar temperature, Scalar pressure)
static Scalar solidHeatCapacity(Scalar temperature)
{
static_assert(AlwaysFalse<C>::value, "Mandatory function not implemented: solidHeatCapacity(t,p)");
DUNE_THROW(Dune::NotImplemented, "solidHeatCapacity(t,p)");
static_assert(AlwaysFalse<C>::value, "Mandatory function not implemented: solidHeatCapacity(t)");
DUNE_THROW(Dune::NotImplemented, "solidHeatCapacity(t)");
}
};
......
......@@ -75,11 +75,9 @@ public:
****************************************/
static constexpr int numPhases = 2; // liquid and gas phases
static constexpr int numComponents = 3; // H2O, Air, NaCl
static constexpr int numSPhases = 1;// precipitated solid phases // TODO: remove
static constexpr int liquidPhaseIdx = 0; // index of the liquid phase
static constexpr int gasPhaseIdx = 1; // index of the gas phase
static constexpr int solidPhaseIdx = 2; // index of the precipitated salt // TODO: remove
static constexpr int phase0Idx = liquidPhaseIdx; // index of the first phase
static constexpr int phase1Idx = gasPhaseIdx; // index of the second phase
......@@ -91,7 +89,7 @@ public:
static constexpr int AirIdx = 1;
static constexpr int comp0Idx = H2OIdx;
static constexpr int comp1Idx = AirIdx;
static constexpr int NaClIdx = 2; // TODO: remove
static constexpr int NaClIdx = 2;
/*!
* \brief Return the human readable name of a fluid phase
......@@ -100,10 +98,10 @@ public:
*/
static std::string phaseName(int phaseIdx)
{
switch (phaseIdx) {
case liquidPhaseIdx: return "liquid";
case gasPhaseIdx: return "gas";
case solidPhaseIdx: return "NaCl";
switch (phaseIdx)
{
case liquidPhaseIdx: return "liquid";
case gasPhaseIdx: return "gas";
}
DUNE_THROW(Dune::InvalidStateException, "Invalid phase index " << phaseIdx);
}
......@@ -220,18 +218,6 @@ public:
DUNE_THROW(Dune::InvalidStateException, "Invalid component index " << compIdx);
}
/*!
* \brief Return the mass density of the precipitate \f$\mathrm{[kg/m^3]}\f$.
*
* \param phaseIdx The index of the precipitated phase to consider
*/
static Scalar precipitateDensity(int phaseIdx)
{
if(phaseIdx != solidPhaseIdx)
DUNE_THROW(Dune::InvalidStateException, "Invalid solid phase index " << solidPhaseIdx);
return NaCl::density();
}
/*!
* \brief Return the saturation vapor pressure of the liquid phase \f$\mathrm{[Pa]}\f$.
*
......@@ -243,25 +229,6 @@ public:
return vaporPressure_(Temperature,salinity);
}
/*!
* \brief Return the salt specific heat capacity \f$\mathrm{[J/molK]}\f$.
*
* \param phaseIdx The index of the precipitated phase to consider
*/
static Scalar precipitateHeatCapacity(int phaseIdx)
{
return NaCl::heatCapacity();
}
/*!
* \brief Return the molar density of the precipitate \f$\mathrm{[mol/m^3]}\f$.
*
* \param phaseIdx The index of the precipitated phase to consider
*/
static Scalar precipitateMolarDensity(int phaseIdx)
{
return precipitateDensity(phaseIdx)/molarMass(phaseIdx); //TODO this only works for this specific case here with phaseIdx=compIdx!
}
/****************************************
* thermodynamic relations
......
......@@ -58,17 +58,6 @@ public:
//! The type of parameter cache objects
using ParameterCache = NullParameterCache;
//! Index of the solid phase
static constexpr int sPhaseIdx = 2;
static std::string phaseName(int phaseIdx)
{
if (phaseIdx == sPhaseIdx)
return "s";
return ParentType::phaseName(phaseIdx);
}
/*!
* \brief Return the enthalpy of a component in a phase.
* \param fluidState A container with the current (physical) state of the fluid
......
......@@ -34,8 +34,6 @@
#include <dumux/material/fluidsystems/base.hh>
#include <dumux/material/components/n2.hh>
#include <dumux/material/components/h2o.hh>
#include <dumux/material/components/cao2h2.hh>
#include <dumux/material/components/cao.hh>
#include <dumux/material/binarycoefficients/h2o_n2.hh>
#include <dumux/material/components/tabulatedcomponent.hh>
......@@ -79,22 +77,15 @@ public:
//! Number of phases in the fluid system
static constexpr int numPhases = 1; //gas phase: N2 and steam
static constexpr int numComponents = 2; // H2O, Air
static constexpr int numSPhases = 2;// solid phases CaO and CaO2H2 TODO: Remove
static constexpr int numSComponents = 2;// CaO2H2, CaO TODO: Remove
static constexpr int gasPhaseIdx = 0; //!< index of the gas phase
static constexpr int phase0Idx = gasPhaseIdx; //!< index of the only phase
static constexpr int cPhaseIdx = 1; // CaO-phaseIdx TODO: Remove
static constexpr int hPhaseIdx = 2; // CaO2H2-phaseIdx TODO: Remove
static constexpr int N2Idx = 0;
static constexpr int H2OIdx = 1;
static constexpr int comp0Idx = N2Idx;
static constexpr int comp1Idx = H2OIdx;
static constexpr int CaOIdx = 2; // CaO-phaseIdx TODO: Remove
static constexpr int CaO2H2Idx = 3; // CaO-phaseIdx TODO: Remove
/****************************************
* Fluid phase related static parameters
......@@ -108,8 +99,6 @@ public:
{
switch (phaseIdx) {
case gasPhaseIdx: return "gas";
case cPhaseIdx : return "CaO";
case hPhaseIdx : return "CaOH2";
}
DUNE_THROW(Dune::InvalidStateException, "Invalid phase index " << phaseIdx);
}
......@@ -197,9 +186,6 @@ public:
{
case H2OIdx: return "H2O";
case N2Idx: return "N2";
case CaOIdx:return "CaO";
case CaO2H2Idx:return "CaO2H2";
}
DUNE_THROW(Dune::InvalidStateException, "Invalid component index " << compIdx);
}
......@@ -215,59 +201,10 @@ public:
{
case H2OIdx: return H2O::molarMass();
case N2Idx: return N2::molarMass();
case CaOIdx: return CaO::molarMass();
case CaO2H2Idx: return CaO2H2::molarMass();
}
DUNE_THROW(Dune::InvalidStateException, "Invalid component index " << compIdx);
}
/*!
* \brief Return the mass density of the solid \f$\mathrm{[kg/m^3]}\f$.
*
* \param phaseIdx The index of the solid phase to consider
*/
static Scalar precipitateDensity(int phaseIdx)
{
if(phaseIdx==cPhaseIdx)
return CaO::density();
if(phaseIdx==hPhaseIdx)
return CaO2H2::density();
else
DUNE_THROW(Dune::InvalidStateException, "Invalid solid index " << phaseIdx);
return 1;
}
/*!
* \brief Return the salt specific heat capacity \f$\mathrm{[J/molK]}\f$.
*
* \param phaseIdx The index of the solid phase to consider
*/
static Scalar precipitateHeatCapacity(int phaseIdx)
{
if(phaseIdx==cPhaseIdx)
return CaO::heatCapacity();
if(phaseIdx==hPhaseIdx)
return CaO2H2::heatCapacity();
else
DUNE_THROW(Dune::InvalidStateException, "Invalid solid index " << phaseIdx);
return 1;
}
/*!
* \brief Return the molar density of the solid \f$\mathrm{[mol/m^3]}\f$.
*
* \param phaseIdx The index of the solid phase to consider
*/
static Scalar precipitateMolarDensity(int phaseIdx)
{
if(phaseIdx==1){
return precipitateDensity(phaseIdx)/ molarMass(CaOIdx);
}
if(phaseIdx==2){
return precipitateDensity(phaseIdx)/molarMass(CaO2H2Idx); }
else
DUNE_THROW(Dune::InvalidStateException, "Invalid solid phase index " << phaseIdx);
}
/****************************************
* thermodynamic relations
......
......@@ -43,7 +43,7 @@ void updateSolidVolumeFractions(const ElemSol &elemSol,
{
for(int sCompIdx = solidState.numComponents- solidState.numInertComponents; sCompIdx < solidState.numComponents; ++sCompIdx)
{
solidState.setVolumeFraction(sCompIdx, problem.spatialParams().inertVolumeFraction(element, scv, sCompIdx));
solidState.setVolumeFraction(sCompIdx, problem.spatialParams().inertVolumeFraction(element, scv, solidState, sCompIdx));
}
if (!(solidState.isInert()))
{
......
......@@ -173,13 +173,10 @@ public:
*
* \param element The current element
* \param scv The sub-control volume inside the element.
* \param elemSol The solution at the dofs connected to the element.
* \return the porosity
*/
template<class ElementSolution>
Scalar porosity(const Element& element,
const SubControlVolume& scv,
const ElementSolution& elemSol) const
const SubControlVolume& scv) const
{
return asImp_().porosityAtPos(scv.center());
}
......@@ -193,97 +190,56 @@ public:
Scalar porosityAtPos(const GlobalPosition& globalPos) const
{
DUNE_THROW(Dune::InvalidStateException,
"The spatial parameters do not provide "
"a porosityAtPos() method.");
}
/*!
* \brief Returns the heat capacity \f$[J / (kg K)]\f$ of the rock matrix.
*
* This is only required for non-isothermal models.
*
* \param element The current element
* \param scv The sub-control volume inside the element.
* \param elemSol The solution at the dofs connected to the element.
*/
template<class ElementSolution>
Scalar solidHeatCapacity(const Element &element,
const SubControlVolume& scv,
const ElementSolution& elemSol) const
{
return asImp_().solidHeatCapacityAtPos(scv.center());
"The spatial parameters do not provide a porosityAtPos() method.");
}
/*!
* \brief Returns the heat capacity \f$[J / (kg K)]\f$ of the rock matrix.
*
* This is only required for non-isothermal models.
*
* \param globalPos The position of the center of the element
*/
Scalar solidHeatCapacityAtPos(const GlobalPosition& globalPos) const
{
DUNE_THROW(Dune::InvalidStateException,
"The spatial parameters do not provide "
"a solidHeatCapacityAtPos() method.");
}
/*!
* \brief Returns the mass density \f$[kg / m^3]\f$ of the rock matrix.
*
* This is only required for non-isothermal models.
*
* \param element The current element
* \param scv The sub-control volume inside the element.
* \param elemSol The solution at the dofs connected to the element.
*/
template<class ElementSolution>
Scalar solidDensity(const Element &element,
const SubControlVolume& scv,
const ElementSolution& elemSol) const
{
return asImp_().solidDensityAtPos(scv.center());
}
/*!
* \brief Returns the mass density \f$[kg / m^3]\f$ of the rock matrix.
*
* This is only required for non-isothermal models.
* \brief Function for defining the porosity.
*
* \param globalPos The position of the center of the element
* \return porosity
* \param globalPos The position of the center of the scv
*/
Scalar solidDensityAtPos(const GlobalPosition& globalPos) const
{
DUNE_THROW(Dune::InvalidStateException,
"The spatial parameters do not provide "
"a solidDensityAtPos() method.");
}
Scalar minimalPorosity(const Element& element,
const SubControlVolume& scv) const
{ return 0.0; }
/*!
* \brief Returns the thermal conductivity \f$\mathrm{[W/(m K)]}\f$ of the porous material.
* \brief Function for defining the porosity.
* That is possibly solution dependent.
*
* \param element The current element
* \param scv The sub-control volume inside the element.
* \param elemSol The solution at the dofs connected to the element.
* \return the porosity
*/
template<class ElementSolution>
Scalar solidThermalConductivity(const Element &element,
const SubControlVolume& scv,
const ElementSolution& elemSol) const
template<class SolidState>
Scalar inertVolumeFraction(const Element& element,
const SubControlVolume& scv,
SolidState& solidState,
int compIdx) const
{
return asImp_().solidThermalConductivityAtPos(scv.center());
return asImp_().inertVolumeFractionAtPos(scv.center(), solidState, compIdx);
}
/*!
* \brief Returns the thermal conductivity \f$\mathrm{[W/(m K)]}\f$ of the porous material.
* \brief Function for defining the porosity.
*
* \param globalPos The position of the center of the element
* \return porosity
* \param globalPos The position of the center of the scv
*/
Scalar solidThermalConductivityAtPos(const GlobalPosition& globalPos) const
template<class SolidState>
Scalar inertVolumeFractionAtPos(const GlobalPosition& globalPos,
SolidState& solidState,
int compIdx) const
{
DUNE_THROW(Dune::InvalidStateException,
"The spatial parameters do not provide "
"a solidThermalConductivityAtPos() method.");
if (solidState.isInert() && solidState.numInertComponents == 1)
{
return 1-asImp_().porosityAtPos(globalPos);
}
else
DUNE_THROW(Dune::InvalidStateException,
"The spatial parameters do not provide "
"a inertVolumeFractionAtPos() method.");
}
/*!
......
......@@ -46,6 +46,9 @@
#include <dumux/material/fluidmatrixinteractions/1p/thermalconductivityaverage.hh>
#include <dumux/material/fluidstates/immiscible.hh>
#include <dumux/material/solidstates/inertsolidstate.hh>
#include <dumux/material/solidsystems/inertsolidphase.hh>
#include <dumux/material/components/constant.hh>
#include <dumux/porousmediumflow/properties.hh>
#include <dumux/porousmediumflow/immiscible/localresidual.hh>
......@@ -89,6 +92,8 @@ struct OnePModelTraits
template<class PV,
class FSY,
class FST,
class SSY,
class SST,
class PT,
class MT>
struct OnePVolumeVariablesTraits
......@@ -96,6 +101,8 @@ struct OnePVolumeVariablesTraits
using PrimaryVariables = PV;
using FluidSystem = FSY;
using FluidState = FST;
using SolidSystem = SSY;
using SolidState = SST;
using PermeabilityType = PT;
using ModelTraits = MT;
};
......@@ -120,10 +127,12 @@ private:
using PV = typename GET_PROP_TYPE(TypeTag, PrimaryVariables);
using FSY = typename GET_PROP_TYPE(TypeTag, FluidSystem);
using FST = typename GET_PROP_TYPE(TypeTag, FluidState);
using SSY = typename GET_PROP_TYPE(TypeTag, SolidSystem);
using SST = typename GET_PROP_TYPE(TypeTag, SolidState);
using MT = typename GET_PROP_TYPE(TypeTag, ModelTraits);
using PT = typename GET_PROP_TYPE(TypeTag, SpatialParams)::PermeabilityType;
using Traits = OnePVolumeVariablesTraits<PV, FSY, FST, PT, MT>;
using Traits = OnePVolumeVariablesTraits<PV, FSY, FST, SSY, SST, PT, MT>;
public:
using type = OnePVolumeVariables<Traits>;
};
......@@ -143,6 +152,24 @@ public:
using type = ImmiscibleFluidState<Scalar, FluidSystem>;
};
//! The two-phase model uses the immiscible fluid state
SET_PROP(OneP, SolidState)
{
private:
using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar);
using SolidSystem = typename GET_PROP_TYPE(TypeTag, SolidSystem);
public:
using type = InertSolidState<Scalar, SolidSystem>;
};
// Set the fluid system
SET_PROP(OneP, SolidSystem)
{
using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar);
using InertComponent = Components::Constant<1, Scalar>;
using type = SolidSystems::InertSolidPhase<Scalar, InertComponent>;
};
///////////////////////////////////////////////////////////////////////////
// properties for the non-isothermal single phase model
///////////////////////////////////////////////////////////////////////////
......
......@@ -26,6 +26,7 @@
#include <dumux/common/properties.hh>
#include <dumux/porousmediumflow/volumevariables.hh>
#include <dumux/porousmediumflow/nonisothermal/volumevariables.hh>
#include <dumux/material/fluidstates/immiscible.hh>
namespace Dumux {
......@@ -38,20 +39,26 @@ namespace Dumux {
* \tparam Traits Class encapsulating types to be used by the vol vars
*/
template<class Traits>
class OnePVolumeVariables : public PorousMediumFlowVolumeVariables< Traits, OnePVolumeVariables<Traits> >
class OnePVolumeVariables
: public PorousMediumFlowVolumeVariables< Traits>
, public EnergyVolumeVariables<Traits, OnePVolumeVariables<Traits> >
{
using ThisType = OnePVolumeVariables<Traits>;
using ParentType = PorousMediumFlowVolumeVariables<Traits, ThisType>;
using ParentType = PorousMediumFlowVolumeVariables<Traits>;
using EnergyVolVars = EnergyVolumeVariables<Traits, OnePVolumeVariables<Traits> >;
using Scalar = typename Traits::PrimaryVariables::value_type;
using Indices = typename Traits::ModelTraits::Indices;
using PermeabilityType = typename Traits::PermeabilityType;