Add New RFT File Writer

This commit implements a new file writing function

    void Opm::RftIO::write()

that is intended to replace the current RFT output functionality
defined in terms of private class 'RFT' in EclipseIO.cpp.  We
support the basic RFT output consisting of

    - Timestamp (elapsed and date)
    - WELLETC metadata including all unit conventions
    - Connection cell (I,J,K), connection cell hostgrid (blank for
      main grid only), connection cell centre depth, connection cell
      pressure, and connection cell water and gas saturations

Connections in inactive cells are omitted.  Note that unit of
measure strings aren't implemented in terms of UnitSystem::name()
due to the strings being padded on the left for centering effect.

Add unit tests to exercise the new writer.
This commit is contained in:
Bård Skaflestad
2019-08-28 10:31:26 +02:00
parent a281a3c6f4
commit 6c47cdb22f
4 changed files with 2492 additions and 0 deletions

View File

@@ -205,6 +205,7 @@ if(ENABLE_ECL_OUTPUT)
src/opm/output/eclipse/RegionCache.cpp
src/opm/output/eclipse/RestartValue.cpp
src/opm/output/eclipse/WriteInit.cpp
src/opm/output/eclipse/WriteRFT.cpp
src/opm/output/data/Solution.cpp
)
endif()
@@ -642,6 +643,7 @@ if(ENABLE_ECL_OUTPUT)
opm/output/eclipse/Tables.hpp
opm/output/eclipse/WindowedArray.hpp
opm/output/eclipse/WriteInit.hpp
opm/output/eclipse/WriteRFT.hpp
opm/output/eclipse/WriteRestartHelpers.hpp
opm/output/OutputWriter.hpp
)

View File

@@ -0,0 +1,89 @@
/*
Copyright (c) 2019 Equinor 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_WRITE_RFT_HPP
#define OPM_WRITE_RFT_HPP
namespace Opm {
class EclipseGrid;
class Schedule;
class UnitSystem;
} // namespace Opm
namespace Opm { namespace data {
class WellRates;
}} // namespace Opm::data
namespace Opm { namespace EclIO { namespace OutputStream {
class RFT;
}}} // namespace Opm::EclIO::OutputStream
namespace Opm { namespace RftIO {
/// Collect RFT data and output to pre-opened output stream.
///
/// RFT data is output for all affected wells at a given timestep,
/// and consists of
///
/// 1) Time stamp (elapsed and date)
/// 2) Metadata about the output (units of measure, well types,
/// data type identifier)
/// 3) Depth, pressure, Sw, Sg, (I,J,K), and hostgrid for each
/// reservoir connection of the affected well.
///
/// If no RFT output is requested at given timestep, then this
/// function returns with no output written to the RFT file.
///
/// \param[in] reportStep Report step time point for which to output
/// RFT data.
///
/// \param[in] elapsed Number of seconds of simulated time until
/// this report step (\p reportStep).
///
/// \param[in] usys Unit system conventions for output. Typically
/// corresponds to run's active unit system (i.e., \code
/// EclipseState::getUnits() \endcode).
///
/// \param[in] grid Run's active grid. Used to determine which
/// reservoir connections are in active grid cells.
///
/// \param[in] schedule Run's SCHEDULE section from which to access
/// the active wells and the RFT configuration.
///
/// \param[in,out] rftFile RFT output stream. On input, appropriately
/// positioned for new content (i.e., at end-of-file). On output,
/// containing new RFT output (if applicable) and positioned after
/// new contents.
void write(const int reportStep,
const double elapsed,
const ::Opm::UnitSystem& usys,
const ::Opm::EclipseGrid& grid,
const ::Opm::Schedule& schedule,
const ::Opm::data::WellRates& wellSol,
::Opm::EclIO::OutputStream::RFT& rftFile);
}} // namespace Opm::RftIO
#endif // OPM_WRITE_RFT_HPP

View File

@@ -0,0 +1,413 @@
/*
Copyright (c) 2019 Equinor ASA
Copyright (c) 2016 Statoil ASA
Copyright (c) 2013-2015 Andreas Lauser
Copyright (c) 2013 SINTEF ICT, Applied Mathematics.
Copyright (c) 2013 Uni Research AS
Copyright (c) 2015 IRIS AS
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/>.
*/
#include <opm/output/eclipse/WriteRFT.hpp>
#include <opm/io/eclipse/OutputStream.hpp>
#include <opm/io/eclipse/PaddedOutputString.hpp>
#include <opm/output/data/Wells.hpp>
#include <opm/output/eclipse/InteHEAD.hpp>
#include <opm/parser/eclipse/EclipseState/Grid/EclipseGrid.hpp>
#include <opm/parser/eclipse/EclipseState/Schedule/Schedule.hpp>
#include <opm/parser/eclipse/EclipseState/Schedule/Well/Connection.hpp>
#include <opm/parser/eclipse/Units/UnitSystem.hpp>
#include <algorithm>
#include <array>
#include <cstddef>
#include <initializer_list>
#include <stdexcept>
#include <string>
#include <vector>
namespace {
enum etcIx : std::vector<Opm::EclIO::PaddedOutputString<8>>::size_type
{
Time = 0,
Well = 1,
LGR = 2,
Depth = 3,
Pressure = 4,
DataType = 5,
WellType = 6,
LiqRate = 7,
GasRate = 8,
ResVRate = 9,
Velocity = 10,
Reserved = 11, // Untouched
Viscosity = 12,
ConcPlyBr = 13,
PlyBrRate = 14,
PlyBrAds = 15,
};
namespace RftUnits {
namespace exceptions {
void metric(std::vector<Opm::EclIO::PaddedOutputString<8>>& welletc)
{
welletc[etcIx::Depth] = " METRES";
welletc[etcIx::Velocity] = " M/SEC";
}
void field(std::vector<Opm::EclIO::PaddedOutputString<8>>& welletc)
{
welletc[etcIx::Depth] = " FEET";
welletc[etcIx::Velocity] = " FT/SEC";
welletc[etcIx::PlyBrRate] = " LB/DAY";
}
void lab(std::vector<Opm::EclIO::PaddedOutputString<8>>& welletc)
{
welletc[etcIx::Time] = " HR";
welletc[etcIx::Pressure] = " ATMA";
welletc[etcIx::Velocity] = " CM/SEC";
welletc[etcIx::ConcPlyBr] = " GM/SCC";
welletc[etcIx::PlyBrRate] = " GM/HR";
welletc[etcIx::PlyBrAds] = " GM/GM";
}
void pvt_m(std::vector<Opm::EclIO::PaddedOutputString<8>>& welletc)
{
// PVT_M is METRIC with pressures in atmospheres.
metric(welletc);
welletc[etcIx::Pressure] = " ATMA";
}
void input(std::vector<Opm::EclIO::PaddedOutputString<8>>& welletc)
{
welletc[etcIx::Time] = " INPUT";
welletc[etcIx::Depth] = " INPUT";
welletc[etcIx::Pressure] = " INPUT";
welletc[etcIx::LiqRate] = " INPUT";
welletc[etcIx::GasRate] = " INPUT";
welletc[etcIx::ResVRate] = " INPUT";
welletc[etcIx::Velocity] = " INPUT";
welletc[etcIx::Viscosity] = " INPUT";
welletc[etcIx::ConcPlyBr] = " INPUT";
welletc[etcIx::PlyBrRate] = " INPUT";
welletc[etcIx::PlyBrAds] = " INPUT";
}
} // namespace exceptions
std::string
centre(std::string s, const std::string::size_type width = 8)
{
if (s.size() > width) { return s.substr(0, width); }
if (s.size() == width) { return s; }
const auto npad = width - s.size();
const auto left = (npad + 1) / 2; // ceil(npad / 2)
return std::string(left, ' ') + s;
}
std::string combine(std::string left, std::string right)
{
return left + '/' + right;
}
void fill(const ::Opm::UnitSystem& usys,
std::vector<Opm::EclIO::PaddedOutputString<8>>& welletc)
{
using M = ::Opm::UnitSystem::measure;
welletc[etcIx::Time] = centre(std::string{usys.name(M::time)} + 'S');
welletc[etcIx::Depth] = centre(usys.name(M::length));
welletc[etcIx::Pressure] = centre(usys.name(M::pressure));
welletc[etcIx::LiqRate] = centre(usys.name(M::liquid_surface_rate));
welletc[etcIx::GasRate] = centre(usys.name(M::gas_surface_rate));
welletc[etcIx::ResVRate] = centre(usys.name(M::rate));
welletc[etcIx::Velocity] = centre(combine(usys.name(M::length),
usys.name(M::time)));
welletc[etcIx::Viscosity] = centre(usys.name(M::viscosity));
welletc[etcIx::ConcPlyBr] =
centre(combine(usys.name(M::mass),
usys.name(M::liquid_surface_volume)));
welletc[etcIx::PlyBrRate] = centre(usys.name(M::mass_rate));
welletc[etcIx::PlyBrAds] = centre(combine(usys.name(M::mass),
usys.name(M::mass)));
}
} // namespace RftUnits
class WellRFT
{
public:
explicit WellRFT(const std::size_t nconn = 0);
void add(const ::Opm::UnitSystem& usys,
const ::Opm::Connection& conn,
const ::Opm::data::Connection& xcon,
const double depth);
std::size_t nConn() const { return this->i_.size(); }
const std::vector<int>& conI() const { return this->i_; }
const std::vector<int>& conJ() const { return this->j_; }
const std::vector<int>& conK() const { return this->k_; }
const std::vector<float>& depth() const { return this->depth_; }
const std::vector<float>& pressure() const { return this->press_; }
const std::vector<float>& swat() const { return this->swat_; }
const std::vector<float>& sgas() const { return this->sgas_; }
const std::vector<Opm::EclIO::PaddedOutputString<8>>&
hostgrid() const
{
return this->host_;
}
private:
std::vector<int> i_;
std::vector<int> j_;
std::vector<int> k_;
std::vector<float> depth_;
std::vector<float> press_;
std::vector<float> swat_;
std::vector<float> sgas_;
std::vector<Opm::EclIO::PaddedOutputString<8>> host_;
};
WellRFT::WellRFT(const std::size_t nconn)
{
if (nconn == 0) { return; }
this->i_.reserve(nconn);
this->j_.reserve(nconn);
this->k_.reserve(nconn);
this->depth_.reserve(nconn);
this->press_.reserve(nconn);
this->swat_ .reserve(nconn);
this->sgas_ .reserve(nconn);
this->host_.reserve(nconn);
}
void WellRFT::add(const ::Opm::UnitSystem& usys,
const ::Opm::Connection& conn,
const ::Opm::data::Connection& xcon,
const double depth)
{
this->i_.push_back(conn.getI() + 1);
this->j_.push_back(conn.getJ() + 1);
this->k_.push_back(conn.getK() + 1);
using M = ::Opm::UnitSystem::measure;
auto cvrt = [&usys](const M meas, const double x) -> float
{
return usys.from_si(meas, x);
};
this->depth_.push_back(cvrt(M::length , depth));
this->press_.push_back(cvrt(M::pressure, xcon.cell_pressure));
this->swat_.push_back(xcon.cell_saturation_water);
this->sgas_.push_back(xcon.cell_saturation_gas);
this->host_.emplace_back();
}
WellRFT
createWellRFT(const int reportStep,
const std::string& wname,
const ::Opm::UnitSystem& usys,
const ::Opm::EclipseGrid& grid,
const ::Opm::Schedule& sched,
const std::vector<Opm::data::Connection>& xcon)
{
auto rft = WellRFT{ xcon.size() };
for (const auto& conn : sched.getWell2(wname, reportStep).getConnections()) {
const auto i = static_cast<std::size_t>(conn.getI());
const auto j = static_cast<std::size_t>(conn.getJ());
const auto k = static_cast<std::size_t>(conn.getK());
if (! grid.cellActive(i, j, k)) {
// Inactive cell. Ignore.
continue;
}
const auto ix = grid.getGlobalIndex(i, j, k);
auto xconPos = std::find_if(xcon.begin(), xcon.end(),
[ix](const ::Opm::data::Connection& c)
{
return c.index == ix;
});
if (xconPos == xcon.end()) {
// RFT data not available for this connection. Unexpected.
continue;
}
rft.add(usys, conn, *xconPos, grid.getCellDepth(ix));
}
return rft;
}
std::vector<Opm::EclIO::PaddedOutputString<8>>
wellETC(const std::string& wellName,
const std::string& dataType,
const ::Opm::UnitSystem& usys)
{
using UT = ::Opm::UnitSystem::UnitType;
auto ret = std::vector<Opm::EclIO::PaddedOutputString<8>>(16);
// Note: ret[etcIx::LGR] is well's LGR. Default constructed
// (i.e., blank) string is sufficient to represent no LGR.
ret[etcIx::Well] = wellName;
// 'P' -> PLT, 'R' -> RFT, 'S' -> Segment
ret[etcIx::DataType] = dataType;
// We support "standard" well type only.
ret[etcIx::WellType] = "STANDARD";
RftUnits::fill(usys, ret);
switch (usys.getType()) {
case UT::UNIT_TYPE_METRIC:
RftUnits::exceptions::metric(ret);
break;
case UT::UNIT_TYPE_FIELD:
RftUnits::exceptions::field(ret);
break;
case UT::UNIT_TYPE_LAB:
RftUnits::exceptions::lab(ret);
break;
case UT::UNIT_TYPE_PVT_M:
RftUnits::exceptions::pvt_m(ret);
break;
case UT::UNIT_TYPE_INPUT:
RftUnits::exceptions::input(ret);
break;
}
return ret;
}
void writeWellHeader(const double elapsed,
const std::string& wellName,
const ::Opm::Schedule& sched,
const ::Opm::UnitSystem& usys,
::Opm::EclIO::OutputStream::RFT& rftFile)
{
{
const auto time =
usys.from_si(::Opm::UnitSystem::measure::time, elapsed);
rftFile.write("TIME", std::vector<float> {
static_cast<float>(time)
});
}
{
const auto timePoint = ::Opm::RestartIO::
getSimulationTimePoint(sched.getStartTime(), elapsed);
rftFile.write("DATE", std::vector<int> {
timePoint.day, // 1..31
timePoint.month, // 1..12
timePoint.year,
});
}
rftFile.write("WELLETC", wellETC(wellName, "R", usys));
}
void write(const WellRFT& rft, ::Opm::EclIO::OutputStream::RFT& rftFile)
{
rftFile.write("CONIPOS", rft.conI());
rftFile.write("CONJPOS", rft.conJ());
rftFile.write("CONKPOS", rft.conK());
rftFile.write("HOSTGRID", rft.hostgrid());
rftFile.write("DEPTH" , rft.depth());
rftFile.write("PRESSURE", rft.pressure());
rftFile.write("SWAT" , rft.swat());
rftFile.write("SGAS" , rft.sgas());
}
void writeWellRFT(const int reportStep,
const double elapsed,
const std::string& wname,
const ::Opm::UnitSystem& usys,
const ::Opm::EclipseGrid& grid,
const ::Opm::Schedule& sched,
const std::vector<Opm::data::Connection>& xcon,
::Opm::EclIO::OutputStream::RFT& rftFile)
{
const auto rft =
createWellRFT(reportStep, wname, usys, grid, sched, xcon);
writeWellHeader(elapsed, wname, sched, usys, rftFile);
if (rft.nConn() > std::size_t{0}) {
write(rft, rftFile);
}
}
} // Anonymous namespace
void Opm::RftIO::write(const int reportStep,
const double elapsed,
const ::Opm::UnitSystem& usys,
const ::Opm::EclipseGrid& grid,
const ::Opm::Schedule& schedule,
const ::Opm::data::WellRates& wellSol,
::Opm::EclIO::OutputStream::RFT& rftFile)
{
const auto& rftCfg = schedule.rftConfig();
if (! rftCfg.active(reportStep)) {
// RFT not yet activated. Nothing to do.
return;
}
for (const auto& wname : schedule.wellNames(reportStep)) {
if (! (rftCfg.rft(wname, reportStep) ||
rftCfg.plt(wname, reportStep)))
{
// RFT output not requested for 'wname' at this time.
continue;
}
// RFT output requested for 'wname' at this time.
writeWellRFT(reportStep, elapsed, wname, usys, grid,
schedule, wellSol.at(wname).connections,
rftFile);
}
}

File diff suppressed because it is too large Load Diff