OilfieldToolsNet/opm-common
4
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

make the fluid states ready for the local-AD framework

Browse Source 
in fact they wouldn't have needed any modification, but returning
constant references instead of copies of the stored values saves quite
a few calls to copy constructors.

besides this, fluid states are now required to export the 'Scalar'
type, which allows to use it outside of the FluidState without
resorting to the c++ 'decltype' construct.
This commit is contained in:
Andreas Lauser 2015-05-21 15:33:14 +02:00
parent 0149901e4f
commit e009e640c5
14 changed files with 109 additions and 62 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

12
opm/material/fluidmatrixinteractions/EffToAbsLaw.hpp
View File

@ -112,7 +112,7 @@ public:
template <class Container, class FluidState>
static void capillaryPressures(Container &values, const Params &params, const FluidState &fs)
{
typedef Opm::SaturationOverlayFluidState<Scalar, FluidState> OverlayFluidState;
typedef Opm::SaturationOverlayFluidState<FluidState> OverlayFluidState;
OverlayFluidState overlayFs(fs);
for (int phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
@ -138,7 +138,7 @@ public:
template <class Container, class FluidState>
static void relativePermeabilities(Container &values, const Params &params, const FluidState &fs)
{
typedef Opm::SaturationOverlayFluidState<Scalar, FluidState> OverlayFluidState;
typedef Opm::SaturationOverlayFluidState<FluidState> OverlayFluidState;
OverlayFluidState overlayFs(fs);
for (int phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
@ -165,7 +165,7 @@ public:
template <class FluidState>
static Scalar pcnw(const Params &params, const FluidState &fs)
{
typedef Opm::SaturationOverlayFluidState<Scalar, FluidState> OverlayFluidState;
typedef Opm::SaturationOverlayFluidState<FluidState> OverlayFluidState;
static_assert(FluidState::numPhases == numPhases,
"The fluid state and the material law must exhibit the same "
@ -294,7 +294,7 @@ public:
template <class FluidState>
static Scalar krw(const Params &params, const FluidState &fs)
{
typedef Opm::SaturationOverlayFluidState<Scalar, FluidState> OverlayFluidState;
typedef Opm::SaturationOverlayFluidState<FluidState> OverlayFluidState;
static_assert(FluidState::numPhases == numPhases,
"The fluid state and the material law must exhibit the same "
@ -322,7 +322,7 @@ public:
template <class FluidState>
static Scalar krn(const Params &params, const FluidState &fs)
{
typedef Opm::SaturationOverlayFluidState<Scalar, FluidState> OverlayFluidState;
typedef Opm::SaturationOverlayFluidState<FluidState> OverlayFluidState;
static_assert(FluidState::numPhases == numPhases,
"The fluid state and the material law must exhibit the same "
@ -353,7 +353,7 @@ public:
static typename std::enable_if< (Traits::numPhases > 2), ScalarT>::type
krg(const Params &params, const FluidState &fs)
{
typedef Opm::SaturationOverlayFluidState<Scalar, FluidState> OverlayFluidState;
typedef Opm::SaturationOverlayFluidState<FluidState> OverlayFluidState;
static_assert(FluidState::numPhases == numPhases,
"The fluid state and the material law must exhibit the same "

24
opm/material/fluidstates/FluidStateCompositionModules.hpp
View File

@ -25,6 +25,7 @@
#define OPM_FLUID_STATE_COMPOSITION_MODULES_HPP
#include <opm/material/common/Valgrind.hpp>
#include <opm/material/common/MathToolbox.hpp>
#include <opm/material/common/ErrorMacros.hpp>
#include <opm/material/common/Exceptions.hpp>
@ -57,7 +58,7 @@ public:
/*!
* \brief The mole fraction of a component in a phase []
*/
Scalar moleFraction(int phaseIdx, int compIdx) const
const Scalar& moleFraction(int phaseIdx, int compIdx) const
{ return moleFraction_[phaseIdx][compIdx]; }
/*!
@ -65,11 +66,13 @@ public:
*/
Scalar massFraction(int phaseIdx, int compIdx) const
{
typedef Opm::MathToolbox<Scalar> Toolbox;
return
std::abs(sumMoleFractions_[phaseIdx])
* moleFraction_[phaseIdx][compIdx]
* FluidSystem::molarMass(compIdx)
/ std::max(1e-40, std::abs(averageMolarMass_[phaseIdx]));
Toolbox::abs(sumMoleFractions_[phaseIdx])
*moleFraction_[phaseIdx][compIdx]
*FluidSystem::molarMass(compIdx)
/ Toolbox::max(1e-40, Toolbox::abs(averageMolarMass_[phaseIdx]));
}
/*!
@ -80,7 +83,7 @@ public:
* component's molar masses weighted by the current mole fraction:
* \f[ \bar M_\alpha = \sum_\kappa M^\kappa x_\alpha^\kappa \f]
*/
Scalar averageMolarMass(int phaseIdx) const
const Scalar& averageMolarMass(int phaseIdx) const
{ return averageMolarMass_[phaseIdx]; }
/*!
@ -100,7 +103,7 @@ public:
* and update the average molar mass [kg/mol] according
* to the current composition of the phase
*/
void setMoleFraction(int phaseIdx, int compIdx, Scalar value)
void setMoleFraction(int phaseIdx, int compIdx, const Scalar& value)
{
Valgrind::CheckDefined(value);
Valgrind::SetUndefined(sumMoleFractions_[phaseIdx]);
@ -125,11 +128,16 @@ public:
template <class FluidState>
void assign(const FluidState& fs)
{
typedef typename FluidState::Scalar FsScalar;
typedef Opm::MathToolbox<FsScalar> FsToolbox;
for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
averageMolarMass_[phaseIdx] = 0;
sumMoleFractions_[phaseIdx] = 0;
for (int compIdx = 0; compIdx < numComponents; ++compIdx) {
moleFraction_[phaseIdx][compIdx] = fs.moleFraction(phaseIdx, compIdx);
moleFraction_[phaseIdx][compIdx] =
FsToolbox::template toLhs<Scalar>(fs.moleFraction(phaseIdx, compIdx));
averageMolarMass_[phaseIdx] += moleFraction_[phaseIdx][compIdx]*FluidSystem::molarMass(compIdx);
sumMoleFractions_[phaseIdx] += moleFraction_[phaseIdx][compIdx];
}

15
opm/material/fluidstates/FluidStateDensityModules.hpp
View File

@ -24,11 +24,12 @@
#ifndef OPM_FLUID_STATE_DENSITY_MODULES_HPP
#define OPM_FLUID_STATE_DENSITY_MODULES_HPP
#include <opm/material/common/Valgrind.hpp>
#include <opm/material/common/ErrorMacros.hpp>
#include <opm/material/common/Exceptions.hpp>
#include <opm/material/common/MathToolbox.hpp>
#include <opm/material/common/Valgrind.hpp>
#include <algorithm>
namespace Opm {
@ -51,7 +52,7 @@ public:
/*!
* \brief The density of a fluid phase [kg/m^3]
*/
Scalar density(int phaseIdx) const
const Scalar& density(int phaseIdx) const
{ return density_[phaseIdx]; }
/*!
@ -69,7 +70,7 @@ public:
/*!
* \brief Set the density of a phase [kg/m^3]
*/
void setDensity(int phaseIdx, Scalar value)
void setDensity(int phaseIdx, const Scalar& value)
{ density_[phaseIdx] = value; }
/*!
@ -79,8 +80,10 @@ public:
template <class FluidState>
void assign(const FluidState& fs)
{
typedef typename FluidState::Scalar FsScalar;
typedef Opm::MathToolbox<FsScalar> FsToolbox;
for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
density_[phaseIdx] = fs.density(phaseIdx);
density_[phaseIdx] = FsToolbox::template toLhs<Scalar>(fs.density(phaseIdx));
}
}
@ -120,7 +123,7 @@ public:
/*!
* \brief The density of a fluid phase [kg/m^3]
*/
Scalar density(int phaseIdx) const
const Scalar& density(int phaseIdx) const
{ OPM_THROW(std::logic_error, "Density is not provided by this fluid state"); }

35
opm/material/fluidstates/FluidStateEnthalpyModules.hpp
View File

@ -24,15 +24,15 @@
#ifndef OPM_FLUID_STATE_ENTHALPY_MODULES_HPP
#define OPM_FLUID_STATE_ENTHALPY_MODULES_HPP
#include <opm/material/common/Valgrind.hpp>
#include <opm/material/common/ErrorMacros.hpp>
#include <opm/material/common/Exceptions.hpp>
#include <opm/material/common/MathToolbox.hpp>
#include <opm/material/common/Valgrind.hpp>
#include <algorithm>
namespace Opm {
/*!
* \brief Module for the modular fluid state which stores the
* enthalpies explicitly.
@ -51,7 +51,7 @@ public:
/*!
* \brief The specific enthalpy of a fluid phase [J/kg]
*/
Scalar enthalpy(int phaseIdx) const
const Scalar& enthalpy(int phaseIdx) const
{ return enthalpy_[phaseIdx]; }
/*!
@ -63,7 +63,7 @@ public:
/*!
* \brief Set the specific enthalpy of a phase [J/kg]
*/
void setEnthalpy(int phaseIdx, Scalar value)
void setEnthalpy(int phaseIdx, const Scalar& value)
{ enthalpy_[phaseIdx] = value; }
/*!
@ -73,8 +73,10 @@ public:
template <class FluidState>
void assign(const FluidState& fs)
{
typedef typename FluidState::Scalar FsScalar;
typedef Opm::MathToolbox<FsScalar> FsToolbox;
for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
enthalpy_[phaseIdx] = fs.enthalpy(phaseIdx);
enthalpy_[phaseIdx] = FsToolbox::template toLhs<Scalar>(fs.enthalpy(phaseIdx));
}
}
@ -100,8 +102,9 @@ protected:
/*!
* \brief Module for the modular fluid state which does not store the
* enthalpies but throws std::logic_error instead.
*/
* enthalpies but returns 0 instead.
*
* Also, the returned values are marked as undefined in Valgrind... */
template <class Scalar,
class FluidSystem,
class Implementation>
@ -114,14 +117,22 @@ public:
/*!
* \brief The specific internal energy of a fluid phase [J/kg]
*/
Scalar internalEnergy(int phaseIdx) const
{ OPM_THROW(std::logic_error, "Internal energy is not provided by this fluid state"); }
const Scalar& internalEnergy(int phaseIdx) const
{
static Scalar tmp = 0;
Valgrind::SetUndefined(tmp);
return tmp;
}
/*!
* \brief The specific enthalpy of a fluid phase [J/kg]
*/
Scalar enthalpy(int phaseIdx) const
{ OPM_THROW(std::logic_error, "Enthalpy is not provided by this fluid state"); }
const Scalar& enthalpy(int phaseIdx) const
{
static Scalar tmp = 0;
Valgrind::SetUndefined(tmp);
return tmp;
}
/*!
* \brief Retrieve all parameters from an arbitrary fluid

10
opm/material/fluidstates/FluidStateFugacityModules.hpp
View File

@ -55,7 +55,7 @@ public:
/*!
* \brief The fugacity coefficient of a component in a phase []
*/
Scalar fugacityCoefficient(int phaseIdx, int compIdx) const
const Scalar& fugacityCoefficient(int phaseIdx, int compIdx) const
{ return fugacityCoefficient_[phaseIdx][compIdx]; }
/*!
@ -67,7 +67,7 @@ public:
/*!
* \brief Set the fugacity of a component in a phase []
*/
void setFugacityCoefficient(int phaseIdx, int compIdx, Scalar value)
void setFugacityCoefficient(int phaseIdx, int compIdx, const Scalar& value)
{ fugacityCoefficient_[phaseIdx][compIdx] = value; }
/*!
@ -137,7 +137,7 @@ public:
/*!
* \brief Set the fugacity of a component in a phase []
*/
void setFugacityCoefficient(int phaseIdx, Scalar value)
void setFugacityCoefficient(int phaseIdx, const Scalar& value)
{ fugacityCoefficient_[phaseIdx] = value; }
/*!
@ -188,13 +188,13 @@ public:
/*!
* \brief The fugacity coefficient of a component in a phase []
*/
Scalar fugacityCoefficient(int phaseIdx, int compIdx) const
const Scalar& fugacityCoefficient(int phaseIdx, int compIdx) const
{ OPM_THROW(std::logic_error, "Fugacity coefficients are not provided by this fluid state"); }
/*!
* \brief The fugacity of a component in a phase [Pa]
*/
Scalar fugacity(int phaseIdx, int compIdx) const
const Scalar& fugacity(int phaseIdx, int compIdx) const
{ OPM_THROW(std::logic_error, "Fugacities coefficients are not provided by this fluid state"); }
/*!

12
opm/material/fluidstates/FluidStatePressureModules.hpp
View File

@ -24,8 +24,8 @@
#ifndef OPM_FLUID_STATE_PRESSURE_MODULES_HPP
#define OPM_FLUID_STATE_PRESSURE_MODULES_HPP
#include <opm/material/common/MathToolbox.hpp>
#include <opm/material/common/Valgrind.hpp>
#include <opm/material/common/Exceptions.hpp>
#include <opm/material/common/ErrorMacros.hpp>
@ -51,14 +51,14 @@ public:
/*!
* \brief The pressure of a fluid phase [Pa]
*/
Scalar pressure(int phaseIdx) const
const Scalar& pressure(int phaseIdx) const
{ return pressure_[phaseIdx]; }
/*!
* \brief Set the pressure of a phase [Pa]
*/
void setPressure(int phaseIdx, Scalar value)
void setPressure(int phaseIdx, const Scalar& value)
{ pressure_[phaseIdx] = value; }
/*!
@ -68,8 +68,10 @@ public:
template <class FluidState>
void assign(const FluidState& fs)
{
typedef typename FluidState::Scalar FsScalar;
typedef Opm::MathToolbox<FsScalar> FsToolbox;
for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
pressure_[phaseIdx] = fs.pressure(phaseIdx);
pressure_[phaseIdx] = FsToolbox::template toLhs<Scalar>(fs.pressure(phaseIdx));
}
}
@ -106,7 +108,7 @@ public:
/*!
* \brief The pressure of a fluid phase [Pa]
*/
Scalar pressure(int phaseIdx) const
const Scalar& pressure(int phaseIdx) const
{ OPM_THROW(std::logic_error, "Pressure is not provided by this fluid state"); }

15
opm/material/fluidstates/FluidStateSaturationModules.hpp
View File

@ -24,11 +24,12 @@
#ifndef OPM_FLUID_STATE_SATURATION_MODULES_HPP
#define OPM_FLUID_STATE_SATURATION_MODULES_HPP
#include <opm/material/common/Valgrind.hpp>
#include <opm/material/common/ErrorMacros.hpp>
#include <opm/material/common/Exceptions.hpp>
#include <opm/material/common/MathToolbox.hpp>
#include <opm/material/common/Valgrind.hpp>
#include <algorithm>
namespace Opm {
@ -51,13 +52,13 @@ public:
/*!
* \brief The saturation of a fluid phase [-]
*/
Scalar saturation(int phaseIdx) const
const Scalar& saturation(int phaseIdx) const
{ return saturation_[phaseIdx]; }
/*!
* \brief Set the saturation of a phase [-]
*/
void setSaturation(int phaseIdx, Scalar value)
void setSaturation(int phaseIdx, const Scalar& value)
{ saturation_[phaseIdx] = value; }
/*!
@ -67,8 +68,10 @@ public:
template <class FluidState>
void assign(const FluidState& fs)
{
typedef typename FluidState::Scalar FsScalar;
typedef Opm::MathToolbox<FsScalar> FsToolbox;
for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
saturation_[phaseIdx] = fs.saturation(phaseIdx);
saturation_[phaseIdx] = FsToolbox::template toLhs<Scalar>(fs.saturation(phaseIdx));
}
}
@ -105,7 +108,7 @@ public:
/*!
* \brief The saturation of a fluid phase [-]
*/
Scalar saturation(int phaseIdx) const
const Scalar& saturation(int phaseIdx) const
{ OPM_THROW(std::runtime_error, "Saturation is not provided by this fluid state"); }
/*!

9
opm/material/fluidstates/FluidStateTemperatureModules.hpp
View File

@ -27,6 +27,7 @@
#include <opm/material/common/Valgrind.hpp>
#include <opm/material/common/MathToolbox.hpp>
#include <opm/material/common/ErrorMacros.hpp>
#include <opm/material/common/Exceptions.hpp>
@ -125,11 +126,15 @@ public:
template <class FluidState>
void assign(const FluidState& fs)
{
temperature_ = fs.temperature(/*phaseIdx=*/0);
typedef Opm::MathToolbox<Scalar> Toolbox;
typedef Opm::MathToolbox<typename FluidState::Scalar> FsToolbox;
temperature_ = FsToolbox::template toLhs<Scalar>(fs.temperature(/*phaseIdx=*/0));
#ifndef NDEBUG
for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
assert(std::abs(fs.temperature(phaseIdx) - temperature_) < 1e-30);
assert(std::abs(FsToolbox::value(fs.temperature(phaseIdx))
- Toolbox::value(temperature_)) < 1e-30);
}
#endif
}

10
opm/material/fluidstates/FluidStateViscosityModules.hpp
View File

@ -25,15 +25,15 @@
#ifndef OPM_FLUID_STATE_VISCOSITY_MODULES_HPP
#define OPM_FLUID_STATE_VISCOSITY_MODULES_HPP
#include <opm/material/common/Valgrind.hpp>
#include <opm/material/common/ErrorMacros.hpp>
#include <opm/material/common/Exceptions.hpp>
#include <opm/material/common/MathToolbox.hpp>
#include <opm/material/common/Valgrind.hpp>
#include <algorithm>
namespace Opm {
/*!
* \brief Module for the modular fluid state which stores the
* viscosities explicitly.
@ -68,8 +68,10 @@ public:
template <class FluidState>
void assign(const FluidState& fs)
{
typedef typename FluidState::Scalar FsScalar;
typedef Opm::MathToolbox<FsScalar> FsToolbox;
for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
viscosity_[phaseIdx] = fs.viscosity(phaseIdx);
viscosity_[phaseIdx] = FsToolbox::template toLhs<Scalar>(fs.viscosity(phaseIdx));
}
}

4
opm/material/fluidstates/PressureOverlayFluidState.hpp
View File

@ -34,10 +34,12 @@ namespace Opm {
* pressures and takes all other quantities from an other
* fluid state.
*/
template <class Scalar, class FluidState>
template <class FluidState>
class PressureOverlayFluidState
{
public:
typedef typename FluidState::Scalar Scalar;
enum { numPhases = FluidState::numPhases };
enum { numComponents = FluidState::numComponents };

4
opm/material/fluidstates/SaturationOverlayFluidState.hpp
View File

@ -34,10 +34,12 @@ namespace Opm {
* saturations and takes all other quantities from an other
* fluid state.
*/
template <class Scalar, class FluidState>
template <class FluidState>
class SaturationOverlayFluidState
{
public:
typedef typename FluidState::Scalar Scalar;
enum { numPhases = FluidState::numPhases };
enum { numComponents = FluidState::numComponents };

4
opm/material/fluidstates/TemperatureOverlayFluidState.hpp
View File

@ -32,10 +32,12 @@ namespace Opm {
* temperatures and takes all other quantities from an other
* fluid state.
*/
template <class Scalar, class FluidState>
template <class FluidState>
class TemperatureOverlayFluidState
{
public:
typedef typename FluidState::Scalar Scalar;
enum { numPhases = FluidState::numPhases };
enum { numComponents = FluidState::numComponents };

11
tests/checkFluidSystem.hpp
View File

@ -51,9 +51,9 @@
* \brief This is a fluid state which makes sure that only the quantities
* allowed are accessed.
*/
template <class Scalar,
template <class ScalarT,
class FluidSystem,
class BaseFluidState = Opm::CompositionalFluidState<Scalar, FluidSystem> >
class BaseFluidState = Opm::CompositionalFluidState<ScalarT, FluidSystem> >
class HairSplittingFluidState
: protected BaseFluidState
{
@ -61,6 +61,8 @@ public:
enum { numPhases = FluidSystem::numPhases };
enum { numComponents = FluidSystem::numComponents };
typedef ScalarT Scalar;
HairSplittingFluidState()
{
// initially, do not allow anything
@ -221,6 +223,11 @@ void checkFluidState(const BaseFluidState &fs)
// a fluid state must provide a checkDefined() method
fs.checkDefined();
// fluid states must export the types which they use as Scalars
typedef typename BaseFluidState::Scalar FsScalar;
static_assert(std::is_same<FsScalar, Scalar>::value,
"Fluid states must export the type they are given as scalar in an unmodified way");
// make sure the fluid state provides all mandatory methods
while (false) {
Scalar OPM_UNUSED val;

6
tests/test_fluidsystems.cpp
View File

@ -108,15 +108,15 @@ int main(int argc, char **argv)
BaseFluidState baseFs;
// TemperatureOverlayFluidState
{ Opm::TemperatureOverlayFluidState<Scalar, BaseFluidState> fs(baseFs);
{ Opm::TemperatureOverlayFluidState<BaseFluidState> fs(baseFs);
checkFluidState<Scalar>(fs); }
// PressureOverlayFluidState
{ Opm::PressureOverlayFluidState<Scalar, BaseFluidState> fs(baseFs);
{ Opm::PressureOverlayFluidState<BaseFluidState> fs(baseFs);
checkFluidState<Scalar>(fs); }
// SaturationOverlayFluidState
{ Opm::SaturationOverlayFluidState<Scalar, BaseFluidState> fs(baseFs);
{ Opm::SaturationOverlayFluidState<BaseFluidState> fs(baseFs);
checkFluidState<Scalar>(fs); }
}