平面に埋め込まれた平面グラフ (平面グラフ) があり、その面を検索したいと考えています。グラフは接続されていませんが、個別にアドレス指定できないいくつかの接続されたグラフで構成されています (たとえば、サブグラフは別のグラフの面に含まれる可能性があります) 特定の 2 次元点を含む多角形 (面) を見つけたいです。多角形は、グラフの面によって形成されます。顔の数がかなり多いので、事前に決定することは避けたいと思います。そのような検索の一般的な複雑さと、それを達成するために使用できる C++ ライブラリ/コーディング アプローチは何ですか。
明確にするために更新: ここでは、xy 平面のグラフを参照しています。
あなたは興味深い挑戦をします。多角形が常に凸面である場合、比較的単純な解決策が可能です。その場合、関心のある点が多角形のすべての側面の同じ側面 (左または右) にあるかどうかを尋ねるだけでよいからです。一般的なケースで特に単純な解決策を私は知りませんが、コーシーの有名な積分公式に間接的に触発されたように、次のコードは任意の多角形で機能するようです。
コードを理解するために Cauchy に精通している必要はありません。コード内のコメントでテクニックが説明されています。
#include <vector>
#include <cstddef>
#include <cstdlib>
#include <cmath>
#include <iostream>
// This program takes its data from the standard input
// stream like this:
//
// 1.2 0.5
// -0.1 -0.2
// 2.7 -0.3
// 2.5 2.9
// 0.1 2.8
//
// Given such input, the program answers whether the
// point (1.2, 0.5) does not lie within the polygon
// whose vertices, in sequence, are (-0.1, -0.2),
// (2.7, -0.3), (2.5, 2.9) and (0.1, 2.8). Naturally,
// the program wants at least three vertices, so it
// requires the input of at least eight numbers (where
// the example has four vertices and thus ten numbers).
//
// This code lacks really robust error handling, which
// could however be added without too much trouble.
// Also, its function angle_swept() could be shortened
// at cost to readability; but this is not done here,
// since the function is already hard enough to grasp as
// it stands.
//
//
const double TWOPI = 8.0 * atan2(1.0, 1.0); // two times pi, or 360 deg
namespace {
struct Point {
double x;
double y;
Point(const double x0 = 0.0, const double y0 = 0.0)
: x(x0), y(y0) {}
};
// As it happens, for the present code's purpose,
// a Point and a Vector want exactly the same
// members and operations; thus, make the one a
// synonym for the other.
typedef Point Vector;
std::istream &operator>>(std::istream &ist, Point &point) {
double x1, y1;
if(ist >> x1 >> y1) point = Point(x1, y1);
return ist;
}
// Calculate the vector from one point to another.
Vector operator-(const Point &point2, const Point &point1) {
return Vector(point2.x - point1.x, point2.y - point1.y);
}
// Calculate the dot product of two Vectors.
// Overload the "*" operator for this purpose.
double operator*(const Vector &vector1, const Vector &vector2) {
return vector1.x*vector2.x + vector1.y*vector2.y;
}
// Calculate the (two-dimensional) cross product of two Vectors.
// Overload the "%" operator for this purpose.
double operator%(const Vector &vector1, const Vector &vector2) {
return vector1.x*vector2.y - vector1.y*vector2.x;
}
// Calculate a Vector's magnitude or length.
double abs(const Vector &vector) {
return std::sqrt(vector.x*vector.x + vector.y*vector.y);
}
// Normalize a vector to unit length.
Vector unit(const Vector &vector) {
const double abs1 = abs(vector);
return Vector(vector.x/abs1, vector.y/abs1);
}
// Imagine standing in the plane at the point of
// interest, facing toward a vertex. Then imagine
// turning to face the next vertex without leaving
// the point. Answer this question: through what
// angle did you just turn, measured in radians?
double angle_swept(
const Point &point, const Point &vertex1, const Point &vertex2
) {
const Vector unit1 = unit(vertex1 - point);
const Vector unit2 = unit(vertex2 - point);
const double dot_product = unit1 * unit2;
const double cross_product = unit1 % unit2;
// (Here we must be careful. Either the dot
// product or the cross product could in theory
// be used to extract the angle but, in
// practice, either the one or the other may be
// numerically problematical. Use whichever
// delivers the better accuracy.)
return (fabs(dot_product) <= fabs(cross_product)) ? (
(cross_product >= 0.0) ? (
// The angle lies between 45 and 135 degrees.
acos(dot_product)
) : (
// The angle lies between -45 and -135 degrees.
-acos(dot_product)
)
) : (
(dot_product >= 0.0) ? (
// The angle lies between -45 and 45 degrees.
asin(cross_product)
) : (
// The angle lies between 135 and 180 degrees
// or between -135 and -180 degrees.
((cross_product >= 0.0) ? TWOPI/2.0 : -TWOPI/2.0)
- asin(cross_product)
)
);
}
}
int main(const int, char **const argv) {
// Read the x and y coordinates of the point of
// interest, followed by the x and y coordinates of
// each vertex in sequence, from std. input.
// Observe that whether the sequence of vertices
// runs clockwise or counterclockwise does
// not matter.
Point point;
std::vector<Point> vertex;
std::cin >> point;
{
Point point1;
while (std::cin >> point1) vertex.push_back(point1);
}
if (vertex.size() < 3) {
std::cerr << argv[0]
<< ": a polygon wants at least three vertices\n";
std::exit(1);
}
// Standing as it were at the point of interest,
// turn to face each vertex in sequence. Keep
// track of the total angle through which you
// have turned.
double cumulative_angle_swept = 0.0;
for (size_t i = 0; i < vertex.size(); ++i) {
// In an N-sided polygon, vertex N is again
// vertex 0. Since j==N is out of range,
// if i==N-1, then let j=0. Otherwise,
// let j=i+1.
const size_t j = (i+1) % vertex.size();
cumulative_angle_swept +=
angle_swept(point, vertex[i], vertex[j]);
}
// Judge the point of interest to lie within the
// polygon if you have turned a complete circuit.
const bool does_the_point_lie_within_the_polygon =
fabs(cumulative_angle_swept) >= TWOPI/2.0;
// Output.
std::cout
<< "The angle swept by the polygon's vertices about the point\n"
<< "of interest is " << cumulative_angle_swept << " radians ("
<< ((360.0/TWOPI)*cumulative_angle_swept) << " degrees).\n"
<< "Therefore, the point lies "
<< (
does_the_point_lie_within_the_polygon
? "within" : "outside of"
)
<< " the polygon.\n";
return !does_the_point_lie_within_the_polygon;
}
もちろん、上記のコードは私が書いたものです。なぜなら、あなたの挑戦は興味深く、私はそれを達成できるかどうかを知りたかったからです。アプリケーションが重要な場合は、コードのテストとレビューの両方を行う必要があります。見つかったバグがあれば、ここで修正してください。私は 2 つまたは 3 つのケースに対してコードをテストしましたが、動作しているように見えますが、重要な任務のためには、より徹底的なテストが必要です。
あなたのアプリケーションで頑張ってください。