mirror of
https://github.com/OPM/ResInsight.git
synced 2025-01-01 03:37:15 -06:00
158 lines
4.6 KiB
C++
158 lines
4.6 KiB
C++
/*
|
|
Copyright 2016 SINTEF ICT, Applied Mathematics.
|
|
Copyright 2016 Statoil ASA.
|
|
|
|
This file is part of the Open Porous Media project (OPM).
|
|
|
|
OPM is free software: you can redistribute it and/or modify
|
|
it under the terms of the GNU General Public License as published by
|
|
the Free Software Foundation, either version 3 of the License, or
|
|
(at your option) any later version.
|
|
|
|
OPM is distributed in the hope that it will be useful,
|
|
but WITHOUT ANY WARRANTY; without even the implied warranty of
|
|
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
|
|
GNU General Public License for more details.
|
|
|
|
You should have received a copy of the GNU General Public License
|
|
along with OPM. If not, see <http://www.gnu.org/licenses/>.
|
|
*/
|
|
|
|
#ifndef OPM_EXAMPLESETUP_HEADER_INCLUDED
|
|
#define OPM_EXAMPLESETUP_HEADER_INCLUDED
|
|
|
|
|
|
|
|
#include <opm/flowdiagnostics/ConnectivityGraph.hpp>
|
|
#include <opm/flowdiagnostics/ConnectionValues.hpp>
|
|
#include <opm/flowdiagnostics/Toolbox.hpp>
|
|
|
|
#include <opm/utility/ECLGraph.hpp>
|
|
#include <opm/utility/ECLWellSolution.hpp>
|
|
|
|
#include <exception>
|
|
#include <iostream>
|
|
#include <sstream>
|
|
#include <stdexcept>
|
|
#include <string>
|
|
#include <utility>
|
|
#include <vector>
|
|
|
|
#include <boost/filesystem.hpp>
|
|
|
|
namespace example {
|
|
|
|
|
|
inline Opm::FlowDiagnostics::ConnectionValues
|
|
extractFluxField(const Opm::ECLGraph& G)
|
|
{
|
|
using ConnVals = Opm::FlowDiagnostics::ConnectionValues;
|
|
|
|
using NConn = ConnVals::NumConnections;
|
|
using NPhas = ConnVals::NumPhases;
|
|
|
|
const auto nconn = NConn{G.numConnections()};
|
|
const auto nphas = NPhas{3};
|
|
|
|
auto flux = ConnVals(nconn, nphas);
|
|
|
|
auto phas = ConnVals::PhaseID{0};
|
|
|
|
for (const auto& p : { Opm::ECLGraph::PhaseIndex::Aqua ,
|
|
Opm::ECLGraph::PhaseIndex::Liquid ,
|
|
Opm::ECLGraph::PhaseIndex::Vapour })
|
|
{
|
|
const std::vector<double> pflux = G.flux(p);
|
|
|
|
if (! pflux.empty()) {
|
|
assert (pflux.size() == nconn.total);
|
|
|
|
auto conn = ConnVals::ConnID{0};
|
|
|
|
for (const auto& v : pflux) {
|
|
flux(conn, phas) = v;
|
|
|
|
conn.id += 1;
|
|
}
|
|
}
|
|
|
|
phas.id += 1;
|
|
}
|
|
|
|
return flux;
|
|
}
|
|
|
|
template <class WellFluxes>
|
|
Opm::FlowDiagnostics::CellSetValues
|
|
extractWellFlows(const Opm::ECLGraph& G,
|
|
const WellFluxes& well_fluxes)
|
|
{
|
|
Opm::FlowDiagnostics::CellSetValues inflow;
|
|
for (const auto& well : well_fluxes) {
|
|
for (const auto& completion : well.completions) {
|
|
const int grid_index = completion.grid_index;
|
|
const auto& ijk = completion.ijk;
|
|
const int cell_index = G.activeCell(ijk, grid_index);
|
|
if (cell_index >= 0) {
|
|
inflow.emplace(cell_index, completion.reservoir_inflow_rate);
|
|
}
|
|
}
|
|
}
|
|
|
|
return inflow;
|
|
}
|
|
|
|
namespace Hack {
|
|
inline Opm::FlowDiagnostics::ConnectionValues
|
|
convert_flux_to_SI(Opm::FlowDiagnostics::ConnectionValues&& fl)
|
|
{
|
|
using Co = Opm::FlowDiagnostics::ConnectionValues::ConnID;
|
|
using Ph = Opm::FlowDiagnostics::ConnectionValues::PhaseID;
|
|
|
|
const auto nconn = fl.numConnections();
|
|
const auto nphas = fl.numPhases();
|
|
|
|
for (auto phas = Ph{0}; phas.id < nphas; ++phas.id) {
|
|
for (auto conn = Co{0}; conn.id < nconn; ++conn.id) {
|
|
fl(conn, phas) /= 86400;
|
|
}
|
|
}
|
|
|
|
return fl;
|
|
}
|
|
}
|
|
|
|
inline Opm::FlowDiagnostics::Toolbox
|
|
initialiseFlowDiagnostics(const Opm::ECLGraph& G)
|
|
{
|
|
const Opm::FlowDiagnostics::ConnectivityGraph connGraph =
|
|
Opm::FlowDiagnostics::ConnectivityGraph{ static_cast<int>(G.numCells()),
|
|
G.neighbours() };
|
|
|
|
// Create the Toolbox.
|
|
|
|
Opm::FlowDiagnostics::Toolbox tool = Opm::FlowDiagnostics::Toolbox{ connGraph };
|
|
|
|
tool.assignPoreVolume(G.poreVolume());
|
|
|
|
Opm::FlowDiagnostics::ConnectionValues connectionsVals = Hack::convert_flux_to_SI(extractFluxField(G));
|
|
tool.assignConnectionFlux(connectionsVals);
|
|
|
|
Opm::ECLWellSolution wsol = Opm::ECLWellSolution{};
|
|
|
|
const std::vector<Opm::ECLWellSolution::WellData> well_fluxes =
|
|
wsol.solution(G.rawResultData(), G.numGrids());
|
|
|
|
tool.assignInflowFlux(extractWellFlows(G, well_fluxes));
|
|
|
|
return tool;
|
|
}
|
|
|
|
|
|
|
|
} // namespace example
|
|
|
|
|
|
|
|
#endif // OPM_EXAMPLESETUP_HEADER_INCLUDED
|