Work in progress

ApplicationLibCode/ProjectDataModel/Faults/CMakeLists_files.cmake

@@ -4,6 +4,7 @@ set(SOURCE_GROUP_HEADER_FILES
     ${CMAKE_CURRENT_LIST_DIR}/RimFaultReactivationModel.h
     ${CMAKE_CURRENT_LIST_DIR}/RimFaultReactivationModelCollection.h
     ${CMAKE_CURRENT_LIST_DIR}/RimFaultReactivationTools.h
+    ${CMAKE_CURRENT_LIST_DIR}/RimFaultReactivationDataAccess.h
 )
 
 set(SOURCE_GROUP_SOURCE_FILES
@@ -12,6 +13,7 @@ set(SOURCE_GROUP_SOURCE_FILES
     ${CMAKE_CURRENT_LIST_DIR}/RimFaultReactivationModel.cpp
     ${CMAKE_CURRENT_LIST_DIR}/RimFaultReactivationModelCollection.cpp
     ${CMAKE_CURRENT_LIST_DIR}/RimFaultReactivationTools.cpp
+    ${CMAKE_CURRENT_LIST_DIR}/RimFaultReactivationDataAccess.cpp
 )
 
 list(APPEND CODE_HEADER_FILES ${SOURCE_GROUP_HEADER_FILES})

ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccess.cpp

@@ -0,0 +1,94 @@
+/////////////////////////////////////////////////////////////////////////////////
+//
+//  Copyright (C) 2021 -    Equinor ASA
+//
+//  ResInsight is free software: you can redistribute it and/or modify
+//  it under the terms of the GNU General Public License as published by
+//  the Free Software Foundation, either version 3 of the License, or
+//  (at your option) any later version.
+//
+//  ResInsight is distributed in the hope that it will be useful, but WITHOUT ANY
+//  WARRANTY; without even the implied warranty of MERCHANTABILITY or
+//  FITNESS FOR A PARTICULAR PURPOSE.
+//
+//  See the GNU General Public License at <http://www.gnu.org/licenses/gpl.html>
+//  for more details.
+//
+/////////////////////////////////////////////////////////////////////////////////
+
+#include "RimFaultReactivationDataAccess.h"
+
+#include "RiaPorosityModel.h"
+
+#include "RigCaseCellResultsData.h"
+#include "RigEclipseCaseData.h"
+#include "RigEclipseResultAddress.h"
+#include "RigMainGrid.h"
+#include "RigResultAccessorFactory.h"
+
+#include "RimEclipseCase.h"
+
+//--------------------------------------------------------------------------------------------------
+///
+//--------------------------------------------------------------------------------------------------
+RimFaultReactivationDataAccess::RimFaultReactivationDataAccess( RimEclipseCase* thecase, size_t timeStepIndex )
+    : m_case( thecase )
+    , m_caseData( nullptr )
+    , m_timeStepIndex( timeStepIndex )
+{
+    if ( m_case ) m_caseData = m_case->eclipseCaseData();
+    if ( m_caseData )
+    {
+        RigEclipseResultAddress resVarAddress( RiaDefines::ResultCatType::DYNAMIC_NATIVE, "PRESSURE" );
+        m_resultAccessor = RigResultAccessorFactory::createFromResultAddress( m_caseData,
+                                                                              0,
+                                                                              RiaDefines::PorosityModelType::MATRIX_MODEL,
+                                                                              timeStepIndex,
+                                                                              resVarAddress );
+    }
+}
+
+//--------------------------------------------------------------------------------------------------
+///
+//--------------------------------------------------------------------------------------------------
+RimFaultReactivationDataAccess::~RimFaultReactivationDataAccess()
+{
+}
+
+//--------------------------------------------------------------------------------------------------
+///
+//--------------------------------------------------------------------------------------------------
+double RimFaultReactivationDataAccess::porePressureAtPosition( cvf::Vec3d position, double defaultPorePressureGradient )
+{
+    double retValue = 0.0;
+
+    size_t cellIdx = cvf::UNDEFINED_SIZE_T;
+    if ( ( m_case != nullptr ) && ( m_caseData != nullptr ) )
+    {
+        auto grid = m_case->mainGrid();
+        if ( grid != nullptr ) cellIdx = grid->findReservoirCellIndexFromPoint( position );
+
+        if ( ( cellIdx != cvf::UNDEFINED_SIZE_T ) && m_resultAccessor.notNull() )
+        {
+            return m_resultAccessor->cellScalar( cellIdx );
+        }
+    }
+
+    return calculatePorePressure( std::abs( position.z() ), defaultPorePressureGradient );
+}
+
+//--------------------------------------------------------------------------------------------------
+///
+//--------------------------------------------------------------------------------------------------
+double RimFaultReactivationDataAccess::calculatePorePressure( double depth, double gradient )
+{
+    return gradient * 9.81 * depth * 1000;
+}
+
+//--------------------------------------------------------------------------------------------------
+///
+//--------------------------------------------------------------------------------------------------
+size_t RimFaultReactivationDataAccess::timeStepIndex() const
+{
+    return m_timeStepIndex;
+}

ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationDataAccess.h

@@ -0,0 +1,50 @@
+/////////////////////////////////////////////////////////////////////////////////
+//
+//  Copyright (C) 2021 -    Equinor ASA
+//
+//  ResInsight is free software: you can redistribute it and/or modify
+//  it under the terms of the GNU General Public License as published by
+//  the Free Software Foundation, either version 3 of the License, or
+//  (at your option) any later version.
+//
+//  ResInsight is distributed in the hope that it will be useful, but WITHOUT ANY
+//  WARRANTY; without even the implied warranty of MERCHANTABILITY or
+//  FITNESS FOR A PARTICULAR PURPOSE.
+//
+//  See the GNU General Public License at <http://www.gnu.org/licenses/gpl.html>
+//  for more details.
+//
+/////////////////////////////////////////////////////////////////////////////////
+
+#pragma once
+
+#include "cvfObject.h"
+#include "cvfVector3.h"
+
+class RimEclipseCase;
+class RigEclipseCaseData;
+class RigResultAccessor;
+
+//==================================================================================================
+///
+///
+//==================================================================================================
+class RimFaultReactivationDataAccess
+{
+public:
+    RimFaultReactivationDataAccess( RimEclipseCase* thecase, size_t timeStepIndex );
+    ~RimFaultReactivationDataAccess();
+
+    double porePressureAtPosition( cvf::Vec3d position, double defaultPorePressureGradient );
+
+    size_t timeStepIndex() const;
+
+protected:
+    double calculatePorePressure( double depth, double gradient );
+
+private:
+    RimEclipseCase*             m_case;
+    RigEclipseCaseData*         m_caseData;
+    size_t                      m_timeStepIndex;
+    cvf::ref<RigResultAccessor> m_resultAccessor;
+};

ApplicationLibCode/ProjectDataModel/Faults/RimFaultReactivationModel.cpp

@@ -41,6 +41,7 @@
 #include "RimEclipseView.h"
 #include "RimFaultInView.h"
 #include "RimFaultInViewCollection.h"
+#include "RimFaultReactivationDataAccess.h"
 #include "RimFaultReactivationTools.h"
 #include "RimPolylineTarget.h"
 #include "RimTimeStepFilter.h"
@@ -615,8 +616,29 @@ bool RimFaultReactivationModel::exportModelSettings()
 //--------------------------------------------------------------------------------------------------
 bool RimFaultReactivationModel::extractAndExportModelData()
 {
+    model()->clearModelData();
+
     exportModelSettings();
 
-    // TODO - get values from eclipse and geomech models here
+    auto eCase = eclipseCase();
+
+    // get the selected time step indexes
+    std::vector<size_t> selectedTimeStepIndexes;
+    for ( auto& timeStep : selectedTimeSteps() )
+    {
+        auto idx = std::find( m_availableTimeSteps.begin(), m_availableTimeSteps.end(), timeStep );
+        if ( idx == m_availableTimeSteps.end() ) return false;
+
+        selectedTimeStepIndexes.push_back( idx - m_availableTimeSteps.begin() );
+    }
+
+    // extract data for each timestep
+    for ( auto timeStepIdx : selectedTimeStepIndexes )
+    {
+        RimFaultReactivationDataAccess dataAccess( eCase, timeStepIdx );
+
+        model()->extractModelData( &dataAccess );
+    }
+
     return true;
 }

ApplicationLibCode/ProjectDataModel/WellPath/RimWellIADataAccess.cpp

@@ -22,6 +22,8 @@
 #include "RigFemPartCollection.h"
 #include "RigFemPartResultsCollection.h"
 #include "RigGeoMechCaseData.h"
+#include "RigHexIntersectionTools.h"
+
 #include "RimGeoMechCase.h"
 
 #include "../cafHexInterpolator/cafHexInterpolator.h" // Use relative path, as this is a header only file not part of a library
@@ -67,11 +69,21 @@ int RimWellIADataAccess::elementIndex( cvf::Vec3d position )
     cvf::BoundingBox bb;
     bb.add( position );
 
+    auto part = m_caseData->femParts()->part( 0 );
+
     std::vector<size_t> closeElements;
-    m_caseData->femParts()->part( 0 )->findIntersectingElementIndices( bb, &closeElements );
+    part->findIntersectingElementIndices( bb, &closeElements );
     if ( closeElements.size() == 0 ) return -1;
 
-    return (int)closeElements[0];
+    for ( auto elmIdx : closeElements )
+    {
+        std::array<cvf::Vec3d, 8> coordinates;
+        if ( !part->fillElementCoordinates( elmIdx, coordinates ) ) continue;
+
+        if ( RigHexIntersectionTools::isPointInCell( position, coordinates.data() ) ) return (int)elmIdx;
+    }
+
+    return -1;
 }
 
 //--------------------------------------------------------------------------------------------------

ApplicationLibCode/ReservoirDataModel/RigFaultReactivationModel.cpp

@@ -101,6 +101,17 @@ void RigFaultReactivationModel::reset()
     }
 }
 
+//--------------------------------------------------------------------------------------------------
+///
+//--------------------------------------------------------------------------------------------------
+void RigFaultReactivationModel::clearModelData()
+{
+    for ( auto part : allGridParts() )
+    {
+        m_3dparts[part]->clearModelData();
+    }
+}
+
 //--------------------------------------------------------------------------------------------------
 ///
 //--------------------------------------------------------------------------------------------------
@@ -328,3 +339,14 @@ std::shared_ptr<RigGriddedPart3d> RigFaultReactivationModel::grid( GridPart part
 {
     return m_3dparts.at( part );
 }
+
+//--------------------------------------------------------------------------------------------------
+///
+//--------------------------------------------------------------------------------------------------
+void RigFaultReactivationModel::extractModelData( RimFaultReactivationDataAccess* dataAccess )
+{
+    for ( auto part : allGridParts() )
+    {
+        m_3dparts[part]->extractModelData( dataAccess );
+    }
+}

ApplicationLibCode/ReservoirDataModel/RigFaultReactivationModel.h

@@ -30,6 +30,7 @@
 #include <vector>
 
 class RigGriddedPart3d;
+class RimFaultReactivationDataAccess;
 
 class RigFRModelPart
 {
@@ -94,6 +95,9 @@ class RigFaultReactivationModel : public cvf::Object
 
     std::shared_ptr<RigGriddedPart3d> grid( GridPart part ) const;
 
+    void clearModelData();
+    void extractModelData( RimFaultReactivationDataAccess* dataAccess );
+
 protected:
     void generateGrids( cvf::Vec3dArray points );
 

ApplicationLibCode/ReservoirDataModel/RigGriddedPart3d.cpp

@@ -18,6 +18,8 @@
 
 #include "RigGriddedPart3d.h"
 
+#include "RimFaultReactivationDataAccess.h"
+
 #include "cvfTextureImage.h"
 
 //--------------------------------------------------------------------------------------------------
@@ -46,6 +48,16 @@ void RigGriddedPart3d::reset()
     m_nodes.clear();
     m_elementIndices.clear();
     m_meshLines.clear();
+
+    clearModelData();
+}
+
+//--------------------------------------------------------------------------------------------------
+///
+//--------------------------------------------------------------------------------------------------
+void RigGriddedPart3d::clearModelData()
+{
+    m_nodePorePressure.clear();
 }
 
 //--------------------------------------------------------------------------------------------------
@@ -344,3 +356,31 @@ const std::map<RigGriddedPart3d::Boundary, std::vector<unsigned int>>& RigGridde
 {
     return m_boundaryNodes;
 }
+
+//--------------------------------------------------------------------------------------------------
+///
+//--------------------------------------------------------------------------------------------------
+const std::vector<double>& RigGriddedPart3d::nodePorePressure( size_t timeStepIndex ) const
+{
+    if ( timeStepIndex >= m_nodePorePressure.size() ) return m_emptyData;
+    return m_nodePorePressure[timeStepIndex];
+}
+
+//--------------------------------------------------------------------------------------------------
+///
+//--------------------------------------------------------------------------------------------------
+void RigGriddedPart3d::extractModelData( RimFaultReactivationDataAccess* dataAccess )
+{
+    const auto timeStepIdx = dataAccess->timeStepIndex();
+
+    if ( m_nodePorePressure.size() <= timeStepIdx )
+    {
+        m_nodePorePressure.resize( timeStepIdx + 1 );
+    }
+
+    for ( auto& node : m_nodes )
+    {
+        double pressure = dataAccess->porePressureAtPosition( node, 1.0 );
+        m_nodePorePressure[timeStepIdx].push_back( pressure );
+    }
+}

ApplicationLibCode/ReservoirDataModel/RigGriddedPart3d.h

@@ -24,6 +24,8 @@
 #include <map>
 #include <vector>
 
+class RimFaultReactivationDataAccess;
+
 //==================================================================================================
 ///
 ///
@@ -51,6 +53,7 @@ class RigGriddedPart3d : public cvf::Object
     ~RigGriddedPart3d() override;
 
     void reset();
+    void clearModelData();
 
     void generateGeometry( std::vector<cvf::Vec3d> inputPoints,
                            int                     nHorzCells,
@@ -59,6 +62,8 @@ class RigGriddedPart3d : public cvf::Object
                            int                     nVertCellsUpper,
                            double                  thickness );
 
+    void extractModelData( RimFaultReactivationDataAccess* dataAccess );
+
     const std::vector<cvf::Vec3d>&                            nodes() const;
     const std::vector<std::vector<unsigned int>>&             elementIndices() const;
     const std::map<BorderSurface, std::vector<unsigned int>>& borderSurfaceElements() const;
@@ -67,6 +72,8 @@ class RigGriddedPart3d : public cvf::Object
     const std::map<Boundary, std::vector<unsigned int>>& boundaryElements() const;
     const std::map<Boundary, std::vector<unsigned int>>& boundaryNodes() const;
 
+    const std::vector<double>& nodePorePressure( size_t timeStepIndex ) const;
+
 protected:
     cvf::Vec3d stepVector( cvf::Vec3d start, cvf::Vec3d stop, int nSteps );
     void       generateMeshlines( std::vector<cvf::Vec3d> cornerPoints, int numHorzCells, int numVertCells );
@@ -80,4 +87,7 @@ class RigGriddedPart3d : public cvf::Object
     std::vector<std::vector<cvf::Vec3d>>               m_meshLines;
     std::map<Boundary, std::vector<unsigned int>>      m_boundaryElements;
     std::map<Boundary, std::vector<unsigned int>>      m_boundaryNodes;
+
+    std::vector<std::vector<double>> m_nodePorePressure;
+    const std::vector<double>        m_emptyData;
 };