From 7d456f1cbb899d29e851d23d6db18be88d0ea2ea Mon Sep 17 00:00:00 2001 From: Kristian Bendiksen Date: Fri, 16 Oct 2020 11:27:46 +0200 Subject: [PATCH] Refactor fracture model result calculation to be independent of display. --- .../RicNewFractureModelPlotFeature.cpp | 65 +- .../RicNewFractureModelPlotFeature.h | 17 +- ...icExportFractureModelPlotToFileFeature.cpp | 2 +- .../RifFractureModelAsymmetricFrkExporter.cpp | 9 +- .../RifFractureModelAsymmetricFrkExporter.h | 4 +- .../RifFractureModelDeviationFrkExporter.cpp | 9 +- .../RifFractureModelDeviationFrkExporter.h | 4 +- .../RifFractureModelGeologicalFrkExporter.cpp | 45 +- .../RifFractureModelGeologicalFrkExporter.h | 4 +- .../RifFractureModelPerfsFrkExporter.cpp | 9 +- .../RifFractureModelPerfsFrkExporter.h | 4 +- .../RifFractureModelPlotExporter.cpp | 16 +- .../RifFractureModelPlotExporter.h | 4 +- .../ProjectDataModel/CMakeLists_files.cmake | 11 + .../Completions/RimFractureModel.cpp | 223 ++++++- .../Completions/RimFractureModel.h | 44 +- .../RimElasticPropertiesCurve.cpp | 300 +-------- .../RimFractureModelCalculator.cpp | 571 ++++++++++++++++++ .../RimFractureModelCalculator.h | 96 +++ .../RimFractureModelCurve.cpp | 404 +------------ .../ProjectDataModel/RimFractureModelCurve.h | 43 -- ...FractureModelElasticPropertyCalculator.cpp | 302 +++++++++ ...imFractureModelElasticPropertyCalculator.h | 66 ++ .../RimFractureModelLayerCalculator.cpp | 157 +++++ .../RimFractureModelLayerCalculator.h | 46 ++ .../ProjectDataModel/RimFractureModelPlot.cpp | 541 +---------------- .../ProjectDataModel/RimFractureModelPlot.h | 54 +- .../RimFractureModelPropertyCalculator.h | 44 ++ .../RimFractureModelStressCalculator.cpp | 150 +++++ .../RimFractureModelStressCalculator.h | 52 ++ .../RimFractureModelStressCurve.cpp | 87 +-- .../RimFractureModelWellLogCalculator.cpp | 401 ++++++++++++ .../RimFractureModelWellLogCalculator.h | 70 +++ .../ProjectDataModel/RimLayerCurve.cpp | 112 +--- .../RimcFractureModelPlot.cpp | 2 +- 35 files changed, 2322 insertions(+), 1646 deletions(-) create mode 100644 ApplicationCode/ProjectDataModel/RimFractureModelCalculator.cpp create mode 100644 ApplicationCode/ProjectDataModel/RimFractureModelCalculator.h create mode 100644 ApplicationCode/ProjectDataModel/RimFractureModelElasticPropertyCalculator.cpp create mode 100644 ApplicationCode/ProjectDataModel/RimFractureModelElasticPropertyCalculator.h create mode 100644 ApplicationCode/ProjectDataModel/RimFractureModelLayerCalculator.cpp create mode 100644 ApplicationCode/ProjectDataModel/RimFractureModelLayerCalculator.h create mode 100644 ApplicationCode/ProjectDataModel/RimFractureModelPropertyCalculator.h create mode 100644 ApplicationCode/ProjectDataModel/RimFractureModelStressCalculator.cpp create mode 100644 ApplicationCode/ProjectDataModel/RimFractureModelStressCalculator.h create mode 100644 ApplicationCode/ProjectDataModel/RimFractureModelWellLogCalculator.cpp create mode 100644 ApplicationCode/ProjectDataModel/RimFractureModelWellLogCalculator.h diff --git a/ApplicationCode/Commands/CompletionCommands/RicNewFractureModelPlotFeature.cpp b/ApplicationCode/Commands/CompletionCommands/RicNewFractureModelPlotFeature.cpp index 1b31aa6cce..223e219a1c 100644 --- a/ApplicationCode/Commands/CompletionCommands/RicNewFractureModelPlotFeature.cpp +++ b/ApplicationCode/Commands/CompletionCommands/RicNewFractureModelPlotFeature.cpp @@ -94,34 +94,10 @@ RimFractureModelPlot* { auto task = progInfo.task( "Creating parameters track", 15 ); - std::map plots; - plots["Porosity"] = {std::make_tuple( "PORO", - RiaDefines::ResultCatType::STATIC_NATIVE, - RimFractureModelCurve::MissingValueStrategy::DEFAULT_VALUE, - false, - RiaDefines::CurveProperty::POROSITY )}; - - plots["Pressure"] = {std::make_tuple( "PRESSURE", - RiaDefines::ResultCatType::DYNAMIC_NATIVE, - RimFractureModelCurve::MissingValueStrategy::LINEAR_INTERPOLATION, - true, - RiaDefines::CurveProperty::INITIAL_PRESSURE ), - std::make_tuple( "PRESSURE", - RiaDefines::ResultCatType::DYNAMIC_NATIVE, - RimFractureModelCurve::MissingValueStrategy::OTHER_CURVE_PROPERTY, - false, - RiaDefines::CurveProperty::PRESSURE )}; - - plots["Permeability"] = {std::make_tuple( "PERMX", - RiaDefines::ResultCatType::STATIC_NATIVE, - RimFractureModelCurve::MissingValueStrategy::DEFAULT_VALUE, - false, - RiaDefines::CurveProperty::PERMEABILITY_X ), - std::make_tuple( "PERMZ", - RiaDefines::ResultCatType::STATIC_NATIVE, - RimFractureModelCurve::MissingValueStrategy::DEFAULT_VALUE, - false, - RiaDefines::CurveProperty::PERMEABILITY_Z )}; + std::map> plots; + plots["Porosity"] = {RiaDefines::CurveProperty::POROSITY}; + plots["Pressure"] = {RiaDefines::CurveProperty::INITIAL_PRESSURE, RiaDefines::CurveProperty::PRESSURE}; + plots["Permeability"] = {RiaDefines::CurveProperty::PERMEABILITY_X, RiaDefines::CurveProperty::PERMEABILITY_Z}; std::set logarithmicPlots; logarithmicPlots.insert( "Permeability" ); @@ -369,13 +345,13 @@ void RicNewFractureModelPlotFeature::createLayersTrack( RimFractureModelPlot* pl //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- -void RicNewFractureModelPlotFeature::createParametersTrack( RimFractureModelPlot* plot, - RimFractureModel* fractureModel, - RimEclipseCase* eclipseCase, - int timeStep, - const QString& trackTitle, - const PlotDefVector& curveConfigurations, - bool isPlotLogarithmic ) +void RicNewFractureModelPlotFeature::createParametersTrack( RimFractureModelPlot* plot, + RimFractureModel* fractureModel, + RimEclipseCase* eclipseCase, + int timeStep, + const QString& trackTitle, + const std::vector& propertyTypes, + bool isPlotLogarithmic ) { RimWellLogTrack* plotTrack = RicNewWellLogPlotFeatureImpl::createWellLogPlotTrack( false, trackTitle, plot ); plotTrack->setFormationCase( eclipseCase ); @@ -387,21 +363,19 @@ void RicNewFractureModelPlotFeature::createParametersTrack( RimFractureModelPlot caf::ColorTable colors = RiaColorTables::wellLogPlotPaletteColors(); int colorIndex = 0; - for ( auto curveConfig : curveConfigurations ) + for ( const RiaDefines::CurveProperty& propertyType : propertyTypes ) { - QString resultVariable = std::get<0>( curveConfig ); - RiaDefines::ResultCatType resultCategoryType = std::get<1>( curveConfig ); - RimFractureModelCurve::MissingValueStrategy missingValueStrategy = std::get<2>( curveConfig ); - bool fixedInitialTimeStep = std::get<3>( curveConfig ); - RiaDefines::CurveProperty curveProperty = std::get<4>( curveConfig ); + QString resultVariable = fractureModel->eclipseResultVariable( propertyType ); + RiaDefines::ResultCatType resultCategoryType = fractureModel->eclipseResultCategory( propertyType ); + // TODO: maybe improve? + bool fixedInitialTimeStep = ( propertyType == RiaDefines::CurveProperty::INITIAL_PRESSURE ); RimFractureModelCurve* curve = new RimFractureModelCurve; - curve->setCurveProperty( curveProperty ); + curve->setCurveProperty( propertyType ); curve->setFractureModel( fractureModel ); curve->setCase( eclipseCase ); curve->setEclipseResultVariable( resultVariable ); curve->setEclipseResultCategory( resultCategoryType ); - curve->setMissingValueStrategy( missingValueStrategy ); curve->setColor( colors.cycledColor3f( colorIndex ) ); curve->setLineStyle( RiuQwtPlotCurve::STYLE_SOLID ); curve->setLineThickness( 2 ); @@ -418,11 +392,6 @@ void RicNewFractureModelPlotFeature::createParametersTrack( RimFractureModelPlot curve->setCurrentTimeStep( timeStep ); } - if ( curveProperty == RiaDefines::CurveProperty::INITIAL_PRESSURE ) - { - curve->setBurdenStrategy( RimFractureModelCurve::BurdenStrategy::GRADIENT ); - } - plotTrack->addCurve( curve ); curve->loadDataAndUpdate( true ); diff --git a/ApplicationCode/Commands/CompletionCommands/RicNewFractureModelPlotFeature.h b/ApplicationCode/Commands/CompletionCommands/RicNewFractureModelPlotFeature.h index f318efdfd4..286c96664f 100644 --- a/ApplicationCode/Commands/CompletionCommands/RicNewFractureModelPlotFeature.h +++ b/ApplicationCode/Commands/CompletionCommands/RicNewFractureModelPlotFeature.h @@ -31,9 +31,6 @@ class RimFractureModelPlot; class RimFractureModelPlotCollection; class RimFractureModel; -typedef std::tuple PlotDef; -typedef std::vector PlotDefVector; - //================================================================================================== /// //================================================================================================== @@ -58,13 +55,13 @@ private: static void createLayersTrack( RimFractureModelPlot* plot, RimFractureModel* fractureModel, RimEclipseCase* eclipseCase ); - static void createParametersTrack( RimFractureModelPlot* plot, - RimFractureModel* fractureModel, - RimEclipseCase* eclipseCase, - int timeStep, - const QString& trackTitle, - const PlotDefVector& curveConfiguration, - bool isPlotLogarithmic ); + static void createParametersTrack( RimFractureModelPlot* plot, + RimFractureModel* fractureModel, + RimEclipseCase* eclipseCase, + int timeStep, + const QString& trackTitle, + const std::vector& propertyTypes, + bool isPlotLogarithmic ); static void createElasticPropertiesTrack( RimFractureModelPlot* plot, RimFractureModel* fractureModel, diff --git a/ApplicationCode/Commands/RicExportFractureModelPlotToFileFeature.cpp b/ApplicationCode/Commands/RicExportFractureModelPlotToFileFeature.cpp index cd22c749fb..3548d9b14c 100644 --- a/ApplicationCode/Commands/RicExportFractureModelPlotToFileFeature.cpp +++ b/ApplicationCode/Commands/RicExportFractureModelPlotToFileFeature.cpp @@ -59,7 +59,7 @@ void RicExportFractureModelPlotToFileFeature::onActionTriggered( bool isChecked if ( directoryPath.isEmpty() ) return; - RifFractureModelPlotExporter::writeToDirectory( fractureModelPlot, + RifFractureModelPlotExporter::writeToDirectory( fractureModelPlot->fractureModel(), fractureModelPlot->fractureModel()->useDetailedFluidLoss(), directoryPath ); diff --git a/ApplicationCode/FileInterface/RifFractureModelAsymmetricFrkExporter.cpp b/ApplicationCode/FileInterface/RifFractureModelAsymmetricFrkExporter.cpp index e3e7150da0..60465fe970 100644 --- a/ApplicationCode/FileInterface/RifFractureModelAsymmetricFrkExporter.cpp +++ b/ApplicationCode/FileInterface/RifFractureModelAsymmetricFrkExporter.cpp @@ -19,7 +19,6 @@ #include "RifFractureModelAsymmetricFrkExporter.h" #include "RimFractureModel.h" -#include "RimFractureModelPlot.h" #include #include @@ -27,14 +26,8 @@ //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- -bool RifFractureModelAsymmetricFrkExporter::writeToFile( RimFractureModelPlot* plot, const QString& filepath ) +bool RifFractureModelAsymmetricFrkExporter::writeToFile( RimFractureModel* fractureModel, const QString& filepath ) { - RimFractureModel* fractureModel = plot->fractureModel(); - if ( !fractureModel ) - { - return false; - } - QFile data( filepath ); if ( !data.open( QFile::WriteOnly | QFile::Truncate ) ) { diff --git a/ApplicationCode/FileInterface/RifFractureModelAsymmetricFrkExporter.h b/ApplicationCode/FileInterface/RifFractureModelAsymmetricFrkExporter.h index 74fa553fcb..75ad194ad6 100644 --- a/ApplicationCode/FileInterface/RifFractureModelAsymmetricFrkExporter.h +++ b/ApplicationCode/FileInterface/RifFractureModelAsymmetricFrkExporter.h @@ -18,7 +18,7 @@ #pragma once -class RimFractureModelPlot; +class RimFractureModel; class QString; class QTextStream; @@ -29,7 +29,7 @@ class QTextStream; class RifFractureModelAsymmetricFrkExporter { public: - static bool writeToFile( RimFractureModelPlot* plot, const QString& filepath ); + static bool writeToFile( RimFractureModel* fractureModel, const QString& filepath ); private: static void appendHeaderToStream( QTextStream& stream ); diff --git a/ApplicationCode/FileInterface/RifFractureModelDeviationFrkExporter.cpp b/ApplicationCode/FileInterface/RifFractureModelDeviationFrkExporter.cpp index 88d644b3c4..a4083399cf 100644 --- a/ApplicationCode/FileInterface/RifFractureModelDeviationFrkExporter.cpp +++ b/ApplicationCode/FileInterface/RifFractureModelDeviationFrkExporter.cpp @@ -19,7 +19,6 @@ #include "RifFractureModelDeviationFrkExporter.h" #include "RimFractureModel.h" -#include "RimFractureModelPlot.h" #include "RimWellPath.h" #include "RigWellPathGeometryExporter.h" @@ -30,14 +29,8 @@ //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- -bool RifFractureModelDeviationFrkExporter::writeToFile( RimFractureModelPlot* plot, const QString& filepath ) +bool RifFractureModelDeviationFrkExporter::writeToFile( RimFractureModel* fractureModel, const QString& filepath ) { - RimFractureModel* fractureModel = plot->fractureModel(); - if ( !fractureModel ) - { - return false; - } - RimWellPath* wellPath = fractureModel->wellPath(); if ( !wellPath ) { diff --git a/ApplicationCode/FileInterface/RifFractureModelDeviationFrkExporter.h b/ApplicationCode/FileInterface/RifFractureModelDeviationFrkExporter.h index 835999786f..ceca7d233b 100644 --- a/ApplicationCode/FileInterface/RifFractureModelDeviationFrkExporter.h +++ b/ApplicationCode/FileInterface/RifFractureModelDeviationFrkExporter.h @@ -20,7 +20,7 @@ #include -class RimFractureModelPlot; +class RimFractureModel; class QString; class QTextStream; @@ -30,7 +30,7 @@ class QTextStream; class RifFractureModelDeviationFrkExporter { public: - static bool writeToFile( RimFractureModelPlot* plot, const QString& filepath ); + static bool writeToFile( RimFractureModel* fractureModel, const QString& filepath ); static void fixupDepthValuesForExport( const std::vector& tvdValues, const std::vector& mdValues, diff --git a/ApplicationCode/FileInterface/RifFractureModelGeologicalFrkExporter.cpp b/ApplicationCode/FileInterface/RifFractureModelGeologicalFrkExporter.cpp index e9e693e6cd..e83eb931db 100644 --- a/ApplicationCode/FileInterface/RifFractureModelGeologicalFrkExporter.cpp +++ b/ApplicationCode/FileInterface/RifFractureModelGeologicalFrkExporter.cpp @@ -18,7 +18,8 @@ #include "RifFractureModelGeologicalFrkExporter.h" -#include "RimFractureModelPlot.h" +#include "RimFractureModel.h" +#include "RimFractureModelCalculator.h" #include #include @@ -26,9 +27,9 @@ //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- -bool RifFractureModelGeologicalFrkExporter::writeToFile( RimFractureModelPlot* plot, - bool useDetailedLoss, - const QString& filepath ) +bool RifFractureModelGeologicalFrkExporter::writeToFile( RimFractureModel* fractureModel, + bool useDetailedLoss, + const QString& filepath ) { std::vector labels; // TVD depth of top of zone (ft) @@ -90,24 +91,24 @@ bool RifFractureModelGeologicalFrkExporter::writeToFile( RimFractureModelPlot* p } std::map> values; - values["dpthlyr"] = plot->calculateTrueVerticalDepth(); - values["strs"] = plot->calculateStress(); - values["strsg"] = plot->calculateStressGradient(); - values["elyr"] = plot->calculateYoungsModulus(); - values["poissonr"] = plot->calculatePoissonsRatio(); - values["tuflyr"] = plot->calculateKIc(); - values["clyrc"] = plot->calculateFluidLossCoefficient(); - values["clyrs"] = plot->calculateSpurtLoss(); - values["pembed"] = plot->calculateProppandEmbedment(); - values["zoneResPres"] = plot->calculateReservoirPressure(); - values["zoneWaterSat"] = plot->calculateImmobileFluidSaturation(); - values["zonePorosity"] = plot->calculatePorosity(); - values["zoneHorizPerm"] = plot->calculateHorizontalPermeability(); - values["zoneVertPerm"] = plot->calculateVerticalPermeability(); - values["zoneTemp"] = plot->calculateTemperature(); - values["zoneRelPerm"] = plot->calculateRelativePermeabilityFactor(); - values["zonePoroElas"] = plot->calculatePoroElasticConstant(); - values["zoneThermalExp"] = plot->calculateThermalExpansionCoefficient(); + values["dpthlyr"] = fractureModel->calculator()->calculateTrueVerticalDepth(); + values["strs"] = fractureModel->calculator()->calculateStress(); + values["strsg"] = fractureModel->calculator()->calculateStressGradient(); + values["elyr"] = fractureModel->calculator()->calculateYoungsModulus(); + values["poissonr"] = fractureModel->calculator()->calculatePoissonsRatio(); + values["tuflyr"] = fractureModel->calculator()->calculateKIc(); + values["clyrc"] = fractureModel->calculator()->calculateFluidLossCoefficient(); + values["clyrs"] = fractureModel->calculator()->calculateSpurtLoss(); + values["pembed"] = fractureModel->calculator()->calculateProppandEmbedment(); + values["zoneResPres"] = fractureModel->calculator()->calculateReservoirPressure(); + values["zoneWaterSat"] = fractureModel->calculator()->calculateImmobileFluidSaturation(); + values["zonePorosity"] = fractureModel->calculator()->calculatePorosity(); + values["zoneHorizPerm"] = fractureModel->calculator()->calculateHorizontalPermeability(); + values["zoneVertPerm"] = fractureModel->calculator()->calculateVerticalPermeability(); + values["zoneTemp"] = fractureModel->calculator()->calculateTemperature(); + values["zoneRelPerm"] = fractureModel->calculator()->calculateRelativePermeabilityFactor(); + values["zonePoroElas"] = fractureModel->calculator()->calculatePoroElasticConstant(); + values["zoneThermalExp"] = fractureModel->calculator()->calculateThermalExpansionCoefficient(); QFile data( filepath ); if ( !data.open( QFile::WriteOnly | QFile::Truncate ) ) diff --git a/ApplicationCode/FileInterface/RifFractureModelGeologicalFrkExporter.h b/ApplicationCode/FileInterface/RifFractureModelGeologicalFrkExporter.h index 75778731d2..b0159fd133 100644 --- a/ApplicationCode/FileInterface/RifFractureModelGeologicalFrkExporter.h +++ b/ApplicationCode/FileInterface/RifFractureModelGeologicalFrkExporter.h @@ -20,7 +20,7 @@ #include -class RimFractureModelPlot; +class RimFractureModel; class QString; class QTextStream; @@ -30,7 +30,7 @@ class QTextStream; class RifFractureModelGeologicalFrkExporter { public: - static bool writeToFile( RimFractureModelPlot* plot, bool useDetailedFluidLoss, const QString& filepath ); + static bool writeToFile( RimFractureModel* plot, bool useDetailedFluidLoss, const QString& filepath ); private: static void appendHeaderToStream( QTextStream& stream ); diff --git a/ApplicationCode/FileInterface/RifFractureModelPerfsFrkExporter.cpp b/ApplicationCode/FileInterface/RifFractureModelPerfsFrkExporter.cpp index 54927d6b41..68ca3d9200 100644 --- a/ApplicationCode/FileInterface/RifFractureModelPerfsFrkExporter.cpp +++ b/ApplicationCode/FileInterface/RifFractureModelPerfsFrkExporter.cpp @@ -21,7 +21,6 @@ #include "RiaLogging.h" #include "RimFractureModel.h" -#include "RimFractureModelPlot.h" #include "RimWellPath.h" #include "RigWellPath.h" @@ -35,14 +34,8 @@ //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- -bool RifFractureModelPerfsFrkExporter::writeToFile( RimFractureModelPlot* plot, const QString& filepath ) +bool RifFractureModelPerfsFrkExporter::writeToFile( RimFractureModel* fractureModel, const QString& filepath ) { - RimFractureModel* fractureModel = plot->fractureModel(); - if ( !fractureModel ) - { - return false; - } - RimWellPath* wellPath = fractureModel->wellPath(); if ( !wellPath ) { diff --git a/ApplicationCode/FileInterface/RifFractureModelPerfsFrkExporter.h b/ApplicationCode/FileInterface/RifFractureModelPerfsFrkExporter.h index 72dce83eba..d04c5f5604 100644 --- a/ApplicationCode/FileInterface/RifFractureModelPerfsFrkExporter.h +++ b/ApplicationCode/FileInterface/RifFractureModelPerfsFrkExporter.h @@ -20,7 +20,7 @@ #include "cvfVector3.h" -class RimFractureModelPlot; +class RimFractureModel; class RimWellPath; class QString; @@ -32,7 +32,7 @@ class QTextStream; class RifFractureModelPerfsFrkExporter { public: - static bool writeToFile( RimFractureModelPlot* plot, const QString& filepath ); + static bool writeToFile( RimFractureModel* fractureModel, const QString& filepath ); private: static void appendHeaderToStream( QTextStream& stream ); diff --git a/ApplicationCode/FileInterface/RifFractureModelPlotExporter.cpp b/ApplicationCode/FileInterface/RifFractureModelPlotExporter.cpp index cd736dbea5..9de4435156 100644 --- a/ApplicationCode/FileInterface/RifFractureModelPlotExporter.cpp +++ b/ApplicationCode/FileInterface/RifFractureModelPlotExporter.cpp @@ -28,12 +28,14 @@ //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- -bool RifFractureModelPlotExporter::writeToDirectory( RimFractureModelPlot* plot, - bool useDetailedFluidLoss, - const QString& directoryPath ) +bool RifFractureModelPlotExporter::writeToDirectory( RimFractureModel* fractureModel, + bool useDetailedFluidLoss, + const QString& directoryPath ) { - return RifFractureModelGeologicalFrkExporter::writeToFile( plot, useDetailedFluidLoss, directoryPath + "/Geological.frk" ) && - RifFractureModelDeviationFrkExporter::writeToFile( plot, directoryPath + "/Deviation.frk" ) && - RifFractureModelPerfsFrkExporter::writeToFile( plot, directoryPath + "/Perfs.frk" ) && - RifFractureModelAsymmetricFrkExporter::writeToFile( plot, directoryPath + "/Asymmetric.frk" ); + return RifFractureModelGeologicalFrkExporter::writeToFile( fractureModel, + useDetailedFluidLoss, + directoryPath + "/Geological.frk" ) && + RifFractureModelDeviationFrkExporter::writeToFile( fractureModel, directoryPath + "/Deviation.frk" ) && + RifFractureModelPerfsFrkExporter::writeToFile( fractureModel, directoryPath + "/Perfs.frk" ) && + RifFractureModelAsymmetricFrkExporter::writeToFile( fractureModel, directoryPath + "/Asymmetric.frk" ); } diff --git a/ApplicationCode/FileInterface/RifFractureModelPlotExporter.h b/ApplicationCode/FileInterface/RifFractureModelPlotExporter.h index ca65caee89..6f5e9f4320 100644 --- a/ApplicationCode/FileInterface/RifFractureModelPlotExporter.h +++ b/ApplicationCode/FileInterface/RifFractureModelPlotExporter.h @@ -20,7 +20,7 @@ #include -class RimFractureModelPlot; +class RimFractureModel; class QString; class QTextStream; @@ -30,5 +30,5 @@ class QTextStream; class RifFractureModelPlotExporter { public: - static bool writeToDirectory( RimFractureModelPlot* plot, bool useDetailedFluidLoss, const QString& directoryPath ); + static bool writeToDirectory( RimFractureModel* fractureModel, bool useDetailedFluidLoss, const QString& directoryPath ); }; diff --git a/ApplicationCode/ProjectDataModel/CMakeLists_files.cmake b/ApplicationCode/ProjectDataModel/CMakeLists_files.cmake index da61be1501..c523183574 100644 --- a/ApplicationCode/ProjectDataModel/CMakeLists_files.cmake +++ b/ApplicationCode/ProjectDataModel/CMakeLists_files.cmake @@ -170,6 +170,12 @@ ${CMAKE_CURRENT_LIST_DIR}/RimAbstractPlotCollection.h ${CMAKE_CURRENT_LIST_DIR}/RimFractureModelPropertyCurve.h ${CMAKE_CURRENT_LIST_DIR}/RimFaciesProperties.h ${CMAKE_CURRENT_LIST_DIR}/RimNonNetLayers.h +${CMAKE_CURRENT_LIST_DIR}/RimFractureModelCalculator.h +${CMAKE_CURRENT_LIST_DIR}/RimFractureModelPropertyCalculator.h +${CMAKE_CURRENT_LIST_DIR}/RimFractureModelWellLogCalculator.h +${CMAKE_CURRENT_LIST_DIR}/RimFractureModelElasticPropertyCalculator.h +${CMAKE_CURRENT_LIST_DIR}/RimFractureModelLayerCalculator.h +${CMAKE_CURRENT_LIST_DIR}/RimFractureModelStressCalculator.h ) @@ -342,6 +348,11 @@ ${CMAKE_CURRENT_LIST_DIR}/RimLayerCurve.cpp ${CMAKE_CURRENT_LIST_DIR}/RimFractureModelStressCurve.cpp ${CMAKE_CURRENT_LIST_DIR}/RimFaciesProperties.cpp ${CMAKE_CURRENT_LIST_DIR}/RimNonNetLayers.cpp +${CMAKE_CURRENT_LIST_DIR}/RimFractureModelCalculator.cpp +${CMAKE_CURRENT_LIST_DIR}/RimFractureModelWellLogCalculator.cpp +${CMAKE_CURRENT_LIST_DIR}/RimFractureModelElasticPropertyCalculator.cpp +${CMAKE_CURRENT_LIST_DIR}/RimFractureModelLayerCalculator.cpp +${CMAKE_CURRENT_LIST_DIR}/RimFractureModelStressCalculator.cpp ) list(APPEND CODE_HEADER_FILES diff --git a/ApplicationCode/ProjectDataModel/Completions/RimFractureModel.cpp b/ApplicationCode/ProjectDataModel/Completions/RimFractureModel.cpp index b787a4e738..f3c6df9d44 100644 --- a/ApplicationCode/ProjectDataModel/Completions/RimFractureModel.cpp +++ b/ApplicationCode/ProjectDataModel/Completions/RimFractureModel.cpp @@ -39,9 +39,12 @@ #include "RimColorLegendItem.h" #include "RimCompletionTemplateCollection.h" #include "RimEclipseCase.h" +#include "RimEclipseResultDefinition.h" #include "RimEclipseView.h" +#include "RimFaciesProperties.h" #include "RimFaultInView.h" #include "RimFaultInViewCollection.h" +#include "RimFractureModelCalculator.h" #include "RimFractureModelPlot.h" #include "RimFractureModelTemplate.h" #include "RimFractureModelTemplateCollection.h" @@ -99,6 +102,25 @@ void caf::AppEnum::setUp() setDefault( RimFractureModel::FractureOrientation::TRANSVERSE_WELL_PATH ); } + +template <> +void caf::AppEnum::setUp() +{ + addItem( RimFractureModel::MissingValueStrategy::DEFAULT_VALUE, "DEFAULT_VALUE", "Default value" ); + addItem( RimFractureModel::MissingValueStrategy::LINEAR_INTERPOLATION, "LINEAR_INTERPOLATION", "Linear interpolation" ); + addItem( RimFractureModel::MissingValueStrategy::OTHER_CURVE_PROPERTY, "OTHER_CURVE_PROPERTY", "Other Curve Property" ); + + setDefault( RimFractureModel::MissingValueStrategy::DEFAULT_VALUE ); +} + +template <> +void caf::AppEnum::setUp() +{ + addItem( RimFractureModel::BurdenStrategy::DEFAULT_VALUE, "DEFAULT_VALUE", "Default value" ); + addItem( RimFractureModel::BurdenStrategy::GRADIENT, "GRADIENT", "Gradient" ); + + setDefault( RimFractureModel::BurdenStrategy::DEFAULT_VALUE ); +} }; // namespace caf //-------------------------------------------------------------------------------------------------- @@ -113,6 +135,8 @@ RimFractureModel::RimFractureModel() m_editFractureModelTemplate.uiCapability()->setUiEditorTypeName( caf::PdmUiToolButtonEditor::uiEditorTypeName() ); m_editFractureModelTemplate.uiCapability()->setUiLabelPosition( caf::PdmUiItemInfo::HIDDEN ); + CAF_PDM_InitScriptableField( &m_timeStep, "TimeStep", 0, "Time Step", "", "", "" ); + CAF_PDM_InitScriptableField( &m_MD, "MeasuredDepth", 0.0, "Measured Depth", "", "", "" ); m_MD.uiCapability()->setUiEditorTypeName( caf::PdmUiDoubleSliderEditor::uiEditorTypeName() ); @@ -199,6 +223,9 @@ RimFractureModel::RimFractureModel() "", "" ); + m_calculator = std::shared_ptr( new RimFractureModelCalculator ); + m_calculator->setFractureModel( this ); + setDeletable( true ); } @@ -782,6 +809,7 @@ void RimFractureModel::defineUiOrdering( QString uiConfigName, caf::PdmUiOrderin uiOrdering.add( &m_fractureModelTemplate, {true, 2, 1} ); uiOrdering.add( &m_editFractureModelTemplate, {false, 1, 0} ); + uiOrdering.add( &m_timeStep ); uiOrdering.add( &m_MD ); uiOrdering.add( &m_extractionType ); uiOrdering.add( &m_anchorPosition ); @@ -976,19 +1004,21 @@ double RimFractureModel::defaultPermeability() const //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- -double RimFractureModel::getDefaultForMissingValue( const QString& keyword ) const +double RimFractureModel::getDefaultForMissingValue( RiaDefines::CurveProperty curveProperty ) const { - if ( keyword == QString( "PORO" ) ) + if ( curveProperty == RiaDefines::CurveProperty::POROSITY ) { return defaultPorosity(); } - else if ( keyword == QString( "PERMX" ) || keyword == QString( "PERMZ" ) ) + else if ( curveProperty == RiaDefines::CurveProperty::PERMEABILITY_X || + curveProperty == RiaDefines::CurveProperty::PERMEABILITY_Z ) { return defaultPermeability(); } else { - RiaLogging::error( QString( "Missing default value for %1." ).arg( keyword ) ); + RiaLogging::error( QString( "Missing default value for %1." ) + .arg( caf::AppEnum( curveProperty ).uiText() ) ); return std::numeric_limits::infinity(); } } @@ -996,9 +1026,9 @@ double RimFractureModel::getDefaultForMissingValue( const QString& keyword ) con //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- -RiaDefines::CurveProperty RimFractureModel::getDefaultPropertyForMissingValues( const QString& keyword ) const +RiaDefines::CurveProperty RimFractureModel::getDefaultPropertyForMissingValues( RiaDefines::CurveProperty curveProperty ) const { - if ( keyword == QString( "PRESSURE" ) ) + if ( curveProperty == RiaDefines::CurveProperty::PRESSURE ) { return RiaDefines::CurveProperty::INITIAL_PRESSURE; } @@ -1009,19 +1039,27 @@ RiaDefines::CurveProperty RimFractureModel::getDefaultPropertyForMissingValues( //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- -double RimFractureModel::getDefaultForMissingOverburdenValue( const QString& keyword ) const +double RimFractureModel::getDefaultForMissingOverburdenValue( RiaDefines::CurveProperty curveProperty ) const { - if ( keyword == QString( "PORO" ) ) + if ( curveProperty == RiaDefines::CurveProperty::POROSITY ) { return defaultOverburdenPorosity(); } - else if ( keyword == QString( "PERMX" ) || keyword == QString( "PERMZ" ) ) + else if ( curveProperty == RiaDefines::CurveProperty::PERMEABILITY_X || + curveProperty == RiaDefines::CurveProperty::PERMEABILITY_Z ) { return defaultOverburdenPermeability(); } + else if ( curveProperty == RiaDefines::CurveProperty::FACIES ) + { + RimColorLegend* faciesColorLegend = getFaciesColorLegend(); + if ( !faciesColorLegend ) return std::numeric_limits::infinity(); + return findFaciesValue( *faciesColorLegend, overburdenFacies() ); + } else { - RiaLogging::error( QString( "Missing default overburden value for %1." ).arg( keyword ) ); + RiaLogging::error( QString( "Missing default overburden value for %1." ) + .arg( caf::AppEnum( curveProperty ).uiText() ) ); return std::numeric_limits::infinity(); } } @@ -1029,19 +1067,27 @@ double RimFractureModel::getDefaultForMissingOverburdenValue( const QString& key //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- -double RimFractureModel::getDefaultForMissingUnderburdenValue( const QString& keyword ) const +double RimFractureModel::getDefaultForMissingUnderburdenValue( RiaDefines::CurveProperty curveProperty ) const { - if ( keyword == QString( "PORO" ) ) + if ( curveProperty == RiaDefines::CurveProperty::POROSITY ) { return defaultUnderburdenPorosity(); } - else if ( keyword == QString( "PERMX" ) || keyword == QString( "PERMZ" ) ) + else if ( curveProperty == RiaDefines::CurveProperty::PERMEABILITY_X || + curveProperty == RiaDefines::CurveProperty::PERMEABILITY_Z ) { return defaultUnderburdenPermeability(); } + else if ( curveProperty == RiaDefines::CurveProperty::FACIES ) + { + RimColorLegend* faciesColorLegend = getFaciesColorLegend(); + if ( !faciesColorLegend ) return std::numeric_limits::infinity(); + return findFaciesValue( *faciesColorLegend, underburdenFacies() ); + } else { - RiaLogging::error( QString( "Missing default underburden value for %1." ).arg( keyword ) ); + RiaLogging::error( QString( "Missing default underburden value for %1." ) + .arg( caf::AppEnum( curveProperty ).uiText() ) ); return std::numeric_limits::infinity(); } } @@ -1049,9 +1095,10 @@ double RimFractureModel::getDefaultForMissingUnderburdenValue( const QString& ke //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- -double RimFractureModel::getOverburdenGradient( const QString& keyword ) const +double RimFractureModel::getOverburdenGradient( RiaDefines::CurveProperty curveProperty ) const { - if ( keyword == QString( "PRESSURE" ) ) + if ( curveProperty == RiaDefines::CurveProperty::PRESSURE || + curveProperty == RiaDefines::CurveProperty::INITIAL_PRESSURE ) { if ( !m_fractureModelTemplate ) { @@ -1061,7 +1108,8 @@ double RimFractureModel::getOverburdenGradient( const QString& keyword ) const } else { - RiaLogging::error( QString( "Missing overburden gradient for %1." ).arg( keyword ) ); + RiaLogging::error( QString( "Missing overburden gradient for %1." ) + .arg( caf::AppEnum( curveProperty ).uiText() ) ); return std::numeric_limits::infinity(); } } @@ -1069,9 +1117,10 @@ double RimFractureModel::getOverburdenGradient( const QString& keyword ) const //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- -double RimFractureModel::getUnderburdenGradient( const QString& keyword ) const +double RimFractureModel::getUnderburdenGradient( RiaDefines::CurveProperty curveProperty ) const { - if ( keyword == QString( "PRESSURE" ) ) + if ( curveProperty == RiaDefines::CurveProperty::PRESSURE || + curveProperty == RiaDefines::CurveProperty::INITIAL_PRESSURE ) { if ( !m_fractureModelTemplate ) { @@ -1082,7 +1131,8 @@ double RimFractureModel::getUnderburdenGradient( const QString& keyword ) const } else { - RiaLogging::error( QString( "Missing underburden gradient for %1." ).arg( keyword ) ); + RiaLogging::error( QString( "Missing underburden gradient for %1." ) + .arg( caf::AppEnum( curveProperty ).uiText() ) ); return std::numeric_limits::infinity(); } } @@ -1112,6 +1162,30 @@ double RimFractureModel::getDefaultValueForProperty( RiaDefines::CurveProperty c } } +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +RimFractureModel::MissingValueStrategy RimFractureModel::missingValueStrategy( RiaDefines::CurveProperty curveProperty ) const +{ + if ( curveProperty == RiaDefines::CurveProperty::INITIAL_PRESSURE ) + return RimFractureModel::MissingValueStrategy::LINEAR_INTERPOLATION; + else if ( curveProperty == RiaDefines::CurveProperty::PRESSURE ) + return RimFractureModel::MissingValueStrategy::OTHER_CURVE_PROPERTY; + else + return RimFractureModel::MissingValueStrategy::DEFAULT_VALUE; +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +RimFractureModel::BurdenStrategy RimFractureModel::burdenStrategy( RiaDefines::CurveProperty curveProperty ) const +{ + if ( curveProperty == RiaDefines::CurveProperty::INITIAL_PRESSURE ) + return RimFractureModel::BurdenStrategy::GRADIENT; + + return RimFractureModel::BurdenStrategy::DEFAULT_VALUE; +} + //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- @@ -1253,6 +1327,14 @@ void RimFractureModel::setMD( double md ) updateBarrierProperties(); } +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +int RimFractureModel::timeStep() const +{ + return m_timeStep; +} + //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- @@ -1445,3 +1527,104 @@ void RimFractureModel::showAllFaults() rimFault->showFault = true; } } + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +std::shared_ptr RimFractureModel::calculator() const +{ + return m_calculator; +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +RimEclipseCase* RimFractureModel::eclipseCase() const +{ + RimProject* proj = RimProject::current(); + if ( proj->eclipseCases().empty() ) return nullptr; + + return proj->eclipseCases()[0]; +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +RiaDefines::ResultCatType RimFractureModel::eclipseResultCategory( RiaDefines::CurveProperty curveProperty ) const +{ + if ( curveProperty == RiaDefines::CurveProperty::PRESSURE || + curveProperty == RiaDefines::CurveProperty::INITIAL_PRESSURE ) + { + return RiaDefines::ResultCatType::DYNAMIC_NATIVE; + } + else if ( curveProperty == RiaDefines::CurveProperty::FACIES ) + { + RimFaciesProperties* faciesProperties = m_fractureModelTemplate->faciesProperties(); + if ( !faciesProperties ) return RiaDefines::ResultCatType::STATIC_NATIVE; + + const RimEclipseResultDefinition* faciesDefinition = faciesProperties->faciesDefinition(); + if ( !faciesDefinition ) return RiaDefines::ResultCatType::STATIC_NATIVE; + + return faciesDefinition->resultType(); + } + else + { + return RiaDefines::ResultCatType::STATIC_NATIVE; + } +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +QString RimFractureModel::eclipseResultVariable( RiaDefines::CurveProperty curveProperty ) const +{ + if ( curveProperty == RiaDefines::CurveProperty::PRESSURE || + curveProperty == RiaDefines::CurveProperty::INITIAL_PRESSURE ) + return "PRESSURE"; + else if ( curveProperty == RiaDefines::CurveProperty::PERMEABILITY_X ) + return "PERMX"; + else if ( curveProperty == RiaDefines::CurveProperty::PERMEABILITY_Z ) + return "PERMZ"; + else if ( curveProperty == RiaDefines::CurveProperty::POROSITY ) + return "PORO"; + else if ( curveProperty == RiaDefines::CurveProperty::FACIES ) + { + if ( !m_fractureModelTemplate ) return ""; + + RimFaciesProperties* faciesProperties = m_fractureModelTemplate->faciesProperties(); + if ( !faciesProperties ) return ""; + + const RimEclipseResultDefinition* faciesDefinition = faciesProperties->faciesDefinition(); + if ( !faciesDefinition ) return ""; + + return faciesDefinition->resultVariable(); + } + else + return ""; +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +RimColorLegend* RimFractureModel::getFaciesColorLegend() const +{ + if ( !m_fractureModelTemplate ) return nullptr; + + RimFaciesProperties* faciesProperties = m_fractureModelTemplate->faciesProperties(); + if ( !faciesProperties ) return nullptr; + + return faciesProperties->colorLegend(); +} + +//------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +double RimFractureModel::findFaciesValue( const RimColorLegend& colorLegend, const QString& name ) +{ + for ( auto item : colorLegend.colorLegendItems() ) + { + if ( item->categoryName() == name ) return item->categoryValue(); + } + + return std::numeric_limits::infinity(); +} diff --git a/ApplicationCode/ProjectDataModel/Completions/RimFractureModel.h b/ApplicationCode/ProjectDataModel/Completions/RimFractureModel.h index 6cc8e5022c..10c4c3ed77 100644 --- a/ApplicationCode/ProjectDataModel/Completions/RimFractureModel.h +++ b/ApplicationCode/ProjectDataModel/Completions/RimFractureModel.h @@ -40,6 +40,8 @@ class RimUserDefinedPolylinesAnnotation; class RimFaciesProperties; class RimFractureModelTemplate; class RimTextAnnotation; +class RimFractureModelCalculator; +class RimColorLegend; //================================================================================================== /// @@ -63,11 +65,26 @@ public: AZIMUTH }; + enum class MissingValueStrategy + { + DEFAULT_VALUE, + LINEAR_INTERPOLATION, + OTHER_CURVE_PROPERTY + }; + + enum class BurdenStrategy + { + DEFAULT_VALUE, + GRADIENT + }; + RimFractureModel( void ); ~RimFractureModel( void ) override; void setMD( double md ); + int timeStep() const; + cvf::Vec3d anchorPosition() const; cvf::Vec3d thicknessDirection() const; @@ -130,18 +147,28 @@ public: double getDefaultValueForProperty( RiaDefines::CurveProperty ) const; bool hasDefaultValueForProperty( RiaDefines::CurveProperty ) const; - RiaDefines::CurveProperty getDefaultPropertyForMissingValues( const QString& keyword ) const; - double getDefaultForMissingOverburdenValue( const QString& keyword ) const; - double getDefaultForMissingUnderburdenValue( const QString& keyword ) const; - double getDefaultForMissingValue( const QString& keyword ) const; - double getOverburdenGradient( const QString& keyword ) const; - double getUnderburdenGradient( const QString& keyword ) const; + RiaDefines::CurveProperty getDefaultPropertyForMissingValues( RiaDefines::CurveProperty curveProperty ) const; + double getDefaultForMissingOverburdenValue( RiaDefines::CurveProperty curveProperty ) const; + double getDefaultForMissingUnderburdenValue( RiaDefines::CurveProperty curveProperty ) const; + double getDefaultForMissingValue( RiaDefines::CurveProperty curveProperty ) const; + double getOverburdenGradient( RiaDefines::CurveProperty curveProperty ) const; + double getUnderburdenGradient( RiaDefines::CurveProperty curveProperty ) const; void setFractureModelTemplate( RimFractureModelTemplate* fractureModelTemplate ); RimFractureModelTemplate* fractureModelTemplate() const; void updateReferringPlots(); + std::shared_ptr calculator() const; + + RimFractureModel::MissingValueStrategy missingValueStrategy( RiaDefines::CurveProperty curveProperty ) const; + RimFractureModel::BurdenStrategy burdenStrategy( RiaDefines::CurveProperty curveProperty ) const; + RimEclipseCase* eclipseCase() const; + RiaDefines::ResultCatType eclipseResultCategory( RiaDefines::CurveProperty curveProperty ) const; + QString eclipseResultVariable( RiaDefines::CurveProperty curveProperty ) const; + + static double findFaciesValue( const RimColorLegend& colorLegend, const QString& name ); + protected: void defineUiOrdering( QString uiConfigName, caf::PdmUiOrdering& uiOrdering ) override; QList calculateValueOptions( const caf::PdmFieldHandle* fieldNeedingOptions, @@ -185,8 +212,11 @@ private: void hideOtherFaults( const QString& targetFaultName ); void showAllFaults(); + RimColorLegend* getFaciesColorLegend() const; + protected: caf::PdmField m_MD; + caf::PdmField m_timeStep; caf::PdmField> m_extractionType; caf::PdmField m_anchorPosition; caf::PdmField m_thicknessDirection; @@ -216,4 +246,6 @@ protected: caf::PdmField m_barrierFaultName; caf::PdmField m_showOnlyBarrierFault; caf::PdmField m_showAllFaults; + + std::shared_ptr m_calculator; }; diff --git a/ApplicationCode/ProjectDataModel/RimElasticPropertiesCurve.cpp b/ApplicationCode/ProjectDataModel/RimElasticPropertiesCurve.cpp index 2912893d0d..afe4098419 100644 --- a/ApplicationCode/ProjectDataModel/RimElasticPropertiesCurve.cpp +++ b/ApplicationCode/ProjectDataModel/RimElasticPropertiesCurve.cpp @@ -19,44 +19,14 @@ #include "RimElasticPropertiesCurve.h" #include "RigEclipseCaseData.h" -#include "RigEclipseWellLogExtractor.h" -#include "RigElasticProperties.h" -#include "RigResultAccessorFactory.h" -#include "RigWellLogCurveData.h" -#include "RigWellPath.h" -#include "RimCase.h" -#include "RimColorLegend.h" -#include "RimColorLegendCollection.h" -#include "RimColorLegendItem.h" #include "RimEclipseCase.h" #include "RimEclipseResultDefinition.h" -#include "RimElasticProperties.h" -#include "RimFaciesProperties.h" #include "RimFractureModel.h" -#include "RimFractureModelPlot.h" -#include "RimFractureModelTemplate.h" +#include "RimFractureModelCalculator.h" #include "RimModeledWellPath.h" -#include "RimProject.h" -#include "RimTools.h" -#include "RimWellLogFile.h" -#include "RimWellLogPlot.h" -#include "RimWellLogTrack.h" -#include "RimWellPath.h" -#include "RimWellPathCollection.h" -#include "RimWellPlotTools.h" -#include "RiuQwtPlotCurve.h" -#include "RiuQwtPlotWidget.h" - -#include "RiaApplication.h" #include "RiaFractureDefines.h" -#include "RiaLogging.h" -#include "RiaPreferences.h" - -#include "cafPdmUiTreeOrdering.h" - -#include CAF_PDM_SOURCE_INIT( RimElasticPropertiesCurve, "ElasticPropertiesCurve" ); @@ -135,174 +105,17 @@ void RimElasticPropertiesCurve::performDataExtraction( bool* isUsingPseudoLength RimEclipseCase* eclipseCase = dynamic_cast( m_case.value() ); if ( eclipseCase && m_fractureModel ) { - RigEclipseWellLogExtractor eclExtractor( eclipseCase->eclipseCaseData(), - m_fractureModel->thicknessDirectionWellPath()->wellPathGeometry(), - "fracture model" ); - - measuredDepthValues = eclExtractor.cellIntersectionMDs(); - tvDepthValues = eclExtractor.cellIntersectionTVDs(); - rkbDiff = eclExtractor.wellPathData()->rkbDiff(); - - // Extract formation data - cvf::ref formationResultAccessor = RigResultAccessorFactory:: - createFromResultAddress( eclipseCase->eclipseCaseData(), - 0, - RiaDefines::PorosityModelType::MATRIX_MODEL, - 0, - RigEclipseResultAddress( RiaDefines::ResultCatType::FORMATION_NAMES, - RiaDefines::activeFormationNamesResultName() ) ); - if ( !formationResultAccessor.notNull() ) + bool isOk = m_fractureModel->calculator()->extractCurveData( curveProperty(), + m_timeStep, + values, + measuredDepthValues, + tvDepthValues, + rkbDiff ); + if ( !isOk ) { - RiaLogging::error( QString( "No formation result found." ) ); return; } - CurveSamplingPointData curveData = - RimWellLogTrack::curveSamplingPointData( &eclExtractor, formationResultAccessor.p() ); - - std::vector formationValues = curveData.data; - - std::vector formationNamesVector = RimWellLogTrack::formationNamesVector( eclipseCase ); - - RimFractureModelTemplate* fractureModelTemplate = m_fractureModel->fractureModelTemplate(); - if ( !fractureModelTemplate ) - { - RiaLogging::error( QString( "No fracture model template found" ) ); - return; - } - - RimFaciesProperties* faciesProperties = fractureModelTemplate->faciesProperties(); - if ( !faciesProperties ) - { - RiaLogging::error( QString( "No facies properties found when extracting elastic properties." ) ); - return; - } - - const RimEclipseResultDefinition* faciesDefinition = faciesProperties->faciesDefinition(); - - // Extract facies data - m_eclipseResultDefinition->setResultVariable( faciesDefinition->resultVariable() ); - m_eclipseResultDefinition->setResultType( faciesDefinition->resultType() ); - m_eclipseResultDefinition->setEclipseCase( eclipseCase ); - m_eclipseResultDefinition->loadResult(); - - cvf::ref faciesResultAccessor = - RigResultAccessorFactory::createFromResultDefinition( eclipseCase->eclipseCaseData(), - 0, - m_timeStep, - m_eclipseResultDefinition ); - - if ( !faciesResultAccessor.notNull() ) - { - RiaLogging::error( QString( "No facies result found." ) ); - return; - } - - std::vector faciesValues; - eclExtractor.curveData( faciesResultAccessor.p(), &faciesValues ); - - // Extract porosity data: get the porosity values from parent - RimFractureModelPlot* fractureModelPlot; - firstAncestorOrThisOfType( fractureModelPlot ); - if ( !fractureModelPlot ) - { - RiaLogging::error( QString( "No porosity data found when extracting elastic properties." ) ); - return; - } - - std::vector poroValues; - fractureModelPlot->getPorosityValues( poroValues ); - if ( poroValues.empty() ) - { - RiaLogging::error( QString( "Empty porosity data found when extracting elastic properties." ) ); - return; - } - - RimColorLegend* colorLegend = faciesProperties->colorLegend(); - if ( !colorLegend ) - { - RiaLogging::error( QString( "No color legend found when extracting elastic properties." ) ); - return; - } - - RimElasticProperties* elasticProperties = fractureModelTemplate->elasticProperties(); - if ( !elasticProperties ) - { - RiaLogging::error( QString( "No elastic properties found" ) ); - return; - } - - double overburdenHeight = m_fractureModel->overburdenHeight(); - if ( overburdenHeight > 0.0 ) - { - double defaultPoroValue = m_fractureModel->defaultOverburdenPorosity(); - QString overburdenFormation = m_fractureModel->overburdenFormation(); - QString overburdenFacies = m_fractureModel->overburdenFacies(); - - addOverburden( formationNamesVector, - formationValues, - faciesValues, - tvDepthValues, - measuredDepthValues, - overburdenHeight, - defaultPoroValue, - overburdenFormation, - findFaciesValue( *colorLegend, overburdenFacies ) ); - } - - double underburdenHeight = m_fractureModel->underburdenHeight(); - if ( underburdenHeight > 0.0 ) - { - double defaultPoroValue = m_fractureModel->defaultUnderburdenPorosity(); - QString underburdenFormation = m_fractureModel->underburdenFormation(); - QString underburdenFacies = m_fractureModel->underburdenFacies(); - - addUnderburden( formationNamesVector, - formationValues, - faciesValues, - tvDepthValues, - measuredDepthValues, - underburdenHeight, - defaultPoroValue, - underburdenFormation, - findFaciesValue( *colorLegend, underburdenFacies ) ); - } - - for ( size_t i = 0; i < tvDepthValues.size(); i++ ) - { - // Avoid using the field name in the match for now - QString fieldName = ""; - QString faciesName = findFaciesName( *colorLegend, faciesValues[i] ); - int idx = static_cast( formationValues[i] ); - QString formationName = formationNamesVector[idx]; - double porosity = poroValues[i]; - - FaciesKey faciesKey = std::make_tuple( fieldName, formationName, faciesName ); - if ( elasticProperties->hasPropertiesForFacies( faciesKey ) ) - { - if ( RimElasticProperties::isScalableProperty( curveProperty() ) ) - { - const RigElasticProperties& rigElasticProperties = elasticProperties->propertiesForFacies( faciesKey ); - double scale = elasticProperties->getPropertyScaling( formationName, faciesName, curveProperty() ); - double val = rigElasticProperties.getValueForPorosity( curveProperty(), porosity, scale ); - values.push_back( val ); - } - else if ( m_fractureModel->hasDefaultValueForProperty( curveProperty() ) ) - { - double val = m_fractureModel->getDefaultValueForProperty( curveProperty() ); - values.push_back( val ); - } - } - else - { - RiaLogging::error( QString( "Missing elastic properties. Field='%1', formation='%2', facies='%3'" ) - .arg( fieldName ) - .arg( formationName ) - .arg( faciesName ) ); - return; - } - } - RiaEclipseUnitTools::UnitSystem eclipseUnitsType = eclipseCase->eclipseCaseData()->unitsType(); if ( eclipseUnitsType == RiaEclipseUnitTools::UnitSystem::UNITS_FIELD ) { @@ -338,32 +151,6 @@ void RimElasticPropertiesCurve::performDataExtraction( bool* isUsingPseudoLength } } -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -QString RimElasticPropertiesCurve::findFaciesName( const RimColorLegend& colorLegend, double value ) -{ - for ( auto item : colorLegend.colorLegendItems() ) - { - if ( item->categoryValue() == static_cast( value ) ) return item->categoryName(); - } - - return "not found"; -} - -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -double RimElasticPropertiesCurve::findFaciesValue( const RimColorLegend& colorLegend, const QString& name ) -{ - for ( auto item : colorLegend.colorLegendItems() ) - { - if ( item->categoryName() == name ) return item->categoryValue(); - } - - return std::numeric_limits::infinity(); -} - //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- @@ -371,74 +158,3 @@ QString RimElasticPropertiesCurve::createCurveAutoName() { return caf::AppEnum::uiText( m_curveProperty() ); } - -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -void RimElasticPropertiesCurve::addOverburden( std::vector& formationNames, - std::vector& formationValues, - std::vector& faciesValues, - std::vector& tvDepthValues, - std::vector& measuredDepthValues, - double overburdenHeight, - double defaultPoroValue, - const QString& formationName, - double faciesValue ) -{ - if ( !faciesValues.empty() ) - { - // Prepend the new "fake" depth for start of overburden - double tvdTop = tvDepthValues[0]; - tvDepthValues.insert( tvDepthValues.begin(), tvdTop ); - tvDepthValues.insert( tvDepthValues.begin(), tvdTop - overburdenHeight ); - - // TODO: this is not always correct - double mdTop = measuredDepthValues[0]; - measuredDepthValues.insert( measuredDepthValues.begin(), mdTop ); - measuredDepthValues.insert( measuredDepthValues.begin(), mdTop - overburdenHeight ); - - formationNames.push_back( formationName ); - - formationValues.insert( formationValues.begin(), formationNames.size() - 1 ); - formationValues.insert( formationValues.begin(), formationNames.size() - 1 ); - - faciesValues.insert( faciesValues.begin(), faciesValue ); - faciesValues.insert( faciesValues.begin(), faciesValue ); - } -} - -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -void RimElasticPropertiesCurve::addUnderburden( std::vector& formationNames, - std::vector& formationValues, - std::vector& faciesValues, - std::vector& tvDepthValues, - std::vector& measuredDepthValues, - double underburdenHeight, - double defaultPoroValue, - const QString& formationName, - double faciesValue ) -{ - if ( !faciesValues.empty() ) - { - size_t lastIndex = tvDepthValues.size() - 1; - - double tvdBottom = tvDepthValues[lastIndex]; - tvDepthValues.push_back( tvdBottom ); - tvDepthValues.push_back( tvdBottom + underburdenHeight ); - - // TODO: this is not always correct - double mdBottom = measuredDepthValues[lastIndex]; - measuredDepthValues.push_back( mdBottom ); - measuredDepthValues.push_back( mdBottom + underburdenHeight ); - - formationNames.push_back( formationName ); - - formationValues.push_back( formationNames.size() - 1 ); - formationValues.push_back( formationNames.size() - 1 ); - - faciesValues.push_back( faciesValue ); - faciesValues.push_back( faciesValue ); - } -} diff --git a/ApplicationCode/ProjectDataModel/RimFractureModelCalculator.cpp b/ApplicationCode/ProjectDataModel/RimFractureModelCalculator.cpp new file mode 100644 index 0000000000..6c2157cece --- /dev/null +++ b/ApplicationCode/ProjectDataModel/RimFractureModelCalculator.cpp @@ -0,0 +1,571 @@ +///////////////////////////////////////////////////////////////////////////////// +// +// Copyright (C) 2020- 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 +// for more details. +// +///////////////////////////////////////////////////////////////////////////////// +#include "RimFractureModelCalculator.h" + +#include "RiaDefines.h" +#include "RiaFractureModelDefines.h" +#include "RiaLogging.h" + +#include "RigEclipseCaseData.h" + +#include "RimEclipseResultDefinition.h" +#include "RimFractureModel.h" +#include "RimFractureModelCalculator.h" +#include "RimFractureModelElasticPropertyCalculator.h" +#include "RimFractureModelLayerCalculator.h" +#include "RimFractureModelPropertyCalculator.h" +#include "RimFractureModelStressCalculator.h" +#include "RimFractureModelWellLogCalculator.h" + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +RimFractureModelCalculator::RimFractureModelCalculator() +{ + m_resultCalculators.push_back( + std::unique_ptr( new RimFractureModelWellLogCalculator( *this ) ) ); + m_resultCalculators.push_back( + std::unique_ptr( new RimFractureModelElasticPropertyCalculator( this ) ) ); + m_resultCalculators.push_back( + std::unique_ptr( new RimFractureModelLayerCalculator( this ) ) ); + m_resultCalculators.push_back( + std::unique_ptr( new RimFractureModelStressCalculator( this ) ) ); +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +void RimFractureModelCalculator::setFractureModel( RimFractureModel* fractureModel ) +{ + m_fractureModel = fractureModel; +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +RimFractureModel* RimFractureModelCalculator::fractureModel() +{ + return m_fractureModel; +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +bool RimFractureModelCalculator::extractCurveData( RiaDefines::CurveProperty curveProperty, + int timeStep, + std::vector& values, + std::vector& measuredDepthValues, + std::vector& tvDepthValues, + double& rkbDiff ) const +{ + for ( const auto& calculator : m_resultCalculators ) + { + if ( calculator->isMatching( curveProperty ) ) + { + return calculator + ->calculate( curveProperty, m_fractureModel, timeStep, values, measuredDepthValues, tvDepthValues, rkbDiff ); + } + } + + return false; +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +std::vector RimFractureModelCalculator::extractValues( RiaDefines::CurveProperty curveProperty, int timeStep ) const +{ + std::vector values; + std::vector measuredDepthValues; + std::vector tvDepthValues; + double rkbDiff; + + extractCurveData( curveProperty, timeStep, values, measuredDepthValues, tvDepthValues, rkbDiff ); + + return values; +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +void RimFractureModelCalculator::calculateLayers( std::vector>& layerBoundaryDepths, + std::vector>& layerBoundaryIndexes ) const +{ + std::vector layerValues; + std::vector measuredDepthValues; + std::vector depths; + double rkbDiff; + + extractCurveData( RiaDefines::CurveProperty::LAYERS, + m_fractureModel->timeStep(), + layerValues, + measuredDepthValues, + depths, + rkbDiff ); + + size_t startIndex = 0; + for ( size_t i = 0; i < depths.size(); i++ ) + { + if ( startIndex != i && ( layerValues[startIndex] != layerValues[i] || i == depths.size() - 1 ) ) + { + layerBoundaryDepths.push_back( std::make_pair( depths[startIndex], depths[i] ) ); + layerBoundaryIndexes.push_back( std::make_pair( startIndex, i ) ); + startIndex = i; + } + } +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +double RimFractureModelCalculator::findValueAtTopOfLayer( const std::vector& values, + const std::vector>& layerBoundaryIndexes, + size_t layerNo ) +{ + size_t index = layerBoundaryIndexes[layerNo].first; + return values.at( index ); +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +double RimFractureModelCalculator::findValueAtBottomOfLayer( const std::vector& values, + const std::vector>& layerBoundaryIndexes, + size_t layerNo ) +{ + size_t index = layerBoundaryIndexes[layerNo].second; + return values.at( index ); +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +void RimFractureModelCalculator::computeAverageByLayer( const std::vector>& layerBoundaryIndexes, + const std::vector& inputVector, + std::vector& result ) +{ + for ( auto boundaryIndex : layerBoundaryIndexes ) + { + double sum = 0.0; + int nValues = 0; + for ( size_t i = boundaryIndex.first; i < boundaryIndex.second; i++ ) + { + sum += inputVector[i]; + nValues++; + } + result.push_back( sum / nValues ); + } +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +void RimFractureModelCalculator::extractTopOfLayerValues( const std::vector>& layerBoundaryIndexes, + const std::vector& inputVector, + std::vector& result ) +{ + for ( size_t i = 0; i < layerBoundaryIndexes.size(); i++ ) + { + result.push_back( findValueAtTopOfLayer( inputVector, layerBoundaryIndexes, i ) ); + } +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +std::vector RimFractureModelCalculator::calculateTrueVerticalDepth() const +{ + std::vector> layerBoundaryDepths; + std::vector> layerBoundaryIndexes; + + calculateLayers( layerBoundaryDepths, layerBoundaryIndexes ); + + std::vector tvdTopZone; + for ( auto p : layerBoundaryDepths ) + { + double depthInFeet = RiaEclipseUnitTools::meterToFeet( p.first ); + tvdTopZone.push_back( depthInFeet ); + } + + return tvdTopZone; +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +std::vector + RimFractureModelCalculator::findCurveAndComputeLayeredAverage( RiaDefines::CurveProperty curveProperty ) const +{ + std::vector> layerBoundaryDepths; + std::vector> layerBoundaryIndexes; + calculateLayers( layerBoundaryDepths, layerBoundaryIndexes ); + + std::vector values = extractValues( curveProperty, m_fractureModel->timeStep() ); + std::vector result; + computeAverageByLayer( layerBoundaryIndexes, values, result ); + return result; +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +std::vector RimFractureModelCalculator::findCurveAndComputeTopOfLayer( RiaDefines::CurveProperty curveProperty ) const +{ + std::vector> layerBoundaryDepths; + std::vector> layerBoundaryIndexes; + calculateLayers( layerBoundaryDepths, layerBoundaryIndexes ); + + std::vector values = extractValues( curveProperty, m_fractureModel->timeStep() ); + std::vector result; + extractTopOfLayerValues( layerBoundaryIndexes, values, result ); + return result; +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +std::vector RimFractureModelCalculator::calculatePorosity() const +{ + return findCurveAndComputeLayeredAverage( RiaDefines::CurveProperty::POROSITY ); +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +std::vector RimFractureModelCalculator::calculateReservoirPressure() const +{ + std::vector pressureBar = findCurveAndComputeTopOfLayer( RiaDefines::CurveProperty::PRESSURE ); + + std::vector pressurePsi; + for ( double p : pressureBar ) + { + pressurePsi.push_back( RiaEclipseUnitTools::barToPsi( p ) ); + } + + return pressurePsi; +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +std::vector RimFractureModelCalculator::calculateHorizontalPermeability() const +{ + return findCurveAndComputeLayeredAverage( RiaDefines::CurveProperty::PERMEABILITY_X ); +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +std::vector RimFractureModelCalculator::calculateVerticalPermeability() const +{ + return findCurveAndComputeLayeredAverage( RiaDefines::CurveProperty::PERMEABILITY_Z ); +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +std::vector RimFractureModelCalculator::calculateStress() const +{ + std::vector stress; + std::vector stressGradients; + std::vector initialStress; + calculateStressWithGradients( stress, stressGradients, initialStress ); + return stress; +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +std::vector RimFractureModelCalculator::calculateInitialStress() const +{ + std::vector stress; + std::vector stressGradients; + std::vector initialStress; + calculateStressWithGradients( stress, stressGradients, initialStress ); + return initialStress; +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +bool RimFractureModelCalculator::calculateStressWithGradients( std::vector& stress, + std::vector& stressGradients, + std::vector& initialStress ) const +{ + // Reference stress + const double verticalStressRef = m_fractureModel->verticalStress(); + const double verticalStressGradientRef = m_fractureModel->verticalStressGradient(); + const double stressDepthRef = m_fractureModel->stressDepth(); + + std::vector> layerBoundaryDepths; + std::vector> layerBoundaryIndexes; + calculateLayers( layerBoundaryDepths, layerBoundaryIndexes ); + + int timeStep = m_fractureModel->timeStep(); + + // Biot coefficient + std::vector biotData = extractValues( RiaDefines::CurveProperty::BIOT_COEFFICIENT, timeStep ); + + // K0 + std::vector k0Data = extractValues( RiaDefines::CurveProperty::K0, timeStep ); + + // Pressure at the give time step + std::vector timeStepPressureData = extractValues( RiaDefines::CurveProperty::PRESSURE, timeStep ); + + // Initial pressure + std::vector initialPressureData = extractValues( RiaDefines::CurveProperty::INITIAL_PRESSURE, timeStep ); + + // Poissons ratio + std::vector poissonsRatioData = extractValues( RiaDefines::CurveProperty::POISSONS_RATIO, timeStep ); + + // Check that we have data from all curves + if ( biotData.empty() || k0Data.empty() || timeStepPressureData.empty() || initialPressureData.empty() || + poissonsRatioData.empty() ) + { + return false; + } + + std::vector stressForGradients; + std::vector pressureForGradients; + std::vector depthForGradients; + + // Calculate the stress + for ( size_t i = 0; i < layerBoundaryDepths.size(); i++ ) + { + double depthTopOfZone = layerBoundaryDepths[i].first; + double depthBottomOfZone = layerBoundaryDepths[i].second; + + // Data from curves at the top zone depth + double k0 = findValueAtTopOfLayer( k0Data, layerBoundaryIndexes, i ); + double biot = findValueAtTopOfLayer( biotData, layerBoundaryIndexes, i ); + double poissonsRatio = findValueAtTopOfLayer( poissonsRatioData, layerBoundaryIndexes, i ); + double initialPressure = findValueAtTopOfLayer( initialPressureData, layerBoundaryIndexes, i ); + double timeStepPressure = findValueAtTopOfLayer( timeStepPressureData, layerBoundaryIndexes, i ); + + // Vertical stress + // Use difference between reference depth and depth of top of zone + double depthDiff = depthTopOfZone - stressDepthRef; + double Sv = verticalStressRef + verticalStressGradientRef * depthDiff; + + double Sh_init = k0 * Sv + initialPressure * ( 1.0 - k0 ); + double pressureDiff = timeStepPressure - initialPressure; + + // Vertical stress diff assumed to be zero + double Sv_diff = 0.0; + double deltaHorizontalStress = poissonsRatio / ( 1.0 - poissonsRatio ) * ( Sv_diff - biot * pressureDiff ) + + ( biot * pressureDiff ); + + double depletionStress = Sh_init + deltaHorizontalStress; + stress.push_back( RiaEclipseUnitTools::barToPsi( depletionStress ) ); + + initialStress.push_back( RiaEclipseUnitTools::barToPsi( Sh_init ) ); + + // Cache some results for the gradients calculation + stressForGradients.push_back( Sv ); + pressureForGradients.push_back( initialPressure ); + depthForGradients.push_back( depthTopOfZone ); + + if ( i == layerBoundaryDepths.size() - 1 ) + { + // Use the bottom of the last layer to compute gradient for last layer + double bottomInitialPressure = findValueAtBottomOfLayer( initialPressureData, layerBoundaryIndexes, i ); + + double bottomDepthDiff = depthBottomOfZone - stressDepthRef; + double bottomSv = verticalStressRef + verticalStressGradientRef * bottomDepthDiff; + stressForGradients.push_back( bottomSv ); + pressureForGradients.push_back( bottomInitialPressure ); + depthForGradients.push_back( depthBottomOfZone ); + } + } + + assert( stressForGradients.size() == layerBoundaryDepths.size() + 1 ); + assert( pressureForGradients.size() == layerBoundaryDepths.size() + 1 ); + assert( depthForGradients.size() == layerBoundaryDepths.size() + 1 ); + + // Second pass to calculate the stress gradients + for ( size_t i = 0; i < layerBoundaryDepths.size(); i++ ) + { + double diffStress = stressForGradients[i + 1] - stressForGradients[i]; + double diffPressure = pressureForGradients[i + 1] - pressureForGradients[i]; + double diffDepth = depthForGradients[i + 1] - depthForGradients[i]; + double k0 = findValueAtTopOfLayer( k0Data, layerBoundaryIndexes, i ); + double stressGradient = ( diffStress * k0 + diffPressure * ( 1.0 - k0 ) ) / diffDepth; + stressGradients.push_back( RiaEclipseUnitTools::barPerMeterToPsiPerFeet( stressGradient ) ); + } + + return true; +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +std::vector RimFractureModelCalculator::calculateStressGradient() const +{ + std::vector stress; + std::vector stressGradients; + std::vector initialStress; + calculateStressWithGradients( stress, stressGradients, initialStress ); + return stressGradients; +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +void RimFractureModelCalculator::calculateTemperature( std::vector& temperatures ) const +{ + // Reference temperature. Unit: degrees celsius + const double referenceTemperature = m_fractureModel->referenceTemperature(); + + // Reference temperature gradient. Unit: degrees Celsius per meter + const double referenceTemperatureGradient = m_fractureModel->referenceTemperatureGradient(); + + // Reference depth for temperature. Unit: meter. + const double referenceTemperatureDepth = m_fractureModel->referenceTemperatureDepth(); + + std::vector> layerBoundaryDepths; + std::vector> layerBoundaryIndexes; + calculateLayers( layerBoundaryDepths, layerBoundaryIndexes ); + + // Calculate the temperatures + for ( size_t i = 0; i < layerBoundaryDepths.size(); i++ ) + { + double depthTopOfZone = layerBoundaryDepths[i].first; + + // Use difference between reference depth and depth of top of zone + double depthDiff = depthTopOfZone - referenceTemperatureDepth; + double temperature = referenceTemperature + referenceTemperatureGradient * depthDiff; + + temperatures.push_back( temperature ); + } +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +std::vector RimFractureModelCalculator::calculateYoungsModulus() const +{ + std::vector valuesGPa = findCurveAndComputeLayeredAverage( RiaDefines::CurveProperty::YOUNGS_MODULUS ); + std::vector valuesMMpsi; + for ( auto value : valuesGPa ) + { + valuesMMpsi.push_back( value * 0.14503773773 ); + } + + return valuesMMpsi; +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +std::vector RimFractureModelCalculator::calculatePoissonsRatio() const +{ + return findCurveAndComputeLayeredAverage( RiaDefines::CurveProperty::POISSONS_RATIO ); +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +std::vector RimFractureModelCalculator::calculateKIc() const +{ + return findCurveAndComputeLayeredAverage( RiaDefines::CurveProperty::K_IC ); +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +std::vector RimFractureModelCalculator::calculateFluidLossCoefficient() const +{ + return findCurveAndComputeLayeredAverage( RiaDefines::CurveProperty::FLUID_LOSS_COEFFICIENT ); +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +std::vector RimFractureModelCalculator::calculateSpurtLoss() const +{ + return findCurveAndComputeLayeredAverage( RiaDefines::CurveProperty::SPURT_LOSS ); +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +std::vector RimFractureModelCalculator::calculateProppandEmbedment() const +{ + return findCurveAndComputeLayeredAverage( RiaDefines::CurveProperty::PROPPANT_EMBEDMENT ); +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +std::vector RimFractureModelCalculator::calculateImmobileFluidSaturation() const +{ + return findCurveAndComputeLayeredAverage( RiaDefines::CurveProperty::IMMOBILE_FLUID_SATURATION ); +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +std::vector RimFractureModelCalculator::calculateTemperature() const +{ + std::vector temperaturesCelsius; + calculateTemperature( temperaturesCelsius ); + + // Convert to Fahrenheit + std::vector temperaturesFahrenheit; + for ( double t : temperaturesCelsius ) + { + temperaturesFahrenheit.push_back( t * 1.8 + 32.0 ); + } + + return temperaturesFahrenheit; +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +std::vector RimFractureModelCalculator::calculateRelativePermeabilityFactor() const +{ + return findCurveAndComputeLayeredAverage( RiaDefines::CurveProperty::RELATIVE_PERMEABILITY_FACTOR ); +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +std::vector RimFractureModelCalculator::calculatePoroElasticConstant() const +{ + return findCurveAndComputeLayeredAverage( RiaDefines::CurveProperty::PORO_ELASTIC_CONSTANT ); +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +std::vector RimFractureModelCalculator::calculateThermalExpansionCoefficient() const +{ + // SI unit is 1/Celsius + std::vector coefficientCelsius = + findCurveAndComputeLayeredAverage( RiaDefines::CurveProperty::THERMAL_EXPANSION_COEFFICIENT ); + + // Field unit is 1/Fahrenheit + std::vector coefficientFahrenheit; + for ( double c : coefficientCelsius ) + { + coefficientFahrenheit.push_back( c / 1.8 ); + } + + return coefficientFahrenheit; +} diff --git a/ApplicationCode/ProjectDataModel/RimFractureModelCalculator.h b/ApplicationCode/ProjectDataModel/RimFractureModelCalculator.h new file mode 100644 index 0000000000..2f7ad6959f --- /dev/null +++ b/ApplicationCode/ProjectDataModel/RimFractureModelCalculator.h @@ -0,0 +1,96 @@ +///////////////////////////////////////////////////////////////////////////////// +// +// Copyright (C) 2020- 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 +// for more details. +// +///////////////////////////////////////////////////////////////////////////////// +#pragma once + +#include "RiaFractureModelDefines.h" + +#include "RimFractureModelPropertyCalculator.h" + +#include +#include + +class RimFractureModel; + +class RimFractureModelCalculator +{ +public: + RimFractureModelCalculator(); + + void setFractureModel( RimFractureModel* fractureModel ); + RimFractureModel* fractureModel(); + + bool extractCurveData( RiaDefines::CurveProperty curveProperty, + int timeStep, + std::vector& values, + std::vector& measuredDepthValues, + std::vector& tvDepthValues, + double& rkbDiff ) const; + + std::vector extractValues( RiaDefines::CurveProperty curveProperty, int timeStep ) const; + + std::vector calculateTrueVerticalDepth() const; + std::vector calculatePorosity() const; + std::vector calculateVerticalPermeability() const; + std::vector calculateHorizontalPermeability() const; + std::vector calculateReservoirPressure() const; + std::vector calculateStress() const; + std::vector calculateInitialStress() const; + std::vector calculateStressGradient() const; + std::vector calculateYoungsModulus() const; + std::vector calculatePoissonsRatio() const; + std::vector calculateKIc() const; + std::vector calculateFluidLossCoefficient() const; + std::vector calculateSpurtLoss() const; + std::vector calculateProppandEmbedment() const; + + std::vector calculateImmobileFluidSaturation() const; + std::vector calculateTemperature() const; + std::vector calculateRelativePermeabilityFactor() const; + std::vector calculatePoroElasticConstant() const; + std::vector calculateThermalExpansionCoefficient() const; + + void calculateTemperature( std::vector& temperatures ) const; + +protected: + std::vector findCurveAndComputeLayeredAverage( RiaDefines::CurveProperty curveProperty ) const; + std::vector findCurveXValuesByProperty( RiaDefines::CurveProperty curveProperty ) const; + std::vector findCurveAndComputeTopOfLayer( RiaDefines::CurveProperty curveProperty ) const; + + void calculateLayers( std::vector>& layerBoundaryDepths, + std::vector>& layerBoundaryIndexes ) const; + bool calculateStressWithGradients( std::vector& stress, + std::vector& stressGradients, + std::vector& initialStress ) const; + + static double findValueAtTopOfLayer( const std::vector& values, + const std::vector>& layerBoundaryIndexes, + size_t layerNo ); + static double findValueAtBottomOfLayer( const std::vector& values, + const std::vector>& layerBoundaryIndexes, + size_t layerNo ); + static void computeAverageByLayer( const std::vector>& layerBoundaryIndexes, + const std::vector& inputVector, + std::vector& result ); + static void extractTopOfLayerValues( const std::vector>& layerBoundaryIndexes, + const std::vector& inputVector, + std::vector& result ); + +private: + RimFractureModel* m_fractureModel; + std::vector> m_resultCalculators; +}; diff --git a/ApplicationCode/ProjectDataModel/RimFractureModelCurve.cpp b/ApplicationCode/ProjectDataModel/RimFractureModelCurve.cpp index 38f706e0cc..cd0d7eb48e 100644 --- a/ApplicationCode/ProjectDataModel/RimFractureModelCurve.cpp +++ b/ApplicationCode/ProjectDataModel/RimFractureModelCurve.cpp @@ -18,67 +18,28 @@ #include "RimFractureModelCurve.h" -#include "RiaFractureModelDefines.h" -#include "RigEclipseCaseData.h" -#include "RigEclipseWellLogExtractor.h" -#include "RigResultAccessorFactory.h" -#include "RigWellLogCurveData.h" -#include "RigWellPath.h" - -#include "RimCase.h" -#include "RimEclipseCase.h" -#include "RimEclipseInputProperty.h" -#include "RimEclipseInputPropertyCollection.h" -#include "RimEclipseResultDefinition.h" -#include "RimFractureModel.h" -#include "RimFractureModelPlot.h" -#include "RimModeledWellPath.h" -#include "RimTools.h" -#include "RimWellLogFile.h" -#include "RimWellLogPlot.h" -#include "RimWellLogTrack.h" -#include "RimWellPath.h" -#include "RimWellPathCollection.h" -#include "RimWellPlotTools.h" - -#include "RiuQwtPlotCurve.h" -#include "RiuQwtPlotWidget.h" - #include "RiaApplication.h" +#include "RiaFractureModelDefines.h" #include "RiaInterpolationTools.h" #include "RiaLogging.h" #include "RiaPreferences.h" +#include "RigEclipseCaseData.h" + +#include "RimEclipseCase.h" +#include "RimEclipseResultDefinition.h" +#include "RimFractureModel.h" +#include "RimFractureModelCalculator.h" +#include "RimFractureModelPlot.h" +#include "RimModeledWellPath.h" + +#include "RiuQwtPlotCurve.h" +#include "RiuQwtPlotWidget.h" + #include "cafPdmUiTreeOrdering.h" CAF_PDM_SOURCE_INIT( RimFractureModelCurve, "FractureModelCurve" ); -namespace caf -{ -template <> -void caf::AppEnum::setUp() -{ - addItem( RimFractureModelCurve::MissingValueStrategy::DEFAULT_VALUE, "DEFAULT_VALUE", "Default value" ); - addItem( RimFractureModelCurve::MissingValueStrategy::LINEAR_INTERPOLATION, - "LINEAR_INTERPOLATION", - "Linear interpolation" ); - addItem( RimFractureModelCurve::MissingValueStrategy::OTHER_CURVE_PROPERTY, - "OTHER_CURVE_PROPERTY", - "Other Curve Property" ); - - setDefault( RimFractureModelCurve::MissingValueStrategy::DEFAULT_VALUE ); -} - -template <> -void caf::AppEnum::setUp() -{ - addItem( RimFractureModelCurve::BurdenStrategy::DEFAULT_VALUE, "DEFAULT_VALUE", "Default value" ); - addItem( RimFractureModelCurve::BurdenStrategy::GRADIENT, "GRADIENT", "Gradient" ); - - setDefault( RimFractureModelCurve::BurdenStrategy::DEFAULT_VALUE ); -} -}; // namespace caf - //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- @@ -90,16 +51,6 @@ RimFractureModelCurve::RimFractureModelCurve() m_fractureModel.uiCapability()->setUiTreeChildrenHidden( true ); m_fractureModel.uiCapability()->setUiHidden( true ); - caf::AppEnum defaultValue = - RimFractureModelCurve::MissingValueStrategy::DEFAULT_VALUE; - CAF_PDM_InitField( &m_missingValueStrategy, "MissingValueStrategy", defaultValue, "Missing Value Strategy", "", "", "" ); - m_missingValueStrategy.uiCapability()->setUiHidden( true ); - - caf::AppEnum defaultBurdenValue = - RimFractureModelCurve::BurdenStrategy::DEFAULT_VALUE; - CAF_PDM_InitField( &m_burdenStrategy, "BurdenStrategy", defaultBurdenValue, "Burden Strategy", "", "", "" ); - m_burdenStrategy.uiCapability()->setUiHidden( true ); - CAF_PDM_InitFieldNoDefault( &m_curveProperty, "CurveProperty", "Curve Property", "", "", "" ); m_curveProperty.uiCapability()->setUiHidden( true ); @@ -181,144 +132,15 @@ void RimFractureModelCurve::performDataExtraction( bool* isUsingPseudoLength ) RimEclipseCase* eclipseCase = dynamic_cast( m_case.value() ); if ( eclipseCase && m_fractureModel ) { - RigEclipseWellLogExtractor eclExtractor( eclipseCase->eclipseCaseData(), - m_fractureModel->thicknessDirectionWellPath()->wellPathGeometry(), - "fracture model" ); - - measuredDepthValues = eclExtractor.cellIntersectionMDs(); - tvDepthValues = eclExtractor.cellIntersectionTVDs(); - rkbDiff = eclExtractor.wellPathData()->rkbDiff(); - - m_eclipseResultDefinition->setEclipseCase( eclipseCase ); - - m_eclipseResultDefinition->loadResult(); - - cvf::ref resAcc = - RigResultAccessorFactory::createFromResultDefinition( eclipseCase->eclipseCaseData(), - 0, - m_timeStep, - m_eclipseResultDefinition ); - - if ( resAcc.notNull() ) + bool isOk = m_fractureModel->calculator()->extractCurveData( curveProperty(), + m_timeStep, + values, + measuredDepthValues, + tvDepthValues, + rkbDiff ); + if ( !isOk ) { - eclExtractor.curveData( resAcc.p(), &values ); - } - else - { - RiaLogging::error( QString( "No result found for %1" ).arg( m_eclipseResultDefinition()->resultVariable() ) ); - } - - double overburdenHeight = m_fractureModel->overburdenHeight(); - if ( overburdenHeight > 0.0 ) - { - addOverburden( tvDepthValues, measuredDepthValues, values ); - } - - double underburdenHeight = m_fractureModel->underburdenHeight(); - if ( underburdenHeight > 0.0 ) - { - addUnderburden( tvDepthValues, measuredDepthValues, values ); - } - - if ( hasMissingValues( values ) ) - { - if ( m_missingValueStrategy() == RimFractureModelCurve::MissingValueStrategy::DEFAULT_VALUE ) - { - // Try to locate a backup accessor (e.g. PORO_1 for PORO) - cvf::ref backupResAcc = - findMissingValuesAccessor( eclipseCase->eclipseCaseData(), - eclipseCase->inputPropertyCollection(), - 0, - m_timeStep, - m_eclipseResultDefinition() ); - - if ( backupResAcc.notNull() ) - { - RiaLogging::info( QString( "Reading missing values from input properties for %1." ) - .arg( m_eclipseResultDefinition()->resultVariable() ) ); - std::vector replacementValues; - eclExtractor.curveData( backupResAcc.p(), &replacementValues ); - - double overburdenHeight = m_fractureModel->overburdenHeight(); - if ( overburdenHeight > 0.0 ) - { - double defaultOverburdenValue = std::numeric_limits::infinity(); - if ( m_burdenStrategy() == RimFractureModelCurve::BurdenStrategy::DEFAULT_VALUE ) - { - defaultOverburdenValue = m_fractureModel->getDefaultForMissingOverburdenValue( - m_eclipseResultDefinition()->resultVariable() ); - } - - replacementValues.insert( replacementValues.begin(), defaultOverburdenValue ); - replacementValues.insert( replacementValues.begin(), defaultOverburdenValue ); - } - - double underburdenHeight = m_fractureModel->underburdenHeight(); - if ( underburdenHeight > 0.0 ) - { - double defaultUnderburdenValue = std::numeric_limits::infinity(); - if ( m_burdenStrategy() == RimFractureModelCurve::BurdenStrategy::DEFAULT_VALUE ) - { - defaultUnderburdenValue = m_fractureModel->getDefaultForMissingUnderburdenValue( - m_eclipseResultDefinition()->resultVariable() ); - } - - replacementValues.push_back( defaultUnderburdenValue ); - replacementValues.push_back( defaultUnderburdenValue ); - } - - replaceMissingValues( values, replacementValues ); - } - - // If the backup accessor is not found, or does not provide all the missing values: - // use default value from the fracture model - if ( !backupResAcc.notNull() || hasMissingValues( values ) ) - { - RiaLogging::info( - QString( "Using default value for %1" ).arg( m_eclipseResultDefinition()->resultVariable() ) ); - - double defaultValue = - m_fractureModel->getDefaultForMissingValue( m_eclipseResultDefinition.value()->resultVariable() ); - - replaceMissingValues( values, defaultValue ); - } - } - else if ( m_missingValueStrategy() == RimFractureModelCurve::MissingValueStrategy::LINEAR_INTERPOLATION ) - { - RiaLogging::info( - QString( "Interpolating missing values for %1" ).arg( m_eclipseResultDefinition()->resultVariable() ) ); - RiaInterpolationTools::interpolateMissingValues( measuredDepthValues, values ); - } - else - { - // Get the missing data from other curve - RimFractureModelPlot* fractureModelPlot; - firstAncestorOrThisOfType( fractureModelPlot ); - if ( !fractureModelPlot ) - { - RiaLogging::error( QString( "No replacement data found for fracture model curve." ) ); - return; - } - - RiaDefines::CurveProperty replacementProperty = m_fractureModel->getDefaultPropertyForMissingValues( - m_eclipseResultDefinition.value()->resultVariable() ); - - std::vector initialValues; - fractureModelPlot->getCurvePropertyValues( replacementProperty, initialValues ); - if ( initialValues.empty() ) - { - RiaLogging::error( QString( "Empty replacement data found for fracture model curve." ) ); - return; - } - - if ( values.size() != initialValues.size() ) - { - RiaLogging::error( QString( "Inconsistent state." ) ); - return; - } - - replaceMissingValues( values, initialValues ); - } + return; } RiaEclipseUnitTools::UnitSystem eclipseUnitsType = eclipseCase->eclipseCaseData()->unitsType(); @@ -355,187 +177,3 @@ void RimFractureModelCurve::performDataExtraction( bool* isUsingPseudoLength ) } } } - -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -void RimFractureModelCurve::setBurdenStrategy( BurdenStrategy strategy ) -{ - m_burdenStrategy = strategy; -} - -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -void RimFractureModelCurve::setMissingValueStrategy( MissingValueStrategy strategy ) -{ - m_missingValueStrategy = strategy; -} - -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -bool RimFractureModelCurve::hasMissingValues( const std::vector& values ) -{ - for ( double v : values ) - { - if ( v == std::numeric_limits::infinity() ) - { - return true; - } - } - - return false; -} - -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -void RimFractureModelCurve::replaceMissingValues( std::vector& values, double defaultValue ) -{ - for ( double& v : values ) - { - if ( v == std::numeric_limits::infinity() ) - { - v = defaultValue; - } - } -} - -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -void RimFractureModelCurve::replaceMissingValues( std::vector& values, const std::vector& replacementValues ) -{ - CVF_ASSERT( values.size() == replacementValues.size() ); - for ( size_t i = 0; i < values.size(); i++ ) - { - if ( values[i] == std::numeric_limits::infinity() ) - { - values[i] = replacementValues[i]; - } - } -} - -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -cvf::ref - RimFractureModelCurve::findMissingValuesAccessor( RigEclipseCaseData* caseData, - RimEclipseInputPropertyCollection* inputPropertyCollection, - int gridIndex, - int timeStepIndex, - RimEclipseResultDefinition* eclipseResultDefinition ) -{ - QString resultName = eclipseResultDefinition->resultVariable(); - - for ( RimEclipseInputProperty* inputProperty : inputPropertyCollection->inputProperties() ) - { - // Look for input properties starting with the same name as result definition - if ( inputProperty && inputProperty->resultName().startsWith( resultName ) ) - { - RiaLogging::info( - QString( "Found missing values result for %1: %2" ).arg( resultName ).arg( inputProperty->resultName() ) ); - - RigEclipseResultAddress resultAddress( RiaDefines::ResultCatType::INPUT_PROPERTY, inputProperty->resultName() ); - caseData->results( eclipseResultDefinition->porosityModel() )->ensureKnownResultLoaded( resultAddress ); - cvf::ref resAcc = - RigResultAccessorFactory::createFromResultAddress( caseData, - gridIndex, - eclipseResultDefinition->porosityModel(), - timeStepIndex, - resultAddress ); - - return resAcc; - } - } - - return nullptr; -} - -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -void RimFractureModelCurve::addOverburden( std::vector& tvDepthValues, - std::vector& measuredDepthValues, - std::vector& values ) const -{ - if ( !values.empty() ) - { - double overburdenHeight = m_fractureModel->overburdenHeight(); - double tvdOverburdenBottom = tvDepthValues[0]; - double tvdOverburdenTop = tvdOverburdenBottom - overburdenHeight; - - double overburdenTopValue = std::numeric_limits::infinity(); - double overburdenBottomValue = std::numeric_limits::infinity(); - if ( m_burdenStrategy() == RimFractureModelCurve::BurdenStrategy::DEFAULT_VALUE ) - { - overburdenTopValue = - m_fractureModel->getDefaultForMissingOverburdenValue( m_eclipseResultDefinition()->resultVariable() ); - overburdenBottomValue = overburdenTopValue; - } - else - { - double gradient = m_fractureModel->getOverburdenGradient( m_eclipseResultDefinition()->resultVariable() ); - overburdenBottomValue = values[0]; - overburdenTopValue = overburdenBottomValue + gradient * -overburdenHeight; - } - - // Prepend the new "fake" depth for start of overburden - tvDepthValues.insert( tvDepthValues.begin(), tvdOverburdenBottom ); - tvDepthValues.insert( tvDepthValues.begin(), tvdOverburdenTop ); - - // TODO: this is not always correct - double mdTop = measuredDepthValues[0]; - measuredDepthValues.insert( measuredDepthValues.begin(), mdTop ); - measuredDepthValues.insert( measuredDepthValues.begin(), mdTop - overburdenHeight ); - - values.insert( values.begin(), overburdenBottomValue ); - values.insert( values.begin(), overburdenTopValue ); - } -} - -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -void RimFractureModelCurve::addUnderburden( std::vector& tvDepthValues, - std::vector& measuredDepthValues, - std::vector& values ) const -{ - if ( !values.empty() ) - { - size_t lastIndex = tvDepthValues.size() - 1; - - double underburdenHeight = m_fractureModel->underburdenHeight(); - double tvdUnderburdenTop = tvDepthValues[lastIndex]; - double tvdUnderburdenBottom = tvdUnderburdenTop + underburdenHeight; - - double underburdenTopValue = std::numeric_limits::infinity(); - double underburdenBottomValue = std::numeric_limits::infinity(); - if ( m_burdenStrategy() == RimFractureModelCurve::BurdenStrategy::DEFAULT_VALUE ) - { - underburdenTopValue = - m_fractureModel->getDefaultForMissingUnderburdenValue( m_eclipseResultDefinition()->resultVariable() ); - underburdenBottomValue = underburdenTopValue; - } - else - { - double gradient = m_fractureModel->getUnderburdenGradient( m_eclipseResultDefinition()->resultVariable() ); - underburdenTopValue = values[lastIndex]; - underburdenBottomValue = underburdenTopValue + gradient * underburdenHeight; - } - - // Append the new "fake" depth for start of underburden - tvDepthValues.push_back( tvdUnderburdenTop ); - tvDepthValues.push_back( tvdUnderburdenBottom ); - - // Append the new "fake" md - // TODO: check if this is correct??? - double mdBottom = measuredDepthValues[lastIndex]; - measuredDepthValues.push_back( mdBottom ); - measuredDepthValues.push_back( mdBottom + underburdenHeight ); - - values.push_back( underburdenTopValue ); - values.push_back( underburdenBottomValue ); - } -} diff --git a/ApplicationCode/ProjectDataModel/RimFractureModelCurve.h b/ApplicationCode/ProjectDataModel/RimFractureModelCurve.h index 78ee37eae9..42f219e4e9 100644 --- a/ApplicationCode/ProjectDataModel/RimFractureModelCurve.h +++ b/ApplicationCode/ProjectDataModel/RimFractureModelCurve.h @@ -21,19 +21,12 @@ #include "RimFractureModelPropertyCurve.h" #include "RimWellLogExtractionCurve.h" -#include "RiuQwtSymbol.h" - #include "cafPdmField.h" #include "cafPdmPtrField.h" #include -class RimWellPath; -class RimWellMeasurement; class RimFractureModel; -class RimEclipseInputPropertyCollection; -class RigEclipseCaseData; -class RigResultAccessor; //================================================================================================== /// @@ -43,19 +36,6 @@ class RimFractureModelCurve : public RimWellLogExtractionCurve, public RimFractu CAF_PDM_HEADER_INIT; public: - enum class MissingValueStrategy - { - DEFAULT_VALUE, - LINEAR_INTERPOLATION, - OTHER_CURVE_PROPERTY - }; - - enum class BurdenStrategy - { - DEFAULT_VALUE, - GRADIENT - }; - RimFractureModelCurve(); ~RimFractureModelCurve() override; @@ -63,10 +43,6 @@ public: void setEclipseResultCategory( RiaDefines::ResultCatType catType ); - void setMissingValueStrategy( MissingValueStrategy strategy ); - - void setBurdenStrategy( BurdenStrategy strategy ); - void setCurveProperty( RiaDefines::CurveProperty ) override; RiaDefines::CurveProperty curveProperty() const override; @@ -75,25 +51,6 @@ protected: void performDataExtraction( bool* isUsingPseudoLength ) override; - static bool hasMissingValues( const std::vector& values ); - static void replaceMissingValues( std::vector& values, double defaultValue ); - static void replaceMissingValues( std::vector& values, const std::vector& replacementValues ); - cvf::ref findMissingValuesAccessor( RigEclipseCaseData* caseData, - RimEclipseInputPropertyCollection* inputPropertyCollection, - int gridIndex, - int timeStepIndex, - RimEclipseResultDefinition* eclipseResultDefinition ); - - void addOverburden( std::vector& tvDepthValues, - std::vector& measuredDepthValues, - std::vector& values ) const; - - void addUnderburden( std::vector& tvDepthValues, - std::vector& measuredDepthValues, - std::vector& values ) const; - caf::PdmPtrField m_fractureModel; - caf::PdmField> m_missingValueStrategy; - caf::PdmField> m_burdenStrategy; caf::PdmField> m_curveProperty; }; diff --git a/ApplicationCode/ProjectDataModel/RimFractureModelElasticPropertyCalculator.cpp b/ApplicationCode/ProjectDataModel/RimFractureModelElasticPropertyCalculator.cpp new file mode 100644 index 0000000000..4064b5bf74 --- /dev/null +++ b/ApplicationCode/ProjectDataModel/RimFractureModelElasticPropertyCalculator.cpp @@ -0,0 +1,302 @@ +///////////////////////////////////////////////////////////////////////////////// +// +// Copyright (C) 2020- 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 +// for more details. +// +///////////////////////////////////////////////////////////////////////////////// +#include "RimFractureModelElasticPropertyCalculator.h" + +#include "RiaDefines.h" +#include "RiaFractureModelDefines.h" +#include "RiaInterpolationTools.h" +#include "RiaLogging.h" + +#include "RigEclipseCaseData.h" +#include "RigEclipseWellLogExtractor.h" +#include "RigElasticProperties.h" +#include "RigResultAccessor.h" +#include "RigResultAccessorFactory.h" +#include "RigWellLogCurveData.h" +#include "RigWellPath.h" + +#include "RimCase.h" +#include "RimColorLegend.h" +#include "RimColorLegendItem.h" +#include "RimEclipseCase.h" +#include "RimEclipseInputProperty.h" +#include "RimEclipseInputPropertyCollection.h" +#include "RimEclipseResultDefinition.h" +#include "RimElasticProperties.h" +#include "RimFaciesProperties.h" +#include "RimFractureModel.h" +#include "RimFractureModelCalculator.h" +#include "RimFractureModelElasticPropertyCalculator.h" +#include "RimFractureModelPlot.h" +#include "RimFractureModelTemplate.h" +#include "RimLayerCurve.h" +#include "RimModeledWellPath.h" +#include "RimTools.h" +#include "RimWellLogFile.h" +#include "RimWellLogPlot.h" +#include "RimWellLogTrack.h" +#include "RimWellPath.h" +#include "RimWellPathCollection.h" +#include "RimWellPlotTools.h" + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +RimFractureModelElasticPropertyCalculator::RimFractureModelElasticPropertyCalculator( + RimFractureModelCalculator* fractureModelCalculator ) + : m_fractureModelCalculator( fractureModelCalculator ) +{ +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +bool RimFractureModelElasticPropertyCalculator::isMatching( RiaDefines::CurveProperty curveProperty ) const +{ + std::vector matching = {RiaDefines::CurveProperty::YOUNGS_MODULUS, + RiaDefines::CurveProperty::POISSONS_RATIO, + RiaDefines::CurveProperty::BIOT_COEFFICIENT, + RiaDefines::CurveProperty::K0, + RiaDefines::CurveProperty::K_IC, + RiaDefines::CurveProperty::PROPPANT_EMBEDMENT, + RiaDefines::CurveProperty::FLUID_LOSS_COEFFICIENT, + RiaDefines::CurveProperty::SPURT_LOSS, + RiaDefines::CurveProperty::RELATIVE_PERMEABILITY_FACTOR, + RiaDefines::CurveProperty::PORO_ELASTIC_CONSTANT, + RiaDefines::CurveProperty::THERMAL_EXPANSION_COEFFICIENT, + RiaDefines::CurveProperty::IMMOBILE_FLUID_SATURATION}; + + return std::find( matching.begin(), matching.end(), curveProperty ) != matching.end(); +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +bool RimFractureModelElasticPropertyCalculator::calculate( RiaDefines::CurveProperty curveProperty, + const RimFractureModel* fractureModel, + int timeStep, + std::vector& values, + std::vector& measuredDepthValues, + std::vector& tvDepthValues, + double& rkbDiff ) const +{ + RimEclipseCase* eclipseCase = fractureModel->eclipseCase(); + if ( !eclipseCase ) + { + return false; + } + + RigEclipseWellLogExtractor eclExtractor( eclipseCase->eclipseCaseData(), + fractureModel->thicknessDirectionWellPath()->wellPathGeometry(), + "fracture model" ); + + measuredDepthValues = eclExtractor.cellIntersectionMDs(); + tvDepthValues = eclExtractor.cellIntersectionTVDs(); + rkbDiff = eclExtractor.wellPathData()->rkbDiff(); + + // Extract formation data + cvf::ref formationResultAccessor = + RigResultAccessorFactory::createFromResultAddress( eclipseCase->eclipseCaseData(), + 0, + RiaDefines::PorosityModelType::MATRIX_MODEL, + 0, + RigEclipseResultAddress( RiaDefines::ResultCatType::FORMATION_NAMES, + RiaDefines::activeFormationNamesResultName() ) ); + if ( !formationResultAccessor.notNull() ) + { + RiaLogging::error( QString( "No formation result found." ) ); + return false; + } + + CurveSamplingPointData curveData = + RimWellLogTrack::curveSamplingPointData( &eclExtractor, formationResultAccessor.p() ); + + std::vector formationValues = curveData.data; + + std::vector formationNamesVector = RimWellLogTrack::formationNamesVector( eclipseCase ); + + RimFractureModelTemplate* fractureModelTemplate = fractureModel->fractureModelTemplate(); + if ( !fractureModelTemplate ) + { + RiaLogging::error( QString( "No fracture model template found" ) ); + return false; + } + + RimFaciesProperties* faciesProperties = fractureModelTemplate->faciesProperties(); + if ( !faciesProperties ) + { + RiaLogging::error( QString( "No facies properties found when extracting elastic properties." ) ); + return false; + } + + RimColorLegend* colorLegend = faciesProperties->colorLegend(); + if ( !colorLegend ) + { + RiaLogging::error( QString( "No color legend found when extracting elastic properties." ) ); + return false; + } + + RimElasticProperties* elasticProperties = fractureModelTemplate->elasticProperties(); + if ( !elasticProperties ) + { + RiaLogging::error( QString( "No elastic properties found" ) ); + return false; + } + + std::vector faciesValues = + m_fractureModelCalculator->extractValues( RiaDefines::CurveProperty::FACIES, timeStep ); + std::vector poroValues = + m_fractureModelCalculator->extractValues( RiaDefines::CurveProperty::POROSITY, timeStep ); + + double overburdenHeight = fractureModel->overburdenHeight(); + if ( overburdenHeight > 0.0 ) + { + QString overburdenFormation = fractureModel->overburdenFormation(); + addOverburden( formationNamesVector, + formationValues, + tvDepthValues, + measuredDepthValues, + overburdenHeight, + overburdenFormation ); + } + + double underburdenHeight = fractureModel->underburdenHeight(); + if ( underburdenHeight > 0.0 ) + { + QString underburdenFormation = fractureModel->underburdenFormation(); + addUnderburden( formationNamesVector, + formationValues, + tvDepthValues, + measuredDepthValues, + underburdenHeight, + underburdenFormation ); + } + + CAF_ASSERT( tvDepthValues.size() == faciesValues.size() ); + CAF_ASSERT( tvDepthValues.size() == poroValues.size() ); + CAF_ASSERT( tvDepthValues.size() == formationValues.size() ); + + for ( size_t i = 0; i < tvDepthValues.size(); i++ ) + { + // Avoid using the field name in the match for now + QString fieldName = ""; + QString faciesName = findFaciesName( *colorLegend, faciesValues[i] ); + int idx = static_cast( formationValues[i] ); + QString formationName = formationNamesVector[idx]; + double porosity = poroValues[i]; + + FaciesKey faciesKey = std::make_tuple( fieldName, formationName, faciesName ); + if ( elasticProperties->hasPropertiesForFacies( faciesKey ) ) + { + if ( RimElasticProperties::isScalableProperty( curveProperty ) ) + { + const RigElasticProperties& rigElasticProperties = elasticProperties->propertiesForFacies( faciesKey ); + double scale = elasticProperties->getPropertyScaling( formationName, faciesName, curveProperty ); + double val = rigElasticProperties.getValueForPorosity( curveProperty, porosity, scale ); + values.push_back( val ); + } + else if ( fractureModel->hasDefaultValueForProperty( curveProperty ) ) + { + double val = fractureModel->getDefaultValueForProperty( curveProperty ); + values.push_back( val ); + } + } + else + { + RiaLogging::error( QString( "Missing elastic properties. Field='%1', formation='%2', facies='%3'" ) + .arg( fieldName ) + .arg( formationName ) + .arg( faciesName ) ); + return false; + } + } + + return true; +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +QString RimFractureModelElasticPropertyCalculator::findFaciesName( const RimColorLegend& colorLegend, double value ) +{ + for ( auto item : colorLegend.colorLegendItems() ) + { + if ( item->categoryValue() == static_cast( value ) ) return item->categoryName(); + } + + return "not found"; +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +void RimFractureModelElasticPropertyCalculator::addOverburden( std::vector& formationNames, + std::vector& formationValues, + std::vector& tvDepthValues, + std::vector& measuredDepthValues, + double overburdenHeight, + const QString& formationName ) +{ + if ( !tvDepthValues.empty() ) + { + // Prepend the new "fake" depth for start of overburden + double tvdTop = tvDepthValues[0]; + tvDepthValues.insert( tvDepthValues.begin(), tvdTop ); + tvDepthValues.insert( tvDepthValues.begin(), tvdTop - overburdenHeight ); + + // TODO: this is not always correct + double mdTop = measuredDepthValues[0]; + measuredDepthValues.insert( measuredDepthValues.begin(), mdTop ); + measuredDepthValues.insert( measuredDepthValues.begin(), mdTop - overburdenHeight ); + + formationNames.push_back( formationName ); + + formationValues.insert( formationValues.begin(), formationNames.size() - 1 ); + formationValues.insert( formationValues.begin(), formationNames.size() - 1 ); + } +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +void RimFractureModelElasticPropertyCalculator::addUnderburden( std::vector& formationNames, + std::vector& formationValues, + std::vector& tvDepthValues, + std::vector& measuredDepthValues, + double underburdenHeight, + const QString& formationName ) +{ + if ( !tvDepthValues.empty() ) + { + size_t lastIndex = tvDepthValues.size() - 1; + + double tvdBottom = tvDepthValues[lastIndex]; + tvDepthValues.push_back( tvdBottom ); + tvDepthValues.push_back( tvdBottom + underburdenHeight ); + + // TODO: this is not always correct + double mdBottom = measuredDepthValues[lastIndex]; + measuredDepthValues.push_back( mdBottom ); + measuredDepthValues.push_back( mdBottom + underburdenHeight ); + + formationNames.push_back( formationName ); + + formationValues.push_back( formationNames.size() - 1 ); + formationValues.push_back( formationNames.size() - 1 ); + } +} diff --git a/ApplicationCode/ProjectDataModel/RimFractureModelElasticPropertyCalculator.h b/ApplicationCode/ProjectDataModel/RimFractureModelElasticPropertyCalculator.h new file mode 100644 index 0000000000..84fd3353fa --- /dev/null +++ b/ApplicationCode/ProjectDataModel/RimFractureModelElasticPropertyCalculator.h @@ -0,0 +1,66 @@ +///////////////////////////////////////////////////////////////////////////////// +// +// Copyright (C) 2020- 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 +// for more details. +// +///////////////////////////////////////////////////////////////////////////////// +#pragma once + +#include "RimFractureModelPropertyCalculator.h" + +#include "RiaFractureModelDefines.h" + +#include + +class RimFractureModelCalculator; +class RimFractureModel; +class RimColorLegend; + +class QString; + +class RimFractureModelElasticPropertyCalculator : public RimFractureModelPropertyCalculator +{ +public: + RimFractureModelElasticPropertyCalculator( RimFractureModelCalculator* calculator ); + + bool calculate( RiaDefines::CurveProperty curveProperty, + const RimFractureModel* fractureModel, + int timeStep, + std::vector& values, + std::vector& measuredDepthValues, + std::vector& tvDepthValues, + double& rkbDiff ) const override; + + bool isMatching( RiaDefines::CurveProperty curveProperty ) const override; + +protected: + static void addOverburden( std::vector& formationNames, + std::vector& formationValues, + std::vector& tvDepthValues, + std::vector& measuredDepthValues, + double overburdenHeight, + const QString& formationName ); + + static void addUnderburden( std::vector& formationNames, + std::vector& formationValues, + std::vector& tvDepthValues, + std::vector& measuredDepthValues, + double underburdenHeight, + const QString& formationName ); + + static QString findFaciesName( const RimColorLegend& colorLegend, double value ); + +private: + RimFractureModelCalculator* m_fractureModelCalculator; +}; diff --git a/ApplicationCode/ProjectDataModel/RimFractureModelLayerCalculator.cpp b/ApplicationCode/ProjectDataModel/RimFractureModelLayerCalculator.cpp new file mode 100644 index 0000000000..4f69e07404 --- /dev/null +++ b/ApplicationCode/ProjectDataModel/RimFractureModelLayerCalculator.cpp @@ -0,0 +1,157 @@ +///////////////////////////////////////////////////////////////////////////////// +// +// Copyright (C) 2020- 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 +// for more details. +// +///////////////////////////////////////////////////////////////////////////////// +#include "RimFractureModelLayerCalculator.h" + +#include "RiaDefines.h" +#include "RiaFractureModelDefines.h" +#include "RiaInterpolationTools.h" +#include "RiaLogging.h" + +#include "RigEclipseCaseData.h" +#include "RigEclipseWellLogExtractor.h" +#include "RigElasticProperties.h" +#include "RigResultAccessor.h" +#include "RigResultAccessorFactory.h" +#include "RigWellLogCurveData.h" +#include "RigWellPath.h" + +#include "RimCase.h" +#include "RimColorLegend.h" +#include "RimColorLegendItem.h" +#include "RimEclipseCase.h" +#include "RimEclipseInputProperty.h" +#include "RimEclipseInputPropertyCollection.h" +#include "RimEclipseResultDefinition.h" +#include "RimElasticProperties.h" +#include "RimFaciesProperties.h" +#include "RimFractureModel.h" +#include "RimFractureModelCalculator.h" +#include "RimFractureModelLayerCalculator.h" +#include "RimFractureModelPlot.h" +#include "RimFractureModelTemplate.h" +#include "RimLayerCurve.h" +#include "RimModeledWellPath.h" +#include "RimTools.h" +#include "RimWellLogFile.h" +#include "RimWellLogPlot.h" +#include "RimWellLogTrack.h" +#include "RimWellPath.h" +#include "RimWellPathCollection.h" +#include "RimWellPlotTools.h" + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +RimFractureModelLayerCalculator::RimFractureModelLayerCalculator( RimFractureModelCalculator* fractureModelCalculator ) + : m_fractureModelCalculator( fractureModelCalculator ) +{ +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +bool RimFractureModelLayerCalculator::isMatching( RiaDefines::CurveProperty curveProperty ) const +{ + return curveProperty == RiaDefines::CurveProperty::LAYERS; +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +bool RimFractureModelLayerCalculator::calculate( RiaDefines::CurveProperty curveProperty, + const RimFractureModel* fractureModel, + int timeStep, + std::vector& values, + std::vector& measuredDepthValues, + std::vector& tvDepthValues, + double& rkbDiff ) const +{ + RimEclipseCase* eclipseCase = fractureModel->eclipseCase(); + if ( !eclipseCase ) + { + return false; + } + + RigEclipseWellLogExtractor eclExtractor( eclipseCase->eclipseCaseData(), + fractureModel->thicknessDirectionWellPath()->wellPathGeometry(), + "fracture model" ); + + rkbDiff = eclExtractor.wellPathData()->rkbDiff(); + + // Extract formation data + cvf::ref formationResultAccessor = + RigResultAccessorFactory::createFromResultAddress( eclipseCase->eclipseCaseData(), + 0, + RiaDefines::PorosityModelType::MATRIX_MODEL, + 0, + RigEclipseResultAddress( RiaDefines::ResultCatType::FORMATION_NAMES, + RiaDefines::activeFormationNamesResultName() ) ); + if ( !formationResultAccessor.notNull() ) + { + RiaLogging::error( QString( "No formation result found." ) ); + return false; + } + + CurveSamplingPointData curveData = + RimWellLogTrack::curveSamplingPointData( &eclExtractor, formationResultAccessor.p() ); + + std::vector formationNamesVector = RimWellLogTrack::formationNamesVector( eclipseCase ); + + double overburdenHeight = fractureModel->overburdenHeight(); + if ( overburdenHeight > 0.0 ) + { + RimWellLogTrack::addOverburden( formationNamesVector, curveData, overburdenHeight ); + } + + double underburdenHeight = fractureModel->underburdenHeight(); + if ( underburdenHeight > 0.0 ) + { + RimWellLogTrack::addUnderburden( formationNamesVector, curveData, underburdenHeight ); + } + + measuredDepthValues = curveData.md; + tvDepthValues = curveData.tvd; + + // Extract facies data + std::vector faciesValues = + m_fractureModelCalculator->extractValues( RiaDefines::CurveProperty::FACIES, timeStep ); + if ( faciesValues.empty() ) + { + RiaLogging::error( QString( "Empty facies data found for layer curve." ) ); + return false; + } + + values.resize( faciesValues.size() ); + + int layerNo = 0; + double previousFormation = -1.0; + double previousFacies = -1.0; + for ( size_t i = 0; i < faciesValues.size(); i++ ) + { + if ( previousFormation != curveData.data[i] || previousFacies != faciesValues[i] ) + { + layerNo++; + } + + values[i] = layerNo; + previousFormation = curveData.data[i]; + previousFacies = faciesValues[i]; + } + + return true; +} diff --git a/ApplicationCode/ProjectDataModel/RimFractureModelLayerCalculator.h b/ApplicationCode/ProjectDataModel/RimFractureModelLayerCalculator.h new file mode 100644 index 0000000000..02e48bcc12 --- /dev/null +++ b/ApplicationCode/ProjectDataModel/RimFractureModelLayerCalculator.h @@ -0,0 +1,46 @@ +///////////////////////////////////////////////////////////////////////////////// +// +// Copyright (C) 2020- 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 +// for more details. +// +///////////////////////////////////////////////////////////////////////////////// +#pragma once + +#include "RimFractureModelPropertyCalculator.h" + +#include "RiaFractureModelDefines.h" + +#include + +class RimFractureModelCalculator; +class RimFractureModel; + +class RimFractureModelLayerCalculator : public RimFractureModelPropertyCalculator +{ +public: + RimFractureModelLayerCalculator( RimFractureModelCalculator* calculator ); + + bool calculate( RiaDefines::CurveProperty curveProperty, + const RimFractureModel* fractureModel, + int timeStep, + std::vector& values, + std::vector& measuredDepthValues, + std::vector& tvDepthValues, + double& rkbDiff ) const override; + + bool isMatching( RiaDefines::CurveProperty curveProperty ) const override; + +private: + RimFractureModelCalculator* m_fractureModelCalculator; +}; diff --git a/ApplicationCode/ProjectDataModel/RimFractureModelPlot.cpp b/ApplicationCode/ProjectDataModel/RimFractureModelPlot.cpp index d5c9d69cd9..05c56cf342 100644 --- a/ApplicationCode/ProjectDataModel/RimFractureModelPlot.cpp +++ b/ApplicationCode/ProjectDataModel/RimFractureModelPlot.cpp @@ -19,22 +19,16 @@ #include "RiaDefines.h" #include "RiaFractureModelDefines.h" -#include "RiaLogging.h" -#include "RicfCommandObject.h" - -#include "RimEclipseCase.h" #include "RimFractureModel.h" #include "RimFractureModelCurve.h" #include "RimFractureModelPropertyCurve.h" -#include "RimLayerCurve.h" #include "RimWellLogTrack.h" -#include "RigWellLogCurveData.h" - #include "cafPdmBase.h" #include "cafPdmFieldScriptingCapability.h" #include "cafPdmObject.h" +#include "cafPdmObjectScriptingCapability.h" #include "cafPdmUiGroup.h" CAF_PDM_SOURCE_INIT( RimFractureModelPlot, "FractureModelPlot" ); @@ -128,124 +122,6 @@ void RimFractureModelPlot::applyDataSource() this->updateConnectedEditors(); } -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -void RimFractureModelPlot::getCurvePropertyValues( RiaDefines::CurveProperty curveProperty, std::vector& values ) const -{ - RimWellLogExtractionCurve* curve = findCurveByProperty( curveProperty ); - if ( curve ) - { - values = curve->curveData()->xValues(); - } -} - -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -void RimFractureModelPlot::getPorosityValues( std::vector& values ) const -{ - getCurvePropertyValues( RiaDefines::CurveProperty::POROSITY, values ); -} - -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -void RimFractureModelPlot::getFaciesValues( std::vector& values ) const -{ - getCurvePropertyValues( RiaDefines::CurveProperty::FACIES, values ); -} - -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -void RimFractureModelPlot::calculateLayers( std::vector>& layerBoundaryDepths, - std::vector>& layerBoundaryIndexes ) const -{ - std::vector curves; - descendantsIncludingThisOfType( curves ); - - if ( curves.empty() ) - { - return; - } - - // Expect to have only one of these - RimLayerCurve* layerCurve = curves[0]; - - const RigWellLogCurveData* curveData = layerCurve->curveData(); - - // Find - std::vector depths = curveData->depths( RiaDefines::DepthTypeEnum::TRUE_VERTICAL_DEPTH ); - std::vector layerValues = curveData->xValues(); - - size_t startIndex = 0; - for ( size_t i = 0; i < depths.size(); i++ ) - { - if ( startIndex != i && ( layerValues[startIndex] != layerValues[i] || i == depths.size() - 1 ) ) - { - layerBoundaryDepths.push_back( std::make_pair( depths[startIndex], depths[i] ) ); - layerBoundaryIndexes.push_back( std::make_pair( startIndex, i ) ); - startIndex = i; - } - } -} - -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -double RimFractureModelPlot::findValueAtTopOfLayer( const std::vector& values, - const std::vector>& layerBoundaryIndexes, - size_t layerNo ) -{ - size_t index = layerBoundaryIndexes[layerNo].first; - return values.at( index ); -} - -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -double RimFractureModelPlot::findValueAtBottomOfLayer( const std::vector& values, - const std::vector>& layerBoundaryIndexes, - size_t layerNo ) -{ - size_t index = layerBoundaryIndexes[layerNo].second; - return values.at( index ); -} - -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -void RimFractureModelPlot::computeAverageByLayer( const std::vector>& layerBoundaryIndexes, - const std::vector& inputVector, - std::vector& result ) -{ - for ( auto boundaryIndex : layerBoundaryIndexes ) - { - double sum = 0.0; - int nValues = 0; - for ( size_t i = boundaryIndex.first; i < boundaryIndex.second; i++ ) - { - sum += inputVector[i]; - nValues++; - } - result.push_back( sum / nValues ); - } -} - -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -void RimFractureModelPlot::extractTopOfLayerValues( const std::vector>& layerBoundaryIndexes, - const std::vector& inputVector, - std::vector& result ) -{ - for ( size_t i = 0; i < layerBoundaryIndexes.size(); i++ ) - { - result.push_back( findValueAtTopOfLayer( inputVector, layerBoundaryIndexes, i ) ); - } -} - //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- @@ -264,418 +140,3 @@ RimWellLogExtractionCurve* RimFractureModelPlot::findCurveByProperty( RiaDefines return nullptr; } - -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -std::vector RimFractureModelPlot::calculateTrueVerticalDepth() const -{ - std::vector> layerBoundaryDepths; - std::vector> layerBoundaryIndexes; - - calculateLayers( layerBoundaryDepths, layerBoundaryIndexes ); - - std::vector tvdTopZone; - for ( auto p : layerBoundaryDepths ) - { - double depthInFeet = RiaEclipseUnitTools::meterToFeet( p.first ); - tvdTopZone.push_back( depthInFeet ); - } - - return tvdTopZone; -} - -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -std::vector RimFractureModelPlot::findCurveAndComputeLayeredAverage( RiaDefines::CurveProperty curveProperty ) const -{ - RimWellLogExtractionCurve* curve = findCurveByProperty( curveProperty ); - if ( !curve ) - { - QString curveName = caf::AppEnum::uiText( curveProperty ); - RiaLogging::error( QString( "No curve for '%1' found" ).arg( curveName ) ); - return std::vector(); - } - - std::vector> layerBoundaryDepths; - std::vector> layerBoundaryIndexes; - calculateLayers( layerBoundaryDepths, layerBoundaryIndexes ); - - const RigWellLogCurveData* curveData = curve->curveData(); - std::vector values = curveData->xValues(); - std::vector result; - computeAverageByLayer( layerBoundaryIndexes, values, result ); - return result; -} - -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -std::vector RimFractureModelPlot::findCurveAndComputeTopOfLayer( RiaDefines::CurveProperty curveProperty ) const -{ - RimWellLogExtractionCurve* curve = findCurveByProperty( curveProperty ); - if ( !curve ) - { - QString curveName = caf::AppEnum::uiText( curveProperty ); - RiaLogging::error( QString( "No curve for '%1' found" ).arg( curveName ) ); - return std::vector(); - } - - std::vector> layerBoundaryDepths; - std::vector> layerBoundaryIndexes; - calculateLayers( layerBoundaryDepths, layerBoundaryIndexes ); - - const RigWellLogCurveData* curveData = curve->curveData(); - std::vector values = curveData->xValues(); - std::vector result; - extractTopOfLayerValues( layerBoundaryIndexes, values, result ); - return result; -} - -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -std::vector RimFractureModelPlot::calculatePorosity() const -{ - return findCurveAndComputeLayeredAverage( RiaDefines::CurveProperty::POROSITY ); -} - -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -std::vector RimFractureModelPlot::calculateReservoirPressure() const -{ - std::vector pressureBar = findCurveAndComputeTopOfLayer( RiaDefines::CurveProperty::PRESSURE ); - - std::vector pressurePsi; - for ( double p : pressureBar ) - { - pressurePsi.push_back( RiaEclipseUnitTools::barToPsi( p ) ); - } - - return pressurePsi; -} - -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -std::vector RimFractureModelPlot::calculateHorizontalPermeability() const -{ - return findCurveAndComputeLayeredAverage( RiaDefines::CurveProperty::PERMEABILITY_X ); -} - -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -std::vector RimFractureModelPlot::calculateVerticalPermeability() const -{ - return findCurveAndComputeLayeredAverage( RiaDefines::CurveProperty::PERMEABILITY_Z ); -} - -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -std::vector RimFractureModelPlot::calculateStress() const -{ - std::vector stress; - std::vector stressGradients; - std::vector initialStress; - calculateStressWithGradients( stress, stressGradients, initialStress ); - return stress; -} - -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -std::vector RimFractureModelPlot::calculateInitialStress() const -{ - std::vector stress; - std::vector stressGradients; - std::vector initialStress; - calculateStressWithGradients( stress, stressGradients, initialStress ); - return initialStress; -} - -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -std::vector RimFractureModelPlot::findCurveXValuesByProperty( RiaDefines::CurveProperty curveProperty ) const -{ - RimWellLogExtractionCurve* curve = findCurveByProperty( curveProperty ); - if ( !curve ) - { - QString curveName = caf::AppEnum::uiText( curveProperty ); - RiaLogging::error( QString( "%1 data not found." ).arg( curveName ) ); - return std::vector(); - } - return curve->curveData()->xValues(); -} - -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -bool RimFractureModelPlot::calculateStressWithGradients( std::vector& stress, - std::vector& stressGradients, - std::vector& initialStress ) const -{ - // Reference stress - const double verticalStressRef = m_fractureModel->verticalStress(); - const double verticalStressGradientRef = m_fractureModel->verticalStressGradient(); - const double stressDepthRef = m_fractureModel->stressDepth(); - - std::vector> layerBoundaryDepths; - std::vector> layerBoundaryIndexes; - calculateLayers( layerBoundaryDepths, layerBoundaryIndexes ); - - // Biot coefficient - std::vector biotData = findCurveXValuesByProperty( RiaDefines::CurveProperty::BIOT_COEFFICIENT ); - - // K0 - std::vector k0Data = findCurveXValuesByProperty( RiaDefines::CurveProperty::K0 ); - - // Pressure at the give time step - std::vector timeStepPressureData = findCurveXValuesByProperty( RiaDefines::CurveProperty::PRESSURE ); - - // Initial pressure - std::vector initialPressureData = findCurveXValuesByProperty( RiaDefines::CurveProperty::INITIAL_PRESSURE ); - - // Poissons ratio - std::vector poissonsRatioData = findCurveXValuesByProperty( RiaDefines::CurveProperty::POISSONS_RATIO ); - - // Check that we have data from all curves - if ( biotData.empty() || k0Data.empty() || timeStepPressureData.empty() || initialPressureData.empty() || - poissonsRatioData.empty() ) - { - return false; - } - - std::vector stressForGradients; - std::vector pressureForGradients; - std::vector depthForGradients; - - // Calculate the stress - for ( size_t i = 0; i < layerBoundaryDepths.size(); i++ ) - { - double depthTopOfZone = layerBoundaryDepths[i].first; - double depthBottomOfZone = layerBoundaryDepths[i].second; - - // Data from curves at the top zone depth - double k0 = findValueAtTopOfLayer( k0Data, layerBoundaryIndexes, i ); - double biot = findValueAtTopOfLayer( biotData, layerBoundaryIndexes, i ); - double poissonsRatio = findValueAtTopOfLayer( poissonsRatioData, layerBoundaryIndexes, i ); - double initialPressure = findValueAtTopOfLayer( initialPressureData, layerBoundaryIndexes, i ); - double timeStepPressure = findValueAtTopOfLayer( timeStepPressureData, layerBoundaryIndexes, i ); - - // Vertical stress - // Use difference between reference depth and depth of top of zone - double depthDiff = depthTopOfZone - stressDepthRef; - double Sv = verticalStressRef + verticalStressGradientRef * depthDiff; - - double Sh_init = k0 * Sv + initialPressure * ( 1.0 - k0 ); - double pressureDiff = timeStepPressure - initialPressure; - - // Vertical stress diff assumed to be zero - double Sv_diff = 0.0; - double deltaHorizontalStress = poissonsRatio / ( 1.0 - poissonsRatio ) * ( Sv_diff - biot * pressureDiff ) + - ( biot * pressureDiff ); - - double depletionStress = Sh_init + deltaHorizontalStress; - stress.push_back( RiaEclipseUnitTools::barToPsi( depletionStress ) ); - - initialStress.push_back( RiaEclipseUnitTools::barToPsi( Sh_init ) ); - - // Cache some results for the gradients calculation - stressForGradients.push_back( Sv ); - pressureForGradients.push_back( initialPressure ); - depthForGradients.push_back( depthTopOfZone ); - - if ( i == layerBoundaryDepths.size() - 1 ) - { - // Use the bottom of the last layer to compute gradient for last layer - double bottomInitialPressure = findValueAtBottomOfLayer( initialPressureData, layerBoundaryIndexes, i ); - - double bottomDepthDiff = depthBottomOfZone - stressDepthRef; - double bottomSv = verticalStressRef + verticalStressGradientRef * bottomDepthDiff; - stressForGradients.push_back( bottomSv ); - pressureForGradients.push_back( bottomInitialPressure ); - depthForGradients.push_back( depthBottomOfZone ); - } - } - - assert( stressForGradients.size() == layerBoundaryDepths.size() + 1 ); - assert( pressureForGradients.size() == layerBoundaryDepths.size() + 1 ); - assert( depthForGradients.size() == layerBoundaryDepths.size() + 1 ); - - // Second pass to calculate the stress gradients - for ( size_t i = 0; i < layerBoundaryDepths.size(); i++ ) - { - double diffStress = stressForGradients[i + 1] - stressForGradients[i]; - double diffPressure = pressureForGradients[i + 1] - pressureForGradients[i]; - double diffDepth = depthForGradients[i + 1] - depthForGradients[i]; - double k0 = findValueAtTopOfLayer( k0Data, layerBoundaryIndexes, i ); - double stressGradient = ( diffStress * k0 + diffPressure * ( 1.0 - k0 ) ) / diffDepth; - stressGradients.push_back( RiaEclipseUnitTools::barPerMeterToPsiPerFeet( stressGradient ) ); - } - - return true; -} - -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -std::vector RimFractureModelPlot::calculateStressGradient() const -{ - std::vector stress; - std::vector stressGradients; - std::vector initialStress; - calculateStressWithGradients( stress, stressGradients, initialStress ); - return stressGradients; -} - -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -void RimFractureModelPlot::calculateTemperature( std::vector& temperatures ) const -{ - // Reference temperature. Unit: degrees celsius - const double referenceTemperature = m_fractureModel->referenceTemperature(); - - // Reference temperature gradient. Unit: degrees Celsius per meter - const double referenceTemperatureGradient = m_fractureModel->referenceTemperatureGradient(); - - // Reference depth for temperature. Unit: meter. - const double referenceTemperatureDepth = m_fractureModel->referenceTemperatureDepth(); - - std::vector> layerBoundaryDepths; - std::vector> layerBoundaryIndexes; - calculateLayers( layerBoundaryDepths, layerBoundaryIndexes ); - - // Calculate the temperatures - for ( size_t i = 0; i < layerBoundaryDepths.size(); i++ ) - { - double depthTopOfZone = layerBoundaryDepths[i].first; - - // Use difference between reference depth and depth of top of zone - double depthDiff = depthTopOfZone - referenceTemperatureDepth; - double temperature = referenceTemperature + referenceTemperatureGradient * depthDiff; - - temperatures.push_back( temperature ); - } -} - -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -std::vector RimFractureModelPlot::calculateYoungsModulus() const -{ - std::vector valuesGPa = findCurveAndComputeLayeredAverage( RiaDefines::CurveProperty::YOUNGS_MODULUS ); - std::vector valuesMMpsi; - for ( auto value : valuesGPa ) - { - valuesMMpsi.push_back( value * 0.14503773773 ); - } - - return valuesMMpsi; -} - -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -std::vector RimFractureModelPlot::calculatePoissonsRatio() const -{ - return findCurveAndComputeLayeredAverage( RiaDefines::CurveProperty::POISSONS_RATIO ); -} - -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -std::vector RimFractureModelPlot::calculateKIc() const -{ - return findCurveAndComputeLayeredAverage( RiaDefines::CurveProperty::K_IC ); -} - -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -std::vector RimFractureModelPlot::calculateFluidLossCoefficient() const -{ - return findCurveAndComputeLayeredAverage( RiaDefines::CurveProperty::FLUID_LOSS_COEFFICIENT ); -} - -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -std::vector RimFractureModelPlot::calculateSpurtLoss() const -{ - return findCurveAndComputeLayeredAverage( RiaDefines::CurveProperty::SPURT_LOSS ); -} - -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -std::vector RimFractureModelPlot::calculateProppandEmbedment() const -{ - return findCurveAndComputeLayeredAverage( RiaDefines::CurveProperty::PROPPANT_EMBEDMENT ); -} - -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -std::vector RimFractureModelPlot::calculateImmobileFluidSaturation() const -{ - return findCurveAndComputeLayeredAverage( RiaDefines::CurveProperty::IMMOBILE_FLUID_SATURATION ); -} - -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -std::vector RimFractureModelPlot::calculateTemperature() const -{ - std::vector temperaturesCelsius; - calculateTemperature( temperaturesCelsius ); - - // Convert to Fahrenheit - std::vector temperaturesFahrenheit; - for ( double t : temperaturesCelsius ) - { - temperaturesFahrenheit.push_back( t * 1.8 + 32.0 ); - } - - return temperaturesFahrenheit; -} - -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -std::vector RimFractureModelPlot::calculateRelativePermeabilityFactor() const -{ - return findCurveAndComputeLayeredAverage( RiaDefines::CurveProperty::RELATIVE_PERMEABILITY_FACTOR ); -} - -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -std::vector RimFractureModelPlot::calculatePoroElasticConstant() const -{ - return findCurveAndComputeLayeredAverage( RiaDefines::CurveProperty::PORO_ELASTIC_CONSTANT ); -} - -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -std::vector RimFractureModelPlot::calculateThermalExpansionCoefficient() const -{ - // SI unit is 1/Celsius - std::vector coefficientCelsius = - findCurveAndComputeLayeredAverage( RiaDefines::CurveProperty::THERMAL_EXPANSION_COEFFICIENT ); - - // Field unit is 1/Fahrenheit - std::vector coefficientFahrenheit; - for ( double c : coefficientCelsius ) - { - coefficientFahrenheit.push_back( c / 1.8 ); - } - - return coefficientFahrenheit; -} diff --git a/ApplicationCode/ProjectDataModel/RimFractureModelPlot.h b/ApplicationCode/ProjectDataModel/RimFractureModelPlot.h index 69bf4d6069..1d74720298 100644 --- a/ApplicationCode/ProjectDataModel/RimFractureModelPlot.h +++ b/ApplicationCode/ProjectDataModel/RimFractureModelPlot.h @@ -21,13 +21,10 @@ #include "RiaFractureModelDefines.h" -#include "cafPdmField.h" #include "cafPdmPtrField.h" -#include - -class RimWellLogExtractionCurve; class RimFractureModel; +class RimWellLogExtractionCurve; class RimFractureModelPlot : public RimDepthTrackPlot { @@ -39,57 +36,8 @@ public: void setFractureModel( RimFractureModel* fractureModel ); RimFractureModel* fractureModel(); - void getPorosityValues( std::vector& values ) const; - void getFaciesValues( std::vector& values ) const; - void getCurvePropertyValues( RiaDefines::CurveProperty curveProperty, std::vector& values ) const; - - std::vector calculateTrueVerticalDepth() const; - std::vector calculatePorosity() const; - std::vector calculateVerticalPermeability() const; - std::vector calculateHorizontalPermeability() const; - std::vector calculateReservoirPressure() const; - std::vector calculateStress() const; - std::vector calculateInitialStress() const; - std::vector calculateStressGradient() const; - std::vector calculateYoungsModulus() const; - std::vector calculatePoissonsRatio() const; - std::vector calculateKIc() const; - std::vector calculateFluidLossCoefficient() const; - std::vector calculateSpurtLoss() const; - std::vector calculateProppandEmbedment() const; - - std::vector calculateImmobileFluidSaturation() const; - std::vector calculateTemperature() const; - std::vector calculateRelativePermeabilityFactor() const; - std::vector calculatePoroElasticConstant() const; - std::vector calculateThermalExpansionCoefficient() const; - - void calculateTemperature( std::vector& temperatures ) const; - protected: - std::vector findCurveAndComputeLayeredAverage( RiaDefines::CurveProperty curveProperty ) const; - std::vector findCurveXValuesByProperty( RiaDefines::CurveProperty curveProperty ) const; - std::vector findCurveAndComputeTopOfLayer( RiaDefines::CurveProperty curveProperty ) const; - - void calculateLayers( std::vector>& layerBoundaryDepths, - std::vector>& layerBoundaryIndexes ) const; RimWellLogExtractionCurve* findCurveByProperty( RiaDefines::CurveProperty curveProperty ) const; - bool calculateStressWithGradients( std::vector& stress, - std::vector& stressGradients, - std::vector& initialStress ) const; - - static double findValueAtTopOfLayer( const std::vector& values, - const std::vector>& layerBoundaryIndexes, - size_t layerNo ); - static double findValueAtBottomOfLayer( const std::vector& values, - const std::vector>& layerBoundaryIndexes, - size_t layerNo ); - static void computeAverageByLayer( const std::vector>& layerBoundaryIndexes, - const std::vector& inputVector, - std::vector& result ); - static void extractTopOfLayerValues( const std::vector>& layerBoundaryIndexes, - const std::vector& inputVector, - std::vector& result ); void defineUiOrdering( QString uiConfigName, caf::PdmUiOrdering& uiOrdering ) override; diff --git a/ApplicationCode/ProjectDataModel/RimFractureModelPropertyCalculator.h b/ApplicationCode/ProjectDataModel/RimFractureModelPropertyCalculator.h new file mode 100644 index 0000000000..b8ff34bcad --- /dev/null +++ b/ApplicationCode/ProjectDataModel/RimFractureModelPropertyCalculator.h @@ -0,0 +1,44 @@ +///////////////////////////////////////////////////////////////////////////////// +// +// Copyright (C) 2020- 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 +// for more details. +// +///////////////////////////////////////////////////////////////////////////////// + +#pragma once + +#include "RiaFractureModelDefines.h" + +#include + +class RimFractureModelCalculator; +class RimFractureModel; + +//================================================================================================== +/// +//================================================================================================== +class RimFractureModelPropertyCalculator +{ +public: + virtual ~RimFractureModelPropertyCalculator(){}; + + virtual bool isMatching( RiaDefines::CurveProperty curveProperty ) const = 0; + virtual bool calculate( RiaDefines::CurveProperty curveProperty, + const RimFractureModel* fractureModel, + int timeStep, + std::vector& values, + std::vector& measuredDepthValues, + std::vector& tvDepthValues, + double& rkbDiff ) const = 0; +}; diff --git a/ApplicationCode/ProjectDataModel/RimFractureModelStressCalculator.cpp b/ApplicationCode/ProjectDataModel/RimFractureModelStressCalculator.cpp new file mode 100644 index 0000000000..5a14fa2d7b --- /dev/null +++ b/ApplicationCode/ProjectDataModel/RimFractureModelStressCalculator.cpp @@ -0,0 +1,150 @@ +///////////////////////////////////////////////////////////////////////////////// +// +// Copyright (C) 2020- 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 +// for more details. +// +///////////////////////////////////////////////////////////////////////////////// +#include "RimFractureModelStressCalculator.h" + +#include "RiaDefines.h" +#include "RiaFractureModelDefines.h" + +#include "RigEclipseCaseData.h" +#include "RigEclipseWellLogExtractor.h" +#include "RigWellPath.h" + +#include "RigWellPathGeometryTools.h" +#include "RimCase.h" +#include "RimEclipseCase.h" +#include "RimFractureModel.h" +#include "RimFractureModelCalculator.h" +#include "RimFractureModelStressCalculator.h" +#include "RimModeledWellPath.h" + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +RimFractureModelStressCalculator::RimFractureModelStressCalculator( RimFractureModelCalculator* fractureModelCalculator ) + : m_fractureModelCalculator( fractureModelCalculator ) +{ +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +bool RimFractureModelStressCalculator::isMatching( RiaDefines::CurveProperty curveProperty ) const +{ + std::vector matching = {RiaDefines::CurveProperty::INITIAL_STRESS, + RiaDefines::CurveProperty::STRESS, + RiaDefines::CurveProperty::STRESS_GRADIENT, + RiaDefines::CurveProperty::TEMPERATURE}; + + return std::find( matching.begin(), matching.end(), curveProperty ) != matching.end(); +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +bool RimFractureModelStressCalculator::calculate( RiaDefines::CurveProperty curveProperty, + const RimFractureModel* fractureModel, + int timeStep, + std::vector& values, + std::vector& measuredDepthValues, + std::vector& tvDepthValues, + double& rkbDiff ) const +{ + RimEclipseCase* eclipseCase = fractureModel->eclipseCase(); + if ( !eclipseCase ) + { + return false; + } + + std::vector tvDepthInFeet = m_fractureModelCalculator->calculateTrueVerticalDepth(); + for ( double f : tvDepthInFeet ) + { + tvDepthValues.push_back( RiaEclipseUnitTools::feetToMeter( f ) ); + } + + if ( curveProperty == RiaDefines::CurveProperty::STRESS ) + { + values = m_fractureModelCalculator->calculateStress(); + std::vector stressGradients = m_fractureModelCalculator->calculateStressGradient(); + addDatapointsForBottomOfLayers( tvDepthValues, values, stressGradients ); + } + else if ( curveProperty == RiaDefines::CurveProperty::INITIAL_STRESS ) + { + values = m_fractureModelCalculator->calculateInitialStress(); + std::vector stressGradients = m_fractureModelCalculator->calculateStressGradient(); + addDatapointsForBottomOfLayers( tvDepthValues, values, stressGradients ); + } + else if ( curveProperty == RiaDefines::CurveProperty::STRESS_GRADIENT ) + { + values = m_fractureModelCalculator->calculateStressGradient(); + } + else if ( curveProperty == RiaDefines::CurveProperty::TEMPERATURE ) + { + m_fractureModelCalculator->calculateTemperature( values ); + } + + if ( eclipseCase ) + { + RigWellPath* wellPathGeometry = fractureModel->thicknessDirectionWellPath()->wellPathGeometry(); + RigEclipseWellLogExtractor eclExtractor( eclipseCase->eclipseCaseData(), wellPathGeometry, "fracture model" ); + + rkbDiff = eclExtractor.wellPathData()->rkbDiff(); + + // Generate MD data by interpolation + const std::vector& mdValuesOfWellPath = wellPathGeometry->measureDepths(); + std::vector tvdValuesOfWellPath = wellPathGeometry->trueVerticalDepths(); + + measuredDepthValues = + RigWellPathGeometryTools::interpolateMdFromTvd( mdValuesOfWellPath, tvdValuesOfWellPath, tvDepthValues ); + CVF_ASSERT( measuredDepthValues.size() == tvDepthValues.size() ); + } + + return true; +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +void RimFractureModelStressCalculator::addDatapointsForBottomOfLayers( std::vector& tvDepthValues, + std::vector& stress, + const std::vector& stressGradients ) +{ + std::vector tvdWithBottomLayers; + std::vector valuesWithBottomLayers; + for ( size_t i = 0; i < stress.size(); i++ ) + { + // Add the data point at top of the layer + double topLayerDepth = tvDepthValues[i]; + double stressValue = stress[i]; + tvdWithBottomLayers.push_back( topLayerDepth ); + valuesWithBottomLayers.push_back( stressValue ); + + // Add extra data points for bottom part of the layer + if ( i < stress.size() - 1 ) + { + double bottomLayerDepth = tvDepthValues[i + 1]; + double diffDepthFeet = RiaEclipseUnitTools::meterToFeet( bottomLayerDepth - topLayerDepth ); + double bottomStress = stressValue + diffDepthFeet * stressGradients[i]; + + tvdWithBottomLayers.push_back( bottomLayerDepth ); + valuesWithBottomLayers.push_back( bottomStress ); + } + } + + stress = valuesWithBottomLayers; + tvDepthValues = tvdWithBottomLayers; +} diff --git a/ApplicationCode/ProjectDataModel/RimFractureModelStressCalculator.h b/ApplicationCode/ProjectDataModel/RimFractureModelStressCalculator.h new file mode 100644 index 0000000000..6692a7850d --- /dev/null +++ b/ApplicationCode/ProjectDataModel/RimFractureModelStressCalculator.h @@ -0,0 +1,52 @@ +///////////////////////////////////////////////////////////////////////////////// +// +// Copyright (C) 2020- 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 +// for more details. +// +///////////////////////////////////////////////////////////////////////////////// +#pragma once + +#include "RimFractureModelPropertyCalculator.h" + +#include "RiaFractureModelDefines.h" + +#include + +class RimFractureModelCalculator; +class RimFractureModel; + +class QString; + +class RimFractureModelStressCalculator : public RimFractureModelPropertyCalculator +{ +public: + RimFractureModelStressCalculator( RimFractureModelCalculator* calculator ); + + bool calculate( RiaDefines::CurveProperty curveProperty, + const RimFractureModel* fractureModel, + int timeStep, + std::vector& values, + std::vector& measuredDepthValues, + std::vector& tvDepthValues, + double& rkbDiff ) const override; + + bool isMatching( RiaDefines::CurveProperty curveProperty ) const override; + +private: + static void addDatapointsForBottomOfLayers( std::vector& tvDepthValues, + std::vector& stress, + const std::vector& stressGradients ); + + RimFractureModelCalculator* m_fractureModelCalculator; +}; diff --git a/ApplicationCode/ProjectDataModel/RimFractureModelStressCurve.cpp b/ApplicationCode/ProjectDataModel/RimFractureModelStressCurve.cpp index 673347cdd2..da977c5f90 100644 --- a/ApplicationCode/ProjectDataModel/RimFractureModelStressCurve.cpp +++ b/ApplicationCode/ProjectDataModel/RimFractureModelStressCurve.cpp @@ -30,6 +30,7 @@ #include "RimCase.h" #include "RimEclipseCase.h" #include "RimFractureModel.h" +#include "RimFractureModelCalculator.h" #include "RimFractureModelPlot.h" #include "RimModeledWellPath.h" #include "RimWellLogFile.h" @@ -107,59 +108,14 @@ void RimFractureModelStressCurve::performDataExtraction( bool* isUsingPseudoLeng *isUsingPseudoLength = false; - // Extract porosity data: get the porosity values from parent - RimFractureModelPlot* fractureModelPlot; - firstAncestorOrThisOfType( fractureModelPlot ); - if ( !fractureModelPlot || !m_fractureModel ) + bool isOk = m_fractureModel->calculator() + ->extractCurveData( curveProperty(), m_timeStep, values, measuredDepthValues, tvDepthValues, rkbDiff ); + if ( !isOk ) { - RiaLogging::error( QString( "No fracture model plot found." ) ); return; } - std::vector tvDepthInFeet = fractureModelPlot->calculateTrueVerticalDepth(); - for ( double f : tvDepthInFeet ) - { - tvDepthValues.push_back( RiaEclipseUnitTools::feetToMeter( f ) ); - } - - if ( m_curveProperty() == RiaDefines::CurveProperty::STRESS ) - { - values = fractureModelPlot->calculateStress(); - std::vector stressGradients = fractureModelPlot->calculateStressGradient(); - addDatapointsForBottomOfLayers( tvDepthValues, values, stressGradients ); - } - else if ( m_curveProperty() == RiaDefines::CurveProperty::INITIAL_STRESS ) - { - values = fractureModelPlot->calculateInitialStress(); - std::vector stressGradients = fractureModelPlot->calculateStressGradient(); - addDatapointsForBottomOfLayers( tvDepthValues, values, stressGradients ); - } - else if ( m_curveProperty() == RiaDefines::CurveProperty::STRESS_GRADIENT ) - { - values = fractureModelPlot->calculateStressGradient(); - } - else if ( m_curveProperty() == RiaDefines::CurveProperty::TEMPERATURE ) - { - fractureModelPlot->calculateTemperature( values ); - } - - RimEclipseCase* eclipseCase = dynamic_cast( m_case.value() ); - if ( eclipseCase ) - { - RigWellPath* wellPathGeometry = m_fractureModel->thicknessDirectionWellPath()->wellPathGeometry(); - RigEclipseWellLogExtractor eclExtractor( eclipseCase->eclipseCaseData(), wellPathGeometry, "fracture model" ); - - rkbDiff = eclExtractor.wellPathData()->rkbDiff(); - - // Generate MD data by interpolation - const std::vector& mdValuesOfWellPath = wellPathGeometry->measureDepths(); - std::vector tvdValuesOfWellPath = wellPathGeometry->trueVerticalDepths(); - - measuredDepthValues = - RigWellPathGeometryTools::interpolateMdFromTvd( mdValuesOfWellPath, tvdValuesOfWellPath, tvDepthValues ); - CVF_ASSERT( measuredDepthValues.size() == tvDepthValues.size() ); - } - + RimEclipseCase* eclipseCase = dynamic_cast( m_case.value() ); RiaEclipseUnitTools::UnitSystem eclipseUnitsType = eclipseCase->eclipseCaseData()->unitsType(); if ( eclipseUnitsType == RiaEclipseUnitTools::UnitSystem::UNITS_FIELD ) { @@ -201,36 +157,3 @@ QString RimFractureModelStressCurve::createCurveAutoName() { return caf::AppEnum::uiText( m_curveProperty() ); } - -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -void RimFractureModelStressCurve::addDatapointsForBottomOfLayers( std::vector& tvDepthValues, - std::vector& stress, - const std::vector& stressGradients ) -{ - std::vector tvdWithBottomLayers; - std::vector valuesWithBottomLayers; - for ( size_t i = 0; i < stress.size(); i++ ) - { - // Add the data point at top of the layer - double topLayerDepth = tvDepthValues[i]; - double stressValue = stress[i]; - tvdWithBottomLayers.push_back( topLayerDepth ); - valuesWithBottomLayers.push_back( stressValue ); - - // Add extra data points for bottom part of the layer - if ( i < stress.size() - 1 ) - { - double bottomLayerDepth = tvDepthValues[i + 1]; - double diffDepthFeet = RiaEclipseUnitTools::meterToFeet( bottomLayerDepth - topLayerDepth ); - double bottomStress = stressValue + diffDepthFeet * stressGradients[i]; - - tvdWithBottomLayers.push_back( bottomLayerDepth ); - valuesWithBottomLayers.push_back( bottomStress ); - } - } - - stress = valuesWithBottomLayers; - tvDepthValues = tvdWithBottomLayers; -} diff --git a/ApplicationCode/ProjectDataModel/RimFractureModelWellLogCalculator.cpp b/ApplicationCode/ProjectDataModel/RimFractureModelWellLogCalculator.cpp new file mode 100644 index 0000000000..b1ef870904 --- /dev/null +++ b/ApplicationCode/ProjectDataModel/RimFractureModelWellLogCalculator.cpp @@ -0,0 +1,401 @@ +///////////////////////////////////////////////////////////////////////////////// +// +// Copyright (C) 2020- 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 +// for more details. +// +///////////////////////////////////////////////////////////////////////////////// +#include "RimFractureModelWellLogCalculator.h" + +#include "RiaDefines.h" +#include "RiaFractureModelDefines.h" +#include "RiaInterpolationTools.h" +#include "RiaLogging.h" + +#include "RigEclipseCaseData.h" +#include "RigEclipseWellLogExtractor.h" +#include "RigResultAccessor.h" +#include "RigResultAccessorFactory.h" +#include "RigWellLogCurveData.h" +#include "RigWellPath.h" + +#include "RimCase.h" +#include "RimEclipseCase.h" +#include "RimEclipseInputProperty.h" +#include "RimEclipseInputPropertyCollection.h" +#include "RimEclipseResultDefinition.h" +#include "RimFractureModel.h" +#include "RimFractureModelCalculator.h" +#include "RimFractureModelPlot.h" +#include "RimLayerCurve.h" +#include "RimModeledWellPath.h" +#include "RimTools.h" +#include "RimWellLogFile.h" +#include "RimWellLogPlot.h" +#include "RimWellLogTrack.h" +#include "RimWellPath.h" +#include "RimWellPathCollection.h" +#include "RimWellPlotTools.h" + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +RimFractureModelWellLogCalculator::RimFractureModelWellLogCalculator( RimFractureModelCalculator& ) +{ +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +bool RimFractureModelWellLogCalculator::isMatching( RiaDefines::CurveProperty curveProperty ) const +{ + std::vector matching = { + RiaDefines::CurveProperty::FACIES, + RiaDefines::CurveProperty::POROSITY, + RiaDefines::CurveProperty::PERMEABILITY_X, + RiaDefines::CurveProperty::PERMEABILITY_Z, + RiaDefines::CurveProperty::INITIAL_PRESSURE, + RiaDefines::CurveProperty::PRESSURE, + }; + + return std::find( matching.begin(), matching.end(), curveProperty ) != matching.end(); +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +bool RimFractureModelWellLogCalculator::calculate( RiaDefines::CurveProperty curveProperty, + const RimFractureModel* fractureModel, + int timeStep, + std::vector& values, + std::vector& measuredDepthValues, + std::vector& tvDepthValues, + double& rkbDiff ) const +{ + RimEclipseCase* eclipseCase = fractureModel->eclipseCase(); + if ( !eclipseCase ) + { + return false; + } + + // TODO: improve this.. + if ( curveProperty == RiaDefines::CurveProperty::INITIAL_PRESSURE ) + { + timeStep = 0; + } + + RigEclipseWellLogExtractor eclExtractor( eclipseCase->eclipseCaseData(), + fractureModel->thicknessDirectionWellPath()->wellPathGeometry(), + "fracture model" ); + + measuredDepthValues = eclExtractor.cellIntersectionMDs(); + tvDepthValues = eclExtractor.cellIntersectionTVDs(); + rkbDiff = eclExtractor.wellPathData()->rkbDiff(); + + RimEclipseResultDefinition eclipseResultDefinition; + eclipseResultDefinition.setEclipseCase( eclipseCase ); + eclipseResultDefinition.setResultType( fractureModel->eclipseResultCategory( curveProperty ) ); + eclipseResultDefinition.setPorosityModel( RiaDefines::PorosityModelType::MATRIX_MODEL ); + eclipseResultDefinition.setResultVariable( fractureModel->eclipseResultVariable( curveProperty ) ); + + eclipseResultDefinition.loadResult(); + + cvf::ref resAcc = + RigResultAccessorFactory::createFromResultDefinition( eclipseCase->eclipseCaseData(), + 0, + timeStep, + &eclipseResultDefinition ); + + if ( resAcc.notNull() ) + { + eclExtractor.curveData( resAcc.p(), &values ); + } + else + { + RiaLogging::error( QString( "No result found for %1" ).arg( eclipseResultDefinition.resultVariable() ) ); + return false; + } + + double overburdenHeight = fractureModel->overburdenHeight(); + if ( overburdenHeight > 0.0 ) + { + addOverburden( curveProperty, fractureModel, tvDepthValues, measuredDepthValues, values ); + } + + double underburdenHeight = fractureModel->underburdenHeight(); + if ( underburdenHeight > 0.0 ) + { + addUnderburden( curveProperty, fractureModel, tvDepthValues, measuredDepthValues, values ); + } + + if ( hasMissingValues( values ) ) + { + if ( fractureModel->missingValueStrategy( curveProperty ) == RimFractureModel::MissingValueStrategy::DEFAULT_VALUE ) + { + // Try to locate a backup accessor (e.g. PORO_1 for PORO) + cvf::ref backupResAcc = findMissingValuesAccessor( eclipseCase->eclipseCaseData(), + eclipseCase->inputPropertyCollection(), + 0, + timeStep, + &eclipseResultDefinition ); + + if ( backupResAcc.notNull() ) + { + RiaLogging::info( QString( "Reading missing values from input properties for %1." ) + .arg( eclipseResultDefinition.resultVariable() ) ); + std::vector replacementValues; + eclExtractor.curveData( backupResAcc.p(), &replacementValues ); + + double overburdenHeight = fractureModel->overburdenHeight(); + if ( overburdenHeight > 0.0 ) + { + double defaultOverburdenValue = std::numeric_limits::infinity(); + if ( fractureModel->burdenStrategy( curveProperty ) == RimFractureModel::BurdenStrategy::DEFAULT_VALUE ) + { + defaultOverburdenValue = fractureModel->getDefaultForMissingOverburdenValue( curveProperty ); + } + + replacementValues.insert( replacementValues.begin(), defaultOverburdenValue ); + replacementValues.insert( replacementValues.begin(), defaultOverburdenValue ); + } + + double underburdenHeight = fractureModel->underburdenHeight(); + if ( underburdenHeight > 0.0 ) + { + double defaultUnderburdenValue = std::numeric_limits::infinity(); + if ( fractureModel->burdenStrategy( curveProperty ) == RimFractureModel::BurdenStrategy::DEFAULT_VALUE ) + { + defaultUnderburdenValue = fractureModel->getDefaultForMissingUnderburdenValue( curveProperty ); + } + + replacementValues.push_back( defaultUnderburdenValue ); + replacementValues.push_back( defaultUnderburdenValue ); + } + + replaceMissingValues( values, replacementValues ); + } + + // If the backup accessor is not found, or does not provide all the missing values: + // use default value from the fracture model + if ( !backupResAcc.notNull() || hasMissingValues( values ) ) + { + RiaLogging::info( QString( "Using default value for %1" ).arg( eclipseResultDefinition.resultVariable() ) ); + + double defaultValue = fractureModel->getDefaultForMissingValue( curveProperty ); + + replaceMissingValues( values, defaultValue ); + } + } + else if ( fractureModel->missingValueStrategy( curveProperty ) == + RimFractureModel::MissingValueStrategy::LINEAR_INTERPOLATION ) + { + RiaLogging::info( + QString( "Interpolating missing values for %1" ).arg( eclipseResultDefinition.resultVariable() ) ); + RiaInterpolationTools::interpolateMissingValues( measuredDepthValues, values ); + } + else + { + // Get the missing data from other curve + RiaDefines::CurveProperty replacementProperty = + fractureModel->getDefaultPropertyForMissingValues( curveProperty ); + + std::vector initialValues; + std::vector initialMds; + std::vector initialTvds; + double initialRkbDiff = -1.0; + calculate( replacementProperty, fractureModel, timeStep, initialValues, initialMds, initialTvds, initialRkbDiff ); + + if ( initialValues.empty() ) + { + RiaLogging::error( QString( "Empty replacement data found for fracture model curve." ) ); + return false; + } + + CVF_ASSERT( values.size() == initialValues.size() ); + replaceMissingValues( values, initialValues ); + } + } + + return true; +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +bool RimFractureModelWellLogCalculator::hasMissingValues( const std::vector& values ) +{ + for ( double v : values ) + { + if ( v == std::numeric_limits::infinity() ) + { + return true; + } + } + + return false; +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +void RimFractureModelWellLogCalculator::replaceMissingValues( std::vector& values, double defaultValue ) +{ + for ( double& v : values ) + { + if ( v == std::numeric_limits::infinity() ) + { + v = defaultValue; + } + } +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +void RimFractureModelWellLogCalculator::replaceMissingValues( std::vector& values, + const std::vector& replacementValues ) +{ + CVF_ASSERT( values.size() == replacementValues.size() ); + for ( size_t i = 0; i < values.size(); i++ ) + { + if ( values[i] == std::numeric_limits::infinity() ) + { + values[i] = replacementValues[i]; + } + } +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +cvf::ref + RimFractureModelWellLogCalculator::findMissingValuesAccessor( RigEclipseCaseData* caseData, + RimEclipseInputPropertyCollection* inputPropertyCollection, + int gridIndex, + int timeStepIndex, + RimEclipseResultDefinition* eclipseResultDefinition ) const +{ + QString resultName = eclipseResultDefinition->resultVariable(); + + for ( RimEclipseInputProperty* inputProperty : inputPropertyCollection->inputProperties() ) + { + // Look for input properties starting with the same name as result definition + if ( inputProperty && inputProperty->resultName().startsWith( resultName ) ) + { + RiaLogging::info( + QString( "Found missing values result for %1: %2" ).arg( resultName ).arg( inputProperty->resultName() ) ); + + RigEclipseResultAddress resultAddress( RiaDefines::ResultCatType::INPUT_PROPERTY, inputProperty->resultName() ); + caseData->results( eclipseResultDefinition->porosityModel() )->ensureKnownResultLoaded( resultAddress ); + cvf::ref resAcc = + RigResultAccessorFactory::createFromResultAddress( caseData, + gridIndex, + eclipseResultDefinition->porosityModel(), + timeStepIndex, + resultAddress ); + + return resAcc; + } + } + + return nullptr; +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +void RimFractureModelWellLogCalculator::addOverburden( RiaDefines::CurveProperty curveProperty, + const RimFractureModel* fractureModel, + std::vector& tvDepthValues, + std::vector& measuredDepthValues, + std::vector& values ) const +{ + if ( !values.empty() ) + { + double overburdenHeight = fractureModel->overburdenHeight(); + double tvdOverburdenBottom = tvDepthValues[0]; + double tvdOverburdenTop = tvdOverburdenBottom - overburdenHeight; + + double overburdenTopValue = std::numeric_limits::infinity(); + double overburdenBottomValue = std::numeric_limits::infinity(); + if ( fractureModel->burdenStrategy( curveProperty ) == RimFractureModel::BurdenStrategy::DEFAULT_VALUE ) + { + overburdenTopValue = fractureModel->getDefaultForMissingOverburdenValue( curveProperty ); + overburdenBottomValue = overburdenTopValue; + } + else + { + double gradient = fractureModel->getOverburdenGradient( curveProperty ); + overburdenBottomValue = values[0]; + overburdenTopValue = overburdenBottomValue + gradient * -overburdenHeight; + } + + // Prepend the new "fake" depth for start of overburden + tvDepthValues.insert( tvDepthValues.begin(), tvdOverburdenBottom ); + tvDepthValues.insert( tvDepthValues.begin(), tvdOverburdenTop ); + + // TODO: this is not always correct + double mdTop = measuredDepthValues[0]; + measuredDepthValues.insert( measuredDepthValues.begin(), mdTop ); + measuredDepthValues.insert( measuredDepthValues.begin(), mdTop - overburdenHeight ); + + values.insert( values.begin(), overburdenBottomValue ); + values.insert( values.begin(), overburdenTopValue ); + } +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +void RimFractureModelWellLogCalculator::addUnderburden( RiaDefines::CurveProperty curveProperty, + const RimFractureModel* fractureModel, + std::vector& tvDepthValues, + std::vector& measuredDepthValues, + std::vector& values ) const +{ + if ( !values.empty() ) + { + size_t lastIndex = tvDepthValues.size() - 1; + + double underburdenHeight = fractureModel->underburdenHeight(); + double tvdUnderburdenTop = tvDepthValues[lastIndex]; + double tvdUnderburdenBottom = tvdUnderburdenTop + underburdenHeight; + + double underburdenTopValue = std::numeric_limits::infinity(); + double underburdenBottomValue = std::numeric_limits::infinity(); + if ( fractureModel->burdenStrategy( curveProperty ) == RimFractureModel::BurdenStrategy::DEFAULT_VALUE ) + { + underburdenTopValue = fractureModel->getDefaultForMissingUnderburdenValue( curveProperty ); + underburdenBottomValue = underburdenTopValue; + } + else + { + double gradient = fractureModel->getUnderburdenGradient( curveProperty ); + underburdenTopValue = values[lastIndex]; + underburdenBottomValue = underburdenTopValue + gradient * underburdenHeight; + } + + // Append the new "fake" depth for start of underburden + tvDepthValues.push_back( tvdUnderburdenTop ); + tvDepthValues.push_back( tvdUnderburdenBottom ); + + // Append the new "fake" md + // TODO: check if this is correct??? + double mdBottom = measuredDepthValues[lastIndex]; + measuredDepthValues.push_back( mdBottom ); + measuredDepthValues.push_back( mdBottom + underburdenHeight ); + + values.push_back( underburdenTopValue ); + values.push_back( underburdenBottomValue ); + } +} diff --git a/ApplicationCode/ProjectDataModel/RimFractureModelWellLogCalculator.h b/ApplicationCode/ProjectDataModel/RimFractureModelWellLogCalculator.h new file mode 100644 index 0000000000..12191649c5 --- /dev/null +++ b/ApplicationCode/ProjectDataModel/RimFractureModelWellLogCalculator.h @@ -0,0 +1,70 @@ +///////////////////////////////////////////////////////////////////////////////// +// +// Copyright (C) 2020- 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 +// for more details. +// +///////////////////////////////////////////////////////////////////////////////// +#pragma once + +#include "RimFractureModelCalculator.h" +#include "RimFractureModelPropertyCalculator.h" + +#include "RiaFractureModelDefines.h" + +#include "cvfObject.h" + +#include + +class RigEclipseCaseData; +class RimEclipseInputPropertyCollection; +class RimEclipseResultDefinition; +class RigResultAccessor; + +class RimFractureModelWellLogCalculator : public RimFractureModelPropertyCalculator +{ +public: + RimFractureModelWellLogCalculator( RimFractureModelCalculator& fractureModelCalculator ); + + bool calculate( RiaDefines::CurveProperty curveProperty, + const RimFractureModel* fractureModel, + int timeStep, + std::vector& values, + std::vector& measuredDepthValues, + std::vector& tvDepthValues, + double& rkbDiff ) const override; + + bool isMatching( RiaDefines::CurveProperty curveProperty ) const override; + +protected: + static bool hasMissingValues( const std::vector& values ); + static void replaceMissingValues( std::vector& values, double defaultValue ); + static void replaceMissingValues( std::vector& values, const std::vector& replacementValues ); + cvf::ref findMissingValuesAccessor( RigEclipseCaseData* caseData, + RimEclipseInputPropertyCollection* inputPropertyCollection, + int gridIndex, + int timeStepIndex, + RimEclipseResultDefinition* eclipseResultDefinition ) const; + + void addOverburden( RiaDefines::CurveProperty curveProperty, + const RimFractureModel* fractureModel, + std::vector& tvDepthValues, + std::vector& measuredDepthValues, + std::vector& values ) const; + + void addUnderburden( RiaDefines::CurveProperty curveProperty, + const RimFractureModel* fractureModel, + std::vector& tvDepthValues, + std::vector& measuredDepthValues, + std::vector& values ) const; +}; diff --git a/ApplicationCode/ProjectDataModel/RimLayerCurve.cpp b/ApplicationCode/ProjectDataModel/RimLayerCurve.cpp index 21020b4096..c5aca6a8e5 100644 --- a/ApplicationCode/ProjectDataModel/RimLayerCurve.cpp +++ b/ApplicationCode/ProjectDataModel/RimLayerCurve.cpp @@ -15,45 +15,14 @@ // for more details. // ///////////////////////////////////////////////////////////////////////////////// - #include "RimLayerCurve.h" #include "RigEclipseCaseData.h" -#include "RigEclipseWellLogExtractor.h" -#include "RigElasticProperties.h" -#include "RigResultAccessorFactory.h" -#include "RigWellLogCurveData.h" -#include "RigWellPath.h" -#include "RimCase.h" -#include "RimColorLegend.h" -#include "RimColorLegendCollection.h" -#include "RimColorLegendItem.h" #include "RimEclipseCase.h" -#include "RimEclipseResultDefinition.h" -#include "RimElasticProperties.h" #include "RimFractureModel.h" -#include "RimFractureModelPlot.h" +#include "RimFractureModelCalculator.h" #include "RimModeledWellPath.h" -#include "RimProject.h" -#include "RimTools.h" -#include "RimWellLogFile.h" -#include "RimWellLogPlot.h" -#include "RimWellLogTrack.h" -#include "RimWellPath.h" -#include "RimWellPathCollection.h" -#include "RimWellPlotTools.h" - -#include "RiuQwtPlotCurve.h" -#include "RiuQwtPlotWidget.h" - -#include "RiaApplication.h" -#include "RiaLogging.h" -#include "RiaPreferences.h" - -#include "cafPdmUiTreeOrdering.h" - -#include CAF_PDM_SOURCE_INIT( RimLayerCurve, "LayerCurve" ); @@ -108,82 +77,17 @@ void RimLayerCurve::performDataExtraction( bool* isUsingPseudoLength ) RimEclipseCase* eclipseCase = dynamic_cast( m_case.value() ); if ( eclipseCase && m_fractureModel ) { - RigEclipseWellLogExtractor eclExtractor( eclipseCase->eclipseCaseData(), - m_fractureModel->thicknessDirectionWellPath()->wellPathGeometry(), - "fracture model" ); - - rkbDiff = eclExtractor.wellPathData()->rkbDiff(); - - // Extract formation data - cvf::ref formationResultAccessor = RigResultAccessorFactory:: - createFromResultAddress( eclipseCase->eclipseCaseData(), - 0, - RiaDefines::PorosityModelType::MATRIX_MODEL, - 0, - RigEclipseResultAddress( RiaDefines::ResultCatType::FORMATION_NAMES, - RiaDefines::activeFormationNamesResultName() ) ); - if ( !formationResultAccessor.notNull() ) + bool isOk = m_fractureModel->calculator()->extractCurveData( curveProperty(), + m_timeStep, + values, + measuredDepthValues, + tvDepthValues, + rkbDiff ); + if ( !isOk ) { - RiaLogging::error( QString( "No formation result found." ) ); return; } - CurveSamplingPointData curveData = - RimWellLogTrack::curveSamplingPointData( &eclExtractor, formationResultAccessor.p() ); - - std::vector formationNamesVector = RimWellLogTrack::formationNamesVector( eclipseCase ); - - double overburdenHeight = m_fractureModel->overburdenHeight(); - if ( overburdenHeight > 0.0 ) - { - RimWellLogTrack::addOverburden( formationNamesVector, curveData, overburdenHeight ); - } - - double underburdenHeight = m_fractureModel->underburdenHeight(); - if ( underburdenHeight > 0.0 ) - { - RimWellLogTrack::addUnderburden( formationNamesVector, curveData, underburdenHeight ); - } - - measuredDepthValues = curveData.md; - tvDepthValues = curveData.tvd; - - // Extract facies data - RimFractureModelPlot* fractureModelPlot; - firstAncestorOrThisOfType( fractureModelPlot ); - if ( !fractureModelPlot ) - { - RiaLogging::error( QString( "No facies data found for layer curve." ) ); - return; - } - - std::vector faciesValues; - fractureModelPlot->getFaciesValues( faciesValues ); - if ( faciesValues.empty() ) - { - RiaLogging::error( QString( "Empty facies data found for layer curve." ) ); - return; - } - - if ( faciesValues.size() != curveData.data.size() || faciesValues.size() != measuredDepthValues.size() ) return; - - values.resize( faciesValues.size() ); - - int layerNo = 0; - double previousFormation = -1.0; - double previousFacies = -1.0; - for ( size_t i = 0; i < faciesValues.size(); i++ ) - { - if ( previousFormation != curveData.data[i] || previousFacies != faciesValues[i] ) - { - layerNo++; - } - - values[i] = layerNo; - previousFormation = curveData.data[i]; - previousFacies = faciesValues[i]; - } - RiaEclipseUnitTools::UnitSystem eclipseUnitsType = eclipseCase->eclipseCaseData()->unitsType(); if ( eclipseUnitsType == RiaEclipseUnitTools::UnitSystem::UNITS_FIELD ) { diff --git a/ApplicationCode/ProjectDataModelCommands/RimcFractureModelPlot.cpp b/ApplicationCode/ProjectDataModelCommands/RimcFractureModelPlot.cpp index 97bb9d5c99..13eb260139 100644 --- a/ApplicationCode/ProjectDataModelCommands/RimcFractureModelPlot.cpp +++ b/ApplicationCode/ProjectDataModelCommands/RimcFractureModelPlot.cpp @@ -44,7 +44,7 @@ caf::PdmObjectHandle* RimcFractureModelPlot_exportToFile::execute() { RimFractureModelPlot* fractureModelPlot = self(); - RifFractureModelPlotExporter::writeToDirectory( fractureModelPlot, + RifFractureModelPlotExporter::writeToDirectory( fractureModelPlot->fractureModel(), fractureModelPlot->fractureModel()->useDetailedFluidLoss(), m_directoryPath() );