Методы метрической геометрии и анализа на неголономных многообразиях в задачах механики сплошных сред (стр. 2 )

Рассмотрим одномерную гиперболическую задачу

(19)

Где в качестве переменной выбирается, скажем , а

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

(20)

где . Для упругих потенциалов общего вида эти матрицы выписаны в работе (Годунов, Пешков, 2008).

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

(21)

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

Построение представления (21) в работе (Годунов, Пешков, 2008) с использованием внутренней энергии вида (16) требует построения сингулярного разложения матрицы . При этом была выписаны в явном виде собственные значения и собственные вектора матрицы , что потребовало использования системы символьных вычислений MAPLE. Итоговые расчетные формулы оказались весьма сложными и громоздкими. Выпишем матрицы для потенциала (18). Поскольку в данном случае в качестве неизвестных можно ограничиться величинами , т. е. первым столбцом матрицы , то сокращения выкладок можно полагать, что размерность этих матриц равна . Таким образом их можно записать в блочном виде

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

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

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

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

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

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

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

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

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

Составим матрицы и

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

Матрица Якоби афинного отображения записывается как .

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

(22)

Здесь - внутренняя энергия деформации.

Численный эксперименты с построением трехмерных отображений проводились с внутренней энергией следующего вида

(23)

Где член вида

(24)

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

Член вида

(25)

представляет собой обезразмеренную газодинамическую часть внутренней энергии. Его можно интерпретировать как меру искажения объема. Весовой коэффициент играет роль модуля всестороннего сжатия.

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

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

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

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

Тогда уравнение стационарноcти функционала в вершине можно записать так:

			(26)
			(27)

Эта система состоит из уравнений, и им соответствуют переменных, из которых составлен вектор .

Пусть обозначает приращение в точке . Линеаризуя уравнение (27), получаем следующее уравнение для

Таким образом, если лежит на границе области , то , и можно представить как линейную комбинацию векторов :

(28)

где - произвольные коэффициенты. Иными словами, равенство (28) означает, что допустимо движение граничных узлов только вдоль касательной плоскости к границе.

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

(29)

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

(30)

где величину можно рассматривать как погрешность задания геометрии.

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

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

	(31)
	(32)

Обозначим через матрицу размера , первыми двумя столбцами которой являются векторы , вычисленные в точке , а последний столбец равен нулю. Если же индекс , то положим .

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

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

	(33)
	(34)

Равенство (34) можно записать как

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

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

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

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

Вариационный метод можно использовать и в том случае, когда алгебраический объем некоторых тетраэдров в начальной сетке равен нулю или отрицателен. При наличии подобных ``вывернутых'' тетраэдров оказался эффективным метод распутывания сеток, предложенный в работе (Гаранжа, Капорин, 1999) . Идея этого метода основана на том, что детерминант в знаменателе формулы (23) заменяется на величину

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

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

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

Рис. 1

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

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

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

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

1.2 Уравнение Бюргерса в механике сплошных сред и неголономных многообразиях

Исследованы свойства поверхностей уровня непрерывно горизонтально дифференцируемых (hc-дифференцируемых) отображений из пространства Карно — Каратеодори M такого, что dim HgM = dim TgM – 1 = N в каждой точке g ÎM, в евклидово пространство размерности N. Данная работа является подготовительной к общей проблеме описания спрямляемых многообразий на пространствах Карно — Каратеодори. Её результаты обобщают некоторые результаты работ [1] и [2], в которых исследованы поверхности уровня hc-дифференцируемых отображений f : H1 → R2.

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

Теорема 1. Пусть отображение f : M → RN hc-дифференцируемо и его дифференциал имеет максимальный ранг. Тогда для каждой точки g Îf-1(0) найдётся окрестность U(g) такая, что f-1(0) Ç U(g) есть образ простой жордановой кривой и имеет топологическую размерность 1.

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

Теорема 2. Пусть отображение f : M → RN hc-дифференцируемо и его дифференциал имеет максимальный ранг. Тогда его поверхность уровня f-1(0) локально является образом простой кривой g: [0,1] → M и в субримановой метрике имеет хаусдорфову размерность dimH g([0,1]) = 2. При этом длина кривой может быть найдена как

H2(g) = lim||s||→0 Sk=1..n d(g(ak-1), g(ak))2,

где s= {0 = a0 < a1 < … < an = 1} и ||s|| = maxk=1..n |ak – ak-1|.

Кроме того, если частные производные отображения f гёльдеровы, то поверхность уровня f-1(0) является 2-регулярной по Альфорсу, её длина конечна, положительна и равна пределу

H2(g) = lim||s||→0 Sk=1..n d(g(ak-1), g(ak))2.

В случае, когда область определения является нильпотентной группой (группой Карно), удаётся выразить длину кривой явно в виде интеграла Римана — Стилтьеса. Пусть группа Карно G = {(x1,…,xN,z) ÎRN+1} моделируется алгеброй Ли, имеющей базис

Xi = ¶xi + Sk=1..N ckixk ¶z,

i=1,…,N,
Z = ¶z, cij = -cji, I,j=1,…,N.

Тогда верна следующая

Теорема 3. Пусть f : G → RN — hc-дифференцируемое отображение с дифференциалом максимального ранга и f(0) = 0. Тогда в некоторой окрестности нуля f-1(0) есть образ простой кривой g: [0,1] → G, и длина кривой может быть найдена, как предел

H2(g) = òg dz + lim||s||→0 Si,j=1..N cij òs xj dxi.

Если к тому же частные производные отображения f гёльдеровы, то длина кривой g положительна, конечна и равна