автореферат диссертации по информатике, вычислительной технике и управлению, 05.13.18, диссертация на тему:Высокопроизводительные методы расчёта дискретных моделей связанных систем тел

кандидата технических наук
Гетманский, Виктор Викторович
город
Волгоград
год
2013
специальность ВАК РФ
05.13.18
цена
450 рублей
Диссертация по информатике, вычислительной технике и управлению на тему «Высокопроизводительные методы расчёта дискретных моделей связанных систем тел»

Автореферат диссертации по теме "Высокопроизводительные методы расчёта дискретных моделей связанных систем тел"

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

Гетманский Виктор Викторович

Высокопроизводительные методы расчёта дискретных моделей связанных систем тел

05.13.18 - Математическое моделирование, численные методы и комплексы программ

АВТОРЕФЕРАТ диссертации на соискание учёной степени кандидата технических наук

? Л /""-! 7/1П

005061807

Волгоград - 2013

005061807

Работа выполнена в Федеральном государственном бюджетном образовательном учреждении высшего профессионального образования «Волгоградский государственный технический университет»

Научный руководитель:

Официальные оппоненты:

Ведущая организация:

доктор технических наук, старший научный сотрудник Горобцов Александр Сергеевич

Крысько Вадим Анатольевич, доктор технических наук, профессор, ФГБОУ ВПО «Саратовский государственный технический университет имени Гагарина Ю.А.», заведующий кафедрой «Математика и моделирование»

Ковалев Роман Васильевич,

кандидат технических наук,

ФГБОУ ВПО «Брянский государственный

технический университет», лаборатория

вычислительной механики,

научный сотрудник

Институт проблем точной механики и управления РАН, г. Саратов

Защита состоится «4» июля 2013 г. в 13— часов на заседании диссертационного совета Д 212.242.08 при ФГБОУ ВПО «Саратовский государственный технический университет имени Гагарина Ю.А.» по адресу: 410054, г. Саратов, ул. Политехническая, 77, корпус 1, ауд.319.

С диссертацией можно ознакомиться в библиотеке ФГБОУ ВПО «Саратовский государственный технический университет имени Гагарина Ю.А.»

Автореферат разослан «30» мая 2013 г.

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

Терентьев А. А.

Общая характеристика работы

Актуальность работы. В практике современных научных и инженерных расчётов динамики конструкций широко используются математические модели динамики связанных систем тел, учитывающие различные физические процессы в отдельных телах, например, тепловое или напряжённо-деформированное состояние. Исследование динамики машин с помощью таких математических моделей проведено в работах К. В. Фролова, А. С. Горобцова, В. Г. Бойкова, Д. Ю. Пого-релова, М.Д. Перминова, A.B. Синева, А.П. Гусенкова, М. Blundell, D. Harty, Т. Gillespie, D. Negrut, Д.Ю. Погорелова, Г.В. Михеева, A.A. Юдакова, R. Craig, М. Bampton, A. Shabana, J. Ambrosio, О. Bachau, C.K. Карцова, Д. А. Ямпольского и др. В случае расчётных областей произвольной формы используется, в частности, метод конечных элементов, в котором динамические модели записываются в виде системы обыкновенных дифференциальных уравнений, редуцируемой для численного решения собственными функциями (метод Крэйга-Бэмптона). Метод реализован в промышленных программных пакетах инженерного анализа (например, ADAMS, DADS, UM). Его недостатками являются трудности определения граничных условий, а также невысокая эффективность распараллеливания вычислений при программной реализации метода.

Таким образом, разработка альтернативных методов представления математических моделей континуальных сред и их решения, позволяющих обойти трудности при задании граничных условий и ориентированных на использование параллельных вычислений, является актуальной. Одним из таких методов является использование дискретноэлементных моделей, рассмотренное в работах A. Tassora, D. Negrut, К. Anderson, P. Fisette, W. Prescott, M. Arnold, O.H. Дмитроченко, Д.Ю. Погорелова, A.M. Кривцова. Результаты максимального ускорения параллельного расчёта в 2-7 раз, полученные в ряде работ на моделях реальных технических объектов, и практически отсутствие публикаций на эту тему в отечественной литературе свидетельствуют об актуальности создания масштабируемых на большое число процессов методов параллельного расчёта описанных моделей. Разработка технологий и программного обеспечения распределённых и высокопроизводительных вычислительных систем соответствует одному их пунктов перечня критических технологий РФ, что также подтверждает актуальность описанной проблемы.

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

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

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

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

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

3. Создание эффективного метода синхронного параллельного расчёта комплексных моделей.

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

5. Расчёт различных моделей с использованием программного комплекса и высокопроизводительных ЭВМ для верификации и апробации разработанной технологии.

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

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

Научная новизна:

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

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

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

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

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

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

Реализация результатов. В рамках НИР по ФЦП "Развитие научного потенциала высшей школы" (2009-2011 гг., проект № 2.1.1/1423) реализовано вычислительное ядро программного комплекса, основанное на расчётном модуле пакета

моделирования динамики систем тел ФРУНД, позволяющее моделировать физические процессы в отдельных телах. Проведён ряд вычислительных экспериментов по анализу динамического состояния элементов подвески машин повышенной проходимости по 4 договорам с промышленными предприятиями (2009-2012 гг.). Получено свидетельство о регистрации программы для ЭВМ, имеется акт внедрения.

На защиту выкосятся следующие основные результаты и положения:

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

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

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

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

5. Разработанное ядро программного комплекса позволяет подключать различные реализации вспомогательных расчётных модулей для расчёта физических процессов теплопроводности и деформации в отдельных телах и проводить параллельный расчёт на кластере или многопроцессорной ЭВМ.

6. При расчёте моделей получен предел ускорения расчёта в 10-12 раз. Метод определения динамических напряжений даёт соответствие результатам расчёта методом конечных элементов. Модель теплового состояния позволяет получить совпадение расчёта динамики нагревания и распределения температуры с экспериментальными данными.

Публикации. По теме диссертации опубликовано 25 печатных работ, из них 6 статей в журналах, рекомендованных ВАК РФ [1-6], 1 свидетельство о регистрации программы [15], 8 статей и тезисов докладов в сборниках трудов международных научных конференций [7-14], 10 статей и тезисов докладов в сборниках трудов прочих научных конференций и семинаров.

Апробация работы. Результаты диссертационной работы представлены на международных конференциях: PARENG 2009 [7], IEEE ICMA 2009 [8], "Информационные технологии в образовании, технике и медицине" (Волгоград, сентябрь 2009) [9], IMSD 2010 [10], СКТ 2010 [11], ПаВТ 2011 [12], Parallel and High-Periormance Computing and Simulation 2012 (Амстердам, апрель 2012) [13], MMTT25 (Волгоград, май 2012) [14].

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

дискретноэлементные модели для расчёта физических процессов в телах, исследована эффективность метода на различных вычислительных системах.

Структура и объем диссертации. Диссертация состоит из введения, 4 глав, заключения и библиографии. Общий объем диссертации 155 страниц текста, включающего 55 рисунков и 13 таблиц. Библиография включает 127 наименовании.

Содержание работы

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

положения. „

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

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

вспомогательными моделями.

В качестве примера моделируемого объекта, на основе которого строятся связанные модели, приведена система подрессоривания автомобиля (рис. 1). Основной

f0-fj - сосредоточенные СИЛЫ Корпус Газ

QîôHicôj [Тпр][Верхний]

f0-f3 - распределенные силы

_..........

/|рама|

2 А Поршень а

чШшипттт v ш™ Дорожный профиль

1 ^Йижнийм/р^

V^vi.esii'.vjU

Фрагмент основной модели системы тел

Вспомогательная модель Вспомогательная модель ГПР нижнего рычага

Рис. 1. Типовая расчётная схема моделируемого объекта

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

находится под действием переменных сил, возникающих в связях этого тела с другими телами. На рис. 1 показано воздействие от колеса (Р) и силы взаимодействия между телами {£■), действующие на нижний рычаг. В результате работы диссипативных сил также возникают тепловые потоки (<3). Силы, действующие на тело, приводят к возникновению в нём напряжённо-деформированного состояния,а тепловые потоки вызывают изменение температуры. В модели НДС использовано допущение о равномерном распределении силы по поверхности контакта в шарнире, связывающем два тела.

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

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

/ ^ V Ьп фв=0 у/0),А„/0) !

<..........- 1 ч \.с У <-.....*

V ^^

Гв_

<р(0 |

)Ао,{1)

Рис. 2. Связь дискретных элементов

Рис. 3. Присоединение упругого тела

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

Граница подвода тепла

С

:

— —

И \контрольный объём

(сЫх,

Рис. 4. Расчётная область для модели нестационарной теплопроводности

Рис, 5, Изменение компоненты х плотности теплового потока в контрольном объём

При моделировании нестационарной теплопроводности в твёрдом изотропном

упруго-демгфи-рующая связь

упругое тело

теле дискретные элементы, формирующие расчётную область, являются контрольными объёмами. Каждый дискретный элемент осуществляет теплообмен с соседними с ним по 6 граням и с окружающей средой. На рис. 4 показано двумерное сечение расчётной области координатной плоскостью. Теплопередача происходит из-за изменения компонент вектора плотности теплового потока (рис. 5).

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

ГМу-Ю£А = е, (1)

\ ЕЭ(у) = О

где М — матрица инерции, g — вектор силовых воздействий, у — вектор координат дискретных элементов, у — вектор координат дискретных элементов, 0 = {й (у),...,йк(у)}т — матрица-столбец, определяющая уравнения связей, Т>у — якобиан, А — вектор множителей Лагранжа, соответствующий реакциям на границах.

Математическая модель НДС описывается уравнением

Му + В^ф^, (2)

где использовано приближённое представление компонент вектора А в виде Ai = = к' = {ш1'-• ■ 1шк}Т' - Функции сил реакций в по-

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

Погрешности по перемещениям каждой к-й связи дискретных элементов описываются системой уравнений

4 = (У{ + АС

•О)

-(У^ + Ащ'д

.(¿Ь

)-Ь0

(3)

где

¿¡к _ подвектор трёх линейных компонент погрешностей вектор-функции податливой связи <1к, ¿1 — подвектор трёх угловых компонент погрешностей, %,] —

индексы связанных дискретных элементов, г^ — радиус-вектор, задающий точку А присоединения связи, верхний индекс (г) обозначает систему координат (СЮ) вектора, А-- — матрица поворота для преобразования из С^ в СЮ, Ьо - начальное расстояние между точками А и В, у< — подвектор координат поступательного движения дискретного элемента с индексом ¿, в(г1,г2) — оператор вычисления трёх углов между проекциями векторов г^ на координатные плоскости уг,хг,ху, ф0 — начальные углы. Для используемых дискретных элементов <р0 = О,Ь0 = 0.

Компоненты симметричного тензора деформаций определяются как шестимерный вектор $ = (£",£") из трёхмерных подвекторов нормальных и сдвиговых деформаций: < е = <

где к — длина ребра кубического дискретного элемента. Переход от деформаций к

напряжениям с использованием параметров Ламе ф Е/(2(1 + и)) выражается формулой

/ о XX \

аУУ Огг Оуг

°хг \ )

/2 ц +

Ф

Ф ■Ф О О О

ф

2^ + ф ф

О О О

ф ф

2 ц + ф О О О

\

о\

о

о

о

о

I £хх

2?»« V У

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

= сткз, (6)

где я = Л.2 — площадь сечения дискретного элемента координатной плоскостью, Ь — шаг дискретизации. Силы демпфирования = (лу<гг,\уА1) определяются через относительные скорости по поступательным и вращательным направлениям:

= ¿(У{ + х г«)) - (у< + х Г£>)) ,

где — вектор угловой скорости дискретного элемента г, са и ё — коэффициенты демпфирования, соответствующие поступательным и вращательным направлениям. Таким образом, функции = + ^ полностью определяются выражениями (3) - (7).

Внешние силы приводятся к системе координат тела (СК1) из глобальной (СКО) и системы координат присоединённого тела (СКЗ) по формуле Г?га), если в СК1, А01А02^-3). если в СК2' (8)

=

ог-^огЧ ,(1)|г(0)

|^|,г(1' = поггп(АТ1А02р^-р51))

где — вектор силы, используемый во вспомогательной модели, /; — вектор силы, используемый в основной модели.

Для векторов поступательных и вращательных ускорений дискретного элемента г в подвекторе а4 = вектора а(7.) с учётом неинерциальности системы отсчёта тела получены выражения в требуемых СК:

= А та№ + х г'1» + (а,!1» х х г'1')) + 2(^11> х V«) ,

^А^ + и^хи,'2'). ^

где г2 - радиус-вектор в системе отсчёта опорного тела, задающий положение центра масс дискретного элемента, СК1 относится к телу, а СК2 к дискретному элементу, ш — угловая скорость, е — угловое ускорение, а — ускорение, V — скорость, кинематические параметры относятся к опорному телу, если нижний индекс равен 1, к дискретному элементу, если 2.

Компонент правой части вектора g уравнения (2) вычисляется по формуле

В = - Ма(4) + зЛу.у), iE.Fi , (10)

где Г.,- — множество индексов дискретных элементов на поверхности, к которой приложена равномерно распределённая реакция — мощность этого множе-

ства б(у,у) — малые силы стабилизации, удерживающие положение дискретных элементов относительно тела.

Математическая модель нестационарной теплопроводности в теле, построенная на основе дискретноэлементной аппроксимации, описывается системой дифференциальных уравнений:

й = 7Д(и) + ^(у, ¿) , (11)

где 7 = {ср)~1 А — коэффициент температуропроводности, с — удельная теплоёмкость, А — коэффициент теплопроводности, ¿}(у, ¿) = (ср)~1ц(у. £) Оператор Лапласа в уравнении (11) в дискретной форме вычисляется методом конечных разностей с использованием семиточечной схемы с поправками на краях расчётной области. В работе получены выражения для учёта граничных условий второго и третьего рода в правой части уравнения (11):

•/ л_ (р{срЬ3п)-\ если У(6Г2, г = 1,2,....я (12)

1> а(и(у4, г) - ие)(сф)-\ если у< е Г3 , где а — коэффициент теплоотдачи, ие — вектор из одинаковых значений температуры окружающей среды, yi — подвектор вектора у с координатами г-го дискретного элемента. Форма границ Г< определяется конструкторской геометрией и набором дискретных элементов, попадающих на граничные поверхности. Рассеиваемая мощность р в демпфере вычисляется как

Р=(Га-Оп, (13)

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

Метод расчёта рассмотренных моделей состоит из следующей последовательности действий:

1. Исходная модель динамики системы тел со связями дополняется вспомогательными моделями НДС и нестационарной теплопроводности для некоторых тел.

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

3. Для модели НДС задаются параметры материала: модуль упругости (Е), коэффициент Пуассона (А), плотность (р). Для модели теплового состояния задаются теплоёмкость (с), коэффициенты теплопроводности (7) и теплоотдачи (а), а также температура окружающей среды (ие). Задаются параметры численного интегрирования: время интегрирования г, шаг интегрирования Дг.

4. Численное решение задачи Коши основной модели, описываемой уравнением вида (1), вспомогательных моделей НДС (2) и теплопередачи (11) методом Рунге-Кутты 4 порядка происходит параллельно с синхронизацией зависимых параметров. В качестве начальных условий приняты исходное расположение и ориентация дискретных элементов для модели НДС и температура окружающей среды для модели теплопередачи.

5. На каждом шаге интегрирования системы вида (1) в основной модели определяются внешние для вспомогательных моделей силы {, кинематические параметры опорного тела А01,г1,у1,а1,е1,а>1.

6. На каждом шаге интегрирования уравнения (2) в модели НДС рассчитывается вектор сил реакций в податливых связях у/^ по формулам (3)-(7). Связываемые параметры учитываются в формулах (8)-(10) при вычислении вектора правых частей g уравнения (2). Исследуемое поле напряжений «т в каждый момент времени определяется приведёнными напряжениями, рассчитанными по гипотезе Мизеса с использованием компонент тензора напряжений, вычисленных по формуле (5).

7. На каждом шаге интегрирования уравнения (11) в модели теплопередачи вычисляются выражения (12), учитывающие граничные условия и включающие связываемый параметр, вычисляемый по формуле (13). Исследуемое поле температур является вектором решений и уравнения (11).

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

Для распараллеливания вычислений используется декомпозиция расчётной области, в связи с чем в работе используются три типа моделей (рис. 6): ос-

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

Декомпозиция 1 2 3..----.2 3 4

Подмодели^.

Основная модель

К

VВспом огательн ая \модель 1 Вспомогательная модель 2

п, п,

.Ш.111 1

\

Подмодели

П3,

131

тш

Рис. 6. Зависимости моделей

Уровень подмоделей требуется для балансировки вычислительной нагрузки при распределении вычислений по вычислительным ресурсам. С точки зрения реализации он наиболее сложен, кроме того, декомпозиция вносит избыточность из-за перекрытия и в общем пространственном случае создаёт большое количество зависимостей между подмоделями. При декомпозиции используются перекрывающие границы для сохранения разрезаемых связей. У элементов, принадлежащих границам разреза Г^- в подобластях П?, отсутствует часть связей, поэтому вектор неизвестных для этих элементов вычисляется в подмодели некорректно. При использовании явного метода численного интегрирования некорректные значения вектора неизвестных можно заменить значениями вычисленными для этих элементов с полным набором связей в смежных подмоделях. Пересылка значений неизвестных между смежными подмоделями необходима перед всеми этапами алгоритма численного интегрирования, использующими вычисления функций правых частей. Для используемого метода Рунге-Кутты 4-го порядка такими этапами являются вычисления 4 коэффициентов.

С учётом особенности алгоритма декомпозиции выведена формула (14), которая показывает условие прекращения балансировки нагрузки:

к

м,

(14)

них - + тах(|Г^.|(1 + «,)) + £ £ |Г^|) < б ■ I , \ 1 1=1 3 з=1 ]

где г — номер подмодели, ] — номер граничной области в подмодели. К — количество подобластей, М, — количество граничных областей в подобласти г, Л'; — количество узлов в подобласти г, Г — множество узлов граничной области, р — количество параллельных процессов, 1,к,а — коэффициенты. Для получения выражения (14) использованы допущения о линейной зависимости вычислительной сложности алгоритмов от размерности модели и времени обмена от размерности границ. Идеально равномерной балансировки вычислений в параллельных процессах добиться нельзя, поэтому в разработанном алгоритме балансировки устанавливается пороговое значение 5, которое ограничивает минимальное приращение времени декомпозицией относительно общего времени параллельного расчёта Экспериментально получено, что <5 е [0,05; 0,1]. При декомпозиции происходит выделение части времени расчёта в параллельные потоки,.что позволяет выполнить больше итераций (рис. 7). Сплошным цветом закрашены на временной шкале выполнения моделирования интервалы времени расчёта основной модели, штриховкой — вспомогательной. Стрелками обозначена пересылка данных при синхронизации параллельных вычислений.

до декомпозиции

Рис. 7. Временная диаграмма до и после декомпозиции

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

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

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

Программный комплекс состоит из ядра, графического интерфейса и постпроцессора для графиков. Автором разработано вычислительное ядро программного комплекса (рис. 8, 9), написанное на языке С++ с использованием объектно-ориентированного подхода.

В методе моделирования применяется динамика систем тел со связями, поэтому программный комплекс основан на модулях пакета ФРУНД для символьной генерации и решения системы уравнений (1) в абсолютных координатах явными численными методами.

Для вспомогательных моделей разработаны отдельные расчётные модули на языке С++, реализующие интерфейсы "1ТЬегша13о1уег" для расчёта нестационарной теплопроводности, "^гевББо^ег" для расчёта динамического напряжённо-деформированного состояния (рис. 8).

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

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

В главе приводятся Методы распараллеливания вычислений с использованием таблицы синхронизации параллельных вычислений — матрицы смежности графа взаимодействия параллельных процессов, в ячейках которой записаны адреса для

9 Инт. командной строки

«component»

Слой парсера командной строки консольного ядра

а

Инт [компиляторов ДУ!нт. генератора Инт. декомпозиции А Инт ^решателя Хйнт"

О V V V V янив,

«component»

Слой генерации

] Обертка над J компиляторами

а

I Обертка

1 генератора ФРУНД

^Препроцессор ^исходных кодов

Препроцессор макросов

Г\

Командная строка компилятора

favmodel

к

fmodel

fcrmodel

pthrewjs MPI

gads ^

«component»

Слой решателя

] Балансировщик 1 нагрузки

1 Модуль 1 декомпозиции

] Модуль синхронизации J параллельного расчета

а

Обертка

решателя ФРУНД

Обертка библиотеки декомпозиции

Обертка решателей для тел

Ч$>

Слой работы с

! Калькулятор

I Транслятор ' модели

К

FModol

«component»

моделью

<

fcrmodel

ISolver

metis

аниматора

АчИнт. создания ^графиков

«component»

Слой обработки результатов

а

Обертка постпроцессора ФРУНД

Провайдер анимации

IStregsSolver ' IThermalSotver

Провайдер формата модели

Провайдер сетки

Провайдер результатов

fanaz

С

Mv

а

Рис. 8. Диаграмма компонентов (UML) ядра программного комплекса

Рис. 9. Диаграмма развёртывания (UML) междисциплинарных расчётных модулей

приёма и отправки данных. Разработанный алгоритм параллельного сопряжённого расчёта рассмотренных моделей обновляет данные во вспомогательных моделях, вычисляемые по формулам (8), (9) и (13) на каждой итерации основного решателя, в соответствии с адресами таблицы синхронизации. Обмен производится при вычислении новых значений в основном решателе, т.е. с шагом интегрирования Дг, обмен при декомпозиции - с шагом вычисления коэффициентов в методе Рунге-Кутты, равным четверти от шага интегрирования (Atsub/4). Таблица динамически перестраивается для двух уровней распараллеливания вспомогательных моделей и подмоделей.

Для работы с распределённой памятью используется MPI (интерфейс обмена сообщениями). Распределённый высокопроизводительный расчёт производится на уровне нескольких решателей. Для организации пересылки зависимых параметров используются операции коллективного обмена данными. Для распараллеливания при декомпозиции расчётной области реализована синхронизация потоков стандарта POSIX с обменом данными через общую память.

В четвёртой главе рассматриваются примеры использования разработанного программного комплекса. Рассмотрены 3 модели транспортных средств: автомобиль повышенной проходимости с колёсной формулой 4x4, грузовой автомобиль с 6 независимыми подвесками и гусеничная платформа. Для первых двух автомобилей заданы вспомогательные модели НДС нижнего и верхнего рычага каждой подвески и нестационарной теплопроводности штока и цилиндра ГПР. Таким образом, для автомобиля повышенной проходимости заданы 16 вспомогательных моделей, для грузового автомобиля — 24. Для гусеничной платформы заданы 28 вспомогательных моделей нестационарной теплопроводности для штоков и цилиндров ГПР.

Рис. 10. Конструкторская геометрия и результаты моделирования динамики автомобиля,

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

Для моделирования использован разработанный программный комплекс, описанный в главе 3. В одном из примеров мультифизического моделирования определены поля температур в штоке и корпусе амортизатора и поле напряжений в нижнем рычаге подвески в каждый момент времени (рис. 10) с использованием расчётной методики, описанной в главе 2 при следующих значениях параметров: шаг пространственной дискретизации (к = 2мм), дающий около 300 тысяч дискретных элементов, модуль упругости (В =2,1е11 Н/м2), коэффициент Пуассона (А =0,28), плотность (р =7700 кг/м3), время интегрирования I = Юс, Д£ = 0,0001с. Нагрузки на рычаг соответствовали силам реакций в статическом положении подвески автомобиля, возникающих в шарнирах крепления рычага к раме, ступице колеса и ГПР. В данном расчёте физическая нелинейность не учитывалась, а геометрическая определялась точными выражениями для деформаций элемента по формуле (4), что необходимо для учёта больших угловых перемещений рычага (повороты на угол до 40°).

Для верификации модели НДС проведено сравнение полученного распределения напряжений в разработанном программном комплексе с результатами статиче-

ских расчётов с помощью МКЭ в промышленном пакете моделирования ЭоПс^ог!« при тех же параметрах модели (рис. 11), при этом использовался вариант закрепления, который дал наименьшее расхождение с результатами расчёта методом дискретных элементов, так как для МКЭ расчёта рассматриваемого рычага задать граничные условия корректно невозможно.

Замеры напряжений вдоль направляющих рычага для обоих методов дали качественное совпадение и максимальное численное расхождение на 30-40% в некоторых точках. Такой результат считается приемлемым для предварительного анализа НДС конструкции, так как в МКЭ погрешности получаемых результатов аналогичные.

Напряжения по направляющей левого и правого плача рычага (МПа)

Рис. 11. Сравнение расчёта НДС рычага подвески

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

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

Для верификации модели теплового состояния проведены сравнения с экспериментом с использованием тепловизора и стойки для исследования работы амортизатора (рис. 12). Моделирование проведено при следующих значениях параметров материала амортизатора: плотность (р = 7700 кг/м3), теплоёмкость (с= 7200 Дж/(кг-К)), коэффициент теплопроводности (7 = 46 Вт/(м2-К)) и коэффициент теплоотдачи (а = 25 Вт/(м2-К)). Температура окружающей среды соответствовала условиям проведения эксперимента (ие =24°С). Результаты моделирования практически полностью соответствуют эксперименту по динамике изменения максимальной температуры и по распределению температуры вдоль направляющей цилиндра амортизатора с отклонением не более, чем на 4%.

Для анализа скорости и масштабируемости метода использованы показатели ускорения S и эффективности Е расчёта при использовании N параллельных потоков: Se[N) — tau/tN, Ee(N) = Se(N)/N, где tau - время последовательного выполнения, t/v - время параллельного выполнения N потоков. В работе описаны следующие результаты исследований быстродействия параллельного метода:

1. Сравнение ускорения расчёта при различной размерности сеток.

2. Сравнение ускорения параллельного расчёта при декомпозиции с использованием технологии ОрепМР и собственной многопоточной реализации.

Динамика нагревания

Т,......-Т,-

Распределение температуры по длине

Т...........Т,-

Результаты моделирования Тепловизор Модель

Рис. 12. Сравнение расчёта нестационарной теплопроводности в амортизаторе

3. Ускорение параллельного расчёта вспомогательных моделей на кластере с MPI и на многопроцессорной ЭВМ с неоднородным доступом к памяти (NUMA)

при использовании ОрепМР И MPI.

В первом случае результат получен в вычислительном эксперименте на кластере с 6 узлами с равномерным распределением вспомогательных решателей по узлам. Параллельное выполнение процессов на узлах реализовано с использованием MPI при распараллеливании решателей с помощью ОрепМР. Эксперимент на малом количестве узлов (до 8) показал достаточно высокую эффективность (выше 0,8) на сетках средней размерности (3 х 105 элементов) и постепенное падение эффективности для сеток малой размерности (3 х 104 элементов). Остальные эксперименты проводились с распараллеливанием по ядрам без внутреннего параллелизма решателей с использованием сеток средней размерности. Для расчёта оценочных значений ускорения и эффективности параллельного расчёта использовался коэффициент деградации а, введённый по аналогии с сетевым законом Амдала. Оценочное ускорение рассчитано по формуле

оШ)__ЬЁ.__(15)

(tSVno(N) + ti/N + tMN -I))' где ii - время последовательного расчёта, tsync - время синхронизации. Коэффициент а определён методом наименьших квадратов для кластера и ЭВМ с NUMA-архитектурой, при этом для различных моделей получено одинаковое значение, что подтверждает корректность формулы оценки ускорения (15). На рис. 13 приведены оценочные и экспериментально полученные ускорения параллельного расчёта.

о 4 8 12 16 20 24

число параллельных процессов

а)

16 20 24 28

число параллельных потоков б)

Рис. 13. Ускорение и эффективность а) на кластере'б) на ЭВМ с ЫиМА На вычислительном кластере при использовании большего количества параллельных процессов ускорение продолжает расти, однако из-за сетевой деградации

наблюдается снижении эффективности, то есть ухудшение масштабируемости. На ЭВМ с NUMA-архитектурой с помощью теста оценки производительности Unpack выявлена высокая деградация параллельных вычислений, из-за которой ускорение при количестве ядер больше 16 начинает падать. Таким образом, для ЭВМ с NUMA архитектурой предел масштабируемости метода составляет 16 ядер. Для снижения накладных расходов при использовании ОрепМР создана собственная многопоточная реализация с использованием потоков POSIX, дающая прирост двукратный ускорения (рис. 13 б). Полученный предел, ускорения в 9-12 раз справедлив для многоядерных систем с NUMA-архитектурой и вычислительных кластеров. Он теоретически преодолим при использовании вычислительных архитектур с массивным параллелизмом, например GPGPU и MIC, для которых время синхронизации и коэффициент деградации существенно меньше.

В заключении приводятся выводы по работе и основные полученные резуль-

ТЭТЫ* *>

'1. Разработаны математические модели связывания расчётных областей для

модели напряжённо-деформированного состояния, модели нестационарной теплопроводности и для модели систем тел.

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

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

4. Разработано ядро программного комплекса с использованием объектно-ориентированного подхода и современных технологий параллельных вычислений (MPI, PThread, ОрепМР). Спроектирована расширяемая архитектура с возможностью присоединения различных реализаций вспомогательных решателей.

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

Основные публикации по теме диссертации

В изданиях, рекомендованных ВАК РФ

1. Гетманский В. В., Горобцов А. С. Решение задач большой размерности в системах моделирования многотельной динамики с использованием параллельных вычислений // Изв. ВолгГТУ. - 2007. - Т. 9, № 3. - С. 10-12.

2. Гетманский В. В., Горобцов А. С., Резников М. В. Параллельное решение систем дифференциально-алгебраических уравнений большой размерности // Информационные технологии. — 2008. - № И. — С. 55-59.

3. Применение файловых операций в MPI программах для распараллеливания вычислительных задач, зависимых по данным / В. В. Гетманский, А. Е. Андреев, Д. Н. Жариков, Е. С. Сергеев // Изв. ВолгГТУ. - 2009. - Т. 6, № 6. - С. 45-47.

4. Гетманский В. В., Сергеев Е. С., Горобцов А. С. Перенос системы моделирования многотельной динамики на вычислительный кластер // Научно-технические ведомости СГПУ. — 2010. — № 3(101).- С. 93-99.

5 Решение систем дифференциально-алгебраических уравнений последовательным исключением множителей Лагранжа / В. В. Гетманский, О. В. Шаповалов, А. Е. Андреев, А. С. Горобцов // Изв. ВолгГТУ. - 2011.- Т. 3, № 10 - С. 31-33.

6 Concurrent simulation of multibody systems coupled with stress-strain and heat transfer solvers / V. V. Getmanski, A. S. Gorobtsov, Sergeev E. S. et al. // Journal of Computational Science. - 2012. - Vol. 3, no. 6. - P. 492 - 497.

В прочих изданиях

7 Getmanski V. V., Gorobtsov A. S., Andreev A. E. Solution of Multibody System Dynamics Problems on Client-Server Architectures with Pipe Data Interchange // Proceedings of PARENG2009, Pecs, Hungary, 6-8 April 2009,- Civil-Comp Press, 2009.

8 Getmanskiy V. V., Gorobtsov A. S., Klimov S. Y. Parallel inverse dynamics method for synthesis of control movement of multidimensional walking robot // Mechatronics and Automation, 2009. (ICMA 2009). International Conference. - 2009. - P. 3168 -3172.

9. Гетманский В. В., Горобцов А. С., Сергеев Е. С. Обработка сверхбольших моделей систем многих тел на параллельных вычислительных комплексах // Информационные технологии в образовании, технике и медицине : матер. Междунар. конф., 21-24 сент. 2009 / ВолгГТУ [и др.].- Волгоград, 2009.- С. 79.

10 Concurrent simulation of large-scale multibody systems using MP1 / V. V. Getmanski, A. S. Gorobtsov, E. S. Sergeev, A. E. Andreev // IMSD10, Lappeenranta, Finland, May 25-27, 2010: book of abstracts. - 2010. - P. 358-359.

11 Гетманский В. В., Сергеев Е. С., Горобцов А. С. Способы ускорения процедуры расчета при междисциплинарном моделировании // СКТ-2010: матер. Междунар. конф,-Т. 1,- Таганрог : ТТИ ЮФУ, 2010,- С. 199-222.

12 Подход к решению междисциплинарных задач на основе системы многотельнои динамики / В. В. Гетманский, Е.С. Сергеев, О. В. Шаповалов, А: С. Горобцов // ПаВТ2011: тр. Междунар. науч. конф. - Челябинск, Изд. центр ЮУрГУ, 2011. - С. 709.

13 Concurrent simulation of multibody system coupled with stress-strain and heat transfer solvers / V. V. Getmanskiy, E. S. Sergeev, A. S. Gorobtsov et al. // Book of Abstracts of YSC 2012, Amsterdam, Netherlands 2-6 April, 2012,- 2012,- P. 43-44,-URL: http://promo.escience.ifmo.ru/content/files/abs.pdf (accessed November 1, 2012).

14 Гетманский В В., Горобцов А. С. Расчет динамики систем тел и физических процессов в отдельных телах на кластере // ММТТ-25: сб. тр. XXV Междунар. науч. конф.-Т. 3,-ВолгГТУ : Саратов, 2012.-С. 171-174.

Авторские документы

15. Свид. о roc. per. программы для ЭВМ №2012611473 от 08.02.2012 г. РФ, МПК(нет). Программная платформа для междисциплинарных расчетов на базе многотельной динамики / В. В. Гетманский [и др.]; ВолгГТУ,—2012.

Гетманский Виктор Викторович

Высокопроизводительные методы расчёта дискретных моделей связанных систем тел

Автореферат

Подписано в печать 39 .ОБ".2013 г. Заказ №594. Тираж 100 экз. Печ. л. 1,0 Формат 60 х 84 1/16. Бумага офсетная. Печать офсетная. Типография ИУНЛ Волгоградского государственного технического университета 400005, Волгоград, ул. Советская, 35

Текст работы Гетманский, Виктор Викторович, диссертация по теме Математическое моделирование, численные методы и комплексы программ

ВОЛГОГРАДСКИЙ ГОСУДАРСТВЕННЫЙ ТЕХНИЧЕСКИЙ

УНИВЕРСИТЕТ

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

0420135^188

ГЕТМАНСКИЙ ВИКТОР ВИКТОРОВИЧ

Высокопроизводительные методы расчёта дискретных моделей связанных систем тел

Специальность 05.13.18 — «Математическое моделирование, численные методы и комплексы

программ»

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

Научный руководитель: доктор технических наук, старший научный сотрудник Горобцов A.C.

Волгоград - 2013

Содержание

Введение..................................................................5

1. Современные технологии сопряжённого моделирования динамики систем тел и физических процессов в отдельных телах ....................................................................13

1.1. Моделирование системы абсолютно твёрдых тел..........13

1.2. Сопряжение физически разнородных моделей с моделью системы абсолютно твёрдых тел..............................19

1.3. Модели упругих тел ..........................................23

1.4. Распараллеливания расчёта сопряжённых моделей .... 27

1.5. Метод подсистем..............................................30

1.6. Методы декомпозиции расчётной области..................31

1.7. Особенности реализации программных комплексов муль-тифизического моделирования................................34

1.8. Выводы по главе 1..............................................36

2. Методы сопряжённого моделирования с использованием параллельных вычислений ..........................................38

2.1. Распараллеливание расчёта динамики системы тел .... 38

2.2. Модель динамического напряжённо-деформированного состояния тела....................................................40

2.3. Модель нестационарного распределения тепла в теле . . 52

2.4. Алгоритм декомпозиции расчётной области................57

2.5. Уровни распараллеливания расчёта сопряжённых моделей 62

2.6. Алгоритм статической балансировки нагрузки..............67

2.7. Алгоритм синхронизации параллельного расчёта..........74

2.8. Выводы по главе 2 ............................................79

3. Программный комплекс для параллельного расчёта мульти-физических моделей, основанных на динамике систем тел . 82

3.1. Архитектура программного комплекса......................82

3.2. Архитектура ядра программного комплекса................86

3.3. Объектное представление предметной области ............93

3.4. Реализация метода параллельного моделирования .... 94

3.5. Выводы по главе 3 ......................106

4. Использование разработанных методов для решения практических задач...........................107

4.1. Область применения разработанных мультифизических моделей ............................................................107

4.2. Верификация вспомогательных моделей..........108

4.3. Модель многокатковой гусеничной платформы............114

4.4. Модель автомобиля с независимыми подвесками..........115

4.5. Модель автомобиля повышенной проходимости............118

4.6. Высокопроизводительные расчёты отдельных деталей в составе модели автомобиля ..................................119

4.7. Методика параллельных расчётов..............120

4.8. Исследование эффективности и ускорения параллельного метода моделирования........................................125

4.9. Выводы по главе 4 ............................................131

Заключение............................... 134

Список использованных источников................ 136

ПРИЛОЖЕНИЕ А. Замеры ускорения и эффективности параллельных вычислений ..................... 151

Список условных сокращений

В работе используются следующие сокращения:

АТТ — абсолютно твёрдое тело

ГПР — гидропневматическая рессора

МКЭ — метод конечных элементов

МСТ — модель системы тел со связями

НДС — напряжённо-деформированное состояние

СК — система координат

СО — система отсчёта

СЛАУ — система линейных алгебраических уравнений

ООП — объектно-ориентированный подход

ОС — операционная система

ТС — транспортное средство

ЭВМ — электронно-вычислительная машина

CAD — Computer-Aided Design - автоматизированное проектирование CAE — Computer-Aided Engineering - инженерный анализ FEM — Finite Element Method - метод конечных элементов DEM — Discrete Element Method - метод дискретных элементов MBS — Multi-body system - система тел со связями MPI — Message Passing Interface - интерфейс обмена сообщениями NUMA — Non-Uniform Memory Access - неоднородный доступ к памяти UML — Unified Modeling Language - унифицированный язык моделирования

Введение

Актуальность работы. В практике современных научных и инженерных расчётов динамики конструкций широко используются математические модели динамики связанных систем тел, учитывающие различные физические процессы в отдельных телах, например, тепловое или напряжённо-деформированное состояние. Исследование динамики машин с помощью таких математических моделей проведено в работах К. В. Фролова, А. С. Горобцова, В. Г. Бойкова, Д. Ю. Погорелова, М. Д. Перминова, A.B. Синева, А. П. Гусенкова, М. Blundell, D. Harty, Т. Gillespie, D. Negrut, Д. Ю. Погорелова, Г. В. Михеева, A.A. Юдакова, R. Craig, М. Bampton, A. Shabana, J. Ambrosio, О. Bachau, C.K. Карцо-ва, Д. А. Ямпольского и др. В случае расчётных областей произвольной формы используется, в частности, метод конечных элементов, в котором динамические модели записываются в виде системы обыкновенных дифференциальных уравнений, редуцируемой для численного решения собственными функциями (метод Крэйга-Бэмптона). Метод реализован в промышленных программных пакетах инженерного анализа (например, ADAMS, DADS, UM). Его недостатками являются трудности определения граничных условий, а также невысокая эффективность распараллеливания вычислений при программной реализации метода.

Таким образом, разработка альтернативных методов представления математических моделей континуальных сред и их решения, позволяющих обойти трудности при задании граничных условий и ориентированных на использование параллельных вычислений, является актуальной. Одним из таких методов является использование дискретноэлементных моделей, рассмотренное в работах A. Tassora, D. Negrut, К. Anderson, P. Fisette, W. Prescott, M. Arnold, O.H. Дмитроченко, Д.Ю. Погорелова, A.M. Кривцова. Результаты максимального ускорения параллельно-

го расчёта в 2-7 раз, полученные в ряде работ на моделях реальных технических объектов, и практически отсутствие публикаций на эту тему в отечественной литературе свидетельствуют об актуальности создания масштабируемых на большое число процессов методов параллельного расчёта описанных моделей. Разработка технологий и программного обеспечения распределённых и высокопроизводительных вычислительных систем соответствует одному их пунктов перечня критических технологий РФ, что также подтверждает актуальность описанной проблемы.

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

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

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

3. Создание эффективного метода синхронного параллельного расчёта комплексных моделей.

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

5. Расчёт различных моделей с использованием программного комплекса и высокопроизводительных ЭВМ для верификации и апробации разработанной технологии.

Методы исследования. Использованы методы математического мо-

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

Научная новизна:

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

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

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

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

Практическая значимость. Предложенная технология синхронного параллельного расчёта мультифизических моделей применима при ре-

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

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

Реализация результатов. В рамках НИР по ФЦП "Развитие научного потенциала высшей школы" (2009-2011 гг., проект № 2.1.1/1423) реализовано вычислительное ядро программного комплекса, основанное на расчётном модуле пакета моделирования динамики систем тел ФРУНД, позволяющее моделировать физические процессы в отдельных телах. Проведён ряд вычислительных экспериментов по анализу динамического состояния элементов подвески машин повышенной проходимости по 4 договорам с промышленными предприятиями (2009-2012 гг.). Получено свидетельство о регистрации программы для ЭВМ, имеется акт внедрения.

На защиту выносятся следующие основные результаты и положения:

1. Разработанная математическая модель обеспечивает связывание расчётных областей и параллельное решение задач определения ди-

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

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

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

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

5. Разработанное ядро программного комплекса позволяет подключать различные реализации вспомогательных расчётных модулей для расчёта физических процессов теплопроводности и деформации в отдельных телах и проводить параллельный расчёт на кластере или многопроцессорной ЭВМ.

6. При расчёте моделей получен предел ускорения расчёта в 10-12 раз. Метод определения динамических напряжений даёт соответствие результатам расчёта методом конечных элементов. Модель теплового состояния позволяет получить совпадение расчёта динамики нагревания и распределения температуры с экспериментальными данными.

Апробация работы. Результаты диссертационной работы представлены на международных конференциях: PARENG 2009 [72], IEEE ICMA 2009 [74], "Информационные технологии в образовании, технике и медицине" (Волгоград, сентябрь 2009) [14], IMSD 2010 [59], СКТ 2010 [10], ПаВТ 2011 [34], Parallel and High-Performance Computing and Simulation 2012 (Амстердам, апрель 2012) [60], MMTT25 (Волгоград, май 2012) [9].

Публикации. По теме диссертации опубликовано 25 печатных работ, из них 6 статей в журналах, рекомендованных ВАК РФ [8, 13, 35, 38, 41, 61], 1 свидетельство о регистрации программы [40], 8 статей и тезисов докладов в сборниках трудов международных научных конференций [9, 10, 14, 34, 59, 60, 72, 74], 10 статей и тезисов докладов в сборниках трудов прочих научных конференций и семинаров.

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

Краткое содержание работы.

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

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

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

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

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

неоднородным доступом к памяти.

Структура и объем диссертации. Диссертация состоит из введения, 4 глав, заключения и библиографии. Общий объем диссертации 155 страниц текста, включающего 55 рисунков и 13 таблиц. Библиография включает 127 наименований.

Глава 1. Современные технологии сопряжённого моделирования динамики систем тел и физических процессов в отдельных телах

1.1. Моделирование системы абсолютно твёрдых тел

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