28#ifndef EWOMS_PVS_MODEL_HH
29#define EWOMS_PVS_MODEL_HH
31#include <opm/common/Exceptions.hpp>
33#include <opm/material/densead/Math.hpp>
34#include <opm/material/fluidmatrixinteractions/MaterialTraits.hpp>
35#include <opm/material/fluidmatrixinteractions/NullMaterial.hpp>
61template <
class TypeTag>
65namespace Opm::Properties {
72 using InheritsFrom = std::tuple<MultiPhaseBaseModel>; };
76template<
class TypeTag>
81template<
class TypeTag>
86template<
class TypeTag>
91template<
class TypeTag>
96template<
class TypeTag>
101template<
class TypeTag>
106template<
class TypeTag>
111template<
class TypeTag>
116template<
class TypeTag>
121template<
class TypeTag>
123{
static constexpr bool value =
false; };
126template<
class TypeTag>
128{
static constexpr bool value =
false; };
131template<
class TypeTag>
135 static constexpr type value = 1.0;
139template<
class TypeTag>
143 static constexpr type value = 1.0;
147template<
class TypeTag>
151 static constexpr type value = 1.0;
156namespace Opm::Parameters {
261template <
class TypeTag>
282 using Element =
typename GridView::template Codim<0>::Entity;
283 using ElementIterator =
typename GridView::template Codim<0>::Iterator;
291 verbosity_ = Parameters::Get<Parameters::PvsVerbosity>();
312 Parameters::Register<Parameters::PvsVerbosity>
313 (
"The verbosity level of the primary variable "
329 if (!(s = EnergyModule::primaryVarName(pvIdx)).empty())
332 std::ostringstream oss;
333 if (pvIdx == Indices::pressure0Idx)
334 oss <<
"pressure_" << FluidSystem::phaseName(0);
335 else if (Indices::switch0Idx <= pvIdx
336 && pvIdx < Indices::switch0Idx + numPhases - 1)
337 oss <<
"switch_" << pvIdx - Indices::switch0Idx;
338 else if (Indices::switch0Idx + numPhases - 1 <= pvIdx
339 && pvIdx < Indices::switch0Idx + numComponents - 1)
340 oss <<
"auxMoleFrac^" << FluidSystem::componentName(pvIdx);
353 if (!(s = EnergyModule::eqName(eqIdx)).empty())
356 std::ostringstream oss;
357 if (Indices::conti0EqIdx <= eqIdx && eqIdx < Indices::conti0EqIdx
359 unsigned compIdx = eqIdx - Indices::conti0EqIdx;
360 oss <<
"continuity^" << FluidSystem::componentName(compIdx);
373 ParentType::updateFailed();
382 ParentType::updateBegin();
387 size_t nDof = this->numTotalDof();
388 for (
unsigned dofIdx = 0; dofIdx < nDof; ++ dofIdx) {
389 if (this->dofTotalVolume(dofIdx) > 0.0) {
391 this->solution(0)[dofIdx][Indices::pressure0Idx];
392 if (referencePressure_ > 0.0)
403 Scalar tmp = EnergyModule::primaryVarWeight(*
this, globalDofIdx, pvIdx);
408 if (Indices::pressure0Idx == pvIdx) {
409 return 10 / referencePressure_;
412 if (Indices::switch0Idx <= pvIdx && pvIdx < Indices::switch0Idx
414 unsigned phaseIdx = pvIdx - Indices::switch0Idx;
416 if (!this->solution(0)[globalDofIdx].phaseIsPresent(phaseIdx))
433 Scalar
eqWeight(
unsigned globalDofIdx,
unsigned eqIdx)
const
435 Scalar tmp = EnergyModule::eqWeight(*
this, globalDofIdx, eqIdx);
440 unsigned compIdx = eqIdx - Indices::conti0EqIdx;
441 assert(compIdx <= numComponents);
444 return FluidSystem::molarMass(compIdx);
452 ParentType::advanceTimeLevel();
461 {
return numSwitched_ > 0; }
466 template <
class DofEntity>
470 ParentType::serializeEntity(outstream, dofEntity);
472 unsigned dofIdx =
static_cast<unsigned>(this->dofMapper().index(dofEntity));
473 if (!outstream.good())
474 throw std::runtime_error(
"Could not serialize DOF "+std::to_string(dofIdx));
476 outstream << this->solution(0)[dofIdx].phasePresence() <<
" ";
482 template <
class DofEntity>
486 ParentType::deserializeEntity(instream, dofEntity);
489 unsigned dofIdx =
static_cast<unsigned>(this->dofMapper().index(dofEntity));
490 if (!instream.good())
491 throw std::runtime_error(
"Could not deserialize DOF "+std::to_string(dofIdx));
495 this->solution(0)[dofIdx].setPhasePresence(tmp);
496 this->solution(1)[dofIdx].setPhasePresence(tmp);
506 void switchPrimaryVars_()
512 std::vector<bool> visited(this->numGridDof(),
false);
513 ElementContext elemCtx(this->simulator_);
515 for (
const auto& elem : elements(this->gridView_, Dune::Partitions::interior)) {
516 elemCtx.updateStencil(elem);
518 size_t numLocalDof = elemCtx.stencil(0).numPrimaryDof();
519 for (
unsigned dofIdx = 0; dofIdx < numLocalDof; ++dofIdx) {
520 unsigned globalIdx = elemCtx.globalSpaceIndex(dofIdx, 0);
522 if (visited[globalIdx])
524 visited[globalIdx] =
true;
527 auto& priVars = this->solution(0)[globalIdx];
528 elemCtx.updateIntensiveQuantities(priVars, dofIdx, 0);
529 const IntensiveQuantities& intQuants = elemCtx.intensiveQuantities(dofIdx, 0);
532 short oldPhasePresence = priVars.phasePresence();
536 priVars.assignNaive(intQuants.fluidState());
538 if (oldPhasePresence != priVars.phasePresence()) {
540 printSwitchedPhases_(elemCtx,
542 intQuants.fluidState(),
554 std::cout <<
"rank " << this->simulator_.gridView().comm().rank()
555 <<
" caught an exception during primary variable switching"
556 <<
"\n" << std::flush;
559 succeeded = this->simulator_.gridView().comm().min(succeeded);
562 throw NumericalProblem(
"A process did not succeed in adapting the primary variables");
567 numSwitched_ = this->gridView_.comm().sum(numSwitched_);
570 this->simulator_.model().newtonMethod().endIterMsg()
571 <<
", num switched=" << numSwitched_;
574 template <
class Flu
idState>
575 void printSwitchedPhases_(
const ElementContext& elemCtx,
577 const FluidState& fs,
578 short oldPhasePresence,
579 const PrimaryVariables& newPv)
const
581 using FsToolbox = Opm::MathToolbox<typename FluidState::Scalar>;
583 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
584 bool oldPhasePresent = (oldPhasePresence& (1 << phaseIdx)) > 0;
585 bool newPhasePresent = newPv.phaseIsPresent(phaseIdx);
586 if (oldPhasePresent == newPhasePresent)
589 const auto& pos = elemCtx.pos(dofIdx, 0);
590 if (oldPhasePresent && !newPhasePresent) {
591 std::cout <<
"'" << FluidSystem::phaseName(phaseIdx)
592 <<
"' phase disappears at position " << pos
593 <<
". saturation=" << fs.saturation(phaseIdx)
598 for (
unsigned compIdx = 0; compIdx < numComponents; ++compIdx)
599 sumx += FsToolbox::value(fs.moleFraction(phaseIdx, compIdx));
601 std::cout <<
"'" << FluidSystem::phaseName(phaseIdx)
602 <<
"' phase appears at position " << pos
603 <<
" sum x = " << sumx << std::flush;
607 std::cout <<
", new primary variables: ";
609 std::cout <<
"\n" << std::flush;
612 void registerOutputModules_()
614 ParentType::registerOutputModules_();
625 mutable Scalar referencePressure_;
629 unsigned numSwitched_;
Provides the auxiliary methods required for consideration of the energy equation.
Definition energymodule.hh:50
A base class for fully-implicit multi-phase porous-media flow models which assume multiple fluid phas...
Definition multiphasebasemodel.hh:153
static void registerParameters()
Register all run-time parameters for the immiscible model.
Definition multiphasebasemodel.hh:179
Implements a rate vector on the boundary for the fully implicit compositional multi-phase primary var...
Definition pvsboundaryratevector.hh:47
Contains all data which is required to calculate all fluxes at a flux integration point for the prima...
Definition pvsextensivequantities.hh:54
The indices for the compositional multi-phase primary variable switching model.
Definition pvsindices.hh:48
Contains the quantities which are are constant within a finite volume in the compositional multi-phas...
Definition pvsintensivequantities.hh:61
Element-wise calculation of the local residual for the compositional multi-phase primary variable swi...
Definition pvslocalresidual.hh:48
A generic compositional multi-phase model using primary-variable switching.
Definition pvsmodel.hh:264
void advanceTimeLevel()
Called by the problem if a time integration was successful, post processing of the solution is done a...
Definition pvsmodel.hh:450
void serializeEntity(std::ostream &outstream, const DofEntity &dofEntity)
Write the current solution for a degree of freedom to a restart file.
Definition pvsmodel.hh:467
std::string eqName(unsigned eqIdx) const
Given an equation index, return a human readable name.
Definition pvsmodel.hh:350
std::string primaryVarName(unsigned pvIdx) const
Given an primary variable index, return a human readable name.
Definition pvsmodel.hh:326
void updateBegin()
Called by the update() method before it tries to apply the newton method.
Definition pvsmodel.hh:380
Scalar eqWeight(unsigned globalDofIdx, unsigned eqIdx) const
Returns the relative weight of an equation.
Definition pvsmodel.hh:433
static std::string name()
Definition pvsmodel.hh:320
bool switched() const
Return true if the primary variables were switched for at least one vertex after the last timestep.
Definition pvsmodel.hh:460
static void registerParameters()
Register all run-time parameters for the PVS compositional model.
Definition pvsmodel.hh:298
Scalar primaryVarWeight(unsigned globalDofIdx, unsigned pvIdx) const
Returns the relative weight of a primary variable for calculating relative errors.
Definition pvsmodel.hh:401
void deserializeEntity(std::istream &instream, const DofEntity &dofEntity)
Reads the current solution variables for a degree of freedom from a restart file.
Definition pvsmodel.hh:483
void updateFailed()
Called by the update() method if it was unsuccessful.
Definition pvsmodel.hh:371
A newton solver which is specific to the compositional multi-phase PVS model.
Definition pvsnewtonmethod.hh:52
Represents the primary variables used in the primary variable switching compositional model.
Definition pvsprimaryvariables.hh:60
Implements a vector representing molar, mass or volumetric rates.
Definition pvsratevector.hh:53
VTK output module for the fluid composition.
Definition vtkcompositionmodule.hpp:57
static void registerParameters()
Register all run-time parameters for the Vtk output module.
Definition vtkcompositionmodule.hpp:86
VTK output module for quantities which make sense for models which incorperate molecular diffusion.
Definition vtkdiffusionmodule.hpp:58
static void registerParameters()
Register all run-time parameters for the Vtk output module.
Definition vtkdiffusionmodule.hpp:87
VTK output module for quantities which make sense for models which assume thermal equilibrium.
Definition vtkenergymodule.hpp:58
static void registerParameters()
Register all run-time parameters for the Vtk output module.
Definition vtkenergymodule.hpp:86
VTK output module for the fluid composition.
Definition vtkphasepresencemodule.hpp:48
static void registerParameters()
Register all run-time parameters for the Vtk output module.
Definition vtkphasepresencemodule.hpp:72
Classes required for molecular diffusion.
Contains the classes required to consider energy as a conservation quantity in a multi-phase module.
A base class for fully-implicit multi-phase porous-media flow models which assume multiple fluid phas...
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition blackoilboundaryratevector.hh:37
typename Properties::Detail::GetPropImpl< TypeTag, Property >::type::type GetPropType
get the type alias defined in the property (equivalent to old macro GET_PROP_TYPE(....
Definition propertysystem.hh:235
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:242
Implements a rate vector on the boundary for the fully implicit compositional multi-phase primary var...
Contains all data which is required to calculate all fluxes at a flux integration point for the prima...
The indices for the compositional multi-phase primary variable switching model.
Contains the quantities which are are constant within a finite volume in the compositional multi-phas...
Element-wise calculation of the local residual for the compositional multi-phase primary variable swi...
A newton solver which is specific to the compositional multi-phase PVS model.
Represents the primary variables used in the primary variable switching compositional model.
Declares the properties required for the compositional multi-phase primary variable switching model.
Implements a vector representing molar, mass or volumetric rates.
The verbosity of the model (0 -> do not print anything, 2 -> spam stdout a lot)
Definition pvsmodel.hh:159
Type of object for specifying boundary conditions.
Definition fvbaseproperties.hh:119
Enable diffusive fluxes?
Definition multiphasebaseproperties.hh:79
Specify whether energy should be considered as a conservation quantity or not.
Definition multiphasebaseproperties.hh:76
Data required to calculate a flux over a face.
Definition fvbaseproperties.hh:149
Enumerations used by the model.
Definition multiphasebaseproperties.hh:48
The secondary variables within a sub-control volume.
Definition fvbaseproperties.hh:133
The type of the local residual function.
Definition fvbaseproperties.hh:94
The type of the model.
Definition basicproperties.hh:88
Specifies the type of the actual Newton method.
Definition nullconvergencewriter.hh:39
A vector of primary variables within a sub-control volume.
Definition fvbaseproperties.hh:130
The basis value for the weight of the mole fraction primary variables.
Definition pvsproperties.hh:51
The basis value for the weight of the pressure primary variable.
Definition pvsproperties.hh:45
The basis value for the weight of the saturation primary variables.
Definition pvsproperties.hh:48
Vector containing volumetric or areal rates of quantities.
Definition fvbaseproperties.hh:116
The type tag for the isothermal single phase problems.
Definition pvsmodel.hh:71
VTK output module for the fluid composition.
VTK output module for quantities which make sense for models which incorperate molecular diffusion.
VTK output module for quantities which make sense for models which assume thermal equilibrium.