My Project
Loading...
Searching...
No Matches
flashprimaryvariables.hh
Go to the documentation of this file.
1// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2// vi: set et ts=4 sw=4 sts=4:
3/*
4 This file is part of the Open Porous Media project (OPM).
5
6 OPM is free software: you can redistribute it and/or modify
7 it under the terms of the GNU General Public License as published by
8 the Free Software Foundation, either version 2 of the License, or
9 (at your option) any later version.
10
11 OPM is distributed in the hope that it will be useful,
12 but WITHOUT ANY WARRANTY; without even the implied warranty of
13 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 GNU General Public License for more details.
15
16 You should have received a copy of the GNU General Public License
17 along with OPM. If not, see <http://www.gnu.org/licenses/>.
18
19 Consult the COPYING file in the top-level source directory of this
20 module for the precise wording of the license and the list of
21 copyright holders.
22*/
28#ifndef EWOMS_FLASH_PRIMARY_VARIABLES_HH
29#define EWOMS_FLASH_PRIMARY_VARIABLES_HH
30
33
34#include <opm/material/common/Valgrind.hpp>
35
36#include <iostream>
37
38namespace Opm {
39
49template <class TypeTag>
51{
53
59
60 // primary variable indices
61 enum { cTot0Idx = Indices::cTot0Idx };
62
64 enum { numComponents = getPropValue<TypeTag, Properties::NumComponents>() };
65
66 using Toolbox = typename Opm::MathToolbox<Evaluation>;
68
69public:
71 { Opm::Valgrind::SetDefined(*this); }
72
82 FlashPrimaryVariables& operator=(const FlashPrimaryVariables& value) = default;
83
84 using ParentType::operator=;
85
89 template <class FluidState>
90 void assignMassConservative(const FluidState& fluidState,
91 const MaterialLawParams&,
92 bool = false)
93 {
94 // there is no difference between naive and mass conservative
95 // assignment in the flash model. (we only need the total
96 // concentrations.)
97 assignNaive(fluidState);
98 }
99
103 template <class FluidState>
104 void assignNaive(const FluidState& fluidState)
105 {
106 // reset everything
107 (*this) = 0.0;
108
109 // assign the phase temperatures. this is out-sourced to
110 // the energy module
111 EnergyModule::setPriVarTemperatures(*this, fluidState);
112
113 // determine the phase presence.
114 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
115 for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
116 this->operator[](cTot0Idx + compIdx) +=
117 fluidState.molarity(phaseIdx, compIdx) * fluidState.saturation(phaseIdx);
118 }
119 }
120 }
121
127 void print(std::ostream& os = std::cout) const
128 {
129 for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
130 os << "(c_tot," << FluidSystem::componentName(compIdx) << " = "
131 << this->operator[](cTot0Idx + compIdx);
132 }
133 os << ")" << std::flush;
134 }
135};
136
137} // namespace Opm
138
139#endif
Provides the auxiliary methods required for consideration of the energy equation.
Definition energymodule.hh:50
Represents the primary variables used by the compositional flow model based on flash calculations.
Definition flashprimaryvariables.hh:51
void print(std::ostream &os=std::cout) const
Prints the names of the primary variables and their values.
Definition flashprimaryvariables.hh:127
void assignNaive(const FluidState &fluidState)
Directly retrieve the primary variables from an arbitrary fluid state.
Definition flashprimaryvariables.hh:104
FlashPrimaryVariables(const FlashPrimaryVariables &value)=default
void assignMassConservative(const FluidState &fluidState, const MaterialLawParams &, bool=false)
< Import base class assignment operators.
Definition flashprimaryvariables.hh:90
Represents the primary variables used by the a model.
Definition fvbaseprimaryvariables.hh:52
Contains the classes required to consider energy as a conservation quantity in a multi-phase module.
Represents the primary variables used by the a model.
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition blackoilboundaryratevector.hh:37
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:242
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