diff --git a/ApplicationCode/ModelVisualization/RivNNCGeometryGenerator.cpp b/ApplicationCode/ModelVisualization/RivNNCGeometryGenerator.cpp index 9629b3c528..003c7e5177 100644 --- a/ApplicationCode/ModelVisualization/RivNNCGeometryGenerator.cpp +++ b/ApplicationCode/ModelVisualization/RivNNCGeometryGenerator.cpp @@ -74,7 +74,7 @@ void RivNNCGeometryGenerator::computeArrays() std::vector vertices; std::vector triangleToNNC; - const cvf::Vec3d offset = m_offset; + const cvf::Vec3f offset( m_offset ); long long numConnections = static_cast( m_nncIndexes.isNull() ? m_nncData->connections().size() : m_nncIndexes->size() ); @@ -122,14 +122,14 @@ void RivNNCGeometryGenerator::computeArrays() if ( isVisible ) { - cvf::Vec3f vx1 = cvf::Vec3f( conn.m_polygon[0] - offset ); + cvf::Vec3f vx1 = conn.m_polygon[0] - offset; cvf::Vec3f vx2; - cvf::Vec3f vx3 = cvf::Vec3f( conn.m_polygon[1] - offset ); + cvf::Vec3f vx3 = conn.m_polygon[1] - offset; for ( size_t vxIdx = 2; vxIdx < conn.m_polygon.size(); ++vxIdx ) { vx2 = vx3; - vx3 = cvf::Vec3f( conn.m_polygon[vxIdx] - offset ); + vx3 = conn.m_polygon[vxIdx] - offset; #pragma omp critical( critical_section_RivNNCGeometryGenerator_computeArrays ) { vertices.push_back( vx1 ); diff --git a/ApplicationCode/ReservoirDataModel/RigCaseCellResultsData.cpp b/ApplicationCode/ReservoirDataModel/RigCaseCellResultsData.cpp index c5ac8a528c..e88004b5a4 100644 --- a/ApplicationCode/ReservoirDataModel/RigCaseCellResultsData.cpp +++ b/ApplicationCode/ReservoirDataModel/RigCaseCellResultsData.cpp @@ -2366,11 +2366,11 @@ void RigCaseCellResultsData::computeNncCombRiTrans() // Connection geometry - cvf::Vec3d faceAreaVec = cvf::Vec3d::ZERO; - cvf::Vec3d faceCenter = cvf::Vec3d::ZERO; + cvf::Vec3f faceAreaVec = cvf::Vec3f::ZERO; + cvf::Vec3f faceCenter = cvf::Vec3f::ZERO; // Polygon center - const std::vector& realPolygon = nncConnections[connIdx].m_polygon; + const std::vector& realPolygon = nncConnections[connIdx].m_polygon; for ( size_t pIdx = 0; pIdx < realPolygon.size(); ++pIdx ) { faceCenter += realPolygon[pIdx]; @@ -2400,7 +2400,7 @@ void RigCaseCellResultsData::computeNncCombRiTrans() ntg = ( *ntgResults )[ntgResIdx]; } - halfCellTrans = halfCellTransmissibility( perm, ntg, centerToFace, faceAreaVec ); + halfCellTrans = halfCellTransmissibility( perm, ntg, centerToFace, cvf::Vec3d( faceAreaVec ) ); } // Neighbor cell half cell transm @@ -2417,7 +2417,7 @@ void RigCaseCellResultsData::computeNncCombRiTrans() ntg = ( *ntgResults )[ntgResIdx]; } - neighborHalfCellTrans = halfCellTransmissibility( perm, ntg, centerToFace, -faceAreaVec ); + neighborHalfCellTrans = halfCellTransmissibility( perm, ntg, centerToFace, -cvf::Vec3d( faceAreaVec ) ); } double newtranTemp = newtran( cdarchy, 1.0, halfCellTrans, neighborHalfCellTrans ); @@ -2686,8 +2686,8 @@ void RigCaseCellResultsData::computeNncCombRiTRANSbyArea() for ( size_t nncConIdx = 0; nncConIdx < riAreaNormTransResults.size(); ++nncConIdx ) { - const std::vector& realPolygon = connections[nncConIdx].m_polygon; - cvf::Vec3d faceAreaVec = cvf::GeometryTools::polygonAreaNormal3D( realPolygon ); + const std::vector& realPolygon = connections[nncConIdx].m_polygon; + cvf::Vec3f faceAreaVec = cvf::GeometryTools::polygonAreaNormal3D( realPolygon ); double areaOfOverlap = faceAreaVec.length(); riAreaNormTransResults[nncConIdx] = ( *transResults )[nncConIdx] / areaOfOverlap; diff --git a/ApplicationCode/ReservoirDataModel/RigCellFaceGeometryTools.cpp b/ApplicationCode/ReservoirDataModel/RigCellFaceGeometryTools.cpp index d74162cd24..d432efdebd 100644 --- a/ApplicationCode/ReservoirDataModel/RigCellFaceGeometryTools.cpp +++ b/ApplicationCode/ReservoirDataModel/RigCellFaceGeometryTools.cpp @@ -119,9 +119,9 @@ cvf::StructGridInterface::FaceType return cvf::StructGridInterface::NO_FACE; } -void assignThreadConnections( std::set>& existingPairs, - RigConnectionContainer& allConnections, - RigConnectionContainer& threadConnections ) +void assignThreadConnections( std::set>& existingPairs, + RigConnectionContainer& allConnections, + RigConnectionContainer& threadConnections ) { for ( size_t i = 0; i < threadConnections.size(); ++i ) { @@ -150,7 +150,7 @@ RigConnectionContainer RigCellFaceGeometryTools::computeOtherNncs( const RigMain // by Eclipse. Use faults as basis for subset of cells to find NNC connection for. The imported connections from // Eclipse are located at the beginning of the connections vector. - std::set> nativeCellPairs; + std::set> nativeCellPairs; for ( size_t i = 0; i < nativeConnections.size(); ++i ) { @@ -168,8 +168,8 @@ RigConnectionContainer RigCellFaceGeometryTools::computeOtherNncs( const RigMain const cvf::Collection& faults = mainGrid->faults(); - std::set> existingPairs; - RigConnectionContainer otherConnections; + std::set> existingPairs; + RigConnectionContainer otherConnections; for ( int faultIdx = 0; faultIdx < (int)faults.size(); faultIdx++ ) { @@ -216,9 +216,9 @@ RigConnectionContainer RigCellFaceGeometryTools::computeOtherNncs( const RigMain } RigConnectionContainer - RigCellFaceGeometryTools::extractConnectionsForFace( const RigFault::FaultFace& face, - const RigMainGrid* mainGrid, - const std::set>& nativeCellPairs ) + RigCellFaceGeometryTools::extractConnectionsForFace( const RigFault::FaultFace& face, + const RigMainGrid* mainGrid, + const std::set>& nativeCellPairs ) { RigConnectionContainer faceConnections; size_t sourceReservoirCellIndex = face.m_nativeReservoirCellIndex; @@ -317,7 +317,7 @@ RigConnectionContainer } } - std::pair candidate( sourceReservoirCellIndex, candidateCellIndex ); + std::pair candidate( sourceReservoirCellIndex, candidateCellIndex ); if ( nativeCellPairs.count( candidate ) > 0 ) { @@ -357,18 +357,18 @@ RigConnectionContainer //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- -std::vector RigCellFaceGeometryTools::extractPolygon( const std::vector& nativeNodes, +std::vector RigCellFaceGeometryTools::extractPolygon( const std::vector& nativeNodes, const std::vector& connectionPolygon, const std::vector& connectionIntersections ) { - std::vector allPolygonNodes; + std::vector allPolygonNodes; for ( size_t polygonIndex : connectionPolygon ) { if ( polygonIndex < nativeNodes.size() ) - allPolygonNodes.push_back( nativeNodes[polygonIndex] ); + allPolygonNodes.push_back( cvf::Vec3f( nativeNodes[polygonIndex] ) ); else - allPolygonNodes.push_back( connectionIntersections[polygonIndex - nativeNodes.size()] ); + allPolygonNodes.push_back( cvf::Vec3f( connectionIntersections[polygonIndex - nativeNodes.size()] ) ); } return allPolygonNodes; diff --git a/ApplicationCode/ReservoirDataModel/RigCellFaceGeometryTools.h b/ApplicationCode/ReservoirDataModel/RigCellFaceGeometryTools.h index 275961d2a4..d506e13969 100644 --- a/ApplicationCode/ReservoirDataModel/RigCellFaceGeometryTools.h +++ b/ApplicationCode/ReservoirDataModel/RigCellFaceGeometryTools.h @@ -51,10 +51,11 @@ public: const RigActiveCellInfo* activeCellInfo, bool includeInactiveCells ); - static RigConnectionContainer extractConnectionsForFace( const RigFault::FaultFace& face, - const RigMainGrid* mainGrid, - const std::set>& nativeCellPairs ); - static std::vector extractPolygon( const std::vector& nativeNodes, + static RigConnectionContainer + extractConnectionsForFace( const RigFault::FaultFace& face, + const RigMainGrid* mainGrid, + const std::set>& nativeCellPairs ); + static std::vector extractPolygon( const std::vector& nativeNodes, const std::vector& connectionPolygon, const std::vector& connectionIntersections ); }; diff --git a/ApplicationCode/ReservoirDataModel/RigNncConnection.cpp b/ApplicationCode/ReservoirDataModel/RigNncConnection.cpp index 48e8d31516..1c7b789e26 100644 --- a/ApplicationCode/ReservoirDataModel/RigNncConnection.cpp +++ b/ApplicationCode/ReservoirDataModel/RigNncConnection.cpp @@ -24,8 +24,8 @@ /// //-------------------------------------------------------------------------------------------------- RigConnection::RigConnection() - : m_c1GlobIdx( cvf::UNDEFINED_SIZE_T ) - , m_c2GlobIdx( cvf::UNDEFINED_SIZE_T ) + : m_c1GlobIdx( cvf::UNDEFINED_UINT ) + , m_c2GlobIdx( cvf::UNDEFINED_UINT ) , m_c1Face( cvf::StructGridInterface::NO_FACE ) { } @@ -33,10 +33,10 @@ RigConnection::RigConnection() //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- -RigConnection::RigConnection( size_t c1GlobIdx, - size_t c2GlobIdx, +RigConnection::RigConnection( unsigned c1GlobIdx, + unsigned c2GlobIdx, cvf::StructGridInterface::FaceType c1Face, - const std::vector& polygon ) + const std::vector& polygon ) : m_c1GlobIdx( c1GlobIdx ) , m_c2GlobIdx( c2GlobIdx ) , m_c1Face( c1Face ) @@ -94,13 +94,16 @@ bool RigConnection::operator<( const RigConnection& other ) const RigConnection RigConnectionContainer::operator[]( size_t i ) const { const auto& globIndices = m_globalIndices[i]; - return RigConnection( globIndices.first, globIndices.second, m_faces[i], m_polygons[i] ); + return RigConnection( globIndices.first, + globIndices.second, + static_cast( m_faces[i] ), + m_polygons[i] ); } //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- -std::pair& RigConnectionContainer::indexPair( size_t i ) +std::pair& RigConnectionContainer::indexPair( size_t i ) { return m_globalIndices[i]; } @@ -108,7 +111,7 @@ std::pair& RigConnectionContainer::indexPair( size_t i ) //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- -cvf::StructGridInterface::FaceType& RigConnectionContainer::face( size_t i ) +unsigned char& RigConnectionContainer::face( size_t i ) { return m_faces[i]; } @@ -116,7 +119,7 @@ cvf::StructGridInterface::FaceType& RigConnectionContainer::face( size_t i ) //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- -std::vector& RigConnectionContainer::polygon( size_t i ) +std::vector& RigConnectionContainer::polygon( size_t i ) { return m_polygons[i]; } diff --git a/ApplicationCode/ReservoirDataModel/RigNncConnection.h b/ApplicationCode/ReservoirDataModel/RigNncConnection.h index 77addb2b79..8cbccb5c28 100644 --- a/ApplicationCode/ReservoirDataModel/RigNncConnection.h +++ b/ApplicationCode/ReservoirDataModel/RigNncConnection.h @@ -31,20 +31,21 @@ class RigConnection { public: RigConnection(); - RigConnection( size_t c1GlobIdx, - size_t c2GlobIdx, + RigConnection( unsigned c1GlobIdx, + unsigned c2GlobIdx, cvf::StructGridInterface::FaceType c1Face, - const std::vector& polygon ); + const std::vector& polygon ); + RigConnection( const RigConnection& rhs ); RigConnection& operator=( RigConnection& rhs ); bool operator<( const RigConnection& other ) const; bool hasCommonArea() const; - size_t m_c1GlobIdx; - size_t m_c2GlobIdx; + unsigned m_c1GlobIdx; + unsigned m_c2GlobIdx; cvf::StructGridInterface::FaceType m_c1Face; - std::vector m_polygon; + std::vector m_polygon; }; class RigConnectionContainer @@ -52,9 +53,9 @@ class RigConnectionContainer public: RigConnection operator[]( size_t i ) const; - std::pair& indexPair( size_t i ); - cvf::StructGridInterface::FaceType& face( size_t i ); - std::vector& polygon( size_t i ); + std::pair& indexPair( size_t i ); + unsigned char& face( size_t i ); + std::vector& polygon( size_t i ); void push_back( const RigConnection& connection ); void insert( const RigConnectionContainer& other ); @@ -63,7 +64,7 @@ public: bool empty() const; private: - std::vector> m_globalIndices; - std::vector m_faces; - std::vector> m_polygons; + std::vector> m_globalIndices; + std::vector m_faces; + std::vector> m_polygons; }; diff --git a/ApplicationCode/ReservoirDataModel/RigReservoirBuilderMock.cpp b/ApplicationCode/ReservoirDataModel/RigReservoirBuilderMock.cpp index 2c116c3904..396b623ef1 100644 --- a/ApplicationCode/ReservoirDataModel/RigReservoirBuilderMock.cpp +++ b/ApplicationCode/ReservoirDataModel/RigReservoirBuilderMock.cpp @@ -592,8 +592,8 @@ void RigReservoirBuilderMock::addNnc( RigMainGrid* grid, size_t c2GlobalIndex = grid->cellIndexFromIJK( i2, j2, k2 ); RigConnection conn; - conn.m_c1GlobIdx = c1GlobalIndex; - conn.m_c2GlobIdx = c2GlobalIndex; + conn.m_c1GlobIdx = static_cast( c1GlobalIndex ); + conn.m_c2GlobIdx = static_cast( c2GlobalIndex ); nncConnections.push_back( conn ); } diff --git a/ApplicationCode/ReservoirDataModel/cvfGeometryTools.cpp b/ApplicationCode/ReservoirDataModel/cvfGeometryTools.cpp index 59d3ccb0cb..0c9a0d47a2 100644 --- a/ApplicationCode/ReservoirDataModel/cvfGeometryTools.cpp +++ b/ApplicationCode/ReservoirDataModel/cvfGeometryTools.cpp @@ -786,6 +786,59 @@ cvf::Vec3d GeometryTools::polygonAreaNormal3D( const std::vector& po } } +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +cvf::Vec3f GeometryTools::polygonAreaNormal3D( const std::vector& polygon ) +{ + size_t pSize = polygon.size(); + switch ( pSize ) + { + case 0: + case 1: + case 2: + { + return cvf::Vec3f::ZERO; + } + break; + case 3: + { + return 0.5f * ( ( polygon[1] - polygon[0] ) ^ ( polygon[2] - polygon[0] ) ); + } + break; + case 4: + { + // Cross product of diagonal = 2*A + return 0.5f * ( ( polygon[2] - polygon[0] ) ^ ( polygon[3] - polygon[1] ) ); + } + break; + default: + { + /// JJS: + // This is possibly not the fastest approach with large polygons, where a scaled projections approach would + // be better, but I suspect this (simpler) approach is faster for small polygons, as long as we do not have + // the polygon normal up front. + // + cvf::Vec3f areaNormal( cvf::Vec3f::ZERO ); + size_t h = ( pSize - 1 ) / 2; + size_t k = ( pSize % 2 ) ? 0 : pSize - 1; + + // First quads + for ( size_t i = 1; i < h; ++i ) + { + areaNormal += ( ( polygon[2 * i] - polygon[0] ) ^ ( polygon[2 * i + 1] - polygon[2 * i - 1] ) ); + } + + // Last triangle or quad + areaNormal += ( ( polygon[2 * h] - polygon[0] ) ^ ( polygon[k] - polygon[2 * h - 1] ) ); + + areaNormal *= 0.5; + + return areaNormal; + } + } +} + //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- diff --git a/ApplicationCode/ReservoirDataModel/cvfGeometryTools.h b/ApplicationCode/ReservoirDataModel/cvfGeometryTools.h index e68a990c3f..8867a32f26 100644 --- a/ApplicationCode/ReservoirDataModel/cvfGeometryTools.h +++ b/ApplicationCode/ReservoirDataModel/cvfGeometryTools.h @@ -76,6 +76,7 @@ public: static double signedAreaPlanarPolygon( const cvf::Vec3d& planeNormal, const std::vector& polygon ); static cvf::Vec3d polygonAreaNormal3D( const std::vector& polygon ); + static cvf::Vec3f polygonAreaNormal3D( const std::vector& polygon ); enum IntersectionStatus {