diff --git a/ApplicationCode/ReservoirDataModel/RigFractureTransCalc.cpp b/ApplicationCode/ReservoirDataModel/RigFractureTransCalc.cpp index 8961d0918e..0f5ab1a76d 100644 --- a/ApplicationCode/ReservoirDataModel/RigFractureTransCalc.cpp +++ b/ApplicationCode/ReservoirDataModel/RigFractureTransCalc.cpp @@ -311,100 +311,6 @@ bool RigFractureTransCalc::planeCellIntersectionPolygons(size_t cellindex, std:: return isCellIntersected; } - - -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -void RigFractureTransCalc::computeUpscaledPropertyFromStimPlanForEclipseCell(double &upscaledAritmStimPlanValue, double &upscaledHarmStimPlanValue, QString resultName, QString resultUnit, size_t timeStepIndex, caf::AppEnum< RimDefines::UnitSystem > unitSystem, size_t cellIndex) -{ - //TODO: A lot of common code with function for calculating transmissibility... - - RimStimPlanFractureTemplate* fracTemplateStimPlan; - if (dynamic_cast(m_fracture->attachedFractureDefinition())) - { - fracTemplateStimPlan = dynamic_cast(m_fracture->attachedFractureDefinition()); - } - else return; - - std::vector > stimPlanCellsAsPolygons; - std::vector stimPlanParameterValues; - fracTemplateStimPlan->getStimPlanDataAsPolygonsAndValues(stimPlanCellsAsPolygons, stimPlanParameterValues, resultName, resultUnit, timeStepIndex); - - //TODO: A lot of common code with function above... Can be cleaned up...? - std::vector fracCells = m_fracture->getPotentiallyFracturedCells(); - - - RigEclipseCaseData* eclipseCaseData = m_case->eclipseCaseData(); - - RifReaderInterface::PorosityModelResultType porosityModel = RifReaderInterface::MATRIX_RESULTS; - RimReservoirCellResultsStorage* gridCellResults = m_case->results(porosityModel); - RigActiveCellInfo* activeCellInfo = eclipseCaseData->activeCellInfo(porosityModel); - - - bool cellIsActive = activeCellInfo->isActive(cellIndex); - - cvf::Vec3d localX; - cvf::Vec3d localY; - cvf::Vec3d localZ; - std::vector > planeCellPolygons; - bool isPlanIntersected = planeCellIntersectionPolygons(cellIndex, planeCellPolygons, localX, localY, localZ); - if (!isPlanIntersected || planeCellPolygons.size() == 0) return; - - //Transform planCell polygon(s) and averageZdirection to x/y coordinate system (where fracturePolygon/stimPlan mesh already is located) - cvf::Mat4f invertedTransMatrix = m_fracture->transformMatrix().getInverted(); - for (std::vector & planeCellPolygon : planeCellPolygons) - { - for (cvf::Vec3d& v : planeCellPolygon) - { - v.transformPoint(static_cast(invertedTransMatrix)); - } - } - - cvf::Vec3d localZinFracPlane; - localZinFracPlane = localZ; - localZinFracPlane.transformVector(static_cast(invertedTransMatrix)); - cvf::Vec3d directionOfLength = cvf::Vec3d::ZERO; - directionOfLength.cross(localZinFracPlane, cvf::Vec3d(0, 0, 1)); - directionOfLength.normalize(); - - std::vector fracPolygon = m_fracture->attachedFractureDefinition()->fracturePolygon(unitSystem); - std::vector > polygonsDescribingFractureInCell; - - double area; - std::vector areaOfFractureParts; - std::vector valuesForFractureParts; - - for (std::vector planeCellPolygon : planeCellPolygons) - { - - for (int i = 0; i < stimPlanParameterValues.size(); i++) - { - double stimPlanParameterValue = stimPlanParameterValues[i]; - if (stimPlanParameterValue != 0) - { - std::vector stimPlanCell = stimPlanCellsAsPolygons[i]; - std::vector >clippedStimPlanPolygons = RigCellGeometryTools::clipPolygons(stimPlanCell, planeCellPolygon); - if (clippedStimPlanPolygons.size() > 0) - { - for (auto clippedStimPlanPolygon : clippedStimPlanPolygons) - { - area = cvf::GeometryTools::polygonAreaNormal3D(clippedStimPlanPolygon).length(); - areaOfFractureParts.push_back(area); - valuesForFractureParts.push_back(stimPlanParameterValue); - } - } - } - } - } - - if (areaOfFractureParts.size() > 0) - { - upscaledAritmStimPlanValue = areaWeightedArithmeticAverage(areaOfFractureParts, valuesForFractureParts); - upscaledHarmStimPlanValue = areaWeightedHarmonicAverage(areaOfFractureParts, valuesForFractureParts); - } -} - //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- @@ -443,11 +349,8 @@ std::pair RigFractureTransCalc::flowAcrossLayersUpscaling(QStrin } } - //Vector for vertical - på tvers av lag... Gitt av orientering av stimPlan-grid... - cvf::Vec3d directionAcrossLayers; cvf::Vec3d directionAlongLayers; - //TODO: Get these vectors properly... directionAcrossLayers = cvf::Vec3d(0, -1, 0); directionAlongLayers = cvf::Vec3d(1, 0, 0); @@ -461,8 +364,6 @@ std::pair RigFractureTransCalc::flowAcrossLayersUpscaling(QStrin std::vector upscaledConductivitiesHA; std::vector upscaledConductivitiesAH; - - //Harmonic weighted mean for (std::vector planeCellPolygon : planeCellPolygons) { @@ -621,50 +522,6 @@ double RigFractureTransCalc::computeAHupscale(RimStimPlanFractureTemplate* fracT } -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -double RigFractureTransCalc::areaWeightedHarmonicAverage(std::vector areaOfFractureParts, std::vector valuesForFractureParts) -{ - //TODO: Unit test? - double fractureCellArea = 0.0; - for (double area : areaOfFractureParts) fractureCellArea += area; - - if (areaOfFractureParts.size() != valuesForFractureParts.size()) return cvf::UNDEFINED_DOUBLE; - - double fractureCellAreaDivvalue = 0.0; - for (int i = 0; i < valuesForFractureParts.size(); i++) - { - fractureCellAreaDivvalue += areaOfFractureParts[i] / valuesForFractureParts[i]; - } - - double upscaledValueHarmonic = fractureCellArea / fractureCellAreaDivvalue; - return upscaledValueHarmonic; -} - -//-------------------------------------------------------------------------------------------------- -/// -//-------------------------------------------------------------------------------------------------- -double RigFractureTransCalc::areaWeightedArithmeticAverage(std::vector areaOfFractureParts, std::vector valuesForFractureParts) -{ - //TODO: Unit test? - double fractureCellArea = 0.0; - for (double area : areaOfFractureParts) fractureCellArea += area; - - if (areaOfFractureParts.size() != valuesForFractureParts.size()) return cvf::UNDEFINED_DOUBLE; - - double fractureCellAreaXvalue = 0.0; - for (int i = 0; i < valuesForFractureParts.size(); i++) - { - fractureCellAreaXvalue += areaOfFractureParts[i] * valuesForFractureParts[i]; - } - - double upscaledValueArithmetic = fractureCellAreaXvalue / fractureCellArea; - return upscaledValueArithmetic; - -} - - //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- @@ -683,8 +540,6 @@ double RigFractureTransCalc::arithmeticAverage(std::vector values) return sumValue / numberOfValues; } - - //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- @@ -719,11 +574,12 @@ void RigFractureTransCalc::computeUpscaledPropertyFromStimPlan( QString resultNa RigFractureData fracData; fracData.reservoirCellIndex = fracCell; - - double upscaledAritmStimPlanValue = cvf::UNDEFINED_DOUBLE; - double upscaledHarmStimPlanValue = cvf::UNDEFINED_DOUBLE; - computeUpscaledPropertyFromStimPlanForEclipseCell(upscaledAritmStimPlanValue, upscaledHarmStimPlanValue, resultName, resultUnit, timeStepIndex, m_unitForCalculation, fracCell); - + + std::pair upscaledCondFlowAcrossLayers = flowAcrossLayersUpscaling(resultName, resultUnit, timeStepIndex, m_unitForCalculation, fracCell); + + double upscaledAritmStimPlanValue = upscaledCondFlowAcrossLayers.first; + double upscaledHarmStimPlanValue = upscaledCondFlowAcrossLayers.second; + if (upscaledAritmStimPlanValue != cvf::UNDEFINED_DOUBLE) { fracData.upscaledAritmStimPlanValue = upscaledAritmStimPlanValue; diff --git a/ApplicationCode/ReservoirDataModel/RigFractureTransCalc.h b/ApplicationCode/ReservoirDataModel/RigFractureTransCalc.h index a5214f4f2d..550101c6b7 100644 --- a/ApplicationCode/ReservoirDataModel/RigFractureTransCalc.h +++ b/ApplicationCode/ReservoirDataModel/RigFractureTransCalc.h @@ -32,7 +32,6 @@ #include #include - class RimFracture; class RimEclipseCase; class RimStimPlanCell; @@ -46,25 +45,19 @@ class RigFractureTransCalc public: explicit RigFractureTransCalc(RimEclipseCase* caseToApply, RimFracture* fracture); - void computeTransmissibility(); - bool planeCellIntersectionPolygons(size_t cellindex, std::vector > & polygons, cvf::Vec3d & localX, cvf::Vec3d & localY, cvf::Vec3d & localZ); + void computeTransmissibility(); + bool planeCellIntersectionPolygons(size_t cellindex, std::vector > & polygons, cvf::Vec3d & localX, cvf::Vec3d & localY, cvf::Vec3d & localZ); - void computeUpscaledPropertyFromStimPlan(QString resultName, QString resultUnit, size_t timeStepIndex); - void computeUpscaledPropertyFromStimPlanForEclipseCell(double &upscaledAritmStimPlanValue, double &upscaledHarmStimPlanValue, QString resultName, QString resultUnit, size_t timeStepIndex, caf::AppEnum< RimDefines::UnitSystem > unitSystem, size_t cellIndex); - - std::pair flowAcrossLayersUpscaling(QString resultName, QString resultUnit, size_t timeStepIndex, RimDefines::UnitSystem unitSystem, size_t eclipseCellIndex); - double computeHAupscale(RimStimPlanFractureTemplate* fracTemplateStimPlan, std::vector stimPlanCells, std::vector planeCellPolygon, cvf::Vec3d directionAlongLayers, cvf::Vec3d directionAcrossLayers); - double computeAHupscale(RimStimPlanFractureTemplate* fracTemplateStimPlan, std::vector stimPlanCells, std::vector planeCellPolygon, cvf::Vec3d directionAlongLayers, cvf::Vec3d directionAcrossLayers); + void computeUpscaledPropertyFromStimPlan(QString resultName, QString resultUnit, size_t timeStepIndex); + std::pair flowAcrossLayersUpscaling(QString resultName, QString resultUnit, size_t timeStepIndex, RimDefines::UnitSystem unitSystem, size_t eclipseCellIndex); + double computeHAupscale(RimStimPlanFractureTemplate* fracTemplateStimPlan, std::vector stimPlanCells, std::vector planeCellPolygon, cvf::Vec3d directionAlongLayers, cvf::Vec3d directionAcrossLayers); + double computeAHupscale(RimStimPlanFractureTemplate* fracTemplateStimPlan, std::vector stimPlanCells, std::vector planeCellPolygon, cvf::Vec3d directionAlongLayers, cvf::Vec3d directionAcrossLayers); + static double arithmeticAverage(std::vector values); - - static double areaWeightedHarmonicAverage(std::vector areaOfFractureParts, std::vector valuesForFractureParts); - static double areaWeightedArithmeticAverage(std::vector areaOfFractureParts, std::vector valuesForFractureParts); - static double arithmeticAverage(std::vector values); - - void computeFlowInFracture(); - void computeFlowIntoTransverseWell(); + void computeFlowInFracture(); + void computeFlowIntoTransverseWell(); static std::vector getRowOfStimPlanCells(std::vector allStimPlanCells, size_t i);