#1364 - pre-proto - Cleanup code related to upscaling of conductivity values for flow across layers

This commit is contained in:
astridkbjorke 2017-03-29 09:43:36 +02:00
parent 2bfb68e17f
commit e9d4acb51a
2 changed files with 15 additions and 166 deletions

View File

@ -311,100 +311,6 @@ bool RigFractureTransCalc::planeCellIntersectionPolygons(size_t cellindex, std::
return isCellIntersected; 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<RimStimPlanFractureTemplate*>(m_fracture->attachedFractureDefinition()))
{
fracTemplateStimPlan = dynamic_cast<RimStimPlanFractureTemplate*>(m_fracture->attachedFractureDefinition());
}
else return;
std::vector<std::vector<cvf::Vec3d> > stimPlanCellsAsPolygons;
std::vector<double> stimPlanParameterValues;
fracTemplateStimPlan->getStimPlanDataAsPolygonsAndValues(stimPlanCellsAsPolygons, stimPlanParameterValues, resultName, resultUnit, timeStepIndex);
//TODO: A lot of common code with function above... Can be cleaned up...?
std::vector<size_t> 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<std::vector<cvf::Vec3d> > 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<cvf::Vec3d> & planeCellPolygon : planeCellPolygons)
{
for (cvf::Vec3d& v : planeCellPolygon)
{
v.transformPoint(static_cast<cvf::Mat4d>(invertedTransMatrix));
}
}
cvf::Vec3d localZinFracPlane;
localZinFracPlane = localZ;
localZinFracPlane.transformVector(static_cast<cvf::Mat4d>(invertedTransMatrix));
cvf::Vec3d directionOfLength = cvf::Vec3d::ZERO;
directionOfLength.cross(localZinFracPlane, cvf::Vec3d(0, 0, 1));
directionOfLength.normalize();
std::vector<cvf::Vec3f> fracPolygon = m_fracture->attachedFractureDefinition()->fracturePolygon(unitSystem);
std::vector<std::vector<cvf::Vec3d> > polygonsDescribingFractureInCell;
double area;
std::vector<double> areaOfFractureParts;
std::vector<double> valuesForFractureParts;
for (std::vector<cvf::Vec3d> planeCellPolygon : planeCellPolygons)
{
for (int i = 0; i < stimPlanParameterValues.size(); i++)
{
double stimPlanParameterValue = stimPlanParameterValues[i];
if (stimPlanParameterValue != 0)
{
std::vector<cvf::Vec3d> stimPlanCell = stimPlanCellsAsPolygons[i];
std::vector<std::vector<cvf::Vec3d> >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<double, double> RigFractureTransCalc::flowAcrossLayersUpscaling(QStrin
} }
} }
//Vector for vertical - på tvers av lag... Gitt av orientering av stimPlan-grid...
cvf::Vec3d directionAcrossLayers; cvf::Vec3d directionAcrossLayers;
cvf::Vec3d directionAlongLayers; cvf::Vec3d directionAlongLayers;
//TODO: Get these vectors properly... //TODO: Get these vectors properly...
directionAcrossLayers = cvf::Vec3d(0, -1, 0); directionAcrossLayers = cvf::Vec3d(0, -1, 0);
directionAlongLayers = cvf::Vec3d(1, 0, 0); directionAlongLayers = cvf::Vec3d(1, 0, 0);
@ -461,8 +364,6 @@ std::pair<double, double> RigFractureTransCalc::flowAcrossLayersUpscaling(QStrin
std::vector<double> upscaledConductivitiesHA; std::vector<double> upscaledConductivitiesHA;
std::vector<double> upscaledConductivitiesAH; std::vector<double> upscaledConductivitiesAH;
//Harmonic weighted mean //Harmonic weighted mean
for (std::vector<cvf::Vec3d> planeCellPolygon : planeCellPolygons) for (std::vector<cvf::Vec3d> planeCellPolygon : planeCellPolygons)
{ {
@ -621,50 +522,6 @@ double RigFractureTransCalc::computeAHupscale(RimStimPlanFractureTemplate* fracT
} }
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
double RigFractureTransCalc::areaWeightedHarmonicAverage(std::vector<double> areaOfFractureParts, std::vector<double> 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<double> areaOfFractureParts, std::vector<double> 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<double> values)
return sumValue / numberOfValues; return sumValue / numberOfValues;
} }
//-------------------------------------------------------------------------------------------------- //--------------------------------------------------------------------------------------------------
/// ///
//-------------------------------------------------------------------------------------------------- //--------------------------------------------------------------------------------------------------
@ -720,9 +575,10 @@ void RigFractureTransCalc::computeUpscaledPropertyFromStimPlan( QString resultNa
RigFractureData fracData; RigFractureData fracData;
fracData.reservoirCellIndex = fracCell; fracData.reservoirCellIndex = fracCell;
double upscaledAritmStimPlanValue = cvf::UNDEFINED_DOUBLE; std::pair<double, double> upscaledCondFlowAcrossLayers = flowAcrossLayersUpscaling(resultName, resultUnit, timeStepIndex, m_unitForCalculation, fracCell);
double upscaledHarmStimPlanValue = cvf::UNDEFINED_DOUBLE;
computeUpscaledPropertyFromStimPlanForEclipseCell(upscaledAritmStimPlanValue, upscaledHarmStimPlanValue, resultName, resultUnit, timeStepIndex, m_unitForCalculation, fracCell); double upscaledAritmStimPlanValue = upscaledCondFlowAcrossLayers.first;
double upscaledHarmStimPlanValue = upscaledCondFlowAcrossLayers.second;
if (upscaledAritmStimPlanValue != cvf::UNDEFINED_DOUBLE) if (upscaledAritmStimPlanValue != cvf::UNDEFINED_DOUBLE)
{ {

View File

@ -32,7 +32,6 @@
#include <vector> #include <vector>
#include <QString> #include <QString>
class RimFracture; class RimFracture;
class RimEclipseCase; class RimEclipseCase;
class RimStimPlanCell; class RimStimPlanCell;
@ -51,16 +50,10 @@ public:
void computeUpscaledPropertyFromStimPlan(QString resultName, QString resultUnit, size_t timeStepIndex); 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<double, double> flowAcrossLayersUpscaling(QString resultName, QString resultUnit, size_t timeStepIndex, RimDefines::UnitSystem unitSystem, size_t eclipseCellIndex); std::pair<double, double> flowAcrossLayersUpscaling(QString resultName, QString resultUnit, size_t timeStepIndex, RimDefines::UnitSystem unitSystem, size_t eclipseCellIndex);
double computeHAupscale(RimStimPlanFractureTemplate* fracTemplateStimPlan, std::vector<RimStimPlanCell *> stimPlanCells, std::vector<cvf::Vec3d> planeCellPolygon, cvf::Vec3d directionAlongLayers, cvf::Vec3d directionAcrossLayers); double computeHAupscale(RimStimPlanFractureTemplate* fracTemplateStimPlan, std::vector<RimStimPlanCell *> stimPlanCells, std::vector<cvf::Vec3d> planeCellPolygon, cvf::Vec3d directionAlongLayers, cvf::Vec3d directionAcrossLayers);
double computeAHupscale(RimStimPlanFractureTemplate* fracTemplateStimPlan, std::vector<RimStimPlanCell *> stimPlanCells, std::vector<cvf::Vec3d> planeCellPolygon, cvf::Vec3d directionAlongLayers, cvf::Vec3d directionAcrossLayers); double computeAHupscale(RimStimPlanFractureTemplate* fracTemplateStimPlan, std::vector<RimStimPlanCell *> stimPlanCells, std::vector<cvf::Vec3d> planeCellPolygon, cvf::Vec3d directionAlongLayers, cvf::Vec3d directionAcrossLayers);
static double areaWeightedHarmonicAverage(std::vector<double> areaOfFractureParts, std::vector<double> valuesForFractureParts);
static double areaWeightedArithmeticAverage(std::vector<double> areaOfFractureParts, std::vector<double> valuesForFractureParts);
static double arithmeticAverage(std::vector<double> values); static double arithmeticAverage(std::vector<double> values);
void computeFlowInFracture(); void computeFlowInFracture();