автореферат диссертации по информатике, вычислительной технике и управлению, 05.13.11, диссертация на тему:Разработка препроцессора подготовки данных для комплексов конечноэлементного моделирования контактных систем

На правах рукописи

ВДОВИЧЕНКО Антон Александрович

РАЗРАБОТКА ПРЕПРОЦЕССОРА ПОДГОТОВКИ ДАННЫХ ДЛЯ КОМПЛЕКСОВ КОНЕЧНОЭЛЕМЕНТНОГО МОДЕЛИРОВАНИЯ КОНТАКТНЫХ СИСТЕМ

Специальности: 05.13.11 - Математическое и программное обеспечение

вычислительных машин, комплексов и компьютерных сетей 05.13.12 - Системы автоматизации проектирования

АВТОРЕФЕРАТ

диссертации на соискание ученой степени кандидата технических наук

Воронеж - 2003

Работа выполнена в Воронежском государственном техническом университете

Научный руководитель Научный консультант Официальные оппоненты:

доктор технических наук, профессор Гольник Эдуард Рувимович!

доктор технических наук, профессор Подвальный Семен Леонидович

доктор физико-математических наук, профессор

Семыкина Татьяна Дмитриевна;

кандидат технических наук, доцент Питолин Андрей Владимирович

Ведущая организация Воронежская государственная -

технологическая академия

Защита состоится "30" октября 2003г. в 1130 в конференц-зале на заседании диссертационного совета Д 212.037.01 Воронежского государственного технического университета по адресу: 394026, г. Воронеж, Московский просп., 14.

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

Автореферат разослан "29" сентября 2003 г.

Ученый секретарь диссертационного совета

Питолин В.М.

ОБЩАЯ ХАРАКТЕРИСТИКА РАБОТЫ

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

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

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

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

Тематика диссертационной работы соответствует научному направлению Воронежского государственного технического университета "Вычислительные системы и программно-аппаратные электротехнические комплексы".

Целью диссертации является разработка методов, алгоритмов и программных средств автоматической подготовки данных для комплексов конечноэлементного моделирования бесконтг согласо-

ванными контактами.

Указанная цель предусматривает решение следующих задач:

1. Произвести анализ существующих методов и алгоритмов генерации сеток конечных элементов, используемых в препроцессорах подготовки данных современных конечноэлементных комплексов;

2. Сформулировать основные требования к разрабатываемому препроцессору, исходя из особенностей моделирования контактных систем;

3. Разработать методы анализа твердотельных геометрических моделей контактных систем и дискретного геометрического моделирования отдельных тел;

4. Разработать программно-алгоритмическое обеспечение, позволяющее в автоматическом режиме проводить выделение согласованных контактных зон, построение сопряженных пар контактирующих узлов и генерацию конечноэлементных сеток в объемах моделируемых тел;

5. Реализовать препроцессор и произвести его интеграцию в комплекс моделирования и анализа контактных систем (МАКС) в качестве модуля подготовки данных.

Методы исследования. В работе использованы методы аналитической и вычислительной геометрии, многопараметрической оптимизации, линейной алгебры, функционального и объектно-ориентированного программирования.

Научная новизна диссертации. Результатами работы, отличающимися научной новизной, являются:

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

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

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

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

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

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

Практическая значимость работы заключается в разработке и внедрении в практику промышленного применения препроцессора подготовки данных, реализованного в качестве программного модуля комплекса конеч-ноэлементного моделирования МАКС, используемого для решения задач промышленного проектирования контактных систем.

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

Апробация работы проведена: на научных конференциях Воронежского государственного технического университета в 2000-2003 годах; на V международной электронной научной конференции «Современные проблемы информатизации в технике и технологиях», (Воронеж, 2000); на П всероссийской научно-технической конференции «Прикладные задачи механики и тепломассообмена в авиастроении» (Воронеж, 2001); на Ш международной конференции «Авиакосмические технологии» (Воронеж, 2002); на международной конференции «Современные сложные системы управления» (Воронеж 2003); на международной научно-технической конференции «Системные проблемы качества, математического моделирования, информационных и электронных технологий» (Сочи, 2003).

Публикации: основные результаты исследований по диссертации отражены в 8 опубликованных работах, в том числе 3 без соавторов. В работах, опубликованных в соавторстве и приведенных в конце автореферата, лично соискателем предложены: в [1,6] - методы автоматизации процедур по обеспечению сопряженности узлов в зонах контакта; [5] - алгоритмы генерации, управления качеством и оптимизации конечноэлементной сетки; [3,4] - методы интеграции препроцессора с комплексами конечноэлементного моделирования.

Структура и объем диссертации: введение, 4 раздела, заключение, список литературы (110 наименований) и 2 приложения; материалы диссертации (без приложений) включают 140 страниц текста, 30 рисунков, 1 таблицу.

ОСНОВНОЕ СОДЕРЖАНИЕ РАБОТЫ

Введение содержит общую характеристику диссертации, кратко изложенную выше в данном автореферате.

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

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

На основе анализа современной литературы, посвященной методам генерации КЭ сеток, рассмотрены преимущества и недостатки существующих подходов. Основное внимание было уделено методам дискретизации на основе различных способов плоской и пространственной триангуляции. Рассмотрены возможности данных методов по построению как регулярных, так и нерегулярных сеток.

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

Наиболее перспективными и универсальными на данный момент являются подходы на основе метода распространяющегося фронта (Advancing Front Method). Исходными данными для этого метода является дискретная модель границы разбиваемой области. Так в двухмерном случае такой моделью служит один или несколько полигонов - это начальный фронт. Основой метода является добавление треугольников или тетраэдров в область, таких, что им принадлежит хотя бы одно ребро или грань из фронта. Такая процедура на каждом шаге изменяет фронт. Когда фронт становится пустым, процедура построения сетки завершается.

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

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

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

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

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

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

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

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

Функция масштаба М = М(р), где р - точка в трехмерном пространстве р еЯ3, вводится с целью характеристики каждой точки поля рассматриваемого пространства над генерируемой сеткой и имеет следующий смысл: для любого узла, находящегося в точке р, значение М(р) определяет «желательную» длину каждого ребра (линейный размер элемента), исходящего из этого узла. В общем случае функциия М(р) определяется множеством масштабирующих точек и имеет вид

1 ш

М(р)^тш{тт№,гг(1-—)-~^У, (;=;..£), (1)

где: Шта - максимальное значение функции масштаба; Р~ коэффициент роста размера элемента, значение которого соответствует конструируемому «потенциальному» полю изменения линейных размеров конечных элементов; т,- коэффициент масштаба г-ой масштабирующей точки; г, - расстояние от г'-й масштабирующей точки до точки р; к ~ число масштабирующих точек.

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

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

Определяется вспомогательная функция масштаба для любой точки каждого ребра е1 тела в виде

Ме, (р) = т т(М0 (р), М,г(р),1,), (2)

где /,• - длина ребра, на котором происходит разбиение; М^ — составляющая функции (2), определяемая особенностямиу'-й поверхности ^ (/ = 1..2), которой принадлежит точка р

М^(р) = тт{ЯтЫ(8^р),Ос}(р),...,Вск(р)),] = 1.2. (3)

В формуле (3) компонента Ят„(Б;р) соответствует минимальному радиусу кривизны поверхности ¿) в точкер; Д/р) - минимальное расстояние от точки р до ребра С„ (]=!..к), где С„ - те ребра грани поверхности которым не принадлежит точка р.

Далее для всех ребер е, проверяется - существует ли точка х, принадлежащая выбранному ребру, такая, что условие

М{х)<МеХх) (4)

не выполняется. При обнаружении такой точки добавляется новая масштабирующая точка с координатами х и коэффициентом масштаба т=Ме[х). Процедура поиска продолжается до тех пор, пока существуют точки, не удовлетворяющие условию (4).

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

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

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

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

Для предложенного метода приводится доказательство инвариантности выбора направления зазора и проецирования узлов в зоне согласованного контакта на примере случая внутреннего контакта двух цилиндров близкого радиуса (рис. 1).

Д <а-М(р),

(5)

Рис. 1. Тестируемая контактная система: «-конфигурация и нагружение; б -геометрия вариантов начальных зазоров

Если принять, что попарно сопряженные узлы располагаются на радиальных направлениях из центра О, цилиндра меньшего радиуса Л, (рис. 16), то начальные зазоры описываются функцией угла ср,.

Если же считать сопряженными узлы, лежащие на радиальных прямых, исходящих из центра Ог окружности большего радиуса Л2, то функция начальных зазоров имеет вид:

Сравнение результатов моделирования зон контактного сжатия, полученных по Герцу и по МКЭ на основе комплекса МАКС в двух вариантах, т.е. с применением функций (6) и (7), показывает, что при согласованном контакте, когда Л2 ~ Л > и начальные зазоры относительно малы, указанные подходы практически тождественны.

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

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

Предложенный и реализованный в данной работе методы дискретного геометрического моделирования основывается на последовательном добавлении элементов в строящуюся КЭ сетку, геометрические параметры которой для каждого элемента определяются с помощью функции масштаба. Главным инструментом, на базе которого конкретизируется порядок добавления элементов, и вторым (после функции масштаба) базисным понятием алгоритма является фронт (граница) генерируемой сетки.

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

(6)

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

Ключевым элементом метода генерации сетки является алгоритм добавления нового узла. Новые узлы добавляются на основе пары смежных ребер фронта в плоском случае (рис. 2) и на основе пары смежных треугольников в трехмерном случае (рис. 3). В обоих случая пара выбирается исходя из условия минимальности угла - плоского или двухгранного соответственно -образуемого ее элементами.

А - .

X X у

Рис. 2. Предусмотренные варианты добавления узлов: а - добавляется 0 узлов; б - добавляется 1 узел; в - добавляется 2 узла. АВ и ВС - выбранная пара ребер

По мере добавления новых элементов область, охватываемая фронтом, уменьшается до тех пор, пока не сожмется до одного треугольного или соответственно тетраэдрального элемента.

Рис. 3. Предусмотренные варианты добавления узлов и новых граней при построении тетраэдральной сетки: а - двухгранный угол «выбранной пары граней ABC и BCD; б ~ добавление одной грани; в - добавление двух граней; г - добавление нового узла X.

Координаты добавляемого узла определяются после нахождения неизвестного угла <р на основе функции масштаба М(р) исходя из соотношений М(А) + М(С) 2

х _ У

ь = -

М(А) M(Q' х2 = а2 + Ъ2 - 2abcos<p [у2 = с2 +b2 -2cb cos (а - <р) в плоском случае (рис. 4а) и х _ у

(8) (9)

М(А)

2

М(Р) а2+Ь2

,У= |XD|,

(10)

(И)

(12)

-2ab-cos ß-cos(p [3'* =с2 +b2 -2cb-cosy-cos(a-gj) в трехмерном случае (рис. 46), где c=|PD|, b определяется на основе функции масштаба из соотношений аналогичных (8), угол /строится аналогично углу ß в плоскости треугольника BCD.

а)

б)

Рис. 4. Добавление нового узла: а - плоский случай; б -трехмерный случай

После завершения этапа генерации элементов к сетке применяется процедура повышения качества на основе последовательного использования методов сглаживания и локальной реконфигурации связей узлов. Критерием качества плоских элементов - треугольников - служит величина минимального угла. Для повышения качества поверхностных сеток используется оптимизационная процедура сглаживания. Целевая функция строится, исходя из принципа максимизации минимального угла, и имеет вид:

в(п) = тт[атт (г,.)] (13)

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

Критерий качества объемных элементов - тетраэдров - вычисляется по формуле:

в(1)=ТГ~ (14)

Л '"шах

где / - тетраэдр, качество которого измеряется; 2 - функция, отражающая качество; V - объем тетраэдра; Я - площадь его поверхности; - длина максимального ребра.

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

Процедура дискретизации производится на основе предложенной иерархической геометрической модели представления данных, построенной для заданной контактной системы или произвольного тела. Внутреннее представление примитивов в такой модели является универсальным и основывается на представлении поверхностей и кривых в виде неоднородных рациональных би-сплайнов (МТЖВЗ-сплайнов). Для уменьшения вычислительных расходов и удобства алгоритмизации процедур геометрического ядра препроцессора производится предварительная обработка данных, которая позволяет перейти от канонического рекуррентного представления МЖВБ-сплайнов к представлению в виде набора дробно-рациональных функций Лр{и,у)

где С,/и,у•) - функция задающая сплайн на конкретном узловом интервале, Ар и В4- полиномы степеней р и д соответственно. Степень и вид этих полиномов определяется на основе степени >ПЖВ8-сплайна к, матрицы контрольных точек Р, матрицы весов \У, и векторов узловых точек X] и х2.

Структура модели представления данных приведена на (рис. 5).

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

= (15)

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

Рис. 5. Структура иерархической модели представления геомет-

рических данных в препроцессоре

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

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

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

Алгоритм добавления элементов в сетку является общим для плоского и пространственного случая. Он основывается на методике, описанной во второй главе, и схематически может быть представлен в виде обобщенной блок-схемы (рис. 6).

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

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

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

Рис. 6. Алгоритм добавления элементов в КЭ сетку

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

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

Четвертая глава в своей первой части содержит описание реализации препроцессора подготовки данных для комплекса МАКС. Приведена структура препроцессора (рис. 7).

Рис. 7 Структура препроцессора подготовки данных

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

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

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

Вторая часть раздела посвящена описанию примера применения препроцессора для построения конечноэлементной модели реальной контактной системы - части кривошипного горячештамповочного пресса (КГШП).

В конструктивную модель (рис. 8), построенную для этой контактной подсистемы с использованием программного комплекса БоНсГУУогкв, включены: двухэксцентриковый вал - 3; шатуны -2,4; опорные втулки станины - 1,5.

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

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

Для сеток, представляющих собой конечноэлементную модель эксцентрикового вала (рис. 9а) и всей подсистемы в сборе (рис. 9е), производится статистический анализ качества формы элементов на основе трех различных критериев, что позволяет оценить корректность результатов работы препроцессора. Наиболее адекватным является критерий на основе соотношения (14), его распределение представлено на (рис. 96), откуда видно, что подавляющее большинство элементов сетки имеют высокое значение параметра качества (0,8-0.9), что соответствует форме тетраэдра близкой к правильной - наиболее оптимальной для МКЭ.

Рис. 8. Конструктивная модель рассматриваемой подсистемы КГШП

Р вел ре двоение параметре мчеета* фориы 0*М««Г

1 8 \ | %Г-Л- , А-% V? "Г* & >м

900-010 о<«-0.20 егмао слс««о о<мн о»«.к> мм?д Знача им параматра о я»«

б) в) Рис. 9. Результаты дискретной аппроксимации твердотельных моделей: а,в) примеры дискретизации эксцентрикового вала и контактной системы в сборе; б) гистограмма распределения параметра качества формы элементов

ОСНОВНЫЕ РЕЗУЛЬТАТЫ РАБОТЫ

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

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

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

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

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

6. С помощью реализованного препроцессора произведена процедура автоматической дискретизации контактной системы, представляющей собой модель части реального промышленного объекта: кривошипного горяче-штамповочного пресса.

ОСНОВНЫЕ РЕЗУЛЬТАТЫ ДИССЕРТАЦИИ ОПУБЛИКОВАНЫ В СЛЕДУЮЩИХ РАБОТАХ:

1. Вдовиченко A.A., Гундорова Н.И., Успехов A.A., Проблемы автоматического построения конечноэлементных сеток в задаче компьютерного моделирования контактных систем // Современные проблемы информатизации в технике и технологиях, Сб. тр. Воронеж: ЦЧКИ, 2000, с. 42-43.

2. Вдовиченко A.A., Алгоритмическое и программное обеспечение генерирования сопряженных ансамблей МКЭ для многотельных контактных систем // Прикладные задачи механики и тепломассообмена в авиастроении, Сб. тр. Воронеж: ЦЧКИ, 2001, с. 191-197.

3. Гольник Э.Р., Успехов A.A., Гундорова Н.И., Вдовиченко A.A., Ляпин С.И., Модернизация программного комплекса МАКС на уровне современных требований механики авиаконструкций // Авиакосмические технологии, Сб. тр. Воронеж: ЦЧКИ, 2002, с. 113-120.

4. Гольник Э.Р., Гундорова Н.И., Ляпин С.И., Вдовиченко A.A., Успехов A.A., Повышение эффективности моделирования МКЭ объектов аэрокосмической техники на основе интегрирования комплексов МАКС и Solidworks // Авиакосмические технологии, Сб. тр. Воронеж, 2002: ЦЧКИ, с. 120-123.

5. Вдовиченко A.A., Подвальный С.Л., Успехов A.A., Алгоритм автоматической генерации управляемых и оптимизируемых конечноэлементных сеток для геометрических тел произвольной конфигурации // Социально-экономическое развитие регионов: реальность и перспективы, Сб. тр. Воронеж 2003: Новый Взгляд, с. 502-507.

6. Вдовиченко A.A., Подвальный С.Л., Успехов A.A., Автоматическое построение системы сопряженных конечноэлементных сеток для произвольной конфигурации контактирующих тел // Современные сложные системы управления, том 2, Сб. тр. Воронеж: ВГАСУ, 2003, с. 38-43.

7. Вдовиченко A.A., Алгоритм генерации конечноэлементной сетки для криволинейных поверхностей заданных в виде NURBS сплайнов // Электротехнические комплексы и системы управления, Сб. тр. Воронеж: ВГТУ, 2003, с. 105-112.

8. Вдовиченко A.A., Обеспечение условия односвязности и ориентированности дискретной модели поверхностной оболочки для метода пространственной дискретизации // Системные проблемы качества, математического моделирования, информационных и электронных технологий. Часть № 8. Москва: Радио и связь, 2003, с. 47-51.

/

ЛР № 066815 от 25.08.99. Подписано в печать 25.09.03.

Формат 60x84/16. Бумага для множительных аппаратов.

Усл. печ. л. 1,0. Тираж 90 экз. Заказ № г

Воронежский государственный технический университет 394026 Воронеж, Московский просп., 14

* 15 1 8 S

2.0С.З-Д

Я

Введение

1. Типовая структура комплекса конечноэлементного моделирования и место препроцессора в ней. Анализ современных методов и алгоритмов генерации конечноэлементных сеток

1.1 Особенности численного моделирования контактных систем

1.2 Основные требования, предъявляемые к сеткам исходя из особенностей методов численного моделирования

1.3 Обзор современных методов генерации КЭ сеток

1.3.1 Методы генерации регулярных сеток

1.3.2 Методы генерации нерегулярных сеток

1.3.3 Методы повышения качества сеток

1.4 Преимущества и недостатки рассмотренных методов

1.5 Обоснование целей и задач диссертации на основе анализа существующих методов генерации сеток и проблем моделирования контактных систем

2. Методы генерации оптимизируемых КЭ сеток с управляемым качеством для моделирования систем контактирующих тел

2.1 Основные особенности предлагаемого подхода

2.2 Функция масштаба, как основное средство управления качеством сетки на этапе ее генерации

2.3 Формулировка условия согласованности контакта. Метод обеспечения сопряженности узлов в областях контакта

2.3.1 Задача о согласованном контакте кругового цилиндра и цилиндрического тела с круговым вырезом

2.3.2 Определение согласованного контакта. Метод построения сопряженных узлов в зонах согласованного контакта

2.4 Метод построения плоской сетки на произвольной поверхности. Учет функции масштаба и зон контакта

2.4.1 Декомпозиция исходной геометрической модели

2.4.2 Подцепление внутренних контуров

2.4.3 Метод построения плоской поверхностной сетки

2.4.4 Построение сеток для криволинейных поверхностей

2.5 Генерация КЭ сетки в объемах контактирующих тел

2.6 Оптимизация геометрической формы полученных конечных элементов. Построение целевой функции и выбор метода оптимизации

2.6.1 Управление качеством триангуляции на поверхности

2.6.2 Управление качеством пространственной сетки

2.7 Выводы

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

3.1 Иерархическая модель представления геометрических данных и ее оптимизация, исходя из особенностей конечноэлементного моделирования

3.2 Построение функции масштаба на основе геометрических особенностей тела и заданных критериев качества

3.3 Выделение зон контакта и построение сопряженных КЭ сеток в этих областях

3.4 Построение КЭ сетки для каждого из контактирующих тел

3.4.1 Алгоритм размещения узлов на ребрах, принадлежащих текущему геометрическому телу

3.4.2 Последовательное построение плоских КЭ сеток на поверхностях обрабатываемого тела. Алгоритм последовательного добавления узлов и элементов

3.4.3 Построение трехмерной сетки в объеме разбиваемого тела. Особенности алгоритма по сравнению с плоским случаем

3.5 Особенности алгоритма оптимизации формы элементов

3.6 Выводы

4. Реализация предложенного метода в качестве подсистемы комплекса МАКС. Примеры построения КЭ сеток систем контактирующих тел

4.1 Структура препроцессора подготовки данных

4.2 Особенности реализации препроцессора в виде динамической библиотеки

4.3 Реализация препроцессора в виде отдельной программы

4.4 Примеры использования препроцессора

4.4.1 Построение геометрической КЭ модели системы контактирующих деталей кривошипного горячештамповочного пресса. Оценка адекватности и качества модели

4.4.2 Примеры использования препроцессора для комплекса МАКС, интегрированного в среду автоматизированного проектирования Solidworks

4.3 Выводы

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

Современный период развития механики, микроэлектроники, электротехники, квантовой физики, биофизики и других областей человеческой деятельности, наступил в результате стремительного перехода (с 1960-ых годов) на компьютерные технологии, позволившие строить и широко использовать дискретные модели и специальные, ориентированные на компьютеры, численные методы. Наиболее универсальным, развитым и эффективным среди них является метод конечных элементов (МКЭ), в теоретическую разработку которого решающий вклад внесли Дж.Аргирис [3,34], А.Ф. Смирнов [32,33], О.Зенкевич [25,26], Р.Галлагер [15], Г.Стренг, Дж.Фикс [35], К. Бате, Е. Вилсон [4], А.В. Александров [2], Н.Н. Шапошников [38] и др. В последние десятилетия как конкурентоспособный с МКЭ следует рассматривать метод граничных элементов (МГЭ). Авторами основополагающих работ по МГЭ являются К. Бреббиа [7], П. Бенерджи, Р. Баттерфилд [6], Ж. Теллес [8], С. Крауч [30] и др.

Основной идеей большинства современных численных методов, и МКЭ в частности, является дискретизация моделируемых объектов в виде множества элементарных примитивов, таких как треугольники и четырехугольники в двухмерном пространстве или тетраэдры и гексаэдры в трехмерном пространстве, - конечных элементов (КЭ) [55,43]. Множество таких элементов, описывающих конкретное геометрическое тело или систему тел, называется конечноэлементной сеткой. Алгоритм такой дискретизации для двухмерных и в особенности трехмерных объектов является отдельной специальной задачей, в большинстве случаев не уступающей, а зачастую превосходящей по сложности сам процесс КЭ моделирования. Эволюция методов генерации сеток конечных элементов идет параллельно с развитием самого МКЭ, и в настоящее время прошла путь от ручных и частично автоматизированных методов [99,109] до полностью автоматических методов [103,105,45,48,64,65], таких как метод распространяющегося фронта [60,83] и подходов, основанных на триангуляции Делоне [60,88,100], обеспечивающих генерацию КЭ сеток для тел произвольной конфигурации.

Таким образом, в настоящее время задачу построения КЭ сеток для отдельных геометрических тел можно считать достаточно изученной и в определенной степени решенной, но все существующие методы обладают теми или иными недостатками, не позволяющими применять их в специальных случаях. Кроме того, не во всех областях, где применяется МКЭ, таких методов достаточно. В частности в механике при решении контактных задач встает проблема построения специального вида КЭ сеток - контактных сеток, требующих обеспечения дополнительного условия сопряженности узлов в зонах контакта. Необходимость построения такого рода сеток появилась достаточно недавно. Она связана с последними достижениями в области численного моделирования контактных систем [18] и в связи с этим в настоящее время не существует полностью детерминированного автоматического метода генерации контактных сеток, так что эту проблему можно считать открытой.

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

Актуальность темы данной диссертации определяется реальными практическими потребностями разработки и реализации методов подготовки данных для задач КЭ моделирования, как отдельных тел, так и контактных систем. В данном исследовании на фоне общей методологии применения метода конечных элементов (МКЭ) [25] обсуждается проблема, первая по необходимости ее неотложного решения, а именно - создания теоретической основы, структуры и управляемых возможностей построения и эксплуатации так называемых препроцессоров подготовки данных.

Преодоление сформулированной Р. Беллманом известной проблемы «проклятия размерности» [5] остается фундаментально важной и чрезвычайно актуальной задачей, поскольку многомерность дискретных моделей МКЭ стремительно возрастает и характеризуется миллионами неизвестных. Подготовка данных для таких массивов информации непосредственно связана с наличием препроцессоров необходимой мощности и имеет свою научную историю. Анализ достижений в этой области, начиная с пионерской системы ASKA [40] вплоть до новейших разработок CAD/CAE/CAM [24], обязывает считать эту проблему механики и информатики открытой.

В настоящее время бурный прогресс в области развития вычислительной техники позволил на порядки увеличить размерность моделей, используемых в задачах промышленного проектирования, но как показывает практика, даже современных вычислительных ресурсов не хватает для моделирования реальных физических моделей в целом. Очевидно, что эта проблема не может быть решена только за счет увеличения вычислительной мощности компьютеров и наращивания объемов памяти. Таким образом, для уменьшения размерности дискретных моделей уже на этапе подготовки данных требуется разработка новых высоких интеллектуальных технологий, в первую очередь, мощного и гибкого программного инструментария — препроцессора, проблемно-ориентированного на подготовку данных и решение современных больших задач МКЭ. Генеральным направлением продвижения к решению проблемы, в том числе, учета алгоритмической специфики и выбора основ построения препроцессора на базе принятой далее триангуляции, является теоретически базирующийся на методе контактных сил и переносных перемещений [18] комплекс МАКС (Моделирование и Анализ Контактных Систем) [21,17].

Для минимизации размерности дискретных моделей на этапе подготовки данных наиболее эффективным является метод варьирования линейных размеров элементов на основе априорных данных о зонах концентрации напряжений, а также геометрических параметров исходной твердотельной модели. Основной недостаток известных технологий построения сеток МКЭ заключается в их чрезмерной априорной жесткости схем, неуправляемости, недостаточной чувствительности к виртуальной вариативности. Уменьшение линейных размеров элементов всей модели приводит к неоправданным вычислительным затратам в областях, где высокая точность моделирования не требуется, при этом процент объема таких областей от общего объема модели в некоторых случаях может приближаться к 100 процентам. В то же время глобальное укрупнение размеров элементов приводит к значительным потерям точности в зонах концентрации напряжений, вплоть до некорректности решения, а также ухудшению качества элементов, в областях, обладающих геометрическими особенностями с размерами значительно меньшими выбранного размера элемента, что также может стать причиной не точного или не корректного решения [79].

Таким образом, наличие метода автоматической подготовки данных, реализующего варьирование размеров элементов в зависимости от требуемой точности и наличия геометрических особенностей и осуществляющего контроль за качеством генерируемых элементов обеспечивает необходимый компромисс между требованиями к точности и корректности решения и допустимыми вычислительными затратами [89,95].

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

В новейшей дискретной механике контактных взаимодействий [18,23,28] доминирует концепция сопряженности узлов конечных элементов смежных тел. Узлы в предположительно контактных зонах называются сопряженными, если они располагаются (назначаются, выбираются) так, что в результате приложения заданных нагрузок узловые координаты становятся попарно совпадающими, а сами узлы контактирующими в соответствующих точках.

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

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

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

В настоящей работе ключом выхода из, казалось бы, теоретически тупикового положения служит предложенная К. Джонсоном [23] идентификация контактных систем на два основных класса взаимодействий, а именно -согласованных и несогласованных по форме тел.

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

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

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

Вследствие большого количества ручных операций процедура подготовки данных является наиболее трудоемким этапом во всем процессе численного моделирования контактных систем в комплексе МАКС на основе метода контактных сил и переносных перемещений [16]. В случае промышленного внедрения комплекса требуются большие затраты на подготовку и обучение специалистов, ориентированных на подготовку данных, что является во многих случаях сдерживающим фактором несмотря на существенный выигрыш в вычислительных ресурсах на этапе численного моделирования. Реализация автоматического метода подготовки данных снимает этот недостаток. Для этого в таком методе должны быть формализованные все этапы построения сопряженных сеток, включая выделение согласованных контактных зон и обеспечения условия сопряженности узлов.

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

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

Указанная цель предусматривает решение следующих задач:

1. Произвести анализ существующих методов и алгоритмов генерации сеток конечных элементов, используемых в препроцессорах подготовки данных современных конечноэлементных комплексов;

2. Сформулировать основные требования к разрабатываемому препроцессору, исходя из особенностей моделирования контактных систем;

3. Разработать методы анализа твердотельных геометрических моделей контактных систем и дискретного геометрического моделирования отдельных тел;

4. Разработать программно-алгоритмическое обеспечение, позволяющее в автоматическом режиме проводить выделение согласованных контактных зон, построение сопряженных пар контактирующих узлов и генерацию конечноэлементных сеток в объемах моделируемых тел;

5. Реализовать препроцессор и произвести его интеграцию в комплекс моделирования и анализа контактных систем (МАКС) в качестве модуля подготовки данных.

Методы исследования. В работе использованы методы аналитической и вычислительной геометрии, многопараметрической оптимизации, линейной алгебры, функционального и объектно-ориентированного программирования.

Научная новизна диссертации. Результатами работы, отличающимися научной новизной, являются:

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

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

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

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

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

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

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

Апробация работы проведена: на научных конференциях Воронежского государственного технического университета в 2000-2003 годах; на V международной электронной научной конференции «Современные проблемы информатизации в технике и технологиях», (Воронеж, 2000); на II всероссийской научно-технической конференции «Прикладные задачи механики и тепломассообмена в авиастроении» (Воронеж, 2001); на III международной конференции «Авиакосмические технологии» (Воронеж, 2002); на международной конференции «Современные сложные системы управления» (Воронеж 2003); на международной научно-технической конференции «Системные проблемы качества, математического моделирования, информационных и электронных технологий» (Сочи, 2003).

Публикации: основные результаты исследований по диссертации отражены в 8 опубликованных работах, в том числе 3 без соавторов.

Структура и объем диссертации: введение, 4 раздела, заключение, список литературы (110 наименований) и 2 приложения; материалы диссертации (без приложений) включают 140 страниц текста, 30 рисунков, 1 таблицу

Основные результаты диссертации:

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

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

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

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

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

6. С помощью реализованного препроцессора произведена процедура автоматической дискретизации контактной системы, представляющей собой модель части реального промышленного объекта: кривошипного горяче-штамповочного пресса.

ЗАКЛЮЧЕНИЕ

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

1. Аведьян А.Б., SolidWorks AP1.- универсальная платформа для разработок пользовательских приложений // САПР и графика. 2002. - № 8.

2. А.В. Александров, Б.Я. Лащеников, Н.Н. Шапошников. Строительная механика. Тонкостенные пространственные системы. М.: Стройиздат, 1983.-488 с.

3. Аргирис Дж.Г. Вычислительные машины и механика // Теоретическая и прикладная механика: Труды XIV Международного конгресса ЮТАМ. М.: Мир, 1979. - С. 15 - 100.

4. Бате К., Вилсон Е. Численные методы анализа и метод конечных элементов. М.: Стройиздат, 1982. - 448 с.

5. Беллман Р., Энджел Э. Динамическое программирование и уравнения в частных производных. М.: Мир, 1974. - 207 с.

6. Бенерджи П., Баттерфилд Р. Методы граничных элементов в прикладных науках. М.: Мир, 1984. - 494 с.

7. Бреббиа К., Уокер С. Применение метода граничных элементов в технике. М.: Мир, 1982. - 248 с.

8. Бреббиа К., Теллес Ж., Вроубел Л. Методы граничных элементов. М.: Мир, 1987.-524 с.

9. Вдовиченко А.А., Алгоритм генерации конечноэлементной сетки для криволинейных поверхностей заданных в виде NURBS сплайнов // Электротехнические комплексы и системы управления, Сб. тр. Воронеж: ВГТУ, 2003, с. 105-112.

10. Вдовиченко А.А., Алгоритмическое и программное обеспечение генерирования сопряженных ансамблей МКЭ для многотельных контактных систем // Прикладные задачи механики и тепломассообмена в авиастроении, Сб. тр. Воронеж: ЦЧКИ, 2001, с. 191-197.

11. Галлагер Р. Метод конечных элементов. Основы. М.: Мир, 1984. - 428 с.

12. Гольник Э.Р. Методы дискретной механики контактных систем и их применение в проектировании прессов: Дис. . докт. техн. наук. Воронеж, 1995. - 377 с.

13. Гольник Э.Р., Радченко И.Г. Дискретное моделирование упругих тел, контактно взаимодействующих при произвольных статических нагрузках, зазорах и натягах // Изв. вуз. Машиностроение. 1987. -№ 12.- С. 11-19.

14. Гольник Э.Р., Радченко И.Г. Обобщение смешанного метода контактных сил и переносных перемещений на класс систем произвольного числа упругих деталей. Изв. вуз. Машиностроение. 1988. - № 10. - С. 17-22.

15. Гольник Э.Р., Радченко И.Г., Лейкин М.А. и др. Программный комплекс МАКС и опыт его применения в САПР тяжелых механических прессов. Кузнечно-штамповочное производство. 1995. - № 3. - С. 19 -22.

16. Гольник Э.Р., Успехов А.А., Гундорова Н.И. и др. Модернизация программного комплекса МАКС на уровне современных требований механики авиаконструкций. Труды третьей международной научно-технической конференции, В. 2002. С. 113-120.

17. Делоне Б. Н. О пустоте сферы. Изв. АН СССР, ОМЕН. 1934, 4. 793800.

18. Джонсон К. Механика контактного взаимодействия. М.: Мир, 1989. -510с.

19. Жук Д.М. CAD/CAE/CAM системы высокого уровня для машиностроения. Информационные технологии, № 1, 1995. С. 22 26.

20. Зенкевич О.С. Метод конечных элементов в технике. М.: Мир, 1975. -541 с.

21. Зенкевич О., Морган К. Конечные элементы и аппроксимация. М.: Мир, 1986.-318 с.

22. Ильман В.М. Экстремальные свойства триангуляции Делоне. Алгоритмы и программы. Вып. 10 (88). М., 1988. 57-66.

23. Иосилевич Г.Б. Концентрация напряжений и деформаций в деталях машин. М.: Машиностроение, 1981, 224 с.

24. Конри Т., Сейрег А. Метод математического программирования в применении к проектированию находящихся в контакте упругих тел // Прикладная механика. 1971.-№2.- С. 96-102.

25. Крауч С., Старфилд А. Методы граничных элементов в механике твердого тела. М.: Мир, 1987. - 328 с.

26. Оден Дж. Конечные элементы в нелинейной механике сплошных сред.- М.: Мир, 1976.-464 с.

27. Смирнов А.Ф., Александров А.В., Шапошников Н.Н., Лащеников Б.Я. Расчет сооружений с применением вычислительных машин. М.: Стройиздат, 1964. - 380 с.

28. Смирнов А.Ф., Александров А.В., Лащеников Б.Я., Шапошников Н.Н. Строительная механика. Стержневые системы. М.: Стройиздат, 1981. -512с.

29. Современные методы расчета сложных статически неопределимых систем: Сб.статей Дж.Аргироса, С.Келси и др. / Под ред. А.П.Филина.- Л.: Судпромгиз, 1961. 876 с.

30. Стренг Г., Фикс Дж. Теория метода конечных элементов. М.: Мир, 1977.-351 с.

31. Хог Э., АрораЯ. Прикладное оптимальное проектирование. М.: Мир, 1983.-478 с.

32. Холл Дж., Уатт Дж. Современные численные методы решения обыкновенных дифференциальных уравнений. М.: Мир, 1979

33. Шапошников Н.Н., Тарабасов Н.Д., Петров В.Б., Мяченков В.И. Расчет машиностроительных конструкций на прочность и жесткость. -М.: Машиностроение, 1981. — 333 с.

34. Aftosmis M.J., Melton J.E., and Berger M.J. Adaptation and Surface Modeling for Cartesian Mesh Methods. In 12th AIAA Computational Fluid Dynamics Conference, number 95-1725-CP, San Diego, USA, 1995.

35. Argyris J.H. u.a. ASKA User's Reference Manual. JSD-Report (Stuttgard).- 1971. № 73. - 43p.

36. Bauer R. Numerische Berechnung von Kapazitaten in dreidimensionalen Verdrahtungsstrukturen. Dissertation, Technische Universitat Wien, 1994.

37. Bern M. Compatible Tetrahedralizations. In Proc. 9th Annual Symposium on Computational Geometry, pages 281-288, San Diego, USA, 1993. ACM

38. Bern M., Plassmann P. E. 2000. Mesh Generation. Handbook of Computational Geometry. Edited by J. Sack and J. Urrutia, Elsevier Science, pp. 291332.

39. Blacker T.D. and Stephenson M.B. Paving: A New Approach to Automated Quadrilateral Mesh Generation. Int.J.Numer.Meth.Eng., pp 811-847, 1991.

40. Buratynski E.K. A Fully Automatic Three-Dimensional Mesh Generator for Complex Geometries. Int.J.Numer.Meth.Eng., pp 931-952, 1990.

41. Canann S. A., Tristano J. R., and Staten M. L., An Approach to Combined Laplaciang and Optimization-Based Smoothing for Triangular, Quadrilateral, and Quad-Dominant Meshes, presented at 7th International Meshing Roundtable, 1998.

42. Castillo J.E. Mathematical Aspects of Numerical Grid Generation. SIAM, Philadelphia, 1991.

43. Cavendish J.C., Field D.A., and Frey W.H. An Approach To Automatic Three-Dimensional Finite Element Mesh Generation. Int.J.Numer.Meth.Eng., pp 329-347, 1985.

44. Chan S.K., Tuba J.S. A finite element method for contact problems of solid bodies. Int. J. Mech. SCI.- 1971. -№ 13,- P. 615-625.

45. Chand R., Haug E.J., Rim K. Analysis of Unbounded Contact Problems by Means of Quadratic Programming. J. Optimization Theory and Applications. 1976. - V. 20 (2). - P. 171 - 190.

46. Conti P. Grid Generation for Three-Dimensional Semiconductor Device Simulation. Hartung-Gorre, 1991.

47. Edelsbrunner H., Preparata F.P., and West D.B. Tetrahedrizing Point Sets in Three Dimensions. J. Symbolic Computation, pp 335-347, 1990.

48. Farm G. NURB Curves and Surfaces from Projective Geometry to Practical Use, A. K. Peters, 1995.

49. Field D.A., Laplacian Smoothing and Delaunay Triangulation, Communications in Applied Numerical Methods, vol. 4, pp. 709-712, 1988.

50. Filipiak M., Mesh Generation, Technology Watch Report, Edinburgh Parallel Computing Centre, 1996.

51. Francavilla A., Zienkiewicz O.C. A note on numerical computation of elastic contact problems // Int. J. Numer. Meth. Eng. 1975. - V. 9. -P. 913-924.

52. Freitag L. A., On Combining Laplacian and Optimization-based Mesh Smoothing Techniques, AMD Trends in Unstructured Mesh Generation, vol. 220, pp. 37-43, 1997.

53. Freitag L.A. and Ollivier-Gooch C. A Comparison of Tetrahedral Mesh Improvement Techniques. In IMRT'96, pages 87-100.

54. Freitag L.A., Ollivier-Gooch C., Tetrahedral mesh improvement using swapping and smoothing, 1997.

55. Frey P.J., Borouchaki H., and George P.L. Delaunay Tetrahedralization Using an Advancing-Front Approach. In IMRT'96, pages 31-43.

56. Gaither, J. A., 1997. A Solid Modeling Topology Data Structure for General Grid Generation, MS Thesis, Mississippi State University.

57. Gartner R. Resolution de problemes de contact elastique aves frottement en atilisant des variales nodales appopices. J. Mech. Appl. 1977. - V. 1, № 3.- p. 247-265.

58. P.L. George. Automatic Mesh Generation and Finite Element Computation. In P.G. Ciarlet and J.L. Lions, editors, Handbook of Numerical Analysis, volume IV. Elsevier, 1996.

59. P.L. George, F. Hecht, and E. Saltel. Automatic Mesh Generator with Specified Boundary. Computer Methods in Applied Mechanics and Engineering, pp 269-288, 1991.

60. George, P. L., Hecht, F., and Saltel, E., 1990, Fully Automatic Mesh Generator for 3D Domains of any Shape, Impact of Computing in Science and Engineering, 2, 187.

61. Haug E.J., Chand R., Pan K.C. Multibody Elastic Contact Analysis by Quadratic Programming. J. Optimization Theory and Application. 1977. -V. 21 (2).- P. 189- 198.

62. J. Hauser and C. Taylor. Numerical Grid Generation in Computational Fluid Dynamics. Pineridge Press, 1986.

63. W.D. Henshaw. Automatic Grid Generation. In A. Iserles, editor, Acta Nu-merica, volume 5, pages 121-148. Cambridge University Press, 1996.

64. L.R. Hermann. Laplacian-Isoparametric Grid Generation Scheme. J. of the Engineering Mechanics Division of the American Society of Civil Engineers, pp 749-756, 1976.

65. Hitschfeld, P. Conti, and W. Fichtner. Mixed Element Trees: A Generalization of Modified Octrees for the Generation of Meshes for the Simulation of Complex 3-D Semiconductor Device Structures. IEEE Trans.Computer-Aided Design, pp 1714-1725, 1993.

66. H. Hoppe, T. DeRose, T. Duchamp, J. McDonald, and W. Stuetzle. Mesh Optimization. In Computer Graphics, SIGGRAPH'93 Proceedings, pages 19-26, 1993.

67. H. Jin and R. I. Tanner, Generation of Unstructured Tetrahedral Meshes by the Advancing Front Technique, Int. J. Num. Meth. Eng. 36, pp 1805-1823, 1993.

68. B. Joe. Construction of Three-Dimensional Delaunay Triangulations Using Local Transformations. Computer Aided Geometric Design, pp 123-142, 1991.

69. В. Joe. Construction of Three-Dimensional Improved-Quality Triangulations Using Local Transformations. SIAM J.Sci.Comput., pp 1292-1307, 1995.

70. В. Joe. Three-Dimensional Triangulations from Local Transformations. SI AM J.Sci.Stat.Comput., pp 718-741,1989.

71. P. Knupp and S. Steinberg. Fundamentals of Grid Generation. CRC Press, 1993.

72. Lattuada R., A triangulation based approach to three dimensional geoscien-tific modelling, Ph.D Thesis, 1998.

73. E. Leitner, W. Bohmayr, P. Fleischmann, E. Strasser, and S. Selberherr. 3D TCAD at TU Vienna. In J. Lorenz, editor, 3-Dimensional Process Simulation, pages 136-161, Wien, 1995. Springer.

74. A. Liu and B. Joe. Quality Local Refinement of Tetrahedral Meshes Based on Bisection. SIAM J.Sci.Comput., pp 1269-1291, 1995.

75. Levy, H. and Lessman, F. Finite Difference Equations. New York: Dover, 1992.

76. R. Lohner, A Parallel Advancing Front Grid Generation Scheme, IMRT'2000, Юр.

77. R. Lohner and P. Parikh. Three-Dimensional Grid Generation by the Advancing Front Method. Int.J.Numer.Meths.Fluids., pp 1135-1149, 1988.

78. D.L. Marcum, 1996b. Unstructured Grid Generation Components for Complete Systems. 5th International conference on Grid Generation in Computational Fluid Simulations, Starkville, MS.

79. D.L. Marcum. Unstructured Grid Generation Using Automatic Point Insertion and Local Reconnection. The Handbook of Grid Generation, edited by J.F. Thompson, B. Soni and N.P. Weatherill, CRC Press, p. 1-18, 1998.

80. D.L. Marcum and N.P. Weatherill. Unstructured Grid Generation Using Iterative Point Insertion and Local Reconnection. In 12th AIAA Applied Aerodynamics Conference, number 94-1926, Colorado Springs, USA, 1994.

81. C.W. Mastin and J.F. Thompson. Quasiconformal Mappings and Grid Generation. SIAM J.Sci.Stat.Comput.,pp 305-316, 1984.

82. D.J. Mavriplis. An Advancing Front Delaunay Triangulation Algorithm Designed for Robustness. J.Comput.Phys., pp 90-101, 1995.

83. S. Mtiller, К. Kells, and W. Fichtner. Automatic Rectangle-Based Adaptive Mesh Generation Without Obtuse Angles. IEEE Trans.Computer-Aided Design, pp 855-863, 1992

84. Murakami S., Ymada Y. Effects of hydrostatic pressure and material anisotropy on the transient creep of thick-walled tubes // Int. J. Mech. Sci. -1974. V. 16, № 3. - P. 145 - 208.

85. NISTIR 4412, Initial Graphics Exchange Specification (IGES), Version 5.0, U. S. National Institute of Standards and Technology, 1990.

86. Owen S. J., Non-Simplicial Unstructured Mesh Generation, Ph.D Thesis, 1999.

87. V. Parthasarathy and S. Kodiyalam, A constrained optimization approach to finite element mesh smoothing, Finite Elements in Analysis and Design, vol. 9, pp. 309-320, 1991.

88. J.W. Peterson. Tessellation of NURB Surfaces. In P.S. Heckbert, editor, Graphics Gems IV, pages 286-320. Academic Press, 1994.

89. M. Radi and S. Selberherr. Three-Dimensional Adaptive Mesh Relaxation. In Meyer and Biesemans. pp 193-196

90. Ruppert J. and Seidel R. On the Difficulty of Tetrahedralizing 3-Dimensional Non-Convex Polyhedra. In Proc. 5th Annual Symposium on Computational Geometry, pages 380-393. ACM, 1989.

91. Ruppert J. and Seidel R. On the Difficulty of Triangulating Three-Dimensional Non-Convex Polyhedra. Discrete & Computational Geometry, pp 227-253, 1992.

92. Z.H. Sahul, R.W. Dutton, and M. Noell. Grid and Geometry Techniques for Multi-Layer Process Simulation. In Selberherr et al., pages 417-420.

93. H. Samet. The Design and Analysis of Spatial Data Structures. Addison-Wesley, 1990.

94. N. Sapidis and R. Perucchio. Delaunay Triangulation of Arbitrarily Shaped Planar Domains. Computer Aided Geometric Design, pp 421-437, 1991.

95. R. Schneiders, R. Schindler, and F. Weiler. Octree-based Generation of Hexahedral Element Meshes. In IMRT'96, pages 205-215.

96. R. Seidel. Constrained Delaunay Triangulations and Voronoi Diagrams with Obstacles. Technical Report 260, Inst, for Information Processing, Graz, Austria, 1988.

97. M.S. Shephard and M.K. Georges. Automatic Three-Dimensional Mesh Generation by the Finite Octree Technique. Int.J.Numer.Meth.Eng., pp 709749, 1991.

98. S.P. Spekreijse. Elliptic Grid Generation Based on Laplace Equations and Algebraic Transformations. J.Comput.Phys., pp 38-61, 1995.

99. V. Srinivasan, L. Nackman, J. Tang, and S. Meshkat. Automatic Mesh Generation Using the Symmetric Axis Transformation of Polygonal Domains.

100. Proc.IEEE, pp 1485-1501, 1992.

101. J. Tautges and S.A. Mitchell. Whisker Weaving: Invalid Connectivity Resolution and Primal Construction Algorithm. In IMRT'95, pages 115-127.

102. Teng S. H., Coarsening, Sampling, and Smoothing: Elements of the Multilevel Method, 1998.

103. J.F. Thompson. A General Three-Dimensional Elliptic Grid Generation System on a Composite Block Structure. Computer Methods in Applied Mechanics and Engineering, pp 377-411, 1987.

104. J.F. Thompson, Z.U.A. Warsi, and C.W. Mastin. Numerical Grid Generation. North Holland, 1985

105. Yamakawa S. and Shimada K., High quality anisotropic tetrahedral mesh geheration via ellipsoidal bubble packing, IMRT'2000, 11 pages.