CC BY
37
УДК 539.3
МЕТОДИКА ЧИСЛЕННОГО РЕШЕНИЯ НЕОСЕСИММЕТРИЧНЫХ ЗАДАЧ НЕЛИНЕЙНОГО ИЗГИБА ПОЛОГИХ ОБОЛОЧЕК ВРАЩЕНИЯ
© А.Н.Шихранов
В статье предлагается методика расчета неосесимметричного напряженно-деформированного состояния (НДС) тонких гибких пологих оболочек вращения и круглых пластинок переменной толщины, находящихся под действием силовых и температурных нагрузок. Методика основана на синтезе методов конечных разностей, матричной прогонки, общей итерации и продолжения решения по параметру. Приводится сравнение результатов ее применения с известными в литературе данными, полученными другими подходами. Показано, что принятый в методике итерационный процесс расчета НДС обладает быстрой сходимостью.
Ключевые слова: оболочки вращения, неосесимметричное деформирование, температурное нагружение, напряженно-деформированное состояние, критические нагрузки, волнообразование.
1. Введение
Тонкостенные оболочки вращения широко используются как конструктивные элементы в машиностроении, авиастроении, в строительной индустрии и во многих других областях техники. Для улучшения технических характеристик и качества проектирования оболочек, оценки их жесткости, прочности и устойчивости с учетом реальных условий работы конструкций необходимо знать их истинное напряженно-деформированное состояние (НДС). Во время эксплуатации оболочек под действием различного рода нагрузок (силовых, температурных и т.д.) возникают большие перемещения, неосесимметричность деформирования. Это приводит к необходимости исследования НДС оболочек с учетом нелинейности и двумерности в математической постановке задачи. Получить решение при такой постановке зачастую возможно лишь численным путем.
В статье изложена численная методика определения нелинейного НДС пологих оболочек вращения переменной толщины и их различных фрагментов (секторов, сегментов, панелей), находящихся под действием несимметричных нагрузок и умеренного температурного поля. Методика основана на синтезе методов конечных разностей, матричной прогонки, общей итерации и продолжения решения по параметру. Показано, что принятый в методике итерационный процесс расчета НДС обладает быстрой сходимостью.
2. Постановка задачи и метод ее решения
Рассматривается НДС тонких пологих оболочек вращения переменной толщины, находящихся под действием несимметричных нагрузок и умеренного температурного поля. Считается, что материал оболочки однородный, термически изотропный, линейно упругий, его физико-
механические характеристики не зависят от температуры. Принимается, что закон распределения температуры по толщине оболочки линейный, начальные напряжения отсутствуют.
Для описания деформирования оболочки используются уравнения геометрически нелинейной теории пологих оболочек среднего изгиба по классификации Х.М.Муштари [1], записанные в перемещениях. При этом температурное поле учитывается введением дополнительных членов в физические соотношения на основе гипотез Дюамеля-Неймана.
Нелинейное поведение конструкции переменной толщины описывается следующей системой дифференциальных уравнений и граничных условий
Ьи=Р в О, (1)
Ш=ц на Г (2)
Здесь и = (ыу,м>)т - вектор перемещений, Р, я - векторы поверхностных (силовых и температурных) и контурных нагрузок, О - область, занимаемая конструкцией, Г - граница этой области, Ь, / - нелинейные дифференциальные операторы. Для решения двумерной краевой задачи (1), (2) используется схема последовательных приближений в форме метода общих итераций М.С. Корнишина [2] с введением коэффициентов релаксации
и(к+1) = и® + т (и* - и®) (3)
Здесь к - номер итерации, т - диагональная матрица коэффициентов релаксации, ненулевые компоненты которой ти,ту,тм, обеспечивают сходимость итерационного процесса (3).
тт* / * * *\ Т
Вектор и =(и V ^ ) является решением следующей линейной краевой задачи
Ьои* =Р - Ь1и(к) (4)
1ои*|г =Я - /Л, (5)
где L = Lo + Ц1, l = 1о + ^, ^, l1 - нелинейные и часть линейных операторов с переменными по угловой координате ф коэффициентами, Lo, 1о _ главные линейные операторы дифференциальных уравнений и граничных условий, описывающие НДС оболочки постоянной толщины, коэффициенты которых не зависят от окружной координаты. Оператор 1о выделен для однотипных по всему контуру Г граничных условий.
При записи соотношений (4), (5) использовалось представление толщины h, жесткостей конструкции на растяжение-сжатие K и изгиб D в виде сумм:
h = Ы + М, K = Ко + Кь D = Do + Dl, (6)
в которых Ы0, К0, D0 являются постоянными, а Ы1, К1, D1 зависят от координат г, ф. Представление жесткостей в виде (6) вводится для расчета оболочек переменной толщины. Величины К0 и D0 вычисляются для оболочки некоторой характерной постоянной толщины ы0.
Матрица дифференциальных операторов L0 имеет вид:
А( W,
г+2
-&*> ^+1 + с^> wi
3 < 0 0
0 А V 0 (7)
0 0 AW.
+D(w>Wi-1 + Е (w)Wг■_2 = К^.
Здесь А(и), Б(и),..., Е(^ - квадратные (раз-
(и)
мерностью пг х пг) матрицы; и,, У,, Wг■, К ’,
1?(Г) !?(w) -
К- 7, К- 7 - векторы искомых перемещений и
Здесь
Ли = А _ а^О/г2 Э2/Эф2, Ду =
=Д _ [1+ (1 - V-1 )]/г" Э2/Эф2,
AW = Л2 + аь Vl = (1 - v)/2, а1 =
=( к1 +2v к к2 + к2) Ко / Do,
Л - оператор Лапласа в полярных координатах, V
- коэффициент Пуассона, к1, к2 _ кривизны оболочки.
Ввиду громоздкости явный вид матрицы дифференциальных операторов L1 и операторов 1о и 11 здесь не приводится.
В случае замкнутых в окружном направлении оболочек система (4), (5) дополняется условиями периодичности искомых функций в окружном направлении
и(г, ф) = и(г, ф + 2п) (8)
Краевая задача (4), (5), (8) решается методом конечных разностей на полярной сетке. Используются центральные разностные схемы второго порядка точности. Система разностных уравнений записывается в виде
А(и )и,+1 + С(и )иг + Б(и )и,_1 = Ки),
г (9)
А(" )У,+| + )У + Б(” )У,_1 - цМ ,
правых частей размерности п^; г = 1, п^ - лучи полярной сетки пг х Пф .
При представлении уравнений равновесия и граничных условий в виде (4), (5) коэффициенты-матрицы системы (9) становятся независящими от угловой координаты ф. В результате этого требуемый объем хранимой на ПЭВМ информации сокращается в пф раз, что приводит к
большой экономии оперативной памяти.
Система разностных уравнений решается методом циклической матричной прогонки в случае оболочек, замкнутых в окружном направлении, и методом матричной прогонки в случае других фрагментов оболочек вращения. При этом прогонка организуется по угловой координате отдельно по каждому из перемещений. Такой выбор направления прогонки обусловлен независимостью матриц в системе (9) от окружной координаты. Данное обстоятельство дает возможность более эффективно реализовать процесс решения двумерной нелинейной краевой задачи.
Известно, что сходимость процесса последовательных приближений при решении систем нелинейных уравнений практически любыми методами существенно зависит от того, насколько корни нулевого приближения близки к истинным. В данной работе для решения нелинейных разностных уравнений используется метод общей итерации в форме (3) в сочетании с экстраполяцией для получения корней нулевого приближения искомых величин. При этом применяются экстраполяционные формулы повышенной точности [2], которые в случае равномерного шага изменения параметра системы имеют вид: w0p+| = 3w^,г - w^>l - w0p, 7=1,2,.. N (10)
],г+1 ^^;,г ^;,м " ],г
w0],1+| = (2w^,г - w^>l) Wj,i/ w0J]1 ,
(11)
7,г "'],г-И "И
7=1,2,.N о о
Здесь, соответственно, w 7,г, w у,г+1- нулевые приближения для X, и Хг+1, а w^■,г■_|, w^■,г■ — точные значения корней при X, и X г-1, где Хг - значения параметра системы, в качестве которого принимается перемещение или нагрузка.
Значения коэффициентов релаксации ти,тМ,xw в (3) выбираются экспериментально. Решение многочисленных практических задач показывает, что сходимость итерационного про-
цесса (3) может быть достигнута за счет надлежащего выбора единого для ти, ту, т у значения т в интервале 0 < т < 1.
Критерием достижения необходимой точности итерационного процесса в вышеописанном алгоритме принимается условие
|и(и) - и(и-1) | / |и(и-1) | < е, (12)
где е - определяющая точность малая величина, и(и)= и(п)(и(п^^(п),у(п)), п- номер итерации.
Для построения кривой деформирования "нагрузка-перемещение" с прохождением предельных точек на ней используется метод продолжения по параметру. В качестве ведущего параметра выбирается либо параметр нагрузки ХР, либо параметр перемещения Ху (соответствующий прогибу в каком-нибудь узле конечноразностной сетки, суммарному прогибу и т.д.), определяемый уравнением
8^ = Ху, (13)
где 8 - линейный оператор.
Если в качестве ведущего параметра выбирается Ху, то искомый параметр нагрузки ХР находится совместным решением третьего уравнения системы (9), записанного в более общем виде Aw = Хр Р ‘, (14)
и уравнения (13):
Хр =(Ху- 8w*) / ^о, (15)
где w0 = А_1Р, w* = А'1¥" .
Последовательно увеличивая перемещение с шагом, определяемым коэффициентом Ху, можно пройти предельные точки и построить всю кривую нагружения. Важным достоинством изложенной процедуры (13)-(15) является то, что она не нарушает квазидиагональную структуру матрицы А.
На завершающем этапе алгоритма по найденным перемещениям в узлах конечноразностной сетки вычисляются другие характеристики НДС оболочки (напряжения, усилия, деформации).
Предлагаемая в настоящей работе численная методика (1)-(15) является развитием подхода, изложенного в [3].
3. Результаты расчетов
При расчетах применяются безразмерные параметры: к = а2/яИ) , р = ра4/ЕИ) , У = у/^о , И = И/И), Г = г/а, Ь = Ъ/а, Мг =Мга2/ЕИ04,
Т =Тв(1+у)а2/И02, V = 0.3, п = ^12(1 -к2) . Здесь
Я , Ь, а - радиус кривизны, внутренний и внешний радиусы оболочки в плане; V, в - коэффициенты Пуассона и температурного расширения; Е - модуль упругости, у - прогиб; Мг- изгибающий радиальный момент, Т - температура.
3.1. Пример 1.
Рассматривается кольцевая пластинка переменной толщины И (г,ф)= И) (1+оо8ф), жестко закрепленная по внешнему контуру г=а под действием на внутреннем контуре г=Ь усилия
Q = Qa3/EИ04 = -8. Пластинка находится в равномерно распределенном температурном поле Т = 3,9. Принимается: Ь = 0,4, И) = 1. В таблице
1 даны значения прогибов и изгибающих моментов на внутреннем контуре пластинки при ф = 0 и ф = п, полученные на трех сетках, и результаты работы [4]. В последних двух колонках приведены число итераций п и соответствующее время ПЭВМ ^ (в секундах), необходимые для получения решения с погрешностью е =10-4. Сходимость итерационного процесса (3) достигается при едином значении т = 0,64 для всех коэффициентов релаксации ти, ^, ту . Из таблицы видно, что уже на сетке 20^20 получается хорошее согласование с результатами работы [4].
Таблица 1
пг X Ф=0 ф=П п Ґ
wmax Мг ,Ъ ^тах Мг ,Ъ
12x10 1,096 3,4752 1,541 1,4624 78 1,0
20x20 1,102 3,5136 1,562 1,5157 92 3,3
30x30 1,103 3,5184 1,570 1,531 100 14,9
[4] 1,104 3,527 1,578 1,541
3.2. Пример 2
Рассматривается задача определения верхних критических нагрузок ("хлопка") жестко закрепленных сферических пологих оболочек под действием неосесимметричного внешнего давления в пределах "полукруга" и "полукольца". На рис.1 для этих задач даются кривые "нагрузка-максимальный прогиб" для трех значений параметра кривизны Л = уЩк = 4; 5; 6 . Сплошные
кривые соответствуют "полукольцу", штриховые
- "полукругу". Результаты получены при е =10-4, И0 =1, ширине "полукольца", равной 0,05 и радиуса 0,5, на сетке 20^20. Ведущим параметром выбран прогиб в точке г=0,5 и ф=п. На рисунке указаны также численные значения критической
нагрузки рсг. На рис.2 приведены формы прогибов в окружном направлении при Х = 4; 6; 8 в момент потери устойчивости оболочки "хлопком" для случая нагрузки в пределах "полукольца".
p -10-
°]
1.2-1
1.0 2.0 ' w
m
Рис.1
О п/2 п 3п/2 2п
О
О
0.7J
w —Я = i
\i6/ i—-—
Рис.2
Видно, что с увеличением высоты подъема оболочки картина волнообразования меняется -происходит перестройка формы прогибов. Значения коэффициентов релаксации ти, ту, тw принимались одинаковыми и равными т в интервале
0,3 < т < 1. Общее число итераций для получения кривой деформирования до достижения прогибом w величины 3Н в случае нагружения оболочки X = 4 в пределах "полукруга" составило 78 итераций.
Заключение
Предложена методика расчета неосесимметричного НДС тонких гибких пологих оболочек вращения и круглых пластинок переменной толщины, находящихся под действием силовых и температурных нагрузок. Методика основана на синтезе методов конечных разностей, матричной прогонки, общей итерации и продолжения решения по параметру. Приводится сравнение результатов ее применения с известными в литературе данными, полученными другими подходами. Показано, что принятый в методике итерационный процесс расчета НДС обладает быстрой сходимостью.
Автор выражает признательность проф. А.А.Аганину за полезное обсуждение результатов статьи.
1. Муштари Х.М., Галимов К.З. Нелинейная теория упругих оболочек. - Казань: Таткнигоиздат, 1957.
- 432 с.
2. Корнишин М.С. Нелинейные задачи теории пластин и пологих оболочек и методы их решения. -М.: Наука, 1964. - 192 с.
3. Шихранов А.Н. Неосесимметричные задачи нелинейного изгиба пологих оболочек вращения и их оптимизация // Образование и наука производству. Сб. трудов Междунар. научно-технической и образовательной конф. Набережные Челны, 2010.
- Ч.1. - Кн.1. - С. 69-71.
4. Григоренко Я.М., Крюков Н.Н., Ахалая Т.Г. Термонапряженное состояние гибких круглых пластин переменной толщины // Прикладная механика. - 1981. - Т.17. - №.7. - С.84-88.
A NUMERICAL TECHNIQUE FOR SOLVING NON-AXIALLY SYMMETRIC PROBLEMS OF NON-LINEARLY BENDING OF SHALLOW SHELLS OF REVOLUTION
0
A.N.Shikhranov
The author offers a numerical technique for evaluating non-axially symmetric stress-stained state of thin flexible shallow shells of revolution and round plates of variable thickness subjected to the force and temperature loads. The technique is based on the synthesis of methods of finite-differences, matrix run, general iteration and extension of solution in a parameter. The comparison of the results of its application with the data obtained by other approaches is presented. The iteration process of computing the stress-stained state used in the technique is shown to have fast convergence.
Key words: rotational shells, non-axially symmetric deformation, temperature loading, stress-strained state, critical loads, wave formation.
Шихранов Анатолий Николаевич - научный сотрудник лаборатории вычислительной динамики сплошной среды Учреждения Российской академии наук Института механики и машиностроения Казанского научного центра РАН.
E-mail: [email protected]
