| Line | Branch | Exec | Source |
|---|---|---|---|
| 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 | // | ||
| 7 | /*! | ||
| 8 | * \file | ||
| 9 | * \ingroup CO2Tests | ||
| 10 | * \brief Definition of a problem, where CO2 is injected into a reservoir. | ||
| 11 | */ | ||
| 12 | |||
| 13 | #ifndef DUMUX_HETEROGENEOUS_PROBLEM_HH | ||
| 14 | #define DUMUX_HETEROGENEOUS_PROBLEM_HH | ||
| 15 | |||
| 16 | #include <dune/grid/common/partitionset.hh> | ||
| 17 | |||
| 18 | #include <dumux/common/properties.hh> | ||
| 19 | #include <dumux/common/parameters.hh> | ||
| 20 | #include <dumux/common/boundarytypes.hh> | ||
| 21 | #include <dumux/common/numeqvector.hh> | ||
| 22 | |||
| 23 | #include <dumux/parallel/vectorcommdatahandle.hh> | ||
| 24 | |||
| 25 | #include <dumux/porousmediumflow/problem.hh> | ||
| 26 | #include <dumux/discretization/box/scvftoscvboundarytypes.hh> | ||
| 27 | #include <dumux/common/gridcapabilities.hh> | ||
| 28 | |||
| 29 | |||
| 30 | namespace Dumux { | ||
| 31 | |||
| 32 | /*! | ||
| 33 | * \ingroup CO2Tests | ||
| 34 | * \brief Definition of a problem, where CO2 is injected into a reservoir. | ||
| 35 | * | ||
| 36 | * The domain is sized 200m times 100m and consists of four layers, a | ||
| 37 | * permeable reservoir layer at the bottom, a barrier rock layer with reduced | ||
| 38 | * permeability, another reservoir layer and at the top a barrier rock layer | ||
| 39 | * with a very low permeablility. | ||
| 40 | * | ||
| 41 | * CO2 is injected at the permeable bottom layer from the left side. | ||
| 42 | * The domain is initially filled with brine. | ||
| 43 | * | ||
| 44 | * The grid is unstructered and permeability and porosity for the elements are | ||
| 45 | * read in from the grid file. The grid file also contains so-called boundary | ||
| 46 | * IDs which can be used during the grid creation in order to differentiate | ||
| 47 | * between different parts of the boundary. | ||
| 48 | * These boundary IDs can be imported into the problem where the boundary | ||
| 49 | * conditions can then be assigned accordingly. | ||
| 50 | * | ||
| 51 | * The model is able to use either mole or mass fractions. The property useMoles | ||
| 52 | * can be set to either true or false in the problem file. Make sure that the | ||
| 53 | * according units are used in the problem setup. | ||
| 54 | * The default setting for useMoles is false. | ||
| 55 | */ | ||
| 56 | template <class TypeTag > | ||
| 57 | class HeterogeneousProblem : public PorousMediumFlowProblem<TypeTag> | ||
| 58 | { | ||
| 59 | using ParentType = PorousMediumFlowProblem<TypeTag>; | ||
| 60 | using GridView = typename GetPropType<TypeTag, Properties::GridGeometry>::GridView; | ||
| 61 | using Scalar = GetPropType<TypeTag, Properties::Scalar>; | ||
| 62 | using FluidSystem = GetPropType<TypeTag, Properties::FluidSystem>; | ||
| 63 | using GridVariables = GetPropType<TypeTag, Properties::GridVariables>; | ||
| 64 | using ElementVolumeVariables = typename GridVariables::GridVolumeVariables::LocalView; | ||
| 65 | using ElementFluxVariablesCache = typename GridVariables::GridFluxVariablesCache::LocalView; | ||
| 66 | using VolumeVariables = typename GridVariables::GridVolumeVariables::VolumeVariables; | ||
| 67 | |||
| 68 | using ModelTraits = GetPropType<TypeTag, Properties::ModelTraits>; | ||
| 69 | using Indices = typename ModelTraits::Indices; | ||
| 70 | |||
| 71 | // copy some indices for convenience | ||
| 72 | enum | ||
| 73 | { | ||
| 74 | // primary variable indices | ||
| 75 | pressureIdx = Indices::pressureIdx, | ||
| 76 | switchIdx = Indices::switchIdx, | ||
| 77 | |||
| 78 | // phase presence index | ||
| 79 | firstPhaseOnly = Indices::firstPhaseOnly, | ||
| 80 | |||
| 81 | // component indices | ||
| 82 | BrineIdx = FluidSystem::BrineIdx, | ||
| 83 | CO2Idx = FluidSystem::CO2Idx, | ||
| 84 | |||
| 85 | // equation indices | ||
| 86 | conti0EqIdx = Indices::conti0EqIdx, | ||
| 87 | contiCO2EqIdx = conti0EqIdx + CO2Idx | ||
| 88 | }; | ||
| 89 | |||
| 90 | #if !ISOTHERMAL | ||
| 91 | enum { | ||
| 92 | temperatureIdx = Indices::temperatureIdx, | ||
| 93 | energyEqIdx = Indices::energyEqIdx, | ||
| 94 | }; | ||
| 95 | #endif | ||
| 96 | |||
| 97 | using PrimaryVariables = GetPropType<TypeTag, Properties::PrimaryVariables>; | ||
| 98 | using NumEqVector = Dumux::NumEqVector<PrimaryVariables>; | ||
| 99 | using BoundaryTypes = Dumux::BoundaryTypes<GetPropType<TypeTag, Properties::ModelTraits>::numEq()>; | ||
| 100 | using Element = typename GridView::template Codim<0>::Entity; | ||
| 101 | using GlobalPosition = typename Element::Geometry::GlobalCoordinate; | ||
| 102 | using GridGeometry = GetPropType<TypeTag, Properties::GridGeometry>; | ||
| 103 | using FVElementGeometry = typename GetPropType<TypeTag, Properties::GridGeometry>::LocalView; | ||
| 104 | using SubControlVolumeFace = typename FVElementGeometry::SubControlVolumeFace; | ||
| 105 | using SubControlVolume = typename FVElementGeometry::SubControlVolume; | ||
| 106 | |||
| 107 | using CO2 = typename FluidSystem::CO2; | ||
| 108 | |||
| 109 | //! Property that defines whether mole or mass fractions are used | ||
| 110 | static constexpr bool useMoles = ModelTraits::useMoles(); | ||
| 111 | |||
| 112 | // the discretization method we are using | ||
| 113 | using DiscretizationMethod = DiscretizationMethods::Box; | ||
| 114 | static constexpr DiscretizationMethod discMethod{}; | ||
| 115 | static constexpr bool isBox = GridGeometry::discMethod == DiscretizationMethods::box; | ||
| 116 | |||
| 117 | // world dimension to access gravity vector | ||
| 118 | static constexpr int dimWorld = GridView::dimensionworld; | ||
| 119 | |||
| 120 | public: | ||
| 121 | template<class SpatialParams> | ||
| 122 | 14 | HeterogeneousProblem(std::shared_ptr<const GridGeometry> gridGeometry, std::shared_ptr<SpatialParams> spatialParams) | |
| 123 | : ParentType(gridGeometry, spatialParams) | ||
| 124 | , injectionTop_(1) | ||
| 125 | , injectionBottom_(2) | ||
| 126 | , dirichletBoundary_(3) | ||
| 127 |
13/44✓ Branch 1 taken 14 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 14 times.
✗ Branch 5 not taken.
✓ Branch 9 taken 14 times.
✗ Branch 10 not taken.
✓ Branch 11 taken 14 times.
✗ Branch 12 not taken.
✓ Branch 13 taken 14 times.
✗ Branch 14 not taken.
✗ Branch 15 not taken.
✓ Branch 16 taken 14 times.
✓ Branch 18 taken 14 times.
✗ Branch 19 not taken.
✓ Branch 21 taken 14 times.
✗ Branch 22 not taken.
✓ Branch 24 taken 14 times.
✗ Branch 25 not taken.
✓ Branch 27 taken 14 times.
✗ Branch 28 not taken.
✓ Branch 30 taken 14 times.
✗ Branch 31 not taken.
✓ Branch 33 taken 14 times.
✗ Branch 34 not taken.
✓ Branch 36 taken 14 times.
✗ Branch 37 not taken.
✗ Branch 38 not taken.
✗ Branch 39 not taken.
✗ Branch 40 not taken.
✗ Branch 41 not taken.
✗ Branch 42 not taken.
✗ Branch 43 not taken.
✗ Branch 46 not taken.
✗ Branch 47 not taken.
✗ Branch 48 not taken.
✗ Branch 49 not taken.
✗ Branch 50 not taken.
✗ Branch 51 not taken.
✗ Branch 52 not taken.
✗ Branch 53 not taken.
✗ Branch 54 not taken.
✗ Branch 55 not taken.
✗ Branch 56 not taken.
✗ Branch 57 not taken.
|
56 | , noFlowBoundary_(4) |
| 128 | { | ||
| 129 |
1/2✓ Branch 1 taken 14 times.
✗ Branch 2 not taken.
|
14 | nTemperature_ = getParam<int>("FluidSystem.NTemperature"); |
| 130 |
1/2✓ Branch 1 taken 14 times.
✗ Branch 2 not taken.
|
14 | nPressure_ = getParam<int>("FluidSystem.NPressure"); |
| 131 |
1/2✓ Branch 1 taken 14 times.
✗ Branch 2 not taken.
|
14 | pressureLow_ = getParam<Scalar>("FluidSystem.PressureLow"); |
| 132 |
1/2✓ Branch 1 taken 14 times.
✗ Branch 2 not taken.
|
14 | pressureHigh_ = getParam<Scalar>("FluidSystem.PressureHigh"); |
| 133 |
1/2✓ Branch 1 taken 14 times.
✗ Branch 2 not taken.
|
14 | temperatureLow_ = getParam<Scalar>("FluidSystem.TemperatureLow"); |
| 134 |
1/2✓ Branch 1 taken 14 times.
✗ Branch 2 not taken.
|
14 | temperatureHigh_ = getParam<Scalar>("FluidSystem.TemperatureHigh"); |
| 135 |
1/2✓ Branch 1 taken 14 times.
✗ Branch 2 not taken.
|
14 | depthBOR_ = getParam<Scalar>("Problem.DepthBOR"); |
| 136 |
2/4✓ Branch 1 taken 14 times.
✗ Branch 2 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 14 times.
|
14 | name_ = getParam<std::string>("Problem.Name"); |
| 137 |
1/2✓ Branch 1 taken 14 times.
✗ Branch 2 not taken.
|
14 | injectionRate_ = getParam<Scalar>("Problem.InjectionRate"); |
| 138 | |||
| 139 | #if !ISOTHERMAL | ||
| 140 |
1/2✓ Branch 1 taken 6 times.
✗ Branch 2 not taken.
|
6 | injectionPressure_ = getParam<Scalar>("Problem.InjectionPressure"); |
| 141 |
1/2✓ Branch 1 taken 6 times.
✗ Branch 2 not taken.
|
6 | injectionTemperature_ = getParam<Scalar>("Problem.InjectionTemperature"); |
| 142 | #endif | ||
| 143 | |||
| 144 | // set the spatial parameters by reading the DGF grid file | ||
| 145 |
2/4✓ Branch 1 taken 14 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 14 times.
✗ Branch 5 not taken.
|
28 | this->spatialParams().getParamsFromGrid(); |
| 146 | |||
| 147 | // initialize the tables of the fluid system | ||
| 148 |
1/2✓ Branch 1 taken 14 times.
✗ Branch 2 not taken.
|
14 | FluidSystem::init(/*Tmin=*/temperatureLow_, |
| 149 | /*Tmax=*/temperatureHigh_, | ||
| 150 | /*nT=*/nTemperature_, | ||
| 151 | /*pmin=*/pressureLow_, | ||
| 152 | /*pmax=*/pressureHigh_, | ||
| 153 | /*np=*/nPressure_); | ||
| 154 | |||
| 155 | // stating in the console whether mole or mass fractions are used | ||
| 156 |
4/4✓ Branch 3 taken 10 times.
✓ Branch 4 taken 4 times.
✓ Branch 5 taken 10 times.
✓ Branch 6 taken 4 times.
|
42 | if (gridGeometry->gridView().comm().rank() == 0) |
| 157 | { | ||
| 158 | if (useMoles) | ||
| 159 | std::cout << "-- Problem uses mole fractions." << std::endl; | ||
| 160 | else | ||
| 161 |
1/2✓ Branch 2 taken 10 times.
✗ Branch 3 not taken.
|
20 | std::cout << "-- Problem uses mass fractions." << std::endl; |
| 162 | } | ||
| 163 | |||
| 164 | // precompute the boundary types for the box method from the cell-centered boundary types | ||
| 165 |
1/2✓ Branch 1 taken 14 times.
✗ Branch 2 not taken.
|
14 | scvfToScvBoundaryTypes_.computeBoundaryTypes(*this); |
| 166 | 14 | } | |
| 167 | |||
| 168 | /*! | ||
| 169 | * \brief Appends all quantities of interest which can be derived | ||
| 170 | * from the solution of the current time step to the VTK | ||
| 171 | * writer. | ||
| 172 | */ | ||
| 173 | template<class VTKWriter> | ||
| 174 | 14 | void addFieldsToWriter(VTKWriter& vtk) | |
| 175 | { | ||
| 176 | 42 | const auto numElements = this->gridGeometry().gridView().size(0); | |
| 177 | 28 | const auto numDofs = this->gridGeometry().numDofs(); | |
| 178 | |||
| 179 | 14 | vtkKxx_.resize(numElements); | |
| 180 | 14 | vtkPorosity_.resize(numElements); | |
| 181 | 14 | vtkBoxVolume_.resize(numDofs, 0.0); | |
| 182 | |||
| 183 |
5/12✓ Branch 1 taken 14 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 14 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 14 times.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✓ Branch 10 taken 14 times.
✓ Branch 12 taken 14 times.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✗ Branch 15 not taken.
|
28 | vtk.addField(vtkKxx_, "Kxx"); |
| 184 |
5/12✓ Branch 1 taken 14 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 14 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 14 times.
✗ Branch 8 not taken.
✓ Branch 9 taken 14 times.
✗ Branch 10 not taken.
✓ Branch 12 taken 14 times.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✗ Branch 15 not taken.
|
42 | vtk.addField(vtkPorosity_, "cellwisePorosity"); |
| 185 |
5/12✓ Branch 1 taken 14 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 14 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 14 times.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✓ Branch 10 taken 14 times.
✓ Branch 12 taken 6 times.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✗ Branch 15 not taken.
|
28 | vtk.addField(vtkBoxVolume_, "boxVolume"); |
| 186 | |||
| 187 | #if !ISOTHERMAL | ||
| 188 |
7/18✓ Branch 1 taken 6 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 6 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 6 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 6 times.
✗ Branch 11 not taken.
✓ Branch 12 taken 6 times.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✓ Branch 15 taken 6 times.
✓ Branch 17 taken 6 times.
✗ Branch 18 not taken.
✗ Branch 19 not taken.
✗ Branch 20 not taken.
✗ Branch 21 not taken.
✗ Branch 22 not taken.
|
18 | vtk.addVolumeVariable([](const VolumeVariables& v){ return v.enthalpy(BrineIdx); }, "enthalpyW"); |
| 189 |
7/18✓ Branch 1 taken 6 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 6 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 6 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 6 times.
✗ Branch 11 not taken.
✓ Branch 12 taken 6 times.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✓ Branch 15 taken 6 times.
✓ Branch 17 taken 6 times.
✗ Branch 18 not taken.
✗ Branch 19 not taken.
✗ Branch 20 not taken.
✗ Branch 21 not taken.
✗ Branch 22 not taken.
|
18 | vtk.addVolumeVariable([](const VolumeVariables& v){ return v.enthalpy(CO2Idx); }, "enthalpyN"); |
| 190 | #else | ||
| 191 | 8 | vtkTemperature_.resize(numDofs, 0.0); | |
| 192 |
5/12✓ Branch 1 taken 8 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 8 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 8 times.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✓ Branch 10 taken 8 times.
✓ Branch 12 taken 8 times.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✗ Branch 15 not taken.
|
16 | vtk.addField(vtkTemperature_, "T"); |
| 193 | #endif | ||
| 194 | |||
| 195 |
2/4✓ Branch 1 taken 14 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 14 times.
✗ Branch 5 not taken.
|
28 | const auto& gridView = this->gridGeometry().gridView(); |
| 196 |
3/8✓ Branch 1 taken 14 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 14 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 3 times.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
|
34 | auto fvGeometry = localView(this->gridGeometry()); |
| 197 |
9/14✓ Branch 1 taken 14 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 14 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 14 times.
✗ Branch 8 not taken.
✓ Branch 9 taken 36670 times.
✓ Branch 10 taken 14 times.
✓ Branch 11 taken 36670 times.
✓ Branch 12 taken 14 times.
✓ Branch 14 taken 36670 times.
✗ Branch 15 not taken.
✓ Branch 19 taken 5 times.
✗ Branch 20 not taken.
|
73378 | for (const auto& element : elements(gridView, Dune::Partitions::interior)) |
| 198 | { | ||
| 199 |
3/6✓ Branch 1 taken 36670 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 36670 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 36670 times.
✗ Branch 8 not taken.
|
110010 | const auto eIdx = this->gridGeometry().elementMapper().index(element); |
| 200 |
1/2✓ Branch 1 taken 36670 times.
✗ Branch 2 not taken.
|
36670 | fvGeometry.bindElement(element); |
| 201 | |||
| 202 |
4/4✓ Branch 0 taken 80674 times.
✓ Branch 1 taken 36670 times.
✓ Branch 2 taken 73340 times.
✓ Branch 3 taken 29336 times.
|
205352 | for (const auto& scv : scvs(fvGeometry)) |
| 203 | { | ||
| 204 | 80674 | const auto dofIdxGlobal = scv.dofIndex(); | |
| 205 | 117344 | vtkBoxVolume_[dofIdxGlobal] += scv.volume(); | |
| 206 | #if ISOTHERMAL | ||
| 207 | 220020 | vtkTemperature_[dofIdxGlobal] = this->spatialParams().temperatureAtPos(scv.dofPosition()); | |
| 208 | #endif | ||
| 209 | } | ||
| 210 | |||
| 211 |
6/6✓ Branch 0 taken 3470 times.
✓ Branch 1 taken 33200 times.
✓ Branch 2 taken 3470 times.
✓ Branch 3 taken 33200 times.
✓ Branch 4 taken 3470 times.
✓ Branch 5 taken 33200 times.
|
110010 | vtkKxx_[eIdx] = this->spatialParams().permeability(eIdx); |
| 212 |
6/8✓ Branch 0 taken 3470 times.
✓ Branch 1 taken 33200 times.
✓ Branch 2 taken 3470 times.
✓ Branch 3 taken 33200 times.
✓ Branch 5 taken 36670 times.
✗ Branch 6 not taken.
✓ Branch 8 taken 36670 times.
✗ Branch 9 not taken.
|
146680 | vtkPorosity_[eIdx] = 1- this->spatialParams().inertVolumeFraction(eIdx); |
| 213 | } | ||
| 214 | |||
| 215 | // communicate box volume at process boundaries for vertex-centered scheme (box) | ||
| 216 | if constexpr (isBox) | ||
| 217 | { | ||
| 218 |
4/4✓ Branch 1 taken 4 times.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 4 times.
✓ Branch 4 taken 2 times.
|
6 | if (gridView.comm().size() > 1) |
| 219 | { | ||
| 220 | VectorCommDataHandleSum<typename GridGeometry::VertexMapper, std::vector<Scalar>, GridView::dimension> | ||
| 221 |
3/6✓ Branch 1 taken 4 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 4 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 4 times.
✗ Branch 8 not taken.
|
12 | sumVolumeHandle(this->gridGeometry().vertexMapper(), vtkBoxVolume_); |
| 222 | if constexpr (Detail::canCommunicate<typename GridView::Traits::Grid, GridView::dimension>) | ||
| 223 |
1/2✓ Branch 1 taken 4 times.
✗ Branch 2 not taken.
|
4 | gridView.communicate(sumVolumeHandle, |
| 224 | Dune::InteriorBorder_InteriorBorder_Interface, | ||
| 225 | Dune::ForwardCommunication); | ||
| 226 | else | ||
| 227 | DUNE_THROW(Dune::InvalidStateException, "Cannot call addFieldsToWriter on multiple processes for a grid that cannot communicate codim-" << GridView::dimension << "-entities."); | ||
| 228 | } | ||
| 229 | } | ||
| 230 | 14 | } | |
| 231 | |||
| 232 | /*! | ||
| 233 | * \name Problem parameters | ||
| 234 | */ | ||
| 235 | // \{ | ||
| 236 | |||
| 237 | /*! | ||
| 238 | * \brief The problem name. | ||
| 239 | * | ||
| 240 | * This is used as a prefix for files generated by the simulation. | ||
| 241 | */ | ||
| 242 | const std::string& name() const | ||
| 243 |
1/2✓ Branch 1 taken 14 times.
✗ Branch 2 not taken.
|
14 | { return name_; } |
| 244 | |||
| 245 | // \} | ||
| 246 | |||
| 247 | /*! | ||
| 248 | * \name Boundary conditions | ||
| 249 | */ | ||
| 250 | // \{ | ||
| 251 | |||
| 252 | /*! | ||
| 253 | * \brief Specifies which kind of boundary condition should be | ||
| 254 | * used for which equation on a given boundary segment. | ||
| 255 | * | ||
| 256 | * \param element The finite element | ||
| 257 | * \param scv The sub-control volume | ||
| 258 | */ | ||
| 259 | ✗ | BoundaryTypes boundaryTypes(const Element &element, | |
| 260 | const SubControlVolume &scv) const | ||
| 261 |
12/12✓ Branch 0 taken 320780 times.
✓ Branch 1 taken 26896 times.
✓ Branch 2 taken 320780 times.
✓ Branch 3 taken 26896 times.
✓ Branch 4 taken 135463 times.
✓ Branch 5 taken 47888 times.
✓ Branch 6 taken 135463 times.
✓ Branch 7 taken 47888 times.
✓ Branch 8 taken 9072 times.
✓ Branch 9 taken 47152 times.
✓ Branch 10 taken 9072 times.
✓ Branch 11 taken 47152 times.
|
1174502 | { return scvfToScvBoundaryTypes_.boundaryTypes(scv); } |
| 262 | |||
| 263 | /*! | ||
| 264 | * \brief Specifies which kind of boundary condition should be | ||
| 265 | * used for which equation on a given boundary segment. | ||
| 266 | * | ||
| 267 | * \param element The finite element | ||
| 268 | * \param scvf The sub-control volume face | ||
| 269 | */ | ||
| 270 | 2846286 | BoundaryTypes boundaryTypes(const Element &element, | |
| 271 | const SubControlVolumeFace &scvf) const | ||
| 272 | { | ||
| 273 | 2846286 | BoundaryTypes bcTypes; | |
| 274 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 2846286 times.
|
2846286 | const auto boundaryId = scvf.boundaryFlag(); |
| 275 | |||
| 276 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 2846286 times.
|
2846286 | if (boundaryId < 1 || boundaryId > 4) |
| 277 | ✗ | DUNE_THROW(Dune::InvalidStateException, "Invalid boundary ID: " << boundaryId); | |
| 278 | |||
| 279 |
2/2✓ Branch 0 taken 468090 times.
✓ Branch 1 taken 2378196 times.
|
2846286 | if (boundaryId == dirichletBoundary_) |
| 280 | bcTypes.setAllDirichlet(); | ||
| 281 | else | ||
| 282 | bcTypes.setAllNeumann(); | ||
| 283 | |||
| 284 | 2846286 | return bcTypes; | |
| 285 | } | ||
| 286 | |||
| 287 | /*! | ||
| 288 | * \brief Evaluates the boundary conditions for a Dirichlet boundary segment. | ||
| 289 | * | ||
| 290 | * \return the Dirichlet values for the conservation equations in | ||
| 291 | * \f$ [ \textnormal{unit of primary variable} ] \f$ | ||
| 292 | * \param globalPos The global position | ||
| 293 | */ | ||
| 294 | PrimaryVariables dirichletAtPos(const GlobalPosition &globalPos) const | ||
| 295 |
1/2✓ Branch 1 taken 90250 times.
✗ Branch 2 not taken.
|
222520 | { return initial_(globalPos); } |
| 296 | |||
| 297 | /*! | ||
| 298 | * \brief Evaluates the boundary conditions for a Neumann boundary segment. | ||
| 299 | * | ||
| 300 | * This is the method for the case where the Neumann condition is | ||
| 301 | * potentially solution dependent and requires some quantities that | ||
| 302 | * are specific to the fully-implicit method. | ||
| 303 | * | ||
| 304 | * \param element The finite element | ||
| 305 | * \param fvGeometry The finite-volume geometry | ||
| 306 | * \param elemVolVars All volume variables for the element | ||
| 307 | * \param elemFluxVarsCache Flux variables caches for all faces in stencil | ||
| 308 | * \param scvf The sub-control volume face | ||
| 309 | * | ||
| 310 | * For this method, the \a values parameter stores the flux | ||
| 311 | * in normal direction of each phase. Negative values mean influx. | ||
| 312 | * E.g. for the mass balance that would the mass flux in \f$ [ kg / (m^2 \cdot s)] \f$. | ||
| 313 | */ | ||
| 314 | ✗ | NumEqVector neumann(const Element& element, | |
| 315 | const FVElementGeometry& fvGeometry, | ||
| 316 | const ElementVolumeVariables& elemVolVars, | ||
| 317 | const ElementFluxVariablesCache& elemFluxVarsCache, | ||
| 318 | const SubControlVolumeFace& scvf) const | ||
| 319 | { | ||
| 320 | 2231280 | const auto boundaryId = scvf.boundaryFlag(); | |
| 321 | |||
| 322 | 2231280 | NumEqVector fluxes(0.0); | |
| 323 | // kg/(m^2*s) or mole/(m^2*s) depending on useMoles | ||
| 324 |
4/6✗ Branch 0 not taken.
✗ Branch 1 not taken.
✓ Branch 2 taken 90432 times.
✓ Branch 3 taken 1102140 times.
✓ Branch 4 taken 80772 times.
✓ Branch 5 taken 957936 times.
|
2231280 | if (boundaryId == injectionBottom_) |
| 325 | { | ||
| 326 | 342408 | fluxes[contiCO2EqIdx] = useMoles ? -injectionRate_/FluidSystem::molarMass(CO2Idx) : -injectionRate_; | |
| 327 | #if !ISOTHERMAL | ||
| 328 | // energy fluxes are always mass specific | ||
| 329 | ✗ | fluxes[energyEqIdx] = -injectionRate_/*kg/(m^2 s)*/*CO2::gasEnthalpy( | |
| 330 | ✗ | injectionTemperature_, injectionPressure_)/*J/kg*/; // W/(m^2) | |
| 331 | #endif | ||
| 332 | } | ||
| 333 | |||
| 334 | 2231280 | return fluxes; | |
| 335 | } | ||
| 336 | |||
| 337 | // \} | ||
| 338 | |||
| 339 | /*! | ||
| 340 | * \name Volume terms | ||
| 341 | */ | ||
| 342 | // \{ | ||
| 343 | |||
| 344 | /*! | ||
| 345 | * \brief Evaluates the initial values at a position. | ||
| 346 | * | ||
| 347 | * \param globalPos The global position | ||
| 348 | * \return The initial values for the conservation equations in | ||
| 349 | * \f$ [ \textnormal{unit of primary variables} ] \f$ | ||
| 350 | */ | ||
| 351 | PrimaryVariables initialAtPos(const GlobalPosition &globalPos) const | ||
| 352 | { | ||
| 353 |
1/2✓ Branch 1 taken 38566 times.
✗ Branch 2 not taken.
|
38566 | return initial_(globalPos); |
| 354 | } | ||
| 355 | |||
| 356 | // \} | ||
| 357 | |||
| 358 | private: | ||
| 359 | /*! | ||
| 360 | * \brief Evaluates the initial values for a control volume. | ||
| 361 | * | ||
| 362 | * The internal method for the initial condition | ||
| 363 | * | ||
| 364 | * \param globalPos The global position | ||
| 365 | */ | ||
| 366 | 261086 | PrimaryVariables initial_(const GlobalPosition &globalPos) const | |
| 367 | { | ||
| 368 | 261086 | PrimaryVariables values(0.0); | |
| 369 | 261086 | values.setState(firstPhaseOnly); | |
| 370 | |||
| 371 | 522172 | const Scalar temp = this->spatialParams().temperatureAtPos(globalPos); | |
| 372 | 261086 | const Scalar densityW = FluidSystem::Brine::liquidDensity(temp, 1e7); | |
| 373 | |||
| 374 | 261086 | const Scalar moleFracLiquidCO2 = 0.00; | |
| 375 | 261086 | const Scalar moleFracLiquidBrine = 1.0 - moleFracLiquidCO2; | |
| 376 | |||
| 377 | 261086 | const Scalar meanM = FluidSystem::molarMass(BrineIdx)*moleFracLiquidBrine | |
| 378 | 261086 | + FluidSystem::molarMass(CO2Idx)*moleFracLiquidCO2; | |
| 379 | |||
| 380 | if(useMoles) // mole-fraction formulation | ||
| 381 | values[switchIdx] = moleFracLiquidCO2; | ||
| 382 | else // mass-fraction formulation | ||
| 383 | 261086 | values[switchIdx] = moleFracLiquidCO2*FluidSystem::molarMass(CO2Idx)/meanM; | |
| 384 | |||
| 385 | 1566516 | values[pressureIdx] = 1.0e5 - densityW*this->spatialParams().gravity(globalPos)[dimWorld-1]*(depthBOR_ - globalPos[dimWorld-1]); | |
| 386 | |||
| 387 | #if !ISOTHERMAL | ||
| 388 | 128584 | values[temperatureIdx] = temp; | |
| 389 | #endif | ||
| 390 | 261086 | return values; | |
| 391 | } | ||
| 392 | |||
| 393 | Scalar depthBOR_; | ||
| 394 | Scalar injectionRate_; | ||
| 395 | |||
| 396 | #if !ISOTHERMAL | ||
| 397 | Scalar injectionPressure_, injectionTemperature_; | ||
| 398 | #endif | ||
| 399 | |||
| 400 | int nTemperature_; | ||
| 401 | int nPressure_; | ||
| 402 | |||
| 403 | std::string name_ ; | ||
| 404 | |||
| 405 | Scalar pressureLow_, pressureHigh_; | ||
| 406 | Scalar temperatureLow_, temperatureHigh_; | ||
| 407 | |||
| 408 | int injectionTop_; | ||
| 409 | int injectionBottom_; | ||
| 410 | int dirichletBoundary_; | ||
| 411 | int noFlowBoundary_; | ||
| 412 | |||
| 413 | // vtk output | ||
| 414 | std::vector<Scalar> vtkKxx_, vtkPorosity_, vtkBoxVolume_, vtkTemperature_; | ||
| 415 | ScvfToScvBoundaryTypes<BoundaryTypes, DiscretizationMethod> scvfToScvBoundaryTypes_; | ||
| 416 | }; | ||
| 417 | |||
| 418 | } // end namespace Dumux | ||
| 419 | |||
| 420 | #endif | ||
| 421 |