Add fault reactivation modeling

This commit is contained in:
jonjenssen
2023-05-03 02:20:05 +02:00
committed by Kristian Bendiksen
parent 14c37c4b2c
commit 32326bfa9b
43 changed files with 2626 additions and 157 deletions

View File

@@ -89,6 +89,9 @@ set(SOURCE_GROUP_HEADER_FILES
${CMAKE_CURRENT_LIST_DIR}/RigTexturedSection.h
${CMAKE_CURRENT_LIST_DIR}/RigPressureDepthData.h
${CMAKE_CURRENT_LIST_DIR}/RigMswCenterLineCalculator.h
${CMAKE_CURRENT_LIST_DIR}/RigBasicPlane.h
${CMAKE_CURRENT_LIST_DIR}/RigGriddedPart3d.h
${CMAKE_CURRENT_LIST_DIR}/RigFaultReactivationModel.h
${CMAKE_CURRENT_LIST_DIR}/RigWellAllocationOverTime.h
${CMAKE_CURRENT_LIST_DIR}/RigWellResultBranch.h
${CMAKE_CURRENT_LIST_DIR}/RigWellResultFrame.h
@@ -178,6 +181,9 @@ set(SOURCE_GROUP_SOURCE_FILES
${CMAKE_CURRENT_LIST_DIR}/RigTexturedSection.cpp
${CMAKE_CURRENT_LIST_DIR}/RigPressureDepthData.cpp
${CMAKE_CURRENT_LIST_DIR}/RigMswCenterLineCalculator.cpp
${CMAKE_CURRENT_LIST_DIR}/RigBasicPlane.cpp
${CMAKE_CURRENT_LIST_DIR}/RigGriddedPart3d.cpp
${CMAKE_CURRENT_LIST_DIR}/RigFaultReactivationModel.cpp
${CMAKE_CURRENT_LIST_DIR}/RigWellAllocationOverTime.cpp
${CMAKE_CURRENT_LIST_DIR}/RigWellResultBranch.cpp
${CMAKE_CURRENT_LIST_DIR}/RigWellResultFrame.cpp

View File

@@ -0,0 +1,167 @@
/////////////////////////////////////////////////////////////////////////////////
//
// Copyright (C) 2023 Equinor ASA
//
// ResInsight is free software: you can redistribute it and/or modify
// it under the terms of the GNU General Public License as published by
// the Free Software Foundation, either version 3 of the License, or
// (at your option) any later version.
//
// ResInsight is distributed in the hope that it will be useful, but WITHOUT ANY
// WARRANTY; without even the implied warranty of MERCHANTABILITY or
// FITNESS FOR A PARTICULAR PURPOSE.
//
// See the GNU General Public License at <http://www.gnu.org/licenses/gpl.html>
// for more details.
//
/////////////////////////////////////////////////////////////////////////////////
#include "RigBasicPlane.h"
#include "cvfTextureImage.h"
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
RigBasicPlane::RigBasicPlane()
: m_isRectValid( false )
, m_maxHorzExtent( 0.0 )
, m_maxVertExtentAbove( 0.0 )
, m_maxVertExtentBelow( 0.0 )
{
m_texture = new cvf::TextureImage();
m_texture->allocate( 1, 1 );
m_texture->fill( cvf::Color4ub( 0, 0, 0, 0 ) );
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
RigBasicPlane::~RigBasicPlane()
{
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
void RigBasicPlane::reset()
{
m_isRectValid = false;
m_rect.clear();
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
bool RigBasicPlane::isValid() const
{
return m_isRectValid;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
void RigBasicPlane::setPlane( cvf::Vec3d anchorPoint, cvf::Vec3d normal )
{
m_planeAnchor = anchorPoint;
m_planeNormal = normal;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
void RigBasicPlane::setColor( cvf::Color3f color )
{
m_color = color;
m_texture->fill( cvf::Color4ub( color.rByte(), color.gByte(), color.bByte(), 255 ) );
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
void RigBasicPlane::setMaxExtentFromAnchor( double maxExtentHorz, double maxExtentVertAbove, double maxExtentVertBelow )
{
m_maxHorzExtent = maxExtentHorz;
m_maxVertExtentAbove = maxExtentVertAbove;
m_maxVertExtentBelow = maxExtentVertBelow;
}
//--------------------------------------------------------------------------------------------------
/// ti
/// 1 ----------- 2
/// | |
/// ml | | mr
/// | |
/// 0 ----------- 3
/// bi
//--------------------------------------------------------------------------------------------------
void RigBasicPlane::updateRect()
{
if ( ( m_maxHorzExtent <= 0.0 ) || ( m_maxVertExtentAbove <= 0.0 ) || ( m_maxVertExtentBelow <= 0.0 ) )
{
m_isRectValid = false;
return;
}
cvf::Vec3d zDirection( 0, 0, 1 );
auto alongPlane = m_planeNormal ^ zDirection;
auto upwards = m_planeNormal ^ alongPlane;
upwards.normalize();
alongPlane.normalize();
const double extHorz = m_maxHorzExtent / 2.0;
auto ml = m_planeAnchor + alongPlane * extHorz;
auto mr = m_planeAnchor - alongPlane * extHorz;
m_rect.resize( 4 );
m_rect[0] = ml - upwards * m_maxVertExtentAbove;
m_rect[1] = ml + upwards * m_maxVertExtentBelow;
m_rect[2] = mr + upwards * m_maxVertExtentBelow;
m_rect[3] = mr - upwards * m_maxVertExtentAbove;
m_isRectValid = true;
m_topIntersect = m_planeAnchor - upwards * m_maxVertExtentAbove;
m_bottomIntersect = m_planeAnchor + upwards * m_maxVertExtentBelow;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
cvf::Vec3dArray RigBasicPlane::rect() const
{
return m_rect;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
cvf::ref<cvf::TextureImage> RigBasicPlane::texture() const
{
return m_texture;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
double RigBasicPlane::maxDepth()
{
if ( !m_isRectValid ) return 0.0;
double maxdepth = 0.0;
for ( auto p : m_rect )
{
maxdepth = std::max( maxdepth, std::abs( p.z() ) );
}
return maxdepth;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
std::pair<cvf::Vec3d, cvf::Vec3d> RigBasicPlane::intersectTopBottomLine()
{
return std::make_pair( m_topIntersect, m_bottomIntersect );
}

View File

@@ -0,0 +1,75 @@
/////////////////////////////////////////////////////////////////////////////////
//
// Copyright (C) 2023 Equinor ASA
//
// ResInsight is free software: you can redistribute it and/or modify
// it under the terms of the GNU General Public License as published by
// the Free Software Foundation, either version 3 of the License, or
// (at your option) any later version.
//
// ResInsight is distributed in the hope that it will be useful, but WITHOUT ANY
// WARRANTY; without even the implied warranty of MERCHANTABILITY or
// FITNESS FOR A PARTICULAR PURPOSE.
//
// See the GNU General Public License at <http://www.gnu.org/licenses/gpl.html>
// for more details.
//
/////////////////////////////////////////////////////////////////////////////////
#pragma once
#include "cvfArray.h"
#include "cvfColor3.h"
#include "cvfObject.h"
#include "cvfPlane.h"
#include "cvfVector3.h"
#include <memory>
#include <vector>
namespace cvf
{
class TextureImage;
}
//==================================================================================================
///
///
//==================================================================================================
class RigBasicPlane : public cvf::Object
{
public:
RigBasicPlane();
~RigBasicPlane() override;
bool isValid() const;
void reset();
void updateRect();
void setPlane( cvf::Vec3d anchorPoint, cvf::Vec3d normal );
void setMaxExtentFromAnchor( double maxExtentHorz, double maxExtentVertAbove, double maxExtentVertBelow );
void setColor( cvf::Color3f color );
double maxDepth();
std::pair<cvf::Vec3d, cvf::Vec3d> intersectTopBottomLine();
cvf::Vec3dArray rect() const;
cvf::ref<cvf::TextureImage> texture() const;
private:
cvf::Vec3d m_planeNormal;
cvf::Vec3d m_planeAnchor;
double m_maxHorzExtent;
double m_maxVertExtentAbove;
double m_maxVertExtentBelow;
cvf::Vec3d m_topIntersect;
cvf::Vec3d m_bottomIntersect;
cvf::Vec3dArray m_rect;
bool m_isRectValid;
cvf::Color3f m_color;
cvf::ref<cvf::TextureImage> m_texture;
};

View File

@@ -0,0 +1,308 @@
/////////////////////////////////////////////////////////////////////////////////
//
// Copyright (C) 2023 Equinor ASA
//
// ResInsight is free software: you can redistribute it and/or modify
// it under the terms of the GNU General Public License as published by
// the Free Software Foundation, either version 3 of the License, or
// (at your option) any later version.
//
// ResInsight is distributed in the hope that it will be useful, but WITHOUT ANY
// WARRANTY; without even the implied warranty of MERCHANTABILITY or
// FITNESS FOR A PARTICULAR PURPOSE.
//
// See the GNU General Public License at <http://www.gnu.org/licenses/gpl.html>
// for more details.
//
/////////////////////////////////////////////////////////////////////////////////
#include "RigFaultReactivationModel.h"
#include "RigGriddedPart3d.h"
#include "RigPolyLinesData.h"
#include "cvfTextureImage.h"
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
RigFaultReactivationModel::RigFaultReactivationModel()
: m_maxZ( 0 )
, m_minZ( 0 )
, m_maxHorzExtent( 0 )
, m_isValid( false )
, m_cellCountHorzPart1( 1 )
, m_cellCountHorzPart2( 1 )
, m_cellCountVertUpper( 1 )
, m_cellCountVertMiddle( 1 )
, m_cellCountVertLower( 1 )
{
for ( auto part : allModelParts() )
{
m_parts[part] = RigFRModelPart();
m_parts[part].texture = new cvf::TextureImage();
m_parts[part].texture->allocate( 1, 1 );
m_parts[part].texture->fill( cvf::Color4ub( 0, 0, 0, 0 ) );
m_parts[part].rect.reserve( 4 );
}
m_cornerIndexes[ModelParts::HiPart1] = { 2, 3, 7, 6 };
m_cornerIndexes[ModelParts::MidPart1] = { 1, 2, 6, 5 };
m_cornerIndexes[ModelParts::LowPart1] = { 0, 1, 5, 4 };
m_cornerIndexes[ModelParts::HiPart2] = { 6, 7, 11, 10 };
m_cornerIndexes[ModelParts::MidPart2] = { 5, 6, 10, 9 };
m_cornerIndexes[ModelParts::LowPart2] = { 4, 5, 9, 8 };
for ( auto part : allGridParts() )
{
m_3dparts[part] = std::make_shared<RigGriddedPart3d>();
}
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
RigFaultReactivationModel::~RigFaultReactivationModel()
{
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
std::vector<RigFaultReactivationModel::ModelParts> RigFaultReactivationModel::allModelParts() const
{
return { ModelParts::HiPart1, ModelParts::MidPart1, ModelParts::LowPart1, ModelParts::HiPart2, ModelParts::MidPart2, ModelParts::LowPart2 };
}
std::vector<RigFaultReactivationModel::GridPart> RigFaultReactivationModel::allGridParts() const
{
return { GridPart::PART1, GridPart::PART2 };
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
void RigFaultReactivationModel::reset()
{
m_isValid = false;
for ( auto part : allModelParts() )
{
m_parts[part].rect.clear();
m_parts[part].rect.reserve( 4 );
}
for ( auto part : allGridParts() )
{
m_3dparts[part]->reset();
}
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
bool RigFaultReactivationModel::isValid() const
{
return m_isValid;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
void RigFaultReactivationModel::setPlane( cvf::Vec3d anchorPoint, cvf::Vec3d normal )
{
m_planeAnchor = anchorPoint;
m_planeNormal = normal;
reset();
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
void RigFaultReactivationModel::setPartColors( cvf::Color3f part1Color, cvf::Color3f part2Color )
{
for ( auto part : { ModelParts::HiPart1, ModelParts::MidPart1, ModelParts::LowPart1 } )
{
m_parts[part].texture->fill( cvf::Color4ub( part1Color.rByte(), part1Color.gByte(), part1Color.bByte(), 255 ) );
}
for ( auto part : { ModelParts::HiPart2, ModelParts::MidPart2, ModelParts::LowPart2 } )
{
m_parts[part].texture->fill( cvf::Color4ub( part2Color.rByte(), part2Color.gByte(), part2Color.bByte(), 255 ) );
}
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
void RigFaultReactivationModel::setMaxExtentFromAnchor( double maxExtentHorz, double minZ, double maxZ )
{
m_maxHorzExtent = maxExtentHorz;
m_minZ = minZ;
m_maxZ = maxZ;
reset();
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
void RigFaultReactivationModel::setFaultPlaneIntersect( cvf::Vec3d faultPlaneTop, cvf::Vec3d faultPlaneBottom )
{
m_faultPlaneIntersectBottom = faultPlaneBottom;
m_faultPlaneIntersectTop = faultPlaneTop;
reset();
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
void RigFaultReactivationModel::setCellCounts( int horzPart1, int horzPart2, int vertUpper, int vertMiddle, int vertLower )
{
m_cellCountHorzPart1 = horzPart1;
m_cellCountHorzPart2 = horzPart2;
m_cellCountVertUpper = vertUpper;
m_cellCountVertMiddle = vertMiddle;
m_cellCountVertLower = vertLower;
reset();
}
//--------------------------------------------------------------------------------------------------
/// 7
/// 3----------|----------- 11
/// | | |
/// | | |
/// | | |
/// 2|---------|----------| 10
/// | \6 |
/// | X Anchor |
/// | \ |
/// 1-------------|------- 9
/// | 5| |
/// | | |
/// | | |
/// | | |
/// 0-------------|-------- 8
/// 4
///
///
///
//--------------------------------------------------------------------------------------------------
void RigFaultReactivationModel::updateRects()
{
reset();
if ( ( m_maxHorzExtent <= 0.0 ) || ( m_minZ == m_maxZ ) )
{
return;
}
cvf::Vec3d zDir( 0, 0, 1 );
auto alongPlane = m_planeNormal ^ zDir;
alongPlane.normalize();
// how far from anchor point we should stop
const double extHorz = m_maxHorzExtent / 2.0;
auto mr = m_planeAnchor + alongPlane * extHorz;
auto ml = m_planeAnchor - alongPlane * extHorz;
cvf::Vec3dArray points;
points.resize( 12 );
points[0] = ml;
points[0].z() = -m_maxZ;
points[1] = ml;
points[1].z() = m_faultPlaneIntersectBottom.z();
points[2] = ml;
points[2].z() = m_faultPlaneIntersectTop.z();
points[3] = ml;
points[3].z() = -m_minZ;
points[4] = m_faultPlaneIntersectBottom;
points[4].z() = -m_maxZ;
points[5] = m_faultPlaneIntersectBottom;
points[6] = m_faultPlaneIntersectTop;
points[7] = m_faultPlaneIntersectTop;
points[7].z() = -m_minZ;
points[8] = mr;
points[8].z() = -m_maxZ;
points[9] = mr;
points[9].z() = m_faultPlaneIntersectBottom.z();
points[10] = mr;
points[10].z() = m_faultPlaneIntersectTop.z();
points[11] = mr;
points[11].z() = -m_minZ;
for ( auto part : allModelParts() )
{
for ( auto i : m_cornerIndexes[part] )
{
m_parts[part].rect.push_back( points[i] );
}
}
m_isValid = true;
generateGrids( points );
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
std::vector<cvf::Vec3d> RigFaultReactivationModel::rect( ModelParts part ) const
{
return m_parts.at( part ).rect;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
cvf::ref<cvf::TextureImage> RigFaultReactivationModel::texture( ModelParts part ) const
{
return m_parts.at( part ).texture;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
const std::vector<std::vector<cvf::Vec3d>>& RigFaultReactivationModel::meshLines( GridPart part ) const
{
return m_3dparts.at( part )->meshLines();
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
void RigFaultReactivationModel::generateGrids( cvf::Vec3dArray points )
{
m_3dparts[GridPart::PART1]->generateGeometry( { points[0], points[1], points[2], points[3], points[4], points[5], points[6], points[7] },
m_cellCountHorzPart1,
m_cellCountVertLower,
m_cellCountVertMiddle,
m_cellCountVertUpper );
m_3dparts[GridPart::PART2]->generateGeometry( { points[8], points[9], points[10], points[11], points[4], points[5], points[6], points[7] },
m_cellCountHorzPart2,
m_cellCountVertLower,
m_cellCountVertMiddle,
m_cellCountVertUpper );
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
std::shared_ptr<RigGriddedPart3d> RigFaultReactivationModel::grid( GridPart part ) const
{
return m_3dparts.at( part );
}

View File

@@ -0,0 +1,124 @@
/////////////////////////////////////////////////////////////////////////////////
//
// Copyright (C) 2023 Equinor ASA
//
// ResInsight is free software: you can redistribute it and/or modify
// it under the terms of the GNU General Public License as published by
// the Free Software Foundation, either version 3 of the License, or
// (at your option) any later version.
//
// ResInsight is distributed in the hope that it will be useful, but WITHOUT ANY
// WARRANTY; without even the implied warranty of MERCHANTABILITY or
// FITNESS FOR A PARTICULAR PURPOSE.
//
// See the GNU General Public License at <http://www.gnu.org/licenses/gpl.html>
// for more details.
//
/////////////////////////////////////////////////////////////////////////////////
#pragma once
#include "cvfArray.h"
#include "cvfColor3.h"
#include "cvfObject.h"
#include "cvfPlane.h"
#include "cvfVector3.h"
#include <map>
#include <memory>
#include <vector>
namespace cvf
{
class TextureImage;
}
class RigGriddedPart3d;
class RigFRModelPart
{
public:
RigFRModelPart(){};
~RigFRModelPart(){};
std::vector<cvf::Vec3d> rect;
cvf::ref<cvf::TextureImage> texture;
};
//==================================================================================================
///
///
//==================================================================================================
class RigFaultReactivationModel : public cvf::Object
{
public:
enum class ModelParts
{
HiPart1 = 0,
MidPart1,
LowPart1,
HiPart2,
MidPart2,
LowPart2
};
enum class GridPart
{
PART1,
PART2
};
public:
RigFaultReactivationModel();
~RigFaultReactivationModel() override;
std::vector<ModelParts> allModelParts() const;
std::vector<GridPart> allGridParts() const;
bool isValid() const;
void reset();
void setPlane( cvf::Vec3d anchorPoint, cvf::Vec3d normal );
void setFaultPlaneIntersect( cvf::Vec3d faultPlaneTop, cvf::Vec3d faultPlaneBottom );
void setMaxExtentFromAnchor( double maxExtentHorz, double minZ, double maxZ );
void setCellCounts( int horzPart1, int horzPart2, int vertUpper, int vertMiddle, int vertLower );
void updateRects();
void setPartColors( cvf::Color3f part1Color, cvf::Color3f part2Color );
std::vector<cvf::Vec3d> rect( ModelParts part ) const;
cvf::ref<cvf::TextureImage> texture( ModelParts part ) const;
const std::vector<std::vector<cvf::Vec3d>>& meshLines( GridPart part ) const;
std::shared_ptr<RigGriddedPart3d> grid( GridPart part ) const;
protected:
void generateGrids( cvf::Vec3dArray points );
private:
cvf::Vec3d m_planeNormal;
cvf::Vec3d m_planeAnchor;
cvf::Vec3d m_faultPlaneIntersectTop;
cvf::Vec3d m_faultPlaneIntersectBottom;
double m_maxHorzExtent;
double m_minZ;
double m_maxZ;
int m_cellCountHorzPart1;
int m_cellCountHorzPart2;
int m_cellCountVertUpper;
int m_cellCountVertMiddle;
int m_cellCountVertLower;
std::map<ModelParts, std::vector<int>> m_cornerIndexes;
std::map<ModelParts, RigFRModelPart> m_parts;
bool m_isValid;
std::map<GridPart, std::shared_ptr<RigGriddedPart3d>> m_3dparts;
};

View File

@@ -0,0 +1,275 @@
/////////////////////////////////////////////////////////////////////////////////
//
// Copyright (C) 2023 Equinor ASA
//
// ResInsight is free software: you can redistribute it and/or modify
// it under the terms of the GNU General Public License as published by
// the Free Software Foundation, either version 3 of the License, or
// (at your option) any later version.
//
// ResInsight is distributed in the hope that it will be useful, but WITHOUT ANY
// WARRANTY; without even the implied warranty of MERCHANTABILITY or
// FITNESS FOR A PARTICULAR PURPOSE.
//
// See the GNU General Public License at <http://www.gnu.org/licenses/gpl.html>
// for more details.
//
/////////////////////////////////////////////////////////////////////////////////
#include "RigGriddedPart3d.h"
#include "cvfTextureImage.h"
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
RigGriddedPart3d::RigGriddedPart3d()
: m_thickness( 10.0 )
{
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
RigGriddedPart3d::~RigGriddedPart3d()
{
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
void RigGriddedPart3d::reset()
{
m_borderSurfaceElements.clear();
m_vertices.clear();
m_elementIndices.clear();
m_meshLines.clear();
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
cvf::Vec3d RigGriddedPart3d::stepVector( cvf::Vec3d start, cvf::Vec3d stop, int nSteps )
{
cvf::Vec3d vec = stop - start;
return vec.getNormalized() * ( vec.length() / nSteps );
}
//--------------------------------------------------------------------------------------------------
/// Point index in input
///
///
/// 3 ----------- 7 *
/// | | *
/// | | *
/// | | *
/// 2 |---------| 6 *
/// | \ *
/// | \ *
/// | \ *
/// 1 -------------| 5 *
/// | | *
/// | | *
/// | | *
/// | | *
/// 0 -------------- 4 *
///
/// Assumes 0->4, 1->5, 2->6 and 3->7 is parallel
///
///
//--------------------------------------------------------------------------------------------------
void RigGriddedPart3d::generateGeometry( std::vector<cvf::Vec3d> inputPoints, int nHorzCells, int nVertCellsLower, int nVertCellsMiddle, int nVertCellsUpper )
{
reset();
cvf::Vec3d step0to1 = stepVector( inputPoints[0], inputPoints[1], nVertCellsLower );
cvf::Vec3d step1to2 = stepVector( inputPoints[1], inputPoints[2], nVertCellsMiddle );
cvf::Vec3d step2to3 = stepVector( inputPoints[2], inputPoints[3], nVertCellsUpper );
cvf::Vec3d step4to5 = stepVector( inputPoints[4], inputPoints[5], nVertCellsLower );
cvf::Vec3d step5to6 = stepVector( inputPoints[5], inputPoints[6], nVertCellsMiddle );
cvf::Vec3d step6to7 = stepVector( inputPoints[6], inputPoints[7], nVertCellsUpper );
cvf::Vec3d step0to4 = stepVector( inputPoints[0], inputPoints[4], nHorzCells );
cvf::Vec3d tVec = step0to4 ^ step0to1;
tVec.normalize();
tVec *= m_thickness;
const std::vector<double> m_thicknessFactors = { -1.0, 0.0, 1.0 };
const int nThicknessCells = 2;
const int nVertCells = nVertCellsLower + nVertCellsMiddle + nVertCellsUpper;
const std::vector<int> vertCells = { nVertCellsLower, nVertCellsMiddle, nVertCellsUpper + 1 };
const std::vector<cvf::Vec3d> firstSteps = { step0to1, step1to2, step2to3 };
const std::vector<cvf::Vec3d> lastSteps = { step4to5, step5to6, step6to7 };
// ** generate vertices
m_vertices.reserve( (size_t)( ( nVertCells + 1 ) * ( nHorzCells + 1 ) ) );
cvf::Vec3d p = inputPoints[0];
cvf::Vec3d pLast = inputPoints[4];
for ( int i = 0; i < (int)vertCells.size(); i++ )
{
for ( int v = 0; v < vertCells[i]; v++ )
{
cvf::Vec3d stepHorz = stepVector( p, pLast, nHorzCells );
cvf::Vec3d p2 = p;
for ( int h = 0; h <= nHorzCells; h++ )
{
for ( int t = 0; t < (int)m_thicknessFactors.size(); t++ )
{
m_vertices.push_back( p2 + m_thicknessFactors[t] * tVec );
}
p2 += stepHorz;
pLast = p2;
}
p += firstSteps[i];
pLast += lastSteps[i];
}
}
// ** generate elements of type hex8
m_elementIndices.resize( (size_t)( nVertCells * nHorzCells * nThicknessCells ) );
m_borderSurfaceElements[BorderSurface::UpperSurface] = {};
m_borderSurfaceElements[BorderSurface::FaultSurface] = {};
m_borderSurfaceElements[BorderSurface::LowerSurface] = {};
int layerIndex = 0;
int elementIdx = 0;
BorderSurface currentRegion = BorderSurface::LowerSurface;
const int nextLayerIdxOff = ( nHorzCells + 1 ) * ( nThicknessCells + 1 );
const int nThicknessOff = nThicknessCells + 1;
for ( int v = 0; v < nVertCells; v++ )
{
if ( v >= nVertCellsLower ) currentRegion = BorderSurface::FaultSurface;
if ( v >= nVertCellsLower + nVertCellsMiddle ) currentRegion = BorderSurface::UpperSurface;
int i = layerIndex;
for ( int h = 0; h < nHorzCells; h++ )
{
for ( int t = 0; t < nThicknessCells; t++ )
{
m_elementIndices[elementIdx].push_back( t + i );
m_elementIndices[elementIdx].push_back( t + i + nThicknessOff );
m_elementIndices[elementIdx].push_back( t + i + nThicknessOff + 1 );
m_elementIndices[elementIdx].push_back( t + i + 1 );
m_elementIndices[elementIdx].push_back( t + nextLayerIdxOff + i );
m_elementIndices[elementIdx].push_back( t + nextLayerIdxOff + i + nThicknessOff );
m_elementIndices[elementIdx].push_back( t + nextLayerIdxOff + i + nThicknessOff + 1 );
m_elementIndices[elementIdx].push_back( t + nextLayerIdxOff + i + 1 );
elementIdx++;
}
i += nThicknessOff;
}
// add elements to border surface in current region
m_borderSurfaceElements[currentRegion].push_back( elementIdx - 2 );
m_borderSurfaceElements[currentRegion].push_back( elementIdx - 1 );
layerIndex += nextLayerIdxOff;
}
// generate meshlines for 2d viz
generateMeshlines( { inputPoints[0], inputPoints[1], inputPoints[5], inputPoints[4] }, nHorzCells, nVertCellsLower );
generateMeshlines( { inputPoints[1], inputPoints[2], inputPoints[6], inputPoints[5] }, nHorzCells, nVertCellsMiddle );
generateMeshlines( { inputPoints[2], inputPoints[3], inputPoints[7], inputPoints[6] }, nHorzCells, nVertCellsUpper );
}
//--------------------------------------------------------------------------------------------------
/// Point index in input
///
/// 1 ____________ 2
/// | /
/// | /
/// | /
/// | /
/// |_______/
/// 0 3
///
/// Assumes 0->3 and 1->2 is parallel
//--------------------------------------------------------------------------------------------------
void RigGriddedPart3d::generateMeshlines( std::vector<cvf::Vec3d> cornerPoints, int numHorzCells, int numVertCells )
{
cvf::Vec3d step0to1 = stepVector( cornerPoints[0], cornerPoints[1], numVertCells );
cvf::Vec3d step0to3 = stepVector( cornerPoints[0], cornerPoints[3], numHorzCells );
cvf::Vec3d step1to2 = stepVector( cornerPoints[1], cornerPoints[2], numHorzCells );
cvf::Vec3d step3to2 = stepVector( cornerPoints[3], cornerPoints[2], numVertCells );
// horizontal lines
cvf::Vec3d startP = cornerPoints[0];
cvf::Vec3d endP = cornerPoints[3];
for ( int v = 0; v <= numVertCells; v++ )
{
m_meshLines.push_back( { startP, endP } );
startP += step0to1;
endP += step3to2;
}
// vertical lines
startP = cornerPoints[0];
endP = cornerPoints[1];
for ( int h = 0; h <= numHorzCells; h++ )
{
m_meshLines.push_back( { startP, endP } );
startP += step0to3;
endP += step1to2;
}
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
const std::vector<cvf::Vec3d>& RigGriddedPart3d::vertices() const
{
return m_vertices;
}
//--------------------------------------------------------------------------------------------------
/// Output elements will be of type HEX8
///
/// 7---------6
/// /| /|
/// / | / |
/// 4---------5 | z
/// | 3------|--2 | y
/// | / | / | /
/// |/ |/ |/
/// 0---------1 ----- x
///
//--------------------------------------------------------------------------------------------------
const std::vector<std::vector<unsigned int>>& RigGriddedPart3d::elementIndices() const
{
return m_elementIndices;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
const std::map<RigGriddedPart3d::BorderSurface, std::vector<unsigned int>>& RigGriddedPart3d::borderSurfaceElements() const
{
return m_borderSurfaceElements;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
const std::vector<std::vector<cvf::Vec3d>>& RigGriddedPart3d::meshLines() const
{
return m_meshLines;
}

View File

@@ -0,0 +1,65 @@
/////////////////////////////////////////////////////////////////////////////////
//
// Copyright (C) 2023 Equinor ASA
//
// ResInsight is free software: you can redistribute it and/or modify
// it under the terms of the GNU General Public License as published by
// the Free Software Foundation, either version 3 of the License, or
// (at your option) any later version.
//
// ResInsight is distributed in the hope that it will be useful, but WITHOUT ANY
// WARRANTY; without even the implied warranty of MERCHANTABILITY or
// FITNESS FOR A PARTICULAR PURPOSE.
//
// See the GNU General Public License at <http://www.gnu.org/licenses/gpl.html>
// for more details.
//
/////////////////////////////////////////////////////////////////////////////////
#pragma once
#include "cvfObject.h"
#include "cvfVector3.h"
#include <map>
#include <vector>
//==================================================================================================
///
///
//==================================================================================================
class RigGriddedPart3d : public cvf::Object
{
public:
enum class BorderSurface
{
UpperSurface = 0,
FaultSurface,
LowerSurface
};
public:
RigGriddedPart3d();
~RigGriddedPart3d() override;
void reset();
void generateGeometry( std::vector<cvf::Vec3d> inputPoints, int nHorzCells, int nVertCellsLower, int nVertCellsMiddle, int nVertCellsUpper );
const std::vector<cvf::Vec3d>& vertices() const;
const std::vector<std::vector<unsigned int>>& elementIndices() const;
const std::map<BorderSurface, std::vector<unsigned int>>& borderSurfaceElements() const;
const std::vector<std::vector<cvf::Vec3d>>& meshLines() const;
protected:
cvf::Vec3d stepVector( cvf::Vec3d start, cvf::Vec3d stop, int nSteps );
void generateMeshlines( std::vector<cvf::Vec3d> cornerPoints, int numHorzCells, int numVertCells );
private:
std::vector<cvf::Vec3d> m_vertices;
std::vector<std::vector<unsigned int>> m_elementIndices;
std::map<BorderSurface, std::vector<unsigned int>> m_borderSurfaceElements;
std::vector<std::vector<cvf::Vec3d>> m_meshLines;
double m_thickness;
};