Бураго Н.Г. Вычислительная механика (1185926), страница 35

Просмтор этого файла доступен только зарегистрированным пользователям. Но у нас супер быстрая регистрация: достаточно только электронной почты!

Текст из файла (страница 35)

формулы (6)-(11)).18.4. Неявные лагранжевы схемыПростейшая неявная схема для задач упругопластичности,основанная на квазилинеаризации исходных нелинейных уравненийи использующая на каждом шаге по времени одну итерацию пометоду Ньютона, имеет вид (обозначения, принятые для явных схем,здесь сохранются):M i ( v in +1 − v in ) / ∆tn = (Fin + [Ft ]in +1 ∆tn )vn +1i=vn +1*i ,( i ∈ Ω \ ΩV )( i ∈ ΩV )[ε p ]kn +1 = [ε p ]nk + (λkn +1[σ nk ] − [ε p ]kn ⋅ [L]kn +1 − ,−[LT ]nk +1 ⋅ [ε p ]kn )∆tn( k ∈ ΩC )[U ]nk +1 = [U ]nk + [σ]nk :[L]nk +1 ∆tn /[ ρ ]kn ( k ∈ ΩC ) (23.22) (19)[a p ]kn +1 = [a p ]nk + [λ ]nk +1[σ]nk :[σ ]nk ∆tn ,( k ∈ ΩC )xin +1 = xin + v in +1∆tn ,(i∈Ω )гдеMC[L]nk +1 = ∑ ∇ nkl vCn +(1k ,l )l =1λn +1kn +1t i[F ]= H ([σ ]nk :[σ]nk − [k p2 ]nk ) H ([σ ]nk :[L]nk )[σ ]nk :[L]nk +1[k p−2 ]nkNC M C= ∑∑ [g ]k =1 l =1n +1C klNB M BHɶ (i − C (k , l )) + ∑∑ [g B ]kln +1 Hɶ (i − B(k , l ))k =1 l =1[g C ]nkl+1 = ∆tn [V ]nk [σ v ]kn +1 ⋅∇ kln(23.23) (20)Cхема (19) имеет почти второй порядок точности и безусловноустойчива.

Шаг по времени, тем не менее, ограничен условиемточности(∆tn ≤ ∆ε max max ([e]nk :[e]kn )1/ 2k∈ΩC)−1(23.24) (21)где ∆ε max - максимально допустимое приращение деформации нашаге по времени. В практических расчетах можно принять203Глава 18. Методы для задач упругопластичности∆ε max = 0.1ε Y, где ε Y- деформация, отвечающая пределутекучести.Порядок вычислений по схеме (19)-(20) почти таков же, каки для явных схем. В случае невной схемы приращения координат∆xin +1 (или скорости v in +1 = ∆xin +1 / ∆tn ) определяются из решениявспомогательной линеаризованной краевой задачи, которая всоответствии с принятыми аппроксимациями представленасистемой линейных алгебраических уравнений (СЛАУ):Mi1  ∆xin +1− v in  = Fin + [Ft ]in +1 ∆tn / 2 ,∆tn  ∆tn∆xin +1 = v in +1∆tni ∈ ΩV )( i ∈ Ω \ ΩV )(22)(23.25)(23)где приращения узловых сил [Ft ]in +1 ∆tn , как видно из формул (19)(23), являются линейными функциями приращений координат∆xin +1 .Разрешающая система уравнений (22)-(23) может бытьзаписана кратко в стандартной формеAy = b(23.26) (24)где y = {∆xin +1}i =V1 , A - матрица "жесткости", b - известный векторправых частей; число неизвестных равно числу узловыхкомпонентов перемещения или скорости ( mNV , где m размерность исходной дифференциальной задачи).

ТрадиционыйМКЭ подразумевает явное формирование СЛАУ (24), а именноопределение матрицы СЛАУ ("конденсация матрицы жесткости") cпоследующей борьбой с известными проблемами минимизациипамяти, требуемой для хранения этой матрицы, оптимизации еезаполнения c целью получениия ленточной структуры матрицыпутем оптимимальной перенумерации узлов сетки, далее следуютпроблема построения алгоритма обращения матрицы, оптимизацииобменов с внешней памятью и так далее.N18.5. Безматричная реализация неявныхсхемРассматрим безматричный способ решения задачи (22)-(23)(Бураго, Кукуджанов, 1988). Матрицу жесткости A специальноформировать и запоминать не потребуется, если воспользоваться204Глава 18. Методы для задач упругопластичностикаким-либо итерационнымневязки уравнений (23.26):методомрешения,использующимNV1  ∆xin +1g = Ay − b =  M i− v in  − Fin − [ Ft ]in +1 ∆tn / 2  ∆tn  ∆tn i =1(23.27)Действительно, вычисление невязок g (y ) для любого заданногоприближенного решения y можно выполнить, не используя явноматрицу жесткости.

Пусть вектор приближенного решенияy = {∆xin +1}iN=V1 задан. Тогда по формулам (19) - (20) дляk = 1, 2,..., N C последовательно вычисляем величины: Lnk+1 , λkn +1 , ...,Fkn ,[Ft ]nk +1 .Окончательно вектор невязки определяетсянепосредственно соотношениями (23.27).Столь же просто определяется и однородная часть невязкиgɶ (y ) = Ay :NV1  ∆xin +1 n +1gɶ = Ay =  M i − [ Ft ]i ∆tn / 2  ∆tn  ∆tn  i =1(23.28)Ясно, что среди множества итерационных методов желательновыбрать наиболее эффективный и простой в реализации. На нашвзгляд таким методом является метод сопряженных градиентов(Хестенес, Штифель, 1952).

Этот итерационный алгоритм стартует снекоторого начального приближения к искомому решению y :g 0 = Ay 0 − bp0 = g0(23.29)и далее подразумевает следующие вычисления для итерацийs = 1, 2,... :α s = (g s −1 ⋅ g s −1 ) /( Ap s −1 ⋅ p s −1 )y s = y s −1 − α s p s −1g s = g s −1 − α s Ap s −1β s = (g s ⋅ g s ) /(g s −1 ⋅ g s −1 )p s = g s − β s p s −1где(23.30)g s и p s - векторы градиента и сопряженного "направленияпоиска".

Метод сопряженных градиентов вырабатывает базис p s205Глава 18. Методы для задач упругопластичности( s = 1, 2,... ) в конечномерном арифметическом пространствевекторов y , поэтому теоретически число итераций, необходимыхдля отыскания решения, не превышает числа искомых компонентоввектора y .Для плохо обусловленных систем уравнений сходимостьпроцесса (23.29) - (23.30) может буть утеряна.

Для предотвращенияэтого применяется продобусловливание путем домножения системыуравнений на приближенную обратную матрицу B −1 . Обобщенныйметод сопряженных градиентов имеет следующий вид:для s = 0 полагается:g 0 = Ay 0 − bh 0 = B −1g 0p0 = h0и далее для s = 1, 2,... имеем:α s = (g s −1 ⋅ h s −1 ) /( Ap s −1 ⋅ p s −1 )y s = y s −1 − α s p s −1g s = g s −1 − α s Ap s −1h s = B −1g sβ s = (g s ⋅ h s ) /(g s −1 ⋅ h s −1 )p s = h s − β s p s −1(23.31)В качестве матрицы B можно использовать диагональную матрицу,составленную из диагональных элементов матрицы А. Такаяматрица В легко обращается один раз перед началом процесса(23.31). Операция B −1g s сводится к покомпонентному умножениюдвух векторов.Итерационный процесс (23.31) останавливается, есливыполнен следующий критерий:(g s ⋅ h s ) < ε *2 ∨ (p s −1 ⋅ p s −1 )α s2 < ε *2(23.32)где величина ε * является максимальным числом, добавлениекоторого к единице дает в результате елиницу для апифметики сограниченным числом разрядов для представления чисел.

Это числоназывается "машинным эпсилон" (Меткалф, 1985) и приближенноравно 10−6 для ЭВМ с четырехбайтовой арифметикой.206Глава 18. Методы для задач упругопластичностиВ случае ( Ap s −1 ⋅ p s −1 ) < ε *2 метод дает отказ из-за того, чтозадача вырождена, а при невыполнении условий (23.32) для s > 2 N( N- число неизвестных) задача является плохо обусловленной.. Вэтих случах метод (23.31) не позволяет определить решение либо всилу неединственности решения (точка ветвления), либо в силуочень плохой обусловленности задачи.

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

Она применима только в случаепостоянной зоны контакта.Модель идеального контакта на несогласованных сетках.В этой модели на контактной границе сетки контактирующих телнесогласованны и имеют несовпадающие узлы. Для каждого узла x +одного тела на контактной границе вводится фиктивный узел x*− вдругом теле, имеющий то же положение ( x + = x*− ) и требуетсявыполнение условий идеального контакта u + = u*− , где векторперемещений u*− в фиктивном узле выражается через узловыеперемещения второго тела интерполяцией. Обычно эти условияучитываются методом штрафа путем добавления в вариационноеуравнение виртуальных работ члена∫ λ (u+− u*− ) ⋅ (δ u + − δ u*− )dSSc207Глава 18.

Методы для задач упругопластичностигде Sc - контактная поверхность, λ >> 1 коэффициент штрафа..Модель идеального контакта на несогласованных сетках предложенаБаженовым (1984, 1994) и также реализована в работе Felippa(2001)..Эта модель нередко используется для проведения расчетовна областях решения сложной формы, состоящих и рядаподобластей, в которых вводятся несогласованные между собойсетки.

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

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

Характеристики

Список файлов книги