Исследования устойчивости упругих пластин и оболочек при помощи конечно-элементного моделирования Текст научной статьи по специальности «Строительство и архитектура»

ИССЛЕДОВАНИЯ УСТОЙЧИВОСТИ УПРУГИХ ПЛАСТИН И ОБОЛОЧЕК ПРИ ПОМОЩИ КОНЕЧНО-ЭЛЕМЕНТНОГО МОДЕЛИРОВАНИЯ

Г.А. МАНУЙЛОВ, канд. техн. наук, доц., С.Б. КОСИЦЫН, д-р техн. наук, профессор, М.М. БЕГИЧЕВ, аспирант

Московский государственный университет путей сообщения (МИИТ) 127994, Москва, ул. Образцова д. 9, стр. 9; [email protected]

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

КЛЮЧЕВЫЕ СЛОВА: пластина, сферическая и сфероидальная оболочки, устойчивость равновесия, конечно-элементное моделирование, геометрически нелинейная постановка задачи, малые несовершенства нагружения.

В настоящей работе сделана попытка использовать промышленные КЭ комплекса NASTRAN и ANSYS для решения весьма трудных классических задач устойчивости тонких упругих (неограниченно) оболочек и пластин в линейной и геометрически нелинейной постановке. Практическое использование каждого комплекса имело свои особенности. В комплексе NASTRAN оболо-чечные 8-ми узловые элементы типа plate parabolic (quad8) во многих задачах не работают, а 4-х узловые плоские элементы «plate» (quad4) не во всех моделях оболочек позволяют получать приемлемые результаты. В частности, трудности возникают в задачах расчета непологих строго выпуклых замкнутых и незамкнутых оболочек. Но зато реализация алгоритма Ньютона - Рафсона в комплексе NASTRAN при построении «хлопковых» участков равновесных диаграмм не встречала существенных затруднений. В комплексе ANSYS оболочечные 8-ми узловые элементы типа shell93, как показала практика, работают вполне надежно. Однако построить диаграммы равновесных состояний, имеющих «хлопковые» участки в ANSYSе при помощи алгоритма Ньютона - Рафсона удавалось только для пологих оболочек (с относительно малыми хлопковыми перемещениями).

В качестве наиболее простой задачи в геометрически нелинейной постановке была исследована устойчивость осесимметричного равновесия круглой шар-нирно опертой пластины (Я=50см, 5=0,5см), нагруженной в центре сосредоточенной силой (рис. 1). КЭ-модель пластины имела порядка 900 элементов (plate в NASTRAN и shell93 в ANSYSе). Использовались модели пластины на осе-симметричных сетках и «свободных» сетках. Потеря устойчивости такой пластины провоцируется действием значительных сжимающих окружных усилий Nt, необходимых для уравновешивания радиальных («цепных») усилий Nr, которые развиваются при больших прогибах пластины (порядка 15-20 толщин пластины). Подобная задача при действии распределенной нагрузки и скользящей заделки на контуре была решена Д.Ю. Пановым и В.И. Феодосьевым [1].

Отметим, что предварительное тестирование модели пластины (элементы plate, NASTRAN) со свободными (т.е. «неправильными») и осесимметричными сетками на задачах изгибного равновесия показало совпадение прогибов в центре с теоретическими линейными решениями в пределах 3-х знаков.

Однако в геометрически нелинейной задаче свободные (неосесимметрич-ные) сетки создавали начальные несовершенства, которые существенно снижали критические нагрузки хлопка из подчиненных предельных точек (max

ДРкр~25,8%). Здесь результаты, полученные на NASTRANe и на ANSYSe почти совпадали (92 т и 94,3 т). Идеальная осесимметричная модель пластины при критической нагрузке (Ркр ~ 121900 кг, Wmax ~ 20,95) прощелкивала в результате симметричной неустойчивой бифуркации в устойчивую неосесимметричную четырехволновую форму равновесия (рис. 1).

Рис. 1

Влияние несовершенств создаваемых неосесимметричной сеткой исследовалось на смешанных моделях, когда центральная часть пластины имела свободную сетку, а приконтурная зона моделировалась осесимметричной сеткой (рис. 2, а). Было установлено, что пока радиус зоны свободной сетки не превышает половины радиуса пластины, критическая нагрузка практически не уменьшается по сравнению в бифуркационной (рис. 2, б). Объясняется это тем, что если основная зона сильного волнообразования и значительных окружных усилий N расположена в пределах осесимметричной сетки, то несовершенство свободной сетки слабо влияет на величину Ркр. Существенное падение критической нагрузки наблюдалось, когда зона свободной сетки захватывала наиболее нагруженную приконтурную область пластины (рис. 2, а). Рассмотренная задача интересна еще тем, что она не имеет линеаризованного аналога (в отличие от большинства задач упругой устойчивости). Здесь производную по Фреше нелинейного оператора (т.е. матрицу линеаризации) можно вычислить для предкри-тического равновесия, но ее нельзя получить заранее, до нагружения пластины.

6)

Падение критических нагрузок при увеличении радиуса зоны _«свободной» сетки_

О 5 10 15 20 25 30 35 40 45 50

Радиус участка «свободной»сетки КЭ,см

Рис. 2

Далее исследовалась устойчивость моментного осесимметричного равновесия в геометрически нелинейной постановке для непологой защемленной оболочки, нагруженной сосредоточенной силой в вершине ^ = 500 см, 5 = 0,25

см, Я/З = 1000, а = 353,5 см, Я = 4/12(1 - и2)а/я5 = 57,5 ). В такой оболочке с

59

начала нагружения под силой образуется правильная осесимметричная луночка (рис.За), которая при достижении первой критической нагрузки Ркр1 бифурка-ционно перестраивается, в неосесимметричную вмятину с тремя окружными полуволнами (рис Зб). Далее при достижении второй критической нагрузки Ркр2 происходит вторичная перестройка формы вмятины (с 3-х полуволновой формы на 4-х полуволновую форму). Подобные перестройки формы вмятины наблюдались в экспериментах А.В. Погорелова с очень тонкими сферическими сегментами, полученными вакуумным напылением [2].

Рис. 3

Задача решалась как при помощи комплекса NASTRAN с элементами plate, так и на КЭ комплекса ANSYS (элементы shell-93). Кривые прогибов оболочки в точке приложения силы имеют примерно одинаковый характер (рис. 3, г, д), однако критические нагрузки из решения на ANSYSе получились несколько меньше, чем на NASTRAN (Ркр1 = 41 и 45, Ркр2 = 61 и 67, соответственно). Это объясняется, по-видимому, большей жесткостью многогранной модели из элементов типа plate по сравнению с ансисовской моделью из конечных элементов shell93.

Исследование устойчивости пологого защемленного сферического сегмента (R = 500 см, 5 = 0,5 см, a = 70,72 см, х = 02(1 -u2)a/RS = 8,129) интересно тем, что эту задачу решал В.Л. Якушев [3] динамическим методом; похожую задачу решал и Н.В. Валишвили (X = 8, [4]). В этой задаче наименьшая критическая нагрузка соответствует неосесимметричной форме потери устойчивости, так как параметр пологости больше граничного значения X* = 5,7. Однако по результатам В.Л. Якушева неосесимметричные формы, провоцируемые малыми силовыми возмущениями ±еР, существуют лишь как промежуточные неустойчивые формы равновесия [3, рис. 4,8]. Конечное равновесие после хлопка - практически осесимметричное. Наименьшая критическая нагрузка по В.Л. Якушеву соответствует четырехволновым возмущениям q^ mm = q^« = 4) = 1,94 кг/см2. То же значение q^ min для близкой по геометрии оболочки получил Н.В. Вилишви-ли. При отсутствии неосесимметричных возмущений оболочка хлопком теряет устойчивость в предельной точке (q^ =2,82 кг/см2 В.Л. Якушев). Решение с помощью NASTRAN (900 КЭ «plate») дало при четырех полуволновых и шести полуволновых возмущениях одинаковую минимальную критическую нагрузку quj, = 1,94 кг/см2. С ростом величины возмущающих сил ±еР критическая на-60

грузка в подчиненных предельных точках заметно падала, согласно зависимости, близкой к «закону 2/3» (рис.4).

Г, 2

S

Ъ 1,8

I 1,6

? 1,4

Влияние несовершенств на потерю устойчивости при 4 полуволнах нагружения

+£ß

0,002 0,004 0,006 0,008 0,01

Величина несовершенства } £

Рис. 4

Однако получить критическую нагрузку в предельной точке q^ при отсутствии возмущений ±гР при помощи NASTRANa не удалось. С помощью КЭ модели в ANSYSе вычисленная критическая нагрузка в предельной точке q^ = 2,85 кг/см2, однако, минимальная нагрузка при неосесимметричных возмущениях оказалась меньше результатов В.Л. Якушева и Н.В. Валишвили: q^ min = 1,81 кг/см2 при трех полуволновых силовых возмущениях. Это подтверждает вывод о том, что КЭ-модели оболочек с 8-ми узловыми элементами plate parabolic оказываются менее жесткими по сравнению с КЭ-моделями из 4х узловых plate элементов.

О ОД 0,2 0,3 0,4 0,5 0,6 0,7 0,8 0,9 1

Рис. 5

Задача устойчивости тонкой сферической оболочки при всестороннем гидростатическом давлении предварительно тестировалась по сходимости напряжений в линейной задаче о равновесии (а = qR/25 = 0,5 кг/см2). Использовалась осесимметричная сетка, КЭ - plate в комплексе NASTRAN с различным числом элементов (от 400 до 4800 элементов на полусфере). На рис. 5 показаны результаты вычислений напряжений вдоль меридиана оболочки, из которых следует, что наибольшая погрешность вычисления напряжений оказывается достаточно малой (~1,5-3%), только при весьма мелких сетках (2400 и 4800 КЭ на полусфере). В остальной части полусферы для вычисления указанных напряжений дос-

таточно 1000-1200 КЭ. Объясняется это тем, что вблизи полюса осесимметрич-ная сетка генерирует весьма вытянутые (т.е. «плохие») четырехугольники неправильной формы.

В задаче устойчивости сферической оболочки, при действии гидростатического давления, линеаризованный анализ в комплексе NASTRAN (осесиммет-ричная сетка с 9481 узлом) дал приближенное значение qKp=2,45 кг/см2 вместо «точного» значения q^= 1,2E7R=2,54 кг/см2. Расчет в геометрически нелинейной постановке на той же КЭ - модели дал очень грубый результат (qKp=1,57 кг/см2). Более того, последовательность моделей с меньшим, но прогрессивно увеличивающимся числом узлов, не показала сходимости вычисленных приближенных критических нагрузок к q^ = 2,54 кг/см2.

Теоретическое решение З. Цолли и Л.С. Лейбензона [5] линеаризованной задачи устойчивости сферической оболочки дает одну и ту же наименьшую критическую нагрузку бифуркации как для осесимметричной формы так и для всех неосесимметричных форм нового моментного равновесия. При использовании осесимметричной почти всюду, кроме полюсов сетки КЭ plate parabolic (ANSYS, shell93, опция eigen buckling), осесимметричная (рис 6б) и несколько первых неосесимметричных (рис. 6а, в, г) форм потери устойчивости развиваются при практически одинаковых критических нагрузках q^= 2,525 кг/см2. Геометрически нелинейный расчет дал на этой же сетке диаграмму развития моментного НДС, которая заканчивалась в предельной точке при qmax= 2,5 кг/см2 кривая 1). В этот момент меридиональное волнообразование и частично широтное захватило почти всю поверхность оболочки (рис. 7).

Рис. 6

Отметим, что при использовании свободной сетки КЭ с элементами shell93 в решении геометрически нелинейной задачи, появление моментного НДС оболочки наблюдались в отдельных малых зонах ее поверхности уже при q = 1,12 кг/см2. Далее эти зоны расширялись, и при qmax=2,48 кг/см2 была достигнута предельная точка (рис. 7, кривая 2). Однако хлопок не был реализован в комплексе ANSYS, ни для модели с осесимметричной сеткой, ни для модели со свободной сеткой.

2,542

1,14

Рис. 7

Известно [5], что форма выпучивания безмоментной замкнутой сферической оболочки при действии внешнего давления описывается при помощи сферических функций. Количество полуволн осесимметричного выпучивания, при наименьшей критической нагрузке,

2Е52 1,21Е52

тт Чкр = „ , = = —^— (у =0,3); (!)

R 2л/э(Г7)

R1

можно определить при помощи выражения:

n(a

v к

J = -0,5 ±J 12(1 - v2) + - - 3v

(2)

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

0,605Е5

°КР =

R

< <т„

(V = 0,3).

Для оболочки из Ст3 (апц~2103 кг/см2) это условие выполняется лишь для весьма тонких оболочек (R/S > 600). Согласно (2) при R/S = 600 количество полуволн выпучивания n ~ 44-45. Для исследуемой оболочки R/S = 1000, и соответствующее число полуволн n ~ 57-58. Однако получит такую форму выпучивания оболочки в эксперименте или в результате вычислений долгое время не удавалось.

С другой стороны А.В. Погореловым [2], в связи с разработкой приближенного вариационно-геометрического метода исследования устойчивости строго выпуклых оболочек вращения, высказывалась мысль о том, что малая податливость срединной поверхности оболочки (е ~ 10-3 при опц=2-103 кг/см2) позволяет считать их в первом приближении практически нерастяжимыми. Но если строго выпуклая замкнутая поверхность (или незамкнутая, но с закрепленным краем)положительной гауссовой кривизны (КГ>0) нерастяжима, то она не допускает (согласно теореме К. Гаусса [2]) бесконечно малых и конечных гладких изгибаний. Такие изгибания допускаются лишь в классе кусочно-регулярных поверхностей [2], что фактически узаконивает образование ребра с зоной сильных деформаций и сосредоточенных кривизн в результате местной хлопковой потери устойчивости в предельной точке части поверхности оболочки. Однако невозможность гладких бесконечно малых и малых конечных изгибаний замкнутой срединной поверхности нерастяжимой оболочки означает невозможность бифуркационной потери устойчивости такой оболочки. Докажем это применительно к замкнутой сферической оболочке, нагруженной внешним давлением. Известно, что при бифуркации работа внешней

нагрузки Tq на бифуркационных перемещениях (т.е. прогибах Ws нового мо-ментного равновесия) равна нулю [6]

Tq = Я q^M = 0

Так как qp = const ф 0, то это означает, что объем оболочки в результате

перехода в новое равновесие не должен измениться (AV=0). Значит сумма объемов «вмятин» (EAV+) и «выступов» (EAV.) в результате бифуркационного волнообразования должна равняться нулю

Tq = qKp (£aV++£AV_) = о, ^ Ау = ^AV++YAV_= 0.

Но в этом случае площадь новой бифуркационной поверхности A's должна быть больше площади сферы As = 4жЯ2, поскольку из изопериметрических представлений [7] среди всех поверхностей, охватывающих одинаковый объем, минимальную поверхность имеет сферическая, для которой справедливо равенство A3 = 36%V2. Для всех других поверхностей A^ > 36nV2. Следовательно, новая «бифуркационная» поверхность моментного равновесия «гофрированной» сферической оболочки действительно имеет большую площадь по сравнению со сферой. Но это возможно только, если допустить растяжимость срединной поверхности сферической оболочки. Получили противоречие, которое и доказывает высказанное выше утверждение. Аналогичные утверждения доказываются при помощи изопериметрических представлений и для любого сферического сегмента, закрепленного вдоль контура, с нерастяжимой срединной поверхностью.

Заметим, что высказанное утверждение о невозможности бифуркационной потери устойчивости требует абсолютной нерастяжимости срединной поверхности строго выпуклой оболочки. При «малой», но конечной растяжимости (е ~ 10-3-10-4) бифуркационная потеря устойчивости безмоментного равновесия строго выпуклой оболочки реализуется, что и было получено в результате экспериментов в комплексе ANSYS.

Для выявления влияния начальных несовершенств формы поверхности сферической оболочки на критические нагрузки внешнего давления, была исследована серия сферических оболочек с вмятинами различного диаметра и глубины (рис. 8 а). Вмятины моделировались зеркальным отражением части поверхности сферы относительно секущей плоскости. Угол раствора вмятины менялся от 2° до 30°. Вблизи края вмятины развивалось осесимметричное мо-ментное НДС типа краевого эффекта, которое провоцировало резкое снижение критической нагрузки хлопка из предельной точки.

Рис. 8

Кривая падения указанных критических нагрузок в предельных точках в зависимости от размера вмятины, полученная с помощью элементов shell93 комплекса ANSYS представлена на рис.8 б. Как видно из этого графика, макси-

мальное падение критической нагрузки по сравнению с бифуркационной достигало 90% для оболочек с наиболее глубокой вмятиной (телесный угол раствора 30°).

На рис. 9 представлены результаты геометрически нелинейного расчета и форма хлопковой потери устойчивости оболочки в виде сфероида (a =100 см, a/b = 4, 5 = 0,5см, R/5 = 200), полученная на NASTRAN (1920 КЭ plate) при действии внутреннего давления.

Рис. 9

Величина внутреннего критического давления по А.В. Погорелову [2] для сильно сплющенных сфероидов оказалась равной

Чрнут = 4E521(43(1 -v>2) = 133,2 кг/сж2.

Критическая нагрузка хлопка в предельной точке (за счет действия машинных возмущений) оказалась несколько меньше (112 кг/см2). Форма потери устойчивости, как и предсказывал А.В. Погорелов [2], имеет по экватору многочисленные вмятины, сильно вытянутые вдоль меридианов.

Л и т е р а т у р а

1. Панов Д.Ю., Феодосьев В.И. О равновесии и потере устойчивости пологих оболочек при больших прогибах - ПММ 1948, т. 12, в. 4, с 398-406.

2. Погорелов А.В. Геометрическая теория устойчивости оболочек - М.: Наука, 1966, 296с.

3. Якушев В.Л. Нелинейные деформации и устойчивость тонких оболочек - М.: Наука, 2004, 275с.

4. Валишвили Н.В. Методы расчета оболочек вращения на ЭВМ - М.: Машиностроение, 1976, 278с.

5. Флюгге В.Статика и динамика оболочек - М.: Стройиздат, 1961, 306с.

6. Мануйлов Г.А., Жуков К.А., Косицын С.Б. Метод неособенных продолжений в задачах устойчивости нелинейно деформируемых систем - Строительная механика и расчет сооружений, 1989, №5, с68-72.

7. Бляшке В. Круг и шар - М.: Наука, 1967, 232с.

THE INVESTIGATION OF STABILITY OF ELASTIC PLATES AND SHELLS WITH THE HELP OF FE MODELING

Manuylov G.A., Kositsyn S.B., Begichev M.M.

The features of the finite element modeling and solution of problems of stability of elastic plates, closed spherical and spheroidal shells, shallow and nonshallow spherical shell segments in linear and geometrically nonlinear formulations, with the influence of small imperfections of the loading.

KEYWORDS: plate, spherical and spheroidal shells, stability of equilibrium, finite element modeling, geometrically nonlinear formulation of the problem, small imperfections of the loading.