Generate model based on an origin in the fault to improve calculation accuracy (#11267)

Generate model based on an origin in the fault to improve calculation accuracy.

ApplicationLibCode/ReservoirDataModel/RigFaultReactivationModelGenerator.cpp

@@ -488,6 +488,20 @@ std::pair<size_t, size_t> RigFaultReactivationModelGenerator::findCellWithInters
     return { i, j };
 }
 
+//--------------------------------------------------------------------------------------------------
+///
+//--------------------------------------------------------------------------------------------------
+std::array<cvf::Vec3d, 12> RigFaultReactivationModelGenerator::shiftOrigin( const std::array<cvf::Vec3d, 12>& points,
+                                                                            const cvf::Vec3d&                 newOrigin )
+{
+    std::array<cvf::Vec3d, 12> retPoints;
+    for ( int i = 0; i < (int)points.size(); i++ )
+    {
+        retPoints[i] = points[i] - newOrigin;
+    }
+    return retPoints;
+}
+
 //--------------------------------------------------------------------------------------------------
 ///
 //--------------------------------------------------------------------------------------------------
@@ -574,6 +588,10 @@ void RigFaultReactivationModelGenerator::generateGeometry( size_t            sta
 
     generatePointsFrontBack();
 
+    // use temp origin in start position, at zero depth
+    cvf::Vec3d origin( m_startPosition );
+    origin.z() = 0.0;
+
     cvf::Vec3d                     tVec = m_modelThickness * m_modelNormal;
     std::vector<cvf::Vec3d>        thicknessVectors;
     std::vector<caf::Line<double>> faultLines;
@@ -581,12 +599,19 @@ void RigFaultReactivationModelGenerator::generateGeometry( size_t            sta
 
     for ( int i = 0; i < 3; i++ )
     {
-        faultLines.push_back( caf::Line<double>( m_topFault + thicknessFactors[i] * tVec, m_bottomFault + thicknessFactors[i] * tVec ) );
+        faultLines.push_back(
+            caf::Line<double>( m_topFault - origin + thicknessFactors[i] * tVec, m_bottomFault - origin + thicknessFactors[i] * tVec ) );
         thicknessVectors.push_back( thicknessFactors[i] * tVec );
     }
 
-    frontPart->generateGeometry( m_frontPoints,
+    std::array<cvf::Vec3d, 12> shiftedFrontPoints = shiftOrigin( m_frontPoints, origin );
-                                 frontReservoirLayers,
+    std::array<cvf::Vec3d, 12> shiftedBackPoints  = shiftOrigin( m_backPoints, origin );
+
+    std::vector<double> frontResZ = extractZValues( frontReservoirLayers );
+    std::vector<double> backResZ  = extractZValues( backReservoirLayers );
+
+    frontPart->generateGeometry( shiftedFrontPoints,
+                                 frontResZ,
                                  m_maxCellHeight,
                                  m_cellSizeHeightFactor,
                                  m_horizontalPartition,
@@ -598,8 +623,8 @@ void RigFaultReactivationModelGenerator::generateGeometry( size_t            sta
     std::reverse( faultLines.begin(), faultLines.end() );
     std::reverse( thicknessVectors.begin(), thicknessVectors.end() );
 
-    backPart->generateGeometry( m_backPoints,
+    backPart->generateGeometry( shiftedBackPoints,
-                                backReservoirLayers,
+                                backResZ,
                                 m_maxCellHeight,
                                 m_cellSizeHeightFactor,
                                 m_horizontalPartition,
@@ -608,6 +633,9 @@ void RigFaultReactivationModelGenerator::generateGeometry( size_t            sta
                                 m_topReservoirBack.z(),
                                 m_faultZoneCells );
 
+    frontPart->shiftNodes( origin );
+    backPart->shiftNodes( origin );
+
     frontPart->generateLocalNodes( m_localCoordTransform );
     backPart->generateLocalNodes( m_localCoordTransform );
 
@@ -792,3 +820,18 @@ std::pair<double, double> RigFaultReactivationModelGenerator::depthTopBottom() c
 {
     return { -m_startDepth, m_bottomDepth };
 }
+
+//--------------------------------------------------------------------------------------------------
+///
+//--------------------------------------------------------------------------------------------------
+std::vector<double> RigFaultReactivationModelGenerator::extractZValues( const std::vector<cvf::Vec3d>& points )
+{
+    std::vector<double> layers;
+
+    for ( auto& p : points )
+    {
+        layers.push_back( p.z() );
+    }
+
+    return layers;
+}

ApplicationLibCode/ReservoirDataModel/RigFaultReactivationModelGenerator.h

@@ -74,9 +74,11 @@ protected:
     static const std::vector<double>     partition( double distance, double startSize, double sizeFactor );
     static std::pair<FaceType, FaceType> sideFacesIJ( FaceType face );
 
-    static cvf::Vec3d extrapolatePoint( cvf::Vec3d startPoint, cvf::Vec3d endPoint, double stopDepth );
+    static cvf::Vec3d                 extrapolatePoint( cvf::Vec3d startPoint, cvf::Vec3d endPoint, double stopDepth );
-    static void       splitLargeLayers( std::map<double, cvf::Vec3d>& layers, double maxHeight );
+    static void                       splitLargeLayers( std::map<double, cvf::Vec3d>& layers, double maxHeight );
-    static void       mergeTinyLayers( std::map<double, cvf::Vec3d>& layers, double minHeight );
+    static void                       mergeTinyLayers( std::map<double, cvf::Vec3d>& layers, double minHeight );
+    static std::vector<double>        extractZValues( const std::vector<cvf::Vec3d>& points );
+    static std::array<cvf::Vec3d, 12> shiftOrigin( const std::array<cvf::Vec3d, 12>& points, const cvf::Vec3d& newOrigin );
 
     std::vector<size_t> buildCellColumn( size_t startCell, FaceType startFace, std::map<double, cvf::Vec3d>& layers );
 

ApplicationLibCode/ReservoirDataModel/RigGriddedPart3d.cpp

@@ -165,21 +165,6 @@ std::vector<double> RigGriddedPart3d::generateGrowingLayers( double zFrom, doubl
     return layers;
 }
 
-//--------------------------------------------------------------------------------------------------
-///
-//--------------------------------------------------------------------------------------------------
-std::vector<double> RigGriddedPart3d::extractZValues( const std::vector<cvf::Vec3d>& points )
-{
-    std::vector<double> layers;
-
-    for ( auto& p : points )
-    {
-        layers.push_back( p.z() );
-    }
-
-    return layers;
-}
-
 //--------------------------------------------------------------------------------------------------
 ///  Point index in input
 ///
@@ -203,7 +188,7 @@ std::vector<double> RigGriddedPart3d::extractZValues( const std::vector<cvf::Vec
 ///
 //--------------------------------------------------------------------------------------------------
 void RigGriddedPart3d::generateGeometry( const std::array<cvf::Vec3d, 12>&    inputPoints,
-                                         const std::vector<cvf::Vec3d>&       reservoirLayers,
+                                         const std::vector<double>&           reservoirZ,
                                          double                               maxCellHeight,
                                          double                               cellSizeFactor,
                                          const std::vector<double>&           horizontalPartition,
@@ -220,8 +205,9 @@ void RigGriddedPart3d::generateGeometry( const std::array<cvf::Vec3d, 12>&    in
 
     layersPerRegion[Regions::LowerUnderburden] = generateGrowingLayers( inputPoints[1].z(), inputPoints[0].z(), maxCellHeight, cellSizeFactor );
     layersPerRegion[Regions::UpperUnderburden] = generateConstantLayers( inputPoints[1].z(), inputPoints[2].z(), maxCellHeight );
-    layersPerRegion[Regions::Reservoir]        = extractZValues( reservoirLayers );
+    layersPerRegion[Regions::Reservoir]        = reservoirZ;
-    layersPerRegion[Regions::LowerOverburden]  = generateConstantLayers( inputPoints[3].z(), inputPoints[4].z(), maxCellHeight );
+
+    layersPerRegion[Regions::LowerOverburden] = generateConstantLayers( inputPoints[3].z(), inputPoints[4].z(), maxCellHeight );
     layersPerRegion[Regions::UpperOverburden] = generateGrowingLayers( inputPoints[4].z(), inputPoints[5].z(), maxCellHeight, cellSizeFactor );
 
     layersPerRegion[Regions::LowerUnderburden].pop_back(); // to avoid overlap with bottom of next region
@@ -304,8 +290,10 @@ void RigGriddedPart3d::generateGeometry( const std::array<cvf::Vec3d, 12>&    in
             }
             else if ( region == Regions::Reservoir )
             {
-                toPos       = reservoirLayers[v];
+                fromPos.z() = reservoirZ[v];
-                fromPos.z() = toPos.z();
+                cvf::Plane zPlane;
+                zPlane.setFromPointAndNormal( fromPos, cvf::Vec3d::Z_AXIS );
+                zPlane.intersect( faultLines[1].start(), faultLines[1].end(), &toPos );
             }
 
             cvf::Vec3d stepHorz = toPos - fromPos;
@@ -318,14 +306,7 @@ void RigGriddedPart3d::generateGeometry( const std::array<cvf::Vec3d, 12>&    in
 
             for ( int h = 0; h <= (int)nHorzCells; h++ )
             {
-                if ( h == (int)nHorzCells )
+                p = toPos - horizontalPartition[h] * stepHorz;
-                {
-                    p = toPos;
-                }
-                else
-                {
-                    p = toPos - horizontalPartition[h] * stepHorz;
-                }
 
                 for ( int t = 0; t <= nThicknessCells; t++, nodeIndex++ )
                 {
@@ -336,7 +317,7 @@ void RigGriddedPart3d::generateGeometry( const std::array<cvf::Vec3d, 12>&    in
                          ( ( region == Regions::Reservoir ) || region == Regions::LowerOverburden || region == Regions::UpperUnderburden ) )
                     {
                         cvf::Plane zPlane;
-                        zPlane.setFromPointAndNormal( nodePoint, cvf::Vec3d::Z_AXIS );
+                        zPlane.setFromPointAndNormal( p, cvf::Vec3d::Z_AXIS );
                         zPlane.intersect( faultLines[t].start(), faultLines[t].end(), &nodePoint );
                     }
 
@@ -687,6 +668,25 @@ void RigGriddedPart3d::generateLocalNodes( const cvf::Mat4d transform )
     }
 }
 
+//--------------------------------------------------------------------------------------------------
+///
+//--------------------------------------------------------------------------------------------------
+void RigGriddedPart3d::shiftNodes( const cvf::Vec3d offset )
+{
+    for ( int i = 0; i < (int)m_nodes.size(); i++ )
+    {
+        m_nodes[i] += offset;
+    }
+
+    for ( int i = 0; i < (int)m_meshLines.size(); i++ )
+    {
+        for ( int j = 0; j < (int)m_meshLines[i].size(); j++ )
+        {
+            m_meshLines[i][j] += offset;
+        }
+    }
+}
+
 //--------------------------------------------------------------------------------------------------
 ///
 //--------------------------------------------------------------------------------------------------

ApplicationLibCode/ReservoirDataModel/RigGriddedPart3d.h

@@ -50,7 +50,7 @@ public:
     void reset();
 
     void generateGeometry( const std::array<cvf::Vec3d, 12>&    inputPoints,
-                           const std::vector<cvf::Vec3d>&       reservoirLayers,
+                           const std::vector<double>&           reservoirZ,
                            double                               maxCellHeight,
                            double                               cellSizeFactor,
                            const std::vector<double>&           horizontalPartition,
@@ -59,6 +59,7 @@ public:
                            double                               topHeight,
                            int                                  nFaultZoneCells );
 
+    void shiftNodes( const cvf::Vec3d offset );
     void generateLocalNodes( const cvf::Mat4d transform );
     void setUseLocalCoordinates( bool useLocalCoordinates );
 
@@ -89,7 +90,6 @@ protected:
     static cvf::Vec3d          stepVector( cvf::Vec3d start, cvf::Vec3d stop, int nSteps );
     static std::vector<double> generateConstantLayers( double zFrom, double zTo, double maxSize );
     static std::vector<double> generateGrowingLayers( double zFrom, double zTo, double maxSize, double growfactor );
-    static std::vector<double> extractZValues( const std::vector<cvf::Vec3d>& points );
 
     void generateVerticalMeshlines( const std::vector<cvf::Vec3d>& cornerPoints, const std::vector<double>& horzPartition );