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

[cleanup] general cleanup

parent 6ef9ad32
......@@ -40,7 +40,7 @@ namespace Components {
*/
template <class Scalar>
class Granite : public Components::Base<Scalar, Granite<Scalar> >
,public Components::Solid<Scalar, Granite<Scalar> >
, public Components::Solid<Scalar, Granite<Scalar> >
{
public:
......
#install headers
install(FILES
basesolidstate.hh
updatesolidvolumefractions.hh
inertsolidstate.hh
compositionalsolidstate.hh
DESTINATION ${CMAKE_INSTALL_INCLUDEDIR}/dumux/material/solidstates)
......@@ -26,10 +26,7 @@
#ifndef DUMUX_SOLID_STATE_COMPOSITIONAL_HH
#define DUMUX_SOLID_STATE_COMPOSITIONAL_HH
#include <dumux/common/valgrind.hh>
#include <limits>
#include "basesolidstate.hh"
#include "updatesolidvolumefractions.hh"
namespace Dumux
{
......@@ -67,7 +64,7 @@ public:
for (int compIdx =0; compIdx < numComponents; ++compIdx)
sumVolumeFraction += volumeFraction(compIdx);
Scalar porosity = 1-sumVolumeFraction;
return std::max(porosity, minPorosity_);
return porosity;
}
/*!
* \brief The mass density \f$\rho_\alpha\f$ of the fluid phase
......@@ -99,8 +96,8 @@ public:
Scalar temperature() const
{ return temperature_; }
Scalar volumeFraction(const int phaseIdx) const
{ return volumeFraction_[phaseIdx]; }
Scalar volumeFraction(const int compIdx) const
{ return volumeFraction_[compIdx]; }
/*****************************************************
* Setter methods. Note that these are not part of the
......@@ -117,11 +114,15 @@ public:
* thermodynamic equilibrium, the result of this method is
* undefined.
*/
template <class CompositionalSolidState>
void assign(const CompositionalSolidState &sst)
template <class SolidState>
void assign(const SolidState &sst)
{
temperature_ = sst.temperature();
density_ = sst.density();
thermalConducivity_ = sst.thermalConductivity();
heatCapacity_ = sst.heatCapacity();
for (int compIdx = 0; compIdx < numComponents; ++compIdx)
volumeFraction_[compIdx] = sst.volumeFraction(compIdx);
}
/*!
......@@ -151,14 +152,10 @@ public:
void setVolumeFraction(const int compIdx, Scalar value)
{ volumeFraction_[compIdx] = value; }
void setMinPorosity(Scalar value)
{ minPorosity_ = value; }
protected:
Scalar density_;
Scalar temperature_;
Scalar volumeFraction_[numComponents];
Scalar minPorosity_;
Scalar heatCapacity_;
Scalar thermalConducivity_;
};
......
......@@ -26,9 +26,7 @@
#ifndef DUMUX_INERT_SOLID_STATE_HH
#define DUMUX_INERT_SOLID_STATE_HH
#include <dumux/common/valgrind.hh>
#include "basesolidstate.hh"
#include <limits>
#include "updatesolidvolumefractions.hh"
namespace Dumux
{
......@@ -47,7 +45,11 @@ public:
numComponents = SolidSystem::numComponents,
numInertComponents = SolidSystem::numInertComponents,
};
static constexpr bool isInert() { return SolidSystem::isInert(); }
static constexpr bool isInert()
{
static_assert(SolidSystem::isInert(), "Only inert solid systems are allowed with the InertSolidState");
return true;
}
/*!
* \brief The average molar mass \f$\overline M_\alpha\f$ of phase \f$\alpha\f$ in \f$\mathrm{[kg/mol]}\f$
......@@ -67,7 +69,7 @@ public:
for (int compIdx =0; compIdx < numComponents; ++compIdx)
sumVolumeFraction += volumeFraction(compIdx);
Scalar porosity = 1-sumVolumeFraction;
return std::max(porosity, minPorosity_);
return porosity;
}
/*!
......@@ -101,8 +103,8 @@ public:
Scalar temperature() const
{ return temperature_; }
Scalar volumeFraction(const int phaseIdx) const
{ return volumeFraction_[phaseIdx]; }
Scalar volumeFraction(const int compIdx) const
{ return volumeFraction_[compIdx]; }
/*****************************************************
* Setter methods. Note that these are not part of the
......@@ -119,11 +121,15 @@ public:
* thermodynamic equilibrium, the result of this method is
* undefined.
*/
template <class InertSolidState>
void assign(const InertSolidState &sst)
template <class SolidState>
void assign(const SolidState &sst)
{
temperature_ = sst.temperature();
density_ = sst.density();
thermalConducivity_ = sst.thermalConductivity();
heatCapacity_ = sst.heatCapacity();
for (int compIdx = 0; compIdx < numComponents; ++compIdx)
volumeFraction_[compIdx] = sst.volumeFraction(compIdx);
}
/*!
......@@ -153,14 +159,10 @@ public:
void setVolumeFraction(const int compIdx, Scalar value)
{ volumeFraction_[compIdx] = value; }
void setMinPorosity(Scalar value)
{ minPorosity_ = value; }
protected:
Scalar density_;
Scalar temperature_;
Scalar volumeFraction_[numComponents];
Scalar minPorosity_;
Scalar heatCapacity_;
Scalar thermalConducivity_;
};
......
......@@ -19,16 +19,10 @@
/*!
* \file
* \ingroup SolidStates
* \brief Represents all relevant thermodynamic quantities of a
* multi-phase fluid system assuming immiscibility and
* thermodynamic equilibrium.
* \brief update the solid volume fractions (inert and reacitve) and set them in the solidstate
*/
#ifndef DUMUX_BASE_SOLID_STATE_HH
#define DUMUX_BASE_SOLID_STATE_HH
#include <dumux/common/valgrind.hh>
#include <limits>
#ifndef DUMUX_UPDATE_SOLID_VOLUME_FRACTION_HH
#define DUMUX_UPDATE_SOLID_VOLUME_FRACTION_HH
namespace Dumux
{
......@@ -39,18 +33,19 @@ void updateSolidVolumeFractions(const ElemSol &elemSol,
const Element &element,
const Scv &scv,
SolidState& solidState,
const int numFluidComponents)
const int solidVolFracOffset)
{
for(int sCompIdx = solidState.numComponents- solidState.numInertComponents; sCompIdx < solidState.numComponents; ++sCompIdx)
{
solidState.setVolumeFraction(sCompIdx, problem.spatialParams().inertVolumeFraction(element, scv, solidState, sCompIdx));
const auto inertVolumeFraction = problem.spatialParams().inertVolumeFraction(element, scv, elemSol, solidState, sCompIdx);
solidState.setVolumeFraction(sCompIdx, inertVolumeFraction);
}
if (!(solidState.isInert()))
{
auto&& priVars = elemSol[scv.localDofIndex()];
for(int sCompIdx = 0; sCompIdx < solidState.numComponents- solidState.numInertComponents; ++sCompIdx)
{
solidState.setVolumeFraction(sCompIdx, priVars[numFluidComponents + sCompIdx]);
solidState.setVolumeFraction(sCompIdx, priVars[solidVolFracOffset + sCompIdx]);
}
}
}
......
......@@ -46,32 +46,23 @@ public:
****************************************/
static constexpr int numComponents = 2;
static constexpr int numInertComponents = isInert1 ? (isInert2 ? 2 : 1) : (isInert2 ? 1 : 0);
static constexpr int componentOneIdx = 0;
static constexpr int componentTwoIdx = 1;
static constexpr int comp0Idx = 0;
static constexpr int comp1Idx = 1;
/*!
* \brief Return the human readable name of a solid phase
*
* \param phaseIdx The index of the solid phase to consider
*/
static std::string phaseName(int phaseIdx)
{
switch (phaseIdx) {
case componentOneIdx: return ComponentOne::name();
case componentTwoIdx: return ComponentTwo::name();
}
DUNE_THROW(Dune::InvalidStateException, "Invalid phase index " << phaseIdx);
}
/*!
* \brief A human readable name for the component.
*
* \param compIdx The index of the component to consider
* \param compIdx The index of the solid phase to consider
*/
static std::string componentName(int compIdx)
{
DUNE_THROW(Dune::InvalidStateException, "ComponentName does not exist " << compIdx);
switch (compIdx)
{
case comp0Idx: return ComponentOne::name();
case comp1Idx: return ComponentTwo::name();
}
DUNE_THROW(Dune::InvalidStateException, "Invalid component index " << compIdx);
}
/*!
......@@ -83,7 +74,7 @@ public:
/*!
* \brief Returns whether the phase is incompressible
*/
static constexpr bool isCompressible(int phaseIdx)
static constexpr bool isCompressible(int compIdx)
{ return false; }
/*!
......@@ -100,13 +91,14 @@ public:
/*!
* \brief The molar mass in \f$\mathrm{[kg/mol]}\f$ of the component.
*/
static Scalar molarMass(int phaseIdx)
static Scalar molarMass(int compIdx)
{
switch (phaseIdx) {
case componentOneIdx: return ComponentOne::molarMass();
case componentTwoIdx: return ComponentTwo::molarMass();
switch (compIdx)
{
case comp0Idx: return ComponentOne::molarMass();
case comp1Idx: return ComponentTwo::molarMass();
}
DUNE_THROW(Dune::InvalidStateException, "Invalid phase index " << phaseIdx);
DUNE_THROW(Dune::InvalidStateException, "Invalid component index " << compIdx);
}
/*!
......@@ -117,8 +109,8 @@ public:
{
Scalar rho1 = ComponentOne::solidDensity(solidState.temperature());
Scalar rho2 = ComponentTwo::solidDensity(solidState.temperature());
Scalar volFrac1 = solidState.volumeFraction(componentOneIdx);
Scalar volFrac2 = solidState.volumeFraction(componentTwoIdx);
Scalar volFrac1 = solidState.volumeFraction(comp0Idx);
Scalar volFrac2 = solidState.volumeFraction(comp1Idx);
return (rho1*volFrac1+
rho2*volFrac2)/(volFrac1+volFrac2);
......@@ -128,28 +120,28 @@ public:
* \brief The density \f$\mathrm{[kg/m^3]}\f$ of the solid phase at a given pressure and temperature.
*/
template <class SolidState>
static Scalar density(const SolidState& solidState, const int phaseIdx)
static Scalar density(const SolidState& solidState, const int compIdx)
{
switch (phaseIdx)
switch (compIdx)
{
case componentOneIdx: return ComponentOne::solidDensity(solidState.temperature());
case componentTwoIdx: return ComponentTwo::solidDensity(solidState.temperature());
case comp0Idx: return ComponentOne::solidDensity(solidState.temperature());
case comp1Idx: return ComponentTwo::solidDensity(solidState.temperature());
}
DUNE_THROW(Dune::InvalidStateException, "Invalid phase index " << phaseIdx);
DUNE_THROW(Dune::InvalidStateException, "Invalid component index " << compIdx);
}
/*!
* \brief The molar density of the solid phase at a given pressure and temperature.
*/
template <class SolidState>
static Scalar molarDensity(const SolidState& solidState, const int phaseIdx)
static Scalar molarDensity(const SolidState& solidState, const int compIdx)
{
switch (phaseIdx)
switch (compIdx)
{
case componentOneIdx: return ComponentOne::solidDensity(solidState.temperature())/ComponentOne::molarMass();
case componentTwoIdx: return ComponentTwo::solidDensity(solidState.temperature())/ComponentTwo::molarMass();
case comp0Idx: return ComponentOne::solidDensity(solidState.temperature())/ComponentOne::molarMass();
case comp1Idx: return ComponentTwo::solidDensity(solidState.temperature())/ComponentTwo::molarMass();
}
DUNE_THROW(Dune::InvalidStateException, "Invalid phase index " << phaseIdx);
DUNE_THROW(Dune::InvalidStateException, "Invalid component index " << compIdx);
}
/*!
......@@ -160,8 +152,8 @@ public:
{
Scalar lambda1 = ComponentOne::solidThermalConductivity(solidState.temperature());
Scalar lambda2 = ComponentTwo::solidThermalConductivity(solidState.temperature());
Scalar volFrac1 = solidState.volumeFraction(componentOneIdx);
Scalar volFrac2 = solidState.volumeFraction(componentTwoIdx);
Scalar volFrac1 = solidState.volumeFraction(comp0Idx);
Scalar volFrac2 = solidState.volumeFraction(comp1Idx);
return (lambda1*volFrac1+
lambda2*volFrac2)/(volFrac1+volFrac2);
......@@ -175,8 +167,8 @@ public:
{
Scalar c1 = ComponentOne::solidHeatCapacity(solidState.temperature());
Scalar c2 = ComponentTwo::solidHeatCapacity(solidState.temperature());
Scalar volFrac1 = solidState.volumeFraction(componentOneIdx);
Scalar volFrac2 = solidState.volumeFraction(componentTwoIdx);
Scalar volFrac1 = solidState.volumeFraction(comp0Idx);
Scalar volFrac2 = solidState.volumeFraction(comp1Idx);
return (c1*volFrac1+
c2*volFrac2)/(volFrac1+volFrac2);
......
......@@ -45,16 +45,6 @@ public:
static constexpr int numComponents = 1;
static constexpr int numInertComponents = 1;
enum { lastInertCompIdx = numInertComponents-1 };
/*!
* \brief Return the human readable name of a solid phase
*
* \param phaseIdx The index of the solid phase to consider
*/
static std::string phaseName(int phaseIdx = 0)
{ return Component::name(); }
/*!
* \brief A human readable name for the component.
*
......@@ -72,7 +62,7 @@ public:
/*!
* \brief Returns whether the phase is incompressible
*/
static constexpr bool isCompressible(int phaseIdx = 0)
static constexpr bool isCompressible(int compIdx = 0)
{ return false; }
/*!
......
......@@ -173,11 +173,15 @@ 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
const SubControlVolume& scv,
const ElementSolution& elemSol) const
{
return asImp_().porosityAtPos(scv.center());
}
......@@ -194,17 +198,7 @@ public:
}
/*!
* \brief Function for defining the porosity.
*
* \return porosity
* \param globalPos The position of the center of the scv
*/
Scalar minimalPorosity(const Element& element,
const SubControlVolume& scv) const
{ return 0.0; }
/*!
* \brief Function for defining the porosity.
* \brief Function for defining the solid volume fraction.
* That is possibly solution dependent.
*
* \param element The current element
......@@ -212,24 +206,25 @@ public:
* \param elemSol The solution at the dofs connected to the element.
* \return the porosity
*/
template<class SolidState>
template<class ElementSolution, class SolidState>
Scalar inertVolumeFraction(const Element& element,
const SubControlVolume& scv,
SolidState& solidState,
const ElementSolution& elemSol,
const SolidState& solidState,
int compIdx) const
{
return asImp_().inertVolumeFractionAtPos(scv.center(), solidState, compIdx);
}
/*!
* \brief Function for defining the porosity.
* \brief Function for defining the solid volume fraction.
*
* \return porosity
* \return solid volume fraction
* \param globalPos The position of the center of the scv
*/
template<class SolidState>
Scalar inertVolumeFractionAtPos(const GlobalPosition& globalPos,
SolidState& solidState,
const SolidState& solidState,
int compIdx) const
{
if (solidState.isInert() && solidState.numInertComponents == 1)
......
......@@ -46,9 +46,6 @@
#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>
......@@ -152,24 +149,6 @@ 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
///////////////////////////////////////////////////////////////////////////
......
......@@ -50,7 +50,7 @@ class OnePVolumeVariables
using Scalar = typename Traits::PrimaryVariables::value_type;
using Indices = typename Traits::ModelTraits::Indices;
using PermeabilityType = typename Traits::PermeabilityType;
static constexpr int numComp = ParentType::numComponents();
static constexpr int numFluidComps = ParentType::numComponents();
public:
//! export the underlying fluid system
using FluidSystem = typename Traits::FluidSystem;
......@@ -81,7 +81,7 @@ public:
completeFluidState(elemSol, problem, element, scv, fluidState_, solidState_);
// porosity and permeability
updateSolidVolumeFractions(elemSol, problem, element, scv, solidState_, numComp);
updateSolidVolumeFractions(elemSol, problem, element, scv, solidState_, numFluidComps);
EnergyVolVars::updateSolidEnergyParams(elemSol, problem, element, scv, solidState_);
permeability_ = problem.spatialParams().permeability(element, scv, elemSol);
};
......@@ -188,7 +188,7 @@ public:
* \brief Return the average porosity \f$\mathrm{[-]}\f$ within the control volume.
*/
Scalar porosity() const
{ return solidState_.porosity();; }
{ return solidState_.porosity(); }
/*!
* \brief Returns the permeability within the control volume in \f$[m^2]\f$.
......
......@@ -61,9 +61,6 @@
#include <dumux/material/fluidmatrixinteractions/1p/thermalconductivityaverage.hh>
#include <dumux/material/fluidstates/compositional.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/1p/model.hh>
......@@ -164,25 +161,6 @@ public:
using type = CompositionalFluidState<Scalar, FluidSystem>;
};
//! The two-phase model uses the immiscible fluid state
SET_PROP(OnePNC, 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(OnePNC, SolidSystem)
{
using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar);
using InertComponent = Components::Constant<1, Scalar>;
using type = SolidSystems::InertSolidPhase<Scalar, InertComponent>;
};
//! Use the model after Millington (1961) for the effective diffusivity
SET_TYPE_PROP(OnePNC, EffectiveDiffusivityModel,
DiffusivityMillingtonQuirk<typename GET_PROP_TYPE(TypeTag, Scalar)>);
......
......@@ -44,14 +44,14 @@ namespace Dumux {
template <class Traits>
class OnePNCVolumeVariables
: public PorousMediumFlowVolumeVariables<Traits>
,public EnergyVolumeVariables<Traits, OnePNCVolumeVariables<Traits> >
, public EnergyVolumeVariables<Traits, OnePNCVolumeVariables<Traits> >
{
using ParentType = PorousMediumFlowVolumeVariables<Traits>;
using EnergyVolVars = EnergyVolumeVariables<Traits, OnePNCVolumeVariables<Traits> >;
using Scalar = typename Traits::PrimaryVariables::value_type;
using PermeabilityType = typename Traits::PermeabilityType;
using Idx = typename Traits::ModelTraits::Indices;
static constexpr int numComp = ParentType::numComponents();
static constexpr int numFluidComps = ParentType::numComponents();
enum
{
......@@ -95,12 +95,8 @@ public:
completeFluidState(elemSol, problem, element, scv, fluidState_, solidState_);
// calculate the remaining quantities
updateSolidVolumeFractions(elemSol, problem, element, scv, solidState_, numComp);
updateSolidVolumeFractions(elemSol, problem, element, scv, solidState_, numFluidComps);
EnergyVolVars::updateSolidEnergyParams(elemSol, problem, element, scv, solidState_);
Scalar minPorosity = problem.spatialParams().minimalPorosity(element, scv);
solidState_.setMinPorosity(minPorosity);
permeability_ = problem.spatialParams().permeability(element, scv, elemSol);
// Second instance of a parameter cache.
......@@ -110,7 +106,7 @@ public:
paramCache.updatePhase(fluidState_, fluidSystemPhaseIdx);
int compIIdx = mainCompMoleOrMassFracIdx;
for (unsigned int compJIdx = 0; compJIdx < numComp; ++compJIdx)
for (unsigned int compJIdx = 0; compJIdx < numFluidComps; ++compJIdx)
{
diffCoeff_[compJIdx] = 0.0;
if(compIIdx != compJIdx)
......@@ -149,10 +145,10 @@ public:
fluidState.setPressure(fluidSystemPhaseIdx, priVars[pressureIdx]);
// calculate the phase composition
Dune::FieldVector<Scalar, numComp> moleFrac;
Dune::FieldVector<Scalar, numFluidComps> moleFrac;
Scalar sumMoleFracNotMainComp = 0;
for (int compIdx = 0; compIdx < numComp; ++compIdx)
for (int compIdx = 0; compIdx < numFluidComps; ++compIdx)
{
if (compIdx != mainCompMoleOrMassFracIdx)
{
......@@ -163,7 +159,7 @@ public:
moleFrac[mainCompMoleOrMassFracIdx] = 1- sumMoleFracNotMainComp;
// Set fluid state mole fractions
for (int compIdx = 0; compIdx < numComp; ++compIdx)
for (int compIdx = 0; compIdx < numFluidComps; ++compIdx)
{
fluidState.setMoleFraction(fluidSystemPhaseIdx, compIdx, moleFrac[compIdx]);
}
......@@ -238,7 +234,7 @@ public:
Scalar moleFraction(int phaseIdx, int compIdx) const
{
// make sure this is only called with admissible indices
assert(compIdx < numComp);
assert(compIdx < numFluidComps