Учебное пособие по курсу лекций (1164067), страница 60
Текст из файла (страница 60)
Типичным примером такой нелинейности является контакт междудвумя телами или между телом и опорой. Размер пятна контакта изменяетсяв процессе деформации, а, следовательно, изменяются внешние нагрузки насистему.Глобальнаяматрица жесткостиопределяетсяпоэлементнымсуммированием матриц жесткости конечных элементов:⎡ k e ⎤ = ∫ [ B ]T [ D ][ B ] dV[ K ] = ∑ ⎡⎣k e ⎤⎦;⎣ ⎦ee[K ] =Vи зависит как от геометрии элементов (матрица градиентов [B] и объем V), так и отфизических свойств среды (матрица упругих коэффициентов [D]).Если упругие деформации малы, то геометрию элементов считаютнеизменной, а матрицу жесткости не зависящей от изменения координатузлов.Если упругие деформации оказываются большими, то изменениякоординат узлов оказывает влияние и на матрицу градиентов и объем392элементов. В этом случае следует учитывать влияние геометрическойнелинейности и на матрицу жесткости.Физическая нелинейность проявляется только во влиянии на матрицужесткости и является следствием нелинейных свойств среды – т.е.нелинейной зависимости напряжений от деформации.
В этом случаенелинейной становится матрица упругих коэффициентов:[ D ] = f D ({U })Типичными примерами нелинейно упругих сред являются резина иполиуретан.Таким образом, для нелинейных задач теории упругости процессорнаястадия процедуры метода конечных элементов усложняется по сравнению срассмотренной в параграфе 4.9.3.
Весь процесс решения разбивают нанесколько шагов (по перемещению какой-то опорной точки или по времени).На каждом шаге решение СнЛУ производят итерационными методами,использующими ее линеаризацию.Один из наиболее часто использующихся методов решения СнЛУ –метод Ньютона. Он предполагает последовательное итерационное уточнениерешения.Введем вектор-функцию:{Φ} ({U }) = [ K ]{U } − {R} = 0Разложим ее в ряд Тейлора и ограничимся линейными членамиразложения:∂ΦUΦU = ΦU 0 +⋅ ({U } − {U 0} ) ≅ 0∂ {U }U =U 0Здесь U0 – начальное приближение вектора перемещений на шаге.Итерационная процедура предполагает линеаризацию исходной системыуравнений на каждом шаге итерации:∂ {Φ}⋅ {∆U }m +1 = − {Φ}m∂ {U }mВ этом уравнении вектор приращений:{∆U }m +1 = {U }m+1 − {U }mМатрица частных производных носит названия якобиана:⎛∂ [ K ] ∂ { R} ⎞∂ {Φ}= ⎜[K ] −−[Y ]m =⎟∂ {U }∂U∂ {U } ⎠{}⎝mmВектор правых частей носит название вектора невязок:{∆Φ}m = −{Φ}m = ({R} − [ K ]{U })mОчевидно, что для линейных задач якобиан равен матрице жесткости, авектор невязок равен нулю.Окончательно на шаге итерации математическая модель представляетсобой СЛАУ относительно приращений узловых перемещений:393[Y ]m ⋅ {∆U }m+1 = {∆Φ}mВ этой системе уравнений якобиан и вектор невязок известны изпредыдущей итерации.
После определения вектора приращений на текущейитерации выполняют итерационное уточнение перемещений по формуле ипереходят к следующей итерации:{U }m+1 = {U }m + {∆U }m+1Итерационный процесс завершают, если норма вектора невязок инорма вектора приращений оказываются меньше некоторого напередзаданного малого числа.{∆U }m ≤ εU ;{∆Φ}m ≤ ε ΦИногда ограничиваются выполнением одного из этих условий.Процесс итерационного решения можно пояснить следующимрисунком, выполненным применительно к одномерному случаю:ΦΦ(U)=0Φ0tan α =Φ1→ΦdΦ=dU −∆UdΦ∆U = −ΦdUΦm∆Um ∆U1Um U1U0UЕсли сходимость итераций не достигается за заданное количествоитераций Mmax, то шаг по перемещению или времени уменьшают иитерационное решение производят заново.
Отсутствие сходимостифизически обычно означает, что на данном временном илипространственном интервале математическая модель существеннонелинейна.После завершения очередного шага вычисляют величины напряженийи деформаций на шаге. (Для работы алгоритма вычисление напряжений идеформаций на шаге не требуется. Его выполняют если необходимопроанализировать историю нагружения.){ε }n = [ B ]n {u e }, {σ }n = [ D ]n {ε }nрассчитываю новуюnПосле этогокоординаты узлов.геометриюдетали-новые394В обобщенном виде процедура процессорной стадии метода конечныхэлементов при решении нелинейных задач теории упругости выглядитследующим образом:Шаг решения n (цикл ПОКА n ≤ N max )Итерация m: определение вектора перемещений {U } на шаге итерационнымn[ ]решением СНлУ K {U } = { R}n(цикл ПОКА ( m ≤ M max ) ∧( ∆U m≤ εU ) ∧ ( ∆Φ m ≤ ε Φ ) )Определение матрицы жесткости (цикл по конечным элементам)eОпределение матрицы жесткости конечного элемента ⎡ k ⎤⎣⎦meАнсамблирование [ K ] = ∑ ⎡ k ⎤m⎣ ⎦emОпределение якобиана YmОпределение вектора невязок {∆Φ}m= ∆Φ}m с целью учетаМодификация СЛАУ Y {∆U }m +1 {mкинематических граничных условийРешение СЛАУ и определение вектора приращений {∆U }m +1Итерационное уточнение вектора перемещений[ ][ ]{U }m+1 = {U }m + {∆U }m+1Вычисление напряжений и деформаций (по запросу)= X + UОпределение новой геометрии детали { X }n +1 { }0 { }n4.9.12Особенностижестко-пластических телпримененияМКЭдляБольшое развитие получила жестко – вязко пластическая постановкаметода конечных элементов, особенно для процессов горячей объемнойштамповки.
В такой постановке пренебрегают упругими деформациями, ноучитывают как деформационное, так и скоростное упрочнение. Инымисловами напряжение текучести σ s считается функцией накопленныхпластических деформаций и скорости деформаций.Физические уравнения для жестко – вязко пластической постановкизадаются в виде соотношений Сен-Венан – Леви – Мизеса.2 σиsij =ε ij3 εигде: sij - девиатор напряжений, ε ij - тензор скоростей деформаций. Чтобы избежатьпутаницы в индексах сокращенной записи тензорных величин несколько изменим395привычнуюиндексациюинтенсивностинапряжений3σи = σi =sij sij2и2ε ijε ij - интенсивность скоростей деформаций.3интенсивности деформаций ε и = ε i =Как было показано в параграфе 4.9.7 в матричной форме взаимосвязьмежду узловыми переменными и напряжениями конечного элемента имеютвид:{s} = ⎡⎣ D p ⎤⎦ [ B ] v(e)Узловыми переменными в этом случае служат не узловые{ }перемещения, а узловые скорости{ v(e)} .
Матрица связи [D ] постоянна иpзависит от вида напряженного состояния.Особенностью физических уравнений для жестко-пластического телаявляется то, что с узловыми переменными (узловыми скоростями)непосредственно связаны компоненты девиатора напряжений, а не тензоранапряжений, как для упругого тела.Вариационную формулировку МКЭ для жестко-пластических телстроят с помощью функционала Маркова, который принимает минимальноезначение на истинном поле скоростейΦ = ∫ σ иε и dV − ∫ pi vi dFVFpПервый интеграл функционала отражает мощность пластическойдеформации, а второй – мощность внешних сил ( pi , vi - проекции удельныхвнешних сил и скоростей на координатные оси).
Можно показать, что вматричной записи интенсивность напряжений принимает вид:22 Tεи =ε ijε ij = {ε } ⎡⎣ D p ⎤⎦ {ε }33С учетом взаимосвязи между скоростями деформации и узловымискоростями {ε } = [B ] v e окончательно получим:{}{ }{ }2 e Tv[ B]T ⎡⎣ D p ⎤⎦ [ B ] ve3Для тела, представленного в виде совокупности конечных элементовиспользуют дискретную формулировку функционала Маркова:εи =Φ = ∑ ∫ σиeVЗдесь{ }{ }{ }{ }2 e TT⎡ee T⎤vBDBvdVvPe−[ ] ⎣ p ⎦[ ]∑3eиспользовано,чтовМКЭсосредоточенным узловым нагрузкамвнешние{Pe} ,силысводяткчто позволяет заменитьинтегрирование суммированием по узлам.Функционал Маркова принимает минимальное значение при:396∂Φ=0∂ veВыполним дифференцирование85:{}∂Φ={ }∂ v=e∂∑ ∫ σи{ }∂ v=∑∫eVeeVσи ⋅{ }{ }2 e T∂v[ B ]T ⎡⎣ D p ⎤⎦ [ B ] ve dV − e3∂ v{ }Pe}{}{e∑ veT{ }112T⋅⋅ ⋅ 2 ⋅ [ B ] ⎡⎣ D p ⎤⎦ [ B ] v e322 e TT⎡e⎤vBDBv[ ] ⎣ p ⎦[ ]3{ }{ }=dV −εи{ }−∑ Pe = 0eПриравняв результат нулю и сгруппировав, получим определяющуюсистему уравнений метода конечных элементов для анализа деформацийжестко-пластических тел:[ K ]{U } = {R}[K ]Здесь:– глобальная матрица жесткости, полученная почленным суммированием матрицжесткостейконечныхэлементов86.конечного элемента с учетомопределяется интегралом:⎡k e ⎤ =⎣ ⎦[ K ] = ∑ ⎡⎣ k e ⎤⎦ .условияМатрицаeпластичностиМизесажесткостиσи = σ s2 σs[ B ]T [ D ][ B ] dVe 3 εи∫V{R} – вектор внешних узловых нагрузок {R} = ∑{Pe} .{U } – вектор неизвестных узловых скоростейeОпределяющая система уравнений нелинейна, поскольку зависит отинтенсивности скоростей деформаций ε и , а та, в свою очередь, от85Правила дифференцирования матричных произведений см.
4.9.8По своей размерности это скорее матрица вязкости, однако такое названиене принято86397неизвестных узловых скоростей. От скоростей деформаций также можетзависеть и напряжение текучести σ s , входящее в интеграл для определенияматрицы жесткости. Кроме того, система уравнений также зависит и отистории нагружения, поскольку напряжение текучести зависит отнакопленных пластических деформаций, мерой которых принято считатькритерий Одквиста q и который в данном случае можно определитьинтегрированием интенсивности скоростей деформации по времени .:σ s = σ s ( q, ε и ) , q = ∫ ε и dttФормулировка МКЭ с помощью функционала Маркова имеетнедостатки:Не учтено условие постоянства объема (условие несжимаемости), поэтомув общем случае решение может этому условию не удовлетворять.Невозможно сразу определить тензор напряжений, поскольку с узловымискоростями связаны компоненты девиатора, а среднее напряжениенеобходимо определять дополнительно.Условие несжимаемости для жестко пластического материаларавносильно равенству нулю первого инварианта тензора скоростейдеформаций:I1 (Tε ) = ε ii = ε x + ε y + ε z = 0Один из наиболее распространенных способов учета несжимаемостиматериала – модификация функционала Маркова с помощью штрафнойфункции.1Φ* = Φ + ∫ Aε v2dV2VЗдесь A - большое положительное число, εV = ε x + ε y + ε z - скоростьобъемной деформации.