From ef90a698b45aced3334ee3aa847a970cd32e40cf Mon Sep 17 00:00:00 2001 From: Kristian Bendiksen Date: Wed, 15 Nov 2023 10:15:16 +0100 Subject: [PATCH] Fault reactivation: compute and export stress. --- .../RifFaultReactivationModelExporter.cpp | 61 ++++-- .../RifFaultReactivationModelExporter.h | 3 +- .../Faults/CMakeLists_files.cmake | 2 + .../Faults/RimFaultReactivationDataAccess.cpp | 39 +++- .../Faults/RimFaultReactivationDataAccessor.h | 5 +- ...imFaultReactivationDataAccessorGeoMech.cpp | 3 +- .../RimFaultReactivationDataAccessorGeoMech.h | 4 +- ...ltReactivationDataAccessorPorePressure.cpp | 2 +- ...aultReactivationDataAccessorPorePressure.h | 4 +- ...RimFaultReactivationDataAccessorStress.cpp | 179 ++++++++++++++++++ .../RimFaultReactivationDataAccessorStress.h | 71 +++++++ ...ultReactivationDataAccessorTemperature.cpp | 2 +- ...FaultReactivationDataAccessorTemperature.h | 4 +- ...FaultReactivationDataAccessorVoidRatio.cpp | 2 +- ...imFaultReactivationDataAccessorVoidRatio.h | 4 +- .../Faults/RimFaultReactivationEnums.h | 8 +- .../Faults/RimFaultReactivationModel.cpp | 10 + .../Faults/RimFaultReactivationModel.h | 2 + .../WellPath/RimWellIADataAccess.cpp | 17 +- .../WellPath/RimWellIADataAccess.h | 6 + 20 files changed, 398 insertions(+), 30 deletions(-) create mode 100644 ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccessorStress.cpp create mode 100644 ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccessorStress.h diff --git a/ApplicationLibCode/FileInterface/RifFaultReactivationModelExporter.cpp b/ApplicationLibCode/FileInterface/RifFaultReactivationModelExporter.cpp index c454d799e8..fbee810856 100644 --- a/ApplicationLibCode/FileInterface/RifFaultReactivationModelExporter.cpp +++ b/ApplicationLibCode/FileInterface/RifFaultReactivationModelExporter.cpp @@ -80,6 +80,7 @@ std::pair RifFaultReactivationModelExporter::exportToStream( bool useGridTemperature = rimModel.useGridTemperature(); bool useGridDensity = rimModel.useGridDensity(); bool useGridElasticProperties = rimModel.useGridElasticProperties(); + bool useGridStress = rimModel.useGridStress(); auto dataAccess = rimModel.dataAccess(); @@ -96,7 +97,7 @@ std::pair RifFaultReactivationModelExporter::exportToStream( }, [&]() { return printInteractionProperties( stream, faultFriction ); }, [&]() { return printBoundaryConditions( stream, *model, partNames, boundaries ); }, - [&]() { return printPredefinedFields( stream, *model, *dataAccess, exportDirectory, partNames, useGridVoidRatio ); }, + [&]() { return printPredefinedFields( stream, *model, *dataAccess, exportDirectory, partNames, useGridVoidRatio, useGridStress ); }, [&]() { return printInteractions( stream, partNames, borders ); }, [&]() { @@ -483,7 +484,8 @@ std::pair const RimFaultReactivationDataAccess& dataAccess, const std::string& exportDirectory, const std::map& partNames, - bool voidRatioFromEclipse ) + bool voidRatioFromEclipse, + bool stressFromGrid ) { // PREDEFINED FIELDS struct PredefinedField @@ -529,6 +531,37 @@ std::pair RifInpExportTools::printHeading( stream, "INCLUDE, input=" + fileName ); } + if ( stressFromGrid ) + { + std::string stressName = "STRESS"; + + // Export the stress to a separate inp file + std::string fileName = stressName + ".inp"; + std::string filePath = createFilePath( exportDirectory, fileName ); + + // Use stress from first time step + size_t timeStep = 0; + bool isOk = writePropertiesToFile( model, + dataAccess, + { RimFaultReactivation::Property::StressTop, + RimFaultReactivation::Property::DepthTop, + RimFaultReactivation::Property::StressBottom, + RimFaultReactivation::Property::DepthBottom, + RimFaultReactivation::Property::LateralStressComponentX, + RimFaultReactivation::Property::LateralStressComponentY }, + {}, + timeStep, + filePath, + partNames, + "", + "" ); + + if ( !isOk ) return { false, "Failed to create " + stressName + " file." }; + + RifInpExportTools::printHeading( stream, "Initial Conditions, TYPE=" + stressName ); + RifInpExportTools::printHeading( stream, "INCLUDE, input=" + fileName ); + } + return { true, "" }; } @@ -654,18 +687,22 @@ bool RifFaultReactivationModelExporter::writePropertiesToFile( const RigFaultRea std::ofstream stream( filePath ); if ( !stream.good() ) return false; - RifInpExportTools::printHeading( stream, "Distribution Table, name=" + tableName + "_Table" ); - std::string propertyNamesLine; - for ( size_t i = 0; i < propertyNames.size(); i++ ) + bool includeHeader = !propertyNames.empty(); + if ( includeHeader ) { - propertyNamesLine += propertyNames[i]; - if ( i != propertyNames.size() - 1 ) propertyNamesLine += ", "; + RifInpExportTools::printHeading( stream, "Distribution Table, name=" + tableName + "_Table" ); + std::string propertyNamesLine; + for ( size_t i = 0; i < propertyNames.size(); i++ ) + { + propertyNamesLine += propertyNames[i]; + if ( i != propertyNames.size() - 1 ) propertyNamesLine += ", "; + } + RifInpExportTools::printLine( stream, propertyNamesLine ); + + RifInpExportTools::printHeading( stream, "Distribution, name=" + tableName + ", location=ELEMENT, Table=" + tableName + "_Table" ); + + RifInpExportTools::printLine( stream, heading ); } - RifInpExportTools::printLine( stream, propertyNamesLine ); - - RifInpExportTools::printHeading( stream, "Distribution, name=" + tableName + ", location=ELEMENT, Table=" + tableName + "_Table" ); - - RifInpExportTools::printLine( stream, heading ); for ( auto [part, partName] : partNames ) { diff --git a/ApplicationLibCode/FileInterface/RifFaultReactivationModelExporter.h b/ApplicationLibCode/FileInterface/RifFaultReactivationModelExporter.h index ff65573fe2..9197d96396 100644 --- a/ApplicationLibCode/FileInterface/RifFaultReactivationModelExporter.h +++ b/ApplicationLibCode/FileInterface/RifFaultReactivationModelExporter.h @@ -73,7 +73,8 @@ private: const RimFaultReactivationDataAccess& dataAccess, const std::string& exportDirectory, const std::map& partNames, - bool useGridVoidRatio ); + bool useGridVoidRatio, + bool useGridStress ); static std::pair printSteps( std::ostream& stream, const RigFaultReactivationModel& model, const RimFaultReactivationDataAccess& dataAccess, diff --git a/ApplicationLibCode/ProjectDataModel/Faults/CMakeLists_files.cmake b/ApplicationLibCode/ProjectDataModel/Faults/CMakeLists_files.cmake index c7024222fa..c904cfb827 100644 --- a/ApplicationLibCode/ProjectDataModel/Faults/CMakeLists_files.cmake +++ b/ApplicationLibCode/ProjectDataModel/Faults/CMakeLists_files.cmake @@ -10,6 +10,7 @@ set(SOURCE_GROUP_HEADER_FILES ${CMAKE_CURRENT_LIST_DIR}/RimFaultReactivationDataAccessorVoidRatio.h ${CMAKE_CURRENT_LIST_DIR}/RimFaultReactivationDataAccessorTemperature.h ${CMAKE_CURRENT_LIST_DIR}/RimFaultReactivationDataAccessorGeoMech.h + ${CMAKE_CURRENT_LIST_DIR}/RimFaultReactivationDataAccessorStress.h ${CMAKE_CURRENT_LIST_DIR}/RimFaultReactivationEnums.h ) @@ -25,6 +26,7 @@ set(SOURCE_GROUP_SOURCE_FILES ${CMAKE_CURRENT_LIST_DIR}/RimFaultReactivationDataAccessorVoidRatio.cpp ${CMAKE_CURRENT_LIST_DIR}/RimFaultReactivationDataAccessorTemperature.cpp ${CMAKE_CURRENT_LIST_DIR}/RimFaultReactivationDataAccessorGeoMech.cpp + ${CMAKE_CURRENT_LIST_DIR}/RimFaultReactivationDataAccessorStress.cpp ) list(APPEND CODE_HEADER_FILES ${SOURCE_GROUP_HEADER_FILES}) diff --git a/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccess.cpp b/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccess.cpp index e3fc59b265..20b961580f 100644 --- a/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccess.cpp +++ b/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccess.cpp @@ -35,6 +35,7 @@ #include "RimFaultReactivationDataAccessor.h" #include "RimFaultReactivationDataAccessorGeoMech.h" #include "RimFaultReactivationDataAccessorPorePressure.h" +#include "RimFaultReactivationDataAccessorStress.h" #include "RimFaultReactivationDataAccessorTemperature.h" #include "RimFaultReactivationDataAccessorVoidRatio.h" #include "RimFaultReactivationEnums.h" @@ -61,6 +62,17 @@ RimFaultReactivationDataAccess::RimFaultReactivationDataAccess( RimEclipseCase* { m_accessors.push_back( std::make_shared( geoMechCase, property ) ); } + + std::vector stressProperties = { RimFaultReactivation::Property::StressTop, + RimFaultReactivation::Property::DepthTop, + RimFaultReactivation::Property::StressBottom, + RimFaultReactivation::Property::DepthBottom, + RimFaultReactivation::Property::LateralStressComponentX, + RimFaultReactivation::Property::LateralStressComponentY }; + for ( auto property : stressProperties ) + { + m_accessors.push_back( std::make_shared( geoMechCase, property ) ); + } } } @@ -104,6 +116,15 @@ std::vector RimFaultReactivationDataAccess::extractModelData( const RigF RimFaultReactivation::Property::VoidRatio, RimFaultReactivation::Property::Temperature }; + auto computeAverageDepth = []( const std::vector& positions, const std::vector& indices ) + { + double sum = 0.0; + for ( size_t idx : indices ) + sum += positions[idx].z(); + + return sum / indices.size(); + }; + std::shared_ptr accessor = getAccessor( property ); if ( accessor ) { @@ -126,9 +147,13 @@ std::vector RimFaultReactivationDataAccess::extractModelData( const RigF size_t numElements = grid->elementIndices().size(); for ( size_t elementIndex = 0; elementIndex < numElements; elementIndex++ ) { - std::vector corners = grid->elementCorners( elementIndex ); - cvf::Vec3d position = RigCaseToCaseCellMapperTools::calculateCellCenter( corners.data() ); - double value = accessor->valueAtPosition( position ); + std::vector corners = grid->elementCorners( elementIndex ); + + double topDepth = computeAverageDepth( corners, { 0, 1, 2, 3 } ); + double bottomDepth = computeAverageDepth( corners, { 4, 5, 6, 7 } ); + + cvf::Vec3d position = RigCaseToCaseCellMapperTools::calculateCellCenter( corners.data() ); + double value = accessor->valueAtPosition( position, topDepth, bottomDepth ); values.push_back( value ); } } @@ -149,7 +174,13 @@ void RimFaultReactivationDataAccess::extractModelData( const RigFaultReactivatio RimFaultReactivation::Property::Temperature, RimFaultReactivation::Property::YoungsModulus, RimFaultReactivation::Property::PoissonsRatio, - RimFaultReactivation::Property::Density }; + RimFaultReactivation::Property::Density, + RimFaultReactivation::Property::StressTop, + RimFaultReactivation::Property::DepthTop, + RimFaultReactivation::Property::StressBottom, + RimFaultReactivation::Property::DepthBottom, + RimFaultReactivation::Property::LateralStressComponentX, + RimFaultReactivation::Property::LateralStressComponentY }; for ( auto property : properties ) { diff --git a/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccessor.h b/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccessor.h index 06a329b54a..2a125fb1b8 100644 --- a/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccessor.h +++ b/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccessor.h @@ -22,6 +22,7 @@ #include "RimFaultReactivationEnums.h" +#include #include class RigMainGrid; @@ -40,7 +41,9 @@ public: virtual bool isMatching( RimFaultReactivation::Property property ) const = 0; - virtual double valueAtPosition( const cvf::Vec3d& position ) const = 0; + virtual double valueAtPosition( const cvf::Vec3d& position, + double topDepth = std::numeric_limits::infinity(), + double bottomDepth = std::numeric_limits::infinity() ) const = 0; virtual bool hasValidDataAtPosition( const cvf::Vec3d& position ) const = 0; diff --git a/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccessorGeoMech.cpp b/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccessorGeoMech.cpp index 2cc9da3406..5465fd67db 100644 --- a/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccessorGeoMech.cpp +++ b/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccessorGeoMech.cpp @@ -97,7 +97,8 @@ bool RimFaultReactivationDataAccessorGeoMech::isMatching( RimFaultReactivation:: //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- -double RimFaultReactivationDataAccessorGeoMech::valueAtPosition( const cvf::Vec3d& position ) const +double RimFaultReactivationDataAccessorGeoMech::valueAtPosition( const cvf::Vec3d& position, double topDepth, double bottomDepth ) const + { if ( !m_resultFrames ) return std::numeric_limits::infinity(); diff --git a/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccessorGeoMech.h b/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccessorGeoMech.h index df9ea703dd..99efa722f2 100644 --- a/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccessorGeoMech.h +++ b/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccessorGeoMech.h @@ -40,7 +40,9 @@ public: bool isMatching( RimFaultReactivation::Property property ) const override; - double valueAtPosition( const cvf::Vec3d& position ) const override; + double valueAtPosition( const cvf::Vec3d& position, + double topDepth = std::numeric_limits::infinity(), + double bottomDepth = std::numeric_limits::infinity() ) const override; bool hasValidDataAtPosition( const cvf::Vec3d& position ) const override; diff --git a/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccessorPorePressure.cpp b/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccessorPorePressure.cpp index c7b5d5e77f..67a137a8ce 100644 --- a/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccessorPorePressure.cpp +++ b/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccessorPorePressure.cpp @@ -86,7 +86,7 @@ bool RimFaultReactivationDataAccessorPorePressure::isMatching( RimFaultReactivat //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- -double RimFaultReactivationDataAccessorPorePressure::valueAtPosition( const cvf::Vec3d& position ) const +double RimFaultReactivationDataAccessorPorePressure::valueAtPosition( const cvf::Vec3d& position, double topDepth, double bottomDepth ) const { if ( ( m_mainGrid != nullptr ) && m_resultAccessor.notNull() ) { diff --git a/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccessorPorePressure.h b/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccessorPorePressure.h index 0e58c4bd4d..bd80a0e34d 100644 --- a/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccessorPorePressure.h +++ b/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccessorPorePressure.h @@ -41,7 +41,9 @@ public: bool isMatching( RimFaultReactivation::Property property ) const override; - double valueAtPosition( const cvf::Vec3d& position ) const override; + double valueAtPosition( const cvf::Vec3d& position, + double topDepth = std::numeric_limits::infinity(), + double bottomDepth = std::numeric_limits::infinity() ) const override; bool hasValidDataAtPosition( const cvf::Vec3d& position ) const override; diff --git a/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccessorStress.cpp b/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccessorStress.cpp new file mode 100644 index 0000000000..c6596816e5 --- /dev/null +++ b/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccessorStress.cpp @@ -0,0 +1,179 @@ +///////////////////////////////////////////////////////////////////////////////// +// +// Copyright (C) 2023 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 "RimFaultReactivationDataAccessorStress.h" + +#include "RiaEclipseUnitTools.h" + +#include "RigFemAddressDefines.h" +#include "RigFemPartCollection.h" +#include "RigFemPartResultsCollection.h" +#include "RigFemResultAddress.h" +#include "RigFemScalarResultFrames.h" +#include "RigGeoMechCaseData.h" +#include "RigResultAccessorFactory.h" + +#include "RimFaultReactivationEnums.h" +#include "RimGeoMechCase.h" +#include "RimWellIADataAccess.h" + +#include +#include + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +RimFaultReactivationDataAccessorStress::RimFaultReactivationDataAccessorStress( RimGeoMechCase* geoMechCase, + RimFaultReactivation::Property property ) + : m_geoMechCase( geoMechCase ) + , m_property( property ) +{ + m_geoMechCaseData = geoMechCase->geoMechData(); +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +RimFaultReactivationDataAccessorStress::~RimFaultReactivationDataAccessorStress() +{ +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +void RimFaultReactivationDataAccessorStress::updateResultAccessor() +{ + const int partIndex = 0; + + auto loadFrameLambda = [&]( auto femParts, RigFemResultAddress addr ) -> RigFemScalarResultFrames* + { + auto result = femParts->findOrLoadScalarResult( partIndex, addr ); + if ( result->frameData( 0, 0 ).empty() ) + { + return nullptr; + } + return result; + }; + + auto femParts = m_geoMechCaseData->femPartResults(); + m_femPart = femParts->parts()->part( partIndex ); + m_s33Frames = loadFrameLambda( femParts, getResultAddress( "ST", "S33" ) ); + m_s11Frames = loadFrameLambda( femParts, getResultAddress( "ST", "S11" ) ); + m_s22Frames = loadFrameLambda( femParts, getResultAddress( "ST", "S22" ) ); + m_porFrames = loadFrameLambda( femParts, RigFemAddressDefines::elementNodalPorBarAddress() ); +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +RigFemResultAddress RimFaultReactivationDataAccessorStress::getResultAddress( const std::string& fieldName, const std::string& componentName ) +{ + return RigFemResultAddress( RIG_ELEMENT_NODAL, fieldName, componentName ); +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +bool RimFaultReactivationDataAccessorStress::isMatching( RimFaultReactivation::Property property ) const +{ + return property == m_property; +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +double RimFaultReactivationDataAccessorStress::valueAtPosition( const cvf::Vec3d& position, double topDepth, double bottomDepth ) const +{ + if ( !m_porFrames || !m_s11Frames || !m_s22Frames || !m_s33Frames || !m_femPart ) return std::numeric_limits::infinity(); + + RimWellIADataAccess iaDataAccess( m_geoMechCase ); + int centerElementIdx = iaDataAccess.elementIndex( position ); + + cvf::Vec3d topPosition( position.x(), position.y(), topDepth ); + int topElementIdx = iaDataAccess.elementIndex( topPosition ); + + cvf::Vec3d bottomPosition( position.x(), position.y(), bottomDepth ); + int bottomElementIdx = iaDataAccess.elementIndex( bottomPosition ); + if ( centerElementIdx != -1 && topElementIdx != -1 && bottomElementIdx != -1 ) + { + int timeStepIndex = 0; + int frameIndex = 0; + + const std::vector& s11Data = m_s11Frames->frameData( timeStepIndex, frameIndex ); + const std::vector& s22Data = m_s22Frames->frameData( timeStepIndex, frameIndex ); + const std::vector& s33Data = m_s33Frames->frameData( timeStepIndex, frameIndex ); + const std::vector& porData = m_porFrames->frameData( timeStepIndex, frameIndex ); + + if ( m_property == RimFaultReactivation::Property::StressTop ) + { + double s33 = interpolatedResultValue( iaDataAccess, m_femPart, topPosition, s33Data ); + double porBar = interpolatedResultValue( iaDataAccess, m_femPart, topPosition, porData ); + return RiaEclipseUnitTools::barToPascal( s33 - porBar ); + } + else if ( m_property == RimFaultReactivation::Property::StressBottom ) + { + double s33 = interpolatedResultValue( iaDataAccess, m_femPart, bottomPosition, s33Data ); + double porBar = interpolatedResultValue( iaDataAccess, m_femPart, bottomPosition, porData ); + return RiaEclipseUnitTools::barToPascal( s33 - porBar ); + } + else if ( m_property == RimFaultReactivation::Property::DepthTop ) + { + return topDepth; + } + else if ( m_property == RimFaultReactivation::Property::DepthBottom ) + { + return bottomDepth; + } + else if ( m_property == RimFaultReactivation::Property::LateralStressComponentX ) + { + double s11 = interpolatedResultValue( iaDataAccess, m_femPart, position, s11Data ); + double s33 = interpolatedResultValue( iaDataAccess, m_femPart, position, s33Data ); + double porBar = interpolatedResultValue( iaDataAccess, m_femPart, position, porData ); + return ( s11 - porBar ) / ( s33 - porBar ); + } + else if ( m_property == RimFaultReactivation::Property::LateralStressComponentY ) + { + double s22 = interpolatedResultValue( iaDataAccess, m_femPart, position, s22Data ); + double s33 = interpolatedResultValue( iaDataAccess, m_femPart, position, s33Data ); + double porBar = interpolatedResultValue( iaDataAccess, m_femPart, position, porData ); + return ( s22 - porBar ) / ( s33 - porBar ); + } + } + + return std::numeric_limits::infinity(); +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +bool RimFaultReactivationDataAccessorStress::hasValidDataAtPosition( const cvf::Vec3d& position ) const +{ + double value = valueAtPosition( position ); + return !std::isinf( value ); +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +double RimFaultReactivationDataAccessorStress::interpolatedResultValue( RimWellIADataAccess& iaDataAccess, + const RigFemPart* femPart, + const cvf::Vec3d& position, + const std::vector& scalarResults ) const +{ + return iaDataAccess.interpolatedResultValue( femPart, scalarResults, RIG_ELEMENT_NODAL, position ); +} diff --git a/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccessorStress.h b/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccessorStress.h new file mode 100644 index 0000000000..53a87f1bb4 --- /dev/null +++ b/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccessorStress.h @@ -0,0 +1,71 @@ +///////////////////////////////////////////////////////////////////////////////// +// +// Copyright (C) 2023 - 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 "RimFaultReactivationDataAccessor.h" +#include "RimFaultReactivationEnums.h" + +#include "RigFemResultAddress.h" + +#include + +class RigFemPartResultsCollection; +class RimGeoMechCase; +class RigGeoMechCaseData; +class RigFemScalarResultFrames; +class RigFemPart; +class RimWellIADataAccess; + +//================================================================================================== +/// +/// +//================================================================================================== +class RimFaultReactivationDataAccessorStress : public RimFaultReactivationDataAccessor +{ +public: + RimFaultReactivationDataAccessorStress( RimGeoMechCase* geoMechCase, RimFaultReactivation::Property property ); + ~RimFaultReactivationDataAccessorStress(); + + bool isMatching( RimFaultReactivation::Property property ) const override; + + double valueAtPosition( const cvf::Vec3d& position, + double topDepth = std::numeric_limits::infinity(), + double bottomDepth = std::numeric_limits::infinity() ) const override; + + bool hasValidDataAtPosition( const cvf::Vec3d& position ) const override; + +private: + void updateResultAccessor() override; + + static RigFemResultAddress getResultAddress( const std::string& fieldName, const std::string& componentName ); + + double interpolatedResultValue( RimWellIADataAccess& iaDataAccess, + const RigFemPart* femPart, + const cvf::Vec3d& position, + const std::vector& scalarResults ) const; + + RimGeoMechCase* m_geoMechCase; + RimFaultReactivation::Property m_property; + RigGeoMechCaseData* m_geoMechCaseData; + RigFemScalarResultFrames* m_s11Frames; + RigFemScalarResultFrames* m_s22Frames; + RigFemScalarResultFrames* m_s33Frames; + RigFemScalarResultFrames* m_porFrames; + const RigFemPart* m_femPart; +}; diff --git a/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccessorTemperature.cpp b/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccessorTemperature.cpp index 625f79aa3f..49364a25d8 100644 --- a/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccessorTemperature.cpp +++ b/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccessorTemperature.cpp @@ -83,7 +83,7 @@ bool RimFaultReactivationDataAccessorTemperature::isMatching( RimFaultReactivati //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- -double RimFaultReactivationDataAccessorTemperature::valueAtPosition( const cvf::Vec3d& position ) const +double RimFaultReactivationDataAccessorTemperature::valueAtPosition( const cvf::Vec3d& position, double topDepth, double bottomDepth ) const { if ( ( m_mainGrid != nullptr ) && m_resultAccessor.notNull() ) { diff --git a/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccessorTemperature.h b/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccessorTemperature.h index 41fd668165..25c828b067 100644 --- a/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccessorTemperature.h +++ b/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccessorTemperature.h @@ -41,7 +41,9 @@ public: bool isMatching( RimFaultReactivation::Property property ) const override; - double valueAtPosition( const cvf::Vec3d& position ) const override; + double valueAtPosition( const cvf::Vec3d& position, + double topDepth = std::numeric_limits::infinity(), + double bottomDepth = std::numeric_limits::infinity() ) const override; bool hasValidDataAtPosition( const cvf::Vec3d& position ) const override; diff --git a/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccessorVoidRatio.cpp b/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccessorVoidRatio.cpp index c72960e5ea..f0c7e91e5a 100644 --- a/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccessorVoidRatio.cpp +++ b/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccessorVoidRatio.cpp @@ -84,7 +84,7 @@ bool RimFaultReactivationDataAccessorVoidRatio::isMatching( RimFaultReactivation //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- -double RimFaultReactivationDataAccessorVoidRatio::valueAtPosition( const cvf::Vec3d& position ) const +double RimFaultReactivationDataAccessorVoidRatio::valueAtPosition( const cvf::Vec3d& position, double topDepth, double bottomDepth ) const { if ( ( m_mainGrid != nullptr ) && m_resultAccessor.notNull() ) { diff --git a/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccessorVoidRatio.h b/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccessorVoidRatio.h index edf6713e33..0b8eafca88 100644 --- a/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccessorVoidRatio.h +++ b/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccessorVoidRatio.h @@ -41,7 +41,9 @@ public: bool isMatching( RimFaultReactivation::Property property ) const override; - double valueAtPosition( const cvf::Vec3d& position ) const override; + double valueAtPosition( const cvf::Vec3d& position, + double topDepth = std::numeric_limits::infinity(), + double bottomDepth = std::numeric_limits::infinity() ) const override; bool hasValidDataAtPosition( const cvf::Vec3d& position ) const override; diff --git a/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationEnums.h b/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationEnums.h index f7c6804eca..d9c140714c 100644 --- a/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationEnums.h +++ b/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationEnums.h @@ -55,7 +55,13 @@ enum class Property Temperature, Density, YoungsModulus, - PoissonsRatio + PoissonsRatio, + StressTop, + StressBottom, + DepthTop, + DepthBottom, + LateralStressComponentX, + LateralStressComponentY }; } // namespace RimFaultReactivation diff --git a/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationModel.cpp b/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationModel.cpp index 4abab61950..32971a7eae 100644 --- a/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationModel.cpp +++ b/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationModel.cpp @@ -123,6 +123,7 @@ RimFaultReactivationModel::RimFaultReactivationModel() CAF_PDM_InitField( &m_useGridTemperature, "UseGridTemperature", true, "Output Grid Temperature" ); CAF_PDM_InitField( &m_useGridDensity, "UseGridDensity", false, "Output Grid Density" ); CAF_PDM_InitField( &m_useGridElasticProperties, "UseGridElasticProperties", false, "Output Grid Elastic Properties" ); + CAF_PDM_InitField( &m_useGridStress, "UseGridStress", false, "Output Grid Stress" ); CAF_PDM_InitFieldNoDefault( &m_targets, "Targets", "Targets" ); m_targets.uiCapability()->setUiEditorTypeName( caf::PdmUiTableViewEditor::uiEditorTypeName() ); @@ -462,6 +463,7 @@ void RimFaultReactivationModel::defineUiOrdering( QString uiConfigName, caf::Pdm propertiesGrp->add( &m_useGridTemperature ); propertiesGrp->add( &m_useGridDensity ); propertiesGrp->add( &m_useGridElasticProperties ); + propertiesGrp->add( &m_useGridStress ); auto trgGroup = uiOrdering.addNewGroup( "Debug" ); trgGroup->setCollapsedByDefault(); @@ -781,3 +783,11 @@ bool RimFaultReactivationModel::useGridElasticProperties() const { return m_useGridElasticProperties(); } + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +bool RimFaultReactivationModel::useGridStress() const +{ + return m_useGridStress(); +} diff --git a/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationModel.h b/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationModel.h index 0760fb3d9d..23c373f78b 100644 --- a/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationModel.h +++ b/ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationModel.h @@ -126,6 +126,7 @@ public: bool useGridTemperature() const; bool useGridDensity() const; bool useGridElasticProperties() const; + bool useGridStress() const; protected: caf::PdmFieldHandle* userDescriptionField() override; @@ -178,6 +179,7 @@ private: caf::PdmField m_useGridTemperature; caf::PdmField m_useGridDensity; caf::PdmField m_useGridElasticProperties; + caf::PdmField m_useGridStress; caf::PdmField m_startCellIndex; caf::PdmField m_startCellFace; diff --git a/ApplicationLibCode/ProjectDataModel/WellPath/RimWellIADataAccess.cpp b/ApplicationLibCode/ProjectDataModel/WellPath/RimWellIADataAccess.cpp index d9cdbe76a4..c2bcf99ee3 100644 --- a/ApplicationLibCode/ProjectDataModel/WellPath/RimWellIADataAccess.cpp +++ b/ApplicationLibCode/ProjectDataModel/WellPath/RimWellIADataAccess.cpp @@ -117,15 +117,26 @@ double RimWellIADataAccess::interpolatedResultValue( QString fieldNa { RigFemResultAddress address( resultType, fieldName.toStdString(), componentName.toStdString() ); + RigFemPart* femPart = m_caseData->femParts()->part( 0 ); + const std::vector& scalarResults = m_caseData->femPartResults()->resultValues( address, 0, timeStep, frameId ); + + return interpolatedResultValue( femPart, scalarResults, resultType, position ); +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +double RimWellIADataAccess::interpolatedResultValue( const RigFemPart* femPart, + const std::vector& scalarResults, + RigFemResultPosEnum resultType, + const cvf::Vec3d& position ) +{ int elmIdx = elementIndex( position ); - RigFemPart* femPart = m_caseData->femParts()->part( 0 ); RigElementType elmType = femPart->elementType( elmIdx ); const int* elementConn = femPart->connectivities( elmIdx ); int elmNodeCount = RigFemTypes::elementNodeCount( elmType ); - const std::vector& scalarResults = m_caseData->femPartResults()->resultValues( address, 0, timeStep, frameId ); - std::array nodeResults; std::array nodeCorners; diff --git a/ApplicationLibCode/ProjectDataModel/WellPath/RimWellIADataAccess.h b/ApplicationLibCode/ProjectDataModel/WellPath/RimWellIADataAccess.h index 71fee0a183..4f240b02a2 100644 --- a/ApplicationLibCode/ProjectDataModel/WellPath/RimWellIADataAccess.h +++ b/ApplicationLibCode/ProjectDataModel/WellPath/RimWellIADataAccess.h @@ -28,6 +28,7 @@ class RimGeoMechCase; class RigGeoMechCaseData; +class RigFemPart; //================================================================================================== /// @@ -49,6 +50,11 @@ public: int timeStep, int frameId ); + double interpolatedResultValue( const RigFemPart* femPart, + const std::vector& scalarResults, + RigFemResultPosEnum resultType, + const cvf::Vec3d& position ); + private: RimGeoMechCase* m_case; RigGeoMechCaseData* m_caseData;