diff --git a/ApplicationCode/ReservoirDataModel/CMakeLists_files.cmake b/ApplicationCode/ReservoirDataModel/CMakeLists_files.cmake index 03fc9e1c83..001434a952 100644 --- a/ApplicationCode/ReservoirDataModel/CMakeLists_files.cmake +++ b/ApplicationCode/ReservoirDataModel/CMakeLists_files.cmake @@ -75,6 +75,7 @@ ${CMAKE_CURRENT_LIST_DIR}/RigEclipseCrossPlotDataExtractor.h ${CMAKE_CURRENT_LIST_DIR}/RigEquil.h ${CMAKE_CURRENT_LIST_DIR}/RigWbsParameter.h ${CMAKE_CURRENT_LIST_DIR}/RigEclipseAllenFaultsStatCalc.h +${CMAKE_CURRENT_LIST_DIR}/RigCellFaceGeometryTools.h ) @@ -147,6 +148,7 @@ ${CMAKE_CURRENT_LIST_DIR}/RigEclipseCrossPlotDataExtractor.cpp ${CMAKE_CURRENT_LIST_DIR}/RigEquil.cpp ${CMAKE_CURRENT_LIST_DIR}/RigWbsParameter.cpp ${CMAKE_CURRENT_LIST_DIR}/RigEclipseAllenFaultsStatCalc.cpp +${CMAKE_CURRENT_LIST_DIR}/RigCellFaceGeometryTools.cpp ) list(APPEND CODE_HEADER_FILES diff --git a/ApplicationCode/ReservoirDataModel/RigCellFaceGeometryTools.cpp b/ApplicationCode/ReservoirDataModel/RigCellFaceGeometryTools.cpp new file mode 100644 index 0000000000..72fd60b40b --- /dev/null +++ b/ApplicationCode/ReservoirDataModel/RigCellFaceGeometryTools.cpp @@ -0,0 +1,112 @@ +///////////////////////////////////////////////////////////////////////////////// +// +// Copyright (C) 2020 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 "RigCellFaceGeometryTools.h" + +#include "RigCell.h" +#include "RigMainGrid.h" + +#include "cvfGeometryTools.h" + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +cvf::StructGridInterface::FaceType + RigCellFaceGeometryTools::calculateCellFaceOverlap( const RigCell& c1, + const RigCell& c2, + const RigMainGrid& mainGrid, + std::vector* connectionPolygon, + std::vector* connectionIntersections ) +{ + // Try to find the shared face + + bool isPossibleNeighborInDirection[6] = {true, true, true, true, true, true}; + + if ( c1.hostGrid() == c2.hostGrid() ) + { + char hasNeighbourInAnyDirection = 0; + + size_t i1, j1, k1; + c1.hostGrid()->ijkFromCellIndex( c1.gridLocalCellIndex(), &i1, &j1, &k1 ); + size_t i2, j2, k2; + c2.hostGrid()->ijkFromCellIndex( c2.gridLocalCellIndex(), &i2, &j2, &k2 ); + + isPossibleNeighborInDirection[cvf::StructGridInterface::POS_I] = ( ( i1 + 1 ) == i2 ); + isPossibleNeighborInDirection[cvf::StructGridInterface::NEG_I] = ( ( i2 + 1 ) == i1 ); + isPossibleNeighborInDirection[cvf::StructGridInterface::POS_J] = ( ( j1 + 1 ) == j2 ); + isPossibleNeighborInDirection[cvf::StructGridInterface::NEG_J] = ( ( j2 + 1 ) == j1 ); + isPossibleNeighborInDirection[cvf::StructGridInterface::POS_K] = ( ( k1 + 1 ) == k2 ); + isPossibleNeighborInDirection[cvf::StructGridInterface::NEG_K] = ( ( k2 + 1 ) == k1 ); + + hasNeighbourInAnyDirection = isPossibleNeighborInDirection[cvf::StructGridInterface::POS_I] + + isPossibleNeighborInDirection[cvf::StructGridInterface::NEG_I] + + isPossibleNeighborInDirection[cvf::StructGridInterface::POS_J] + + isPossibleNeighborInDirection[cvf::StructGridInterface::NEG_J] + + isPossibleNeighborInDirection[cvf::StructGridInterface::POS_K] + + isPossibleNeighborInDirection[cvf::StructGridInterface::NEG_K]; + + // If cell 2 is not adjancent with respect to any of the six ijk directions, + // assume that we have no overlapping area. + + if ( !hasNeighbourInAnyDirection ) + { + // Add to search map + // m_cellIdxToFaceToConnectionIdxMap[m_connections[cnIdx].m_c1GlobIdx][cvf::StructGridInterface::NO_FACE].push_back(cnIdx); + // m_cellIdxToFaceToConnectionIdxMap[m_connections[cnIdx].m_c2GlobIdx][cvf::StructGridInterface::NO_FACE].push_back(cnIdx); + + // cvf::Trace::show("NNC: No direct neighbors : C1: " + cvf::String((int)m_connections[cnIdx].m_c1GlobIdx) + + // " C2: " + cvf::String((int)m_connections[cnIdx].m_c2GlobIdx)); + return cvf::StructGridInterface::NO_FACE; + } + } + + for ( unsigned char fIdx = 0; fIdx < 6; ++fIdx ) + { + if ( !isPossibleNeighborInDirection[fIdx] ) + { + continue; + } + + // Calculate connection polygon + + std::vector polygon; + std::vector intersections; + std::array face1; + std::array face2; + c1.faceIndices( ( cvf::StructGridInterface::FaceType )( fIdx ), &face1 ); + c2.faceIndices( cvf::StructGridInterface::oppositeFace( ( cvf::StructGridInterface::FaceType )( fIdx ) ), &face2 ); + + bool foundOverlap = + cvf::GeometryTools::calculateOverlapPolygonOfTwoQuads( &polygon, + &intersections, + (cvf::EdgeIntersectStorage*)nullptr, + cvf::wrapArrayConst( &mainGrid.nodes() ), + face1.data(), + face2.data(), + 1e-6 ); + + if ( foundOverlap ) + { + if ( connectionPolygon ) ( *connectionPolygon ) = polygon; + if ( connectionIntersections ) ( *connectionIntersections ) = intersections; + return ( cvf::StructGridInterface::FaceType )( fIdx ); + } + } + + return cvf::StructGridInterface::NO_FACE; +} diff --git a/ApplicationCode/ReservoirDataModel/RigCellFaceGeometryTools.h b/ApplicationCode/ReservoirDataModel/RigCellFaceGeometryTools.h new file mode 100644 index 0000000000..0c3eac5d56 --- /dev/null +++ b/ApplicationCode/ReservoirDataModel/RigCellFaceGeometryTools.h @@ -0,0 +1,41 @@ +///////////////////////////////////////////////////////////////////////////////// +// +// Copyright (C) 2020 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 "cvfCollection.h" +#include "cvfStructGrid.h" +#include "cvfVector3.h" + +#include + +class RigCell; +class RigMainGrid; + +//================================================================================================== +/// +//================================================================================================== +class RigCellFaceGeometryTools +{ +public: + static cvf::StructGridInterface::FaceType calculateCellFaceOverlap( const RigCell& c1, + const RigCell& c2, + const RigMainGrid& mainGrid, + std::vector* connectionPolygon, + std::vector* connectionIntersections ); +}; diff --git a/ApplicationCode/ReservoirDataModel/RigNNCData.cpp b/ApplicationCode/ReservoirDataModel/RigNNCData.cpp index 3294caa092..0d3d8ef0bd 100644 --- a/ApplicationCode/ReservoirDataModel/RigNNCData.cpp +++ b/ApplicationCode/ReservoirDataModel/RigNNCData.cpp @@ -18,8 +18,11 @@ ///////////////////////////////////////////////////////////////////////////////// #include "RigNNCData.h" + +#include "RigCellFaceGeometryTools.h" #include "RigEclipseResultAddress.h" #include "RigMainGrid.h" + #include "cvfGeometryTools.h" //-------------------------------------------------------------------------------------------------- @@ -43,7 +46,11 @@ void RigNNCData::processConnections( const RigMainGrid& mainGrid ) std::vector connectionIntersections; cvf::StructGridInterface::FaceType connectionFace = cvf::StructGridInterface::NO_FACE; - connectionFace = calculateCellFaceOverlap( c1, c2, mainGrid, &connectionPolygon, &connectionIntersections ); + connectionFace = RigCellFaceGeometryTools::calculateCellFaceOverlap( c1, + c2, + mainGrid, + &connectionPolygon, + &connectionIntersections ); if ( connectionFace != cvf::StructGridInterface::NO_FACE ) { @@ -71,114 +78,6 @@ void RigNNCData::processConnections( const RigMainGrid& mainGrid ) } } -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -cvf::StructGridInterface::FaceType RigNNCData::calculateCellFaceOverlap( const RigCell& c1, - const RigCell& c2, - const RigMainGrid& mainGrid, - std::vector* connectionPolygon, - std::vector* connectionIntersections ) -{ - // Try to find the shared face - - bool isPossibleNeighborInDirection[6] = {true, true, true, true, true, true}; - - if ( c1.hostGrid() == c2.hostGrid() ) - { - char hasNeighbourInAnyDirection = 0; - - size_t i1, j1, k1; - c1.hostGrid()->ijkFromCellIndex( c1.gridLocalCellIndex(), &i1, &j1, &k1 ); - size_t i2, j2, k2; - c2.hostGrid()->ijkFromCellIndex( c2.gridLocalCellIndex(), &i2, &j2, &k2 ); - - isPossibleNeighborInDirection[cvf::StructGridInterface::POS_I] = ( ( i1 + 1 ) == i2 ); - isPossibleNeighborInDirection[cvf::StructGridInterface::NEG_I] = ( ( i2 + 1 ) == i1 ); - isPossibleNeighborInDirection[cvf::StructGridInterface::POS_J] = ( ( j1 + 1 ) == j2 ); - isPossibleNeighborInDirection[cvf::StructGridInterface::NEG_J] = ( ( j2 + 1 ) == j1 ); - isPossibleNeighborInDirection[cvf::StructGridInterface::POS_K] = ( ( k1 + 1 ) == k2 ); - isPossibleNeighborInDirection[cvf::StructGridInterface::NEG_K] = ( ( k2 + 1 ) == k1 ); - - hasNeighbourInAnyDirection = isPossibleNeighborInDirection[cvf::StructGridInterface::POS_I] + - isPossibleNeighborInDirection[cvf::StructGridInterface::NEG_I] + - isPossibleNeighborInDirection[cvf::StructGridInterface::POS_J] + - isPossibleNeighborInDirection[cvf::StructGridInterface::NEG_J] + - isPossibleNeighborInDirection[cvf::StructGridInterface::POS_K] + - isPossibleNeighborInDirection[cvf::StructGridInterface::NEG_K]; - - // If cell 2 is not adjancent with respect to any of the six ijk directions, - // assume that we have no overlapping area. - - if ( !hasNeighbourInAnyDirection ) - { - // Add to search map - // m_cellIdxToFaceToConnectionIdxMap[m_connections[cnIdx].m_c1GlobIdx][cvf::StructGridInterface::NO_FACE].push_back(cnIdx); - // m_cellIdxToFaceToConnectionIdxMap[m_connections[cnIdx].m_c2GlobIdx][cvf::StructGridInterface::NO_FACE].push_back(cnIdx); - - // cvf::Trace::show("NNC: No direct neighbors : C1: " + cvf::String((int)m_connections[cnIdx].m_c1GlobIdx) + - // " C2: " + cvf::String((int)m_connections[cnIdx].m_c2GlobIdx)); - return cvf::StructGridInterface::NO_FACE; - } - } - -#if 0 - // Possibly do some testing to avoid unneccesary overlap calculations - cvf::Vec3d normal; - for ( char fIdx = 0; fIdx < 6; ++fIdx ) - { - if ( isPossibleNeighborInDirection[fIdx] ) - { - cvf::Vec3d fc1 = c1.faceCenter((cvf::StructGridInterface::FaceType)(fIdx)); - cvf::Vec3d fc2 = c2.faceCenter(cvf::StructGridInterface::oppositeFace((cvf::StructGridInterface::FaceType)(fIdx))); - cvf::Vec3d fc1ToFc2 = fc2 - fc1; - normal = c1.faceNormalWithAreaLenght((cvf::StructGridInterface::FaceType)(fIdx)); - normal.normalize(); - // Check that face centers are approx in the face plane - if ( normal.dot(fc1ToFc2) < 0.01*fc1ToFc2.length() ) - { - - } - } - } -#endif - - for ( unsigned char fIdx = 0; fIdx < 6; ++fIdx ) - { - if ( !isPossibleNeighborInDirection[fIdx] ) - { - continue; - } - - // Calculate connection polygon - - std::vector polygon; - std::vector intersections; - std::array face1; - std::array face2; - c1.faceIndices( ( cvf::StructGridInterface::FaceType )( fIdx ), &face1 ); - c2.faceIndices( cvf::StructGridInterface::oppositeFace( ( cvf::StructGridInterface::FaceType )( fIdx ) ), &face2 ); - - bool foundOverlap = - cvf::GeometryTools::calculateOverlapPolygonOfTwoQuads( &polygon, - &intersections, - (cvf::EdgeIntersectStorage*)nullptr, - cvf::wrapArrayConst( &mainGrid.nodes() ), - face1.data(), - face2.data(), - 1e-6 ); - - if ( foundOverlap ) - { - if ( connectionPolygon ) ( *connectionPolygon ) = polygon; - if ( connectionIntersections ) ( *connectionIntersections ) = intersections; - return ( cvf::StructGridInterface::FaceType )( fIdx ); - } - } - - return cvf::StructGridInterface::NO_FACE; -} - //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- diff --git a/ApplicationCode/ReservoirDataModel/RigNNCData.h b/ApplicationCode/ReservoirDataModel/RigNNCData.h index 13059ea2e8..93132e8237 100644 --- a/ApplicationCode/ReservoirDataModel/RigNNCData.h +++ b/ApplicationCode/ReservoirDataModel/RigNNCData.h @@ -69,12 +69,6 @@ public: void processConnections( const RigMainGrid& mainGrid ); - static cvf::StructGridInterface::FaceType calculateCellFaceOverlap( const RigCell& c1, - const RigCell& c2, - const RigMainGrid& mainGrid, - std::vector* connectionPolygon, - std::vector* connectionIntersections ); - void setConnections( std::vector& connections ); const std::vector& connections() const; diff --git a/ApplicationCode/ReservoirDataModel/RigSimulationWellCenterLineCalculator.cpp b/ApplicationCode/ReservoirDataModel/RigSimulationWellCenterLineCalculator.cpp index 26249621b1..94b4c275d5 100644 --- a/ApplicationCode/ReservoirDataModel/RigSimulationWellCenterLineCalculator.cpp +++ b/ApplicationCode/ReservoirDataModel/RigSimulationWellCenterLineCalculator.cpp @@ -20,6 +20,7 @@ #include "RigSimulationWellCenterLineCalculator.h" #include "RigCell.h" +#include "RigCellFaceGeometryTools.h" #include "RigEclipseCaseData.h" #include "RigMainGrid.h" @@ -800,11 +801,12 @@ private: std::vector poygonIndices; std::vector intersections; - auto contactFace = RigNNCData::calculateCellFaceOverlap( c1, - c2, - *( m_eclipseCaseData->mainGrid() ), - &poygonIndices, - &intersections ); + auto contactFace = + RigCellFaceGeometryTools::calculateCellFaceOverlap( c1, + c2, + *( m_eclipseCaseData->mainGrid() ), + &poygonIndices, + &intersections ); if ( contactFace != cvf::StructGridInterface::NO_FACE ) { diff --git a/ApplicationCode/ReservoirDataModel/RigSimulationWellCenterLineCalculator.h b/ApplicationCode/ReservoirDataModel/RigSimulationWellCenterLineCalculator.h index bd1d757b35..b99d8ba2e2 100644 --- a/ApplicationCode/ReservoirDataModel/RigSimulationWellCenterLineCalculator.h +++ b/ApplicationCode/ReservoirDataModel/RigSimulationWellCenterLineCalculator.h @@ -16,6 +16,7 @@ // for more details. // ///////////////////////////////////////////////////////////////////////////////// + #pragma once #include "RigSimWellData.h"