mirror of
https://github.com/OPM/ResInsight.git
synced 2024-12-28 18:01:08 -06:00
334 lines
16 KiB
C++
334 lines
16 KiB
C++
/////////////////////////////////////////////////////////////////////////////////
|
|
//
|
|
// Copyright (C) 2018- 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 "RigEclipseContourMapProjection.h"
|
|
|
|
#include "RigActiveCellInfo.h"
|
|
#include "RigCaseCellResultsData.h"
|
|
#include "RigCell.h"
|
|
#include "RigCellGeometryTools.h"
|
|
#include "RigContourMapCalculator.h"
|
|
#include "RigContourMapGrid.h"
|
|
#include "RigEclipseCaseData.h"
|
|
#include "RigEclipseResultAddress.h"
|
|
#include "RigHexIntersectionTools.h"
|
|
#include "RigMainGrid.h"
|
|
|
|
#include <algorithm>
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
///
|
|
//--------------------------------------------------------------------------------------------------
|
|
RigEclipseContourMapProjection::RigEclipseContourMapProjection( const RigContourMapGrid& contourMapGrid,
|
|
RigEclipseCaseData& eclipseCaseData,
|
|
RigCaseCellResultsData& resultData )
|
|
: RigContourMapProjection( contourMapGrid )
|
|
, m_eclipseCaseData( eclipseCaseData )
|
|
, m_resultData( resultData )
|
|
, m_kLayers( 0u )
|
|
, m_useActiveCellInfo( true )
|
|
{
|
|
m_mainGrid = m_eclipseCaseData.mainGrid();
|
|
m_activeCellInfo = m_eclipseCaseData.activeCellInfo( RiaDefines::PorosityModelType::MATRIX_MODEL );
|
|
m_kLayers = m_mainGrid->cellCountK();
|
|
}
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
///
|
|
//--------------------------------------------------------------------------------------------------
|
|
RigEclipseContourMapProjection::~RigEclipseContourMapProjection()
|
|
{
|
|
}
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
///
|
|
//--------------------------------------------------------------------------------------------------
|
|
void RigEclipseContourMapProjection::generateAndSaveResults( const RigEclipseResultAddress& resultAddress,
|
|
RigContourMapCalculator::ResultAggregationType resultAggregation,
|
|
int timeStep )
|
|
{
|
|
std::tie( m_useActiveCellInfo, m_aggregatedResults ) =
|
|
generateResults( *this, m_contourMapGrid, m_resultData, resultAddress, resultAggregation, timeStep );
|
|
}
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
///
|
|
//--------------------------------------------------------------------------------------------------
|
|
std::vector<double> RigEclipseContourMapProjection::generateResults( const RigEclipseResultAddress& resultAddress,
|
|
RigContourMapCalculator::ResultAggregationType resultAggregation,
|
|
int timeStep ) const
|
|
{
|
|
std::pair<bool, std::vector<double>> result =
|
|
generateResults( *this, m_contourMapGrid, m_resultData, resultAddress, resultAggregation, timeStep );
|
|
return result.second;
|
|
}
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
///
|
|
//--------------------------------------------------------------------------------------------------
|
|
std::pair<bool, std::vector<double>>
|
|
RigEclipseContourMapProjection::generateResults( const RigEclipseContourMapProjection& contourMapProjection,
|
|
const RigContourMapGrid& contourMapGrid,
|
|
RigCaseCellResultsData& resultData,
|
|
const RigEclipseResultAddress& resultAddress,
|
|
RigContourMapCalculator::ResultAggregationType resultAggregation,
|
|
int timeStep )
|
|
{
|
|
size_t nCells = contourMapProjection.numberOfCells();
|
|
|
|
std::vector<double> aggregatedResults = std::vector<double>( nCells, std::numeric_limits<double>::infinity() );
|
|
|
|
bool useActiveCellInfo = resultAddress.isValid() && resultData.hasResultEntry( resultAddress ) &&
|
|
resultData.isUsingGlobalActiveIndex( resultAddress );
|
|
|
|
auto isTernaryResult = []( const RigEclipseResultAddress& address ) -> bool
|
|
{
|
|
return address.resultCatType() == RiaDefines::ResultCatType::DYNAMIC_NATIVE &&
|
|
( address.resultName().compare( RiaResultNames::ternarySaturationResultName(), Qt::CaseInsensitive ) == 0 );
|
|
};
|
|
|
|
if ( !isTernaryResult( resultAddress ) )
|
|
{
|
|
std::vector<double> gridResultValues;
|
|
if ( RigContourMapCalculator::isColumnResult( resultAggregation ) )
|
|
{
|
|
resultData.ensureKnownResultLoaded( RigEclipseResultAddress( RiaDefines::ResultCatType::STATIC_NATIVE, "PORO" ) );
|
|
resultData.ensureKnownResultLoaded( RigEclipseResultAddress( RiaDefines::ResultCatType::STATIC_NATIVE, "NTG" ) );
|
|
resultData.ensureKnownResultLoaded( RigEclipseResultAddress( RiaDefines::ResultCatType::STATIC_NATIVE, "DZ" ) );
|
|
if ( resultAggregation == RigContourMapCalculator::OIL_COLUMN || resultAggregation == RigContourMapCalculator::HYDROCARBON_COLUMN )
|
|
{
|
|
resultData.ensureKnownResultLoaded(
|
|
RigEclipseResultAddress( RiaDefines::ResultCatType::DYNAMIC_NATIVE, RiaResultNames::soil() ) );
|
|
}
|
|
if ( resultAggregation == RigContourMapCalculator::GAS_COLUMN || resultAggregation == RigContourMapCalculator::HYDROCARBON_COLUMN )
|
|
{
|
|
resultData.ensureKnownResultLoaded(
|
|
RigEclipseResultAddress( RiaDefines::ResultCatType::DYNAMIC_NATIVE, RiaResultNames::sgas() ) );
|
|
}
|
|
gridResultValues = calculateColumnResult( resultData, resultAggregation, timeStep );
|
|
}
|
|
else
|
|
{
|
|
// TODO: this was RimEclipseCellColors->hasStaticResult()
|
|
if ( resultAddress.resultCatType() == RiaDefines::ResultCatType::STATIC_NATIVE && timeStep > 0 ) timeStep = 0;
|
|
|
|
resultData.ensureKnownResultLoaded( resultAddress );
|
|
// When loading a project file, grid calculator results are not computed the first time this function is
|
|
// called. Must check if result is loaded. See RimReloadCaseTools::updateAll3dViews()
|
|
if ( resultAddress.isValid() && resultData.hasResultEntry( resultAddress ) && resultData.isResultLoaded( resultAddress ) )
|
|
{
|
|
gridResultValues = resultData.cellScalarResults( resultAddress, timeStep );
|
|
}
|
|
}
|
|
|
|
if ( !gridResultValues.empty() )
|
|
{
|
|
#pragma omp parallel for
|
|
for ( int index = 0; index < static_cast<int>( nCells ); ++index )
|
|
{
|
|
cvf::Vec2ui ij = contourMapGrid.ijFromCellIndex( index );
|
|
const std::vector<std::pair<size_t, double>>& matchingCells = contourMapProjection.cellsAtIJ( ij.x(), ij.y() );
|
|
aggregatedResults[index] =
|
|
RigContourMapCalculator::calculateValueInMapCell( contourMapProjection, matchingCells, gridResultValues, resultAggregation );
|
|
}
|
|
}
|
|
}
|
|
|
|
return { useActiveCellInfo, aggregatedResults };
|
|
}
|
|
|
|
std::vector<double> RigEclipseContourMapProjection::calculateColumnResult( RigContourMapCalculator::ResultAggregationType resultAggregation,
|
|
int timeStep ) const
|
|
{
|
|
return calculateColumnResult( m_resultData, resultAggregation, timeStep );
|
|
}
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
///
|
|
//--------------------------------------------------------------------------------------------------
|
|
std::vector<double> RigEclipseContourMapProjection::calculateColumnResult( RigCaseCellResultsData& resultData,
|
|
RigContourMapCalculator::ResultAggregationType resultAggregation,
|
|
int timeStep )
|
|
{
|
|
bool hasPoroResult = resultData.hasResultEntry( RigEclipseResultAddress( RiaDefines::ResultCatType::STATIC_NATIVE, "PORO" ) );
|
|
bool hasNtgResult = resultData.hasResultEntry( RigEclipseResultAddress( RiaDefines::ResultCatType::STATIC_NATIVE, "NTG" ) );
|
|
bool hasDzResult = resultData.hasResultEntry( RigEclipseResultAddress( RiaDefines::ResultCatType::STATIC_NATIVE, "DZ" ) );
|
|
|
|
if ( !( hasPoroResult && hasNtgResult && hasDzResult ) )
|
|
{
|
|
return std::vector<double>();
|
|
}
|
|
|
|
const std::vector<double>& poroResults =
|
|
resultData.cellScalarResults( RigEclipseResultAddress( RiaDefines::ResultCatType::STATIC_NATIVE, "PORO" ), 0 );
|
|
const std::vector<double>& ntgResults =
|
|
resultData.cellScalarResults( RigEclipseResultAddress( RiaDefines::ResultCatType::STATIC_NATIVE, "NTG" ), 0 );
|
|
const std::vector<double>& dzResults =
|
|
resultData.cellScalarResults( RigEclipseResultAddress( RiaDefines::ResultCatType::STATIC_NATIVE, "DZ" ), 0 );
|
|
|
|
CVF_ASSERT( poroResults.size() == ntgResults.size() && ntgResults.size() == dzResults.size() );
|
|
|
|
std::vector<double> resultValues( poroResults.size(), 0.0 );
|
|
|
|
if ( resultAggregation == RigContourMapCalculator::OIL_COLUMN || resultAggregation == RigContourMapCalculator::HYDROCARBON_COLUMN )
|
|
{
|
|
const std::vector<double>& soilResults =
|
|
resultData.cellScalarResults( RigEclipseResultAddress( RiaDefines::ResultCatType::DYNAMIC_NATIVE, RiaResultNames::soil() ),
|
|
timeStep );
|
|
for ( size_t cellResultIdx = 0; cellResultIdx < resultValues.size(); ++cellResultIdx )
|
|
{
|
|
resultValues[cellResultIdx] = soilResults[cellResultIdx];
|
|
}
|
|
}
|
|
|
|
if ( resultAggregation == RigContourMapCalculator::GAS_COLUMN || resultAggregation == RigContourMapCalculator::HYDROCARBON_COLUMN )
|
|
{
|
|
bool hasGasResult =
|
|
resultData.hasResultEntry( RigEclipseResultAddress( RiaDefines::ResultCatType::DYNAMIC_NATIVE, RiaResultNames::sgas() ) );
|
|
if ( hasGasResult )
|
|
{
|
|
const std::vector<double>& sgasResults =
|
|
resultData.cellScalarResults( RigEclipseResultAddress( RiaDefines::ResultCatType::DYNAMIC_NATIVE, RiaResultNames::sgas() ),
|
|
timeStep );
|
|
for ( size_t cellResultIdx = 0; cellResultIdx < resultValues.size(); ++cellResultIdx )
|
|
{
|
|
resultValues[cellResultIdx] += sgasResults[cellResultIdx];
|
|
}
|
|
}
|
|
}
|
|
|
|
for ( size_t cellResultIdx = 0; cellResultIdx < resultValues.size(); ++cellResultIdx )
|
|
{
|
|
resultValues[cellResultIdx] *= poroResults[cellResultIdx] * ntgResults[cellResultIdx] * dzResults[cellResultIdx];
|
|
}
|
|
|
|
return resultValues;
|
|
}
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
///
|
|
//--------------------------------------------------------------------------------------------------
|
|
std::vector<size_t> RigEclipseContourMapProjection::findIntersectingCells( const cvf::BoundingBox& bbox ) const
|
|
{
|
|
return m_mainGrid->findIntersectingCells( bbox );
|
|
}
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
///
|
|
//--------------------------------------------------------------------------------------------------
|
|
size_t RigEclipseContourMapProjection::kLayer( size_t globalCellIdx ) const
|
|
{
|
|
const RigCell& cell = m_mainGrid->cell( globalCellIdx );
|
|
size_t mainGridCellIdx = cell.mainGridCellIndex();
|
|
size_t i, j, k;
|
|
m_mainGrid->ijkFromCellIndex( mainGridCellIdx, &i, &j, &k );
|
|
return k;
|
|
}
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
///
|
|
//--------------------------------------------------------------------------------------------------
|
|
size_t RigEclipseContourMapProjection::kLayers() const
|
|
{
|
|
return m_kLayers;
|
|
}
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
///
|
|
//--------------------------------------------------------------------------------------------------
|
|
double RigEclipseContourMapProjection::calculateOverlapVolume( size_t globalCellIdx, const cvf::BoundingBox& bbox ) const
|
|
{
|
|
std::array<cvf::Vec3d, 8> hexCorners;
|
|
|
|
const RigCell& cell = m_mainGrid->cell( globalCellIdx );
|
|
|
|
size_t localCellIdx = cell.gridLocalCellIndex();
|
|
RigGridBase* localGrid = cell.hostGrid();
|
|
|
|
localGrid->cellCornerVertices( localCellIdx, hexCorners.data() );
|
|
|
|
cvf::BoundingBox overlapBBox;
|
|
std::array<cvf::Vec3d, 8> overlapCorners;
|
|
if ( RigCellGeometryTools::estimateHexOverlapWithBoundingBox( hexCorners, bbox, &overlapCorners, &overlapBBox ) )
|
|
{
|
|
double overlapVolume = RigCellGeometryTools::calculateCellVolume( overlapCorners );
|
|
return overlapVolume;
|
|
}
|
|
return 0.0;
|
|
}
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
///
|
|
//--------------------------------------------------------------------------------------------------
|
|
double RigEclipseContourMapProjection::calculateRayLengthInCell( size_t globalCellIdx,
|
|
const cvf::Vec3d& highestPoint,
|
|
const cvf::Vec3d& lowestPoint ) const
|
|
{
|
|
std::array<cvf::Vec3d, 8> hexCorners;
|
|
|
|
const RigCell& cell = m_mainGrid->cell( globalCellIdx );
|
|
|
|
size_t localCellIdx = cell.gridLocalCellIndex();
|
|
RigGridBase* localGrid = cell.hostGrid();
|
|
|
|
localGrid->cellCornerVertices( localCellIdx, hexCorners.data() );
|
|
std::vector<HexIntersectionInfo> intersections;
|
|
|
|
if ( RigHexIntersectionTools::lineHexCellIntersection( highestPoint, lowestPoint, hexCorners.data(), 0, &intersections ) )
|
|
{
|
|
double lengthInCell = ( intersections.back().m_intersectionPoint - intersections.front().m_intersectionPoint ).length();
|
|
return lengthInCell;
|
|
}
|
|
return 0.0;
|
|
}
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
///
|
|
//--------------------------------------------------------------------------------------------------
|
|
double RigEclipseContourMapProjection::getParameterWeightForCell( size_t cellResultIdx, const std::vector<double>& cellWeights ) const
|
|
{
|
|
if ( cellWeights.empty() ) return 1.0;
|
|
|
|
double result = std::max( cellWeights[cellResultIdx], 0.0 );
|
|
if ( result < 1.0e-6 )
|
|
{
|
|
result = 0.0;
|
|
}
|
|
return result;
|
|
}
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
///
|
|
//--------------------------------------------------------------------------------------------------
|
|
size_t RigEclipseContourMapProjection::gridResultIndex( size_t globalCellIdx ) const
|
|
{
|
|
if ( m_useActiveCellInfo ) return m_activeCellInfo->cellResultIndex( globalCellIdx );
|
|
|
|
return globalCellIdx;
|
|
}
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
///
|
|
//--------------------------------------------------------------------------------------------------
|
|
std::vector<bool> RigEclipseContourMapProjection::getMapCellVisibility( int viewStepIndex,
|
|
RigContourMapCalculator::ResultAggregationType resultAggregation )
|
|
|
|
{
|
|
return std::vector<bool>( numberOfCells(), true );
|
|
}
|