From 7369ce56ba681874ebece44746ef646f63eb5b50 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Bj=C3=B8rn=20Erik=20Jensen?= Date: Mon, 20 Aug 2018 11:47:05 +0200 Subject: [PATCH] #3114 compdat export. Improve performance --- .../RicExportFractureCompletionsImpl.cpp | 17 +++++-- .../RigCellGeometryTools.cpp | 25 +++++----- .../ReservoirDataModel/RigCellGeometryTools.h | 10 ++-- .../RigHexIntersectionTools.cpp | 16 +++---- .../RigHexIntersectionTools.h | 2 +- .../RigHexIntersectionTools-Test.cpp | 47 +++++++++++++++++++ 6 files changed, 88 insertions(+), 29 deletions(-) diff --git a/ApplicationCode/Commands/CompletionExportCommands/RicExportFractureCompletionsImpl.cpp b/ApplicationCode/Commands/CompletionExportCommands/RicExportFractureCompletionsImpl.cpp index 643541af42..78a30a6c2f 100644 --- a/ApplicationCode/Commands/CompletionExportCommands/RicExportFractureCompletionsImpl.cpp +++ b/ApplicationCode/Commands/CompletionExportCommands/RicExportFractureCompletionsImpl.cpp @@ -142,10 +142,14 @@ std::vector // To handle several fractures in the same eclipse cell we need to keep track of the transmissibility // to the well from each fracture intersecting the cell and sum these transmissibilities at the end. // std::map > - std::map> eclCellIdxToTransPrFractureMap; + //std::map> eclCellIdxToTransPrFractureMap; - for (const RimFracture* fracture : fractures) + std::vector> sharedComplForFracture(fractures.size()); + +#pragma omp parallel for + for (int i = 0; i < (int)fractures.size(); i++) { + RimFracture* fracture = fractures[i]; RimFractureTemplate* fracTemplate = fracture->fractureTemplate(); if (!fracTemplate) continue; @@ -320,7 +324,7 @@ std::vector double trans = transCondenser.condensedTransmissibility( externalCell, {true, RigTransmissibilityCondenser::CellAddress::WELL, 1}); - eclCellIdxToTransPrFractureMap[externalCell.m_globalCellIdx][fracture] = trans; + //eclCellIdxToTransPrFractureMap[externalCell.m_globalCellIdx][fracture] = trans; RigCompletionData compDat(wellPathName, RigCompletionDataGridCell(externalCell.m_globalCellIdx, caseToApply->mainGrid()), @@ -459,7 +463,8 @@ std::vector } std::copy( - allCompletionsForOneFracture.begin(), allCompletionsForOneFracture.end(), std::back_inserter(fractureCompletions)); + allCompletionsForOneFracture.begin(), allCompletionsForOneFracture.end(), std::back_inserter(sharedComplForFracture[i])); + if (outputStreamForIntermediateResultsText) { @@ -483,5 +488,9 @@ std::vector } } + for (const auto& completions : sharedComplForFracture) + { + std::copy(completions.begin(), completions.end(), std::back_inserter(fractureCompletions)); + } return fractureCompletions; } diff --git a/ApplicationCode/ReservoirDataModel/RigCellGeometryTools.cpp b/ApplicationCode/ReservoirDataModel/RigCellGeometryTools.cpp index f1cc1478e1..506bae0ac9 100644 --- a/ApplicationCode/ReservoirDataModel/RigCellGeometryTools.cpp +++ b/ApplicationCode/ReservoirDataModel/RigCellGeometryTools.cpp @@ -151,7 +151,7 @@ void RigCellGeometryTools::findCellLocalXYZ(const std::array& hex //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- -double RigCellGeometryTools::polygonLengthInLocalXdirWeightedByArea(std::vector polygonToCalcLengthOf) +double RigCellGeometryTools::polygonLengthInLocalXdirWeightedByArea(const std::vector& polygonToCalcLengthOf) { //Find bounding box cvf::BoundingBox polygonBBox; @@ -250,19 +250,20 @@ cvf::Vec3d fromClipperPoint(const ClipperLib::IntPoint& clipPoint) //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- -std::vector > RigCellGeometryTools::intersectPolygons(std::vector polygon1, std::vector polygon2) +std::vector > RigCellGeometryTools::intersectPolygons(const std::vector& polygon1, + const std::vector& polygon2) { std::vector > clippedPolygons; // Convert to int for clipper library and store as clipper "path" ClipperLib::Path polygon1path; - for (cvf::Vec3d& v : polygon1) + for (const cvf::Vec3d& v : polygon1) { polygon1path.push_back(toClipperPoint(v)); } ClipperLib::Path polygon2path; - for (cvf::Vec3d& v : polygon2) + for (const cvf::Vec3d& v : polygon2) { polygon2path.push_back(toClipperPoint(v)); } @@ -458,28 +459,30 @@ std::vector > RigCellGeometryTools::clipPolylineByPolygo //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- -std::pair RigCellGeometryTools::getLineThroughBoundingBox(cvf::Vec3d lineDirection, cvf::BoundingBox polygonBBox, cvf::Vec3d pointOnLine) +std::pair RigCellGeometryTools::getLineThroughBoundingBox(const cvf::Vec3d& lineDirection, + const cvf::BoundingBox& polygonBBox, + const cvf::Vec3d& pointOnLine) { cvf::Vec3d bboxCorners[8]; polygonBBox.cornerVertices(bboxCorners); cvf::Vec3d startPoint = pointOnLine; cvf::Vec3d endPoint = pointOnLine; - + cvf::Vec3d lineDir = lineDirection; //To avoid doing many iterations in loops below linedirection should be quite large. - lineDirection.normalize(); - lineDirection = lineDirection * polygonBBox.extent().length() / 5; + lineDir.normalize(); + lineDir = lineDir * polygonBBox.extent().length() / 5; //Extend line in positive direction while (polygonBBox.contains(startPoint)) { - startPoint = startPoint + lineDirection; + startPoint = startPoint + lineDir; } //Extend line in negative direction while (polygonBBox.contains(endPoint)) { - endPoint = endPoint - lineDirection; + endPoint = endPoint - lineDir; } std::pair line; @@ -508,7 +511,7 @@ double RigCellGeometryTools::getLengthOfPolygonAlongLine(const std::pair RigCellGeometryTools::unionOfPolygons(std::vector> polygons) +std::vector RigCellGeometryTools::unionOfPolygons(const std::vector>& polygons) { // Convert to int for clipper library and store as clipper "path" std::vector polygonPaths; diff --git a/ApplicationCode/ReservoirDataModel/RigCellGeometryTools.h b/ApplicationCode/ReservoirDataModel/RigCellGeometryTools.h index ded82413a0..49b4488e16 100644 --- a/ApplicationCode/ReservoirDataModel/RigCellGeometryTools.h +++ b/ApplicationCode/ReservoirDataModel/RigCellGeometryTools.h @@ -38,10 +38,10 @@ public: static void findCellLocalXYZ(const std::array& hexCorners, cvf::Vec3d &localXdirection, cvf::Vec3d &localYdirection, cvf::Vec3d &localZdirection); - static double polygonLengthInLocalXdirWeightedByArea(std::vector polygon2d); + static double polygonLengthInLocalXdirWeightedByArea(const std::vector& polygon2d); - static std::vector> intersectPolygons(std::vector polygon1, - std::vector polygon2); + static std::vector> intersectPolygons(const std::vector& polygon1, + const std::vector& polygon2); static std::vector> subtractPolygons(const std::vector& sourcePolygon, const std::vector>& polygonToSubtract); @@ -51,11 +51,11 @@ public: const std::vector& polygon, ZInterpolationType interpolType = USE_ZERO); - static std::pair getLineThroughBoundingBox(cvf::Vec3d lineDirection, cvf::BoundingBox polygonBBox, cvf::Vec3d pointOnLine); + static std::pair getLineThroughBoundingBox(const cvf::Vec3d& lineDirection, const cvf::BoundingBox& polygonBBox, const cvf::Vec3d& pointOnLine); static double getLengthOfPolygonAlongLine(const std::pair& line, const std::vector& polygon); - static std::vector unionOfPolygons(std::vector> polygons); + static std::vector unionOfPolygons(const std::vector>& polygons); private: static std::vector ajustPolygonToAvoidIntersectionsAtVertex(const std::vector& polyLine, diff --git a/ApplicationCode/ReservoirDataModel/RigHexIntersectionTools.cpp b/ApplicationCode/ReservoirDataModel/RigHexIntersectionTools.cpp index 3d51738c39..db40714f13 100644 --- a/ApplicationCode/ReservoirDataModel/RigHexIntersectionTools.cpp +++ b/ApplicationCode/ReservoirDataModel/RigHexIntersectionTools.cpp @@ -114,12 +114,16 @@ bool RigHexIntersectionTools::isPointInCell(const cvf::Vec3d point, //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- -bool RigHexIntersectionTools::planeHexCellIntersection(cvf::Vec3d* hexCorners, cvf::Plane fracturePlane, std::list >& intersectionLineSegments) +bool RigHexIntersectionTools::planeHexCellIntersection(cvf::Vec3d* hexCorners, const cvf::Plane& fracturePlane, std::list >& intersectionLineSegments) { - bool isCellIntersected = false; + bool isCellIntersected = false; + cvf::ubyte faceVertexIndices[4]; + caf::HexGridIntersectionTools::ClipVx triangleIntersectionPoint1; + caf::HexGridIntersectionTools::ClipVx triangleIntersectionPoint2; + bool isMostVxesOnPositiveSideOfP1 = false; + for (int face = 0; face < 6; ++face) { - cvf::ubyte faceVertexIndices[4]; cvf::StructGridInterface::cellFaceVertexIndices(static_cast(face), faceVertexIndices); cvf::Vec3d faceCenter = cvf::GeometryTools::computeFaceCenter(hexCorners[faceVertexIndices[0]], hexCorners[faceVertexIndices[1]], hexCorners[faceVertexIndices[2]], hexCorners[faceVertexIndices[3]]); @@ -127,10 +131,6 @@ bool RigHexIntersectionTools::planeHexCellIntersection(cvf::Vec3d* hexCorners, c for (int i = 0; i < 4; i++) { int next = i < 3 ? i + 1 : 0; - caf::HexGridIntersectionTools::ClipVx triangleIntersectionPoint1; - caf::HexGridIntersectionTools::ClipVx triangleIntersectionPoint2; - - bool isMostVxesOnPositiveSideOfP1 = false; bool isIntersectingPlane = caf::HexGridIntersectionTools::planeTriangleIntersection(fracturePlane, hexCorners[faceVertexIndices[i]], 0, @@ -141,7 +141,7 @@ bool RigHexIntersectionTools::planeHexCellIntersection(cvf::Vec3d* hexCorners, c if (isIntersectingPlane) { isCellIntersected = true; - intersectionLineSegments.push_back({ triangleIntersectionPoint1.vx, triangleIntersectionPoint2.vx }); + intersectionLineSegments.emplace_back( triangleIntersectionPoint1.vx, triangleIntersectionPoint2.vx ); } } } return isCellIntersected; diff --git a/ApplicationCode/ReservoirDataModel/RigHexIntersectionTools.h b/ApplicationCode/ReservoirDataModel/RigHexIntersectionTools.h index e40bd6507c..5265d43989 100644 --- a/ApplicationCode/ReservoirDataModel/RigHexIntersectionTools.h +++ b/ApplicationCode/ReservoirDataModel/RigHexIntersectionTools.h @@ -64,7 +64,7 @@ struct RigHexIntersectionTools static bool isPointInCell(const cvf::Vec3d point, const cvf::Vec3d hexCorners[8]); static bool planeHexCellIntersection(cvf::Vec3d* hexCorners, - cvf::Plane fracturePlane, + const cvf::Plane& fracturePlane, std::list >& intersectionLineSegments); static bool planeHexIntersectionPolygons(std::array hexCorners, diff --git a/ApplicationCode/UnitTests/RigHexIntersectionTools-Test.cpp b/ApplicationCode/UnitTests/RigHexIntersectionTools-Test.cpp index 3a0175df77..426af26bc4 100644 --- a/ApplicationCode/UnitTests/RigHexIntersectionTools-Test.cpp +++ b/ApplicationCode/UnitTests/RigHexIntersectionTools-Test.cpp @@ -19,6 +19,8 @@ #include "gtest/gtest.h" #include "RigHexIntersectionTools.h" +#include "QDateTime" +#include "QDebug" //-------------------------------------------------------------------------------------------------- /// @@ -48,3 +50,48 @@ TEST(RigHexIntersectionTools, planeHexCellIntersectionTest) EXPECT_FALSE(isCellIntersected); } + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +TEST(RigHexIntersectionTools, DISABLED_planeHexCellIntersectionPerformanceTest) +{ + 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[3] = cvf::Vec3d(0, 0, 1); + hexCorners[4] = cvf::Vec3d(0, 1, 1); + hexCorners[5] = cvf::Vec3d(1, 1, 0); + hexCorners[6] = cvf::Vec3d(1, 0, 1); + hexCorners[7] = cvf::Vec3d(1, 1, 1); + + bool isCellIntersected = false; + cvf::Plane fracturePlaneIntersecting; + cvf::Plane fracturePlaneNotIntersecting; + + fracturePlaneNotIntersecting.setFromPointAndNormal(cvf::Vec3d(1.5, 1.5, 1.5), cvf::Vec3d(1, 0, 0)); + fracturePlaneIntersecting.setFromPointAndNormal(cvf::Vec3d(0.5, 0.5, 0.5), cvf::Vec3d(1, 0, 0)); + + QTime timeTotal; + timeTotal.start(); + + for (int run = 0; run < 5; run++) + { + QTime timeLocal; + timeLocal.start(); + + for (int i = 0; i < 2000000; i++) + { + std::list > intersectionLineSegments; + + isCellIntersected = RigHexIntersectionTools::planeHexCellIntersection(hexCorners, fracturePlaneIntersecting, intersectionLineSegments); + isCellIntersected = RigHexIntersectionTools::planeHexCellIntersection(hexCorners, fracturePlaneNotIntersecting, intersectionLineSegments); + } + + qDebug() << "Time rim elapsed: " << timeLocal.elapsed(); + } + + qDebug() << "Time total elapsed: " << timeTotal.elapsed(); +} +