Moved Faluts to MainGrid. New fault calculation

Now finding the uncovered geometrical faults and adds them to a separate generated fault.
2025-02-25 18:55:39 -06:00 · 2013-12-12 22:15:14 +01:00
parent b23490106e
commit bf6fba191e
10 changed files with 189 additions and 133 deletions
--- a/ApplicationCode/ReservoirDataModel/RigMainGrid.cpp
+++ b/ApplicationCode/ReservoirDataModel/RigMainGrid.cpp
@@ -192,3 +192,122 @@ RigNNCData* RigMainGrid::nncData()
    return m_nncData.p();
 }

+//--------------------------------------------------------------------------------------------------
+/// 
+//--------------------------------------------------------------------------------------------------
+void RigMainGrid::setFaults(const cvf::Collection<RigFault>& faults)
+{
+    m_faults = faults;
+
+#pragma omp parallel for 
+    for (int i = 0; i < m_faults.size(); i++)
+    {
+        m_faults[i]->computeFaultFacesFromCellRanges(this->mainGrid());
+    }
+}
+
+
+//--------------------------------------------------------------------------------------------------
+/// 
+//--------------------------------------------------------------------------------------------------
+void RigMainGrid::calculateFaults()
+{
+    //RigFault::initFaultsPrCellAccumulator(m_cells.size());
+    cvf::ref<RigFaultsPrCellAccumulator> faultsPrCellAcc = new RigFaultsPrCellAccumulator(m_cells.size());
+
+    // Spread fault idx'es on the cells from the faults
+    for (size_t fIdx = 0 ; fIdx < m_faults.size(); ++fIdx)
+    {
+        m_faults[fIdx]->accumulateFaultsPrCell(faultsPrCellAcc.p(), static_cast<int>(fIdx));
+    }
+
+    // Find the geometrical faults that is in addition
+
+    RigFault * unNamedFault = new RigFault;
+    int unNamedFaultIdx = static_cast<int>(m_faults.size());
+
+    for (size_t gcIdx = 0 ; gcIdx < m_cells.size(); ++gcIdx)
+    {
+        if ( m_cells[gcIdx].isInvalid())
+        {
+            continue;
+        }
+
+        size_t neighborGlobalCellIdx;
+        size_t neighborGridCellIdx;
+        size_t i, j, k;
+        RigGridBase* hostGrid; ;
+        bool firstNO_FAULTFaceForCell = true;
+
+        for (char faceIdx = 0; faceIdx < 6; ++faceIdx)
+        {
+            cvf::StructGridInterface::FaceType face = cvf::StructGridInterface::FaceType(faceIdx);
+
+            if (faultsPrCellAcc->faultIdx(gcIdx, face) == RigFaultsPrCellAccumulator::NO_FAULT)
+            {
+                // Find neighbor cell
+                if (firstNO_FAULTFaceForCell) // To avoid doing this for every face, and only when detecting a NO_FAULT
+                {
+                    hostGrid = m_cells[gcIdx].hostGrid();
+                    hostGrid->ijkFromCellIndex(m_cells[gcIdx].cellIndex(), &i,&j, &k);
+                    firstNO_FAULTFaceForCell = false;
+                }
+
+                if(!hostGrid->cellIJKNeighbor(i, j, k, face, &neighborGridCellIdx))
+                {
+                    continue;
+                }
+
+                neighborGlobalCellIdx = hostGrid->globalGridCellIndex(neighborGridCellIdx);
+                if (m_cells[neighborGlobalCellIdx].isInvalid())
+                {
+                    continue;
+                }
+
+                double tolerance = 1e-6;
+               
+
+                caf::SizeTArray4 faceIdxs;
+                m_cells[gcIdx].faceIndices(face, &faceIdxs);
+                caf::SizeTArray4 nbFaceIdxs;
+                m_cells[neighborGlobalCellIdx].faceIndices(StructGridInterface::oppositeFace(face), &nbFaceIdxs);
+
+                const std::vector<cvf::Vec3d>& vxs = m_mainGrid->nodes();
+
+                bool sharedFaceVertices = true;
+                if (sharedFaceVertices && vxs[faceIdxs[0]].pointDistance(vxs[nbFaceIdxs[0]]) > tolerance ) sharedFaceVertices = false;
+                if (sharedFaceVertices && vxs[faceIdxs[1]].pointDistance(vxs[nbFaceIdxs[3]]) > tolerance ) sharedFaceVertices = false;
+                if (sharedFaceVertices && vxs[faceIdxs[2]].pointDistance(vxs[nbFaceIdxs[2]]) > tolerance ) sharedFaceVertices = false;
+                if (sharedFaceVertices && vxs[faceIdxs[3]].pointDistance(vxs[nbFaceIdxs[1]]) > tolerance ) sharedFaceVertices = false;
+
+                if (sharedFaceVertices)
+                {
+                    continue;
+                }
+
+                // To avoid doing this calculation for the opposite face
+
+                faultsPrCellAcc->setFaultIdx(gcIdx, face, unNamedFaultIdx);
+                faultsPrCellAcc->setFaultIdx(neighborGlobalCellIdx, StructGridInterface::oppositeFace(face), unNamedFaultIdx);
+
+                // Add as fault face only if the grid index is less than the neighbors
+
+                if (gcIdx < neighborGlobalCellIdx)
+                {
+                    RigFault::FaultFace ff(gcIdx, cvf::StructGridInterface::FaceType(faceIdx), neighborGlobalCellIdx);
+                    unNamedFault->faultFaces().push_back(ff);
+                }
+                else
+                {
+                    CVF_ASSERT(false); // Should never occur. because we flag the opposite face in the faultsPrCellAcc
+                }
+            }
+        }
+    }
+
+    if (unNamedFault->faultFaces().size())
+    {
+        unNamedFault->setName("Unnamed grid faults");
+        m_faults.push_back(unNamedFault);
+    }
+}