книги / Решение нелинейных задач теории оболочек на ЭВМ
..pdfМетод Ньютона. Как известно (см. гл. 3 или 8), метод Ньютона приводит к следующей итерационной схеме:
{[К0]+ |
[К*([<**])]) [Ad*] = [Qu (№*1)1, |
|
[dk+l] = |
\dk] -f [àdk], k = 0, 1, 2, . . . |
(9.148) |
Для существенно нелинейных задач метод Ньютона состоит в |
||
последовательном приложении приращений нагрузки и |
нахожде |
нии приращений перемещений на всех промежуточных шагах нагру жения. Процесс продолжается до тех пор, пока величина [Дdk\ не станет достаточно малой или пока неуравновешенность сил [Qu] не станет достаточно малой.
Модифицированный метод Ньютона. Применение метода Ньюто на приводит при вычислении матрицы [Ко + К* ([41)1 из (9.148) к значительным затратам машинного времени* Расчет матрицы [/(*] длителен сам по себе, кроме того, каждый раз матрица [К*1 кор ректируется, и новая матрица коэффициентов подвергается обра щению. Поэтому часто применяют модифицированный метод Нью тона, в котором матрица [К*] остается неизменной для нескольких итераций и приращений нагрузки, а после того, как скорость схо димости начнет уменьшаться, матрица корректируется на основе перемещений, которые система имеет в данный момент.
Кроме того, для более точной оценки начальных перемещений при переходе от одного приращения к следующему метод часто до полняется различными способами экстраполяции. Для ускорения сходимости перемещения часто пересчитываются посредством ре лаксации. В вычислительную программу могут быть введены разли чные схемы автоматического уменьшения или увеличения величины шага по нагрузке.
Метод возмущений. Уравнение рановесия (9.143) запишем в виде
[/Col la i = P ÎQI - [Q* ([<*!)], |
(9.149) |
где перемещения [d] считаются функциями параметра нагрузки р.
В качестве такого параметра нагрузки р может быть выбран коэф фициент приведения или безразмерный параметр нагрузки. При известной точке равновесия на кривой нагрузка—перемещение про изводные по этой кривой можно использовать для расчета следу ющей точки. Разложим_перемещение [d] в ряд Тейлора около извест
ного положения (fdj,p):
dj (р -j- Ар) dj (р) + d/ (р) Ар + - j dj (р) (Ар)2 +
+ |
*£• dj (р) (Ар)3 + . . |
(9.150) |
|||
где |
|
|
|
|
|
ddj \р) |
* /« |
d2dj(p) |
/ = |
1, 2, |
|
dp ' |
■dP2 |
||||
’ |
|
функции можно записать в виде |
|
|
|
||
N |
|
ÔQ |
т |
= [/с*] ш . |
|
V? (М)1 = 23 |
|
dd |
(9.153) |
||
/=1 |
1J |
|
|
||
Тогда уравнение (9.152) с учетом (9.153) перепишется в виде |
|||||
([Ко] + |
[К*]) [d] = |
[<7j- |
(9.154) |
Продифференцировав уравнение (9.152) еще раз по параметру ру получим
[Ко] [d] = — [ ф ([d])]. |
(9.155) |
В уравнения (9.152), (9.155) не входит матрица [/С*]. Поскольку на вычисление матрицы [/<*] затрачивается достаточно много ма шинного времени, то можно предположить, что уравнения (9.152) и (9.155) более предпочительны. Однако в них входят производные от нелинейной вектор1-функции [Q*], которые могут привести к за труднениям, связанным с устойчивостью вычислительной схемы, так как производные должны быть выражены через конечные раз ности.
Так, в некоторых задачах пластичности, решаемых МКЭ, ис пользовались левая разностная формула для [ф] и правая для [d]:
И = ^ fld i+1J — Id1]). |
<9.15в> |
Подстановка этих формул в уравнение (9.152) приводит к рекуррент
ному соотношению для |
[d * 1] на шаге |
t -f-1 : |
|
|
I d * 1] = |
[d1] + |
[ К о ] - 1(Д~р ÎQ] - |
[Q*£] + [Q*(l‘-1 )]). |
(9.157) |
Оказывается, что |
это |
конечно-разностное уравнение становится |
неустойчивым для существенно нелинейных задач.
Метод самокорректирующихся начальных значений. Решение, получаемое с помощью уравнений (9.152), (9.154), (9.155), может отклоняться от точного. Поэтому разработано несколько вариантов модификации метода начальных значений [71, 79]. Чтобы скорректи ровать отклонение полученного решения по отношению к точному, необходимо в уравнении (9.152) или (9.154) учесть неуравновешен
ность усилий |
(9.146). |
|
|
Z и добавляя к |
|
Умножая |
(9.146) на скалярный коэффициент |
||||
правой части уравнений (9.152), (9.154), придем |
к уравнениям |
||||
[/(о] ([d] + Z [d]) = |
(1 + ~pZ) [QJ - |
Z («5*1 + |
{ IQ*I); |
||
(I/fol + (K*l) [d] + |
2 № 1 Id] = (1 + |
pZ) IQI - |
Z IQ*I. (!).i58> |
Член, добавляемый к каждому из уравнений, равен нулю для истинного равновесного положения, т. е. в случае, Когда решение
оказывается на истинной кривой нагрузка—перемещение. При Z =
1 *= ^второе уравнение (9.158) эквивалентно уравнению самокоррек
тирующего метода приращений (9.145). Очевидно, что поскольку уравнения (9.152) и (9.154) можно записать в виде
[Qt/] + Z[Qu] = [0], |
(9.159) |
неуравновешенность [Qo] будет экспоненциально уменьшаться, при
ближаясь к нулю с возрастанием параметра р (именно .поэтому метод назван самокорректирующимся). Уравнение (9.155) приводит ся к самокорректирующемуся виду путем добавления в правую
часть неуравновешенности Z[Q(/] и ее первой производной C[Qy] по параметру р
[Ко] [à] - - [Q*] + Z [Qt/] -J-С [Qcl. |
(9.160) |
Скалярные коэффициенты Z u C иллюстрируют для уравнения (9.160) влияние дополнительных усилий. С помощью (9.J46) последнее урав нение приводится к виду
[Ко] аЙ + C[d] + Z [d]) = (С + pZ) [Q] -
- (Z LQ*]-+ C [Q*] + [Q*]). |
(9.161) |
Это уравнение допускает простую физическую интерпретацию, ко торая оказывается полезной при его решении. Левая часть его
эквивалентна уравнениям простого гармонического колебания с коэффициентом затухания С и недемпфируемой собственной частотой’
J/Z . Правая часть трактуется как вынуждающая сила [79, 71] (без указания размерности). От выбора коэффициентов С и Z за висит быстрота сходимости итерационного процесса.
Заметим, что при пользовании уравнениями (9.157) или (9.160) требуется только один раз во всем решении обращать матрицу жесткости [Ко]. Воспользовавшись соотношением (9Л53), можно получить разрешающие уравнения в ином виде, однако при этом будет введена матрица [К*] и вычислительные преимущества урав нений будут утрачены. Для всех значений нагрузок не было, найдено устойчивых методов интегрирования уравнения (9.158).
Независимо от вида [к] метод приводит к расходящемуся реше нию, начиная с некоторых значений нагрузок.
Уравнение (9.161) решалось с помощью левой четырехточечной разностной формулы для [к] и левой трехточечной разностной фор
мулы для [к] (интерполяции назад). Член [Q*] отбрасывался как малый, а для оставшихся нелинейных членов с псевдосилой было получено приближение с помощью линейной экстраполяции по предыдущим приращениям нагрузки.
Рекуррентное соотношение для [d1] на /-ом шаге приращения
нагрузки |
(соответствующем параметру Нагрузки |
р1)\ |
|
|
|
о т = ( V [« о !-1( с + m m - |
|
|
|
- |
+ i(Q * « -i> _ |
Q»((-S)] (tTp + |- ) z ) + (5 + |
2ДрС) Id1-'] - |
|
- |
(4 + - ^ ) [d ‘- |
2! + Id'-3j) (2 + I ДрС + |
ZAp^. |
(9.162) |
После нескольких численных экспериментов было установлено что параметры Z и С по соображениям устойчивости должны быть приняты в виде функций, уменьшающихся при возрастании на грузки.
В результате решения большого числа задач для оболочек вращения обнаружено, что обычный метод .приращений часто при водит к существенным отклонениям от истинного решения. До биться точности можно разбиением пути нагружения на мелкие шаги. Но при этом для реализации метода требуется больше вре мени, чем для реализации по методу Ньютона. Если необходима лишь приближенная оценка решения, то можно использовать ме тод приращений. Добавление поправочного члена в метод прира щений повышает точность решения. Самокорректирующийся метод приращений при той же простоте, что и для обычного метода приращений, достигает точности метода Ньютона.
Методы с применением итераций или последовательных подста новок для решения уравнений равновесия могут использоваться лишь в ограниченном числе случаев. Итерационные методы при менимы лишь в случаях, когда нелинейное' решение для переме щений отличается от линейного не более чем в полтора-два раза. При большей нелинейности метод просто не сходится. Метод достаточно прост и требует минимальных затрат машинного време ни. Его можно применять для получения начальных значений перемещений в других методах. Метод Ньютона наиболее точен, но требует наибольших затрат машинного врёмени. Модифициро ванные варианты метода Ньютона снижают расход машинного времени, но требуют от программиста достаточного опыта.
Метод возмущения требует больших затрат машинного времени. Различные варианты метода начальных значений заслуживают дальнейшего изучения. Самокорректирующиеся варианты этого метода Aarof такую же точность, как и метод Ньютона, и требуют
значительно |
меньших затрат машинного времени, чем для |
послед |
него. В них матрица жесткости обращается всего один |
раз за |
|
весь расчет. |
|
|
Для задач со слабой нелинейностью можно рекомендовать применять самокорректирующиеся методы приращений или простой итерации как наиболее экономичные. В случае существенно не линейных задач можно применять метод Ньютона или самокор ректирующийся метод начальных значений.
Процедура метода итераций начинается с определения линейного решения при максимальной нагрузке. Это решение используется далее в качестве пер вого приближения при вычислении нелинейных членов. Итерационная про цедура продолжается до момента наступления сходимости или до заданного числа итераций, если процесс расходится. Метод хорошо зарекомендовал себя для слабо нелинейных систем. Для сильно нелинейных систем решение может приобрести колебательный характер, хотя и будет медленно сходиться.
Комбинированный метод приращений и итераций включает следующий •момеьт. После первого приращения нагрузки определяют перемещения. Путем применения итераций на этом же этапе нагружения уточняют решения. Ите рации продолжаются до наступления сходимости. Затем перемещения экстра полируются на следующие нагружения и итерации вновь продолжаются до наступления сходимости. Продолжаем шаговые нагружения и итерационную процедуру до тех пор, пока нагрузка не Достигнет своего максимального или еаданного значения. Этод метод позволяет исследовать сильно нелинейные системы и может быть использован для оценки критической нагрузки. Пара метр экстраполяции может быть принят равным единице, пока критерий схо димости выполняется точно.
Расчет рассматриваемой задачи методом приращения нагрузки и комбини рованным методом привел к одним и тем же результатам. На рис. 9.12 и 9.13
дано |
сравнение |
меридиональных напряжений, |
полученных |
этими методами |
|||||
в линейной и нелинейной |
постановках для двух значений |
радиуса |
тора R = |
||||||
= 3,2 |
см и R = |
6,35 см |
соответственно. |
Решение |
задачи |
в |
линейной поста |
||
новке |
получено также и тремя другими |
методами: |
наблюдается полное сов |
||||||
падение результатов. На |
рисунках по оси ординат отложена величина мериди. |
||||||||
онального напряжения с*, деленная на 700 кГ/см2, так что Oi = о} |
700 кГ/см2, |
||||||||
вдоль |
меридиана |
s, деленное на 2,54, т, |
е. s = |
s'.24,5 см. |
Кривой |
1 обозна |
чено напряжение на внешней поверхности оболочки, кривой 2 — на внутренней.
Пунктирной линией показано решение линейной задачи. |
|
|
120 |
мм, |
||||||||||||
|
Пример |
2. |
Пусть |
пологая сферическая оболочка радиуса R = |
||||||||||||
толщиной 0,4 мм со стрелой подъема 2,2 мм защемлена |
по |
контуру |
и нахо |
|||||||||||||
дится |
под |
действием |
сосредоточенной |
|
|
|
|
|
||||||||
силы Р, приложенной в полюсе. Модуль |
|
|
|
|
|
|||||||||||
упругости |
Е = |
6,9 • |
104 М Н,м2, коэффи |
|
|
|
|
|
||||||||
циент Пуассона |
|
v = |
0,3. |
одной гармо |
|
|
|
|
|
|||||||
|
Методом |
п. |
4 с учетом |
|
|
|
|
|
||||||||
ники задача приводится к нелинейной |
|
|
|
|
|
|||||||||||
алгебраической |
системе |
уравнений |
вида |
|
|
|
|
|
||||||||
{9.109), |
которая |
решается |
несколькими |
|
|
|
|
|
||||||||
методами. [79]. |
Решение, полученное ме |
|
|
|
|
|
||||||||||
тодом |
Ньютона, |
|
хорошо |
согласуется с |
|
|
|
|
|
|||||||
конечно-разностным решением при при |
|
|
|
|
|
|||||||||||
ращениях нагрузки Др — 4,45 Н |
и |
Др — |
|
|
|
|
|
|||||||||
= 22,25 Н и принимается за истинное. |
|
|
|
|
|
|||||||||||
Вычисления |
показывают, |
что |
поведение |
|
|
|
|
|
||||||||
оболочки существенно нелинейно, причем |
|
|
|
|
|
|||||||||||
имеется начальный участок мягкой не |
|
|
|
|
|
|||||||||||
линейности и последующий участок жест |
|
|
|
|
|
|||||||||||
кой |
нелинейности |
(рис. |
9.14, кривая 1). |
|
|
|
|
|
||||||||
Пунктиром |
показано |
решение по линей |
|
|
|
|
|
|||||||||
ной |
теорий. |
|
|
|
|
|
|
|
|
результаты вычислений нагрузка — прогиб |
||||||
|
На рис. 9.14 приведены также |
|||||||||||||||
методом Ньютона и модифицированным методом приращений при Др = |
4,45 Н |
|||||||||||||||
(кривая 1), обычном методом приращений при Др = |
4,45 |
Н |
(кривая 2) и |
при |
||||||||||||
Др = |
0,56 |
Н |
(кривая |
3) |
и |
модифицированным |
методом |
приращений |
при |
|||||||
Ар = |
22,25 |
Н |
(кривая 4). |
|
|
|
|
малого приращения |
||||||||
|
Из |
рисунка |
|
видно, |
что если не принять достаточно |
|||||||||||
нагрузки, то |
решение будет отклоняться от истинного при использовании метода |
|||||||||||||||
приращений. |
Модифицированный |
самокорректирующийся |
вариант |
метода |
приращений (9.145) обеспечивает значительно лучшую сходимость решения по сравнению с обычным вариантом. При приращении нагрузки -Др — 4,45 Н решение хорошо согласуется с точным. При Др = 22,25 Н решение по моди фицированному методу прирашений несколько отклоняется от точного реше ния. Отметим, что прямые итерации при решении уравнений равновесия не дают сходимости по нагрузке, превышающей 90 Н, а с применением релакса ционной схемы процесс сходится только до 134' Н.
Интегрирование уравнения (9.152) в методе начальных значений выпол нено способом Адамса второго порядка с предсказанием и исправлением. Хотя
|
|
|
Рис. 9.16 |
при Др = |
1,11 Н |
этот |
способ дает очень точные результаты, все же при |
нагрузке |
выше |
55,5 Н |
решение перестанет быть сходящимся. По-видимому, |
здесь необходимо воспользоваться методом Адамса второго порядка для инте грирования уравнения (9.154), в которое входит матрица [К*1- Однако при
этом требуется значительно больше |
машинного времени. |
|
|
|
|||
Самокорректирующийся метод начальных значений (9.1581 позволил полу |
|||||||
чить интересные результаты (рис. 9.15). Рассматривалось два случая,' |
|
||||||
I) Z = Др ^р |
Apy/i *G = |
• А/7 >Z, Др — 0,89 Н; |
|
|
|||
2) 2 ~ Др [P Ар]1/*’ С = |
~2 Z0,2> |
Лр = 0,22 Н' |
|
|
|
||
Результаты вычислений весьма точно показывают поведение системы с |
суще |
||||||
ственной нелинейностью. В обоих случаях решение колеблется |
около |
истин |
|||||
ного (пунктир — первый |
случай, |
штрих-пунктир — второй). |
Эти |
колебания |
|||
затухают после нескольких циклов. |
напряжения |
в наружном |
слое |
в |
центре |
||
На рис. 9.16 дана зависимость |
оболочки от нагрузки. Напряжения, вычисленные самокорректирующимся мето дом начальных значений, также колеблются около истинного, но эти Колебания менее заметны. В ряде случаев фаза колебаний напряжений не совпадает
строчки (9.164), относятся как 8 : 3 : 6. Учитывая точность третьего случая,’
отмечается, что метод начальных значений заслуживает рассмотрения в качестве
практического метода |
анализа. |
самокорректирую |
|||
Пунктирной кривой на рис. 9.17 нанесено решение по |
|||||
щемуся методу начальных значений при Ар = 0,445 Н, сплошной — по |
методу |
||||
Ньютона |
при Др = 4,45 Н, прямая линия — решение по |
линейной |
теории. |
||
Результаты счета с помощью самокорректирующегося метода начальных |
зна |
||||
чений при Ар = 0,89 |
Н также дают колеблющиеся решения относительно |
||||
истинных |
равновесных |
состояний, но с несколько большей |
амплитудой, |
чем |
|
при Ар = |
0,445 Н (соответствующая кривая на десунке не приведена). |
|
|
Пример 4. Рассмотрим конструкцию, состоящую из цилиндрической, торо идальной оболочек и оболочки отрицательной гауссовой кривизны [70].
Пусть оболочка находится под действием внутреннего давления, края
оболочки жестко защемлены. |
геометрические и механические харак |
||||||||||
При расчетах приняты следующие |
|||||||||||
теристики: |
|
|
|
|
|
|
|
|
|
|
|
|
# ц = |
252 мм; RT —65,2 мм; |
/?0 = |
138 мм; 7/ц = |
62 |
мм; |
|
|
|||
|
h = |
3,2 мм: Е = 7,3 Ю4^ г ; |
v = |
0,333. |
|
|
|
|
|||
Задача приводится к алгебраической системе уравнений при 50 элементах |
|||||||||||
и одной гармонике в ряде Фурье. |
|
|
|
|
|
с |
приращениями |
||||
Для получения решения применялся метод Ньютона |
|||||||||||
нагрузки |
Ар = 22,25 Н. Проводилось сравнение с решением, |
полученным с по |
|||||||||
мощью самокорректирующегося метода начальных значений |
(9.158) |
при |
Ар = |
||||||||
= 0,1, Z |
Ар (Р |
Др)1/2 = 1 , С —0,05]Л£. |
fîpu |
расчете по методу |
начальных |
||||||
значений полная нагрузка была сведена |
к такому |
масштабу, |
что параметр |
||||||||
нагрузки Р изменялся от 0 до 100. Это |
позволило |
придать |
некоторое |
едино |
|||||||
образие решению, так что если найдены |
|
Âp, |
Z, |
С, |
приводящие к |
реальным |
результатам, этими значениями можно пользоваться для большинства задач независимо от характера их нелинейности. Таким образом, значение Д/? = 0,1
- |
о с |
соответствует величине приращения нагрузки |
при полной нагрузке |
3,5 |
Ар = 0 ,3 4 ^
Закрнтическое |
поведение конструкции |
легко рассчитывается по |
методу |
Ньютона при Др = |
Н |
через неустойчивую область |
метод |
34,5-р . При переходе |
„кН
сходится к значению 2,00 -~а , а затем совершает скачок в закрнтическое состо-
кИ
яние при 2,05 |
(сплошная |
кривая на рис. |
9.18). |
|
|
|||
Метод начальных значений при Ар = 3,46 |
Н |
Дает очень |
точное |
значение |
||||
^ |
||||||||
в докритическом |
состоянии, |
а затем |
переходит |
к закритической форме при |
||||
мерно при Р = 2,07— Хотя |
решение несколько |
колеблется, |
ио это |
колебание |
||||
происходит около |
истинной |
кривой, |
н колебания затухают (пунктирная кривая |
|||||
на рнс. 9.18). |
|
быстрого |
затухания |
колебаний решений, |
если при |
|||
Можно добиться более |
||||||||
|
кН |
|
|
|
|
|
|
|
нагрузке около 222,5 начать новый счет с меньшим приращением нагрузки. С уменьшением колебаний приращение нагрузки можно вновь увеличить. На