#1091 - pre-proto - Expanding the findCellAverageZdirection function to find local X, Y and Z for a cell.

ApplicationCode/ProjectDataModel/RimFracture.cpp

@@ -269,10 +269,15 @@ void RimFracture::computeTransmissibility()
 
 
         std::vector<std::vector<cvf::Vec3d> > planeCellPolygons;
-        bool isPlanIntersected = planeCellIntersectionPolygons(fracCell, planeCellPolygons);
+        cvf::Vec3d localX;
+        cvf::Vec3d localY;
+        cvf::Vec3d localZ;
+
+        bool isPlanIntersected = planeCellIntersectionPolygons(fracCell, planeCellPolygons, localX, localY, localZ);
+        
         if (!isPlanIntersected || planeCellPolygons.size()==0) continue;
 
-        //Transform planCell polygon(s) to x/y coordinate system (where fracturePolygon already is located)
+        //Transform planCell polygon(s) and averageZdirection to x/y coordinate system (where fracturePolygon already is located)
         cvf::Mat4f invertedTransMatrix = transformMatrix().getInverted();
         for (std::vector<cvf::Vec3d> & planeCellPolygon : planeCellPolygons)
         {
@@ -283,6 +288,9 @@ void RimFracture::computeTransmissibility()
         }
 
 
+        //TODO: Make copy of z dir vector, we need it both in fracture coords and domain cords
+        localZ.transformVector(static_cast<cvf::Mat4d>(invertedTransMatrix));
+
         RigFractureData fracData; 
         fracData.reservoirCellIndex = fracCell;
 
@@ -323,7 +331,8 @@ void RimFracture::computeTransmissibility()
 //--------------------------------------------------------------------------------------------------
 /// 
 //--------------------------------------------------------------------------------------------------
-bool RimFracture::planeCellIntersectionPolygons(size_t cellindex, std::vector<std::vector<cvf::Vec3d> > & polygons)
+bool RimFracture::planeCellIntersectionPolygons(size_t cellindex, std::vector<std::vector<cvf::Vec3d> > & polygons, 
+                                                cvf::Vec3d & localX, cvf::Vec3d & localY, cvf::Vec3d & localZ)
 {
     
     cvf::Plane fracturePlane;
@@ -372,6 +381,8 @@ bool RimFracture::planeCellIntersectionPolygons(size_t cellindex, std::vector<st
        
     RigCellGeometryTools::createPolygonFromLineSegments(intersectionLineSegments, polygons);
 
+    RigCellGeometryTools::findCellLocalXYZ(hexCorners, localX, localY, localZ);
+
     return isCellIntersected;
 }
 

ApplicationCode/ProjectDataModel/RimFracture.h

@@ -92,7 +92,7 @@ private:
     QString                         createOneBasedIJK() const;
 
     //Functions for area calculations - should these be in separate class
-    bool planeCellIntersectionPolygons(size_t cellindex, std::vector<std::vector<cvf::Vec3d> > & polygons);
+    bool planeCellIntersectionPolygons(size_t cellindex, std::vector<std::vector<cvf::Vec3d> > & polygons, cvf::Vec3d & localX, cvf::Vec3d & localY, cvf::Vec3d & localZ);
 
     void calculateFracturePlaneCellPolygonOverlapArea(std::vector<std::vector<cvf::Vec3d> > planeCellPolygons,
         std::vector<cvf::Vec3f> fracturePolygon, double & area);

ApplicationCode/ReservoirDataModel/RigCellGeometryTools.cpp

@@ -139,10 +139,9 @@ void RigCellGeometryTools::createPolygonFromLineSegments(std::list<std::pair<cvf
 }
 
 //--------------------------------------------------------------------------------------------------
-/// Returns the a vector normal to the vectors between face centers in I and J direction 
 ///  
 //--------------------------------------------------------------------------------------------------
-void RigCellGeometryTools::findCellAverageZdirection(cvf::Vec3d * hexCorners, cvf::Vec3d &averageZdirection)
+void RigCellGeometryTools::findCellLocalXYZ(cvf::Vec3d * hexCorners, cvf::Vec3d &localXdirection, cvf::Vec3d &localYdirection, cvf::Vec3d &localZdirection)
 {
     cvf::ubyte faceVertexIndices[4];
     cvf::StructGridInterface::FaceEnum face;
@@ -150,6 +149,7 @@ void RigCellGeometryTools::findCellAverageZdirection(cvf::Vec3d * hexCorners, cv
     face = cvf::StructGridInterface::NEG_I;
     cvf::StructGridInterface::cellFaceVertexIndices(face, faceVertexIndices);
     cvf::Vec3d faceCenterNegI = cvf::GeometryTools::computeFaceCenter(hexCorners[faceVertexIndices[0]], hexCorners[faceVertexIndices[1]], hexCorners[faceVertexIndices[2]], hexCorners[faceVertexIndices[3]]);
+    //TODO: Should we use face centroids instead of face centers?
 
     face = cvf::StructGridInterface::POS_I;
     cvf::StructGridInterface::cellFaceVertexIndices(face, faceVertexIndices);
@@ -166,8 +166,18 @@ void RigCellGeometryTools::findCellAverageZdirection(cvf::Vec3d * hexCorners, cv
     cvf::Vec3d faceCenterCenterVectorI = faceCenterPosI - faceCenterNegI;
     cvf::Vec3d faceCenterCenterVectorJ = faceCenterPosJ - faceCenterNegJ;
 
-    averageZdirection.cross(faceCenterCenterVectorI, faceCenterCenterVectorJ);
+    localZdirection.cross(faceCenterCenterVectorI, faceCenterCenterVectorJ);
-    averageZdirection.normalize();
+    localZdirection.normalize();
+
+    cvf::Vec3d crossPoductJZ;
+    crossPoductJZ.cross(faceCenterCenterVectorJ, localZdirection);
+    localXdirection = faceCenterCenterVectorI + crossPoductJZ;
+    localXdirection.normalize();
+
+    cvf::Vec3d crossPoductIZ;
+    crossPoductIZ.cross(faceCenterCenterVectorI, localZdirection);
+    localYdirection = faceCenterCenterVectorJ - crossPoductIZ;
+    localYdirection.normalize();
+    
 }
 
-

ApplicationCode/ReservoirDataModel/RigCellGeometryTools.h

@@ -37,6 +37,7 @@ public:
 
     static void createPolygonFromLineSegments(std::list<std::pair<cvf::Vec3d, cvf::Vec3d>> &intersectionLineSegments, std::vector<std::vector<cvf::Vec3d>> &polygons);
 
-    static void findCellAverageZdirection(cvf::Vec3d * hexCorners, cvf::Vec3d &averageZdirection);
+    static void findCellLocalXYZ(cvf::Vec3d * hexCorners, cvf::Vec3d &localXdirection, cvf::Vec3d &localYdirection, cvf::Vec3d &localZdirection);
+
 
 };

ApplicationCode/UnitTests/RigCellGeometryTools-Test.cpp

@@ -111,24 +111,35 @@ TEST(RigCellGeometryTools, planeHexCellIntersectionTest)
 TEST(RigCellGeometryTools, findCellAverageZTest)
 {
     cvf::Vec3d hexCorners[8];
+ 
     hexCorners[0] = cvf::Vec3d(0, 0, 0);
     hexCorners[1] = cvf::Vec3d(1, 0, 0);
-    hexCorners[2] = cvf::Vec3d(0, 1, 0);
+    hexCorners[2] = cvf::Vec3d(1, 1, 0);
-    hexCorners[3] = cvf::Vec3d(1, 1, 0);
+    hexCorners[3] = cvf::Vec3d(0, 1, 0);
 
     hexCorners[4] = cvf::Vec3d(0, 0, 1);
     hexCorners[5] = cvf::Vec3d(1, 0, 1);
     hexCorners[6] = cvf::Vec3d(1, 1, 1);
-    hexCorners[7] = cvf::Vec3d(1, 0, 1);
+    hexCorners[7] = cvf::Vec3d(0, 1, 1);
 
-    cvf::Vec3d averageZdirection;
 
-    RigCellGeometryTools::findCellAverageZdirection(hexCorners, averageZdirection);
+    cvf::Vec3d localX;
+    cvf::Vec3d localY;
+    cvf::Vec3d localZ;
 
-    EXPECT_DOUBLE_EQ(averageZdirection[0], 0);
+    RigCellGeometryTools::findCellLocalXYZ(hexCorners, localX, localY, localZ);
-    EXPECT_DOUBLE_EQ(averageZdirection[1], 0);
-    EXPECT_DOUBLE_EQ(averageZdirection[2], 1);
 
+    EXPECT_DOUBLE_EQ(localX[0], 1);
+    EXPECT_DOUBLE_EQ(localX[1], 0);
+    EXPECT_DOUBLE_EQ(localX[2], 0);
+
+    EXPECT_DOUBLE_EQ(localY[0], 0);
+    EXPECT_DOUBLE_EQ(localY[1], 1);
+    EXPECT_DOUBLE_EQ(localY[2], 0);
+
+    EXPECT_DOUBLE_EQ(localZ[0], 0);
+    EXPECT_DOUBLE_EQ(localZ[1], 0);
+    EXPECT_DOUBLE_EQ(localZ[2], 1);
 
 }