diff --git a/ApplicationCode/ProjectDataModel/RimStimPlanCell.cpp b/ApplicationCode/ProjectDataModel/RimStimPlanCell.cpp index 1817aecad1..4901f1958a 100644 --- a/ApplicationCode/ProjectDataModel/RimStimPlanCell.cpp +++ b/ApplicationCode/ProjectDataModel/RimStimPlanCell.cpp @@ -18,9 +18,6 @@ #include "RimStimPlanCell.h" - -CAF_PDM_SOURCE_INIT(RimStimPlanCell, "RimStimPlanCell"); - //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- diff --git a/ApplicationCode/ProjectDataModel/RimStimPlanCell.h b/ApplicationCode/ProjectDataModel/RimStimPlanCell.h index 92ec869e19..795fccc053 100644 --- a/ApplicationCode/ProjectDataModel/RimStimPlanCell.h +++ b/ApplicationCode/ProjectDataModel/RimStimPlanCell.h @@ -18,8 +18,6 @@ #pragma once -#include "cafPdmObject.h" - #include "cvfBase.h" #include "cvfVector3.h" @@ -29,9 +27,8 @@ /// /// //================================================================================================== -class RimStimPlanCell : public caf::PdmObject +class RimStimPlanCell { - CAF_PDM_HEADER_INIT; public: RimStimPlanCell(); diff --git a/ApplicationCode/ReservoirDataModel/RigFractureTransCalc.cpp b/ApplicationCode/ReservoirDataModel/RigFractureTransCalc.cpp index e8f9172677..8961d0918e 100644 --- a/ApplicationCode/ReservoirDataModel/RigFractureTransCalc.cpp +++ b/ApplicationCode/ReservoirDataModel/RigFractureTransCalc.cpp @@ -408,7 +408,7 @@ void RigFractureTransCalc::computeUpscaledPropertyFromStimPlanForEclipseCell(dou //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- -double RigFractureTransCalc::HAflowAcrossLayersUpscale(QString resultName, QString resultUnit, size_t timeStepIndex, RimDefines::UnitSystem unitSystem, size_t eclipseCellIndex) +std::pair RigFractureTransCalc::flowAcrossLayersUpscaling(QString resultName, QString resultUnit, size_t timeStepIndex, RimDefines::UnitSystem unitSystem, size_t eclipseCellIndex) { //TODO: A lot of common code with function for calculating transmissibility... @@ -418,7 +418,7 @@ double RigFractureTransCalc::HAflowAcrossLayersUpscale(QString resultName, QStri { fracTemplateStimPlan = dynamic_cast(m_fracture->attachedFractureDefinition()); } - else return cvf::UNDEFINED_DOUBLE; + else return std::make_pair(cvf::UNDEFINED_DOUBLE, cvf::UNDEFINED_DOUBLE); std::vector stimPlanCells = fracTemplateStimPlan->getStimPlanCells(resultName, resultUnit, timeStepIndex); @@ -431,7 +431,7 @@ double RigFractureTransCalc::HAflowAcrossLayersUpscale(QString resultName, QStri cvf::Vec3d localZ; std::vector > planeCellPolygons; bool isPlanIntersected = planeCellIntersectionPolygons(eclipseCellIndex, planeCellPolygons, localX, localY, localZ); - if (!isPlanIntersected || planeCellPolygons.size() == 0) return cvf::UNDEFINED_DOUBLE; + if (!isPlanIntersected || planeCellPolygons.size() == 0) return std::make_pair(cvf::UNDEFINED_DOUBLE, cvf::UNDEFINED_DOUBLE); //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(); @@ -458,79 +458,167 @@ double RigFractureTransCalc::HAflowAcrossLayersUpscale(QString resultName, QStri std::vector fracPolygon = m_fracture->attachedFractureDefinition()->fracturePolygon(unitSystem); std::vector > polygonsDescribingFractureInCell; - std::vector upscaledConductivities; + std::vector upscaledConductivitiesHA; + std::vector upscaledConductivitiesAH; + //Harmonic weighted mean for (std::vector planeCellPolygon : planeCellPolygons) { - std::vector DcolSum; - std::vector lavgCol; - std::vector CondHarmCol; - - for (size_t j = 0; j < fracTemplateStimPlan->stimPlanGridNumberOfColums(); j++) - { - std::vector conductivitiesInStimPlanCells; - std::vector lengthsLiOfStimPlanCol; - std::vector heightsDioFStimPlanCells; - - std::vector stimPlanCellsCol = getColOfStimPlanCells(stimPlanCells, j); - for (RimStimPlanCell* stimPlanCell : stimPlanCellsCol) - { - if (stimPlanCell->getValue() > 1e-7) - { - std::vector >clippedStimPlanPolygons = RigCellGeometryTools::clipPolygons(stimPlanCell->getPolygon(), planeCellPolygon); - if (clippedStimPlanPolygons.size() > 0) - { - for (auto clippedStimPlanPolygon : clippedStimPlanPolygons) - { - conductivitiesInStimPlanCells.push_back(stimPlanCell->getValue()); - lengthsLiOfStimPlanCol.push_back(RigCellGeometryTools::polygonAreaWeightedLength(directionAlongLayers, clippedStimPlanPolygon)); - heightsDioFStimPlanCells.push_back(RigCellGeometryTools::polygonAreaWeightedLength(directionAcrossLayers, clippedStimPlanPolygon)); - } - } - } - } - //Regne ut average - double sumDiDivCondLi = 0.0; - double sumDi = 0.0; - double sumLiDi = 0.0; - for (int i = 0; i < conductivitiesInStimPlanCells.size(); i++) - { - sumDiDivCondLi += heightsDioFStimPlanCells[i] / (conductivitiesInStimPlanCells[i] * lengthsLiOfStimPlanCol[i]); - sumDi += heightsDioFStimPlanCells[i]; - sumLiDi += heightsDioFStimPlanCells[i] * lengthsLiOfStimPlanCol[i]; - } - - if (sumDiDivCondLi != 0) - { - DcolSum.push_back(sumDi); - double lAvgValue = sumLiDi / sumDi; - lavgCol.push_back(lAvgValue); - CondHarmCol.push_back(1 / (sumLiDi*sumDiDivCondLi)); - } - } - - //Do arithmetic upscaling based on harmonic upscaled values for coloums - double sumCondHLiDivDi = 0.0; - double sumLi = 0.0; - double sumDiLi = 0.0; - for (int i = 0; i < CondHarmCol.size(); i++) - { - sumLi += lavgCol[i]; - sumDiLi += DcolSum[i] * lavgCol[i]; - sumCondHLiDivDi += CondHarmCol[i] * lavgCol[i] / DcolSum[i]; - } - double Davg = sumDiLi / sumLi; - double condHA = (Davg / sumLi) * sumCondHLiDivDi; - - upscaledConductivities.push_back(condHA); + double condHA = computeHAupscale(fracTemplateStimPlan, stimPlanCells, planeCellPolygon, directionAlongLayers, directionAcrossLayers); + upscaledConductivitiesHA.push_back(condHA); + double condAH = computeAHupscale(fracTemplateStimPlan, stimPlanCells, planeCellPolygon, directionAlongLayers, directionAcrossLayers); + upscaledConductivitiesAH.push_back(condAH); } //TODO: Is this the right way of handling getting several values for each cell? - return arithmeticAverage(upscaledConductivities); + return std::make_pair(arithmeticAverage(upscaledConductivitiesHA), arithmeticAverage(upscaledConductivitiesAH)); +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +double RigFractureTransCalc::computeHAupscale(RimStimPlanFractureTemplate* fracTemplateStimPlan, std::vector stimPlanCells, std::vector planeCellPolygon, cvf::Vec3d directionAlongLayers, cvf::Vec3d directionAcrossLayers) +{ + std::vector DcolSum; + std::vector lavgCol; + std::vector CondHarmCol; + + for (size_t j = 0; j < fracTemplateStimPlan->stimPlanGridNumberOfColums(); j++) + { + std::vector conductivitiesInStimPlanCells; + std::vector lengthsLiOfStimPlanCol; + std::vector heightsDiOfStimPlanCells; + + std::vector stimPlanCellsCol = getColOfStimPlanCells(stimPlanCells, j); + for (RimStimPlanCell* stimPlanCell : stimPlanCellsCol) + { + if (stimPlanCell->getValue() > 1e-7) + { + std::vector >clippedStimPlanPolygons = RigCellGeometryTools::clipPolygons(stimPlanCell->getPolygon(), planeCellPolygon); + if (clippedStimPlanPolygons.size() > 0) + { + for (auto clippedStimPlanPolygon : clippedStimPlanPolygons) + { + conductivitiesInStimPlanCells.push_back(stimPlanCell->getValue()); + lengthsLiOfStimPlanCol.push_back(RigCellGeometryTools::polygonAreaWeightedLength(directionAlongLayers, clippedStimPlanPolygon)); + heightsDiOfStimPlanCells.push_back(RigCellGeometryTools::polygonAreaWeightedLength(directionAcrossLayers, clippedStimPlanPolygon)); + } + } + } + } + //Regne ut average + double sumDiDivCondLi = 0.0; + double sumDi = 0.0; + double sumLiDi = 0.0; + for (int i = 0; i < conductivitiesInStimPlanCells.size(); i++) + { + sumDiDivCondLi += heightsDiOfStimPlanCells[i] / (conductivitiesInStimPlanCells[i] * lengthsLiOfStimPlanCol[i]); + sumDi += heightsDiOfStimPlanCells[i]; + sumLiDi += heightsDiOfStimPlanCells[i] * lengthsLiOfStimPlanCol[i]; + } + + if (sumDiDivCondLi != 0) + { + DcolSum.push_back(sumDi); + double lAvgValue = sumLiDi / sumDi; + lavgCol.push_back(lAvgValue); + CondHarmCol.push_back(1 / (sumLiDi*sumDiDivCondLi)); + } + } + + //Do arithmetic upscaling based on harmonic upscaled values for coloums + double sumCondHLiDivDi = 0.0; + double sumLi = 0.0; + double sumDiLi = 0.0; + for (int i = 0; i < CondHarmCol.size(); i++) + { + sumLi += lavgCol[i]; + sumDiLi += DcolSum[i] * lavgCol[i]; + sumCondHLiDivDi += CondHarmCol[i] * lavgCol[i] / DcolSum[i]; + } + double Davg = sumDiLi / sumLi; + double condHA = (Davg / sumLi) * sumCondHLiDivDi; + + return condHA; +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +double RigFractureTransCalc::computeAHupscale(RimStimPlanFractureTemplate* fracTemplateStimPlan, std::vector stimPlanCells, std::vector planeCellPolygon, cvf::Vec3d directionAlongLayers, cvf::Vec3d directionAcrossLayers) +{ + std::vector DcolAvg; + std::vector liColSum; + std::vector CondAritCol; + + for (size_t j = 0; j < fracTemplateStimPlan->stimPlanGridNumberOfColums(); j++) + { + std::vector conductivitiesInStimPlanCells; + std::vector lengthsLiOfStimPlanCol; + std::vector heightsDiOfStimPlanCells; + + std::vector stimPlanCellsCol = getColOfStimPlanCells(stimPlanCells, j); + for (RimStimPlanCell* stimPlanCell : stimPlanCellsCol) + { + if (stimPlanCell->getValue() > 1e-7) + { + std::vector >clippedStimPlanPolygons = RigCellGeometryTools::clipPolygons(stimPlanCell->getPolygon(), planeCellPolygon); + if (clippedStimPlanPolygons.size() > 0) + { + for (auto clippedStimPlanPolygon : clippedStimPlanPolygons) + { + conductivitiesInStimPlanCells.push_back(stimPlanCell->getValue()); + lengthsLiOfStimPlanCol.push_back(RigCellGeometryTools::polygonAreaWeightedLength(directionAlongLayers, clippedStimPlanPolygon)); + heightsDiOfStimPlanCells.push_back(RigCellGeometryTools::polygonAreaWeightedLength(directionAcrossLayers, clippedStimPlanPolygon)); + } + } + } + } + //Calculate sums needed for (arithmetic) average for coloum + double sumCondLiDivDi = 0.0; + double sumDi = 0.0; + double sumLiDi = 0.0; + double sumLi = 0.0; + for (int i = 0; i < conductivitiesInStimPlanCells.size(); i++) + { + sumCondLiDivDi += (conductivitiesInStimPlanCells[i] * lengthsLiOfStimPlanCol[i]) / heightsDiOfStimPlanCells[i]; + sumDi += heightsDiOfStimPlanCells[i]; + sumLiDi += heightsDiOfStimPlanCells[i] * lengthsLiOfStimPlanCol[i]; + sumLi += lengthsLiOfStimPlanCol[i]; + } + + if (sumCondLiDivDi != 0) + { + //Calculating art avg + double dAvg = sumLiDi / sumLi; + DcolAvg.push_back(dAvg); + liColSum.push_back(sumLi); + CondAritCol.push_back(dAvg / sumLi * sumCondLiDivDi); + } + } + + //Do harmonic upscaling based on arithmetric upscaled values for coloums + double sumDiDivCondALi = 0.0; + double sumDi = 0.0; + double sumDiLi = 0.0; + for (int i = 0; i < CondAritCol.size(); i++) + { + sumDi += DcolAvg[i]; + sumDiLi += DcolAvg[i] * liColSum[i]; + sumDiDivCondALi += DcolAvg[i] / (CondAritCol[i] * liColSum[i]); + + } + double Lavg = sumDiLi / sumDi; + double condAH = (sumDi / Lavg) * (1 / sumDiDivCondALi); + + return condAH; + + + } //-------------------------------------------------------------------------------------------------- @@ -803,4 +891,3 @@ double RigFractureTransCalc::convertConductivtyValue(double Kw, RimDefines::Unit return cvf::UNDEFINED_DOUBLE; } - diff --git a/ApplicationCode/ReservoirDataModel/RigFractureTransCalc.h b/ApplicationCode/ReservoirDataModel/RigFractureTransCalc.h index 5212dcda7b..a5214f4f2d 100644 --- a/ApplicationCode/ReservoirDataModel/RigFractureTransCalc.h +++ b/ApplicationCode/ReservoirDataModel/RigFractureTransCalc.h @@ -36,6 +36,7 @@ class RimFracture; class RimEclipseCase; class RimStimPlanCell; +class RimStimPlanFractureTemplate; //================================================================================================== /// @@ -52,7 +53,9 @@ public: 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); - double HAflowAcrossLayersUpscale(QString resultName, QString resultUnit, size_t timeStepIndex, RimDefines::UnitSystem unitSystem, size_t eclipseCellIndex); + 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);