changed: ewoms/models/flash -> opm/models/flash

This commit is contained in:
Arne Morten Kvarving 2019-09-16 11:58:51 +02:00
parent 16f4bdcf02
commit 72b5e42016
14 changed files with 1464 additions and 5 deletions

View File

@ -33,7 +33,7 @@
#endif
#include <opm/models/utils/start.hh>
#include <ewoms/models/flash/flashmodel.hh>
#include <opm/models/flash/flashmodel.hh>
#include <opm/models/discretization/ecfv/ecfvdiscretization.hh>
#include "problems/co2injectionflash.hh"
#include "problems/co2injectionproblem.hh"

View File

@ -30,7 +30,7 @@
#include <opm/material/common/quad.hpp>
#include <opm/models/utils/start.hh>
#include <ewoms/models/flash/flashmodel.hh>
#include <opm/models/flash/flashmodel.hh>
#include <opm/models/discretization/ecfv/ecfvdiscretization.hh>
#include "problems/co2injectionflash.hh"
#include "problems/co2injectionproblem.hh"

View File

@ -30,7 +30,7 @@
#include <opm/material/common/quad.hpp>
#include <opm/models/utils/start.hh>
#include <ewoms/models/flash/flashmodel.hh>
#include <opm/models/flash/flashmodel.hh>
#include <opm/models/discretization/vcfv/vcfvdiscretization.hh>
#include "problems/co2injectionflash.hh"
#include "problems/co2injectionproblem.hh"

View File

@ -33,7 +33,7 @@
#endif
#include <opm/models/utils/start.hh>
#include <ewoms/models/flash/flashmodel.hh>
#include <opm/models/flash/flashmodel.hh>
#include <opm/models/discretization/vcfv/vcfvdiscretization.hh>
#include "problems/co2injectionflash.hh"
#include "problems/co2injectionproblem.hh"

View File

@ -28,7 +28,7 @@
#include "config.h"
#include <opm/models/utils/start.hh>
#include <ewoms/models/flash/flashmodel.hh>
#include <opm/models/flash/flashmodel.hh>
#include "problems/diffusionproblem.hh"
BEGIN_PROPERTIES

View File

@ -0,0 +1,205 @@
// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
// vi: set et ts=4 sw=4 sts=4:
/*
This file is part of the Open Porous Media project (OPM).
OPM 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.
OPM 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 OPM. If not, see <http://www.gnu.org/licenses/>.
Consult the COPYING file in the top-level source directory of this
module for the precise wording of the license and the list of
copyright holders.
*/
/*!
* \file
*
* \copydoc Opm::FlashBoundaryRateVector
*/
#ifndef EWOMS_FLASH_BOUNDARY_RATE_VECTOR_HH
#define EWOMS_FLASH_BOUNDARY_RATE_VECTOR_HH
#include "flashproperties.hh"
#include <opm/models/common/energymodule.hh>
#include <opm/material/common/Valgrind.hpp>
namespace Opm {
/*!
* \ingroup FlashModel
*
* \brief Implements a boundary vector for the fully implicit
* compositional multi-phase model which is based on flash
* calculations.
*/
template <class TypeTag>
class FlashBoundaryRateVector : public GET_PROP_TYPE(TypeTag, RateVector)
{
typedef typename GET_PROP_TYPE(TypeTag, RateVector) ParentType;
typedef typename GET_PROP_TYPE(TypeTag, ExtensiveQuantities) ExtensiveQuantities;
typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
enum { numEq = GET_PROP_VALUE(TypeTag, NumEq) };
enum { numPhases = GET_PROP_VALUE(TypeTag, NumPhases) };
enum { numComponents = GET_PROP_VALUE(TypeTag, NumComponents) };
enum { conti0EqIdx = Indices::conti0EqIdx };
enum { enableEnergy = GET_PROP_VALUE(TypeTag, EnableEnergy) };
typedef Opm::EnergyModule<TypeTag, enableEnergy> EnergyModule;
typedef Opm::MathToolbox<Evaluation> Toolbox;
public:
FlashBoundaryRateVector() : ParentType()
{}
/*!
* \copydoc
* ImmiscibleBoundaryRateVector::ImmiscibleBoundaryRateVector(Scalar)
*/
FlashBoundaryRateVector(const Evaluation& value) : ParentType(value)
{}
/*!
* \copydoc ImmiscibleBoundaryRateVector::ImmiscibleBoundaryRateVector(const
* ImmiscibleBoundaryRateVector& )
*/
FlashBoundaryRateVector(const FlashBoundaryRateVector& value) = default;
FlashBoundaryRateVector& operator=(const FlashBoundaryRateVector& value) = default;
/*!
* \copydoc ImmiscibleBoundaryRateVector::setFreeFlow
*/
template <class Context, class FluidState>
void setFreeFlow(const Context& context,
unsigned bfIdx,
unsigned timeIdx,
const FluidState& fluidState)
{
ExtensiveQuantities extQuants;
extQuants.updateBoundary(context, bfIdx, timeIdx, fluidState);
const auto& insideIntQuants = context.intensiveQuantities(bfIdx, timeIdx);
unsigned focusDofIdx = context.focusDofIndex();
unsigned interiorDofIdx = context.interiorScvIndex(bfIdx, timeIdx);
////////
// advective fluxes of all components in all phases
////////
(*this) = Evaluation(0.0);
for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
Evaluation density;
if (fluidState.pressure(phaseIdx) > insideIntQuants.fluidState().pressure(phaseIdx)) {
if (focusDofIdx == interiorDofIdx)
density = fluidState.density(phaseIdx);
else
density = Opm::getValue(fluidState.density(phaseIdx));
}
else if (focusDofIdx == interiorDofIdx)
density = insideIntQuants.fluidState().density(phaseIdx);
else
density = Opm::getValue(insideIntQuants.fluidState().density(phaseIdx));
for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
Evaluation molarity;
if (fluidState.pressure(phaseIdx) > insideIntQuants.fluidState().pressure(phaseIdx)) {
if (focusDofIdx == interiorDofIdx)
molarity = fluidState.molarity(phaseIdx, compIdx);
else
molarity = Opm::getValue(fluidState.molarity(phaseIdx, compIdx));
}
else if (focusDofIdx == interiorDofIdx)
molarity = insideIntQuants.fluidState().molarity(phaseIdx, compIdx);
else
molarity = Opm::getValue(insideIntQuants.fluidState().molarity(phaseIdx, compIdx));
// add advective flux of current component in current
// phase
(*this)[conti0EqIdx + compIdx] += extQuants.volumeFlux(phaseIdx)*molarity;
}
if (enableEnergy) {
Evaluation specificEnthalpy;
if (fluidState.pressure(phaseIdx) > insideIntQuants.fluidState().pressure(phaseIdx)) {
if (focusDofIdx == interiorDofIdx)
specificEnthalpy = fluidState.enthalpy(phaseIdx);
else
specificEnthalpy = Opm::getValue(fluidState.enthalpy(phaseIdx));
}
else if (focusDofIdx == interiorDofIdx)
specificEnthalpy = insideIntQuants.fluidState().enthalpy(phaseIdx);
else
specificEnthalpy = Opm::getValue(insideIntQuants.fluidState().enthalpy(phaseIdx));
Evaluation enthalpyRate = density*extQuants.volumeFlux(phaseIdx)*specificEnthalpy;
EnergyModule::addToEnthalpyRate(*this, enthalpyRate);
}
}
// thermal conduction
EnergyModule::addToEnthalpyRate(*this, EnergyModule::thermalConductionRate(extQuants));
#ifndef NDEBUG
for (unsigned i = 0; i < numEq; ++i) {
Opm::Valgrind::CheckDefined((*this)[i]);
}
#endif
}
/*!
* \copydoc ImmiscibleBoundaryRateVector::setInFlow
*/
template <class Context, class FluidState>
void setInFlow(const Context& context,
unsigned bfIdx,
unsigned timeIdx,
const FluidState& fluidState)
{
this->setFreeFlow(context, bfIdx, timeIdx, fluidState);
// we only allow fluxes in the direction opposite to the outer unit normal
for (unsigned eqIdx = 0; eqIdx < numEq; ++eqIdx) {
Evaluation& val = this->operator[](eqIdx);
val = Toolbox::min(0.0, val);
}
}
/*!
* \copydoc ImmiscibleBoundaryRateVector::setOutFlow
*/
template <class Context, class FluidState>
void setOutFlow(const Context& context,
unsigned bfIdx,
unsigned timeIdx,
const FluidState& fluidState)
{
this->setFreeFlow(context, bfIdx, timeIdx, fluidState);
// we only allow fluxes in the same direction as the outer unit normal
for (unsigned eqIdx = 0; eqIdx < numEq; ++eqIdx) {
Evaluation& val = this->operator[](eqIdx);
val = Toolbox::max(0.0, val);
}
}
/*!
* \copydoc ImmiscibleBoundaryRateVector::setNoFlow
*/
void setNoFlow()
{ (*this) = Evaluation(0.0); }
};
} // namespace Opm
#endif

View File

@ -0,0 +1,96 @@
// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
// vi: set et ts=4 sw=4 sts=4:
/*
This file is part of the Open Porous Media project (OPM).
OPM 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.
OPM 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 OPM. If not, see <http://www.gnu.org/licenses/>.
Consult the COPYING file in the top-level source directory of this
module for the precise wording of the license and the list of
copyright holders.
*/
/*!
* \file
*
* \copydoc Opm::FlashExtensiveQuantities
*/
#ifndef EWOMS_FLASH_EXTENSIVE_QUANTITIES_HH
#define EWOMS_FLASH_EXTENSIVE_QUANTITIES_HH
#include "flashproperties.hh"
#include <opm/models/common/multiphasebaseextensivequantities.hh>
#include <opm/models/common/energymodule.hh>
#include <opm/models/common/diffusionmodule.hh>
namespace Opm {
/*!
* \ingroup FlashModel
* \ingroup ExtensiveQuantities
*
* \brief This template class contains the data which is required to
* calculate all fluxes of components over a face of a finite
* volume for the compositional multi-phase model based on
* flash calculations.
*
* This means pressure and concentration gradients, phase densities at
* the integration point, etc.
*/
template <class TypeTag>
class FlashExtensiveQuantities
: public MultiPhaseBaseExtensiveQuantities<TypeTag>
, public EnergyExtensiveQuantities<TypeTag, GET_PROP_VALUE(TypeTag, EnableEnergy)>
, public DiffusionExtensiveQuantities<TypeTag, GET_PROP_VALUE(TypeTag, EnableDiffusion)>
{
typedef MultiPhaseBaseExtensiveQuantities<TypeTag> ParentType;
typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
enum { enableDiffusion = GET_PROP_VALUE(TypeTag, EnableDiffusion) };
typedef Opm::DiffusionExtensiveQuantities<TypeTag, enableDiffusion> DiffusionExtensiveQuantities;
enum { enableEnergy = GET_PROP_VALUE(TypeTag, EnableEnergy) };
typedef Opm::EnergyExtensiveQuantities<TypeTag, enableEnergy> EnergyExtensiveQuantities;
public:
/*!
* \copydoc MultiPhaseBaseExtensiveQuantities::update
*/
void update(const ElementContext& elemCtx, unsigned scvfIdx, unsigned timeIdx)
{
ParentType::update(elemCtx, scvfIdx, timeIdx);
DiffusionExtensiveQuantities::update_(elemCtx, scvfIdx, timeIdx);
EnergyExtensiveQuantities::update_(elemCtx, scvfIdx, timeIdx);
}
/*!
* \copydoc MultiPhaseBaseExtensiveQuantities::updateBoundary
*/
template <class Context, class FluidState>
void updateBoundary(const Context& context,
unsigned bfIdx,
unsigned timeIdx,
const FluidState& fluidState)
{
ParentType::updateBoundary(context, bfIdx, timeIdx, fluidState);
DiffusionExtensiveQuantities::updateBoundary_(context, bfIdx, timeIdx, fluidState);
EnergyExtensiveQuantities::updateBoundary_(context, bfIdx, timeIdx, fluidState);
}
};
} // namespace Opm
#endif

View File

@ -0,0 +1,72 @@
// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
// vi: set et ts=4 sw=4 sts=4:
/*
This file is part of the Open Porous Media project (OPM).
OPM 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.
OPM 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 OPM. If not, see <http://www.gnu.org/licenses/>.
Consult the COPYING file in the top-level source directory of this
module for the precise wording of the license and the list of
copyright holders.
*/
/*!
* \file
*
* \copydoc Opm::FlashIndices
*/
#ifndef EWOMS_FLASH_INDICES_HH
#define EWOMS_FLASH_INDICES_HH
#include "flashproperties.hh"
#include <opm/models/common/energymodule.hh>
namespace Opm {
/*!
* \ingroup FlashModel
*
* \brief Defines the primary variable and equation indices for the
* compositional multi-phase model based on flash calculations.
*
* \tparam PVOffset The first index in a primary variable vector.
*/
template <class TypeTag, int PVOffset>
class FlashIndices
: public EnergyIndices<PVOffset + GET_PROP_VALUE(TypeTag, NumComponents),
GET_PROP_VALUE(TypeTag, EnableEnergy)>
{
enum { numComponents = GET_PROP_VALUE(TypeTag, NumComponents) };
enum { enableEnergy = GET_PROP_VALUE(TypeTag, EnableEnergy) };
typedef Opm::EnergyIndices<PVOffset + numComponents, enableEnergy> EnergyIndices;
public:
//! number of equations/primary variables
static const int numEq = numComponents + EnergyIndices::numEq_;
// Primary variable indices
//! Index of the total concentration of the first component in the
//! pore space.
static const int cTot0Idx = PVOffset;
// equation indices
//! Index of the mass conservation equation for the first
//! component.
static const int conti0EqIdx = PVOffset;
};
} // namespace Opm
#endif

View File

@ -0,0 +1,217 @@
// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
// vi: set et ts=4 sw=4 sts=4:
/*
This file is part of the Open Porous Media project (OPM).
OPM 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.
OPM 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 OPM. If not, see <http://www.gnu.org/licenses/>.
Consult the COPYING file in the top-level source directory of this
module for the precise wording of the license and the list of
copyright holders.
*/
/*!
* \file
*
* \copydoc Opm::FlashIntensiveQuantities
*/
#ifndef EWOMS_FLASH_INTENSIVE_QUANTITIES_HH
#define EWOMS_FLASH_INTENSIVE_QUANTITIES_HH
#include "flashproperties.hh"
#include "flashindices.hh"
#include <opm/models/common/energymodule.hh>
#include <opm/models/common/diffusionmodule.hh>
#include <opm/material/fluidstates/CompositionalFluidState.hpp>
#include <opm/material/common/Valgrind.hpp>
#include <dune/common/fvector.hh>
#include <dune/common/fmatrix.hh>
namespace Opm {
/*!
* \ingroup FlashModel
* \ingroup IntensiveQuantities
*
* \brief Contains the intensive quantities of the flash-based compositional multi-phase model
*/
template <class TypeTag>
class FlashIntensiveQuantities
: public GET_PROP_TYPE(TypeTag, DiscIntensiveQuantities)
, public DiffusionIntensiveQuantities<TypeTag, GET_PROP_VALUE(TypeTag, EnableDiffusion) >
, public EnergyIntensiveQuantities<TypeTag, GET_PROP_VALUE(TypeTag, EnableEnergy) >
, public GET_PROP_TYPE(TypeTag, FluxModule)::FluxIntensiveQuantities
{
typedef typename GET_PROP_TYPE(TypeTag, DiscIntensiveQuantities) ParentType;
typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
typedef typename GET_PROP_TYPE(TypeTag, MaterialLaw) MaterialLaw;
typedef typename GET_PROP_TYPE(TypeTag, MaterialLawParams) MaterialLawParams;
typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
typedef typename GET_PROP_TYPE(TypeTag, FluxModule) FluxModule;
typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
typedef typename GET_PROP_TYPE(TypeTag, ThreadManager) ThreadManager;
// primary variable indices
enum { cTot0Idx = Indices::cTot0Idx };
enum { numPhases = GET_PROP_VALUE(TypeTag, NumPhases) };
enum { numComponents = GET_PROP_VALUE(TypeTag, NumComponents) };
enum { enableDiffusion = GET_PROP_VALUE(TypeTag, EnableDiffusion) };
enum { enableEnergy = GET_PROP_VALUE(TypeTag, EnableEnergy) };
enum { dimWorld = GridView::dimensionworld };
typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
typedef typename GET_PROP_TYPE(TypeTag, FlashSolver) FlashSolver;
typedef Dune::FieldVector<Evaluation, numComponents> ComponentVector;
typedef Dune::FieldMatrix<Scalar, dimWorld, dimWorld> DimMatrix;
typedef typename FluxModule::FluxIntensiveQuantities FluxIntensiveQuantities;
typedef Opm::DiffusionIntensiveQuantities<TypeTag, enableDiffusion> DiffusionIntensiveQuantities;
typedef Opm::EnergyIntensiveQuantities<TypeTag, enableEnergy> EnergyIntensiveQuantities;
public:
//! The type of the object returned by the fluidState() method
typedef Opm::CompositionalFluidState<Evaluation, FluidSystem, enableEnergy> FluidState;
FlashIntensiveQuantities()
{ }
FlashIntensiveQuantities(const FlashIntensiveQuantities& other) = default;
FlashIntensiveQuantities& operator=(const FlashIntensiveQuantities& other) = default;
/*!
* \copydoc IntensiveQuantities::update
*/
void update(const ElementContext& elemCtx, unsigned dofIdx, unsigned timeIdx)
{
ParentType::update(elemCtx, dofIdx, timeIdx);
EnergyIntensiveQuantities::updateTemperatures_(fluidState_, elemCtx, dofIdx, timeIdx);
const auto& priVars = elemCtx.primaryVars(dofIdx, timeIdx);
const auto& problem = elemCtx.problem();
Scalar flashTolerance = EWOMS_GET_PARAM(TypeTag, Scalar, FlashTolerance);
// extract the total molar densities of the components
ComponentVector cTotal;
for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx)
cTotal[compIdx] = priVars.makeEvaluation(cTot0Idx + compIdx, timeIdx);
const auto *hint = elemCtx.thermodynamicHint(dofIdx, timeIdx);
if (hint) {
// use the same fluid state as the one of the hint, but
// make sure that we don't overwrite the temperature
// specified by the primary variables
Evaluation T = fluidState_.temperature(/*phaseIdx=*/0);
fluidState_.assign(hint->fluidState());
fluidState_.setTemperature(T);
}
else
FlashSolver::guessInitial(fluidState_, cTotal);
// compute the phase compositions, densities and pressures
typename FluidSystem::template ParameterCache<Evaluation> paramCache;
const MaterialLawParams& materialParams =
problem.materialLawParams(elemCtx, dofIdx, timeIdx);
FlashSolver::template solve<MaterialLaw>(fluidState_,
materialParams,
paramCache,
cTotal,
flashTolerance);
// calculate relative permeabilities
MaterialLaw::relativePermeabilities(relativePermeability_,
materialParams, fluidState_);
Opm::Valgrind::CheckDefined(relativePermeability_);
// set the phase viscosities
for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
paramCache.updatePhase(fluidState_, phaseIdx);
const Evaluation& mu = FluidSystem::viscosity(fluidState_, paramCache, phaseIdx);
fluidState_.setViscosity(phaseIdx, mu);
mobility_[phaseIdx] = relativePermeability_[phaseIdx] / mu;
Opm::Valgrind::CheckDefined(mobility_[phaseIdx]);
}
/////////////
// calculate the remaining quantities
/////////////
// porosity
porosity_ = problem.porosity(elemCtx, dofIdx, timeIdx);
Opm::Valgrind::CheckDefined(porosity_);
// intrinsic permeability
intrinsicPerm_ = problem.intrinsicPermeability(elemCtx, dofIdx, timeIdx);
// update the quantities specific for the velocity model
FluxIntensiveQuantities::update_(elemCtx, dofIdx, timeIdx);
// energy related quantities
EnergyIntensiveQuantities::update_(fluidState_, paramCache, elemCtx, dofIdx, timeIdx);
// update the diffusion specific quantities of the intensive quantities
DiffusionIntensiveQuantities::update_(fluidState_, paramCache, elemCtx, dofIdx, timeIdx);
}
/*!
* \copydoc ImmiscibleIntensiveQuantities::fluidState
*/
const FluidState& fluidState() const
{ return fluidState_; }
/*!
* \copydoc ImmiscibleIntensiveQuantities::intrinsicPermeability
*/
const DimMatrix& intrinsicPermeability() const
{ return intrinsicPerm_; }
/*!
* \copydoc ImmiscibleIntensiveQuantities::relativePermeability
*/
const Evaluation& relativePermeability(unsigned phaseIdx) const
{ return relativePermeability_[phaseIdx]; }
/*!
* \copydoc ImmiscibleIntensiveQuantities::mobility
*/
const Evaluation& mobility(unsigned phaseIdx) const
{
return mobility_[phaseIdx];
}
/*!
* \copydoc ImmiscibleIntensiveQuantities::porosity
*/
const Evaluation& porosity() const
{ return porosity_; }
private:
DimMatrix intrinsicPerm_;
FluidState fluidState_;
Evaluation porosity_;
Evaluation relativePermeability_[numPhases];
Evaluation mobility_[numPhases];
};
} // namespace Opm
#endif

View File

@ -0,0 +1,198 @@
// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
// vi: set et ts=4 sw=4 sts=4:
/*
This file is part of the Open Porous Media project (OPM).
OPM 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.
OPM 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 OPM. If not, see <http://www.gnu.org/licenses/>.
Consult the COPYING file in the top-level source directory of this
module for the precise wording of the license and the list of
copyright holders.
*/
/*!
* \file
*
* \copydoc Opm::FlashLocalResidual
*/
#ifndef EWOMS_FLASH_LOCAL_RESIDUAL_HH
#define EWOMS_FLASH_LOCAL_RESIDUAL_HH
#include "flashproperties.hh"
#include <opm/models/common/diffusionmodule.hh>
#include <opm/models/common/energymodule.hh>
#include <opm/material/common/Valgrind.hpp>
namespace Opm {
/*!
* \ingroup FlashModel
*
* \brief Calculates the local residual of the compositional multi-phase
* model based on flash calculations.
*/
template <class TypeTag>
class FlashLocalResidual: public GET_PROP_TYPE(TypeTag, DiscLocalResidual)
{
typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
typedef typename GET_PROP_TYPE(TypeTag, EqVector) EqVector;
typedef typename GET_PROP_TYPE(TypeTag, RateVector) RateVector;
typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
typedef typename GET_PROP_TYPE(TypeTag, IntensiveQuantities) IntensiveQuantities;
typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
enum { numEq = GET_PROP_VALUE(TypeTag, NumEq) };
enum { numPhases = GET_PROP_VALUE(TypeTag, NumPhases) };
enum { numComponents = GET_PROP_VALUE(TypeTag, NumComponents) };
enum { conti0EqIdx = Indices::conti0EqIdx };
enum { enableDiffusion = GET_PROP_VALUE(TypeTag, EnableDiffusion) };
typedef Opm::DiffusionModule<TypeTag, enableDiffusion> DiffusionModule;
enum { enableEnergy = GET_PROP_VALUE(TypeTag, EnableEnergy) };
typedef Opm::EnergyModule<TypeTag, enableEnergy> EnergyModule;
typedef Opm::MathToolbox<Evaluation> Toolbox;
public:
/*!
* \copydoc ImmiscibleLocalResidual::addPhaseStorage
*/
template <class LhsEval>
void addPhaseStorage(Dune::FieldVector<LhsEval, numEq>& storage,
const ElementContext& elemCtx,
unsigned dofIdx,
unsigned timeIdx,
unsigned phaseIdx) const
{
const IntensiveQuantities& intQuants = elemCtx.intensiveQuantities(dofIdx, timeIdx);
const auto& fs = intQuants.fluidState();
// compute storage term of all components within all phases
for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
unsigned eqIdx = conti0EqIdx + compIdx;
storage[eqIdx] +=
Toolbox::template decay<LhsEval>(fs.molarity(phaseIdx, compIdx))
* Toolbox::template decay<LhsEval>(fs.saturation(phaseIdx))
* Toolbox::template decay<LhsEval>(intQuants.porosity());
}
EnergyModule::addPhaseStorage(storage, elemCtx.intensiveQuantities(dofIdx, timeIdx), phaseIdx);
}
/*!
* \copydoc FvBaseLocalResidual::computeStorage
*/
template <class LhsEval>
void computeStorage(Dune::FieldVector<LhsEval, numEq>& storage,
const ElementContext& elemCtx,
unsigned dofIdx,
unsigned timeIdx) const
{
storage = 0;
for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
addPhaseStorage(storage, elemCtx, dofIdx, timeIdx, phaseIdx);
EnergyModule::addSolidEnergyStorage(storage, elemCtx.intensiveQuantities(dofIdx, timeIdx));
}
/*!
* \copydoc FvBaseLocalResidual::computeFlux
*/
void computeFlux(RateVector& flux,
const ElementContext& elemCtx,
unsigned scvfIdx,
unsigned timeIdx) const
{
flux = 0.0;
addAdvectiveFlux(flux, elemCtx, scvfIdx, timeIdx);
Opm::Valgrind::CheckDefined(flux);
addDiffusiveFlux(flux, elemCtx, scvfIdx, timeIdx);
Opm::Valgrind::CheckDefined(flux);
}
/*!
* \copydoc ImmiscibleLocalResidual::addAdvectiveFlux
*/
void addAdvectiveFlux(RateVector& flux,
const ElementContext& elemCtx,
unsigned scvfIdx,
unsigned timeIdx) const
{
const auto& extQuants = elemCtx.extensiveQuantities(scvfIdx, timeIdx);
unsigned focusDofIdx = elemCtx.focusDofIndex();
for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
// data attached to upstream and the finite volume of the current phase
unsigned upIdx = static_cast<unsigned>(extQuants.upstreamIndex(phaseIdx));
const IntensiveQuantities& up = elemCtx.intensiveQuantities(upIdx, timeIdx);
// this is a bit hacky because it is specific to the element-centered
// finite volume scheme. (N.B. that if finite differences are used to
// linearize the system of equations, it does not matter.)
if (upIdx == focusDofIdx) {
Evaluation tmp =
up.fluidState().molarDensity(phaseIdx)
* extQuants.volumeFlux(phaseIdx);
for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
flux[conti0EqIdx + compIdx] +=
tmp*up.fluidState().moleFraction(phaseIdx, compIdx);
}
}
else {
Evaluation tmp =
Toolbox::value(up.fluidState().molarDensity(phaseIdx))
* extQuants.volumeFlux(phaseIdx);
for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
flux[conti0EqIdx + compIdx] +=
tmp*Toolbox::value(up.fluidState().moleFraction(phaseIdx, compIdx));
}
}
}
EnergyModule::addAdvectiveFlux(flux, elemCtx, scvfIdx, timeIdx);
}
/*!
* \copydoc ImmiscibleLocalResidual::addDiffusiveFlux
*/
void addDiffusiveFlux(RateVector& flux,
const ElementContext& elemCtx,
unsigned scvfIdx,
unsigned timeIdx) const
{
DiffusionModule::addDiffusiveFlux(flux, elemCtx, scvfIdx, timeIdx);
EnergyModule::addDiffusiveFlux(flux, elemCtx, scvfIdx, timeIdx);
}
/*!
* \copydoc ImmiscibleLocalResidual::computeSource
*/
void computeSource(RateVector& source,
const ElementContext& elemCtx,
unsigned dofIdx,
unsigned timeIdx) const
{
Opm::Valgrind::SetUndefined(source);
elemCtx.problem().source(source, elemCtx, dofIdx, timeIdx);
Opm::Valgrind::CheckDefined(source);
}
};
} // namespace Opm
#endif

View File

@ -0,0 +1,318 @@
// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
// vi: set et ts=4 sw=4 sts=4:
/*
This file is part of the Open Porous Media project (OPM).
OPM 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.
OPM 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 OPM. If not, see <http://www.gnu.org/licenses/>.
Consult the COPYING file in the top-level source directory of this
module for the precise wording of the license and the list of
copyright holders.
*/
/*!
* \file
*
* \copydoc Opm::FlashModel
*/
#ifndef EWOMS_FLASH_MODEL_HH
#define EWOMS_FLASH_MODEL_HH
#include <opm/material/densead/Math.hpp>
#include "flashproperties.hh"
#include "flashprimaryvariables.hh"
#include "flashlocalresidual.hh"
#include "flashratevector.hh"
#include "flashboundaryratevector.hh"
#include "flashintensivequantities.hh"
#include "flashextensivequantities.hh"
#include "flashindices.hh"
#include <opm/models/common/multiphasebasemodel.hh>
#include <opm/models/common/energymodule.hh>
#include <ewoms/io/vtkcompositionmodule.hh>
#include <ewoms/io/vtkenergymodule.hh>
#include <ewoms/io/vtkdiffusionmodule.hh>
#include <opm/material/fluidmatrixinteractions/NullMaterial.hpp>
#include <opm/material/fluidmatrixinteractions/MaterialTraits.hpp>
#include <opm/material/constraintsolvers/NcpFlash.hpp>
#include <sstream>
#include <string>
namespace Opm {
template <class TypeTag>
class FlashModel;
}
BEGIN_PROPERTIES
//! The type tag for the isothermal single phase problems
NEW_TYPE_TAG(FlashModel, INHERITS_FROM(MultiPhaseBaseModel,
VtkComposition,
VtkEnergy,
VtkDiffusion));
//! Use the FlashLocalResidual function for the flash model
SET_TYPE_PROP(FlashModel, LocalResidual,
Opm::FlashLocalResidual<TypeTag>);
//! Use the NCP flash solver by default
SET_TYPE_PROP(FlashModel, FlashSolver,
Opm::NcpFlash<typename GET_PROP_TYPE(TypeTag, Scalar),
typename GET_PROP_TYPE(TypeTag, FluidSystem)>);
//! Let the flash solver choose its tolerance by default
SET_SCALAR_PROP(FlashModel, FlashTolerance, -1.0);
//! the Model property
SET_TYPE_PROP(FlashModel, Model, Opm::FlashModel<TypeTag>);
//! the PrimaryVariables property
SET_TYPE_PROP(FlashModel, PrimaryVariables, Opm::FlashPrimaryVariables<TypeTag>);
//! the RateVector property
SET_TYPE_PROP(FlashModel, RateVector, Opm::FlashRateVector<TypeTag>);
//! the BoundaryRateVector property
SET_TYPE_PROP(FlashModel, BoundaryRateVector, Opm::FlashBoundaryRateVector<TypeTag>);
//! the IntensiveQuantities property
SET_TYPE_PROP(FlashModel, IntensiveQuantities, Opm::FlashIntensiveQuantities<TypeTag>);
//! the ExtensiveQuantities property
SET_TYPE_PROP(FlashModel, ExtensiveQuantities, Opm::FlashExtensiveQuantities<TypeTag>);
//! The indices required by the flash-baseed isothermal compositional model
SET_TYPE_PROP(FlashModel, Indices, Opm::FlashIndices<TypeTag, /*PVIdx=*/0>);
// The updates of intensive quantities tend to be _very_ expensive for this
// model, so let's try to minimize the number of required ones
SET_BOOL_PROP(FlashModel, EnableIntensiveQuantityCache, true);
// since thermodynamic hints are basically free if the cache for intensive quantities is
// enabled, and this model usually shows quite a performance improvment if they are
// enabled, let's enable them by default.
SET_BOOL_PROP(FlashModel, EnableThermodynamicHints, true);
// disable molecular diffusion by default
SET_BOOL_PROP(FlashModel, EnableDiffusion, false);
//! Disable the energy equation by default
SET_BOOL_PROP(FlashModel, EnableEnergy, false);
END_PROPERTIES
namespace Opm {
/*!
* \ingroup FlashModel
*
* \brief A compositional multi-phase model based on flash-calculations
*
* This model assumes a flow of \f$M \geq 1\f$ fluid phases
* \f$\alpha\f$, each of which is assumed to be a mixture \f$N \geq
* M\f$ chemical species (denoted by the upper index \f$\kappa\f$).
*
* By default, the standard multi-phase Darcy approach is used to determine
* the velocity, i.e.
* \f[
* \mathbf{v}_\alpha =
* - \frac{k_{r\alpha}}{\mu_\alpha} \mathbf{K}
* \left(\mathbf{grad}\, p_\alpha
* - \varrho_{\alpha} \mathbf{g} \right) \;,
* \f]
* although the actual approach which is used can be specified via the
* \c FluxModule property. For example, the velocity model can by
* changed to the Forchheimer approach by
* \code
* SET_TYPE_PROP(MyProblemTypeTag, FluxModule, Opm::ForchheimerFluxModule<TypeTag>);
* \endcode
*
* The core of the model is the conservation mass of each component by
* means of the equation
* \f[
* \sum_\alpha \frac{\partial\;\phi c_\alpha^\kappa S_\alpha }{\partial t}
* - \sum_\alpha \mathrm{div} \left\{ c_\alpha^\kappa \mathbf{v}_\alpha \right\}
* - q^\kappa = 0 \;.
* \f]
*
* To determine the quanties that occur in the equations above, this
* model uses <i>flash calculations</i>. A flash solver starts with
* the total mass or molar mass per volume for each component and,
* calculates the compositions, saturations and pressures of all
* phases at a given temperature. For this the flash solver has to use
* some model assumptions internally. (Often these are the same
* primary variable switching or NCP assumptions as used by the other
* fully implicit compositional multi-phase models provided by eWoms.)
*
* Using flash calculations for the flow model has some disadvantages:
* - The accuracy of the flash solver needs to be sufficient to
* calculate the parital derivatives using numerical differentiation
* which are required for the Newton scheme.
* - Flash calculations tend to be quite computationally expensive and
* are often numerically unstable.
*
* It is thus adviced to increase the target tolerance of the Newton
* scheme or a to use type for scalar values which exhibits higher
* precision than the standard \c double (e.g. \c quad) if this model
* ought to be used.
*
* The model uses the following primary variables:
* - The total molar concentration of each component:
* \f$c^\kappa = \sum_\alpha S_\alpha x_\alpha^\kappa \rho_{mol, \alpha}\f$
* - The absolute temperature $T$ in Kelvins if the energy equation enabled.
*/
template <class TypeTag>
class FlashModel
: public MultiPhaseBaseModel<TypeTag>
{
typedef MultiPhaseBaseModel<TypeTag> ParentType;
typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
typedef typename GET_PROP_TYPE(TypeTag, Simulator) Simulator;
typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
enum { numComponents = GET_PROP_VALUE(TypeTag, NumComponents) };
enum { enableDiffusion = GET_PROP_VALUE(TypeTag, EnableDiffusion) };
enum { enableEnergy = GET_PROP_VALUE(TypeTag, EnableEnergy) };
typedef Opm::EnergyModule<TypeTag, enableEnergy> EnergyModule;
public:
FlashModel(Simulator& simulator)
: ParentType(simulator)
{}
/*!
* \brief Register all run-time parameters for the immiscible model.
*/
static void registerParameters()
{
ParentType::registerParameters();
// register runtime parameters of the VTK output modules
Opm::VtkCompositionModule<TypeTag>::registerParameters();
if (enableDiffusion)
Opm::VtkDiffusionModule<TypeTag>::registerParameters();
if (enableEnergy)
Opm::VtkEnergyModule<TypeTag>::registerParameters();
EWOMS_REGISTER_PARAM(TypeTag, Scalar, FlashTolerance,
"The maximum tolerance for the flash solver to "
"consider the solution converged");
}
/*!
* \copydoc FvBaseDiscretization::name
*/
static std::string name()
{ return "flash"; }
/*!
* \copydoc FvBaseDiscretization::primaryVarName
*/
std::string primaryVarName(unsigned pvIdx) const
{
const std::string& tmp = EnergyModule::primaryVarName(pvIdx);
if (tmp != "")
return tmp;
std::ostringstream oss;
if (Indices::cTot0Idx <= pvIdx && pvIdx < Indices::cTot0Idx
+ numComponents)
oss << "c_tot," << FluidSystem::componentName(/*compIdx=*/pvIdx
- Indices::cTot0Idx);
else
assert(false);
return oss.str();
}
/*!
* \copydoc FvBaseDiscretization::eqName
*/
std::string eqName(unsigned eqIdx) const
{
const std::string& tmp = EnergyModule::eqName(eqIdx);
if (tmp != "")
return tmp;
std::ostringstream oss;
if (Indices::conti0EqIdx <= eqIdx && eqIdx < Indices::conti0EqIdx
+ numComponents) {
unsigned compIdx = eqIdx - Indices::conti0EqIdx;
oss << "continuity^" << FluidSystem::componentName(compIdx);
}
else
assert(false);
return oss.str();
}
/*!
* \copydoc FvBaseDiscretization::primaryVarWeight
*/
Scalar primaryVarWeight(unsigned globalDofIdx, unsigned pvIdx) const
{
Scalar tmp = EnergyModule::primaryVarWeight(*this, globalDofIdx, pvIdx);
if (tmp > 0)
return tmp;
unsigned compIdx = pvIdx - Indices::cTot0Idx;
// make all kg equal. also, divide the weight of all total
// compositions by 100 to make the relative errors more
// comparable to the ones of the other models (at 10% porosity
// the medium is fully saturated with water at atmospheric
// conditions if 100 kg/m^3 are present!)
return FluidSystem::molarMass(compIdx) / 100.0;
}
/*!
* \copydoc FvBaseDiscretization::eqWeight
*/
Scalar eqWeight(unsigned globalDofIdx, unsigned eqIdx) const
{
Scalar tmp = EnergyModule::eqWeight(*this, globalDofIdx, eqIdx);
if (tmp > 0)
return tmp;
unsigned compIdx = eqIdx - Indices::conti0EqIdx;
// make all kg equal
return FluidSystem::molarMass(compIdx);
}
void registerOutputModules_()
{
ParentType::registerOutputModules_();
// add the VTK output modules which are meaningful for the model
this->addOutputModule(new Opm::VtkCompositionModule<TypeTag>(this->simulator_));
if (enableDiffusion)
this->addOutputModule(new Opm::VtkDiffusionModule<TypeTag>(this->simulator_));
if (enableEnergy)
this->addOutputModule(new Opm::VtkEnergyModule<TypeTag>(this->simulator_));
}
};
} // namespace Opm
#endif

View File

@ -0,0 +1,150 @@
// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
// vi: set et ts=4 sw=4 sts=4:
/*
This file is part of the Open Porous Media project (OPM).
OPM 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.
OPM 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 OPM. If not, see <http://www.gnu.org/licenses/>.
Consult the COPYING file in the top-level source directory of this
module for the precise wording of the license and the list of
copyright holders.
*/
/*!
* \file
*
* \copydoc Opm::FlashPrimaryVariables
*/
#ifndef EWOMS_FLASH_PRIMARY_VARIABLES_HH
#define EWOMS_FLASH_PRIMARY_VARIABLES_HH
#include "flashindices.hh"
#include "flashproperties.hh"
#include <opm/models/discretization/common/fvbaseprimaryvariables.hh>
#include <opm/models/common/energymodule.hh>
#include <opm/material/constraintsolvers/NcpFlash.hpp>
#include <opm/material/fluidstates/CompositionalFluidState.hpp>
#include <opm/material/common/Valgrind.hpp>
#include <opm/material/common/Unused.hpp>
#include <dune/common/fvector.hh>
#include <iostream>
namespace Opm {
/*!
* \ingroup FlashModel
*
* \brief Represents the primary variables used by the compositional
* flow model based on flash calculations.
*
* This class is basically a Dune::FieldVector which can retrieve its
* contents from an aribitatry fluid state.
*/
template <class TypeTag>
class FlashPrimaryVariables : public FvBasePrimaryVariables<TypeTag>
{
typedef FvBasePrimaryVariables<TypeTag> ParentType;
typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
typedef typename GET_PROP_TYPE(TypeTag, MaterialLawParams) MaterialLawParams;
typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
// primary variable indices
enum { cTot0Idx = Indices::cTot0Idx };
enum { numPhases = GET_PROP_VALUE(TypeTag, NumPhases) };
enum { numComponents = GET_PROP_VALUE(TypeTag, NumComponents) };
typedef typename Opm::MathToolbox<Evaluation> Toolbox;
typedef Opm::EnergyModule<TypeTag, GET_PROP_VALUE(TypeTag, EnableEnergy)> EnergyModule;
public:
FlashPrimaryVariables() : ParentType()
{ Opm::Valgrind::SetDefined(*this); }
/*!
* \copydoc ImmisciblePrimaryVariables::ImmisciblePrimaryVariables(Scalar)
*/
FlashPrimaryVariables(Scalar value) : ParentType(value)
{
Opm::Valgrind::CheckDefined(value);
Opm::Valgrind::SetDefined(*this);
}
/*!
* \copydoc ImmisciblePrimaryVariables::ImmisciblePrimaryVariables(const
* ImmisciblePrimaryVariables& )
*/
FlashPrimaryVariables(const FlashPrimaryVariables& value) = default;
FlashPrimaryVariables& operator=(const FlashPrimaryVariables& value) = default;
/*!
* \copydoc ImmisciblePrimaryVariables::assignMassConservative
*/
template <class FluidState>
void assignMassConservative(const FluidState& fluidState,
const MaterialLawParams& matParams OPM_UNUSED,
bool isInEquilibrium OPM_UNUSED= false)
{
// there is no difference between naive and mass conservative
// assignment in the flash model. (we only need the total
// concentrations.)
assignNaive(fluidState);
}
/*!
* \copydoc ImmisciblePrimaryVariables::assignNaive
*/
template <class FluidState>
void assignNaive(const FluidState& fluidState)
{
// reset everything
(*this) = 0.0;
// assign the phase temperatures. this is out-sourced to
// the energy module
EnergyModule::setPriVarTemperatures(*this, fluidState);
// determine the phase presence.
for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
this->operator[](cTot0Idx + compIdx) +=
fluidState.molarity(phaseIdx, compIdx) * fluidState.saturation(phaseIdx);
}
}
}
/*!
* \brief Prints the names of the primary variables and their values.
*
* \param os The \c std::ostream which should be used for the output.
*/
void print(std::ostream& os = std::cout) const
{
for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
os << "(c_tot," << FluidSystem::componentName(compIdx) << " = "
<< this->operator[](cTot0Idx + compIdx);
}
os << ")" << std::flush;
}
};
} // namespace Opm
#endif

View File

@ -0,0 +1,59 @@
// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
// vi: set et ts=4 sw=4 sts=4:
/*
This file is part of the Open Porous Media project (OPM).
OPM 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.
OPM 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 OPM. If not, see <http://www.gnu.org/licenses/>.
Consult the COPYING file in the top-level source directory of this
module for the precise wording of the license and the list of
copyright holders.
*/
/*!
* \file
* \ingroup FlashModel
*
* \brief Declares the properties required by the compositional
* multi-phase model based on flash calculations.
*/
#ifndef EWOMS_FLASH_PROPERTIES_HH
#define EWOMS_FLASH_PROPERTIES_HH
#include <opm/models/common/multiphasebaseproperties.hh>
#include <ewoms/io/vtkcompositionmodule.hh>
#include <ewoms/io/vtkenergymodule.hh>
#include <ewoms/io/vtkdiffusionmodule.hh>
BEGIN_PROPERTIES
//! Provides the thermodynamic relations
NEW_PROP_TAG(FluidSystem);
//! The type of the flash constraint solver
NEW_PROP_TAG(FlashSolver);
//! The maximum accepted error of the flash solver
NEW_PROP_TAG(FlashTolerance);
//! The thermal conduction law which ought to be used
NEW_PROP_TAG(ThermalConductionLaw);
//! The parameters of the thermal conduction law
NEW_PROP_TAG(ThermalConductionLawParams);
//! Specifies whether energy should be considered as a conservation quantity or not
NEW_PROP_TAG(EnableEnergy);
//! Enable diffusive fluxes?
NEW_PROP_TAG(EnableDiffusion);
END_PROPERTIES
#endif

View File

@ -0,0 +1,144 @@
// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
// vi: set et ts=4 sw=4 sts=4:
/*
This file is part of the Open Porous Media project (OPM).
OPM 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.
OPM 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 OPM. If not, see <http://www.gnu.org/licenses/>.
Consult the COPYING file in the top-level source directory of this
module for the precise wording of the license and the list of
copyright holders.
*/
/*!
* \file
*
* \copydoc Opm::FlashRateVector
*/
#ifndef EWOMS_FLASH_RATE_VECTOR_HH
#define EWOMS_FLASH_RATE_VECTOR_HH
#include <dune/common/fvector.hh>
#include <opm/models/common/energymodule.hh>
#include <opm/material/constraintsolvers/NcpFlash.hpp>
#include <opm/material/common/Valgrind.hpp>
#include "flashintensivequantities.hh"
namespace Opm {
/*!
* \ingroup FlashModel
*
* \copydoc Opm::ImmiscibleRateVector
*/
template <class TypeTag>
class FlashRateVector
: public Dune::FieldVector<typename GET_PROP_TYPE(TypeTag, Evaluation),
GET_PROP_VALUE(TypeTag, NumEq)>
{
typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
enum { conti0EqIdx = Indices::conti0EqIdx };
enum { numComponents = GET_PROP_VALUE(TypeTag, NumComponents) };
enum { numEq = GET_PROP_VALUE(TypeTag, NumEq) };
typedef Dune::FieldVector<Evaluation, numEq> ParentType;
typedef Opm::EnergyModule<TypeTag, GET_PROP_VALUE(TypeTag, EnableEnergy)> EnergyModule;
public:
FlashRateVector() : ParentType()
{ Opm::Valgrind::SetUndefined(*this); }
/*!
* \copydoc ImmiscibleRateVector::ImmiscibleRateVector(Scalar)
*/
FlashRateVector(const Evaluation& value) : ParentType(value)
{}
/*!
* \copydoc ImmiscibleRateVector::ImmiscibleRateVector(const
* ImmiscibleRateVector& )
*/
FlashRateVector(const FlashRateVector& value) : ParentType(value)
{}
/*!
* \copydoc ImmiscibleRateVector::setMassRate
*/
void setMassRate(const ParentType& value)
{
// convert to molar rates
ParentType molarRate(value);
for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx)
molarRate[conti0EqIdx + compIdx] /= FluidSystem::molarMass(compIdx);
setMolarRate(molarRate);
}
/*!
* \copydoc ImmiscibleRateVector::setMolarRate
*/
void setMolarRate(const ParentType& value)
{ ParentType::operator=(value); }
/*!
* \copydoc ImmiscibleRateVector::setEnthalpyRate
*/
void setEnthalpyRate(const Evaluation& rate)
{ EnergyModule::setEnthalpyRate(*this, rate); }
/*!
* \copydoc ImmiscibleRateVector::setVolumetricRate
*/
template <class FluidState, class RhsEval>
void setVolumetricRate(const FluidState& fluidState, unsigned phaseIdx, const RhsEval& volume)
{
for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx)
(*this)[conti0EqIdx + compIdx] =
fluidState.density(phaseIdx, compIdx)
* fluidState.moleFraction(phaseIdx, compIdx)
* volume;
EnergyModule::setEnthalpyRate(*this, fluidState, phaseIdx, volume);
}
/*!
* \brief Assignment operator from a scalar or a function evaluation
*/
template <class RhsEval>
FlashRateVector& operator=(const RhsEval& value)
{
for (unsigned i=0; i < this->size(); ++i)
(*this)[i] = value;
return *this;
}
/*!
* \brief Assignment operator from another rate vector
*/
FlashRateVector& operator=(const FlashRateVector& other)
{
for (unsigned i=0; i < this->size(); ++i)
(*this)[i] = other[i];
return *this;
}
};
} // namespace Opm
#endif