электронное научно-техническое издание
НАУКА и ОБРАЗОВАНИЕ
Эя № ФС 77 - 305БЭ. Государствен над регистрация №0421100025. ISSN 1994-04OS_
Алгоритм автоматического построения сеток из четырёхузловых конечных элементов
77-30569/332595
# 02, февраль 2012
Станкевич И. В., Яковлев М. Е., Си Ту Хтет
УДК 519.3
МГТУ им. Н.Э. Баумана [email protected], [email protected] nyinthamg355@gmail .com
Широкое внедрение конечно-элементной технологии, как составной части комплексной автоматизации сквозного цикла: проектирование - конструирование -изготовление, определило появление большого количества комплексов и пакетов прикладных программ (КПП и ИЛИ). Условно всё их многообразие можно разделить на две большие группы: исследовательские и профессиональные. Первые имеют весьма узкую специализированную ориентацию. Определяющим здесь является быстрота разработки, отладки и оперативное проведение численных исследований. Как правило, исследовательские КПП сопровождаются самим авторским коллективом. После доработки пользовательского интерфейса, а иногда, и функционального ядра с учётом опыта эксплуатации, исследовательские комплексы программ могут быть доведены до уровня профессиональных ПИП, имеющих высокий сервисно-диагностический уровень. Например, к известным профессиональным ППП можно отнести "средние" и "тяжёлые" CAD/CAE/CAM системы, созданные в нашей стране такие, как МАРС, ЛИРА, МЕГРЭ-3Д, FEMHCA, МАК, АСТРА, КАСКАД-2, ASTA, APM WinMachine, и за рубежом - I-DEAS, CATIA, Pro/ENGINEER, MSC/NASTRAN, ANSYS, MARC+Mentat II, MATRA, Unigraphics, COSMOS/M, STAADIII, BEASY, CADdy и некоторые другие. Эти пакеты обладают очень высокой универсальностью и обеспечивают решение практически любой задачи, не содержащей особых сложностей. Для пользователей основная трудность состоит в овладении навыками сопровождения, причём процесс полного освоения такого пакета может оказаться достаточно длительным и трудоёмким, а иногда, и дорогостоящим.
В настоящее время сохраняется необходимость в дальнейшем развитии существующих и создании новых исследовательских КПП, что вызвано, во-первых, ростом требований к уровню проводимых численных исследований и, во-вторых, динамически расширяющимися возможностями современных вычислительных средств. Практика численных исследований убедительно показывает, что наряду с созданием и развитием программных комплексов общего назначения необходимо вести разработку целевых программ для решения задач в рамках одной или нескольких идейно близких математических моделей, поскольку такие программы значительно повышают эффективность вычислительного эксперимента в соответствующей предметной области.
Адаптация МКЭ (метод конечных элементов) к новым задачам требует пересмотра концептуальных взглядов на отдельные этапы его реализации и, в частности, на подходы к построению самих математических моделей, описывающих сложные теплофизические и физико-механические процессы, где и закладывается начальная посылка получения достоверных результатов. При этом остаётся проблема организации вычислительных процедур с максимальной степенью экономичности при соблюдении достаточно высокой точности.
Первым этапом при проведении численных исследований является построение сетки конечноэлементной модели. Автоматическое построение сеток конечноэлементных моделей с помощью симплекс-элементов для двумерных задач рассмотрено во многих работах, например, [1-4]. В данной работе рассматривается построение двумерных конечноэлементных моделей с использованием четырёхузловых конечных элементов, которые находят широкое применение, например, при математическом моделировании кинетики напряженно-деформированного состояния при исследовании процессов обработки металлов давлением или контактного взаимодействия вязкоупругопластических тел.
Общая последовательность построения состоит в следующем. В начале с помощью CAD-системы создаётся геометрическая 2Б-модель исследуемой конструкции. Для этого на данном этапе можно использовать "средние" системы типа AutoCAD 2011, SolidEdge, SolidWorks, ADEM, T-FLEX CAD, KOMnAC-3D, оснащённые вполне развитым аппаратом двухмерного параметрического моделирования. В отдельных случаях исходную геометрическую модель необходимо перестроить в модель, удобную для построения конечноэлементной сетки с учётом особенностей термомеханического нагружения. Затем, уже адаптированная геометрическая модель разбивается на конечные элементы.
В простейшем случае разбивку на четырёхузловые конечные элементы можно выполнить, поступив следующим образом. Исследуемая конструкция с помощью пакета AutoCAD 2011 "вручную" разбивается на достаточно крупные геометрические примитивы, которые в плане можно рассматривать как 8- и 6-узловые квадратичные конечные элементы, а затем эти геометрические примитивы уже в автоматическом режиме разбиваются на конечные элементы. После этого выполняется процедура автоматической сшивки конечноэлементных моделей выделенных геометрических примитивов.
В тех случаях, когда по ходу проведения расчёта необходимо перестраивать сетку, например, при исследовании кинетики напряжённо-деформированного состояния в процессе обработки материала давлением, всякие "ручные" процедуры, за исключением этапа подготовки начальной геометрии, должны быть исключены. В рамках данной работы предлагается алгоритм, позволяющий автоматически строить конечноэлементные модели, состоящие в основном из четырехузловых конечных элементов и некоторого весьма ограниченного количества трехузловых конечных элементов. Алгоритм реализуется в два этапа: на первом - триангулируется область анализа, а на втором -треугольная сетка перестраивается в четырехугольную.
Начальной информацией является информация о контуре области. Контур области задаётся в виде линий, каждая из которых рассматривается как одномерный квадратичный конечный элемент (рисунок 1,а). Нумеруются начальный, конечный и некоторый средний узел каждой линии. В соответствии с конечноэлементной технологией в пределах каждой
контурной линии узлы имеют локальную нумерацию j, k), а в пределах всего
контура те же узлы имеют глобальную нумерацию (таблица 1, рисунок 1,а). Кроме того,
каждому узлу п = 1, в описании контура ставятся в соответствие координаты (Хп, Уп ) и размер шага сетки кп (таблица 2). Для контура, изображенного на рисунке
1,а, N0 = 8.
Рисунок 1 - Начальное описание контура и его разбивка
Таблица 1
Описание контурных линий (в соответствии с рисунком 2)
Номер контурной линии Локальные узлы
1 ] к
1 1 2 3
2 3 4 5
3 5 6 7
4 7 8 1
Таблица 2 Описание узлов контурных линий
Номер узла х, мм у, мм И, мм
1 У1 К
2 Х2 У2 И2
п хп Уп К
N о ХМо Ум0 К
Далее, все контурные линии последовательно разбиваются на отрезки прямых с учётом размеров заданных шагов Ип . Таким образом, формируется массив контурных
узлов исходной области к = 1, К0 (рисунок 1,б).
Далее выполняется триангуляция области. Алгоритм триангуляции состоит в следующем. Рассматривается текущий контур, на котором размещены контурные узлы. На первом шаге текущий контур, естественно, совпадает с исходным контуром, который ограничивает триангулируемую область. На контуре определяются три последовательных узла, например, р, Ц и г, имеющие минимальный угол а между построенными на них
векторами цр и ЦГ (рисунок 2). Далее, анализируется геометрическое расположение оставшихся контурных узлов относительно трех выделенных узлов. При этом устанавливается: попадает ли рассматриваемый контурный узел, например, ? в сектор анализа, характеризуемый величиной угла а, и если да, то вычисляются расстояния между узлами р, г и ?: р? и Г. Если хотя бы одно из этих расстояний меньше
максимальной длины векторов цр и цг, то проверяются две возможности построения треугольников, показанные на рисунке 3. Выбирается та, которая обеспечивает минимум максимального угла в треугольниках. Кроме того, проверяется пересекаемость ребрами треугольников контура. Если построение проходит, то текущий контур разбивается на два подконтура, а затем каждый анализируется отдельно.
г
Рисунок 2 - Поиск минимального угла
а б
Рисунок 3 - Построение элементов: случай, когда на контуре найден узел
В том случае, когда не находится узла, который попадал бы в сектор, характеризуемый величиной угла а, или расстояния, соединяющие узлы, больше
максимальной длины векторов др и дг, то выполняется построение на основе выделенных узлов р, д и г. В зависимости от величины угла а возможны два
п
варианта построения, если угол а < —, то треугольник строится на выделенных узлах р ,
д, г и узел д из описания контура исключается (рисунок 4а). В противном случае, когда п
а>— определяется новый узел g путем деления угла а пополам и строятся два
треугольника (рисунок 4,б), при этом проверяется пересекаемость ребрами треугольников контура. Если построение не проходит, то узел g размещается на середине отрезка рг. Далее узел g вводится в описание текущего контура, а узел д во всех случаях из описания контура исключается.
Рисунок 4 - Построение элементов
После построения конечноэлементной модели выполняется ее оптимизация. Оптимизация состоит в коррекции координат внутренних узлов, то есть тех узлов, которые были построены внутри триангулируемых подобластей. Координаты узлов, расположенных на исходном контуре, не изменяются. Оптимизация координат внутренних узлов выполняется следующим образом. Определяются номера конечных элементов, в состав узлов которых входит данный узел, и их число, то есть, фактически определяются все конечные элементы, которые окружают данный узел. Затем, вычисляются координаты центров тяжести всех найденных конечных элементов. В качестве новых координат узла принимаются координаты, вычисленные как среднее арифметическое координат центров тяжести конечных элементов, в состав узлов которых входит данный узел. Эта процедура выполняется несколько раз. Как
показывают исследования, для оптимизации необходимо выполнить 5 - 7 итераций. Кроме того, при оптимизации необходимо оценивать соотношение длин сторон треугольников. Если максимальное отношение длин сторон треугольника превышает некоторую заданную величину х (обычно X = 2), то наибольшая по длине сторона треугольника делится пополам и соседние треугольники разбиваются на два.
Построение треугольной сетки и ее оптимизация завершают первый этап. На втором этапе треугольная сетка перестраивается в четырехугольную. Перестроение реализуется следующим образом. Последовательно рассматриваются треугольники исходной сетки и их соседи, которые с данным треугольником имеют общие стороны. В
качестве примера на рисунке 5 показан рассматриваемый треугольник (е), имеющий узлы р , ц, г и три соседних треугольника (а), (в) и (/). Четырехузловой конечный элемент можно построить, объединяя треугольник (е), например, с треугольником (а)
или (в) или (у). Выбирается комбинация, у которой минимален максимальный угол в
четырехугольнике. Треугольники, формирующие четырехугольник из дальнейшего геометрического анализа исключаются.
Рисунок 5 - Перестроение сетки
После процедуры перестроения треугольной сетки в четырехугольную, состоящей в объединении выбранных пар треугольников, может остаться некоторое количество треугольников, которые не удалось объединить. Как показывают исследования, количество треугольных конечных элементов составляет не более 5 - 7 % от общего числа конечных элементов, формирующих конечноэлементную модель, и не оказывает заметного влияния на результаты численных исследований.
Второй этап построения сетки также завершается итерационной процедурой коррекции координат внутренних узлов. Определяются номера всех конечных элементов, в состав узлов которых входит рассматриваемый узел, и их общее число. Затем, вычисляются координаты центров тяжести всех найденных конечных элементов. В качестве новых координат узла принимаются координаты, вычисленные как среднее арифметическое координат центров тяжести конечных элементов, содержащих данный узел. Эта процедура выполняется несколько раз. Как показывают исследования, для оптимизации необходимо выполнить 3 - 5 итераций.
Если процесс решения прикладной задачи предусматривает перестроение конечноэлементной модели после нескольких этапов нагружения, то геометрической основой построения новой сетки являются узлы тех конечных элементов, стороны которых образуют внешний контур области.
Изложенный алгоритм был реализован в виде комплекса прикладных программ "8ЕТКА-4К-2Б". На рисунках 6 и 7 показаны увеличенные фрагменты конечноэлементных моделей, построенных с помощью данного комплекса. В канонических областях сетки имеют выраженный структурированный характер, а в областях с криволинейной границей сетки отличаются некоторой нерегулярностью.
ш щщш
Рисунок 6 - Сетка в канонической области
Рисунок 7 - Сетка в области с криволинейной границей
Таким образом, разработанный алгоритм позволяет строить сетки в геометрически сложных двухмерных областях и может быть использован при создании алгоритмов для построения сеток в трехмерных областях.
Список литературы
1. Сабоннадьер Ж.-К., Кулон Ж.-Л. Метод конечных элементов и САПР. - М.: Мир,1989.-192 с.
2. Математика и САПР: В 2 кн. / П. Шенен [и др.] - М.:Мир, 1988. - Кн. 1. - 204 с.
3. Математика и САПР: В 2 кн. / П. Жермен-Лакур [и др.] - М.:Мир, 1989. - Кн. 2. - 264 с.
4. Скворцов А.В. Обзор алгоритмов построения триангуляции Делоне // Вычислительные методы и программирование, 2002, №3.-с. 14-39.
electronic scientific and technical periodical
SCIENCE and EDUCATION
_EL № KS 77 - 3Ü56'». .V;II421100025, ISSN 1994-jMOg_
Algorithm of automatic construction of meshes from quadrangular finite elements
77-30569/332595 # 02, February 2012
Stankevich I.V., Yakovlev M.E., Si Tu Htet
Bauman Moscow Technical University [email protected] [email protected] nyintham g355@ gmail .com
An algorithm for construction of finite element models consisting of quadrangular finite elements for two-dimensional domains with complex geometric design was developed. The main distinctive feature of the proposed algorithm by comparison with the existing ones was a possibility of construction of irregular meshes in multicoherent areas without their preliminary splitting into the one-coherent. It was established that meshes in areas of canonical form had expressed structured character, but meshes in areas with curved boundary were irregular. The developed algorithm can be used in the general finite element technology for solving the thermo-mechanical problems, which arise at creation and operational development of objects of the aviation and space technics.
Publications with keywords: algorithm, finite element technology, finite element model, quadrangular finite element, triangulation, mesh, geometric analysis
Publications with words: algorithm, finite element technology, finite element model, quadrangular finite element, triangulation, mesh, geometric analysis
Reference
1. Sabonnad'er Zh.-K., Kulon Zh.-L., The finite element method and CAD, Moscow, Mir, 1989, 192 p.
2. Shenen P., et al., Mathematics and CAD: In 2 books, Book 1, Moscow, Mir, 1988, 204 p.
3. Shenen P., et al., Mathematics and CAD: In 2 books, Book 2, Moscow, Mir, 1989, 264 p.
4. Skvortsov A.V., Overview of algorithms for constructing of Delaunay triangulation, Vychislitel'nye metody i programmirovanie 3 (2002) 14-39.