Геометрия — 1

Вычислительная геометрия — раздел информатики, изучающий алгоритмы для решения геометрических задач. В этой и следующей лекциях мы познакомимся с основными понятиями и научимся работать с точками, векторами, прямыми и отрезками на плоскости.

Основные понятия

Плоскость задаётся двумя перпендикулярными осями координат \(Ox\) и \(Oy\). Каждая точка на плоскости однозначно определяется парой чисел \((x, y)\) — её координатами.

Вектор — направленный отрезок, то есть отрезок, для которого указано, какой конец является началом, а какой — концом. Вектор на плоскости задаётся двумя числами — координатами по горизонтали и вертикали.

Вектор \(\overrightarrow{AB}\)

Отрезок — часть прямой, ограниченная двумя точками (концами отрезка). В отличие от вектора, он неориентирован.

Луч — часть прямой, ограниченная с одной стороны точкой (началом луча) и бесконечная в другом направлении.

Точка \(\simeq\) вектор

Поскольку и точка, и вектор задаются парой чисел, удобно считать их одним объектом. Каждой точке можно сопоставить её радиус-вектор — вектор из начала координат в эту точку. В дальнейшем мы будем хранить и точки, и векторы одной структурой.

Напоминание: линейная функция

Функция \(y = kx + b\) задаёт прямую на плоскости. Коэффициент \(k\) — это угловой коэффициент (наклон прямой). Он равен тангенсу угла наклона прямой к оси \(Ox\):

\[ k = \tan \alpha = \frac{\Delta y}{\Delta x} \]

Это понятие нам пригодится при обсуждении способов задания прямой.

Хранение точек и векторов в C++

Создадим структуру для представления точки (и вектора) на плоскости:

struct Point {
    int x = 0;
    int y = 0;

    Point() = default;
    Point(int x, int y) : x(x), y(y) {}
};

Функция с тем же именем, что и структура, вызывается при создании объекта — это конструктор. Здесь Point() создаёт точку \((0, 0)\), а Point(3, 5) — точку с координатами \((3, 5)\).

Мы выбрали целочисленный тип int для хранения координат. Если все входные координаты целые, то многие операции (скалярное и векторное произведения, проверки на коллинеарность) тоже можно проводить в целых числах, избегая проблем с точностью чисел с плавающей точкой.

Операторы

В C++ можно перегружать стандартные операторы, например +, -, *. Определим сложение и вычитание векторов, а также умножение вектора на скаляр:

Point operator+(Point a, Point b) { return {a.x + b.x, a.y + b.y}; }
Point operator-(Point a, Point b) { return {a.x - b.x, a.y - b.y}; }
Point operator*(int k, Point a)   { return {k * a.x, k * a.y}; }
Point operator*(Point a, int k)   { return {k * a.x, k * a.y}; }

Теперь можно писать Point c = a + b; или Point d = a - b;, что делает код более читаемым.

Длина вектора

Длина (модуль) вектора вычисляется по теореме Пифагора: \(|\vec{a}| = \sqrt{x^2 + y^2}\).

double length(Point a) {
    return std::hypot(a.x, a.y);
}

Функция std::hypot из <cmath> вычисляет \(\sqrt{x^2 + y^2}\) точнее и безопаснее, чем std::sqrt(a.x * a.x + a.y * a.y), избегая промежуточного переполнения.

Угол вектора

Для вычисления угла вектора относительно оси \(Ox\) удобно использовать функцию atan2:

Тригонометрический круг
double angle(Point a) {
    return std::atan2(a.y, a.x);
}

Функция std::atan2(y, x) возвращает угол в радианах из промежутка \([-\pi, +\pi]\). Она работает корректнее, чем std::atan(y / x), поскольку учитывает знаки обоих аргументов и не требует деления (что было бы проблемой при \(x = 0\)).

Для перевода в градусы: \(\alpha_{\circ} = \alpha_{\text{рад}} \cdot \frac{180}{\pi}\).

Помимо сложения, вычитания и умножения на скаляр, у векторов есть две особенно важные операции — скалярное и векторное произведения. Они являются основным инструментом вычислительной геометрии.

Скалярное произведение

Скалярное произведение (англ. dot product) двух векторов равно произведению их длин и косинуса угла между ними:

\[ \vec{a} \cdot \vec{b} = |\vec{a}| \cdot |\vec{b}| \cdot \cos \theta \]

Его можно вычислить через координаты:

\[ \vec{a} \cdot \vec{b} = x_a x_b + y_a y_b \]

Геометрически скалярное произведение равно длине проекции вектора \(\vec{b}\) на направление \(\vec{a}\), умноженной на \(|\vec{a}|\):

Скалярное произведение как проекция

Свойства:

  • Скалярное произведение симметрично: \(\vec{a} \cdot \vec{b} = \vec{b} \cdot \vec{a}\).

  • Если угол между векторами острый, произведение положительное.

  • Если угол тупой — отрицательное.

  • Перпендикулярные вектора имеют нулевое скалярное произведение (так как \(\cos 90° = 0\)).

int dot(Point a, Point b) {
    return a.x * b.x + a.y * b.y;
}

Векторное произведение

Векторное произведение (англ. cross product, также называется косым или псевдоскалярным) для двух векторов на плоскости равно произведению их длин и синуса угла между ними. Знак синуса зависит от порядка операндов:

\[ \vec{a} \times \vec{b} = |\vec{a}| \cdot |\vec{b}| \cdot \sin \theta \]

Координатная формула:

\[ \vec{a} \times \vec{b} = x_a y_b - y_a x_b \]

Геометрически, модуль векторного произведения — это площадь параллелограмма, натянутого на вектора \(\vec{a}\) и \(\vec{b}\):

Площадь параллелограмма

Свойства:

  • Векторное произведение антисимметрично: \(\vec{a} \times \vec{b} = -(\vec{b} \times \vec{a})\).

  • Коллинеарные вектора имеют нулевое векторное произведение (так как \(\sin 0° = 0\)).

  • Если \(\vec{b}\) находится «слева» от \(\vec{a}\) (поворот против часовой стрелки), произведение положительное.

  • Если «справа» (по часовой стрелке) — отрицательное.

int cross(Point a, Point b) {
    return a.x * b.y - a.y * b.x;
}
Note

Формально, в пространстве векторное произведение определяется как вектор, перпендикулярный обоим исходным, длиной \(|\vec{a}| \cdot |\vec{b}| \cdot \sin \theta\). Это определение работает в трёхмерной геометрии и требует знания матриц. На плоскости мы используем его скалярный аналог — проекцию результирующего вектора на ось \(Oz\).

Tip

Также часто перегружают operator* для скалярного произведения и operator% (или operator^) для векторного, чтобы записывать выражения компактнее. Мы используем именованные функции dot и cross для ясности, но при необходимости можно определить:

int operator*(Point a, Point b) { return a.x * b.x + a.y * b.y; }
int operator%(Point a, Point b) { return a.x * b.y - a.y * b.x; }

Применение произведений

Скалярное и векторное произведения тесно связаны с углами и позволяют вычислять ориентированные углы и площади.

Например, угол между двумя векторами можно найти, передав в atan2 векторное произведение (пропорционально синусу) и скалярное (пропорционально косинусу):

double angleBetween(Point a, Point b) {
    return std::atan2(cross(a, b), dot(a, b));
}

В дальнейшем свойства произведений помогут нам определять взаимное расположение точек, прямых и отрезков, поэтому важно их хорошо понять и запомнить.

Прямые на плоскости

Способ 1: общее уравнение прямой

Любую прямую на плоскости можно задать уравнением:

\[ Ax + By + C = 0 \]

где \(A\), \(B\), \(C\) — коэффициенты, причём \(A\) и \(B\) не равны нулю одновременно.

При этом:

  • Вектор \((A, B)\) является нормалью прямой — он ей перпендикулярен.

  • Вектор \((-B, A)\) является направляющим вектором — он параллелен прямой.

Кроме того, неравенство \(Ax + By + C > 0\) задаёт полуплоскость — множество точек по одну сторону от прямой.

Способ 2: точка и направляющий вектор

Прямую можно задать точкой \(P_0\), лежащей на ней, и направляющим вектором \(\vec{d}\). Тогда любая точка \(P\) на прямой выражается как:

\[ P = P_0 + t \cdot \vec{d}, \quad t \in \mathbb{R} \]

В координатной форме:

\[ \begin{cases} x = x_0 + t \cdot d_x \\ y = y_0 + t \cdot d_y \end{cases} \]

Способ 3: две точки

Если даны две различные точки \(A\) и \(B\), то прямая через них определяется однозначно. Направляющий вектор — \(\vec{d} = B - A\), а за начальную точку берём \(A\):

\[ P = A + t \cdot (B - A), \quad t \in \mathbb{R} \]

Тем самым задача сводится к способу 2.

Вывод знакомого уравнения. Из параметрического представления можно получить уравнение \(y = kx + b\). Если \(B_x \ne A_x\), выразим \(t\) из первого координатного уравнения:

\[ t = \frac{x - A_x}{B_x - A_x} \]

Подставив во второе:

\[ y = A_y + \frac{B_y - A_y}{B_x - A_x} \cdot (x - A_x) \]

Здесь \(k = \frac{B_y - A_y}{B_x - A_x}\) — это угловой коэффициент, тот самый \(\tan \alpha\).

Warning

Представление \(y = kx + b\) не работает для вертикальных прямых (\(B_x = A_x\)). Поэтому в программировании обычно используют общее уравнение \(Ax + By + C = 0\) или параметрическую форму.

Переходы между представлениями

На практике прямая часто задана в одном виде, а работать удобнее с другим:

  • Две точки \(\to\) уравнение. По двум точкам \(A\) и \(B\) коэффициенты уравнения: \(A_l = A_y - B_y\), \(\ B_l = B_x - A_x\), \(\ C_l = -(A_l \cdot A_x + B_l \cdot A_y)\).

  • Уравнение \(\to\) направляющий вектор. \(\vec{d} = (-B, A)\).

  • Уравнение \(\to\) нормаль. \(\vec{n} = (A, B)\).

  • Уравнение \(\to\) две точки. Достаточно подобрать любые две точки, удовлетворяющие уравнению:

Point p1, p2;
if (A == 0) {
    p1 = {0, -C / B};
    p2 = {1, -C / B};
} else {
    p1 = {-C / A, 0};
    p2 = {-(C + B) / A, 1};
}

Расстояния

Расстояние между двумя точками

Расстояние между точками \(A\) и \(B\) — длина вектора \(\overrightarrow{AB}\):

\[ \rho(A, B) = |\overrightarrow{AB}| = \sqrt{(B_x - A_x)^2 + (B_y - A_y)^2} \]

double dist(Point a, Point b) {
    return length(b - a);
}

Расстояние от точки до прямой

Через уравнение прямой. Если прямая задана как \(Ax + By + C = 0\), а точка имеет координаты \((x_0, y_0)\):

\[ \rho = \frac{|Ax_0 + By_0 + C|}{\sqrt{A^2 + B^2}} \]

Об этой формуле можно думать как о скалярном произведении вектора-точки на нормированный (\(\frac{1}{\sqrt{A^2+B^2}}\)) вектор нормали, то есть как о длине проекции точки на нормаль прямой.

Выражение \(Ax_0 + By_0 + C\) (без модуля) определяет знаковое расстояние: его знак показывает, по какую сторону от прямой находится точка.

Через векторное произведение. Если прямая задана двумя точками \(A\) и \(B\), можно выразить расстояние через площадь треугольника. Площадь треугольника \(\triangle PAB\) равна \(S = \frac{1}{2} \cdot |\overrightarrow{AB}| \cdot h\), откуда:

\[ h = \frac{2S}{|\overrightarrow{AB}|} = \frac{|\overrightarrow{PA} \times \overrightarrow{PB}|}{|\overrightarrow{AB}|} \]

В числителе стоит модуль векторного произведения — он равен удвоенной площади треугольника \(\triangle PAB\).

double distToLine(Point p, Point a, Point b) {
    return std::abs(cross(a - p, b - p)) / length(b - a);
}

Расстояние от точки до луча

Луч начинается в точке \(A\) и направлен в сторону \(B\). Расстояние от точки \(P\) до луча зависит от того, куда проецируется \(P\) на прямую \(AB\):

  • Если проекция попадает на луч (вперёд от \(A\)), расстояние равно расстоянию от \(P\) до прямой \(AB\).

  • Если проекция попадает за начало луча, ближайшая точка луча — \(A\), и расстояние равно \(\rho(P, A)\).

Скалярное произведение позволяет различить эти случаи: \(\overrightarrow{AP} \cdot \overrightarrow{AB} \geq 0\) означает, что проекция лежит на луче.

double distToRay(Point p, Point a, Point b) {
    if (dot(p - a, b - a) < 0) {
        return length(p - a);
    }
    return distToLine(p, a, b);
}

Расстояние от точки до отрезка

Для отрезка нужно проверить обе границы. Если проекция точки лежит между концами отрезка — берём расстояние до прямой, иначе — до ближайшего конца:

double distToSegment(Point p, Point a, Point b) {
    if (dot(p - a, b - a) < 0) {
        return length(p - a);
    }
    if (dot(p - b, a - b) < 0) {
        return length(p - b);
    }
    return distToLine(p, a, b);
}

Принадлежность точки геометрическим объектам

Точка на прямой

Точка \(P\) лежит на прямой через \(A\) и \(B\) тогда и только тогда, когда вектора \(\overrightarrow{AP}\) и \(\overrightarrow{AB}\) коллинеарны:

\[ \overrightarrow{AP} \times \overrightarrow{AB} = 0 \]

bool onLine(Point p, Point a, Point b) {
    return cross(p - a, b - a) == 0;
}

Точка на луче

Точка \(P\) лежит на луче из \(A\) в направлении \(B\), если она лежит на прямой \(AB\) и при этом расположена в том же направлении от \(A\), что и \(B\):

\[ \overrightarrow{AP} \times \overrightarrow{AB} = 0 \quad \text{и} \quad \overrightarrow{AP} \cdot \overrightarrow{AB} \geq 0 \]

bool onRay(Point p, Point a, Point b) {
    return cross(p - a, b - a) == 0
        && dot(p - a, b - a) >= 0;
}

Точка на отрезке

Точка \(P\) лежит на отрезке \(AB\), если она лежит на прямой \(AB\) и находится между \(A\) и \(B\):

\[ \overrightarrow{AP} \times \overrightarrow{AB} = 0 \quad \text{и} \quad 0 \leq \overrightarrow{AP} \cdot \overrightarrow{AB} \leq \overrightarrow{AB} \cdot \overrightarrow{AB} \]

Второе условие означает, что проекция \(P\) на направление \(\overrightarrow{AB}\) лежит между \(0\) и \(|\overrightarrow{AB}|^2\).

bool onSegment(Point p, Point a, Point b) {
    return cross(p - a, b - a) == 0
        && dot(p - a, b - a) >= 0
        && dot(p - a, b - a) <= dot(b - a, b - a);
}
Note

Проверка принадлежности точки выпуклому многоугольнику — тема следующей лекции.

Пересечения

Пересечение двух прямых

Найти точку пересечения двух прямых — значит найти точку, удовлетворяющую обоим уравнениям:

\[ \begin{cases} A_1 x + B_1 y + C_1 = 0 \\ A_2 x + B_2 y + C_2 = 0 \end{cases} \]

Решая систему, получаем:

\[ x = \frac{B_1 C_2 - B_2 C_1}{A_1 B_2 - A_2 B_1}, \quad y = \frac{A_2 C_1 - A_1 C_2}{A_1 B_2 - A_2 B_1} \]

Знаменатель \(D = A_1 B_2 - A_2 B_1\) — это векторное произведение нормалей двух прямых. Если \(D = 0\), прямые параллельны (или совпадают), и точки пересечения в обычном смысле нет. Этот случай нужно обрабатывать отдельно.

struct Line {
    int a, b, c;
};

bool intersectLines(Line l1, Line l2, double& x, double& y) {
    int d = l1.a * l2.b - l2.a * l1.b;
    if (d == 0) {
        return false;
    }
    x = (double)(l1.b * l2.c - l2.b * l1.c) / d;
    y = (double)(l2.a * l1.c - l1.a * l2.c) / d;
    return true;
}

Пересечение двух отрезков

Для проверки пересечения отрезков используется следующее наблюдение: отрезки \(P_1 P_2\) и \(Q_1 Q_2\) пересекаются, если концы каждого отрезка лежат по разные стороны от прямой, содержащей другой отрезок.

Пересечение отрезков

В терминах векторного произведения, «\(Q_1\) и \(Q_2\) лежат по разные стороны от прямой \(P_1 P_2\)» записывается как:

\[ (\overrightarrow{P_1 P_2} \times \overrightarrow{P_1 Q_1}) \cdot (\overrightarrow{P_1 P_2} \times \overrightarrow{P_1 Q_2}) < 0 \]

Записав аналогичное условие для \(P_1, P_2\) относительно прямой \(Q_1 Q_2\), получаем полный критерий пересечения (в невырожденном случае). Отдельно нужно обработать вырожденные случаи, когда конец одного отрезка лежит на другом:

bool segmentsIntersect(Point p1, Point p2, Point q1, Point q2) {
    long long d1 = cross(p2 - p1, q1 - p1);
    long long d2 = cross(p2 - p1, q2 - p1);
    long long d3 = cross(q2 - q1, p1 - q1);
    long long d4 = cross(q2 - q1, p2 - q1);

    if (d1 * d2 < 0 && d3 * d4 < 0) {
        return true;
    }

    if (d1 == 0 && onSegment(q1, p1, p2)) return true;
    if (d2 == 0 && onSegment(q2, p1, p2)) return true;
    if (d3 == 0 && onSegment(p1, q1, q2)) return true;
    if (d4 == 0 && onSegment(p2, q1, q2)) return true;

    return false;
}

Если отрезки пересекаются и нужно найти саму точку пересечения, достаточно найти пересечение прямых, содержащих эти отрезки (с помощью intersectLines).

Проекция точки на прямую

Проекция точки \(P\) на прямую — это ближайшая к \(P\) точка прямой. Если прямая задана точкой \(A\) и направляющим вектором \(\vec{d}\), то проекцию можно выразить через скалярное произведение:

\[ \text{Пр}_{\vec{d}} \overrightarrow{AP} = \frac{\overrightarrow{AP} \cdot \vec{d}}{|\vec{d}|^2} \cdot \vec{d} \]

Геометрический смысл: находим параметр \(t = \frac{\overrightarrow{AP} \cdot \vec{d}}{|\vec{d}|^2}\), который показывает, какая доля направляющего вектора соответствует проекции, и прибавляем \(t \cdot \vec{d}\) к начальной точке \(A\).

Вместо того чтобы раскрывать формулы покоординатно и получать громоздкий ответ, вычислим всё по частям — так логика алгоритма остаётся наглядной и код проще отлаживать:

std::pair<double, double> project(Point p, Point a, Point d) {
    Point v = p - a;
    double t = (double)dot(v, d) / dot(d, d);
    return {a.x + t * d.x, a.y + t * d.y};
}