66

Подсистема алгоритма Рапперта

1. Постановка прикладной задачи.

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

Многочисленные методы построения разбиений для дву- и трехмерных областей с геометрической точки зрения подразделяются на три основных класса:

1. построение разбиения, осуществляемого преобразованием отображения разбиения области с геометрически простой формой;

2. построение разбиения, осуществляемое преобразованием уже существующего разбиения;

3. прямое, элемент за элементом, построение разбиения, начиная с задания распределения точек в области или на ее границе.

Для построения разбиения чаще всего применяются следующие методы:

• построение разбиения “вручную” с представлением всей необходимой информации в виде структуры данных;

• построение покрытия области делением нескольких крупных элементов на более мелкие;

• построение покрытия области элементами, начиная с задания распределения точек на ее границе;

• построение покрытия области на основе облака точек, распределенных внутри области;

• построение разбиения с помощью геометрического (симметрия, локальное или глобальное деление, т.д.) и/или топологического преобразования уже существующего разбиения;

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

Проблема оптимального разбиения пространства задачи является подчас очень сложной. На каждый элемент разбиения могут накладываться довольно жесткие ограничения. К тому же, в пространстве задачи могут быть некие характерные области, где параметры меняются довольно резко. Например, в задачах механики твердого деформируемого тела, такие области образуются вблизи концентраторов напряжений и в них напряжения меняются чрезвычайно резко. Такие области требуют более частого разбиения.

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

Среди наиболее же последних методов стоит выделить метод генерации треугольных иррациональных конечно-элементных сеток Рапперта.

Этот алгоритм довольно новый, он был опубликован в 1994 году. Его разработал Джим Рапперт ещё, будучи студентом, в рамках проекта, поддерживаемым NASA. Алгоритм адаптирован для метода конечных элементов и основные требования, предъявляемые к конечно-элементным сеткам удовлетворяются в нём автоматически. К его преимуществам можно отнести и то, что алгоритм легко параметризуется и управляется.

Если говорить точнее, этот алгоритм не является алгоритмом генерации, а лишь алгоритмом улучшения качеств сетки. Тем не менее, этот метод в настоящее время считается наиболее перспективным и заслуживающим внимания.

Именно этот метод является предметом рассмотрения в дипломной работе и именно на его основе реализована модификация системы Sigma.

Задание на модификацию системы:

• На основе алгоритма Рапперта разработать учебную подсистему триангуляции для CAE системы Sigma с возможностью гибкого формирования алгоритма генерации иррациональных конечно-элементных сеток

2. Описание алгоритма Рапперта.

2. 1. Основные понятия

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

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

Узел (node) — точка в которой сходятся грани и соприкасаются несколько элементов. В отличие от вершин элементов узел общий для всех соприкасающихся в одной точке элементов.

Сегмент (segment) — это отрезок, соединяющий две соседние точки лежащие на границе области разбиения, другими словами этот отрезок принадлежит контурам области.

Грань (edge) — это отрезок по которому граничат два соприкасающихся треугольника. Используя предыдущее определение можно сказать, что все множество отрезков, соединяющих вершины элементов, в любой момент работы алгоритма складывается из сегментов и граней.

Включенная точка (encroached point) — любая вершина текущей сетки, которая находится внутри окружности, радиусом которой является любой сегмент.

"Неправильный треугольник" (bad triangle) — это треугольник неудовлетворяющий глобальным условиям, которые накладываются на генерируемую сетку.

Вытянутый треугольник (skinny triangle) — треугольник, имеющий одну сторону, намного отличающуюся от других двух. То есть треугольник слишком вытянут вдоль некоторой прямой. На критериях "вытянутости" мы остановимся ниже.

2.2. Базовые процедуры алгоритма

Алгоритм опирается на две базовые процедуры: a) разбиение неправильного треугольника с введением нового узла и b) разбиение сегмента, принадлежащего границе области разбиения с введением нового узла.

При разбиении треугольника происходит следующее:

1. Вычисляются координаты центра окружности, описанной вокруг треугольника, подлежащего разбиению.

2. В эту точку добавляется новый узел.

3. Удаляется треугольник, подлежащий разбиению и прилегающие и добавляются новые, как показано на рисунке 23.

В
торая базовая процедура состоит в следующем:

1. Если в окружность, диаметром которой является сегмент, принадлежащий границе области разбиения, попадает точка, не принадлежащая этому сегменту, то сегмент делится на две части.

2
. Удаляется треугольник, которому принадлежал первоначальный сегмент, и добавляются два новых треугольника (см. рис. 24).

2.3. Принцип работы

Принцип работы алгоритма состоит в цикле определения "неправильных" треугольников, т.е. треугольников, не удовлетворяющих определенному критерию, и произведение над ним одной из базовых процедур. В оригинальном алгоритме, предложенным Джимом Раппертом используется один критерий — минимальный угол в треугольниках сетки. В каждом треугольнике сетки определяется минимальный угол и, затем, сравниваются минимальные углы для всех треугольников с целью нахождения минимального угла сетки. Другими словами, минимальный угол в сетке равен минимально возможному углу, образованному двумя сегментами, либо двумя гранями, либо сегментом и гранью, имеющими общий узел.

Для задачи генерации конечно-элементной сетки этого критерия недостаточно, дело в том, что функции, значения которых мы хотим получить в узлах сетки, нелинейны, при этом они могут иметь области с большими значениями своих производных (области концентрации напряжений). При больших размерах конечного элемента, велик риск попросту упустить максимальные значения напряжений, поэтому был введен еще один критерий — условие максимальной площади треугольника. Сочетание этих двух ограничений гарантирует, что в сетке углы всех треугольников не менее значения, заданного пользователем, а площади всех треугольников не более другого значения, также задаваемого пользователем.

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

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

Все описанные процедуры выполняются в определенной последовательности, они имеют разный приоритет. Общая схема работы алгоритма следующая:

1. Производится поиск включенной (encroached) точки перебором всех точек и сегментов (только сегментов, но не граней).

2. Если такая точка найдена то над сегментом, который включает ее выполняется процедура деления сегмента пополам и производится возврат на шаг 1. Если включенных точек не было найдено, то алгоритм переходит на следующий шаг.

3. Производится поиск треугольника с минимальным углом.

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

5. Производится поиск треугольника максимальной площадью.

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

7. Все условия удовлетворены. Конец работы алгоритма

Автор в своей работе приводит строгие математические доказательства, гарантирующие правильную работу алгоритма при минимальном угле α < 20°.

2.4. Модификация алгоритма Рапперта

В работе для Journal of Algorithms раскрыт только лишь чисто математический аспект проблемы, хотя вероятней всего алгоритм разрабатывался для решения прикладных задач и для генерации конечно-элементной сетки в том числе. К тому же, и математическая сторона была раскрыта не полностью, а только лишь с точки зрения сходимости и некоторых особенностей. Ряд самых неочевидных подводных камней не поднимается на обсуждение. Видимо это объясняется тем, что это проект NASA и, вероятно, имеет какую-то степень закрытости.

В том варианте, в котором он был предложен автором, алгоритм не обладал хорошей сходимостью и устойчивостью. Опытным путем было установлено, что углы α > 25° были недоступны для него. Для повышения качества работы алгоритма, помимо двух вышеперечисленных процедур была введена еще одна — процедура поворота диагонали. К сожалению, автор опустил этот момент и не упомянул о нем, несмотря на то что, ее роль в работе алгоритма очень высока. Эта процедура, на самом деле, неявно присутствует в оригинальном алгоритме Рапперта (см.рис. 23, 24), но все же эту процедуру лучше выделить в отдельную и остановиться на ней подробно. Это упрощает и понимание работы алгоритма, и его программирование. На рисунках видно, что помимо непосредственного добавления новых треугольников, производятся и еще некоторые операции. Эти операции есть ничто иное как смена диагонали.

Суть ее заключается в следующем. Два треугольника, имеющие общую грань, образуют четырехугольник, разделенный диагональю. Внутри такого образования диагональ имеет два возможных положения. На рисунке 25 показаны эти положения — отрезки ac и bd.

2.5. Критерии оценки улучшения свойств сетки

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

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

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

К каждому возможному варианту необходимо предъявлять ряд требований. Во-первых, это адекватность этого критерия, т.е. улучшает ли он свойства сетки на самом деле, и если да, то с какими оговорками. Может случится, что способ хорош внешне, прост, не трудоемок, но работает при ряде ограничений. Эти ограничения и могут свести на нет все преимущества предлагаемого метода. Во-вторых, необходимо учитывать машинную трудоемкость выбираемого метода. И в-третьих, эффективность его будет зависеть еще и от сложности алгоритма, реализующего предлагаемый метод.

2.5.1. Длина диагонали

Первый, возможно, самый очевидный критерий — это длина диагонали (см. рис. 26). Действительно, зрительно, диагональ ac кажется гораздо более удачной диагонали bd для минимального угла двух треугольников и поменяв эту диагональ мы, вероятно, добьемся улучшения свойств сетки. У этого способа, к тому же, есть преимущества: он не трудоемок в вычислительном плане и его просто запрограммировать. Но есть у него и критический недостаток — он неадекватен. На рисунке 26 показан случай, где использование критерия минимальной длины диагонали может привести к серьезному ухудшению свойств сетки.

Д
иагональ bd короче диагонали ac, поэтому, по критерию минимальной длины диагонали, диагональ ac надо поменять на bd. Но это приведет к появлению очень тонкого треугольника, а их появление крайне нежелательно, а часто и вовсе становится критическим. Этот случай уже не так просто формализуется, и запрограммировав его, легко убедиться, что это дополнительное условие сильно "утяжеляет" алгоритм, сводя на нет все его преимущества.

2.5.2. Минимальный угол в кластере

С
ледующий очевидный критерий можно получить сравнивая минимальные углы в двух треугольниках до и после смены диагонали, для предсказания эффективности этого действия (см. рис. 27). Такой критерий является абсолютно адекватным, т.к. он совпадает с основным глобальным критерием сетки — критерием минимального угла и это его существенное преимущество.

Применение такого критерия действительно кажется заманчивым в силу его объективности, но у него есть один специфический недостаток — трудоемкость. Для него необходимо вычислить 6 углов, а по объему операций трудоемкость этого превосходит более чем в 10 раз вычисление двух длин отрезков.