version 3.10-dev
co2/primaryvariableswitch.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// SPDX-FileCopyrightInfo: Copyright © DuMux Project contributors, see AUTHORS.md in root folder
5// SPDX-License-Identifier: GPL-3.0-or-later
6//
13#ifndef DUMUX_2P2C_CO2_PRIMARY_VARIABLE_SWITCH_HH
14#define DUMUX_2P2C_CO2_PRIMARY_VARIABLE_SWITCH_HH
15
16#include <iostream>
17
20
21namespace Dumux
22{
32: public PrimaryVariableSwitch< TwoPTwoCCO2PrimaryVariableSwitch >
33{
35 friend ParentType;
36
37public:
38 using ParentType::ParentType;
39
40protected:
41 // perform variable switch at a degree of freedom location
42 template<class VolumeVariables, class IndexType, class GlobalPosition>
43 bool update_(typename VolumeVariables::PrimaryVariables& priVars,
44 const VolumeVariables& volVars,
45 IndexType dofIdxGlobal,
46 const GlobalPosition& globalPos)
47 {
48 using Scalar = typename VolumeVariables::PrimaryVariables::value_type;
49
50 using FluidSystem = typename VolumeVariables::FluidSystem;
51 static constexpr int phase0Idx = FluidSystem::phase0Idx;
52 static constexpr int phase1Idx = FluidSystem::phase1Idx;
53 static constexpr int comp0Idx = FluidSystem::comp0Idx;
54 static constexpr int comp1Idx = FluidSystem::comp1Idx;
55
56 static constexpr bool useMoles = VolumeVariables::useMoles();
57 static constexpr auto formulation = VolumeVariables::priVarFormulation();
58 static_assert( (formulation == TwoPFormulation::p0s1 || formulation == TwoPFormulation::p1s0),
59 "Chosen TwoPFormulation not supported!");
60
61 using Indices = typename VolumeVariables::Indices;
62 static constexpr int switchIdx = Indices::switchIdx;
63
64 // evaluate primary variable switch
65 bool wouldSwitch = false;
66 int phasePresence = priVars.state();
67 int newPhasePresence = phasePresence;
68
69 // the param cache to evaluate the equilibrium mole fraction
70 typename FluidSystem::ParameterCache paramCache;
71
72 // check if a primary var switch is necessary
73 if (phasePresence == Indices::secondPhaseOnly)
74 {
75 // calculate wetting component mole fraction in the second phase
76 Scalar xnw = volVars.moleFraction(phase1Idx, comp0Idx);
77 Scalar xnwMax = FluidSystem::equilibriumMoleFraction(volVars.fluidState(), paramCache, phase1Idx);
78
79 // if it is larger than the equilibirum mole fraction switch
80 if(xnw > xnwMax)
81 wouldSwitch = true;
82
83 if (this->wasSwitched_[dofIdxGlobal])
84 xnwMax *= 1.02;
85
86 // if it is larger than the equilibirum mole fraction switch: first phase appears
87 if (xnw > xnwMax)
88 {
89 // wetting phase appears
90 if (this->verbosity() > 1)
91 std::cout << "First phase (" << FluidSystem::phaseName(phase0Idx) << ") appears at dof " << dofIdxGlobal
92 << ", coordinates: " << globalPos
93 << ", x^" << FluidSystem::componentName(comp0Idx) << "_" << FluidSystem::phaseName(phase1Idx) << " > x_equilibrium: "
94 << xnw << " > " << xnwMax << std::endl;
95 newPhasePresence = Indices::bothPhases;
96 if (formulation == TwoPFormulation::p1s0)
97 priVars[switchIdx] = 0.0;
98 else
99 priVars[switchIdx] = 1.0;
100 }
101 }
102 else if (phasePresence == Indices::firstPhaseOnly)
103 {
104 // calculate second component mole fraction in the wetting phase
105 Scalar xwn = volVars.moleFraction(phase0Idx, comp1Idx);
106 Scalar xwnMax = FluidSystem::equilibriumMoleFraction(volVars.fluidState(), paramCache, phase0Idx);
107
108 // if it is larger than the equilibirum mole fraction switch
109 if(xwn > xwnMax)
110 wouldSwitch = true;
111
112 if (this->wasSwitched_[dofIdxGlobal])
113 xwnMax *= 1.02;
114
115 // if it is larger than the equilibirum mole fraction switch second phase appears
116 if(xwn > xwnMax)
117 {
118 // Second phase appears
119 if (this->verbosity() > 1)
120 std::cout << "Second phase (" << FluidSystem::phaseName(phase1Idx) << ") appears at dof " << dofIdxGlobal
121 << ", coordinates: " << globalPos
122 << ", x^" << FluidSystem::componentName(comp1Idx) << "_" << FluidSystem::phaseName(phase0Idx) << " > x_equilibrium: "
123 << xwn << " > " << xwnMax << std::endl;
124 newPhasePresence = Indices::bothPhases;
125 if (formulation == TwoPFormulation::p1s0)
126 priVars[switchIdx] = 0.999;
127 else
128 priVars[switchIdx] = 0.001;
129 }
130 }
131 // TODO: this is the same as for the 2p2c model maybe factor out
132 else if (phasePresence == Indices::bothPhases)
133 {
134 Scalar Smin = 0.0;
135 if (this->wasSwitched_[dofIdxGlobal])
136 Smin = -0.01;
137
138 if (volVars.saturation(phase1Idx) <= Smin)
139 {
140 wouldSwitch = true;
141 // nonwetting phase disappears
142 if (this->verbosity() > 1)
143 std::cout << "Second phase (" << FluidSystem::phaseName(phase1Idx) << ") disappears at dof " << dofIdxGlobal
144 << ", coordinates: " << globalPos
145 << ", S_" << FluidSystem::phaseName(phase1Idx) << ": " << volVars.saturation(phase1Idx)
146 << std::endl;
147 newPhasePresence = Indices::firstPhaseOnly;
148
149 if(useMoles) // mole-fraction formulation
150 priVars[switchIdx] = volVars.moleFraction(phase0Idx, comp1Idx);
151 else // mass-fraction formulation
152 priVars[switchIdx] = volVars.massFraction(phase0Idx, comp1Idx);
153 }
154 else if (volVars.saturation(phase0Idx) <= Smin)
155 {
156 wouldSwitch = true;
157 // wetting phase disappears
158 if (this->verbosity() > 1)
159 std::cout << "First phase (" << FluidSystem::phaseName(phase0Idx) << ") disappears at dof " << dofIdxGlobal
160 << ", coordinates: " << globalPos
161 << ", S_" << FluidSystem::phaseName(phase0Idx) << ": " << volVars.saturation(phase0Idx)
162 << std::endl;
163 newPhasePresence = Indices::secondPhaseOnly;
164
165 if(useMoles) // mole-fraction formulation
166 priVars[switchIdx] = volVars.moleFraction(phase1Idx, comp0Idx);
167 else // mass-fraction formulation
168 priVars[switchIdx] = volVars.massFraction(phase1Idx, comp0Idx);
169 }
170 }
171
172 priVars.setState(newPhasePresence);
173 this->wasSwitched_[dofIdxGlobal] = wouldSwitch;
174 return phasePresence != newPhasePresence;
175 }
176};
177
178} // end namespace Dumux
179
180#endif
The primary variable switch controlling the phase presence state variable.
Definition: compositional/primaryvariableswitch.hh:49
std::vector< bool > wasSwitched_
Definition: compositional/primaryvariableswitch.hh:433
int verbosity() const
The verbosity level.
Definition: compositional/primaryvariableswitch.hh:271
The primary variable switch for the 2p2c-CO2 model controlling the phase presence state variable.
Definition: co2/primaryvariableswitch.hh:33
bool update_(typename VolumeVariables::PrimaryVariables &priVars, const VolumeVariables &volVars, IndexType dofIdxGlobal, const GlobalPosition &globalPos)
Definition: co2/primaryvariableswitch.hh:43
The primary variable switch base class for compositional models.
Defines an enumeration for the formulations accepted by the two-phase model.
@ p1s0
first phase saturation and second phase pressure as primary variables
@ p0s1
first phase pressure and second phase saturation as primary variables
std::string phasePresence() noexcept
I/O name of phase presence.
Definition: name.hh:135
Definition: adapt.hh:17