#3695 Well path creation : Crash when creating a vertical well

Caused by two identical reference points next to each other
This commit is contained in:
Magne Sjaastad 2018-11-16 08:43:19 +01:00
parent e2cc7fb423
commit 9cd695278e
3 changed files with 137 additions and 74 deletions

View File

@ -1,32 +1,31 @@
/////////////////////////////////////////////////////////////////////////////////
//
// Copyright (C) 2018- Statoil 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>
//
// See the GNU General Public License at <http://www.gnu.org/licenses/gpl.html>
// for more details.
//
/////////////////////////////////////////////////////////////////////////////////
#include "RiaPolyArcLineSampler.h"
#include "cvfGeometryTools.h"
#include "cvfMatrix4.h"
#include "RiaArcCurveCalculator.h"
#include "cvfGeometryTools.h"
#include "cvfMatrix4.h"
//--------------------------------------------------------------------------------------------------
///
///
//--------------------------------------------------------------------------------------------------
RiaPolyArcLineSampler::RiaPolyArcLineSampler(const cvf::Vec3d& startTangent,
const std::vector<cvf::Vec3d>& lineArcEndPoints)
RiaPolyArcLineSampler::RiaPolyArcLineSampler(const cvf::Vec3d& startTangent, const std::vector<cvf::Vec3d>& lineArcEndPoints)
: m_startTangent(startTangent)
, m_lineArcEndPoints(lineArcEndPoints)
, m_samplingsInterval(0.15)
@ -35,17 +34,16 @@ RiaPolyArcLineSampler::RiaPolyArcLineSampler(const cvf::Vec3d& startTangent,
, m_points(nullptr)
, m_meshDs(nullptr)
{
}
//--------------------------------------------------------------------------------------------------
///
///
//--------------------------------------------------------------------------------------------------
void RiaPolyArcLineSampler::sampledPointsAndMDs(double sampleInterval,
bool isResamplingLines,
void RiaPolyArcLineSampler::sampledPointsAndMDs(double sampleInterval,
bool isResamplingLines,
std::vector<cvf::Vec3d>* points,
std::vector<double>* mds)
std::vector<double>* mds)
{
CVF_ASSERT(sampleInterval > 0.0);
@ -56,41 +54,43 @@ void RiaPolyArcLineSampler::sampledPointsAndMDs(double sampleInterval,
points->clear();
mds->clear();
if (m_lineArcEndPoints.size() < 2) return ;
std::vector<cvf::Vec3d> pointsNoDuplicates = RiaPolyArcLineSampler::pointsWithoutDuplicates(m_lineArcEndPoints);
m_points = points;
m_meshDs = mds;
if (pointsNoDuplicates.size() < 2) return;
m_points = points;
m_meshDs = mds;
m_totalMD = startMD;
cvf::Vec3d p1 = m_lineArcEndPoints[0];
cvf::Vec3d p2 = m_lineArcEndPoints[1];
cvf::Vec3d p1 = pointsNoDuplicates[0];
cvf::Vec3d p2 = pointsNoDuplicates[1];
m_points->push_back(p1);
m_meshDs->push_back(m_totalMD);
cvf::Vec3d t2 = m_startTangent;
cvf::Vec3d t2 = m_startTangent;
for (size_t pIdx = 0; pIdx < m_lineArcEndPoints.size() - 1 ; ++pIdx)
for (size_t pIdx = 0; pIdx < pointsNoDuplicates.size() - 1; ++pIdx)
{
sampleSegment(t2, m_lineArcEndPoints[pIdx], m_lineArcEndPoints[pIdx + 1] , &t2);
sampleSegment(t2, pointsNoDuplicates[pIdx], pointsNoDuplicates[pIdx + 1], &t2);
}
return ;
return;
}
//--------------------------------------------------------------------------------------------------
///
///
//--------------------------------------------------------------------------------------------------
void RiaPolyArcLineSampler::sampleSegment(cvf::Vec3d t1, cvf::Vec3d p1, cvf::Vec3d p2, cvf::Vec3d* endTangent)
{
cvf::Vec3d p1p2 = p2 - p1;
CVF_ASSERT (p1p2.lengthSquared() > 1e-20);
CVF_ASSERT(p1p2.lengthSquared() > 1e-20);
if (cvf::GeometryTools::getAngle(t1, p1p2) < 1e-5)
{
sampleLine(p1, p2, endTangent);
}
}
else // resample arc
{
sampleArc(t1, p1, p2, endTangent);
@ -98,18 +98,39 @@ void RiaPolyArcLineSampler::sampleSegment(cvf::Vec3d t1, cvf::Vec3d p1, cvf::Vec
}
//--------------------------------------------------------------------------------------------------
///
///
//--------------------------------------------------------------------------------------------------
void RiaPolyArcLineSampler::sampleLine(cvf::Vec3d p1, cvf::Vec3d p2, cvf::Vec3d* endTangent )
std::vector<cvf::Vec3d> RiaPolyArcLineSampler::pointsWithoutDuplicates(const std::vector<cvf::Vec3d>& points)
{
std::vector<cvf::Vec3d> outputPoints;
cvf::Vec3d previousPoint = cvf::Vec3d::UNDEFINED;
const double threshold = 1e-6;
for (const auto& p : points)
{
if (previousPoint.isUndefined() || ((previousPoint - p).lengthSquared()) > threshold)
{
outputPoints.push_back(p);
previousPoint = p;
}
}
return outputPoints;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
void RiaPolyArcLineSampler::sampleLine(cvf::Vec3d p1, cvf::Vec3d p2, cvf::Vec3d* endTangent)
{
cvf::Vec3d p1p2 = p2 - p1;
double p1p2Length = p1p2.length();
if ( p1p2Length > m_samplingsInterval && m_isResamplingLines )
if (p1p2Length > m_samplingsInterval && m_isResamplingLines)
{
cvf::Vec3d tp1p2 = p1p2 / p1p2Length;
double mdInc = m_samplingsInterval;
while ( mdInc < p1p2Length )
double mdInc = m_samplingsInterval;
while (mdInc < p1p2Length)
{
cvf::Vec3d ps = p1 + mdInc * tp1p2;
m_points->push_back(ps);
@ -119,37 +140,37 @@ void RiaPolyArcLineSampler::sampleLine(cvf::Vec3d p1, cvf::Vec3d p2, cvf::Vec3d*
}
m_totalMD += p1p2Length;
m_points->push_back(p2);
m_meshDs->push_back(m_totalMD);
m_meshDs->push_back(m_totalMD);
(*endTangent) = p1p2.getNormalized();
}
//--------------------------------------------------------------------------------------------------
///
///
//--------------------------------------------------------------------------------------------------
void RiaPolyArcLineSampler::sampleArc(cvf::Vec3d t1, cvf::Vec3d p1, cvf::Vec3d p2, cvf::Vec3d* endTangent)
{
// Find arc CS
RiaArcCurveCalculator CS_rad(p1, t1, p2);
double radius = CS_rad.radius();
cvf::Mat4d arcCS = CS_rad.arcCS();
double radius = CS_rad.radius();
cvf::Mat4d arcCS = CS_rad.arcCS();
double angleInc = m_samplingsInterval/ radius;
double angleInc = m_samplingsInterval / radius;
cvf::Vec3d C = CS_rad.center();
cvf::Vec3d N = CS_rad.normal();
// Sample arc by
// Sample arc by
// Rotate vector an increment, and transform to arc CS
double arcAngle = cvf::GeometryTools::getAngle(N, p1 - C, p2 - C);
if (arcAngle/angleInc > 5000)
if (arcAngle / angleInc > 5000)
{
angleInc = arcAngle/5000;
angleInc = arcAngle / 5000;
}
for ( double angle = angleInc; angle < arcAngle; angle += angleInc )
for (double angle = angleInc; angle < arcAngle; angle += angleInc)
{
cvf::Vec3d C_to_incP = cvf::Vec3d::X_AXIS;
C_to_incP *= radius;
@ -159,9 +180,8 @@ void RiaPolyArcLineSampler::sampleArc(cvf::Vec3d t1, cvf::Vec3d p1, cvf::Vec3d p
m_points->push_back(C_to_incP);
m_meshDs->push_back(m_totalMD + angle * radius);
}
m_totalMD += arcAngle*radius;
m_totalMD += arcAngle * radius;
m_points->push_back(p2);
m_meshDs->push_back(m_totalMD);

View File

@ -1,17 +1,17 @@
/////////////////////////////////////////////////////////////////////////////////
//
// Copyright (C) 2018- Statoil 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>
//
// See the GNU General Public License at <http://www.gnu.org/licenses/gpl.html>
// for more details.
//
/////////////////////////////////////////////////////////////////////////////////
@ -25,14 +25,14 @@
class RiaPolyArcLineSampler
{
public:
RiaPolyArcLineSampler(const cvf::Vec3d& startTangent,
const std::vector<cvf::Vec3d>& lineArcEndPoints);
RiaPolyArcLineSampler(const cvf::Vec3d& startTangent, const std::vector<cvf::Vec3d>& lineArcEndPoints);
void sampledPointsAndMDs(double sampleInterval,
bool isResamplingLines,
std::vector<cvf::Vec3d>* points,
std::vector<double>* mds);
void sampledPointsAndMDs(double sampleInterval,
bool isResamplingLines,
std::vector<cvf::Vec3d>* points,
std::vector<double>* mds);
static std::vector<cvf::Vec3d> pointsWithoutDuplicates(const std::vector<cvf::Vec3d>& points);
private:
void sampleLine(cvf::Vec3d p1, cvf::Vec3d p2, cvf::Vec3d* endTangent);
@ -42,11 +42,10 @@ private:
std::vector<cvf::Vec3d>* m_points; // Internal temporary pointers to collections beeing filled.
std::vector<double>* m_meshDs;
double m_samplingsInterval;
bool m_isResamplingLines;
double m_totalMD;
cvf::Vec3d m_startTangent;
std::vector<cvf::Vec3d> m_lineArcEndPoints;
};
double m_samplingsInterval;
bool m_isResamplingLines;
double m_totalMD;
cvf::Vec3d m_startTangent;
std::vector<cvf::Vec3d> m_lineArcEndPoints;
};

View File

@ -6,28 +6,72 @@
//--------------------------------------------------------------------------------------------------
TEST(RiaPolyArcLineSampler, Basic)
{
std::vector<cvf::Vec3d> points { {0,0,0}, {0,0,-10}, {0,10,-20},{0,20,-20}, {0,30, -30} };
std::vector<cvf::Vec3d> points{{0, 0, 0}, {0, 0, -10}, {0, 10, -20}, {0, 20, -20}, {0, 30, -30}};
RiaPolyArcLineSampler sampler({0,0,-1}, points);
RiaPolyArcLineSampler sampler({0, 0, -1}, points);
std::vector<cvf::Vec3d> sampledPoints;
std::vector<double> mds;
std::vector<cvf::Vec3d> sampledPoints;
std::vector<double> mds;
sampler.sampledPointsAndMDs(2, true, &sampledPoints, &mds);
#if 0
#if 0
for (size_t pIdx = 0; pIdx < sampledPoints.size(); ++pIdx)
{
std::cout << sampledPoints[pIdx].x() << " "
<< sampledPoints[pIdx].y() << " "
<< sampledPoints[pIdx].z() << " md: " << mds[pIdx] << std::endl;
}
#endif
}
#endif
EXPECT_EQ(27, (int)sampledPoints.size());
EXPECT_NEAR(51.4159, mds.back(), 1e-4);
EXPECT_NEAR(51.4159, mds.back(), 1e-4);
EXPECT_NEAR(points[2].y(), sampledPoints[12].y(), 2);
EXPECT_NEAR(points[2].z(), sampledPoints[12].z(), 2);
EXPECT_NEAR(points[2].y(), sampledPoints[12].y(), 2);
EXPECT_NEAR(points[2].z(), sampledPoints[12].z(), 2);
EXPECT_NEAR(points[4].y(), sampledPoints[25].y(), 2);
EXPECT_NEAR(points[4].z(), sampledPoints[25].z(), 2);
EXPECT_NEAR(points[4].y(), sampledPoints[25].y(), 2);
EXPECT_NEAR(points[4].z(), sampledPoints[25].z(), 2);
}
//--------------------------------------------------------------------------------------------------
TEST(RiaPolyArcLineSampler, TestInvalidInput)
{
{
// Two identical points after each other
std::vector<cvf::Vec3d> points{{0, 0, -20}, {0, 0, -20}};
RiaPolyArcLineSampler sampler({0, 0, -1}, points);
std::vector<cvf::Vec3d> sampledPoints;
std::vector<double> mds;
sampler.sampledPointsAndMDs(2, true, &sampledPoints, &mds);
EXPECT_EQ(0, (int)sampledPoints.size());
}
{
std::vector<cvf::Vec3d> points;
RiaPolyArcLineSampler sampler({0, 0, -1}, points);
std::vector<cvf::Vec3d> sampledPoints;
std::vector<double> mds;
sampler.sampledPointsAndMDs(2, true, &sampledPoints, &mds);
EXPECT_EQ(0, (int)sampledPoints.size());
}
{
std::vector<cvf::Vec3d> points{{0, 0, 0}};
RiaPolyArcLineSampler sampler({0, 0, -1}, points);
std::vector<cvf::Vec3d> sampledPoints;
std::vector<double> mds;
sampler.sampledPointsAndMDs(2, true, &sampledPoints, &mds);
EXPECT_EQ(0, (int)sampledPoints.size());
}
}