Содержание к диссертации
Введение
Глава 1. Классификация точки относительно контура 23
1.1. Введение 23
1.2. Суммирование угла 26
1.3. Несколько контуров 27
1.4. Угол вращения вдоль ребра 28
1.5. Число обусловленности 29
1.6. Ускорение
1.6.1. Предподсчёт 31
1.6.2. Запрос 31
1.6.3. Временная сложность 33
1.6.4. Сравнение 36
1.6.5. Следствия 38
1.6.6. Абсолютный угол поворота 40
1.6.7. Результаты
1.7. NURBS-кривые 44
1.8. Заключение 46
Глава 2 . Вычисление NURBS 48
2.1. Введение 48
2.2. Предварительные сведения 50
2.3. Прямое вычисление: поиск спа на 51
2.4. Базисные функции
2.4.1. Формула Де Бура 52
2.4.2. Степенной базис 52
2.4.3. Смешивание контрольных точек 53
2.5. Подход в 3D 54
2.5.1. Алгоритм Де Бура 54
2.5.2. Кусочно-полиномиальное представление 54
2.6. Вычисление производных 55
2.6.1. Формула Де Бура 56
2.6.2. Многочлены 56
2.6.3. Смешивание контрольных точек з
2.7. Предлагаемый подход 58
2.7.1. Векторизация 58
2.7.2. Поиск спана 60
2.7.3. Базисные функции 61
2.7.4. Смешивание контрольных точек 63
2.7.5. Производные в 3D 64
2.7.6. Последние замечания
2.8. Сравнение производительности 65
2.9. Заключение 69
Глава 3. Метод Ньютона-Рафсона 70
3.1. Введение 70
3.1.1. Форсирование ранга 72
3.2. Теория 73
3.2.1. Предположения 73
3.2.2. Сходимость метода Ньютона 74
3.2.3. Трассировка кривой 81
3.3. Пересечение поверхностей 86
3.3.1. Общее описание 86
3.3.2. Трассировка и сходимость 87
3.4. Экспериментальные данные 89
3.4.1. Метод Ньютона-Рафсона 89
3.4.2. Трассировка кривой 90
3.5. Заключение 93
Глава 4. Пересечение кривых и поверхностей 95
4.1. Постановка задачи 95
4.2. Аналитические случаи 98
4.3. Общая схема алгоритма 99
4.4. Обзор литературы 100
4.5. Глобальный анализ 101
4.6. Локальный анализ 104
4.7. Режимы и системы уравнений
4.7.1. Несингулярная система 105
4.7.2. Сингулярная система 106
4.7.3. Масштабирование 107
4.7.4. Реализация
4.8. Метод Ньютона 110
4.9. Трассировка
4.9.1. Аппроксимация кривой 113
4.9.2. Определение скорости 113
4.9.3. Предиктор 115
4.9.4. Корректор 116
4.9.5. Построение и проверка звена 117
4.9.6. Адаптивность величины шага 119
4.9.7. Экспериментальные данные 119
4.10. Недопущение дубликатов 122
4.10.1. Проверка столкновения 122
4.10.2. Вычисление ширины 123
4.10.3. Определение точки или кривой
4.11. Поиск наложения поверхностей 125
4.12. Якорные точки
4.12.1. Замыкание кривой в контур 129
4.12.2. Завершение в якорной точке 130
4.12.3. Завершение на границе 131
4.13. Обработка границ 133
4.13.1. Пересечения в области определения 133
4.13.2. Метод Ньютона на границе 133
4.13.3. Кривая или точка 134
4.13.4. Трассировка вдоль границы 135
4.13.5. Гладкое сочленение поверхностей 1 4.14. Экспериментальные данные 138
4.15. Заключение 141
Заключение 143
Список литературы
- Число обусловленности
- Прямое вычисление: поиск спа
- Кусочно-полиномиальное представление
- Сходимость метода Ньютона
Число обусловленности
Определение положения точки относительно контура — фундаментальная задача вычислительной геометрии. В ней требуется определить для заданной точки на плоскости, лежит ли она внутри заданного простого многоугольника. Иногда требуется выбрать один из трёх вариантов: внутри, снаружи или на границе. В классической постановке многоугольник определяется последовательностью своих вершин, а все рёбра являются отрезками. В обобщённой постановке многоугольник определяется последовательностью своих рёбер, причём рёбра могут быть криволинейными.
Задача определения положения точки относительно контура, как и более общая задача локализации точки на планарной карте, часто возникает в различных приложениях компьютерной графики и геометрического моделирования, а также в контексте геоинформационных систем, систем автоматизации проектирования и пространственных баз данных. Алгоритмы для решения этой задачи применяются при прямой визуализации обрезанных NURBS поверхностей методом трассировки лучей [46, 47], при отсечении многоугольников и выполнении двумерных булевых операций [48, 49]. При работе с картой в ГИС эти алгоритмы могут определять, например, в какую область попал курсором мыши пользователь.
Следует заметить, что в области САПР имеется своя специфика: что геометрические данные часто задаются неточно. Рёбра контура могут быть немного смещены относительно истинного положения, а также между соседними рёбрами многоугольника могут быть небольшие промежутки. По этой причине важным критерием выбора алгоритма становится его устойчивость к неточностям входных данных и численным погрешностям.
В статье [50] содержится хороший обзор простых методов решения задачи в классической постановке. Другой обзор простых методов решения задачи есть в [51]. Надёжность существующих реализаций, использующих числа с плавающей точкой, изучается в [52].
Метод пересечения луча. На практике для решения задачи чаще всего применяется метод пересечения луча. В этом методе используется тот факт, что любой луч, выпущенный из точки, пересекает многоугольник в нечётном количестве точек тогда и только тогда, когда начало луча лежит внутри многоугольника. Это правило неверно, если луч проходит через вершину1.
Реализация, использующая числа с плавающей точкой, также может вернуть неверный результат, если луч проходит близко к вершине. В классической постановке задачи этой проблемы можно избежать, если выбрать направление луча вдоль оси координат и правильно рассмотреть все случаи. Однако, такое решение проблемы не работает в обобщённой постановке. Другая идея заключается в том, чтобы выпустить несколько лучей в случайных направлениях, а потом установить конечный вердикт методом голосования. Этот вариант суще-1 или касается криволинейного ребра ственно уменьшает вероятность ошибки, однако полностью исключить возможность ошибки не получится.
Если проверять много точек относительно одного многоугольника, то можно ускорить точечные запросы за счёт предподсчёта. Иерархия ограничивающих объёмов (BVH) позволяет быстро отсечь рёбра, которые точно не пересекаются с лучом. Это сильно ускоряет запросы, хотя нетривиальные верхние оценки на временную сложность запроса автору неизвестны.
Метод суммирования угла вращения. Другой хорошо известный подход к задаче определения положения точки относительно контура — суммирование угла вращения. Общая идея такова: поставим начало вектора в проверяемую точку, а конец вектора проведём вдоль контура. После одного оборота нужно определить, сколько вращений вокруг точки совершил вектор. Это число равно нулю в том и только в том случае, когда проверяемая точка лежит снаружи контура.
Проклятие метода суммирования угла — массовое использование обратных тригонометрических функций. По этой причине данный метод часто считают самым медленным и плохим решением задачи2. В качестве компенсации, метод тривиального суммирования угла является обратно устойчивым в классическом смысле, а также устойчивым к возмущению рёбер контура [5]. По этой важнейшей причине было решено использовать метод суммирования угла в контексте геометрического ядра САПР, и ему по большей части посвящена данная глава работы.
В работе [53] была сделана попытка вычислить число вращения без обратных тригонометрических функций. Фактически, авторы предлагают запустить четыре луча параллельно координатным осям и обработать точки пересечения некоторым образом. Эта модификация метода очень похожа на метод пересечения луча, а значит в ней теряется устойчивость метода суммирования угла.
В данной главе показано, что запросы на суммирование угла вращения могут быть ускорены при помощи предподсчитанной иерархии ограничивающих объёмов. Идея состоит в том, что можно заменить связную подпоследовательность рёбер на один отрезок, если выполняются некоторые условия. Это позволяет быстро отсекать наборы рёбер. В результате запрос выполняется за ( log ) времени в худшем случае. Эта оценка выполняется только в том случае, когда все ограничивающие брусы в BVH точные. Через обозначено число обусловленности задачи, которое вводится в данной главе.
Обзор актуальной литературы. Заметим, что задача определения положения точки в классической постановке может быть решена за оптимальное время (log). Известно несколько оптимальных решений более общей задачи локализации точки в подразбиении (см. [54–56]). Эти алгоритмы существенно более сложные, чем те простые алгоритмы, которые описаны в данной работе. Некоторые из них могут быть обобщены на криволинейный случай. Анализ их устойчивости выходит за рамки данной работы. Однако весьма вероятно, что у них ещё больше проблем с устойчивостью, чем у простых алгоритмов.
2 из [50]: “The worst algorithm in the world for testing points is the angle summation method”. Несмотря на существование оптимальных алгоритмов, задача определения точки относительно контура продолжает активно изучаться. В работах [57, 58] запросы ускоряются при помощи равномерной сетки. В [57] многоугольник растеризуется на равномерной сетке. Затем множество пустых клеток, то есть клеток, не пересекающихся с многоугольником, разбивается на компоненты связности, классификация которых на внутренние и внешние выполняется методом пересечения луча по центру одной из клеток. Если проверяемая точка попадает в пустую клетку, то ответ определяется сразу, а в противном случае точка классифицируется по ближайшей точке многоугольника, которую предлагается искать среди предварительно записанных рёбер, пересекающих клетку. Это приводит к неверной работе алгоритма, что было исправлено в [58].
В работе [59] для ускорения алгоритма используется квадродерево. Идея заключается в том, чтобы адаптивным образом разбить рабочий брус на мелкие брусы, внутри каждого из которых лежит очень простая часть многоугольника. Брусы, не пересекающие многоугольник, явно помечаются как внутренние и внешние, а при попадании проверяемой точки в остальные брусы предлагается решить задачу разбором случаев. Следует заметить, что в работе разбираются лишь случаи пустой клетки, клетки, рассечённой ребром, и клетки с вершиной и двумя инцидентными рёбрами, хотя в результате предлагаемого построения могут получаться и другие варианты клеток.
В [60] используется структура данных с отсечением по полярному углу. Начало координат ставится примерно в центр многоугольника, и вся плоскость разбивается на области лучами, исходящими из начала координат. Легко понять, что это эквивалентно разбиению отрезка [0;2) полярных координат. Разбиение выполняется адаптивным образом, аналогично построению дерева отрезков по полярному углу. Для каждой области сохраняются все рёбра, через неё проходящие. Это позволяет классифицировать точку методом пересечения луча или методом заметающих треугольников, просматривая лишь те рёбра, которые сохранены для области, содержащей проверяемую точку.
Задачу можно решать, разбивая область на подобласти простого вида. Например, в [61] многоугольник разбивается на выпуклые многоугольники. При проверке точки используется дерево бинарного разбиения пространства (BSP) для того, чтобы быстро определить выпуклый многоугольник, содержащий точку.
В работе [62] используется метод заметающих треугольников. Для многоугольника предлагается просмотреть каждый треугольник, образованный ребром многоугольника и началом координат. Ответ определяется по количеству треугольников, содержащих точку, причём с учётом знака, зависящего от ориентации соответствующего ребра многоугольника. Подход расширяется на случай криволинейных рёбер двух типов. Для участка квадратичной кривой используется её уравнение для определения попадания точки в образуемый ею сегмент. Кубическая кривая Безье сначала при необходимости разбивается на части, а для частей попадание точки в образуемый кубической кривой сегмент определяется по неявному представлению кривой в виде уравнения.
В [49, Раздел 2] предлагается алгоритм локализации точек в случае, когда рёбра представлены кубическими кривыми. Выполняется поиск ближайшей точки пересечения луча, выпущенного из проверяемой точки, с рёбрами контура. Затем проверяемая точка класси фицируется на основе локальных дифференциальных свойств контура в окрестности полученной точки пересечения. Операции над кубическими кривыми выполняются посредством решения кубических уравнений.
При внимательном рассмотрении вышеперечисленных работ можно заметить, что в них отсутствуют строго доказанные нетривиальные оценки быстродействия запроса в худшем случае (разумеется, исключая классические работы [54-56] с доказанно оптимальными алгоритмами). Если асимптотические оценки сложности запроса и присутствуют, то они опираются на эвристические доводы и показывают время работы «в среднем», причём смысл термина «в среднем» не определяется. Доказанная в данной главе оценка 0{К\оє%) вы-полняется в худшем случае. Из неё легко следует время работы 0(logn) для выпуклых и звёздных многоугольников, что косвенно указывает на её нетривиальность.
Содержание главы. Данная глава организована следующим образом. Раздел 1.2 содержит общие определения и тривиальный алгоритм суммирования угла. Элементарное обобщение метода на случай многосвязной области, заданной несколькими контурами, приводится в разделе 1.3. В разделе 1.4 показано, как вычислить угол вращения вдоль отрезка прямой и дуги окружности. В разделе 1.5 определяется число обусловленности К, играющее важную роль в доказанной оценке быстродействия ускоренного алгоритма. В разделе 1.6 представлен алгоритм, ускоренный при помощи BVH, вместе с теоремой о временной сложности и экспериментальными данными. Адаптация алгоритма к NURBS кривым и смешанным контурам описана в разделе 1.7.
Прямое вычисление: поиск спа
Допустим, нужно найти все производные вплоть до второго порядка, т.е. = 2. Тогда для кубической кривой будет затрачено времени. Заметим, что это близко к объёму вычислений, требуемому для одного смешивания контрольных точек. Однако в реальности выполняется по одной операции смешивания на каждую вычисляемую производную. В частности, выполняется 3 операции смешивания для кривой и 6 для поверхности (при = 2). Это наблюдение ещё раз подтверждает, что смешивание контрольных точек — наиболее трудоёмкая часть вычисления NURBS.
Выбор алгоритма вычисления NURBS зависел от нужд разрабатываемого геометрического ядра.
Прежде всего, был необходим хорошо известный и надёжный алгоритм как твёрдая основа для дальнейших экспериментов. Поэтому был реализован стандартный подход, описанный в [35]. Базисные функции вычисляются через формулу Де Бура (см. раздел 2.4.1). Производные базисных функций также вычисляются известным образом (см. раздел 2.6.1). Далее выполняется смешивание контрольных точек (как в разделе 2.4.3). Наконец, в случае рационального NURBS находятся трёхмерные производные (как описано в разделе 2.6.3).
Потом возникла идея упростить вычисление базисных функций и избавиться от медленных операций деления. Был реализован предподсчёт базисных функций в степенном базисе, и они использовались для вычисления. Изначально использовались многочлены в виде , и вскоре была обнаружена серьёзная проблема с точностью для NURBS с большим количеством контрольных точек. Из-за этого было решено переключиться на многочлены в сдвинутом базисе (-), как описано в разделе 2.4.2. В результате перехода с формулы Де Бура на многочлены вычисление NURBS стало быстрее. Однако, некоторые операции геометрического ядра (например, пересечение поверхностей) всё ещё тратили существенное время на вычисление NURBS поверхностей.
Ситуация с NURBS кривыми была удовлетворительной: одно вычисление NURBS кривой не занимает много времени, особенно по сравнению со временем, затрачиваемым на обработку результата в алгоритме на более высоком уровне. Также можно заметить, что NURBS кривые создаются очень часто, в том числе в результате кубической аппроксимации некоторых кривых. Такие кривые обычно содержат сотни точек, а количество выполняемых вычислений для них довольно мало. Поэтому не было смысла предподсчитывать что-либо для дальнейшей оптимизации вычисления кривых.
В тот момент было ясно, что алгоритмический потенциал по ускорению вычисления NURBS поверхностей уже исчерпан. Мы не пробовали реализовать алгоритм Де Бура в 3D пространстве (описанный в разделе 2.5.1), потому что он концептуально медленнее обоих алгоритмов, использующих базисные функции. Единственный известный нам алгоритм, который выполняет меньше операций — вычисление предподсчитанных многочленов от двух переменных в 3D пространстве (описано в разделе 2.5.2). Однако, эта идея с предподсчи-тыванием кусочно-полиномиального представления была отвергнута из-за больших затрат памяти.
В конце концов было решено улучшить производительность на низком уровне. Количество операций при одном вычислении NURBS поверхности довольно велико, так что векторизация может дать положительный эффект. Как замечено в разделе 2.2, SSE2 сегодня повсеместно поддерживается в процессорах x86, так что можно использовать 16-байтные регистры и обрабатывать по два вещественных числа с двойной точностью одновременно.
Есть два принципиально различных способа векторизовать вычислительную процедуру. (вычисляется слева направо)
Первый способ — ускорить один запуск процедуры (т.е. вычисление в одной точке). Второй — векторизовать по нескольким вызовам процедуры (т.е. вычислять по 2 или 4 точки сразу). Пример использования второго способа для вычисления кубического многочлена представлен на рисунке 2.1. Этот второй способ не следует применять к вычислению NURBS по двум причинам:
1. Изменяется интерфейс процедуры: пользователь должен вычислять по несколько точек за раз, а не по одной. Это не представляется возможным для некоторых алгоритмов вроде метода Ньютона или трассировки кривой, которые существенно последовательные по своей природе. Кроме того, это просто неудобно — адаптировать алгоритмы к новому интерфейсу.
2. На каждом спане свой набор данных, используемых в вычислении. Если на вход подаётся 4 точки, они наверняка попадут в разные спаны. В результате совершенно различный набор данных (т.е. узлов, весов и контрольных точек) будет необходим для вычисления каждой точки. Непонятно, как эффективно векторизовать такой алгоритм.
Если все спаны совпадают, то для каждой точки используются одни и те же данные. Однако этого сложно добиться. На самом деле, это возможно только в том случае, когда пользователь запрашивает вычисления большого количества точек сразу (т.е. существенно больше, чем прямоугольных спанов в NURBS поверхности). В таком случае точки можно сперва распределить по спанам, а затем запустить на каждом спане серию векторизованных вычислений группами по 4 точки. В результате алгоритм должен использовать массовые вычисления, чтобы выиграть от такой процедуры.
Обе причины приводят к одному заключению: векторизованная таким образом процедура пригодна только для массового вычисления NURBS. Однако в массовом случае скорее всего можно использовать GPU.
В первом способе векторизации точки вычисляются по одной, но векторизуется внутренняя логика. Преимущество этого способа в том, что интерфейс не меняется, а значит каждый алгоритм, который вычисляет точки на NURBS, может сразу выиграть от векторизации без каких-либо изменений. К счастью, наиболее затратная по времени часть (т.е. смешивание контрольных точек) может быть легко векторизована. Более того, можно вычислять все базисные функции одновременно, если правильно разместить коэффициенты полиномов в памяти. 12 10
Кусочно-полиномиальное представление
Если ничего не предпринимать в связи с отсутствием аффинной инвариантности, то на практических сценариях численные методы будут работать существенно хуже, чем могли бы. Параметры кривых и поверхностей не имеют никакого физического и геометрического смысла, и пользователь может установить их области определения произвольно: как очень большими, так и очень маленькими. При существенно различных характерных размерах переменные системы получаются плохо масштабированными, и минимизация нормы вектора сдвига, которая выполняется в методе Ньютона с псевдообратной матрицей для недоопреде-лённой системы, теряет смысл, потому что некоторые переменные почти не будут изменяться.
Подробнее проблема масштабирования для метода Ньютона с псевдообратной матрицей рассмотрена в [121]. Простой практический приём решения этой проблемы заключается в масштабировании переменных и уравнений [122]. Предлагается выполнить замену переменных на новые, отличающиеся коэффициентом, и умножить каждое уравнение на некоторую константу. Конкретные коэффициенты можно подобрать так, чтобы в результате переменные и уравнения были примерно равной величины. Данный подход похож на предобуславливание непосредственно нелинейной системы уравнений [123].
Для масштабирования переменных заметим, что трёхмерное пространство имеет физический смысл, и можно использовать его как эталон для выбора масштаба. Предлагается выбрать масштаб переменных так, чтобы новые параметризации многообразий были примерно натуральными, то есть чтобы единичное изменение параметра многообразия приводило к примерно единичному смещению точки в пространстве. Для этого достаточно вычислить в небольшом количестве пробных точек частные производные многообразий по всем параметрам и для каждого параметра найти среднюю норму производной. Эти величины предлагается использовать в качестве масштабов параметров.
При масштабировании уравнений важно учитывать тот факт, что это влияет на понятие є-решения. В несингулярном уравнении невязка является в точности смещением в пространстве. Невязка по норме меньше є ровно тогда, когда соответствующие точки находятся на расстоянии є в пространстве. Поскольку трёхмерное пространство считается изотропным, никакое масштабирование для уравнений несингулярной системы выполнять не нужно.
Для сингулярных систем невязку можно разложить на две компоненты: d — невязку уравнений совпадения из — невязку уравнений сингулярности. Первые уравнения не нуждаются в масштабировании, однако вторые имеют кардинально отличающийся масштаб. Если полагать, что параметры измеряются в секундах, а координаты — в метрах, то невязка уравнений совпадения будет иметь размерность метров, а уравнений сингулярности — метр2/сек2 для пересечения кривых и метр3/сек3 для остальных типов пересечений.
Благодаря масштабированию переменных можно считать метр и секунду соразмерными, поскольку параметризации примерно натуральны. Однако предлагается другой подход, который никак не зависит от масштабирования переменных. Можно изменить уравнения сингулярности так, чтобы участвующие в них частные производные нормировались перед использованием в векторном или смешанном произведении. Например, если обозначить опе 109 ратор нормирования вектора через Л/, то уравнение сингулярности для кривой и поверхности примет вид ,fdS\ ,fdS\ ,(dC\\ Л/ р— х Л/ —— ЛМ— = 0. аи ov at
После добавления нормирования в сингулярные уравнения их невязка автоматически станет безразмерной. Более того, при малых значениях невязка имеет смысл некоторого угла. Для пересечения кривых норма невязки равна синусу угла между касательными. Для пересечения кривой и поверхности невязка примерно равна телесному углу трёхгранного угла между касательной кривой и изопараметрическими касательными поверхности, что в силу регулярности поверхности пропорционально с небольшим коэффициентом углу между касательной кривой и касательной плоскостью поверхности. Аналогичный геометрический смысл можно указать и для пересечения поверхностей.
Предлагается умножить безразмерные сингулярные уравнения на характерный масштаб задачи, определяемый по ограничивающим брусам пересекаемых многообразий. Тогда примерный геометрический смысл их невязки станет расстоянием между касательными пространствами на расстоянии характерного размера задачи (при малом её значении и с точностью до коэффициента для поверхностей с неортогональными частными производными).
Подробный анализ ранга систем с нормированными уравнениями сингулярности не проводился, однако эксперименты подтверждают интуитивное предположение о том, что результаты анализа останутся такими же.
Для практического использования системы уравнений нужно уметь вычислять её невязку и матрицу Якоби (то есть значения и первые производные). Ручное выписывание производных занимает много времени и чревато ошибками, которые сложно обнаружить.
Поэтому для вычисления производных используется простейшее автоматическое дифференцирование. Любое промежуточное значение хранится вместе со своими производными по всем переменным системы в небольшом классе C++ с перегруженными арифметическими операторами. В результате достаточно записать формулы для вычисления левых частей уравнений, и матрица Якоби автоматически посчитается в процессе. Использование шаблонов языка C++ позволяет, не дублируя текст программы, дать возможность компилятору полноценно оптимизировать код вычисления для каждого типа пересечения.
В процессе оптимизации в некоторых операциях автоматическое дифференцирование было заменено на ручные формулы. Явные формулы позволяют выполнить нормирование одного вектора производной поверхности за 12 сложений, 24 умножения и один обратный квадратный корень. Для пересечения поверхностей вычисление уравнений сингулярности после нормирований сводится к вычислению 8-ми смешанных произведений (см. (В.1)), в которых есть общие пары векторов, которые можно переиспользовать. Ручная обработка производных позволяет вычислить уравнения сингулярности по нормированным частным производным за 51 сложение и 72 умножения.
Каждое многообразие при вычислении системы в точке вычисляется ровно один раз, при этом находятся все производные вплоть до 1-ого порядка в случае несингулярной системы и до 2-ого — в случае сингулярной.
Было замечено, что при использовании системы уравнений численные методы локального анализа иногда вычисляют значение в одной и той же точке несколько раз. При этом такие два одинаковых вычисления чаще всего происходят либо подряд, либо между ними есть очень мало других вычислений. Отследить это переиспользование на внешнем уровне довольно сложно, поэтому в функцию вычисления системы уравнений было добавлено кэширование. Дополнительно затрачиваемое на кэширование время мало, зато близкие кратные вычисления системы устраняются.
В кэше хранятся 32 последние уникальные точки вычисления системы. При новом вычислении в первую очередь проверяется совпадение входных данных с предыдущим запросом, а затем со всеми значениями в кэше. Если значение в кэше уже есть, новое вычисление не выполняется, а если нет, то оно выполняется, и новая точка сохраняется в кэше, возможно вытесняя при этом самую старую точку. Для ускорения проверки точек на равенство используется простой хэш от бинарного представления входных данных.
Сходимость метода Ньютона
Пусть имеется произвольная поверхность, и её разрезали на две части изопараметри-ческой кривой. Если теперь запустить алгоритм пересечения двух полученных частей, в результате должна получиться одна кривая разреза. Под гладким сочленением поверхностей имеется ввиду именно этот случай.
Такой случай вполне может возникнуть в реальной практике. Несмотря на то, что его легко обработать специальной проверкой, для общего алгоритма трассировки он представляет большую проблему, потому что вдоль кривой пересечения нельзя определить направление трассировки исходя из дифференциальных свойств поверхностей (все производные на обеих поверхностях совпадают).
Раз поверхности вдоль общей границы локально совпадают, то ядро матрицы Якоби () системы уравнений будет иметь размерность 2, как при совпадении поверхностей. Фактически, это ядро определяет плоскость направлений в 4D пространстве переменных, при движении в сторону которых обе поверхности будут выдавать одинаковую точку в 3D.
В качестве начальной скорости трассировки можно выбрать любое из этих направлений. Благодаря сдвигу точки на границу области в методе предиктор-корректор точка после первого шага будет действительно на общей граничной кривой поверхностей. Величина первого шага скорее всего получится очень маленькой, потому что из-за исходно неверно определённого вектора скорости порядок аппроксимации будет ниже четвёртого.
При выполнении последующих шагов трассировки направление в предикторе будет определяться как проекция смещения между точками на ядро матрицы Якоби, что указано в разделе 4.9.2. Поскольку смещение направлено вдоль общей граничной кривой, то на всех дальнейших шагах направление будет определяться верно, и можно будет выполнять шаги трассировки полноценной величины.
Примеры сочленения представлены на рисунке 4.16. Проблемный случай представлен на правой половине: поверхности вместе составляют одну гладкую поверхность. В этом случае трассировка успешно находит всю кривую с требуемой точностью (10-5 при единичном размере теста), которая состоит из 203 кубических звеньев. Аналогичный пример приводится в [129].
На левой половине представлен модифицированный пример, на котором сочленение негладкое. В этом случае успешно находятся три несингулярные части кривой сочленения, разделённые сингулярными точками (в которых нормали поверхностей совпадают). В сумме получается 76 звеньев на все три кривые пересечения.
Результаты работы всего алгоритма поиска пересечений представлены на страницах 139 и 140. Все примеры относятся к наиболее сложному и интересному случаю пересечения поверхностей. На второй странице приведена краткая легенда для подписей. Для каждого примера подписано: сколько пересечений и каких типов было найдено («результат»), сколько времени потрачено на работу всего алгоритма, а также его частей («время»), и сколько рабочих брусов получено в результате глобального анализа («брусы»). Все примеры были решены общим алгоритмом пересечения. Для этого в некоторых случаях аналитические поверхности были представлены в виде NURBS, а в некоторых случаях аналитические алгоритмы пересечения были специально отключены. Каждый пример был решён программой 100 раз подряд, время приводится в среднем на один запуск. Во всех запусках использовалась толерантность = 10-5.
В правой половине первой страницы с примерами показаны случаи полностью несингулярного пересечения. В этих случаях глобальный анализ завершается досрочно, определив отсутствие сингулярностей при помощи критериев из приложения Б. Количество рабочих брусов составляет около тысячи, благодаря чему алгоритм пересечения успешно находит все пересечения примерно за десятую долю секунды. Стоит отметить, что в одном из примеров участвует сильно искривлённая поверхность: геликоид с малым периодом. В других примерах используется «волновая поверхность», представленная в виде NURBS с 100 100 контрольными точками.
В левой половине первой страницы с примерами показаны случаи с точками ветвления, например, классический случай крестового пересечения цилиндров, а также окружности Вилларсо при рассечении тора плоскостью. Два теста из этого набора решаются довольно быстро (за десятки миллисекунд) благодаря тому, что количество рабочих брусов на итерациях глобального анализа растёт медленно (это возможно благодаря использованию точных критериев непересечения лоскутов поверхностей, см. приложение А). В примере со сферой и цилиндром полюс сферы попадает точно в точку ветвления, вызывая рост количества брусов в глобальном анализе, из-за чего время решения примера оказывается довольно большим.
На второй странице с примерами приводятся наиболее сложные для общего алгоритма случаи, в которых имеется либо сингулярная кривая пересечения, либо частичное совпадение поверхностей. Практически во всех этих случаях глобальный анализ работает плохо и завершается по жёсткому ограничению на количество рабочих брусов. Все эти случаи решаются корректно, хотя крайне медленно. Приводятся случаи поверхностей, касающихся по кривой нетривиального вида, и граничного касания поверхностей. Один из примеров показывает гладкое сочленение поверхностей по границе, успешно трассирующееся общим алгоритмом. Наконец, имеется пример с частичным совпадением поверхностей, в котором находится нетривиальная кривая пересечения.