#12957 Sector Model Export: add non-uniform cell refinement

Add a new "Non-Uniform" refinement mode that coexists with the existing
uniform refinement. Users can specify a cell range, direction, and
fractional subdivision widths to create variable-resolution grids.

Core infrastructure:
- RigNonUniformRefinement: per-dimension refinement with cumulative
  fractions and offset tables for O(1)/O(log n) lookups
- RiaCellDividingTools: non-uniform createHexCornerCoords overload
- RigGridExportAdapter: non-uniform constructor, mapping, and
  coordinate transformation overloads

Pipeline integration:
- RigSimulationInputSettings: non-uniform refinement storage
- RigSimulationInputTool: non-uniform paths for grid construction,
  data cropping, BOX context, well keywords, and WELLDIMS
- RifEclipseInputFileTools: non-uniform extractKeywordData overload

UI:
- RicExportSectorModelUi: Uniform/Non-Uniform mode selector with
  direction, cell range, and fractional widths fields
This commit is contained in:
Kristian Bendiksen committed 2026-03-30 11:24:56 +02:00
1 parent 871db2e143
commit 06f985b1ed
17 files changed
+1496 -63

No files matched your search

@@ -32,6 +32,11 @@
//--------------------------------------------------------------------------------------------------
std::vector<cvf::Vec3d> splitLine( cvf::Vec3d ptStart, cvf::Vec3d ptEnd, size_t partCount );
//--------------------------------------------------------------------------------------------------
/// Splits a line using cumulative fractions instead of equal parts
//--------------------------------------------------------------------------------------------------
std::vector<cvf::Vec3d> splitLineNonUniform( cvf::Vec3d ptStart, cvf::Vec3d ptEnd, const std::vector<double>& cumulativeFractions );
//--------------------------------------------------------------------------------------------------
/// Calculates all points on a face described by edge points from all four edges.
/// The result is a grid of points including the given edge points
@@ -50,6 +55,15 @@ std::vector<std::vector<cvf::Vec3d>> calcFacePoints( const std::vector<cvf::Vec3
const std::vector<cvf::Vec3d> edgeYPtsLow,
const std::vector<cvf::Vec3d> edgeYPtsHigh );
//--------------------------------------------------------------------------------------------------
/// Calculates face points using non-uniform X-direction fractions for interior rows
//--------------------------------------------------------------------------------------------------
std::vector<std::vector<cvf::Vec3d>> calcFacePointsNonUniform( const std::vector<cvf::Vec3d>& edgeXPtsLow,
const std::vector<cvf::Vec3d>& edgeXPtsHigh,
const std::vector<cvf::Vec3d>& edgeYPtsLow,
const std::vector<cvf::Vec3d>& edgeYPtsHigh,
const std::vector<double>& xCumulativeFractions );
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
@@ -204,3 +218,135 @@ std::vector<std::vector<cvf::Vec3d>> calcFacePoints( const std::vector<cvf::Vec3
pts.push_back( edgeXPtsHigh );
return pts;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
std::vector<cvf::Vec3d> splitLineNonUniform( cvf::Vec3d ptStart, cvf::Vec3d ptEnd, const std::vector<double>& cumulativeFractions )
{
std::vector<cvf::Vec3d> pts = { ptStart };
for ( size_t i = 0; i < cumulativeFractions.size() - 1; ++i )
{
double f = cumulativeFractions[i];
pts.push_back( cvf::Vec3d( ptStart.x() + ( ptEnd.x() - ptStart.x() ) * f,
ptStart.y() + ( ptEnd.y() - ptStart.y() ) * f,
ptStart.z() + ( ptEnd.z() - ptStart.z() ) * f ) );
}
pts.push_back( ptEnd );
return pts;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
std::vector<std::vector<cvf::Vec3d>> calcFacePointsNonUniform( const std::vector<cvf::Vec3d>& edgeXPtsLow,
const std::vector<cvf::Vec3d>& edgeXPtsHigh,
const std::vector<cvf::Vec3d>& edgeYPtsLow,
const std::vector<cvf::Vec3d>& edgeYPtsHigh,
const std::vector<double>& xCumulativeFractions )
{
CVF_ASSERT( edgeYPtsLow.size() == edgeYPtsHigh.size() );
size_t ySize = edgeYPtsLow.size();
std::vector<std::vector<cvf::Vec3d>> pts;
pts.push_back( edgeXPtsLow );
for ( size_t y = 1; y < ySize - 1; y++ )
{
auto interiorPts = splitLineNonUniform( edgeYPtsLow[y], edgeYPtsHigh[y], xCumulativeFractions );
pts.push_back( interiorPts );
}
pts.push_back( edgeXPtsHigh );
return pts;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
std::vector<cvf::Vec3d> RiaCellDividingTools::createHexCornerCoords( const std::array<cvf::Vec3d, 8>& mainCellCorners,
const std::vector<double>& cumulativeFractionsX,
const std::vector<double>& cumulativeFractionsY,
const std::vector<double>& cumulativeFractionsZ )
{
size_t nx = cumulativeFractionsX.size();
size_t ny = cumulativeFractionsY.size();
size_t nz = cumulativeFractionsZ.size();
std::array<std::pair<size_t, size_t>, 12> edgeCorners = {
std::make_pair( 0, 1 ),
std::make_pair( 3, 2 ),
std::make_pair( 4, 5 ),
std::make_pair( 7, 6 ), // X
std::make_pair( 0, 3 ),
std::make_pair( 4, 7 ),
std::make_pair( 1, 2 ),
std::make_pair( 5, 6 ), // Y
std::make_pair( 0, 4 ),
std::make_pair( 1, 5 ),
std::make_pair( 3, 7 ),
std::make_pair( 2, 6 ), // Z
};
// Split edges using appropriate cumulative fractions per dimension
const std::vector<double>* dimFractions[3] = { &cumulativeFractionsX, &cumulativeFractionsY, &cumulativeFractionsZ };
std::array<std::vector<cvf::Vec3d>, 12> edgePoints;
for ( int i = 0; i < 12; i++ )
{
int dimIndex = i / 4;
edgePoints[i] =
splitLineNonUniform( mainCellCorners[edgeCorners[i].first], mainCellCorners[edgeCorners[i].second], *dimFractions[dimIndex] );
}
std::vector<std::vector<std::vector<cvf::Vec3d>>> nodes;
nodes.reserve( ( nx + 1 ) * ( ny + 1 ) * ( nz + 1 ) );
auto xyFacePtsLow = calcFacePointsNonUniform( edgePoints[0], edgePoints[1], edgePoints[4], edgePoints[6], cumulativeFractionsX );
auto xyFacePtsHigh = calcFacePointsNonUniform( edgePoints[2], edgePoints[3], edgePoints[5], edgePoints[7], cumulativeFractionsX );
auto yzFacePtsLow = calcFacePointsNonUniform( edgePoints[4], edgePoints[5], edgePoints[8], edgePoints[10], cumulativeFractionsY );
auto yzFacePtsHigh = calcFacePointsNonUniform( edgePoints[6], edgePoints[7], edgePoints[9], edgePoints[11], cumulativeFractionsY );
auto xzFacePtsLow = calcFacePointsNonUniform( edgePoints[0], edgePoints[2], edgePoints[8], edgePoints[9], cumulativeFractionsX );
auto xzFacePtsHigh = calcFacePointsNonUniform( edgePoints[1], edgePoints[3], edgePoints[10], edgePoints[11], cumulativeFractionsX );
nodes.push_back( xyFacePtsLow );
for ( size_t z = 1; z < nz; z++ )
{
auto xyFacePoints =
calcFacePointsNonUniform( xzFacePtsLow[z], xzFacePtsHigh[z], yzFacePtsLow[z], yzFacePtsHigh[z], cumulativeFractionsX );
nodes.push_back( xyFacePoints );
}
nodes.push_back( xyFacePtsHigh );
std::vector<cvf::Vec3d> coords;
coords.reserve( nx * ny * nz * 8 );
for ( size_t z = 1; z < nz + 1; z++ )
{
for ( size_t y = 1; y < ny + 1; y++ )
{
for ( size_t x = 1; x < nx + 1; x++ )
{
std::array<cvf::Vec3d, 8> cs;
cs[0] = nodes[z - 1][y - 1][x - 1];
cs[1] = nodes[z - 1][y - 1][x];
cs[2] = nodes[z - 1][y][x];
cs[3] = nodes[z - 1][y][x - 1];
cs[4] = nodes[z][y - 1][x - 1];
cs[5] = nodes[z][y - 1][x];
cs[6] = nodes[z][y][x];
cs[7] = nodes[z][y][x - 1];
coords.insert( coords.end(), cs.begin(), cs.end() );
}
}
}
return coords;
}
@@ -31,5 +31,11 @@ class RiaCellDividingTools
public:
static std::vector<cvf::Vec3d> createHexCornerCoords( std::array<cvf::Vec3d, 8> mainCellCorners, size_t nx, size_t ny, size_t nz );
// Non-uniform overload: takes cumulative fractions per dimension instead of integer counts
static std::vector<cvf::Vec3d> createHexCornerCoords( const std::array<cvf::Vec3d, 8>& mainCellCorners,
const std::vector<double>& cumulativeFractionsX,
const std::vector<double>& cumulativeFractionsY,
const std::vector<double>& cumulativeFractionsZ );
static double computeFlowDistance( const std::array<cvf::Vec3d, 8>& cellVertices, const cvf::Vec3d& areaCenter );
};
@@ -101,6 +101,11 @@ void RicExportSectorModelFeature::doExport( RicExportSectorModelUi* exportSettin
settings.setMax( exportSettings->max() );
settings.setRefinement( exportSettings->refinement() );
if ( exportSettings->hasNonUniformRefinement() )
{
settings.setNonUniformRefinement( exportSettings->nonUniformRefinement() );
}
std::vector<std::string> keywordsToRemove;
for ( auto& s : exportSettings->keywordsToRemove() )
{
@@ -44,6 +44,17 @@
#include <utility>
namespace caf
{
template <>
void AppEnum<RicExportSectorModelUi::RefinementMode>::setUp()
{
addItem( RicExportSectorModelUi::UNIFORM, "UNIFORM", "Uniform" );
addItem( RicExportSectorModelUi::NON_UNIFORM, "NON_UNIFORM", "Non-Uniform" );
setDefault( RicExportSectorModelUi::UNIFORM );
}
} // namespace caf
CAF_PDM_SOURCE_INIT( RicExportSectorModelUi, "RicExportSectorModelUi" );
//--------------------------------------------------------------------------------------------------
@@ -99,6 +110,20 @@ RicExportSectorModelUi::RicExportSectorModelUi()
m_refinementCountJ.setRange( 1, 10 );
m_refinementCountK.setRange( 1, 10 );
CAF_PDM_InitFieldNoDefault( &m_refinementMode, "RefinementMode", "Refinement Mode" );
m_refinementMode.uiCapability()->setUiEditorTypeName( caf::PdmUiRadioButtonEditor::uiEditorTypeName() );
CAF_PDM_InitField( &m_nonUniformDirection, "NonUniformDirection", 0, "Direction" );
m_nonUniformDirection.setRange( 0, 2 );
CAF_PDM_InitField( &m_nonUniformRangeStart, "NonUniformRangeStart", 1, "Cell Range Start" );
m_nonUniformRangeStart.setMinValue( 1 );
CAF_PDM_InitField( &m_nonUniformRangeEnd, "NonUniformRangeEnd", 1, "Cell Range End" );
m_nonUniformRangeEnd.setMinValue( 1 );
CAF_PDM_InitField( &m_nonUniformIntervals, "NonUniformIntervals", QString( "0.5, 0.5" ), "Fractional Widths" );
CAF_PDM_InitFieldNoDefault( &m_bcpropKeywords, "BcpropKeywords", "BCPROP Keywords" );
m_bcpropKeywords.uiCapability()->setUiEditorTypeName( caf::PdmUiTableViewEditor::uiEditorTypeName() );
m_bcpropKeywords.uiCapability()->setUiLabelPosition( caf::PdmUiItemInfo::LabelPosition::HIDDEN );
@@ -241,13 +266,35 @@ void RicExportSectorModelUi::defineUiOrdering( QString uiConfigName, caf::PdmUiO
uiOrdering.add( &m_refineGrid );
uiOrdering.addNewLabel( "" );
uiOrdering.add( &m_refinementCountI, { .newRow = true, .totalColumnSpan = 2, .leftLabelColumnSpan = 1 } );
uiOrdering.appendToRow( &m_refinementCountJ );
uiOrdering.appendToRow( &m_refinementCountK );
uiOrdering.add( &m_refinementMode );
uiOrdering.addNewLabel( "" );
m_refinementCountI.uiCapability()->setUiReadOnly( !m_refineGrid() );
m_refinementCountJ.uiCapability()->setUiReadOnly( !m_refineGrid() );
m_refinementCountK.uiCapability()->setUiReadOnly( !m_refineGrid() );
bool isEnabled = m_refineGrid();
if ( m_refinementMode() == UNIFORM )
{
uiOrdering.add( &m_refinementCountI, { .newRow = true, .totalColumnSpan = 2, .leftLabelColumnSpan = 1 } );
uiOrdering.appendToRow( &m_refinementCountJ );
uiOrdering.appendToRow( &m_refinementCountK );
m_refinementCountI.uiCapability()->setUiReadOnly( !isEnabled );
m_refinementCountJ.uiCapability()->setUiReadOnly( !isEnabled );
m_refinementCountK.uiCapability()->setUiReadOnly( !isEnabled );
}
else
{
uiOrdering.add( &m_nonUniformDirection );
uiOrdering.add( &m_nonUniformRangeStart );
uiOrdering.add( &m_nonUniformRangeEnd );
uiOrdering.add( &m_nonUniformIntervals );
m_nonUniformDirection.uiCapability()->setUiReadOnly( !isEnabled );
m_nonUniformRangeStart.uiCapability()->setUiReadOnly( !isEnabled );
m_nonUniformRangeEnd.uiCapability()->setUiReadOnly( !isEnabled );
m_nonUniformIntervals.uiCapability()->setUiReadOnly( !isEnabled );
}
m_refinementMode.uiCapability()->setUiReadOnly( !isEnabled );
}
else if ( uiConfigName == m_pageNames[WizardPageEnum::BoundaryConditions] )
{
@@ -466,13 +513,65 @@ void RicExportSectorModelUi::setMax( const caf::VecIjk0& max )
//--------------------------------------------------------------------------------------------------
cvf::Vec3st RicExportSectorModelUi::refinement() const
{
if ( !m_refineGrid() )
if ( !m_refineGrid() || m_refinementMode() == NON_UNIFORM )
{
return cvf::Vec3st( 1, 1, 1 );
}
return cvf::Vec3st( m_refinementCountI(), m_refinementCountJ(), m_refinementCountK() );
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
RicExportSectorModelUi::RefinementMode RicExportSectorModelUi::refinementMode() const
{
if ( !m_refineGrid() ) return UNIFORM;
return m_refinementMode();
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
RigNonUniformRefinement RicExportSectorModelUi::nonUniformRefinement() const
{
cvf::Vec3st sectorSize( max().x() - min().x() + 1, max().y() - min().y() + 1, max().z() - min().z() + 1 );
RigNonUniformRefinement result( sectorSize );
if ( !m_refineGrid() || m_refinementMode() != NON_UNIFORM ) return result;
// Parse the fractional widths from the text field
QStringList parts = m_nonUniformIntervals().split( ",", Qt::SkipEmptyParts );
std::vector<double> widths;
for ( const auto& part : parts )
{
bool ok = false;
double value = part.trimmed().toDouble( &ok );
if ( ok && value > 0.0 ) widths.push_back( value );
}
if ( widths.empty() ) return result;
auto dim = static_cast<RigNonUniformRefinement::Dimension>( m_nonUniformDirection() );
int start = std::max( 0, m_nonUniformRangeStart() - 1 ); // Convert to 0-based, clamp
int end = std::min( static_cast<int>( result.sectorSize( dim ) ) - 1, m_nonUniformRangeEnd() - 1 );
if ( start > end ) return result;
// Distribute the widths across the combined range of cells
result.distributeWidthsAcrossCells( dim, start, end, widths );
return result;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
bool RicExportSectorModelUi::hasNonUniformRefinement() const
{
return m_refineGrid() && m_refinementMode() == NON_UNIFORM;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
@@ -25,6 +25,7 @@
#include "cafVecIjk.h"
#include "RiaModelExportDefines.h"
#include "RigNonUniformRefinement.h"
#include <QString>
#include <QStringList>
@@ -43,8 +44,15 @@ class RicExportSectorModelUi : public caf::PdmObject
{
CAF_PDM_HEADER_INIT;
enum RefinementMode
{
UNIFORM,
NON_UNIFORM
};
using GridBoxSelectionEnum = caf::AppEnum<RiaModelExportDefines::GridBoxSelection>;
using BoundaryConditionEnum = caf::AppEnum<RiaModelExportDefines::BoundaryCondition>;
using RefinementModeEnum = caf::AppEnum<RefinementMode>;
public:
RicExportSectorModelUi();
@@ -63,6 +71,10 @@ public:
cvf::Vec3st refinement() const;
RefinementMode refinementMode() const;
RigNonUniformRefinement nonUniformRefinement() const;
bool hasNonUniformRefinement() const;
std::vector<QString> keywordsToRemove() const;
std::vector<RimKeywordBcprop*> bcpropKeywords() const;
@@ -126,6 +138,12 @@ private:
caf::PdmField<int> m_refinementCountJ;
caf::PdmField<int> m_refinementCountK;
caf::PdmField<RefinementModeEnum> m_refinementMode;
caf::PdmField<int> m_nonUniformDirection; // 0=I, 1=J, 2=K
caf::PdmField<int> m_nonUniformRangeStart;
caf::PdmField<int> m_nonUniformRangeEnd;
caf::PdmField<QString> m_nonUniformIntervals;
caf::PdmField<BoundaryConditionEnum> m_boundaryCondition;
caf::PdmChildArrayField<RimKeywordBcprop*> m_bcpropKeywords;
caf::PdmField<double> m_porvMultiplier;
@@ -335,6 +335,93 @@ std::expected<std::vector<double>, std::string> RifEclipseInputFileTools::extrac
return filteredResults;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
std::expected<std::vector<double>, std::string> RifEclipseInputFileTools::extractKeywordData( RigEclipseCaseData* eclipseCase,
const QString& keyword,
const cvf::Vec3st& min,
const cvf::Vec3st& maxIn,
const RigNonUniformRefinement& nonUniformRefinement )
{
RigCaseCellResultsData* cellResultsData = eclipseCase->results( RiaDefines::PorosityModelType::MATRIX_MODEL );
RigActiveCellInfo* activeCells = cellResultsData->activeCellInfo();
RigMainGrid* mainGrid = eclipseCase->mainGrid();
cvf::Vec3st max = maxIn;
if ( max.isUndefined() )
{
max = cvf::Vec3st( mainGrid->cellCountI() - 1, mainGrid->cellCountJ() - 1, mainGrid->cellCountK() - 1 );
}
auto allResultAddresses = cellResultsData->existingResults();
auto findResultAddress = [&allResultAddresses]( const QString& keyword ) -> RigEclipseResultAddress
{
for ( auto it = allResultAddresses.rbegin(); it != allResultAddresses.rend(); ++it )
{
if ( it->resultName() == keyword ) return *it;
}
return {};
};
RigEclipseResultAddress resAddr = findResultAddress( keyword );
if ( !cellResultsData->hasResultEntry( resAddr ) )
return std::unexpected( "Keyword '" + keyword.toStdString() + "' not found in results" );
if ( !cellResultsData->ensureKnownResultLoaded( resAddr ) )
return std::unexpected( "Keyword '" + keyword.toStdString() + "' result not loaded." );
if ( cellResultsData->cellScalarResults( resAddr ).empty() )
return std::unexpected( "Keyword '" + keyword.toStdString() + "' result empty." );
std::vector<double> resultValues = cellResultsData->cellScalarResults( resAddr )[0];
if ( resultValues.empty() ) return std::unexpected( "Keyword '" + keyword.toStdString() + "' has no data" );
double defaultExportValue = RiaResultNames::isCategoryResult( keyword ) ? 1.0 : 0.0;
RiaDefines::PorosityModelType porosityModel = RiaDefines::PorosityModelType::MATRIX_MODEL;
auto resultAcc = RigResultAccessorFactory::createFromResultAddress( eclipseCase, 0, porosityModel, 0, resAddr );
std::vector<double> filteredResults;
size_t sectorNi = max.x() - min.x() + 1;
size_t sectorNj = max.y() - min.y() + 1;
size_t sectorNk = max.z() - min.z() + 1;
for ( size_t origK = 0; origK < sectorNk; ++origK )
{
size_t mainK = min.z() + origK;
size_t subCountK = nonUniformRefinement.subcellCount( RigNonUniformRefinement::DimK, origK );
for ( size_t rk = 0; rk < subCountK; ++rk )
{
for ( size_t origJ = 0; origJ < sectorNj; ++origJ )
{
size_t mainJ = min.y() + origJ;
size_t subCountJ = nonUniformRefinement.subcellCount( RigNonUniformRefinement::DimJ, origJ );
for ( size_t rj = 0; rj < subCountJ; ++rj )
{
for ( size_t origI = 0; origI < sectorNi; ++origI )
{
size_t mainI = min.x() + origI;
size_t reservoirCellIndex = mainGrid->cellIndexFromIJK( mainI, mainJ, mainK );
size_t resIndex = activeCells->cellResultIndex( reservoirCellIndex );
size_t subCountI = nonUniformRefinement.subcellCount( RigNonUniformRefinement::DimI, origI );
double value = ( resIndex != cvf::UNDEFINED_SIZE_T ) ? resultAcc->cellScalarGlobIdx( reservoirCellIndex )
: defaultExportValue;
for ( size_t ri = 0; ri < subCountI; ++ri )
{
filteredResults.push_back( value );
}
}
}
}
}
}
return filteredResults;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
@@ -25,6 +25,7 @@
#include "RiaDefines.h"
#include "RigFault.h"
#include "RigNonUniformRefinement.h"
#include "ert/ecl/ecl_kw.h"
@@ -82,6 +83,12 @@ public:
const cvf::Vec3st& max = cvf::Vec3st::UNDEFINED,
const cvf::Vec3st& refinement = cvf::Vec3st( 1, 1, 1 ) );
static std::expected<std::vector<double>, std::string> extractKeywordData( RigEclipseCaseData* eclipseCase,
const QString& keyword,
const cvf::Vec3st& min,
const cvf::Vec3st& max,
const RigNonUniformRefinement& nonUniformRefinement );
static void saveFault( QString completeFilename,
const RigMainGrid* mainGrid,
const std::vector<RigFault::FaultFace>& faultFaces,
@@ -88,6 +88,7 @@ set(SOURCE_GROUP_HEADER_FILES
${CMAKE_CURRENT_LIST_DIR}/RigDoglegTools.h
${CMAKE_CURRENT_LIST_DIR}/RigResdataGridConverter.h
${CMAKE_CURRENT_LIST_DIR}/RigGridExportAdapter.h
${CMAKE_CURRENT_LIST_DIR}/RigNonUniformRefinement.h
)
set(SOURCE_GROUP_SOURCE_FILES
@@ -175,6 +176,7 @@ set(SOURCE_GROUP_SOURCE_FILES
${CMAKE_CURRENT_LIST_DIR}/RigDoglegTools.cpp
${CMAKE_CURRENT_LIST_DIR}/RigResdataGridConverter.cpp
${CMAKE_CURRENT_LIST_DIR}/RigGridExportAdapter.cpp
${CMAKE_CURRENT_LIST_DIR}/RigNonUniformRefinement.cpp
)
list(APPEND CODE_HEADER_FILES ${SOURCE_GROUP_HEADER_FILES})
@@ -26,6 +26,7 @@
#include "RigCell.h"
#include "RigEclipseCaseData.h"
#include "RigMainGrid.h"
#include "RigNonUniformRefinement.h"
#include "cafVecIjk.h"
#include "cvfAssert.h"
@@ -64,6 +65,35 @@ static std::array<cvf::Vec3d, 8> generateRefinedCellCorners( const std::array<cv
return refinedCorners;
}
//--------------------------------------------------------------------------------------------------
/// Generate refined cell corners using non-uniform cumulative fractions
//--------------------------------------------------------------------------------------------------
static std::array<cvf::Vec3d, 8> generateRefinedCellCornersNonUniform( const std::array<cvf::Vec3d, 8>& originalCorners,
const std::vector<double>& cumulativeFractionsX,
const std::vector<double>& cumulativeFractionsY,
const std::vector<double>& cumulativeFractionsZ,
size_t subI,
size_t subJ,
size_t subK )
{
auto allRefinedCorners =
RiaCellDividingTools::createHexCornerCoords( originalCorners, cumulativeFractionsX, cumulativeFractionsY, cumulativeFractionsZ );
size_t nx = cumulativeFractionsX.size();
size_t ny = cumulativeFractionsY.size();
size_t subcellIndex = subK * ( nx * ny ) + subJ * nx + subI;
size_t cornerStartIndex = subcellIndex * 8;
std::array<cvf::Vec3d, 8> refinedCorners;
for ( size_t i = 0; i < 8; ++i )
{
refinedCorners[i] = allRefinedCorners[cornerStartIndex + i];
}
return refinedCorners;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
@@ -99,6 +129,46 @@ RigGridExportAdapter::RigGridExportAdapter( RigEclipseCaseData* eclipseCase,
m_refinedNI = ( maxActual.x() - min.x() + 1 ) * refinement.x();
m_refinedNJ = ( maxActual.y() - min.y() + 1 ) * refinement.y();
m_refinedNK = ( maxActual.z() - min.z() + 1 ) * refinement.z();
// Create equivalent non-uniform refinement for internal use
cvf::Vec3st sectorSize( maxActual.x() - min.x() + 1, maxActual.y() - min.y() + 1, maxActual.z() - min.z() + 1 );
m_nonUniformRefinement = RigNonUniformRefinement::fromUniform( refinement, sectorSize );
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
RigGridExportAdapter::RigGridExportAdapter( RigEclipseCaseData* eclipseCase,
const cvf::Vec3st& min,
const cvf::Vec3st& max,
const RigNonUniformRefinement& nonUniformRefinement,
const cvf::UByteArray* cellVisibilityOverrideForActnum )
: m_mainGrid( nullptr )
, m_activeCellInfo( nullptr )
, m_cellVisibilityOverride( cellVisibilityOverrideForActnum )
, m_min( min )
, m_max( max )
, m_refinement( 1, 1, 1 )
, m_nonUniformRefinement( nonUniformRefinement )
, m_refinedNI( 0 )
, m_refinedNJ( 0 )
, m_refinedNK( 0 )
{
CVF_ASSERT( eclipseCase );
m_mainGrid = eclipseCase->mainGrid();
m_activeCellInfo = eclipseCase->activeCellInfo( RiaDefines::PorosityModelType::MATRIX_MODEL );
CVF_ASSERT( m_mainGrid );
CVF_ASSERT( m_activeCellInfo );
cvf::Vec3st maxActual =
max.isUndefined() ? cvf::Vec3st( m_mainGrid->cellCountI() - 1, m_mainGrid->cellCountJ() - 1, m_mainGrid->cellCountK() - 1 ) : max;
m_max = maxActual;
m_refinedNI = m_nonUniformRefinement.totalRefinedCount( RigNonUniformRefinement::DimI );
m_refinedNJ = m_nonUniformRefinement.totalRefinedCount( RigNonUniformRefinement::DimJ );
m_refinedNK = m_nonUniformRefinement.totalRefinedCount( RigNonUniformRefinement::DimK );
}
//--------------------------------------------------------------------------------------------------
@@ -157,14 +227,33 @@ std::array<cvf::Vec3d, 4> RigGridExportAdapter::getFaceCorners( size_t i, size_t
// Apply coordinate transformations if needed
applyCoordinateTransformation( originalCorners );
// Generate refined cell corners using ResInsight's refinement algorithm
auto refinedCorners = generateRefinedCellCorners( originalCorners,
m_refinement.x(),
m_refinement.y(),
m_refinement.z(),
mapping.subI,
mapping.subJ,
mapping.subK );
// Generate refined cell corners using the appropriate refinement algorithm
std::array<cvf::Vec3d, 8> refinedCorners;
if ( m_nonUniformRefinement.hasNonUniformRefinement() )
{
size_t sectorI = mapping.originalI - m_min.x();
size_t sectorJ = mapping.originalJ - m_min.y();
size_t sectorK = mapping.originalK - m_min.z();
refinedCorners =
generateRefinedCellCornersNonUniform( originalCorners,
m_nonUniformRefinement.cumulativeFractions( RigNonUniformRefinement::DimI, sectorI ),
m_nonUniformRefinement.cumulativeFractions( RigNonUniformRefinement::DimJ, sectorJ ),
m_nonUniformRefinement.cumulativeFractions( RigNonUniformRefinement::DimK, sectorK ),
mapping.subI,
mapping.subJ,
mapping.subK );
}
else
{
refinedCorners = generateRefinedCellCorners( originalCorners,
m_refinement.x(),
m_refinement.y(),
m_refinement.z(),
mapping.subI,
mapping.subJ,
mapping.subK );
}
// Extract the 4 corners for the requested face from the refined corners
std::array<cvf::Vec3d, 4> faceCorners;
@@ -299,10 +388,23 @@ std::array<cvf::Vec3d, 8> RigGridExportAdapter::getOriginalCellCorners( size_t o
std::array<cvf::Vec3d, 8>
RigGridExportAdapter::getRefinedCellCorners( size_t origI, size_t origJ, size_t origK, size_t subI, size_t subJ, size_t subK ) const
{
// Get original cell corners first
std::array<cvf::Vec3d, 8> originalCorners = getOriginalCellCorners( origI, origJ, origK );
// Generate refined subcell corners
if ( m_nonUniformRefinement.hasNonUniformRefinement() )
{
size_t sectorI = origI - m_min.x();
size_t sectorJ = origJ - m_min.y();
size_t sectorK = origK - m_min.z();
return generateRefinedCellCornersNonUniform( originalCorners,
m_nonUniformRefinement.cumulativeFractions( RigNonUniformRefinement::DimI, sectorI ),
m_nonUniformRefinement.cumulativeFractions( RigNonUniformRefinement::DimJ, sectorJ ),
m_nonUniformRefinement.cumulativeFractions( RigNonUniformRefinement::DimK, sectorK ),
subI,
subJ,
subK );
}
return generateRefinedCellCorners( originalCorners, m_refinement.x(), m_refinement.y(), m_refinement.z(), subI, subJ, subK );
}
@@ -343,15 +445,29 @@ RigGridExportAdapter::CellMapping RigGridExportAdapter::mapRefinedToOriginal( si
{
CellMapping mapping;
// Calculate which original cell this refined cell belongs to
mapping.originalI = m_min.x() + refinedI / m_refinement.x();
mapping.originalJ = m_min.y() + refinedJ / m_refinement.y();
mapping.originalK = m_min.z() + refinedK / m_refinement.z();
if ( m_nonUniformRefinement.hasNonUniformRefinement() )
{
auto [origI, subI] = m_nonUniformRefinement.mapRefinedToOriginal( RigNonUniformRefinement::DimI, refinedI );
auto [origJ, subJ] = m_nonUniformRefinement.mapRefinedToOriginal( RigNonUniformRefinement::DimJ, refinedJ );
auto [origK, subK] = m_nonUniformRefinement.mapRefinedToOriginal( RigNonUniformRefinement::DimK, refinedK );
// Calculate subcell indices within the original cell
mapping.subI = refinedI % m_refinement.x();
mapping.subJ = refinedJ % m_refinement.y();
mapping.subK = refinedK % m_refinement.z();
mapping.originalI = m_min.x() + origI;
mapping.originalJ = m_min.y() + origJ;
mapping.originalK = m_min.z() + origK;
mapping.subI = subI;
mapping.subJ = subJ;
mapping.subK = subK;
}
else
{
mapping.originalI = m_min.x() + refinedI / m_refinement.x();
mapping.originalJ = m_min.y() + refinedJ / m_refinement.y();
mapping.originalK = m_min.z() + refinedK / m_refinement.z();
mapping.subI = refinedI % m_refinement.x();
mapping.subJ = refinedJ % m_refinement.y();
mapping.subK = refinedK % m_refinement.z();
}
return mapping;
}
@@ -385,7 +501,7 @@ cvf::Vec3st RigGridExportAdapter::refinement() const
//--------------------------------------------------------------------------------------------------
bool RigGridExportAdapter::hasRefinement() const
{
return m_refinement.x() > 1 || m_refinement.y() > 1 || m_refinement.z() > 1;
return m_refinement.x() > 1 || m_refinement.y() > 1 || m_refinement.z() > 1 || m_nonUniformRefinement.hasNonUniformRefinement();
}
//--------------------------------------------------------------------------------------------------
@@ -467,3 +583,57 @@ std::expected<caf::VecIjk0, QString> RigGridExportAdapter::transformIjkToSectorC
( ijk.z() - min.z() ) * refinement.z() );
}
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
const RigNonUniformRefinement& RigGridExportAdapter::nonUniformRefinement() const
{
return m_nonUniformRefinement;
}
//--------------------------------------------------------------------------------------------------
/// Transform IJK coordinates from global grid space to sector-relative space with non-uniform refinement
//--------------------------------------------------------------------------------------------------
std::expected<caf::VecIjk0, QString> RigGridExportAdapter::transformIjkToSectorCoordinates( const caf::VecIjk0& ijk,
const caf::VecIjk0& min,
const caf::VecIjk0& max,
const RigNonUniformRefinement& nonUniformRefinement,
bool applyRefinementCentering,
bool isBoxMaxCoordinate )
{
RigBoundingBoxIjk<caf::VecIjk0> sectorBox( min, max );
if ( !sectorBox.contains( ijk ) )
{
return std::unexpected( QString( "IJK coordinates (%1) are outside sector bounds [(%2), (%3)]" )
.arg( QString::fromStdString( ijk.toString() ) )
.arg( QString::fromStdString( min.toString() ) )
.arg( QString::fromStdString( max.toString() ) ) );
}
size_t sectorI = ijk.x() - min.x();
size_t sectorJ = ijk.y() - min.y();
size_t sectorK = ijk.z() - min.z();
if ( applyRefinementCentering )
{
return caf::VecIjk0( nonUniformRefinement.cumulativeOffset( RigNonUniformRefinement::DimI, sectorI ) +
nonUniformRefinement.subcellCount( RigNonUniformRefinement::DimI, sectorI ) / 2,
nonUniformRefinement.cumulativeOffset( RigNonUniformRefinement::DimJ, sectorJ ) +
nonUniformRefinement.subcellCount( RigNonUniformRefinement::DimJ, sectorJ ) / 2,
nonUniformRefinement.cumulativeOffset( RigNonUniformRefinement::DimK, sectorK ) +
nonUniformRefinement.subcellCount( RigNonUniformRefinement::DimK, sectorK ) / 2 );
}
else if ( isBoxMaxCoordinate )
{
return caf::VecIjk0( nonUniformRefinement.cumulativeOffset( RigNonUniformRefinement::DimI, sectorI + 1 ) - 1,
nonUniformRefinement.cumulativeOffset( RigNonUniformRefinement::DimJ, sectorJ + 1 ) - 1,
nonUniformRefinement.cumulativeOffset( RigNonUniformRefinement::DimK, sectorK + 1 ) - 1 );
}
else
{
return caf::VecIjk0( nonUniformRefinement.cumulativeOffset( RigNonUniformRefinement::DimI, sectorI ),
nonUniformRefinement.cumulativeOffset( RigNonUniformRefinement::DimJ, sectorJ ),
nonUniformRefinement.cumulativeOffset( RigNonUniformRefinement::DimK, sectorK ) );
}
}
@@ -18,6 +18,8 @@
#pragma once
#include "RigNonUniformRefinement.h"
#include "cafVecIjk.h"
#include "cvfArray.h"
@@ -66,6 +68,12 @@ public:
const cvf::Vec3st& refinement,
const cvf::UByteArray* cellVisibilityOverrideForActnum = nullptr );
RigGridExportAdapter( RigEclipseCaseData* eclipseCase,
const cvf::Vec3st& min,
const cvf::Vec3st& max,
const RigNonUniformRefinement& nonUniformRefinement,
const cvf::UByteArray* cellVisibilityOverrideForActnum = nullptr );
// Grid dimensions (after refinement)
size_t cellCountI() const;
size_t cellCountJ() const;
@@ -88,6 +96,9 @@ public:
cvf::Vec3st refinement() const;
bool hasRefinement() const;
// Non-uniform refinement accessor
const RigNonUniformRefinement& nonUniformRefinement() const;
// Coordinate transformation utilities
static std::expected<caf::VecIjk0, QString> transformIjkToSectorCoordinates( const caf::VecIjk0& originalIjk,
const caf::VecIjk0& min,
@@ -96,6 +107,13 @@ public:
bool applyRefinementCentering = false,
bool isBoxMaxCoordinate = false );
static std::expected<caf::VecIjk0, QString> transformIjkToSectorCoordinates( const caf::VecIjk0& originalIjk,
const caf::VecIjk0& min,
const caf::VecIjk0& max,
const RigNonUniformRefinement& nonUniformRefinement,
bool applyRefinementCentering = false,
bool isBoxMaxCoordinate = false );
private:
// Internal methods to handle original vs refined cell access
std::array<cvf::Vec3d, 8> getOriginalCellCorners( size_t origI, size_t origJ, size_t origK ) const;
@@ -119,6 +137,8 @@ private:
cvf::Vec3st m_max;
cvf::Vec3st m_refinement;
RigNonUniformRefinement m_nonUniformRefinement;
// Cached refined grid dimensions
size_t m_refinedNI;
size_t m_refinedNJ;
@@ -0,0 +1,248 @@
/////////////////////////////////////////////////////////////////////////////////
//
// Copyright (C) 2026- 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 "RigNonUniformRefinement.h"
#include <algorithm>
#include <numeric>
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
RigNonUniformRefinement::RigNonUniformRefinement()
: m_sectorSize( 0, 0, 0 )
{
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
RigNonUniformRefinement::RigNonUniformRefinement( const cvf::Vec3st& sectorSize )
: m_sectorSize( sectorSize )
{
// Initialize each cell with default {1.0} (no refinement)
for ( auto dim : { DimI, DimJ, DimK } )
{
m_fractions[dim].resize( m_sectorSize[dim], std::vector<double>{ 1.0 } );
rebuildOffsets( dim );
}
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
void RigNonUniformRefinement::setCumulativeFractions( Dimension dim, size_t origIndex, const std::vector<double>& cumulativeFractions )
{
if ( origIndex >= m_sectorSize[dim] ) return;
m_fractions[dim][origIndex] = cumulativeFractions;
rebuildOffsets( dim );
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
std::vector<double> RigNonUniformRefinement::widthsToCumulativeFractions( const std::vector<double>& widths )
{
if ( widths.empty() ) return { 1.0 };
double sum = std::accumulate( widths.begin(), widths.end(), 0.0 );
if ( sum <= 0.0 ) return { 1.0 };
std::vector<double> cumulative;
cumulative.reserve( widths.size() );
double runningSum = 0.0;
for ( size_t i = 0; i < widths.size(); ++i )
{
runningSum += widths[i] / sum;
cumulative.push_back( runningSum );
}
// Ensure the last value is exactly 1.0
cumulative.back() = 1.0;
return cumulative;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
size_t RigNonUniformRefinement::subcellCount( Dimension dim, size_t origIndex ) const
{
if ( origIndex >= m_sectorSize[dim] ) return 1;
return m_fractions[dim][origIndex].size();
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
size_t RigNonUniformRefinement::cumulativeOffset( Dimension dim, size_t origIndex ) const
{
// Allow origIndex == m_sectorSize[dim] to get the total count (end sentinel)
if ( origIndex > m_sectorSize[dim] ) return 0;
return m_offsets[dim][origIndex];
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
size_t RigNonUniformRefinement::totalRefinedCount( Dimension dim ) const
{
if ( m_sectorSize[dim] == 0 ) return 0;
return m_offsets[dim][m_sectorSize[dim]];
}
//--------------------------------------------------------------------------------------------------
/// Maps a refined index back to {originalIndex, subIndex} using binary search on the offset table
//--------------------------------------------------------------------------------------------------
std::pair<size_t, size_t> RigNonUniformRefinement::mapRefinedToOriginal( Dimension dim, size_t refinedIndex ) const
{
if ( m_sectorSize[dim] == 0 ) return { 0, 0 };
// Binary search: find the largest origIndex such that m_offsets[dim][origIndex] <= refinedIndex
// m_offsets[dim] has size m_sectorSize[dim]+1
// We want upper_bound - 1
auto it = std::upper_bound( m_offsets[dim].begin(), m_offsets[dim].end(), refinedIndex );
size_t origIndex = static_cast<size_t>( std::distance( m_offsets[dim].begin(), it ) ) - 1;
size_t subIndex = refinedIndex - m_offsets[dim][origIndex];
return { origIndex, subIndex };
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
const std::vector<double>& RigNonUniformRefinement::cumulativeFractions( Dimension dim, size_t origIndex ) const
{
static const std::vector<double> defaultFraction = { 1.0 };
if ( origIndex >= m_sectorSize[dim] ) return defaultFraction;
return m_fractions[dim][origIndex];
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
size_t RigNonUniformRefinement::sectorSize( Dimension dim ) const
{
return m_sectorSize[dim];
}
//--------------------------------------------------------------------------------------------------
/// Distribute fractional widths across a range of cells.
///
/// The widths define subdivision boundaries in the combined range of cells [startIndex, endIndex].
/// Each original cell occupies 1/numCells of the global range (in index space).
/// Global subdivision boundaries that fall within a cell are converted to cell-local fractions.
///
/// Example: 3 cells, widths [0.2, 0.3, 0.5] -> cumulative [0.2, 0.5, 1.0]
/// Cell 0 [0, 0.333]: boundary 0.2 falls here -> local fracs [0.6, 1.0] -> 2 subcells
/// Cell 1 [0.333, 0.667]: boundary 0.5 falls here -> local fracs [0.5, 1.0] -> 2 subcells
/// Cell 2 [0.667, 1.0]: no interior boundaries -> local fracs [1.0] -> 1 subcell
//--------------------------------------------------------------------------------------------------
void RigNonUniformRefinement::distributeWidthsAcrossCells( Dimension dim, size_t startIndex, size_t endIndex, const std::vector<double>& widths )
{
if ( startIndex > endIndex || endIndex >= m_sectorSize[dim] || widths.empty() ) return;
// Convert widths to global cumulative fractions
auto globalFractions = widthsToCumulativeFractions( widths );
size_t numCells = endIndex - startIndex + 1;
for ( size_t c = 0; c < numCells; ++c )
{
double cellStart = static_cast<double>( c ) / static_cast<double>( numCells );
double cellEnd = static_cast<double>( c + 1 ) / static_cast<double>( numCells );
double cellWidth = cellEnd - cellStart;
std::vector<double> localFractions;
// Find global boundaries that fall strictly within this cell
for ( double gf : globalFractions )
{
if ( gf > cellStart && gf < cellEnd )
{
// Convert to cell-local fraction
localFractions.push_back( ( gf - cellStart ) / cellWidth );
}
}
// Always add 1.0 as the final boundary (the cell end)
localFractions.push_back( 1.0 );
setCumulativeFractions( dim, startIndex + c, localFractions );
}
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
bool RigNonUniformRefinement::hasNonUniformRefinement() const
{
for ( auto dim : { DimI, DimJ, DimK } )
{
if ( totalRefinedCount( dim ) != m_sectorSize[dim] ) return true;
}
return false;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
RigNonUniformRefinement RigNonUniformRefinement::fromUniform( const cvf::Vec3st& refinement, const cvf::Vec3st& sectorSize )
{
RigNonUniformRefinement result( sectorSize );
size_t refCount[3] = { refinement.x(), refinement.y(), refinement.z() };
for ( auto dim : { DimI, DimJ, DimK } )
{
if ( refCount[dim] <= 1 ) continue;
// Create uniform cumulative fractions
std::vector<double> uniformFractions;
uniformFractions.reserve( refCount[dim] );
for ( size_t i = 1; i <= refCount[dim]; ++i )
{
uniformFractions.push_back( static_cast<double>( i ) / static_cast<double>( refCount[dim] ) );
}
uniformFractions.back() = 1.0;
for ( size_t idx = 0; idx < result.m_sectorSize[dim]; ++idx )
{
result.m_fractions[dim][idx] = uniformFractions;
}
result.rebuildOffsets( dim );
}
return result;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
void RigNonUniformRefinement::rebuildOffsets( Dimension dim )
{
m_offsets[dim].resize( m_sectorSize[dim] + 1 );
m_offsets[dim][0] = 0;
for ( size_t i = 0; i < m_sectorSize[dim]; ++i )
{
m_offsets[dim][i + 1] = m_offsets[dim][i] + m_fractions[dim][i].size();
}
}
@@ -0,0 +1,93 @@
/////////////////////////////////////////////////////////////////////////////////
//
// Copyright (C) 2026- 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.
//
/////////////////////////////////////////////////////////////////////////////////
#pragma once
#include "cvfVector3.h"
#include <cstddef>
#include <utility>
#include <vector>
//==================================================================================================
//
// Non-uniform per-dimension refinement specification for sector model export.
//
// Since refinement "crosses over" the model (all cells at a given I index share the same
// I-refinement), each dimension is independent.
//
// For each dimension: a vector-of-vectors indexed by original cell index (relative to sector),
// containing cumulative fractions in (0, 1]. Cells with no explicit refinement default to {1.0}
// (1 subcell = no refinement).
//
// A cumulative offset table per dimension enables O(1) original->refined lookup
// and O(log n) reverse lookup via binary search.
//
//==================================================================================================
class RigNonUniformRefinement
{
public:
enum Dimension : size_t
{
DimI = 0,
DimJ = 1,
DimK = 2
};
RigNonUniformRefinement();
explicit RigNonUniformRefinement( const cvf::Vec3st& sectorSize );
// Set cumulative fractions for a specific cell index in a dimension.
// The fractions must be strictly increasing and end at 1.0.
void setCumulativeFractions( Dimension dim, size_t origIndex, const std::vector<double>& cumulativeFractions );
// Set fractional widths for a specific cell index in a dimension.
// Widths are normalized so they sum to 1.0, then converted to cumulative fractions.
static std::vector<double> widthsToCumulativeFractions( const std::vector<double>& widths );
// Query methods
size_t subcellCount( Dimension dim, size_t origIndex ) const;
size_t cumulativeOffset( Dimension dim, size_t origIndex ) const;
size_t totalRefinedCount( Dimension dim ) const;
std::pair<size_t, size_t> mapRefinedToOriginal( Dimension dim, size_t refinedIndex ) const;
const std::vector<double>& cumulativeFractions( Dimension dim, size_t origIndex ) const;
size_t sectorSize( Dimension dim ) const;
bool hasNonUniformRefinement() const;
// Distribute fractional widths across a range of cells in a given dimension.
// The widths span the combined range [startIndex, endIndex] — subdivision boundaries
// are placed in the global fraction space and then assigned to individual cells.
// Each cell in the range occupies an equal share (1/numCells) of the global range.
void distributeWidthsAcrossCells( Dimension dim, size_t startIndex, size_t endIndex, const std::vector<double>& widths );
// Create a non-uniform refinement that is equivalent to uniform refinement
static RigNonUniformRefinement fromUniform( const cvf::Vec3st& refinement, const cvf::Vec3st& sectorSize );
private:
void rebuildOffsets( Dimension dim );
// Per-dimension data: m_fractions[dim][origIndex] = cumulative fractions
std::vector<std::vector<double>> m_fractions[3];
// Per-dimension cumulative offset table: m_offsets[dim][n] = total refined cells for cells 0..n-1
std::vector<size_t> m_offsets[3];
// Sector dimensions (number of original cells per dimension)
cvf::Vec3st m_sectorSize;
};
@@ -93,6 +93,30 @@ void RigSimulationInputSettings::setRefinement( const cvf::Vec3st& refinement )
m_refinement = refinement;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
const RigNonUniformRefinement& RigSimulationInputSettings::nonUniformRefinement() const
{
return m_nonUniformRefinement;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
void RigSimulationInputSettings::setNonUniformRefinement( const RigNonUniformRefinement& refinement )
{
m_nonUniformRefinement = refinement;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
bool RigSimulationInputSettings::hasNonUniformRefinement() const
{
return m_nonUniformRefinement.hasNonUniformRefinement();
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
@@ -20,6 +20,7 @@
#include "RiaModelExportDefines.h"
#include "RigModelPaddingSettings.h"
#include "RigNonUniformRefinement.h"
#include "cafAppEnum.h"
#include "cafVecIjk.h"
@@ -60,6 +61,11 @@ public:
cvf::Vec3st refinement() const;
void setRefinement( const cvf::Vec3st& refinement );
// Non-uniform refinement
const RigNonUniformRefinement& nonUniformRefinement() const;
void setNonUniformRefinement( const RigNonUniformRefinement& refinement );
bool hasNonUniformRefinement() const;
// Keyword handling
const std::vector<std::string>& keywordsToRemove() const;
void setKeywordsToRemove( const std::vector<std::string>& keywords );
@@ -89,6 +95,7 @@ private:
caf::VecIjk0 m_min;
caf::VecIjk0 m_max;
cvf::Vec3st m_refinement;
RigNonUniformRefinement m_nonUniformRefinement;
std::vector<std::string> m_keywordsToRemove;
std::vector<Opm::DeckRecord> m_bcpropKeywords;
BoundaryConditionEnum m_boundaryCondition;
@@ -34,6 +34,7 @@
#include "RigGridExportAdapter.h"
#include "RigMainGrid.h"
#include "RigModelPaddingSettings.h"
#include "RigNonUniformRefinement.h"
#include "RigPadModel.h"
#include "RigResdataGridConverter.h"
#include "RigSimulationInputSettings.h"
@@ -193,8 +194,16 @@ std::expected<void, QString> RigSimulationInputTool::updateCornerPointGridInDeck
RifOpmFlowDeckFile& deckFile )
{
// Get grid bounds for extraction
auto createGridAdapter = [&]()
{
if ( settings.hasNonUniformRefinement() )
{
return RigGridExportAdapter( eclipseCase->eclipseCaseData(), settings.min(), settings.max(), settings.nonUniformRefinement() );
}
return RigGridExportAdapter( eclipseCase->eclipseCaseData(), settings.min(), settings.max(), settings.refinement() );
};
RigGridExportAdapter gridAdapter( eclipseCase->eclipseCaseData(), settings.min(), settings.max(), settings.refinement() );
RigGridExportAdapter gridAdapter = createGridAdapter();
std::vector<float> coordArray;
std::vector<float> zcornArray;
@@ -293,11 +302,16 @@ std::expected<void, QString> RigSimulationInputTool::replaceKeywordValuesInDeckF
}
// Try to extract keyword data
auto result = RifEclipseInputFileTools::extractKeywordData( eclipseCase->eclipseCaseData(),
keyword,
settings.min(),
settings.max(),
settings.refinement() );
auto result = settings.hasNonUniformRefinement() ? RifEclipseInputFileTools::extractKeywordData( eclipseCase->eclipseCaseData(),
keyword,
settings.min(),
settings.max(),
settings.nonUniformRefinement() )
: RifEclipseInputFileTools::extractKeywordData( eclipseCase->eclipseCaseData(),
keyword,
settings.min(),
settings.max(),
settings.refinement() );
if ( result )
{
// Replace keyword values in deck with extracted data
@@ -337,6 +351,14 @@ std::set<std::string> RigSimulationInputTool::cropDataKeywordsInDeckFile( RimEcl
const auto max = settings.max();
const auto refinement = settings.refinement();
// Build non-uniform refinement for cropping (handles both uniform and non-uniform cases)
const bool hasNonUniform = settings.hasNonUniformRefinement();
const auto nonUniformRef =
hasNonUniform
? settings.nonUniformRefinement()
: RigNonUniformRefinement::fromUniform( refinement,
cvf::Vec3st( max.x() - min.x() + 1, max.y() - min.y() + 1, max.z() - min.z() + 1 ) );
auto keywords = deckFile.keywords( false );
for ( const auto& name : keywords )
@@ -357,20 +379,35 @@ std::set<std::string> RigSimulationInputTool::cropDataKeywordsInDeckFile( RimEcl
auto cropAndReplace = [&]<typename T>( const std::vector<T>& fullData )
{
std::vector<T> croppedData;
cvf::Vec3st refinedMin( min.x() * refinement.x(), min.y() * refinement.y(), min.z() * refinement.z() );
cvf::Vec3st refinedMax( ( max.x() + 1 ) * refinement.x(), ( max.y() + 1 ) * refinement.y(), ( max.z() + 1 ) * refinement.z() );
for ( size_t k = refinedMin.z(); k < refinedMax.z(); ++k )
size_t sectorNk = max.z() - min.z() + 1;
size_t sectorNj = max.y() - min.y() + 1;
size_t sectorNi = max.x() - min.x() + 1;
for ( size_t origK = 0; origK < sectorNk; ++origK )
{
size_t mainK = k / refinement.z();
for ( size_t j = refinedMin.y(); j < refinedMax.y(); ++j )
size_t mainK = min.z() + origK;
size_t subCountK = nonUniformRef.subcellCount( RigNonUniformRefinement::DimK, origK );
for ( size_t rk = 0; rk < subCountK; ++rk )
{
size_t mainJ = j / refinement.y();
for ( size_t i = refinedMin.x(); i < refinedMax.x(); ++i )
for ( size_t origJ = 0; origJ < sectorNj; ++origJ )
{
size_t mainI = i / refinement.x();
size_t sourceIdx = mainI + mainJ * fullNx + mainK * fullNx * fullNy;
croppedData.push_back( fullData[sourceIdx] );
size_t mainJ = min.y() + origJ;
size_t subCountJ = nonUniformRef.subcellCount( RigNonUniformRefinement::DimJ, origJ );
for ( size_t rj = 0; rj < subCountJ; ++rj )
{
for ( size_t origI = 0; origI < sectorNi; ++origI )
{
size_t mainI = min.x() + origI;
size_t sourceIdx = mainI + mainJ * fullNx + mainK * fullNx * fullNy;
size_t subCountI = nonUniformRef.subcellCount( RigNonUniformRefinement::DimI, origI );
for ( size_t ri = 0; ri < subCountI; ++ri )
{
croppedData.push_back( fullData[sourceIdx] );
}
}
}
}
}
}
@@ -419,7 +456,19 @@ std::expected<void, QString> RigSimulationInputTool::addBorderBoundaryConditions
}
// Create grid adapter for refined grid operations
RigGridExportAdapter gridAdapter( eclipseCase->eclipseCaseData(), settings.min(), settings.max(), settings.refinement(), visibility.p() );
auto createBorderGridAdapter = [&]()
{
if ( settings.hasNonUniformRefinement() )
{
return RigGridExportAdapter( eclipseCase->eclipseCaseData(),
settings.min(),
settings.max(),
settings.nonUniformRefinement(),
visibility.p() );
}
return RigGridExportAdapter( eclipseCase->eclipseCaseData(), settings.min(), settings.max(), settings.refinement(), visibility.p() );
};
RigGridExportAdapter gridAdapter = createBorderGridAdapter();
// Create refined visibility array matching refined grid dimensions
std::vector<int> refinedVisibility = createRefinedVisibility( gridAdapter );
@@ -657,13 +706,14 @@ static std::expected<Modification, QString> processBoxKeywordEntry( const Opm::D
/// Crop a data keyword to the intersection of its BOX context and the sector.
/// Returns nullopt if no modification is needed (e.g. keyword has unexpected structure).
//--------------------------------------------------------------------------------------------------
static std::optional<Modification> cropDataKeywordInBoxContext( const Opm::DeckKeyword& kw,
const Opm::FileDeck::Index& index,
const caf::VecIjk0& boxMin,
const caf::VecIjk0& boxMax,
const caf::VecIjk0& sectorMin,
const caf::VecIjk0& sectorMax,
const cvf::Vec3st& refinement )
static std::optional<Modification> cropDataKeywordInBoxContext( const Opm::DeckKeyword& kw,
const Opm::FileDeck::Index& index,
const caf::VecIjk0& boxMin,
const caf::VecIjk0& boxMax,
const caf::VecIjk0& sectorMin,
const caf::VecIjk0& sectorMax,
const cvf::Vec3st& refinement,
const RigNonUniformRefinement* nonUniformRef = nullptr )
{
const std::string& name = kw.name();
@@ -691,11 +741,15 @@ static std::optional<Modification> cropDataKeywordInBoxContext( const Opm::DeckK
for ( size_t k = intMin.z(); k <= static_cast<size_t>( intMax.z() ); ++k )
{
for ( size_t rk = 0; rk < refinement.z(); ++rk )
size_t subCountK = nonUniformRef ? nonUniformRef->subcellCount( RigNonUniformRefinement::DimK, k - sectorMin.z() )
: refinement.z();
for ( size_t rk = 0; rk < subCountK; ++rk )
{
for ( size_t j = intMin.y(); j <= static_cast<size_t>( intMax.y() ); ++j )
{
for ( size_t rj = 0; rj < refinement.y(); ++rj )
size_t subCountJ = nonUniformRef ? nonUniformRef->subcellCount( RigNonUniformRefinement::DimJ, j - sectorMin.y() )
: refinement.y();
for ( size_t rj = 0; rj < subCountJ; ++rj )
{
for ( size_t i = intMin.x(); i <= static_cast<size_t>( intMax.x() ); ++i )
{
@@ -704,7 +758,9 @@ static std::optional<Modification> cropDataKeywordInBoxContext( const Opm::DeckK
size_t boxRelK = k - boxMin.z();
size_t srcIdx = boxRelI + boxRelJ * boxNi + boxRelK * boxNi * boxNj;
for ( size_t ri = 0; ri < refinement.x(); ++ri )
size_t subCountI = nonUniformRef ? nonUniformRef->subcellCount( RigNonUniformRefinement::DimI, i - sectorMin.x() )
: refinement.x();
for ( size_t ri = 0; ri < subCountI; ++ri )
{
croppedData.push_back( sourceData[srcIdx] );
}
@@ -843,6 +899,9 @@ std::expected<void, QString> RigSimulationInputTool::transformKeywordsInDeckFile
const auto sectorMax = settings.max();
const auto refinement = settings.refinement();
// Build non-uniform refinement pointer for functions that support it
const RigNonUniformRefinement* nonUniformRefPtr = settings.hasNonUniformRefinement() ? &settings.nonUniformRefinement() : nullptr;
// Map keyword names to their record processor functions
const std::map<std::string, RecordProcessorFunc> recordProcessors = {
{ "EQUALS", processEqualsRecord },
@@ -878,7 +937,7 @@ std::expected<void, QString> RigSimulationInputTool::transformKeywordsInDeckFile
}
else if ( insideBox && kw.isDataKeyword() )
{
auto result = cropDataKeywordInBoxContext( kw, it, boxMin, boxMax, sectorMin, sectorMax, refinement );
auto result = cropDataKeywordInBoxContext( kw, it, boxMin, boxMax, sectorMin, sectorMax, refinement, nonUniformRefPtr );
if ( result ) modifications.push_back( std::move( result.value() ) );
}
else if ( auto procIt = recordProcessors.find( name ); procIt != recordProcessors.end() )
@@ -1013,7 +1072,28 @@ std::expected<void, QString> RigSimulationInputTool::updateWelldimsKeyword( RimE
if ( wellDims.empty() ) return std::unexpected( QString( "Missing WELLDIMS keyword" ) );
// Scale the max connections with the refinement
int newMaxConnections = wellDims[1] * static_cast<int>( settings.refinement().x() * settings.refinement().y() * settings.refinement().z() );
int refinementProduct = 1;
if ( settings.hasNonUniformRefinement() )
{
const auto& nuRef = settings.nonUniformRefinement();
// Use the maximum subcell product across all cell combinations
size_t maxProduct = 1;
for ( auto dim : { RigNonUniformRefinement::DimI, RigNonUniformRefinement::DimJ, RigNonUniformRefinement::DimK } )
{
size_t maxSubcells = 1;
for ( size_t idx = 0; idx < nuRef.sectorSize( dim ); ++idx )
{
maxSubcells = std::max( maxSubcells, nuRef.subcellCount( dim, idx ) );
}
maxProduct *= maxSubcells;
}
refinementProduct = static_cast<int>( maxProduct );
}
else
{
refinementProduct = static_cast<int>( settings.refinement().x() * settings.refinement().y() * settings.refinement().z() );
}
int newMaxConnections = wellDims[1] * refinementProduct;
deckFile.setWelldims( wellDims[0], newMaxConnections, wellDims[2], wellDims[3] );
@@ -1133,8 +1213,20 @@ std::expected<Opm::DeckRecord, QString> RigSimulationInputTool::processWelspecsR
// Transform with clamped coordinates (this will always succeed)
// Center the coordinate in the refined cell block
size_t sectorI = ( clampedI - settings.min().x() ) * settings.refinement().x() + ( settings.refinement().x() + 1 ) / 2;
size_t sectorJ = ( clampedJ - settings.min().y() ) * settings.refinement().y() + ( settings.refinement().y() + 1 ) / 2;
size_t sectorI, sectorJ;
if ( settings.hasNonUniformRefinement() )
{
const auto& nuRef = settings.nonUniformRefinement();
size_t relI = clampedI - settings.min().x();
size_t relJ = clampedJ - settings.min().y();
sectorI = nuRef.cumulativeOffset( RigNonUniformRefinement::DimI, relI ) + nuRef.subcellCount( RigNonUniformRefinement::DimI, relI ) / 2;
sectorJ = nuRef.cumulativeOffset( RigNonUniformRefinement::DimJ, relJ ) + nuRef.subcellCount( RigNonUniformRefinement::DimJ, relJ ) / 2;
}
else
{
sectorI = ( clampedI - settings.min().x() ) * settings.refinement().x() + ( settings.refinement().x() + 1 ) / 2;
sectorJ = ( clampedJ - settings.min().y() ) * settings.refinement().y() + ( settings.refinement().y() + 1 ) / 2;
}
if ( origI < static_cast<int>( settings.min().x() ) || origI > static_cast<int>( settings.max().x() ) ||
origJ < static_cast<int>( settings.min().y() ) || origJ > static_cast<int>( settings.max().y() ) )
@@ -1190,12 +1282,16 @@ std::expected<Opm::DeckRecord, QString> RigSimulationInputTool::processCompdatRe
// Transform K1
caf::VecIjk0 origIjkK1( origI, origJ, origK1 );
auto transformResultK1 =
RigGridExportAdapter::transformIjkToSectorCoordinates( origIjkK1, settings.min(), settings.max(), settings.refinement(), true );
settings.hasNonUniformRefinement()
? RigGridExportAdapter::transformIjkToSectorCoordinates( origIjkK1, settings.min(), settings.max(), settings.nonUniformRefinement(), true )
: RigGridExportAdapter::transformIjkToSectorCoordinates( origIjkK1, settings.min(), settings.max(), settings.refinement(), true );
// Transform K2
caf::VecIjk0 origIjkK2( origI, origJ, origK2 );
auto transformResultK2 =
RigGridExportAdapter::transformIjkToSectorCoordinates( origIjkK2, settings.min(), settings.max(), settings.refinement(), true );
settings.hasNonUniformRefinement()
? RigGridExportAdapter::transformIjkToSectorCoordinates( origIjkK2, settings.min(), settings.max(), settings.nonUniformRefinement(), true )
: RigGridExportAdapter::transformIjkToSectorCoordinates( origIjkK2, settings.min(), settings.max(), settings.refinement(), true );
if ( !transformResultK1 )
{
@@ -1253,8 +1349,17 @@ std::expected<Opm::DeckRecord, QString> RigSimulationInputTool::processCompsegsR
// Transform I, J, K (first three items)
caf::VecIjk1 origIjk = extractIjk( record, 0, 1, 2 );
auto transformResult =
RigGridExportAdapter::transformIjkToSectorCoordinates( origIjk.toZeroBased(), settings.min(), settings.max(), settings.refinement(), true );
auto transformResult = settings.hasNonUniformRefinement()
? RigGridExportAdapter::transformIjkToSectorCoordinates( origIjk.toZeroBased(),
settings.min(),
settings.max(),
settings.nonUniformRefinement(),
true )
: RigGridExportAdapter::transformIjkToSectorCoordinates( origIjk.toZeroBased(),
settings.min(),
settings.max(),
settings.refinement(),
true );
if ( !transformResult )
{
@@ -2440,11 +2545,12 @@ std::expected<void, QString> RigSimulationInputTool::exportEditNncKeyword( RimEc
QString( "Found %1 internal EDITNNC connections in sector (out of %2 total)" ).arg( sectorConnections.size() ).arg( allConnections.size() ) );
// Step 3: Transform and refine
bool hasRefinement = ( settings.refinement().x() != 1 || settings.refinement().y() != 1 || settings.refinement().z() != 1 );
bool hasRef = ( settings.refinement().x() != 1 || settings.refinement().y() != 1 || settings.refinement().z() != 1 ) ||
settings.hasNonUniformRefinement();
std::vector<TransformedNNCConnection> transformedConnections;
if ( hasRefinement )
if ( hasRef )
{
// Refine to subcells with corresponding positions
for ( const auto& conn : sectorConnections )
@@ -23,6 +23,7 @@ set(SOURCE_UNITTEST_FILES
${CMAKE_CURRENT_LIST_DIR}/RigReservoir-Test.cpp
${CMAKE_CURRENT_LIST_DIR}/RigResdataGridConverter-Test.cpp
${CMAKE_CURRENT_LIST_DIR}/RigGridExportAdapter-Test.cpp
${CMAKE_CURRENT_LIST_DIR}/RigNonUniformRefinement-Test.cpp
${CMAKE_CURRENT_LIST_DIR}/RigStatisticsMath-Test.cpp
${CMAKE_CURRENT_LIST_DIR}/RigWellPathIntersectionTools-Test.cpp
${CMAKE_CURRENT_LIST_DIR}/RimWellLogExtractionCurveImpl-Test.cpp
@@ -0,0 +1,394 @@
/////////////////////////////////////////////////////////////////////////////////
//
// Copyright (C) 2026- 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 "gtest/gtest.h"
#include "RigNonUniformRefinement.h"
#include "RiaCellDividingTools.h"
#include "RigGridExportAdapter.h"
#include "cafVecIjk.h"
//--------------------------------------------------------------------------------------------------
/// Test fromUniform() matches expected uniform counts
//--------------------------------------------------------------------------------------------------
TEST( RigNonUniformRefinement, FromUniform )
{
auto nuRef = RigNonUniformRefinement::fromUniform( cvf::Vec3st( 3, 2, 1 ), cvf::Vec3st( 4, 5, 6 ) );
// Total refined count should match uniform: sectorSize * refinement
EXPECT_EQ( nuRef.totalRefinedCount( RigNonUniformRefinement::DimI ), 4u * 3 );
EXPECT_EQ( nuRef.totalRefinedCount( RigNonUniformRefinement::DimJ ), 5u * 2 );
EXPECT_EQ( nuRef.totalRefinedCount( RigNonUniformRefinement::DimK ), 6u * 1 );
// Each cell in dim 0 should have 3 subcells
for ( size_t i = 0; i < 4; ++i )
{
EXPECT_EQ( nuRef.subcellCount( RigNonUniformRefinement::DimI, i ), 3u );
}
// Each cell in dim 1 should have 2 subcells
for ( size_t j = 0; j < 5; ++j )
{
EXPECT_EQ( nuRef.subcellCount( RigNonUniformRefinement::DimJ, j ), 2u );
}
// Each cell in dim 2 should have 1 subcell
for ( size_t k = 0; k < 6; ++k )
{
EXPECT_EQ( nuRef.subcellCount( RigNonUniformRefinement::DimK, k ), 1u );
}
EXPECT_TRUE( nuRef.hasNonUniformRefinement() );
}
//--------------------------------------------------------------------------------------------------
/// Test no refinement case
//--------------------------------------------------------------------------------------------------
TEST( RigNonUniformRefinement, NoRefinement )
{
auto nuRef = RigNonUniformRefinement::fromUniform( cvf::Vec3st( 1, 1, 1 ), cvf::Vec3st( 3, 3, 3 ) );
EXPECT_EQ( nuRef.totalRefinedCount( RigNonUniformRefinement::DimI ), 3u );
EXPECT_EQ( nuRef.totalRefinedCount( RigNonUniformRefinement::DimJ ), 3u );
EXPECT_EQ( nuRef.totalRefinedCount( RigNonUniformRefinement::DimK ), 3u );
EXPECT_FALSE( nuRef.hasNonUniformRefinement() );
}
//--------------------------------------------------------------------------------------------------
/// Test subcellCount and cumulativeOffset for mixed refinement
//--------------------------------------------------------------------------------------------------
TEST( RigNonUniformRefinement, MixedRefinement )
{
RigNonUniformRefinement nuRef( cvf::Vec3st( 3, 1, 1 ) );
// Cell 0: 2 subcells, Cell 1: 5 subcells, Cell 2: 1 subcell (default)
auto fractions2 = RigNonUniformRefinement::widthsToCumulativeFractions( { 0.5, 0.5 } );
auto fractions5 = RigNonUniformRefinement::widthsToCumulativeFractions( { 0.2, 0.4, 0.5, 0.6, 0.8 } );
nuRef.setCumulativeFractions( RigNonUniformRefinement::DimI, 0, fractions2 );
nuRef.setCumulativeFractions( RigNonUniformRefinement::DimI, 1, fractions5 );
EXPECT_EQ( nuRef.subcellCount( RigNonUniformRefinement::DimI, 0 ), 2u );
EXPECT_EQ( nuRef.subcellCount( RigNonUniformRefinement::DimI, 1 ), 5u );
EXPECT_EQ( nuRef.subcellCount( RigNonUniformRefinement::DimI, 2 ), 1u );
EXPECT_EQ( nuRef.cumulativeOffset( RigNonUniformRefinement::DimI, 0 ), 0u );
EXPECT_EQ( nuRef.cumulativeOffset( RigNonUniformRefinement::DimI, 1 ), 2u );
EXPECT_EQ( nuRef.cumulativeOffset( RigNonUniformRefinement::DimI, 2 ), 7u );
EXPECT_EQ( nuRef.totalRefinedCount( RigNonUniformRefinement::DimI ), 8u ); // 2 + 5 + 1
EXPECT_TRUE( nuRef.hasNonUniformRefinement() );
}
//--------------------------------------------------------------------------------------------------
/// Test mapRefinedToOriginal round-trips with cumulativeOffset
//--------------------------------------------------------------------------------------------------
TEST( RigNonUniformRefinement, MapRefinedToOriginalRoundTrip )
{
RigNonUniformRefinement nuRef( cvf::Vec3st( 4, 1, 1 ) );
// Cell 0: 3 subcells, Cell 1: 1 subcell, Cell 2: 2 subcells, Cell 3: 4 subcells
nuRef.setCumulativeFractions( RigNonUniformRefinement::DimI, 0, RigNonUniformRefinement::widthsToCumulativeFractions( { 1.0, 1.0, 1.0 } ) );
// Cell 1: default (1 subcell)
nuRef.setCumulativeFractions( RigNonUniformRefinement::DimI, 2, RigNonUniformRefinement::widthsToCumulativeFractions( { 1.0, 1.0 } ) );
nuRef.setCumulativeFractions( RigNonUniformRefinement::DimI,
3,
RigNonUniformRefinement::widthsToCumulativeFractions( { 1.0, 1.0, 1.0, 1.0 } ) );
EXPECT_EQ( nuRef.totalRefinedCount( RigNonUniformRefinement::DimI ), 10u ); // 3 + 1 + 2 + 4
// Verify round-trip: for each original cell, cumulativeOffset should map back correctly
for ( size_t origIdx = 0; origIdx < 4; ++origIdx )
{
size_t offset = nuRef.cumulativeOffset( RigNonUniformRefinement::DimI, origIdx );
auto [mappedOrig, mappedSub] = nuRef.mapRefinedToOriginal( RigNonUniformRefinement::DimI, offset );
EXPECT_EQ( mappedOrig, origIdx );
EXPECT_EQ( mappedSub, 0u );
}
// Check specific refined indices
auto [orig0, sub0] = nuRef.mapRefinedToOriginal( RigNonUniformRefinement::DimI, 0 ); // First subcell of cell 0
EXPECT_EQ( orig0, 0u );
EXPECT_EQ( sub0, 0u );
auto [orig1, sub1] = nuRef.mapRefinedToOriginal( RigNonUniformRefinement::DimI, 2 ); // Last subcell of cell 0
EXPECT_EQ( orig1, 0u );
EXPECT_EQ( sub1, 2u );
auto [orig2, sub2] = nuRef.mapRefinedToOriginal( RigNonUniformRefinement::DimI, 3 ); // Cell 1 (only subcell)
EXPECT_EQ( orig2, 1u );
EXPECT_EQ( sub2, 0u );
auto [orig3, sub3] = nuRef.mapRefinedToOriginal( RigNonUniformRefinement::DimI, 4 ); // First subcell of cell 2
EXPECT_EQ( orig3, 2u );
EXPECT_EQ( sub3, 0u );
auto [orig4, sub4] = nuRef.mapRefinedToOriginal( RigNonUniformRefinement::DimI, 9 ); // Last subcell of cell 3
EXPECT_EQ( orig4, 3u );
EXPECT_EQ( sub4, 3u );
}
//--------------------------------------------------------------------------------------------------
/// Test widthsToCumulativeFractions normalization
//--------------------------------------------------------------------------------------------------
TEST( RigNonUniformRefinement, WidthsToCumulativeFractions )
{
// Equal widths
auto fractions = RigNonUniformRefinement::widthsToCumulativeFractions( { 1.0, 1.0, 1.0 } );
EXPECT_EQ( fractions.size(), 3u );
EXPECT_NEAR( fractions[0], 1.0 / 3.0, 1e-10 );
EXPECT_NEAR( fractions[1], 2.0 / 3.0, 1e-10 );
EXPECT_DOUBLE_EQ( fractions[2], 1.0 );
// Non-equal widths
auto fractions2 = RigNonUniformRefinement::widthsToCumulativeFractions( { 1.0, 3.0 } );
EXPECT_EQ( fractions2.size(), 2u );
EXPECT_NEAR( fractions2[0], 0.25, 1e-10 );
EXPECT_DOUBLE_EQ( fractions2[1], 1.0 );
// Empty input
auto fractions3 = RigNonUniformRefinement::widthsToCumulativeFractions( {} );
EXPECT_EQ( fractions3.size(), 1u );
EXPECT_DOUBLE_EQ( fractions3[0], 1.0 );
}
//--------------------------------------------------------------------------------------------------
/// Test single cell edge case
//--------------------------------------------------------------------------------------------------
TEST( RigNonUniformRefinement, SingleCell )
{
RigNonUniformRefinement nuRef( cvf::Vec3st( 1, 1, 1 ) );
EXPECT_EQ( nuRef.totalRefinedCount( RigNonUniformRefinement::DimI ), 1u );
EXPECT_EQ( nuRef.subcellCount( RigNonUniformRefinement::DimI, 0 ), 1u );
EXPECT_EQ( nuRef.cumulativeOffset( RigNonUniformRefinement::DimI, 0 ), 0u );
auto [orig, sub] = nuRef.mapRefinedToOriginal( RigNonUniformRefinement::DimI, 0 );
EXPECT_EQ( orig, 0u );
EXPECT_EQ( sub, 0u );
}
//--------------------------------------------------------------------------------------------------
/// Test distributeWidthsAcrossCells: 3 equal widths across 3 cells → 1 subcell each
//--------------------------------------------------------------------------------------------------
TEST( RigNonUniformRefinement, DistributeEqualWidthsAcrossEqualCells )
{
RigNonUniformRefinement nuRef( cvf::Vec3st( 3, 1, 1 ) );
// 3 equal widths across 3 cells: each boundary lands on a cell boundary
// so each cell gets exactly 1 subcell
nuRef.distributeWidthsAcrossCells( RigNonUniformRefinement::DimI, 0, 2, { 1.0, 1.0, 1.0 } );
EXPECT_EQ( nuRef.subcellCount( RigNonUniformRefinement::DimI, 0 ), 1u );
EXPECT_EQ( nuRef.subcellCount( RigNonUniformRefinement::DimI, 1 ), 1u );
EXPECT_EQ( nuRef.subcellCount( RigNonUniformRefinement::DimI, 2 ), 1u );
EXPECT_EQ( nuRef.totalRefinedCount( RigNonUniformRefinement::DimI ), 3u );
}
//--------------------------------------------------------------------------------------------------
/// Test distributeWidthsAcrossCells: 6 equal widths across 3 cells → 2 subcells each
//--------------------------------------------------------------------------------------------------
TEST( RigNonUniformRefinement, DistributeDoubleWidthsAcrossCells )
{
RigNonUniformRefinement nuRef( cvf::Vec3st( 3, 1, 1 ) );
// 6 equal widths across 3 cells: 2 subcell boundaries per cell
nuRef.distributeWidthsAcrossCells( RigNonUniformRefinement::DimI, 0, 2, { 1.0, 1.0, 1.0, 1.0, 1.0, 1.0 } );
EXPECT_EQ( nuRef.subcellCount( RigNonUniformRefinement::DimI, 0 ), 2u );
EXPECT_EQ( nuRef.subcellCount( RigNonUniformRefinement::DimI, 1 ), 2u );
EXPECT_EQ( nuRef.subcellCount( RigNonUniformRefinement::DimI, 2 ), 2u );
EXPECT_EQ( nuRef.totalRefinedCount( RigNonUniformRefinement::DimI ), 6u );
}
//--------------------------------------------------------------------------------------------------
/// Test distributeWidthsAcrossCells: uneven widths distributed across cells
//--------------------------------------------------------------------------------------------------
TEST( RigNonUniformRefinement, DistributeUnevenWidthsAcrossCells )
{
RigNonUniformRefinement nuRef( cvf::Vec3st( 3, 1, 1 ) );
// 5 widths [0.1, 0.1, 0.1, 0.1, 0.6] across 3 cells
// Cumulative fractions: [0.1, 0.2, 0.3, 0.4, 1.0]
// Cell 0 [0, 0.333]: boundaries 0.1, 0.2, 0.3 fall inside → 4 subcells
// Cell 1 [0.333, 0.667]: boundary 0.4 falls inside → 2 subcells
// Cell 2 [0.667, 1.0]: no interior boundaries → 1 subcell
nuRef.distributeWidthsAcrossCells( RigNonUniformRefinement::DimI, 0, 2, { 0.1, 0.1, 0.1, 0.1, 0.6 } );
EXPECT_EQ( nuRef.subcellCount( RigNonUniformRefinement::DimI, 0 ), 4u );
EXPECT_EQ( nuRef.subcellCount( RigNonUniformRefinement::DimI, 1 ), 2u );
EXPECT_EQ( nuRef.subcellCount( RigNonUniformRefinement::DimI, 2 ), 1u );
EXPECT_EQ( nuRef.totalRefinedCount( RigNonUniformRefinement::DimI ), 7u );
}
//--------------------------------------------------------------------------------------------------
/// Test distributeWidthsAcrossCells: single cell range
//--------------------------------------------------------------------------------------------------
TEST( RigNonUniformRefinement, DistributeWidthsSingleCell )
{
RigNonUniformRefinement nuRef( cvf::Vec3st( 5, 1, 1 ) );
// Apply 3 widths to just cell 2
nuRef.distributeWidthsAcrossCells( RigNonUniformRefinement::DimI, 2, 2, { 0.2, 0.3, 0.5 } );
EXPECT_EQ( nuRef.subcellCount( RigNonUniformRefinement::DimI, 0 ), 1u ); // Unaffected
EXPECT_EQ( nuRef.subcellCount( RigNonUniformRefinement::DimI, 1 ), 1u ); // Unaffected
EXPECT_EQ( nuRef.subcellCount( RigNonUniformRefinement::DimI, 2 ), 3u ); // 3 subcells
EXPECT_EQ( nuRef.subcellCount( RigNonUniformRefinement::DimI, 3 ), 1u ); // Unaffected
EXPECT_EQ( nuRef.subcellCount( RigNonUniformRefinement::DimI, 4 ), 1u ); // Unaffected
EXPECT_EQ( nuRef.totalRefinedCount( RigNonUniformRefinement::DimI ), 7u ); // 1+1+3+1+1
}
//--------------------------------------------------------------------------------------------------
/// Test distributeWidthsAcrossCells: many widths concentrate in first cell
//--------------------------------------------------------------------------------------------------
TEST( RigNonUniformRefinement, DistributeConcentratedWidths )
{
RigNonUniformRefinement nuRef( cvf::Vec3st( 2, 1, 1 ) );
// 4 widths [0.1, 0.1, 0.1, 0.7] across 2 cells
// Cumulative fractions: [0.1, 0.2, 0.3, 1.0]
// Cell 0 [0, 0.5]: boundaries 0.1, 0.2, 0.3 fall inside → 4 subcells
// Cell 1 [0.5, 1.0]: no interior boundaries → 1 subcell
nuRef.distributeWidthsAcrossCells( RigNonUniformRefinement::DimI, 0, 1, { 0.1, 0.1, 0.1, 0.7 } );
EXPECT_EQ( nuRef.subcellCount( RigNonUniformRefinement::DimI, 0 ), 4u );
EXPECT_EQ( nuRef.subcellCount( RigNonUniformRefinement::DimI, 1 ), 1u );
EXPECT_EQ( nuRef.totalRefinedCount( RigNonUniformRefinement::DimI ), 5u );
}
//--------------------------------------------------------------------------------------------------
/// Test non-uniform geometry: equal fractions should match uniform createHexCornerCoords
//--------------------------------------------------------------------------------------------------
TEST( RiaCellDividingTools, NonUniformEqualFractionsMatchUniform )
{
std::array<cvf::Vec3d, 8> corners = { cvf::Vec3d( 0, 0, 0 ),
cvf::Vec3d( 1, 0, 0 ),
cvf::Vec3d( 1, 1, 0 ),
cvf::Vec3d( 0, 1, 0 ),
cvf::Vec3d( 0, 0, 1 ),
cvf::Vec3d( 1, 0, 1 ),
cvf::Vec3d( 1, 1, 1 ),
cvf::Vec3d( 0, 1, 1 ) };
size_t nx = 3, ny = 2, nz = 2;
// Uniform result
auto uniformResult = RiaCellDividingTools::createHexCornerCoords( corners, nx, ny, nz );
// Non-uniform with equal fractions
std::vector<double> fracX = { 1.0 / 3, 2.0 / 3, 1.0 };
std::vector<double> fracY = { 0.5, 1.0 };
std::vector<double> fracZ = { 0.5, 1.0 };
auto nonUniformResult = RiaCellDividingTools::createHexCornerCoords( corners, fracX, fracY, fracZ );
ASSERT_EQ( uniformResult.size(), nonUniformResult.size() );
for ( size_t i = 0; i < uniformResult.size(); ++i )
{
EXPECT_NEAR( uniformResult[i].x(), nonUniformResult[i].x(), 1e-10 ) << "Mismatch at index " << i << " (x)";
EXPECT_NEAR( uniformResult[i].y(), nonUniformResult[i].y(), 1e-10 ) << "Mismatch at index " << i << " (y)";
EXPECT_NEAR( uniformResult[i].z(), nonUniformResult[i].z(), 1e-10 ) << "Mismatch at index " << i << " (z)";
}
}
//--------------------------------------------------------------------------------------------------
/// Test non-uniform fractions produce geometrically correct subcells
//--------------------------------------------------------------------------------------------------
TEST( RiaCellDividingTools, NonUniformFractionsGeometry )
{
std::array<cvf::Vec3d, 8> corners = { cvf::Vec3d( 0, 0, 0 ),
cvf::Vec3d( 10, 0, 0 ),
cvf::Vec3d( 10, 10, 0 ),
cvf::Vec3d( 0, 10, 0 ),
cvf::Vec3d( 0, 0, 10 ),
cvf::Vec3d( 10, 0, 10 ),
cvf::Vec3d( 10, 10, 10 ),
cvf::Vec3d( 0, 10, 10 ) };
// Non-uniform X: 20% and 80%
std::vector<double> fracX = { 0.2, 1.0 };
std::vector<double> fracY = { 1.0 }; // No refinement
std::vector<double> fracZ = { 1.0 }; // No refinement
auto result = RiaCellDividingTools::createHexCornerCoords( corners, fracX, fracY, fracZ );
// Should produce 2*1*1 = 2 subcells, each with 8 corners
ASSERT_EQ( result.size(), 16u );
// First subcell: x range [0, 2]
EXPECT_NEAR( result[0].x(), 0.0, 1e-10 ); // corner 0 x
EXPECT_NEAR( result[1].x(), 2.0, 1e-10 ); // corner 1 x
// Second subcell: x range [2, 10]
EXPECT_NEAR( result[8].x(), 2.0, 1e-10 ); // corner 0 x
EXPECT_NEAR( result[9].x(), 10.0, 1e-10 ); // corner 1 x
}
//--------------------------------------------------------------------------------------------------
/// Test transformIjkToSectorCoordinates with non-uniform refinement
//--------------------------------------------------------------------------------------------------
TEST( RigGridExportAdapter, TransformIjkNonUniform )
{
RigNonUniformRefinement nuRef( cvf::Vec3st( 3, 2, 1 ) );
// I: cell 0 -> 2 subcells, cell 1 -> 3 subcells, cell 2 -> 1 subcell
nuRef.setCumulativeFractions( RigNonUniformRefinement::DimI, 0, RigNonUniformRefinement::widthsToCumulativeFractions( { 1.0, 1.0 } ) );
nuRef.setCumulativeFractions( RigNonUniformRefinement::DimI, 1, RigNonUniformRefinement::widthsToCumulativeFractions( { 1.0, 1.0, 1.0 } ) );
// Cell 2: default (1 subcell)
// J: cell 0 -> 1 subcell, cell 1 -> 2 subcells
nuRef.setCumulativeFractions( RigNonUniformRefinement::DimJ, 1, RigNonUniformRefinement::widthsToCumulativeFractions( { 1.0, 1.0 } ) );
caf::VecIjk0 sectorMin( 10, 20, 30 );
caf::VecIjk0 sectorMax( 12, 21, 30 ); // 3 cells I, 2 cells J, 1 cell K
// Test min coordinate (default): cell (10, 20, 30) -> offset(0), offset(0), offset(0) = (0, 0, 0)
auto result = RigGridExportAdapter::transformIjkToSectorCoordinates( caf::VecIjk0( 10, 20, 30 ), sectorMin, sectorMax, nuRef );
ASSERT_TRUE( result.has_value() );
EXPECT_EQ( result->x(), 0u );
EXPECT_EQ( result->y(), 0u );
EXPECT_EQ( result->z(), 0u );
// Test cell (11, 20, 30) -> offset(1) = 2, offset(0) = 0, offset(0) = 0
auto result2 = RigGridExportAdapter::transformIjkToSectorCoordinates( caf::VecIjk0( 11, 20, 30 ), sectorMin, sectorMax, nuRef );
ASSERT_TRUE( result2.has_value() );
EXPECT_EQ( result2->x(), 2u );
EXPECT_EQ( result2->y(), 0u );
// Test max coordinate: cell (12, 21, 30) -> offset(3)-1 = 5, offset(2)-1 = 2, offset(1)-1 = 0
auto result3 =
RigGridExportAdapter::transformIjkToSectorCoordinates( caf::VecIjk0( 12, 21, 30 ), sectorMin, sectorMax, nuRef, false, true );
ASSERT_TRUE( result3.has_value() );
EXPECT_EQ( result3->x(), 5u ); // 2 + 3 + 1 - 1 = 5
EXPECT_EQ( result3->y(), 2u ); // 1 + 2 - 1 = 2
EXPECT_EQ( result3->z(), 0u );
// Test centering: cell (11, 20, 30) -> offset(1) + subcellCount(1)/2 = 2 + 1 = 3
auto result4 = RigGridExportAdapter::transformIjkToSectorCoordinates( caf::VecIjk0( 11, 20, 30 ), sectorMin, sectorMax, nuRef, true );
ASSERT_TRUE( result4.has_value() );
EXPECT_EQ( result4->x(), 3u ); // offset(1)=2, subcellCount(1)=3, center=2+3/2=3
EXPECT_EQ( result4->y(), 0u ); // offset(0)=0, subcellCount(0)=1, center=0+1/2=0
// Test out of bounds
auto result5 = RigGridExportAdapter::transformIjkToSectorCoordinates( caf::VecIjk0( 9, 20, 30 ), sectorMin, sectorMax, nuRef );
EXPECT_FALSE( result5.has_value() );
}