My Project
WellConnectionAuxiliaryModule.hpp
1/*
2 Copyright 2017 Dr. Blatt - HPC-Simulation-Software & Services
3 Copyright 2017 Statoil ASA.
4
5 This file is part of the Open Porous Media project (OPM).
6
7 OPM is free software: you can redistribute it and/or modify
8 it under the terms of the GNU General Public License as published by
9 the Free Software Foundation, either version 3 of the License, or
10 (at your option) any later version.
11
12 OPM is distributed in the hope that it will be useful,
13 but WITHOUT ANY WARRANTY; without even the implied warranty of
14 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15 GNU General Public License for more details.
16
17 You should have received a copy of the GNU General Public License
18 along with OPM. If not, see <http://www.gnu.org/licenses/>.
19*/
20
21#ifndef OPM_WELLCONNECTIONAUXILIARYMODULE_HEADER_INCLUDED
22#define OPM_WELLCONNECTIONAUXILIARYMODULE_HEADER_INCLUDED
23
24#include <opm/models/discretization/common/baseauxiliarymodule.hh>
25
26#include <opm/grid/CpGrid.hpp>
27
28#include <opm/input/eclipse/EclipseState/EclipseState.hpp>
29
30namespace Opm
31{
32template<class TypeTag>
34 : public BaseAuxiliaryModule<TypeTag>
35{
36 using GlobalEqVector = GetPropType<TypeTag, Properties::GlobalEqVector>;
37 using SparseMatrixAdapter = GetPropType<TypeTag, Properties::SparseMatrixAdapter>;
38
39public:
40
41 using NeighborSet = typename
42 ::Opm::BaseAuxiliaryModule<TypeTag>::NeighborSet;
43
44 WellConnectionAuxiliaryModule(const Schedule& schedule,
45 const Dune::CpGrid& grid)
46 {
47 // Create cartesian to compressed mapping
48 const auto& globalCell = grid.globalCell();
49 const auto& cartesianSize = grid.logicalCartesianSize();
50
51 auto size = cartesianSize[0]*cartesianSize[1]*cartesianSize[2];
52
53 std::vector<int> cartesianToCompressed(size, -1);
54 auto begin = globalCell.begin();
55
56 for ( auto cell = begin, end= globalCell.end(); cell != end; ++cell )
57 {
58 cartesianToCompressed[ *cell ] = cell - begin;
59 }
60
61 const auto& schedule_wells = schedule.getWellsatEnd();
62 wells_.reserve(schedule_wells.size());
63
64 // initialize the additional cell connections introduced by wells.
65 for ( const auto& well : schedule_wells )
66 {
67 std::vector<int> compressed_well_perforations;
68 // All possible completions of the well
69 const auto& completionSet = well.getConnections();
70 compressed_well_perforations.reserve(completionSet.size());
71
72 for ( size_t c=0; c < completionSet.size(); c++ )
73 {
74 const auto& completion = completionSet.get(c);
75 int i = completion.getI();
76 int j = completion.getJ();
77 int k = completion.getK();
78 int cart_grid_idx = i + cartesianSize[0]*(j + cartesianSize[1]*k);
79 int compressed_idx = cartesianToCompressed[cart_grid_idx];
80
81 if ( compressed_idx >= 0 ) // Ignore completions in inactive/remote cells.
82 {
83 compressed_well_perforations.push_back(compressed_idx);
84 }
85 }
86
87 if ( ! compressed_well_perforations.empty() )
88 {
89 std::sort(compressed_well_perforations.begin(),
90 compressed_well_perforations.end());
91
92 wells_.push_back(compressed_well_perforations);
93 }
94 }
95 }
96
97 unsigned numDofs() const
98 {
99 // No extra dofs are inserted for wells.
100 return 0;
101 }
102
103 void addNeighbors(std::vector<NeighborSet>& neighbors) const
104 {
105 for(const auto& well_perforations : wells_)
106 {
107 for(const auto& perforation : well_perforations)
108 neighbors[perforation].insert(well_perforations.begin(),
109 well_perforations.end());
110 }
111 }
112
113 void applyInitial()
114 {}
115
116 void linearize(SparseMatrixAdapter& , GlobalEqVector&)
117 {
118 // Linearization is done in StandardDenseWells
119 }
120
121 private:
122
123 std::vector<std::vector<int> > wells_;
124};
125
126} // end namespace OPM
127#endif
Definition: WellConnectionAuxiliaryModule.hpp:35
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: BlackoilPhases.hpp:27