#12950 Generate BCCON/BCPROP for simulation input.

This commit is contained in:
Kristian Bendiksen
2025-11-12 08:39:20 +01:00
parent 4291ebd91c
commit b1646d9a1f
6 changed files with 303 additions and 5 deletions
@@ -64,7 +64,9 @@
#include "opm/input/eclipse/Deck/DeckItem.hpp"
#include "opm/input/eclipse/Deck/DeckKeyword.hpp"
#include "opm/input/eclipse/Deck/DeckRecord.hpp"
#include "opm/input/eclipse/Parser/ParserKeyword.hpp"
#include "opm/input/eclipse/Parser/ParserKeywords/C.hpp"
#include "opm/input/eclipse/Parser/ParserKeywords/S.hpp"
#include "opm/input/eclipse/Parser/ParserKeywords/W.hpp"
#include "opm/input/eclipse/Utility/Typetools.hpp"
@@ -319,6 +321,9 @@ std::expected<void, QString> RicExportEclipseSectorModelFeature::exportSimulatio
{
RigEclipseResultTools::generateBorderResult( &eclipseCase, bordnumVisibility, RiaResultNames::bordnum() );
// Generate BCCON result to assign values 1-6 based on which face of the box the border cells are on
RigEclipseResultTools::generateBcconResult( &eclipseCase, exportSettings.min(), exportSettings.max() );
// Generate OPERNUM result based on BORDNUM (border cells get max existing OPERNUM + 1)
RigEclipseResultTools::generateOperNumResult( &eclipseCase );
}
@@ -488,6 +493,28 @@ std::expected<void, QString> RicExportEclipseSectorModelFeature::addBorderBounda
if ( !borderCellFaces.empty() )
{
// Transform border cell face coordinates to sector-relative coordinates
for ( auto& face : borderCellFaces )
{
auto transformResult =
transformIjkToSectorCoordinates( face.ijk, exportSettings.min(), exportSettings.max(), exportSettings.refinement() );
if ( !transformResult )
{
RiaLogging::warning( QString( "Failed to transform border cell face at (%1, %2, %3): %4" )
.arg( face.ijk.x() )
.arg( face.ijk.y() )
.arg( face.ijk.z() )
.arg( transformResult.error() ) );
continue;
}
// Update the IJK coordinates to sector-relative (1-based Eclipse coordinates)
// Note: transformIjkToSectorCoordinates returns 1-based coordinates, but we need to convert back to 0-based
// for the BorderCellFace struct
face.ijk = cvf::Vec3st( transformResult->x() - 1, transformResult->y() - 1, transformResult->z() - 1 );
}
// Create BCCON keyword using the factory
Opm::DeckKeyword bcconKw = RimKeywordFactory::bcconKeyword( borderCellFaces );
@@ -513,7 +540,7 @@ std::expected<void, QString> RicExportEclipseSectorModelFeature::addBorderBounda
Opm::DeckKeyword bcpropKw = RimKeywordFactory::bcpropKeyword( borderCellFaces, bcpropRecords );
// Replace BCPROP keyword in GRID section
if ( !deckFile.replaceKeyword( "GRID", bcpropKw ) )
if ( !deckFile.replaceKeyword( Opm::ParserKeywords::SCHEDULE::keywordName, bcpropKw ) )
{
return std::unexpected( "Failed to replace BCPROP keyword in deck file" );
}
@@ -204,6 +204,19 @@ static std::optional<Opm::FileDeck::Index> findSectionInsertionPoint( std::uniqu
insertIdx++;
}
// If insertIdx is at stop() (end of deck), we can't insert there as it's beyond valid blocks
// We need to insert at the end of the last block instead
if ( insertIdx == fileDeck->stop() )
{
// Decrement to get the last valid keyword
auto lastKeyword = fileDeck->stop();
--lastKeyword;
// We want to insert after the last keyword, so increment the keyword_index
insertIdx.file_index = lastKeyword.file_index;
insertIdx.keyword_index = lastKeyword.keyword_index + 1;
}
return insertIdx;
}
@@ -479,6 +479,7 @@ Opm::DeckKeyword bcconKeyword( const std::vector<RigEclipseResultTools::BorderCe
Opm::DeckKeyword kw{ Opm::ParserKeywords::BCCON() };
int bcconIndex = 1;
for ( const auto& borderFace : borderCellFaces )
{
// Convert from 0-based to 1-based Eclipse indexing
@@ -491,7 +492,7 @@ Opm::DeckKeyword bcconKeyword( const std::vector<RigEclipseResultTools::BorderCe
// Create items for the record
std::vector<Opm::DeckItem> items;
items.push_back( RifOpmDeckTools::item( B::INDEX::itemName, borderFace.boundaryCondition ) );
items.push_back( RifOpmDeckTools::item( B::INDEX::itemName, bcconIndex ) );
items.push_back( RifOpmDeckTools::item( B::I1::itemName, i1 ) );
items.push_back( RifOpmDeckTools::item( B::I2::itemName, i1 ) );
items.push_back( RifOpmDeckTools::item( B::J1::itemName, j1 ) );
@@ -501,6 +502,7 @@ Opm::DeckKeyword bcconKeyword( const std::vector<RigEclipseResultTools::BorderCe
items.push_back( RifOpmDeckTools::item( B::DIRECTION::itemName, faceStr ) );
kw.addRecord( Opm::DeckRecord{ std::move( items ) } );
bcconIndex++;
}
return kw;
@@ -521,6 +523,7 @@ Opm::DeckKeyword bcpropKeyword( const std::vector<RigEclipseResultTools::BorderC
Opm::DeckKeyword kw{ Opm::ParserKeywords::BCPROP() };
int bcIndex = 1;
// Add one entry per boundary condition
for ( const auto& bc : boundaryConditions )
{
@@ -528,7 +531,7 @@ Opm::DeckKeyword bcpropKeyword( const std::vector<RigEclipseResultTools::BorderC
// Find the corresponding property record
// The properties vector should be indexed by boundaryCondition - 1
size_t propIndex = static_cast<size_t>( bc.boundaryCondition - 1 );
size_t propIndex = static_cast<size_t>( bc.boundaryCondition );
if ( propIndex < boundaryConditionProperties.size() )
{
const auto& propRecord = boundaryConditionProperties[propIndex];
@@ -537,16 +540,20 @@ Opm::DeckKeyword bcpropKeyword( const std::vector<RigEclipseResultTools::BorderC
std::vector<Opm::DeckItem> items;
// Add INDEX field
items.push_back( RifOpmDeckTools::item( B::INDEX::itemName, bc.boundaryCondition ) );
items.push_back( RifOpmDeckTools::item( B::INDEX::itemName, bcIndex ) );
// Copy all items from the property record (which doesn't include INDEX)
for ( size_t i = 0; i < propRecord.size(); ++i )
{
items.push_back( propRecord.getItem( i ) );
if ( propRecord.getItem( i ).name() != B::INDEX::itemName )
{
items.push_back( propRecord.getItem( i ) );
}
}
kw.addRecord( Opm::DeckRecord{ std::move( items ) } );
}
bcIndex++;
}
return kw;
@@ -237,6 +237,86 @@ int findMaxOperNumValue( RimEclipseCase* eclipseCase )
return maxValue;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
void generateBcconResult( RimEclipseCase* eclipseCase, const cvf::Vec3st& min, const cvf::Vec3st& max )
{
if ( eclipseCase == nullptr ) return;
auto resultsData = eclipseCase->results( RiaDefines::PorosityModelType::MATRIX_MODEL );
if ( !resultsData ) return;
auto grid = eclipseCase->eclipseCaseData()->mainGrid();
if ( !grid ) return;
// Check if BORDNUM result exists
RigEclipseResultAddress bordNumAddr( RiaDefines::ResultCatType::GENERATED, RiaDefines::ResultDataType::INTEGER, RiaResultNames::bordnum() );
if ( !resultsData->hasResultEntry( bordNumAddr ) )
{
RiaLogging::warning( "BORDNUM result not found - cannot generate BCCON result" );
return;
}
resultsData->ensureKnownResultLoaded( bordNumAddr );
auto bordNumValues = resultsData->cellScalarResults( bordNumAddr, 0 );
if ( bordNumValues.empty() ) return;
auto activeReservoirCellIdxs =
eclipseCase->eclipseCaseData()->activeCellInfo( RiaDefines::PorosityModelType::MATRIX_MODEL )->activeReservoirCellIndices();
size_t reservoirCellCount =
eclipseCase->eclipseCaseData()->activeCellInfo( RiaDefines::PorosityModelType::MATRIX_MODEL )->reservoirCellCount();
std::vector<int> result( reservoirCellCount, 0 );
// Iterate through all active cells
for ( auto activeCellIdx : activeReservoirCellIdxs )
{
// Check if this cell is a border cell
int borderValue = static_cast<int>( bordNumValues[activeCellIdx.value()] );
if ( borderValue != BorderType::BORDER_CELL ) continue;
// Get IJK indices for this cell
size_t i, j, k;
if ( !grid->ijkFromCellIndex( activeCellIdx.value(), &i, &j, &k ) ) continue;
// Determine which face of the box this cell is on
// Priority: I faces, then J faces, then K faces (for corner/edge cells)
int bcconValue = 0;
if ( i == min.x() )
{
bcconValue = 1; // I- face
}
else if ( i == max.x() )
{
bcconValue = 2; // I+ face
}
else if ( j == min.y() )
{
bcconValue = 3; // J- face
}
else if ( j == max.y() )
{
bcconValue = 4; // J+ face
}
else if ( k == min.z() )
{
bcconValue = 5; // K- face
}
else if ( k == max.z() )
{
bcconValue = 6; // K+ face
}
result[activeCellIdx.value()] = bcconValue;
}
RigEclipseResultTools::createResultVector( *eclipseCase, "BCCON", result );
eclipseCase->updateConnectedEditors();
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
@@ -51,6 +51,8 @@ int generateOperNumResult( RimEclipseCase* eclipseCase, int borderCellValue = -1
int findMaxOperNumValue( RimEclipseCase* eclipseCase );
void generateBcconResult( RimEclipseCase* eclipseCase, const cvf::Vec3st& min, const cvf::Vec3st& max );
std::vector<BorderCellFace> generateBorderCellFaces( RimEclipseCase* eclipseCase );
} // namespace RigEclipseResultTools
@@ -206,3 +206,172 @@ TEST( RigEclipseResultToolsTest, BorderCellBcconGeneration )
content.contains( " Y " ) || content.contains( " Z " );
EXPECT_TRUE( hasFaceDirection ) << "BCCON entries don't contain face directions";
}
//--------------------------------------------------------------------------------------------------
/// Test BCCON result generation with face numbering 1-6
///
/// This test verifies that:
/// 1. We can generate BORDNUM result identifying border and interior cells
/// 2. generateBcconResult() assigns values 1-6 to border cells based on which face of the box they're on
/// 3. Values are assigned correctly: 1=I-, 2=I+, 3=J-, 4=J+, 5=K-, 6=K+
///
/// Test process:
/// 1. Load test model
/// 2. Create custom visibility for a box region
/// 3. Generate BORDNUM result
/// 4. Generate BCCON result with box bounds
/// 5. Verify border cells have correct BCCON values (1-6) based on their position
//--------------------------------------------------------------------------------------------------
TEST( RigEclipseResultToolsTest, BcconResultWithFaceNumbering )
{
// Setup test data directory
QDir baseFolder( TEST_MODEL_DIR );
bool subFolderExists = baseFolder.cd( "Case_with_10_timesteps/Real0" );
ASSERT_TRUE( subFolderExists ) << "Test model directory not found";
QString inputFilename( "BRUGGE_0000.EGRID" );
QString inputFilePath = baseFolder.absoluteFilePath( inputFilename );
ASSERT_TRUE( QFile::exists( inputFilePath ) ) << "Test model file not found: " << inputFilePath.toStdString();
// Step 1: Load grid
std::unique_ptr<RimEclipseResultCase> testCase( new RimEclipseResultCase );
cvf::ref<RigEclipseCaseData> caseData = new RigEclipseCaseData( testCase.get() );
cvf::ref<RifReaderEclipseOutput> readerInterfaceEcl = new RifReaderEclipseOutput;
bool success = readerInterfaceEcl->open( inputFilePath, caseData.p() );
ASSERT_TRUE( success ) << "Failed to load grid";
testCase->setReservoirData( caseData.p() );
const RigMainGrid* grid = caseData->mainGrid();
ASSERT_NE( grid, nullptr ) << "Grid is null";
// Step 2: Create custom visibility for a box region
size_t cellCount = grid->cellCount();
cvf::ref<cvf::UByteArray> customVisibility = new cvf::UByteArray( cellCount );
customVisibility->setAll( 0 ); // Start with all invisible
// Define a box in the middle of the grid
size_t startI = grid->cellCountI() / 4;
size_t endI = 3 * grid->cellCountI() / 4;
size_t startJ = grid->cellCountJ() / 4;
size_t endJ = 3 * grid->cellCountJ() / 4;
size_t startK = grid->cellCountK() / 4;
size_t endK = 3 * grid->cellCountK() / 4;
cvf::Vec3st min( startI, startJ, startK );
cvf::Vec3st max( endI - 1, endJ - 1, endK - 1 );
// Make the box visible
for ( size_t i = startI; i < endI; ++i )
{
for ( size_t j = startJ; j < endJ; ++j )
{
for ( size_t k = startK; k < endK; ++k )
{
size_t cellIndex = grid->cellIndexFromIJK( i, j, k );
( *customVisibility )[cellIndex] = 1;
}
}
}
// Step 3: Generate BORDNUM result
RigEclipseResultTools::generateBorderResult( testCase.get(), customVisibility, "BORDNUM" );
// Verify BORDNUM was created
auto resultsData = testCase->results( RiaDefines::PorosityModelType::MATRIX_MODEL );
ASSERT_NE( resultsData, nullptr );
RigEclipseResultAddress bordNumAddr( RiaDefines::ResultCatType::GENERATED, RiaDefines::ResultDataType::INTEGER, "BORDNUM" );
ASSERT_TRUE( resultsData->hasResultEntry( bordNumAddr ) ) << "BORDNUM result not created";
// Step 4: Generate BCCON result
RigEclipseResultTools::generateBcconResult( testCase.get(), min, max );
// Verify BCCON was created
RigEclipseResultAddress bcconAddr( RiaDefines::ResultCatType::GENERATED, RiaDefines::ResultDataType::INTEGER, "BCCON" );
ASSERT_TRUE( resultsData->hasResultEntry( bcconAddr ) ) << "BCCON result not created";
// Get BCCON and BORDNUM values
resultsData->ensureKnownResultLoaded( bcconAddr );
resultsData->ensureKnownResultLoaded( bordNumAddr );
auto bcconValues = resultsData->cellScalarResults( bcconAddr, 0 );
auto bordNumValues = resultsData->cellScalarResults( bordNumAddr, 0 );
ASSERT_FALSE( bcconValues.empty() ) << "BCCON values are empty";
ASSERT_FALSE( bordNumValues.empty() ) << "BORDNUM values are empty";
// Step 5: Verify BCCON values are correct
int countI_minus = 0, countI_plus = 0;
int countJ_minus = 0, countJ_plus = 0;
int countK_minus = 0, countK_plus = 0;
for ( size_t i = startI; i < endI; ++i )
{
for ( size_t j = startJ; j < endJ; ++j )
{
for ( size_t k = startK; k < endK; ++k )
{
size_t cellIndex = grid->cellIndexFromIJK( i, j, k );
int borderValue = static_cast<int>( bordNumValues[cellIndex] );
int bcconValue = static_cast<int>( bcconValues[cellIndex] );
// Only check border cells
if ( borderValue == RigEclipseResultTools::BorderType::BORDER_CELL )
{
// Verify BCCON value is in valid range 1-6
EXPECT_GE( bcconValue, 1 ) << "BCCON value out of range at (" << i << "," << j << "," << k << ")";
EXPECT_LE( bcconValue, 6 ) << "BCCON value out of range at (" << i << "," << j << "," << k << ")";
// Check specific face values
if ( i == min.x() )
{
EXPECT_EQ( bcconValue, 1 ) << "I- face should have BCCON=1 at (" << i << "," << j << "," << k << ")";
countI_minus++;
}
else if ( i == max.x() )
{
EXPECT_EQ( bcconValue, 2 ) << "I+ face should have BCCON=2 at (" << i << "," << j << "," << k << ")";
countI_plus++;
}
else if ( j == min.y() )
{
EXPECT_EQ( bcconValue, 3 ) << "J- face should have BCCON=3 at (" << i << "," << j << "," << k << ")";
countJ_minus++;
}
else if ( j == max.y() )
{
EXPECT_EQ( bcconValue, 4 ) << "J+ face should have BCCON=4 at (" << i << "," << j << "," << k << ")";
countJ_plus++;
}
else if ( k == min.z() )
{
EXPECT_EQ( bcconValue, 5 ) << "K- face should have BCCON=5 at (" << i << "," << j << "," << k << ")";
countK_minus++;
}
else if ( k == max.z() )
{
EXPECT_EQ( bcconValue, 6 ) << "K+ face should have BCCON=6 at (" << i << "," << j << "," << k << ")";
countK_plus++;
}
}
else
{
// Interior cells should have BCCON=0
EXPECT_EQ( bcconValue, 0 ) << "Interior cell should have BCCON=0 at (" << i << "," << j << "," << k << ")";
}
}
}
}
// Verify we found border cells on all 6 faces
EXPECT_GT( countI_minus, 0 ) << "No border cells found on I- face";
EXPECT_GT( countI_plus, 0 ) << "No border cells found on I+ face";
EXPECT_GT( countJ_minus, 0 ) << "No border cells found on J- face";
EXPECT_GT( countJ_plus, 0 ) << "No border cells found on J+ face";
EXPECT_GT( countK_minus, 0 ) << "No border cells found on K- face";
EXPECT_GT( countK_plus, 0 ) << "No border cells found on K+ face";
}