Added unit conversion for MSW data pluss correctin for complum

This commit is contained in:
Jostein Alvestad 2018-07-12 12:25:24 +02:00
parent 9d0604ad99
commit 17966d098b
4 changed files with 48 additions and 21 deletions

View File

@ -30,6 +30,7 @@
namespace Opm {
class Phases;
class Schedule;
class UnitSystem;
} // Opm
namespace Opm { namespace RestartIO { namespace Helpers {
@ -49,6 +50,7 @@ namespace Opm { namespace RestartIO { namespace Helpers {
void captureDeclaredMSWData(const Opm::Schedule& sched,
const std::size_t rptStep,
const Opm::UnitSystem& units,
const std::vector<int>& inteHead);
/// Retrieve Integer Multisegment well data Array.

View File

@ -30,6 +30,8 @@
#include <cassert>
#include <cstddef>
#include <iostream>
#include <cmath>
namespace VI = Opm::RestartIO::Helpers::VectorItems;
@ -52,9 +54,12 @@ namespace {
{
// make seqIndex to Connection map
std::map <std::size_t, const Opm::Connection*> seqIndConnMap;
for (const auto conn : conns) {
for (const auto & conn : conns) {
//for (auto it = conns.begin(); it < conns.end(); it++) {
std::size_t sI = conn.getSeqIndex();
seqIndConnMap.insert(std::make_pair(sI, &conn));
std::cout << "seqIndex: " << sI << " i: " << conn.getI() << " j: " << conn.getJ() << " i: " << conn.getK() <<
" &conn->getSeqIndex(): " << (&conn)->getSeqIndex() << std::endl;
}
return seqIndConnMap;
}
@ -63,9 +68,11 @@ namespace {
{
// make CompSegSeqIndex to Connection map
std::map <std::size_t, const Opm::Connection*> cs_seqIndConnMap;
for (const auto conn : conns) {
for (const auto & conn : conns) {
std::size_t sI = conn.getCompSegSeqIndex();
cs_seqIndConnMap.insert(std::make_pair(sI, &conn));
std::cout << "getCompSegSeqIndex: " << sI << " i: " << conn.getI() << " j: " << conn.getJ() << " i: " << conn.getK() <<
" &conn->getCompSegSeqIndex(): " << (&conn)->getCompSegSeqIndex() << std::endl;
}
return cs_seqIndConnMap;
}
@ -97,19 +104,29 @@ namespace {
//sort connections according to input sequence in COMPDAT
sIToConn = mapSeqIndexToConnection(conns);
}
for (auto con : conns) {
std::cout << "loop act conns - seqIndex: " << con.getSeqIndex() << " i: " << con.getI() << " j: " << con.getJ() << " i: " << con.getK() << std::endl;
}
for (auto conp : sIToConn) {
std::cout << "loop sIToConn - seqIndex: " << conp.first << " conp.second->getSeqIndex(): " << conp.second->getSeqIndex() << std::endl;
}
std::vector<const Opm::Connection*> connSI;
int niSI = well->getConnections(sim_step).size();
std::cout << "No of Connections: " << niSI << std::endl;
for (int iSI = 0; iSI < niSI; iSI++) {
const auto searchSI = sIToConn.find(static_cast<std::size_t>(iSI));
std::cout << "Loop over connections, iSI: " << iSI << " searchSI->first: " << searchSI->first
<< " searchSI->second->getSeqIndex(): " << searchSI->second->getSeqIndex() << std::endl;
if (searchSI != sIToConn.end()) {
connSI.push_back(searchSI->second);
}
std::cout << "Store connection: " << searchSI->second->getSeqIndex() << " searchSI->first: " << searchSI->first << std::endl;
}
}
std::cout << "connSI-size " << connSI.size() << std::endl;
for (auto nConn = connSI.size(), connID = 0*nConn;
connID < nConn; ++connID)
{
connOp(*well, wellID, *connSI[connID], connID);
connOp(*well, wellID, *(connSI[connID]), connID);
}
}
}
@ -154,7 +171,7 @@ namespace {
// draining and imbibition curves at connections.
iConn[9] = conn.sat_tableId;
iConn[12] = conn.complnum;
iConn[12] = std::abs(conn.complnum);
iConn[13] = conn.dir;
iConn[14] = conn.attachedToSegment()
? conn.segment_number : 0;
@ -243,9 +260,10 @@ captureDeclaredConnData(const Schedule& sched,
(const Well& /* well */, const std::size_t wellID,
const Connection& conn , const std::size_t connID) -> void
{
std::cout << "Connection loop - wellID: " << wellID << " connID: " << connID << std::endl;
auto ic = this->iConn_(wellID, connID);
auto sc = this->sConn_(wellID, connID);
std::cout << "Conn loop - conn- seqIndex: " << conn.getSeqIndex() << std::endl;
IConn::staticContrib(conn, connID, ic);
SConn::staticContrib(conn, units, sc);
});

View File

@ -29,6 +29,8 @@
#include <opm/parser/eclipse/EclipseState/Schedule/Connection.hpp>
#include <opm/parser/eclipse/EclipseState/Schedule/WellConnections.hpp>
#include <opm/parser/eclipse/Units/UnitSystem.hpp>
//#include <opm/output/data/Wells.hpp>
#include <algorithm>
@ -351,9 +353,11 @@ namespace {
void staticContrib(const Opm::Well& well,
const std::size_t rptStep,
const std::vector<int>& inteHead,
const Opm::UnitSystem& units,
RSegArray& rSeg)
{
if (well.isMultiSegment(rptStep)) {
using M = ::Opm::UnitSystem::measure;
//loop over segment set and print out information
auto welSegSet = well.getWellSegments(rptStep);
auto completionSet = well.getCompletions(rptStep);
@ -361,9 +365,9 @@ namespace {
//treat the top segment individually
int ind_seg = 1;
auto ind = welSegSet.segmentNumberToIndex(ind_seg);
rSeg[0] = welSegSet.lengthTopSegment();
rSeg[1] = welSegSet.depthTopSegment();
rSeg[5] = welSegSet.volumeTopSegment();
rSeg[0] = units.from_si(M::length, welSegSet.lengthTopSegment());
rSeg[1] = units.from_si(M::length, welSegSet.depthTopSegment());
rSeg[5] = units.from_si(M::volume, welSegSet.volumeTopSegment());
rSeg[6] = rSeg[0];
rSeg[7] = rSeg[1];
@ -387,14 +391,16 @@ namespace {
auto outSeg = welSegSet[ind].outletSegment();
auto ind_ofs = welSegSet.segmentNumberToIndex(outSeg);
auto iS = (ind_seg-1)*noElmSeg;
rSeg[iS + 0] = welSegSet[ind].totalLength() - welSegSet[ind_ofs].totalLength();
rSeg[iS + 1] = welSegSet[ind].depth() - welSegSet[ind_ofs].depth();
rSeg[iS + 2] = welSegSet[ind].internalDiameter();
rSeg[iS + 3] = welSegSet[ind].roughness();
rSeg[iS + 4] = welSegSet[ind].crossArea();
rSeg[iS + 5] = welSegSet[ind].volume();
rSeg[iS + 6] = welSegSet[ind].totalLength();
rSeg[iS + 7] = welSegSet[ind].depth();
rSeg[iS + 0] = units.from_si(M::length, (welSegSet[ind].totalLength() - welSegSet[ind_ofs].totalLength()));
rSeg[iS + 1] = units.from_si(M::length, (welSegSet[ind].depth() - welSegSet[ind_ofs].depth()));
rSeg[iS + 2] = units.from_si(M::length, (welSegSet[ind].internalDiameter()));
rSeg[iS + 3] = units.from_si(M::length, (welSegSet[ind].roughness()));
rSeg[iS + 4] = units.from_si(M::length, (welSegSet[ind].crossArea()));
//repeat unit conversion to convert for area instead of length
rSeg[iS + 4] = units.from_si(M::length, rSeg[iS + 4]);
rSeg[iS + 5] = units.from_si(M::volume, (welSegSet[ind].volume()));
rSeg[iS + 6] = units.from_si(M::length, (welSegSet[ind].totalLength()));
rSeg[iS + 7] = units.from_si(M::length, (welSegSet[ind].depth()));
// segment pressure (to be added!!)
rSeg[iS + 11] = 0;
@ -517,6 +523,7 @@ void
Opm::RestartIO::Helpers::AggregateMSWData::
captureDeclaredMSWData(const Schedule& sched,
const std::size_t rptStep,
const Opm::UnitSystem& units,
const std::vector<int>& inteHead)
{
const auto& wells = sched.getWells(rptStep);
@ -540,12 +547,12 @@ captureDeclaredMSWData(const Schedule& sched,
// Extract Contributions to RSeg Array
{
MSWLoop(msw, [rptStep, inteHead, this]
MSWLoop(msw, [&units, rptStep, inteHead, this]
(const Well& well, const std::size_t mswID) -> void
{
auto rmsw = this->rSeg_[mswID];
RSeg::staticContrib(well, rptStep, inteHead, rmsw);
RSeg::staticContrib(well, rptStep, inteHead, units, rmsw);
});
}

View File

@ -510,7 +510,7 @@ void writeMSWData(::Opm::RestartIO::ecl_rst_file_type * rst_file,
// write ISEG, RSEG, ILBS and ILBR to restart file
const size_t simStep = static_cast<size_t> (sim_step);
auto MSWData = Helpers::AggregateMSWData(ih);
MSWData.captureDeclaredMSWData(schedule, simStep, ih);
MSWData.captureDeclaredMSWData(schedule, simStep, units, ih);
write_kw(rst_file, EclKW<int> ("ISEG", MSWData.getISeg()));
write_kw(rst_file, EclKW<int> ("ILBS", MSWData.getILBs()));
write_kw(rst_file, EclKW<int> ("ILBR", MSWData.getILBr()));