mirror of
https://github.com/OPM/ResInsight.git
synced 2026-10-27 23:50:19 -05:00
Magne Sjaastad
d52f595c41
Replaced all references to RiaStatisticsTools with RigStatisticsTools across the codebase, including header includes, static method calls, and test files. Updated CMake configuration to remove RiaStatisticsTools and add RigStatisticsTools. Renamed and updated related unit tests. No changes to logic or algorithms; this is a namespace and organizational refactor for consistency. Also added ApplicationLibCode include directories as PRIVATE to relevant targets.
1251 lines
63 KiB
C++
1251 lines
63 KiB
C++
/////////////////////////////////////////////////////////////////////////////////
|
|
//
|
|
// Copyright (C) 2024- Equinor ASA
|
|
//
|
|
// ResInsight 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.
|
|
//
|
|
// ResInsight 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 at <http://www.gnu.org/licenses/gpl.html>
|
|
// for more details.
|
|
//
|
|
/////////////////////////////////////////////////////////////////////////////////
|
|
|
|
#include "RigWellTargetMapping.h"
|
|
|
|
#include "RiaDefines.h"
|
|
#include "RiaEclipseUnitTools.h"
|
|
#include "RiaLogging.h"
|
|
#include "RiaPorosityModel.h"
|
|
#include "RiaResultNames.h"
|
|
#include "RiaWeightedMeanCalculator.h"
|
|
|
|
#include "RigActiveCellInfo.h"
|
|
#include "RigCaseCellResultsData.h"
|
|
#include "RigEclipseCaseData.h"
|
|
#include "RigEclipseResultAddress.h"
|
|
#include "RigEclipseResultInfo.h"
|
|
#include "RigFloodingSettings.h"
|
|
#include "RigHydrocarbonFlowTools.h"
|
|
#include "RigMainGrid.h"
|
|
#include "RigNNCData.h"
|
|
#include "RigStatisticsMath.h"
|
|
|
|
#include "RimEclipseCase.h"
|
|
#include "RimEclipseCaseEnsemble.h"
|
|
#include "RimEclipseView.h"
|
|
#include "RimProject.h"
|
|
#include "RimPropertyFilterCollection.h"
|
|
#include "RimRegularGridCase.h"
|
|
#include "RimTools.h"
|
|
|
|
#include "cafProgressInfo.h"
|
|
#include "cafVecIjk.h"
|
|
|
|
#include "cvfBoundingBox.h"
|
|
#include "cvfMath.h"
|
|
#include "cvfStructGrid.h"
|
|
|
|
#include <cmath>
|
|
#include <limits>
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
///
|
|
//--------------------------------------------------------------------------------------------------
|
|
void RigWellTargetMapping::generateCandidates( RimEclipseCase* eclipseCase,
|
|
size_t timeStepIdx,
|
|
VolumeType volumeType,
|
|
VolumesType volumesType,
|
|
VolumeResultType volumeResultType,
|
|
const RigFloodingSettings& floodingSettings,
|
|
const ClusteringLimits& limits,
|
|
bool skipUndefinedResults,
|
|
bool setTimeStepInView = true )
|
|
{
|
|
if ( !eclipseCase->ensureReservoirCaseIsOpen() ) return;
|
|
|
|
auto activeCellCount = getActiveCellCount( eclipseCase );
|
|
if ( !activeCellCount )
|
|
{
|
|
RiaLogging::error( "No active cells found" );
|
|
return;
|
|
}
|
|
|
|
auto resultsData = eclipseCase->results( RiaDefines::PorosityModelType::MATRIX_MODEL );
|
|
if ( !resultsData ) return;
|
|
|
|
auto caseData = eclipseCase->eclipseCaseData();
|
|
if ( !caseData ) return;
|
|
|
|
RiaDefines::EclipseUnitSystem unitsType = caseData->unitsType();
|
|
|
|
DataContainer data;
|
|
data.volume = getVolumeVector( *resultsData, unitsType, volumeType, volumesType, volumeResultType, timeStepIdx, floodingSettings );
|
|
if ( data.volume.empty() )
|
|
{
|
|
RiaLogging::error( "Unable to produce volume vector." );
|
|
return;
|
|
}
|
|
|
|
RigEclipseResultAddress soilAddress( RiaDefines::ResultCatType::DYNAMIC_NATIVE, RiaResultNames::soil() );
|
|
if ( resultsData->ensureKnownResultLoaded( soilAddress ) )
|
|
{
|
|
data.saturationOil = resultsData->cellScalarResults( soilAddress, timeStepIdx );
|
|
}
|
|
|
|
RigEclipseResultAddress sgasAddress( RiaDefines::ResultCatType::DYNAMIC_NATIVE, RiaResultNames::sgas() );
|
|
if ( resultsData->ensureKnownResultLoaded( sgasAddress ) )
|
|
{
|
|
data.saturationGas = resultsData->cellScalarResults( sgasAddress, timeStepIdx );
|
|
}
|
|
else
|
|
{
|
|
if ( volumeType == VolumeType::GAS || volumeType == VolumeType::HYDROCARBON )
|
|
{
|
|
RiaLogging::error( "Missing SGAS result needed for well target mapping volume type: " +
|
|
caf::AppEnum<VolumeType>::uiText( volumeType ) );
|
|
return;
|
|
}
|
|
}
|
|
|
|
RigEclipseResultAddress pressureAddress( RiaDefines::ResultCatType::DYNAMIC_NATIVE, "PRESSURE" );
|
|
resultsData->ensureKnownResultLoaded( pressureAddress );
|
|
data.pressure = resultsData->cellScalarResults( pressureAddress, timeStepIdx );
|
|
|
|
RigEclipseResultAddress permeabilityXAddress( RiaDefines::ResultCatType::STATIC_NATIVE, "PERMX" );
|
|
resultsData->ensureKnownResultLoaded( permeabilityXAddress );
|
|
data.permeabilityX = resultsData->cellScalarResults( permeabilityXAddress, 0 );
|
|
|
|
RigEclipseResultAddress transmissibilityXAddress( RiaDefines::ResultCatType::STATIC_NATIVE, "TRANX" );
|
|
resultsData->ensureKnownResultLoaded( transmissibilityXAddress );
|
|
data.transmissibilityX = resultsData->cellScalarResults( transmissibilityXAddress, 0 );
|
|
|
|
RigEclipseResultAddress transmissibilityYAddress( RiaDefines::ResultCatType::STATIC_NATIVE, "TRANY" );
|
|
resultsData->ensureKnownResultLoaded( transmissibilityYAddress );
|
|
data.transmissibilityY = resultsData->cellScalarResults( transmissibilityYAddress, 0 );
|
|
|
|
RigEclipseResultAddress transmissibilityZAddress( RiaDefines::ResultCatType::STATIC_NATIVE, "TRANZ" );
|
|
resultsData->ensureKnownResultLoaded( transmissibilityZAddress );
|
|
data.transmissibilityZ = resultsData->cellScalarResults( transmissibilityZAddress, 0 );
|
|
|
|
auto mainGrid = eclipseCase->eclipseCaseData()->mainGrid();
|
|
if ( mainGrid != nullptr )
|
|
{
|
|
data.transmissibilityNNC = mainGrid->nncData()->staticConnectionScalarResultByName( RiaDefines::propertyNameCombTrans() );
|
|
}
|
|
else
|
|
{
|
|
data.transmissibilityNNC = nullptr;
|
|
}
|
|
|
|
std::vector<double> filterVector;
|
|
if ( !limits.filter.empty() )
|
|
{
|
|
filterVector = limits.filter;
|
|
}
|
|
else if ( limits.filterAddress.isValid() )
|
|
{
|
|
resultsData->ensureKnownResultLoaded( limits.filterAddress );
|
|
filterVector = resultsData->cellScalarResults( limits.filterAddress, timeStepIdx );
|
|
}
|
|
|
|
// Disable filtering if not configured: set all items to 1 (which includes all cells).
|
|
if ( filterVector.empty() || filterVector.size() != data.pressure.size() )
|
|
{
|
|
filterVector = std::vector<double>( data.pressure.size(), 1.0 );
|
|
}
|
|
|
|
const QString logKeyword = "RigWellTargetMapping";
|
|
|
|
std::vector<int> clusters( activeCellCount.value(), 0 );
|
|
auto start = std::chrono::high_resolution_clock::now();
|
|
int numClusters = limits.maxNumTargets;
|
|
int maxIterations = limits.maxIterations;
|
|
int numClustersFound = 0;
|
|
for ( int clusterId = 1; clusterId <= numClusters; clusterId++ )
|
|
{
|
|
std::optional<caf::VecIjk0> startCell = findStartCell( eclipseCase, timeStepIdx, volumeType, limits, data, filterVector, clusters );
|
|
|
|
if ( startCell.has_value() )
|
|
{
|
|
RiaLogging::info( QString( "Cluster %1 start cell: [%2 %3 %4] " )
|
|
.arg( clusterId )
|
|
.arg( startCell->i() + 1 )
|
|
.arg( startCell->j() + 1 )
|
|
.arg( startCell->k() + 1 ),
|
|
logKeyword );
|
|
|
|
growCluster( eclipseCase, startCell.value(), volumeType, limits, data, filterVector, clusters, clusterId, timeStepIdx, maxIterations );
|
|
numClustersFound++;
|
|
}
|
|
else
|
|
{
|
|
RiaLogging::error( "No suitable starting cell found" );
|
|
break;
|
|
}
|
|
}
|
|
|
|
RiaLogging::info( QString( "Found %1 clusters." ).arg( numClustersFound ) );
|
|
|
|
auto finish = std::chrono::high_resolution_clock::now();
|
|
|
|
auto milliseconds = std::chrono::duration_cast<std::chrono::milliseconds>( finish - start );
|
|
RiaLogging::info( QString( "Time spent: %1 ms" ).arg( milliseconds.count() ), logKeyword );
|
|
|
|
QString resultName = RigWellTargetMapping::wellTargetResultName();
|
|
createResultVector( *eclipseCase, resultName, clusters, timeStepIdx );
|
|
|
|
std::vector<ClusterStatistics> statistics =
|
|
generateStatistics( eclipseCase, data.pressure, data.permeabilityX, numClustersFound, timeStepIdx, resultName );
|
|
std::vector<double> totalPorvSoil( clusters.size(), std::numeric_limits<double>::infinity() );
|
|
std::vector<double> totalPorvSgas( clusters.size(), std::numeric_limits<double>::infinity() );
|
|
std::vector<double> totalPorvSoilAndSgas( clusters.size(), std::numeric_limits<double>::infinity() );
|
|
std::vector<double> totalFipOil( clusters.size(), std::numeric_limits<double>::infinity() );
|
|
std::vector<double> totalFipGas( clusters.size(), std::numeric_limits<double>::infinity() );
|
|
std::vector<double> totalRfipOil( clusters.size(), std::numeric_limits<double>::infinity() );
|
|
std::vector<double> totalRfipGas( clusters.size(), std::numeric_limits<double>::infinity() );
|
|
std::vector<double> totalSfipOil( clusters.size(), std::numeric_limits<double>::infinity() );
|
|
std::vector<double> totalSfipGas( clusters.size(), std::numeric_limits<double>::infinity() );
|
|
|
|
auto addValuesForClusterId = []( std::vector<double>& values, const std::vector<int>& clusters, int clusterId, double value )
|
|
{
|
|
#pragma omp parallel for
|
|
for ( int i = 0; i < static_cast<int>( clusters.size() ); i++ )
|
|
{
|
|
if ( clusters[i] == clusterId ) values[i] = value;
|
|
}
|
|
};
|
|
|
|
int clusterId = 1;
|
|
for ( const auto& s : statistics )
|
|
{
|
|
auto logIfValid = [&logKeyword]( const QString& pattern, double value )
|
|
{
|
|
if ( !std::isinf( value ) && !std::isnan( value ) ) RiaLogging::info( pattern.arg( value ), logKeyword );
|
|
};
|
|
|
|
RiaLogging::info( QString( "Cluster #%1 Statistics" ).arg( s.id ), logKeyword );
|
|
RiaLogging::info( QString( "Number of cells: %1" ).arg( s.numCells ), logKeyword );
|
|
logIfValid( "Total PORV*SOIL: %1", s.totalPorvSoil );
|
|
logIfValid( "Total PORV*SOIL: %1", s.totalPorvSoil );
|
|
logIfValid( "Total PORV*SGAS: %1", s.totalPorvSgas );
|
|
logIfValid( "Total PORV*(SOIL+SGAS): %1", s.totalPorvSoilAndSgas );
|
|
logIfValid( "Total FIPOIL: %1", s.totalFipOil );
|
|
logIfValid( "Total FIPGAS: %1", s.totalFipGas );
|
|
logIfValid( "Total RFIPOIL: %1", s.totalRfipOil );
|
|
logIfValid( "Total RFIPGAS: %1", s.totalRfipGas );
|
|
logIfValid( "Total SFIPOIL: %1", s.totalSfipOil );
|
|
logIfValid( "Total SFIPGAS: %1", s.totalSfipGas );
|
|
logIfValid( "Average Permeability: %1", s.permeability );
|
|
logIfValid( "Average Pressure: %1", s.pressure );
|
|
|
|
addValuesForClusterId( totalPorvSoil, clusters, clusterId, s.totalPorvSoil );
|
|
addValuesForClusterId( totalPorvSgas, clusters, clusterId, s.totalPorvSgas );
|
|
addValuesForClusterId( totalPorvSoilAndSgas, clusters, clusterId, s.totalPorvSoilAndSgas );
|
|
addValuesForClusterId( totalFipOil, clusters, clusterId, s.totalFipOil );
|
|
addValuesForClusterId( totalFipGas, clusters, clusterId, s.totalFipGas );
|
|
addValuesForClusterId( totalRfipOil, clusters, clusterId, s.totalRfipOil );
|
|
addValuesForClusterId( totalRfipGas, clusters, clusterId, s.totalRfipGas );
|
|
addValuesForClusterId( totalSfipOil, clusters, clusterId, s.totalSfipOil );
|
|
addValuesForClusterId( totalSfipGas, clusters, clusterId, s.totalSfipGas );
|
|
|
|
clusterId++;
|
|
}
|
|
|
|
if ( skipUndefinedResults )
|
|
{
|
|
const int ts = (int)timeStepIdx;
|
|
createResultVectorIfDefined( *eclipseCase, "TOTAL_PORV_SOIL", totalPorvSoil, ts );
|
|
createResultVectorIfDefined( *eclipseCase, "TOTAL_PORV_SGAS", totalPorvSgas, ts );
|
|
createResultVectorIfDefined( *eclipseCase, "TOTAL_PORV_SOIL_SGAS", totalPorvSoilAndSgas, ts );
|
|
createResultVectorIfDefined( *eclipseCase, "TOTAL_FIPOIL", totalFipOil, ts );
|
|
createResultVectorIfDefined( *eclipseCase, "TOTAL_FIPGAS", totalFipGas, ts );
|
|
createResultVectorIfDefined( *eclipseCase, "TOTAL_RFIPOIL", totalRfipOil, ts );
|
|
createResultVectorIfDefined( *eclipseCase, "TOTAL_RFIPGAS", totalRfipGas, ts );
|
|
createResultVectorIfDefined( *eclipseCase, "TOTAL_SFIPOIL", totalSfipOil, ts );
|
|
createResultVectorIfDefined( *eclipseCase, "TOTAL_SFIPGAS", totalSfipGas, ts );
|
|
}
|
|
else
|
|
{
|
|
createResultVector( *eclipseCase, "TOTAL_PORV_SOIL", totalPorvSoil, timeStepIdx );
|
|
createResultVector( *eclipseCase, "TOTAL_PORV_SGAS", totalPorvSgas, timeStepIdx );
|
|
createResultVector( *eclipseCase, "TOTAL_PORV_SOIL_SGAS", totalPorvSoilAndSgas, timeStepIdx );
|
|
createResultVector( *eclipseCase, "TOTAL_FIPOIL", totalFipOil, timeStepIdx );
|
|
createResultVector( *eclipseCase, "TOTAL_FIPGAS", totalFipGas, timeStepIdx );
|
|
createResultVector( *eclipseCase, "TOTAL_RFIPOIL", totalRfipOil, timeStepIdx );
|
|
createResultVector( *eclipseCase, "TOTAL_RFIPGAS", totalRfipGas, timeStepIdx );
|
|
createResultVector( *eclipseCase, "TOTAL_SFIPOIL", totalSfipOil, timeStepIdx );
|
|
createResultVector( *eclipseCase, "TOTAL_SFIPGAS", totalSfipGas, timeStepIdx );
|
|
}
|
|
eclipseCase->updateResultAddressCollection();
|
|
|
|
// Update views and property filters
|
|
RimProject* proj = RimProject::current();
|
|
proj->scheduleCreateDisplayModelAndRedrawAllViews();
|
|
for ( auto view : eclipseCase->reservoirViews() )
|
|
{
|
|
if ( auto eclipseView = dynamic_cast<RimEclipseView*>( view ) )
|
|
{
|
|
if ( setTimeStepInView )
|
|
{
|
|
eclipseView->setCurrentTimeStep( (int)timeStepIdx );
|
|
}
|
|
eclipseView->scheduleReservoirGridGeometryRegen();
|
|
eclipseView->propertyFilterCollection()->updateConnectedEditors();
|
|
}
|
|
}
|
|
}
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
///
|
|
//--------------------------------------------------------------------------------------------------
|
|
std::optional<caf::VecIjk0> RigWellTargetMapping::findStartCell( RimEclipseCase* eclipseCase,
|
|
size_t timeStepIdx,
|
|
const VolumeType volumeType,
|
|
const ClusteringLimits& limits,
|
|
const DataContainer& data,
|
|
const std::vector<double>& filterVector,
|
|
const std::vector<int>& clusters )
|
|
{
|
|
auto resultsData = eclipseCase->results( RiaDefines::PorosityModelType::MATRIX_MODEL );
|
|
if ( !resultsData )
|
|
{
|
|
RiaLogging::error( "No results data found for eclipse case" );
|
|
return {};
|
|
}
|
|
|
|
size_t startCell = std::numeric_limits<size_t>::max();
|
|
double maxVolume = -std::numeric_limits<double>::max();
|
|
const size_t numReservoirCells = resultsData->activeCellInfo()->reservoirCellCount();
|
|
for ( size_t reservoirCellIdx = 0; reservoirCellIdx < numReservoirCells; reservoirCellIdx++ )
|
|
{
|
|
size_t resultIndex = resultsData->activeCellInfo()->cellResultIndex( reservoirCellIdx );
|
|
if ( resultIndex != cvf::UNDEFINED_SIZE_T && clusters[resultIndex] == 0 )
|
|
{
|
|
const double cellVolume = data.volume[resultIndex];
|
|
const double cellPressure = data.pressure[resultIndex];
|
|
|
|
const bool isSaturationValid = isSaturationSufficient( volumeType, data, limits, resultIndex );
|
|
|
|
const double cellPermeabiltyX = data.permeabilityX[resultIndex];
|
|
const bool permeabilityValid = ( cellPermeabiltyX >= limits.permeability );
|
|
|
|
const bool filterValue = !std::isinf( filterVector[resultIndex] ) && filterVector[resultIndex] > 0.0;
|
|
|
|
if ( cellVolume > maxVolume && cellPressure >= limits.pressure && permeabilityValid && filterValue && isSaturationValid )
|
|
{
|
|
maxVolume = cellVolume;
|
|
startCell = reservoirCellIdx;
|
|
}
|
|
}
|
|
}
|
|
|
|
if ( startCell == std::numeric_limits<size_t>::max() ) return {};
|
|
|
|
return eclipseCase->mainGrid()->ijkFromCellIndex( startCell );
|
|
}
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
///
|
|
//--------------------------------------------------------------------------------------------------
|
|
void RigWellTargetMapping::growCluster( RimEclipseCase* eclipseCase,
|
|
const caf::VecIjk0& startCell,
|
|
const VolumeType volumeType,
|
|
const ClusteringLimits& limits,
|
|
const DataContainer& data,
|
|
const std::vector<double>& filterVector,
|
|
std::vector<int>& clusters,
|
|
int clusterId,
|
|
size_t timeStepIdx,
|
|
int maxIterations )
|
|
{
|
|
auto resultsData = eclipseCase->results( RiaDefines::PorosityModelType::MATRIX_MODEL );
|
|
|
|
// Initially only the start cell is found
|
|
size_t reservoirCellIdx = eclipseCase->mainGrid()->cellIndexFromIJK( startCell.i(), startCell.j(), startCell.k() );
|
|
std::vector<size_t> foundCells = { reservoirCellIdx };
|
|
assignClusterIdToCells( *resultsData->activeCellInfo(), foundCells, clusters, clusterId );
|
|
|
|
for ( int i = 0; i < maxIterations; i++ )
|
|
{
|
|
foundCells = findCandidates( eclipseCase, foundCells, volumeType, limits, data, filterVector, clusters );
|
|
if ( foundCells.empty() ) break;
|
|
assignClusterIdToCells( *resultsData->activeCellInfo(), foundCells, clusters, clusterId );
|
|
}
|
|
}
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
///
|
|
//--------------------------------------------------------------------------------------------------
|
|
std::vector<size_t> RigWellTargetMapping::findCandidates( RimEclipseCase* eclipseCase,
|
|
const std::vector<size_t>& previousCells,
|
|
const VolumeType volumeType,
|
|
const ClusteringLimits& limits,
|
|
const DataContainer& data,
|
|
const std::vector<double>& filterVector,
|
|
std::vector<int>& clusters )
|
|
{
|
|
std::vector<size_t> candidates;
|
|
auto resultsData = eclipseCase->results( RiaDefines::PorosityModelType::MATRIX_MODEL );
|
|
auto mainGrid = eclipseCase->eclipseCaseData()->mainGrid();
|
|
|
|
const std::vector<CellFaceType> faces = {
|
|
cvf::StructGridInterface::FaceType::POS_I,
|
|
cvf::StructGridInterface::FaceType::NEG_I,
|
|
cvf::StructGridInterface::FaceType::POS_J,
|
|
cvf::StructGridInterface::FaceType::NEG_J,
|
|
cvf::StructGridInterface::FaceType::POS_K,
|
|
cvf::StructGridInterface::FaceType::NEG_K,
|
|
};
|
|
|
|
for ( size_t cellIdx : previousCells )
|
|
{
|
|
const RigCell& cell = mainGrid->cell( cellIdx );
|
|
if ( cell.isInvalid() ) continue;
|
|
|
|
RigGridBase* grid = cell.hostGrid();
|
|
size_t gridLocalCellIndex = cell.gridLocalCellIndex();
|
|
size_t resultIndex = resultsData->activeCellInfo()->cellResultIndex( cellIdx );
|
|
|
|
size_t i, j, k;
|
|
|
|
grid->ijkFromCellIndex( gridLocalCellIndex, &i, &j, &k );
|
|
|
|
for ( CellFaceType face : faces )
|
|
{
|
|
size_t gridLocalNeighborCellIdx;
|
|
if ( grid->cellIJKNeighbor( i, j, k, face, &gridLocalNeighborCellIdx ) )
|
|
{
|
|
size_t neighborResvCellIdx = grid->reservoirCellIndex( gridLocalNeighborCellIdx );
|
|
size_t neighborResultIndex = resultsData->activeCellInfo()->cellResultIndex( neighborResvCellIdx );
|
|
if ( neighborResultIndex != cvf::UNDEFINED_SIZE_T && clusters[neighborResultIndex] == 0 )
|
|
{
|
|
double permeability = data.permeabilityX[neighborResultIndex];
|
|
double transmissibility = getTransmissibilityValueForFace( data.transmissibilityX,
|
|
data.transmissibilityY,
|
|
data.transmissibilityZ,
|
|
face,
|
|
resultIndex,
|
|
neighborResultIndex );
|
|
bool filterValue = !std::isinf( filterVector[neighborResultIndex] ) && filterVector[neighborResultIndex] > 0.0;
|
|
|
|
const bool isSaturationValid = isSaturationSufficient( volumeType, data, limits, neighborResultIndex );
|
|
|
|
if ( data.pressure[neighborResultIndex] > limits.pressure && permeability > limits.permeability &&
|
|
transmissibility > limits.transmissibility && filterValue && isSaturationValid )
|
|
{
|
|
candidates.push_back( neighborResvCellIdx );
|
|
clusters[neighborResultIndex] = -1;
|
|
}
|
|
}
|
|
}
|
|
}
|
|
|
|
if ( data.transmissibilityNNC != nullptr )
|
|
{
|
|
auto nncCells = nncConnectionCellAndResult( cellIdx, mainGrid );
|
|
for ( auto& [cellInfo, nncResultIdx] : nncCells )
|
|
{
|
|
auto& [otherCellIdx, face] = cellInfo;
|
|
double transmissibility = data.transmissibilityNNC->at( nncResultIdx );
|
|
|
|
size_t otherResultIndex = resultsData->activeCellInfo()->cellResultIndex( otherCellIdx );
|
|
|
|
double permeability = data.permeabilityX[otherResultIndex];
|
|
|
|
bool filterValue = !std::isinf( filterVector[otherResultIndex] ) && filterVector[otherResultIndex] > 0.0;
|
|
|
|
const bool isSaturationValid = isSaturationSufficient( volumeType, data, limits, otherResultIndex );
|
|
|
|
if ( data.pressure[otherResultIndex] > limits.pressure && permeability > limits.permeability &&
|
|
transmissibility > limits.transmissibility && filterValue && isSaturationValid )
|
|
{
|
|
candidates.push_back( otherCellIdx );
|
|
clusters[otherResultIndex] = -1;
|
|
}
|
|
}
|
|
}
|
|
}
|
|
|
|
return candidates;
|
|
}
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
///
|
|
//--------------------------------------------------------------------------------------------------
|
|
void RigWellTargetMapping::assignClusterIdToCells( const RigActiveCellInfo& activeCellInfo,
|
|
const std::vector<size_t>& cells,
|
|
std::vector<int>& clusters,
|
|
int clusterId )
|
|
{
|
|
for ( size_t reservoirCellIdx : cells )
|
|
{
|
|
size_t resultIndex = activeCellInfo.cellResultIndex( reservoirCellIdx );
|
|
if ( resultIndex != cvf::UNDEFINED_SIZE_T ) clusters[resultIndex] = clusterId;
|
|
}
|
|
}
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
///
|
|
//--------------------------------------------------------------------------------------------------
|
|
void RigWellTargetMapping::createDynamicResultEntry( RigCaseCellResultsData* resultsData, const RigEclipseResultAddress& address )
|
|
{
|
|
if ( !resultsData->hasResultEntry( address ) )
|
|
{
|
|
resultsData->createResultEntry( address, false );
|
|
|
|
RigEclipseResultAddress addrToMaxTimeStepCountResult;
|
|
resultsData->maxTimeStepCount( &addrToMaxTimeStepCountResult );
|
|
const std::vector<RigEclipseTimeStepInfo> timeStepInfos = resultsData->timeStepInfos( addrToMaxTimeStepCountResult );
|
|
resultsData->setTimeStepInfos( address, timeStepInfos );
|
|
}
|
|
}
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
///
|
|
//--------------------------------------------------------------------------------------------------
|
|
void RigWellTargetMapping::createResultVector( RimEclipseCase& eclipseCase,
|
|
const QString& resultName,
|
|
const std::vector<int>& clusterIds,
|
|
size_t timeStepIdx )
|
|
{
|
|
RigEclipseResultAddress resultAddress( RiaDefines::ResultCatType::GENERATED, RiaDefines::ResultDataType::INTEGER, resultName );
|
|
|
|
auto resultsData = eclipseCase.results( RiaDefines::PorosityModelType::MATRIX_MODEL );
|
|
|
|
createDynamicResultEntry( resultsData, resultAddress );
|
|
|
|
std::vector<double>* resultVector = resultsData->modifiableCellScalarResult( resultAddress, timeStepIdx );
|
|
resultVector->resize( clusterIds.size(), std::numeric_limits<double>::infinity() );
|
|
|
|
std::fill( resultVector->begin(), resultVector->end(), std::numeric_limits<double>::infinity() );
|
|
|
|
for ( size_t idx = 0; idx < clusterIds.size(); idx++ )
|
|
{
|
|
if ( clusterIds[idx] > 0 )
|
|
{
|
|
resultVector->at( idx ) = 1.0 * clusterIds[idx];
|
|
}
|
|
}
|
|
|
|
std::set<int> uniqueClusterIds( clusterIds.begin(), clusterIds.end() );
|
|
uniqueClusterIds.erase( 0 ); // Remove unassigned cluster id
|
|
resultsData->recalculateStatistics( resultAddress, uniqueClusterIds );
|
|
}
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
///
|
|
//--------------------------------------------------------------------------------------------------
|
|
void RigWellTargetMapping::createResultVectorIfDefined( RimEclipseCase& eclipseCase,
|
|
const QString& resultName,
|
|
const std::vector<double>& values,
|
|
int timeStepIdx )
|
|
{
|
|
// Avoid creating the result vector if all values are inf/nan
|
|
if ( std::all_of( values.begin(), values.end(), []( auto v ) { return std::isinf( v ) || std::isnan( v ); } ) ) return;
|
|
|
|
if ( timeStepIdx < 0 )
|
|
{
|
|
createStaticResultVector( eclipseCase, resultName, values );
|
|
}
|
|
else
|
|
{
|
|
createResultVector( eclipseCase, resultName, values, timeStepIdx );
|
|
}
|
|
}
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
///
|
|
//--------------------------------------------------------------------------------------------------
|
|
void RigWellTargetMapping::createResultVector( RimEclipseCase& eclipseCase,
|
|
const QString& resultName,
|
|
const std::vector<double>& values,
|
|
size_t timeStepIdx )
|
|
{
|
|
RigEclipseResultAddress resultAddress( RiaDefines::ResultCatType::GENERATED, resultName );
|
|
|
|
auto resultsData = eclipseCase.results( RiaDefines::PorosityModelType::MATRIX_MODEL );
|
|
|
|
createDynamicResultEntry( resultsData, resultAddress );
|
|
|
|
auto resultVector = resultsData->modifiableCellScalarResult( resultAddress, timeStepIdx );
|
|
resultVector->resize( values.size(), std::numeric_limits<double>::infinity() );
|
|
std::copy( values.begin(), values.end(), resultVector->begin() );
|
|
|
|
resultsData->recalculateStatistics( resultAddress );
|
|
}
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
///
|
|
//--------------------------------------------------------------------------------------------------
|
|
void RigWellTargetMapping::createStaticResultVector( RimEclipseCase& eclipseCase, const QString& resultName, const std::vector<int>& intValues )
|
|
{
|
|
RigEclipseResultAddress resultAddress( RiaDefines::ResultCatType::GENERATED, resultName );
|
|
|
|
auto resultsData = eclipseCase.results( RiaDefines::PorosityModelType::MATRIX_MODEL );
|
|
|
|
resultsData->addStaticScalarResult( RiaDefines::ResultCatType::GENERATED, resultName, false, intValues.size() );
|
|
|
|
std::vector<double>* resultVector = resultsData->modifiableCellScalarResult( resultAddress, 0 );
|
|
resultVector->resize( intValues.size(), std::numeric_limits<double>::infinity() );
|
|
|
|
std::fill( resultVector->begin(), resultVector->end(), std::numeric_limits<double>::infinity() );
|
|
|
|
for ( size_t idx = 0; idx < intValues.size(); idx++ )
|
|
{
|
|
if ( intValues[idx] > 0 )
|
|
{
|
|
resultVector->at( idx ) = intValues[idx];
|
|
}
|
|
}
|
|
|
|
resultsData->recalculateStatistics( resultAddress );
|
|
}
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
///
|
|
//--------------------------------------------------------------------------------------------------
|
|
void RigWellTargetMapping::createStaticResultVector( RimEclipseCase& eclipseCase, const QString& resultName, const std::vector<double>& values )
|
|
{
|
|
RigEclipseResultAddress resultAddress( RiaDefines::ResultCatType::GENERATED, resultName );
|
|
|
|
auto resultsData = eclipseCase.results( RiaDefines::PorosityModelType::MATRIX_MODEL );
|
|
|
|
resultsData->addStaticScalarResult( RiaDefines::ResultCatType::GENERATED, resultName, false, values.size() );
|
|
|
|
std::vector<double>* resultVector = resultsData->modifiableCellScalarResult( resultAddress, 0 );
|
|
resultVector->resize( values.size(), std::numeric_limits<double>::infinity() );
|
|
|
|
std::copy( values.begin(), values.end(), resultVector->begin() );
|
|
|
|
resultsData->recalculateStatistics( resultAddress );
|
|
}
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
///
|
|
//--------------------------------------------------------------------------------------------------
|
|
std::optional<size_t> RigWellTargetMapping::getActiveCellCount( RimEclipseCase* eclipseCase )
|
|
{
|
|
auto resultsData = eclipseCase->results( RiaDefines::PorosityModelType::MATRIX_MODEL );
|
|
if ( !resultsData ) return {};
|
|
|
|
return resultsData->activeCellInfo()->reservoirActiveCellCount();
|
|
}
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
///
|
|
//--------------------------------------------------------------------------------------------------
|
|
double RigWellTargetMapping::getValueForFace( const std::vector<double>& x,
|
|
const std::vector<double>& y,
|
|
const std::vector<double>& z,
|
|
CellFaceType face,
|
|
size_t resultIndex )
|
|
{
|
|
if ( face == cvf::StructGridInterface::FaceType::POS_I || face == cvf::StructGridInterface::FaceType::NEG_I ) return x[resultIndex];
|
|
if ( face == cvf::StructGridInterface::FaceType::POS_J || face == cvf::StructGridInterface::FaceType::NEG_J ) return y[resultIndex];
|
|
if ( face == cvf::StructGridInterface::FaceType::POS_K || face == cvf::StructGridInterface::FaceType::NEG_K ) return z[resultIndex];
|
|
return std::numeric_limits<double>::infinity();
|
|
}
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
///
|
|
//--------------------------------------------------------------------------------------------------
|
|
double RigWellTargetMapping::getTransmissibilityValueForFace( const std::vector<double>& x,
|
|
const std::vector<double>& y,
|
|
const std::vector<double>& z,
|
|
CellFaceType face,
|
|
size_t resultIndex,
|
|
size_t neighborResultIndex )
|
|
{
|
|
// For negative directions (NEG_I, NEG_J, NEG_K) use the value from the neighbor cell
|
|
bool isPos = face == cvf::StructGridInterface::FaceType::POS_I || face == cvf::StructGridInterface::FaceType::POS_J ||
|
|
face == cvf::StructGridInterface::FaceType::POS_K;
|
|
size_t index = isPos ? resultIndex : neighborResultIndex;
|
|
return getValueForFace( x, y, z, face, index );
|
|
}
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
///
|
|
//--------------------------------------------------------------------------------------------------
|
|
std::vector<double> RigWellTargetMapping::loadVectorByName( RigCaseCellResultsData& resultsData, const QString& resultName, size_t timeStepIdx )
|
|
{
|
|
RigEclipseResultAddress address( RiaDefines::ResultCatType::DYNAMIC_NATIVE, resultName );
|
|
if ( !resultsData.ensureKnownResultLoaded( address ) ) return {};
|
|
return resultsData.cellScalarResults( address, timeStepIdx );
|
|
}
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
///
|
|
//--------------------------------------------------------------------------------------------------
|
|
QString RigWellTargetMapping::getOilVectorName( VolumesType volumesType )
|
|
{
|
|
switch ( volumesType )
|
|
{
|
|
case VolumesType::RESERVOIR_VOLUMES_COMPUTED:
|
|
return RiaResultNames::riPorvSoil();
|
|
case VolumesType::RESERVOIR_VOLUMES:
|
|
return "RFIPOIL";
|
|
case VolumesType::SURFACE_VOLUMES_SFIP:
|
|
return "SFIPOIL";
|
|
case VolumesType::SURFACE_VOLUMES_FIP:
|
|
return "FIPOIL";
|
|
default:
|
|
{
|
|
CAF_ASSERT( false );
|
|
return "";
|
|
}
|
|
}
|
|
}
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
///
|
|
//--------------------------------------------------------------------------------------------------
|
|
QString RigWellTargetMapping::getGasVectorName( VolumesType volumesType )
|
|
{
|
|
switch ( volumesType )
|
|
{
|
|
case VolumesType::RESERVOIR_VOLUMES_COMPUTED:
|
|
return RiaResultNames::riPorvSgas();
|
|
case VolumesType::RESERVOIR_VOLUMES:
|
|
return "RFIPGAS";
|
|
case VolumesType::SURFACE_VOLUMES_SFIP:
|
|
return "SFIPGAS";
|
|
case VolumesType::SURFACE_VOLUMES_FIP:
|
|
return "FIPGAS";
|
|
default:
|
|
{
|
|
CAF_ASSERT( false );
|
|
return "";
|
|
}
|
|
}
|
|
}
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
///
|
|
//--------------------------------------------------------------------------------------------------
|
|
std::vector<double> RigWellTargetMapping::loadOilVectorByName( RigCaseCellResultsData& resultsData,
|
|
VolumesType volumesType,
|
|
VolumeResultType volumeResultType,
|
|
size_t timeStepIdx,
|
|
const RigFloodingSettings& floodingSettings )
|
|
{
|
|
std::vector<double> volume = loadVectorByName( resultsData, getOilVectorName( volumesType ), timeStepIdx );
|
|
if ( volumeResultType == RigWellTargetMapping::VolumeResultType::MOBILE )
|
|
{
|
|
std::vector<double> residualOil = RigHydrocarbonFlowTools::residualOilData( resultsData,
|
|
RigHydrocarbonFlowTools::ResultType::MOBILE_OIL,
|
|
floodingSettings,
|
|
volume.size() );
|
|
if ( volumesType == RigWellTargetMapping::VolumesType::RESERVOIR_VOLUMES_COMPUTED )
|
|
{
|
|
const std::vector<double>& porvResults =
|
|
resultsData.cellScalarResults( RigEclipseResultAddress( RiaDefines::ResultCatType::STATIC_NATIVE, RiaResultNames::porv() ), 0 );
|
|
|
|
for ( size_t i = 0; i < volume.size(); i++ )
|
|
{
|
|
volume[i] = std::max( volume[i] - ( porvResults[i] * residualOil[i] ), 0.0 );
|
|
}
|
|
}
|
|
else
|
|
{
|
|
const std::vector<double>& soilResults =
|
|
resultsData.cellScalarResults( RigEclipseResultAddress( RiaDefines::ResultCatType::DYNAMIC_NATIVE, RiaResultNames::soil() ),
|
|
timeStepIdx );
|
|
|
|
for ( size_t i = 0; i < volume.size(); i++ )
|
|
{
|
|
if ( soilResults[i] != 0.0 )
|
|
{
|
|
volume[i] = std::max( volume[i] * ( soilResults[i] - residualOil[i] ) / soilResults[i], 0.0 );
|
|
}
|
|
else
|
|
{
|
|
volume[i] = 0.0;
|
|
}
|
|
}
|
|
}
|
|
}
|
|
|
|
return volume;
|
|
}
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
///
|
|
//--------------------------------------------------------------------------------------------------
|
|
std::vector<double> RigWellTargetMapping::loadGasVectorByName( RigCaseCellResultsData& resultsData,
|
|
RiaDefines::EclipseUnitSystem unitsType,
|
|
VolumesType volumesType,
|
|
VolumeResultType volumeResultType,
|
|
size_t timeStepIdx,
|
|
const RigFloodingSettings& floodingSettings )
|
|
{
|
|
std::vector<double> volume = loadVectorByName( resultsData, getGasVectorName( volumesType ), timeStepIdx );
|
|
|
|
if ( volumeResultType == RigWellTargetMapping::VolumeResultType::MOBILE )
|
|
{
|
|
std::vector<double> residualGas = RigHydrocarbonFlowTools::residualGasData( resultsData,
|
|
RigHydrocarbonFlowTools::ResultType::MOBILE_GAS,
|
|
floodingSettings,
|
|
volume.size() );
|
|
if ( volumesType == RigWellTargetMapping::VolumesType::RESERVOIR_VOLUMES_COMPUTED )
|
|
{
|
|
const std::vector<double>& porvResults =
|
|
resultsData.cellScalarResults( RigEclipseResultAddress( RiaDefines::ResultCatType::STATIC_NATIVE, RiaResultNames::porv() ), 0 );
|
|
|
|
for ( size_t i = 0; i < volume.size(); i++ )
|
|
{
|
|
volume[i] = std::max( volume[i] - ( porvResults[i] * residualGas[i] ), 0.0 );
|
|
}
|
|
}
|
|
else
|
|
{
|
|
const std::vector<double>& sgasResults =
|
|
resultsData.cellScalarResults( RigEclipseResultAddress( RiaDefines::ResultCatType::DYNAMIC_NATIVE, RiaResultNames::sgas() ),
|
|
timeStepIdx );
|
|
|
|
for ( size_t i = 0; i < volume.size(); i++ )
|
|
{
|
|
if ( sgasResults[i] != 0.0 )
|
|
{
|
|
volume[i] = std::max( volume[i] * ( sgasResults[i] - residualGas[i] ) / sgasResults[i], 0.0 );
|
|
}
|
|
else
|
|
{
|
|
volume[i] = 0.0;
|
|
}
|
|
}
|
|
}
|
|
}
|
|
|
|
// Convert to gas volumes to oil equivalents
|
|
for ( size_t i = 0; i < volume.size(); i++ )
|
|
{
|
|
volume[i] = RiaEclipseUnitTools::convertSurfaceGasFlowRateToOilEquivalents( unitsType, volume[i] );
|
|
}
|
|
|
|
return volume;
|
|
}
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
///
|
|
//--------------------------------------------------------------------------------------------------
|
|
std::vector<double> RigWellTargetMapping::getVolumeVector( RigCaseCellResultsData& resultsData,
|
|
RiaDefines::EclipseUnitSystem unitsType,
|
|
VolumeType volumeType,
|
|
VolumesType volumesType,
|
|
VolumeResultType volumeResultType,
|
|
size_t timeStepIdx,
|
|
const RigFloodingSettings& floodingSettings )
|
|
{
|
|
if ( volumeType == VolumeType::OIL )
|
|
{
|
|
return loadOilVectorByName( resultsData, volumesType, volumeResultType, timeStepIdx, floodingSettings );
|
|
}
|
|
else if ( volumeType == VolumeType::GAS )
|
|
{
|
|
return loadGasVectorByName( resultsData, unitsType, volumesType, volumeResultType, timeStepIdx, floodingSettings );
|
|
}
|
|
else if ( volumeType == VolumeType::HYDROCARBON )
|
|
{
|
|
std::vector<double> oilVolume = loadOilVectorByName( resultsData, volumesType, volumeResultType, timeStepIdx, floodingSettings );
|
|
std::vector<double> gasVolume =
|
|
loadGasVectorByName( resultsData, unitsType, volumesType, volumeResultType, timeStepIdx, floodingSettings );
|
|
if ( oilVolume.empty() || gasVolume.empty() || oilVolume.size() != gasVolume.size() ) return {};
|
|
|
|
std::vector<double> volume;
|
|
volume.resize( oilVolume.size(), std::numeric_limits<double>::infinity() );
|
|
for ( size_t i = 0; i < oilVolume.size(); i++ )
|
|
{
|
|
volume[i] = oilVolume[i] + gasVolume[i];
|
|
}
|
|
return volume;
|
|
}
|
|
|
|
CAF_ASSERT( false && "Unknown volume type" );
|
|
return {};
|
|
}
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
///
|
|
//--------------------------------------------------------------------------------------------------
|
|
std::vector<RigWellTargetMapping::ClusterStatistics> RigWellTargetMapping::generateStatistics( RimEclipseCase* eclipseCase,
|
|
const std::vector<double>& pressure,
|
|
const std::vector<double>& permeabilityX,
|
|
int numClustersFound,
|
|
size_t timeStepIdx,
|
|
const QString& clusterResultName )
|
|
{
|
|
std::vector<ClusterStatistics> statistics( numClustersFound );
|
|
|
|
auto resultsData = eclipseCase->results( RiaDefines::PorosityModelType::MATRIX_MODEL );
|
|
if ( !resultsData ) return statistics;
|
|
|
|
auto loadData = []( RigCaseCellResultsData* resultsData, RiaDefines::ResultCatType categoryType, const QString& name, size_t timeStepIdx )
|
|
{
|
|
RigEclipseResultAddress address( categoryType, name );
|
|
std::vector<double> values;
|
|
if ( resultsData->ensureKnownResultLoaded( address ) )
|
|
{
|
|
values = resultsData->cellScalarResults( address, timeStepIdx );
|
|
}
|
|
|
|
return values;
|
|
};
|
|
|
|
const std::vector<double> porv = loadData( resultsData, RiaDefines::ResultCatType::STATIC_NATIVE, RiaResultNames::porv(), 0 );
|
|
const std::vector<double> porvSoil =
|
|
loadData( resultsData, RiaDefines::ResultCatType::DYNAMIC_NATIVE, RiaResultNames::riPorvSoil(), timeStepIdx );
|
|
const std::vector<double> porvSgas =
|
|
loadData( resultsData, RiaDefines::ResultCatType::DYNAMIC_NATIVE, RiaResultNames::riPorvSgas(), timeStepIdx );
|
|
const std::vector<double> porvSoilAndSgas =
|
|
loadData( resultsData, RiaDefines::ResultCatType::DYNAMIC_NATIVE, RiaResultNames::riPorvSoilSgas(), timeStepIdx );
|
|
const std::vector<double> fipOil = loadData( resultsData, RiaDefines::ResultCatType::DYNAMIC_NATIVE, "FIPOIL", timeStepIdx );
|
|
const std::vector<double> fipGas = loadData( resultsData, RiaDefines::ResultCatType::DYNAMIC_NATIVE, "FIPGAS", timeStepIdx );
|
|
const std::vector<double> sfipOil = loadData( resultsData, RiaDefines::ResultCatType::DYNAMIC_NATIVE, "SFIPOIL", timeStepIdx );
|
|
const std::vector<double> sfipGas = loadData( resultsData, RiaDefines::ResultCatType::DYNAMIC_NATIVE, "SFIPGAS", timeStepIdx );
|
|
const std::vector<double> rfipOil = loadData( resultsData, RiaDefines::ResultCatType::DYNAMIC_NATIVE, "RFIPOIL", timeStepIdx );
|
|
const std::vector<double> rfipGas = loadData( resultsData, RiaDefines::ResultCatType::DYNAMIC_NATIVE, "RFIPGAS", timeStepIdx );
|
|
|
|
RigEclipseResultAddress clusterAddress( RiaDefines::ResultCatType::GENERATED, clusterResultName );
|
|
resultsData->ensureKnownResultLoaded( clusterAddress );
|
|
const std::vector<double>& clusterIds = resultsData->cellScalarResults( clusterAddress, timeStepIdx );
|
|
|
|
std::vector<RiaWeightedMeanCalculator<double>> permeabilityCalculators( numClustersFound );
|
|
std::vector<RiaWeightedMeanCalculator<double>> pressureCalculators( numClustersFound );
|
|
|
|
for ( size_t idx = 0; idx < clusterIds.size(); idx++ )
|
|
{
|
|
if ( !std::isinf( clusterIds[idx] ) && static_cast<int>( clusterIds[idx] ) > 0 )
|
|
{
|
|
size_t i = clusterIds[idx] - 1;
|
|
if ( i < static_cast<size_t>( numClustersFound ) )
|
|
{
|
|
statistics[i].id = clusterIds[idx];
|
|
statistics[i].numCells++;
|
|
if ( idx < porvSoil.size() ) statistics[i].totalPorvSoil += porvSoil[idx];
|
|
if ( idx < porvSgas.size() ) statistics[i].totalPorvSgas += porvSgas[idx];
|
|
if ( idx < porvSoilAndSgas.size() ) statistics[i].totalPorvSoilAndSgas += porvSoilAndSgas[idx];
|
|
|
|
if ( idx < fipOil.size() ) statistics[i].totalFipOil += fipOil[idx];
|
|
if ( idx < fipGas.size() ) statistics[i].totalFipGas += fipGas[idx];
|
|
|
|
if ( idx < rfipOil.size() ) statistics[i].totalRfipOil += rfipOil[idx];
|
|
if ( idx < rfipGas.size() ) statistics[i].totalRfipGas += rfipGas[idx];
|
|
|
|
if ( idx < sfipOil.size() ) statistics[i].totalSfipOil += sfipOil[idx];
|
|
if ( idx < sfipGas.size() ) statistics[i].totalSfipGas += sfipGas[idx];
|
|
|
|
permeabilityCalculators[i].addValueAndWeight( permeabilityX[idx], porv[idx] );
|
|
|
|
pressureCalculators[i].addValueAndWeight( pressure[idx], porv[idx] );
|
|
}
|
|
}
|
|
}
|
|
|
|
for ( int i = 0; i < numClustersFound; i++ )
|
|
{
|
|
statistics[i].permeability = permeabilityCalculators[i].weightedMean();
|
|
statistics[i].pressure = pressureCalculators[i].weightedMean();
|
|
|
|
// Invalidate results for empty vectors
|
|
if ( porvSoil.empty() ) statistics[i].totalPorvSoil = std::numeric_limits<double>::infinity();
|
|
if ( porvSgas.empty() ) statistics[i].totalPorvSgas = std::numeric_limits<double>::infinity();
|
|
if ( porvSoilAndSgas.empty() ) statistics[i].totalPorvSoilAndSgas = std::numeric_limits<double>::infinity();
|
|
if ( fipOil.empty() ) statistics[i].totalFipOil = std::numeric_limits<double>::infinity();
|
|
if ( fipGas.empty() ) statistics[i].totalFipGas = std::numeric_limits<double>::infinity();
|
|
if ( rfipOil.empty() ) statistics[i].totalRfipOil = std::numeric_limits<double>::infinity();
|
|
if ( rfipGas.empty() ) statistics[i].totalRfipGas = std::numeric_limits<double>::infinity();
|
|
if ( sfipOil.empty() ) statistics[i].totalSfipOil = std::numeric_limits<double>::infinity();
|
|
if ( sfipGas.empty() ) statistics[i].totalSfipGas = std::numeric_limits<double>::infinity();
|
|
}
|
|
|
|
return statistics;
|
|
}
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
///
|
|
//--------------------------------------------------------------------------------------------------
|
|
RimRegularGridCase* RigWellTargetMapping::generateEnsembleCandidates( RimEclipseCaseEnsemble& ensemble,
|
|
size_t timeStepIdx,
|
|
const cvf::Vec3st& resultGridCellCount,
|
|
VolumeType volumeType,
|
|
VolumesType volumesType,
|
|
VolumeResultType volumeResultType,
|
|
const RigFloodingSettings& floodingSettings,
|
|
const ClusteringLimits& limits )
|
|
{
|
|
RiaLogging::debug( "Generating ensemble statistics" );
|
|
|
|
caf::ProgressInfo progInfo( ensemble.cases().size() * 2, "Generating ensemble statistics" );
|
|
|
|
for ( auto eclipseCase : ensemble.cases() )
|
|
{
|
|
auto task = progInfo.task( "Generating realization statistics.", 1 );
|
|
|
|
generateCandidates( eclipseCase, timeStepIdx, volumeType, volumesType, volumeResultType, floodingSettings, limits, false );
|
|
}
|
|
|
|
cvf::BoundingBox boundingBox;
|
|
for ( auto eclipseCase : ensemble.cases() )
|
|
{
|
|
cvf::BoundingBox bb = computeBoundingBoxForResult( *eclipseCase, RigWellTargetMapping::wellTargetResultName(), timeStepIdx );
|
|
boundingBox.add( bb );
|
|
}
|
|
|
|
RiaLogging::debug(
|
|
QString( "Clusters bounding box min: [%1 %2 %3]" ).arg( boundingBox.min().x() ).arg( boundingBox.min().y() ).arg( boundingBox.min().z() ) );
|
|
RiaLogging::debug(
|
|
QString( "Clusters bounding box max: [%1 %2 %3]" ).arg( boundingBox.max().x() ).arg( boundingBox.max().y() ).arg( boundingBox.max().z() ) );
|
|
|
|
RimRegularGridCase* targetCase = new RimRegularGridCase;
|
|
targetCase->setBoundingBox( boundingBox );
|
|
targetCase->setCellCount( resultGridCellCount );
|
|
targetCase->createModel();
|
|
|
|
std::vector<int> occurrence;
|
|
|
|
std::map<QString, std::vector<std::vector<double>>> resultNamesAndSamples;
|
|
resultNamesAndSamples["TOTAL_PORV_SOIL"] = {};
|
|
resultNamesAndSamples["TOTAL_PORV_SGAS"] = {};
|
|
resultNamesAndSamples["TOTAL_PORV_SOIL_SGAS"] = {};
|
|
resultNamesAndSamples["TOTAL_FIPOIL"] = {};
|
|
resultNamesAndSamples["TOTAL_FIPGAS"] = {};
|
|
resultNamesAndSamples["TOTAL_RFIPOIL"] = {};
|
|
resultNamesAndSamples["TOTAL_RFIPGAS"] = {};
|
|
resultNamesAndSamples["TOTAL_SFIPOIL"] = {};
|
|
resultNamesAndSamples["TOTAL_SFIPGAS"] = {};
|
|
|
|
for ( auto eclipseCase : ensemble.cases() )
|
|
{
|
|
auto task = progInfo.task( "Accumulating results.", 1 );
|
|
|
|
accumulateResultsForSingleCase( *eclipseCase, *targetCase, resultNamesAndSamples, occurrence, timeStepIdx );
|
|
}
|
|
|
|
auto createFractionVector = []( const std::vector<int>& occurrence, int maxRealizationCount ) -> std::vector<double>
|
|
{
|
|
std::vector<double> fractions( occurrence.size() );
|
|
std::transform( occurrence.begin(),
|
|
occurrence.end(),
|
|
fractions.begin(),
|
|
[maxRealizationCount]( int value ) { return static_cast<double>( value ) / maxRealizationCount; } );
|
|
|
|
return fractions;
|
|
};
|
|
|
|
createStaticResultVector( *targetCase, "OCCURRENCE", occurrence );
|
|
std::vector<double> probability = createFractionVector( occurrence, static_cast<int>( ensemble.cases().size() ) );
|
|
createStaticResultVector( *targetCase, "PROBABILITY", probability );
|
|
|
|
for ( auto [resultName, vec] : resultNamesAndSamples )
|
|
{
|
|
computeStatisticsAndCreateVectors( *targetCase, resultName, vec );
|
|
}
|
|
|
|
return targetCase;
|
|
}
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
///
|
|
//--------------------------------------------------------------------------------------------------
|
|
void RigWellTargetMapping::computeStatisticsAndCreateVectors( RimEclipseCase& targetCase,
|
|
const QString& resultName,
|
|
const std::vector<std::vector<double>>& vec )
|
|
{
|
|
const RigCaseCellResultsData* targetResultsData = targetCase.results( RiaDefines::PorosityModelType::MATRIX_MODEL );
|
|
if ( !targetResultsData ) return;
|
|
|
|
const RigActiveCellInfo* targetActiveCellInfo = targetResultsData->activeCellInfo();
|
|
if ( !targetActiveCellInfo ) return;
|
|
|
|
const size_t targetNumActiveCells = targetActiveCellInfo->reservoirActiveCellCount();
|
|
|
|
int nCells = static_cast<int>( targetNumActiveCells );
|
|
std::vector<double> p10Results( nCells, std::numeric_limits<double>::infinity() );
|
|
std::vector<double> p50Results( nCells, std::numeric_limits<double>::infinity() );
|
|
std::vector<double> p90Results( nCells, std::numeric_limits<double>::infinity() );
|
|
std::vector<double> meanResults( nCells, std::numeric_limits<double>::infinity() );
|
|
std::vector<double> minResults( nCells, std::numeric_limits<double>::infinity() );
|
|
std::vector<double> maxResults( nCells, std::numeric_limits<double>::infinity() );
|
|
|
|
#pragma omp parallel for
|
|
for ( int i = 0; i < nCells; i++ )
|
|
{
|
|
size_t numSamples = vec.size();
|
|
std::vector<double> samples( numSamples, 0.0 );
|
|
for ( size_t s = 0; s < numSamples; s++ )
|
|
samples[s] = vec[s][i];
|
|
|
|
double p10, p50, p90, mean;
|
|
RigStatisticsMath::calculateStatisticsCurves( samples, &p10, &p50, &p90, &mean, RigStatisticsMath::PercentileStyle::SWITCHED );
|
|
|
|
if ( RigStatisticsTools::isValidNumber( p10 ) ) p10Results[i] = p10;
|
|
if ( RigStatisticsTools::isValidNumber( p50 ) ) p50Results[i] = p50;
|
|
if ( RigStatisticsTools::isValidNumber( p90 ) ) p90Results[i] = p90;
|
|
if ( RigStatisticsTools::isValidNumber( mean ) ) meanResults[i] = mean;
|
|
|
|
double minValue = RigStatisticsTools::minimumValue( samples );
|
|
if ( RigStatisticsTools::isValidNumber( minValue ) && minValue < std::numeric_limits<double>::max() ) minResults[i] = minValue;
|
|
|
|
double maxValue = RigStatisticsTools::maximumValue( samples );
|
|
if ( RigStatisticsTools::isValidNumber( maxValue ) && maxValue > -std::numeric_limits<double>::max() ) maxResults[i] = maxValue;
|
|
}
|
|
|
|
createResultVectorIfDefined( targetCase, resultName + "_P10", p10Results );
|
|
createResultVectorIfDefined( targetCase, resultName + "_P50", p50Results );
|
|
createResultVectorIfDefined( targetCase, resultName + "_P90", p90Results );
|
|
createResultVectorIfDefined( targetCase, resultName + "_MEAN", meanResults );
|
|
createResultVectorIfDefined( targetCase, resultName + "_MIN", minResults );
|
|
createResultVectorIfDefined( targetCase, resultName + "_MAX", maxResults );
|
|
}
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
///
|
|
//--------------------------------------------------------------------------------------------------
|
|
void RigWellTargetMapping::accumulateResultsForSingleCase( RimEclipseCase& eclipseCase,
|
|
RimEclipseCase& targetCase,
|
|
std::map<QString, std::vector<std::vector<double>>>& resultNamesAndSamples,
|
|
std::vector<int>& occupancy,
|
|
size_t timeStepIdx )
|
|
{
|
|
RigCaseCellResultsData* resultsData = eclipseCase.results( RiaDefines::PorosityModelType::MATRIX_MODEL );
|
|
if ( !resultsData ) return;
|
|
const RigMainGrid* mainGrid = eclipseCase.mainGrid();
|
|
if ( !mainGrid ) return;
|
|
const RigActiveCellInfo* activeCellInfo = resultsData->activeCellInfo();
|
|
if ( !activeCellInfo ) return;
|
|
|
|
const RigCaseCellResultsData* targetResultsData = targetCase.results( RiaDefines::PorosityModelType::MATRIX_MODEL );
|
|
|
|
const RigActiveCellInfo* targetActiveCellInfo = targetResultsData->activeCellInfo();
|
|
|
|
const size_t targetNumReservoirCells = targetActiveCellInfo->reservoirCellCount();
|
|
const size_t targetNumActiveCells = targetActiveCellInfo->reservoirActiveCellCount();
|
|
|
|
occupancy.resize( targetNumActiveCells, 0 );
|
|
|
|
RigEclipseResultAddress clustersNumAddress( RiaDefines::ResultCatType::GENERATED, RigWellTargetMapping::wellTargetResultName() );
|
|
resultsData->ensureKnownResultLoaded( clustersNumAddress );
|
|
const std::vector<double>& clusterNum = resultsData->cellScalarResults( clustersNumAddress, timeStepIdx );
|
|
|
|
std::map<QString, const std::vector<double>*> namedInputVector;
|
|
|
|
for ( const auto& [resultName, vec] : resultNamesAndSamples )
|
|
{
|
|
RigEclipseResultAddress resultAddress( RiaDefines::ResultCatType::GENERATED, resultName );
|
|
resultsData->ensureKnownResultLoaded( resultAddress );
|
|
const std::vector<double>& resultVector = resultsData->cellScalarResults( resultAddress, timeStepIdx );
|
|
namedInputVector[resultName] = &resultVector;
|
|
}
|
|
|
|
std::map<QString, std::vector<double>> namedOutputVector;
|
|
for ( const auto& [resultName, vec] : resultNamesAndSamples )
|
|
{
|
|
namedOutputVector[resultName] = std::vector( targetNumActiveCells, std::numeric_limits<double>::infinity() );
|
|
}
|
|
|
|
for ( size_t targetCellIdx = 0; targetCellIdx < targetNumReservoirCells; targetCellIdx++ )
|
|
{
|
|
const RigCell& nativeCell = targetCase.mainGrid()->cell( targetCellIdx );
|
|
cvf::Vec3d cellCenter = nativeCell.center();
|
|
|
|
size_t targetResultIndex = targetActiveCellInfo->cellResultIndex( targetCellIdx );
|
|
|
|
size_t cellIdx = mainGrid->findReservoirCellIndexFromPoint( cellCenter );
|
|
if ( cellIdx != cvf::UNDEFINED_SIZE_T && activeCellInfo->isActive( cellIdx ) && targetResultIndex != cvf::UNDEFINED_SIZE_T )
|
|
{
|
|
size_t resultIndex = resultsData->activeCellInfo()->cellResultIndex( cellIdx );
|
|
if ( !std::isinf( clusterNum[resultIndex] ) && clusterNum[resultIndex] > 0 )
|
|
{
|
|
occupancy[targetResultIndex]++;
|
|
for ( const auto& [resultName, vec] : resultNamesAndSamples )
|
|
{
|
|
namedOutputVector[resultName][targetResultIndex] = namedInputVector[resultName]->at( resultIndex );
|
|
}
|
|
}
|
|
}
|
|
}
|
|
|
|
for ( const auto& [resultName, vec] : resultNamesAndSamples )
|
|
{
|
|
resultNamesAndSamples[resultName].push_back( namedOutputVector[resultName] );
|
|
}
|
|
}
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
///
|
|
//--------------------------------------------------------------------------------------------------
|
|
cvf::BoundingBox RigWellTargetMapping::computeBoundingBoxForResult( RimEclipseCase& eclipseCase, const QString& resultName, size_t timeStepIndex )
|
|
{
|
|
RigCaseCellResultsData* resultsData = eclipseCase.results( RiaDefines::PorosityModelType::MATRIX_MODEL );
|
|
const RigMainGrid* mainGrid = eclipseCase.mainGrid();
|
|
const RigActiveCellInfo* activeCellInfo = resultsData->activeCellInfo();
|
|
const size_t numReservoirCells = activeCellInfo->reservoirCellCount();
|
|
|
|
RigEclipseResultAddress clustersNumAddress( RiaDefines::ResultCatType::GENERATED, resultName );
|
|
resultsData->ensureKnownResultLoaded( clustersNumAddress );
|
|
const std::vector<double>& clusterNum = resultsData->cellScalarResults( clustersNumAddress, timeStepIndex );
|
|
|
|
cvf::BoundingBox boundingBox;
|
|
for ( size_t reservoirCellIndex = 0; reservoirCellIndex < numReservoirCells; reservoirCellIndex++ )
|
|
{
|
|
size_t targetResultIndex = activeCellInfo->cellResultIndex( reservoirCellIndex );
|
|
if ( reservoirCellIndex != cvf::UNDEFINED_SIZE_T && activeCellInfo->isActive( reservoirCellIndex ) &&
|
|
targetResultIndex != cvf::UNDEFINED_SIZE_T && !std::isinf( clusterNum[targetResultIndex] ) && clusterNum[targetResultIndex] > 0 )
|
|
{
|
|
const RigCell& nativeCell = mainGrid->cell( reservoirCellIndex );
|
|
boundingBox.add( nativeCell.boundingBox() );
|
|
}
|
|
}
|
|
|
|
return boundingBox;
|
|
}
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
///
|
|
//--------------------------------------------------------------------------------------------------
|
|
std::list<std::pair<std::pair<size_t, RigWellTargetMapping::CellFaceType>, size_t>>
|
|
RigWellTargetMapping::nncConnectionCellAndResult( size_t cellIdx, RigMainGrid* mainGrid )
|
|
{
|
|
std::list<std::pair<std::pair<size_t, CellFaceType>, size_t>> foundCells;
|
|
|
|
if ( mainGrid->nncData() == nullptr ) return foundCells;
|
|
|
|
auto& connections = mainGrid->nncData()->allConnections();
|
|
for ( size_t i = 0; i < connections.size(); i++ )
|
|
{
|
|
if ( connections[i].c1GlobIdx() == cellIdx )
|
|
{
|
|
foundCells.push_back( std::make_pair( std::make_pair( connections[i].c2GlobIdx(), connections[i].face() ), i ) );
|
|
}
|
|
}
|
|
|
|
return foundCells;
|
|
}
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
///
|
|
//--------------------------------------------------------------------------------------------------
|
|
QString RigWellTargetMapping::wellTargetResultName()
|
|
{
|
|
return "WELL_TARGET_NUM";
|
|
}
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
///
|
|
//--------------------------------------------------------------------------------------------------
|
|
bool RigWellTargetMapping::isSaturationSufficient( const VolumeType volumeType, const DataContainer& data, const ClusteringLimits& limits, size_t idx )
|
|
{
|
|
bool needsValidOil = volumeType == VolumeType::OIL || volumeType == VolumeType::HYDROCARBON;
|
|
bool needsValidGas = volumeType == VolumeType::GAS || volumeType == VolumeType::HYDROCARBON;
|
|
// For hydrocarbon it is enough that one of the saturations is above the limit.
|
|
if ( needsValidOil && data.saturationOil[idx] >= limits.saturationOil ) return true;
|
|
if ( needsValidGas && data.saturationGas[idx] >= limits.saturationGas ) return true;
|
|
return false;
|
|
}
|