Разработка и применение метода коллокаций и наименьших невязок к решению задач механики анизотропных слоистых пластин Текст научной статьи по специальности «Физика»

Вычислительные технологии

Том 19, № 5, 2014

Разработка и применение метода коллокаций

и наименьших невязок к решению задач механики

>к

анизотропных слоистых пластин*

С. К. Голушко1,2, С. В. Идимешев1, В. П. ШАПЕЕВ3'4

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

Ключевые слова: метод коллокаций и наименьших невязок, полиномы Чебышё-ва и Лагранжа, композит, краевая задача, теория упругости, теория пластин.

Введение

Настоящая статья является продолжением работы [1] по адаптации метода коллокаций и наименьших невязок (КНН) к задачам механики многослойных композитных пластин и значительно расширяет класс исследуемых задач. Ранее метод КНН был описан и реализован для задач изгиба прямоугольных изотропных пластин, в том числе для пластин, лежащих на упругом основании. Ниже будут рассмотрены прямоугольные пластины, состоящие из произвольного числа ортотропных слоев, направления ортотропии которых являются произвольными в плоскости пластин, что приводит к анизотропии пластин в целом.

Анизотропия и слоистая (неоднородная) структура приводят к сложному распределению напряжений и деформаций в пластинах. Расчёт напряжённо-деформированного состояния (НДС) многослойных анизотропных пластин можно попытаться осуществить, опираясь на пространственную теорию упругости с заданием условий контактов между слоями. При таком подходе при построении вычислительных моделей возникнет ряд трудностей. С одной стороны, малая толщина слоев по сравнению с другими геометрическими размерами пластин потребует применения очень мелких сеток, что ухудшит обусловленность вычислительной задачи. При использовании ячеек большего размера для проекционно-сеточных методов ячейки сетки будут слишком вытянуты в плоскости

1 Конструкторско-технологический институт вычислительной техники СО РАН, Новосибирск, Россия

2Институт вычислительных технологий СО РАН, Новосибирск, Россия

3Институт теоретической и прикладной механики СО РАН, Новосибирск, Россия

4Новосибирский национальный исследовательский государственный университет, Россия Контактный e-mail: idimeshev@gmail.com

*Работа выполнена при финансовой поддержке РФФИ (гранты № 13-01-12032-офи_м и 13-01-00227), гранта Президента РФ по поддержке ведущих научных школ (№ НШ-5006.2014.9).

пластины, что приведёт к ухудшению качества численного решения. С другой стороны, наличие малого параметра в ряде случаев позволяет заранее сформулировать некоторые допущения о характере НДС конструкции и понизить размерность исходной задачи, исключив из рассмотрения направление вдоль толщины пластины. К настоящему времени разработано большое количество уточнённых теорий пластин, использующих различные гипотезы о характере распределения напряжений, деформаций и перемещений по толщине оболочки [2]. В настоящей работе при расчёте анизотропных слоистых пластин помимо классической теории Кирхгофа — Лява [3-5] рассмотрены широко распространённая теория Тимошенко [6-8] и теория ломаной линии [9, 10].

Усложнение математических моделей повышает требования к используемым численным методам. При решении задач механики пространственных конструкций наиболее востребованными остаются конечно-элементный [11-13] и гранично-элементный [14-16] подходы. В представленном исследовании развивается другой проекционно-сеточный метод для численного решения указанного класса задач. Ниже описан модифицированный метод КНН, основанный на специальном расположении узлов кол-локаций и выборе базиса в виде комбинации полиномов Лагранжа, который оказался весьма эффективным при решении краевых задач механики композитных конструкций.

1. Постановка задачи

Рассматривается задача изгиба прямоугольной пластины, состоящей из N изотропных или ортотропных произвольно ориентированных в плоскости пластины слоев (рис. 1). Предполагается, что пластина находится под действием поперечной нагрузки произвольного вида д(ж,у), слои жёстко скреплены между собой, а на торцах пластины заданы соответствующие краевые условия. Необходимо рассчитать напряжённо-деформированное состояние такой пластины.

Рис. 1. Прямоугольная слоистая пластина под действием поперечной нагрузки; параметры а, Ь, Н означают размеры пластины в направлениях х, у, г соответственно, ги — координата верхней границы к-го слоя, к = 1,..., N

В рамках пространственной теории упругости физические соотношения для каждого ортотропного слоя можно записать следующим образом [17]:

пх х хх 01 1 012 0

пу у уу 01 2 022 0

пкх 01 3 023 0

пх у 01 6 026 0

к 13 к 23 к 33 к 36

0Ъ\

026

066/ к

/ек \

XX к

'УУ к 'XX

\еху/ к

024 04!

025

55

где п п п п п п е е е е е е

^ "хх уу' хх) ху ? хх з ух' хх уу ? ^хх) ху ? хх ? ух

жений и деформаций; 0

к

компоненты тензоров напря-коэффициенты матрицы жёсткости к-го слоя в главной системе координат (х,у) (г = 1,..., 6, = 1,..., 6). Коэффициенты выражаются через коэффициенты матрицы жёсткости Ск в материальной системе координат в виде

01 1 012 01 3 016

01 2 022 02 3 026

01 3 023 03 3 036

01 6 026 03 6 066

^к ск (я? )т,

044 04 5

045 05 5

я2к с? (як )т

где

щ

( сое2 вк ят 0

2 в к

ят2 в к

сов

2в

0 - 8т2вк\ 0 вт 2в к 10

\(яп2в к)/2 — (ят 2вк)/2 0 соя 2вк /

Ск

С1 1 С12 С13 0

С12 С22 С23 0

С13 С23 С33 0

0 0 0 С6 6

як

соя в к вт в к — вт в к соя в к

С2

чк С44

0

С5

55

При этом Ск выражаются через технические постоянные следующим образом:

С1 1

С

23

1 ^23^23

е22 Е3к А ;

. ^32 +

" ЕкЕкА , А

-(к _ ^21 + ^31 ^23

С12 =

е2к е2 а

^ к _ 1 — ^кз ^21

С33 =

ЕЕа '

/-<к С13

С44

4 + ^32 Е2кЕ3к А ;

^ к _ 1 — ^к3^3к1 С22 =

ск

23

С

55

ск

31

ее а

с66 = ^12,

1 ^к2^21 ^23^32 ^13^31 2^21^32^23

е2Е2 Е32

Здесь Ек, Ск- , ^к — модули Юнга , модули сдвига , коэффициенты Пуассона к-го слоя

пластины соответственно.

Систему уравнений , описывающую НДС к-го слоя , можно выписать относительно

Iк(ж,у,г) , ик(

функций перемещений ик(х,у, г) , ик(х,у,г) , (х,у, г) в виде

0

д2и к + ок5д2ик + (о к? + о к6) ^и2 + (о к3 + о к5)д2™ к

к дк к д и1 , к д и1 , (глк I к ^ д и2 , ^ к , к

11 дХ + '

(окз + ок6)^ + 066

'66 ду2 д2и1

дг 2 д2и

д2и

дхду д2и

дхду

к дХ22 + 022 ду22 + 0к4 дг22 + (0к3 + 0к4)

дхдг д2и>к дудг

0,

0,

(ок3 + 055) + (ок3 + о4к4) -Й22 + ок5 + 0к4 ^ + ок^ = 0.

дхдг

дудг

дх2

дУ2

дг2

0

0

На границах между слоями потребуем непрерывности функций щ, щ, пк, а^, а^, а^, что отвечает условию "жёсткого" контакта. На верхней и нижней гранях пластины справедливы зависимости

аМ (х,У,Н/2) = (х,У,Н/2) = ауг (х,У,Н/2) = О,

(Х - V2) = (Х - V2) = (Х - V2) = 0.

Для последующего сравнения с результатами работы [18] на торцевых гранях примем краевые условия, соответствующие шарнирному закреплению:

аЖж(0,У,г) = и2 (0,У,г) = П (0,У,г) = аЖж (а,У,г) = и2 (а,У,г) = П (а,У,г) = ° а^у(х, 0, г) = и2(х, 0, г) = и)к(х, 0, г) = а^(ж, Ь, г) = ик(х, Ь, г) = и)к(ж, Ь, г) = 0.

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

Сформулируем постановку исходной задачи в рамках уточнённой теории пластин [9, 10]. Далее воспользуемся следующими обозначениями: м^х, у), и2(х,у), п(х, у) — неизвестные функции перемещений точек срединной поверхности (г = 0), вк(х,у), в!(ж, у) — неизвестные функции угла поворота к-го слоя вдоль направлений х и у соответственно.

Рассмотрим случай, когда слои имеют постоянную толщину. Пусть N * — номер слоя, через который проходит срединная плоскость г = 0. Определим вспомогательные величины, необходимые для выбора срединной плоскости в качестве поверхности отсчёта:

в

где

(м * -1)

= Е рвп П=1

00

Н 0

Р = Н (О

Г*

1, 2,

0 0 0

\Л,1 Л-2 Н

0 0 0

1

0 0

0 0

В рамках теории ломаной линии [9, 10] распределение перемещений по толщине зависит от перемещений срединной плоскости следующим образом:

N

к к п = п, и = и,- —

в* + £ РвП + (г - гк-1)вк, г = 1, 2

П=1

Выражения для деформаций имеют вид

N

Й, = + £ Р^ + (* - гк-1)

А;

х

П=1

х

жу

5и1 дм2 у х

N

+ £ Р

п=1

х

дв^ дУ

N

уу

ди2 + £ рдвП + -гк-,)-двк

П=1

дУ

дУ

дУ

+

дв: дх

^ + (* - *-.)(М + двк

у

х

= вк +

уг

вк + |У,

ду

ек =0

*

к

к

к

к

е

жг

Физические соотношения (1) в рассматриваемых ниже теориях пластин будут иметь более простой вид вследствие принятия гипотезы об отсутствии обжатия е^ = 0. Уравнения равновесия запишем в виде

дТп + =0 + дГк = 0 + _

-/

<9x <9y ' <9x <9y ' <9x <9y

дФ1 дФ1 дФ1 дФ1

дф11 + ^ _ Q = 0, —Фк.2 + _ Q =

<9x <9y ' <9x <9y

где

Здесь Т-, ^ — удельные нормальные касательные и поперечные усилия, Мк, фк — удельные изгибающие и крутящие моменты и их обобщённые представления.

Из соотношений (1)-(4) получена разрешающая система уравнений в кинематических переменных м^х, у), и2(х,у), т(х,у), вк(х,у), вг(х, у), которая из-за громоздкости не приводится.

Порядок получившейся системы дифференциальных уравнений в частных производных будет равен 4Ж + 6. Выпишем краевые условия для теории ломаной линии на примере условий, моделирующих шарнирное закрепление:

ди двк

х = 0, а : т = 0, и = 0, —1 = 0, в2к = 0, = 0, к = 1,...,Ж,

х 2 х

у = 0,6 : т = 0, и = 0, ^ = 0, = 0, = 0, к = 1,...,Ж

у 1 у

Соотношения теории Тимошенко [6] могут быть получены как частный случай теории ломаной линии [9, 10], когда

в1 = вк, в2 = в2, для всех к.

Количество уравнений в этом случае уменьшается, а порядок разрешающей системы дифференциальных уравнений понижается до 10. Переход от теории Тимошенко к классической теории Кирхгофа — Лява [3] можно осуществить, если принять

дт дт

в1 = - "я"", в2 = а . х у

Порядок разрешающей системы дифференциальных уравнений теперь будет равен 8.

0

2. Модифицированный метод коллокаций и наименьших невязок

Идеи и принципы метода КНН берут начало в работах Ашера (Ascher U.), Кристи-ансена (Christiansen J.) и Рассела (Russel R. D.) [19-21], с помощью которых был создан известный пакет прикладных программ COLSYS для решения жёстких систем

нелинейных ОДУ. Далее этому методу были посвящены исследования А.Г. Слепцова, В.П. Шапеева, Ю.И. Шокина и других авторов. В работах [22-25] впервые реализован вариант метода для систем уравнений в частных производных. Варианты с адаптивными сетками в областях с криволинейными границами предложены в работе [26]. В [27] разработан вариант метода КНН повышенного порядка точности численного решения. В методе КНН, описанном в [1], исходная прямоугольная область П=[0, а] х [0, Ь] покрывалась прямоугольной сеткой и решение в каждой ячейке представлялось в виде полинома. Решение задачи во всей области П аппроксимировалось кусочно-полиномиальной функцией на этой сетке. В случае большого количества подобластей при декомпозиции исходной области может происходить замедление сходимости итерационного процесса. Основная идея модифицированного метода КНН заключается в уменьшении количества ячеек сетки за счёт значительно более высокой точности аппроксимации решения в отдельной ячейке. Так, для исходной прямоугольной области достаточно использовать лишь одну ячейку, совпадающую с самой областью. Более высокая точность аппроксимации достигается специальным выбором точек коллокаций в узлах полинома Чебышёва и базисных функций в виде произведений полиномов Лагранжа.

Как и ранее в [1], определим в ячейке локальную систему координат (а, в) такую, что а е [—1,1] и в е [—1,1]:

= х — а/2 = у — Ь/2 а а/2 , в Ь/2 .

Пусть и (а, в) — вектор неизвестных функций ит(а,в). Решение в ячейке представим в виде [28]

¿к (а, в ) = ¿А(а)ЬВ (в) =

А В

(а, в) = ^^ с™¿к(а, в)

г=1 к=1

( А ( Л ( В .. Л

-р-т (а — а к)

А=1 (а — а к)

V к=г /

Д (в — в!)

(вк — в!) \!=к )

где с™ — неизвестные коэффициенты для т-й функции, ¿А(а), ¿В(в) — полиномы Лагранжа степени А и В, а^, вк — координаты точек коллокаций (г = 1,...,А, 3 = 1,...,В).

В модифицированном методе КНН в качестве координат а^, вк точек коллокаций выбираются корни многочленов Чебышёва степени А и В соответственно.

Для определения неизвестных коэффициентов с™ сформируем систему линейных алгебраических уравнений (СЛАУ)

Яст = Ь (6)

и назовем её локальной СЛАУ. Если исходную задачу записать в виде

¿(х,у)[и(х,у)] = /(х,у) для (х,у) е П У)[и(Х у)] = у) для (Х у) е дП

то локальная СЛАУ будет состоять из следующих частей:

Рис. 2. Расположение точек записи уравнений коллокаций и граничных условий в ячейке

— уравнения коллокаций — требования удовлетворения исходным уравнениям задачи в точках коллокации авз

Ь(а,в)[и(а,в)]= I(а, в);

— граничные условия, которые задаются в некоторых точках на границах ячейки (аг,0), (аг, Ь), (0, вз) и (а, вз)

1(а,в )[и(а, в)] = д(а,в).

Схема расположения точек записи уравнений для случая A=B=6 приведена на рис. 2.

Получившаяся СЛАУ (6) является переопределённой. В силу относительно небольшого размера СЛАУ для её решения применяется прямой метод отражений Хаусхол-дера с выбором главного элемента [29].

Покажем преимущества модифицированного метода КНН на примере решения задач изгиба композитных слоистых пластин с известными аналитическими решениями.

3. Ортотропная пластина

Рассмотрим квадратную пластину, выполненную из ортотропного материала: Е1 = 25Е2, 012 = ^13 = 0.5Е2, 023 = 0.2Е2, V = 0.25, Е2=5ГПа. Пластина шарнирно закреплена по всему контуру и находится под действием поперечной равномерной нагрузки q0=0.1 МПа. В рамках теории Тимошенко эта задача решена методом Фурье в работе [17]. Рассчитаем НДС такой пластины модифицированным методом КНН и проведём сравнение результатов численных экспериментов с аналитическим решением из [17].

Для оценки погрешности численного решения будем использовать величину отклонения в точке .

т(х,у) - т*(х,у)

Ег(х, у) =

100%,

т(х,у)

где гш*(х,у) — численное решение, гш(х,у) — точное решение.

В табл. 1 представлены результаты расчётов НДС ортотропной пластины модифицированным методом КНН и сравнение с данными из работы [17] при а = Ь = 1 м, Н = 1/20 м. Из приведённых данных видно, что для совпадения пяти значащих цифр

Таблица 1. Значения и погрешности функций прогиба "(а/2, Ь/ 2) и компонент тензора напряжений ахх(а/2, Ь/2, Н/2), ауу(а/2, Ь/2, Н/2), аху(а, Ь, — Н/2) для ортотропной прямоугольной пластины на разных сетках

Сетка А х В ■ш Ег Охх Ег ауу Ег ^ху Ег

Работа [17] 1.1619-3 - 3.1312+7 - 1.0880+6 - 1.9480+6 -

4 х 4 6.0910-3 >100 9.9815+7 >100 9.8848+6 >100 2.3601+7 >100

8 х 8 1.1080-3 4.64 3.1401+7 0.28 8.6890+5 20.14 1.9159+6 1.65

16 х 16 1.1619-3 0.00 3.1300+7 0.04 1.0908+6 0.26 1.9431+6 0.25

функции прогиба w необходимо 16 х 16 точек коллокаций, расположенных в узлах Чебы-шёва. Компоненты тензора напряжений ахх, ауу, аху не входят в разрешающую систему уравнений и определяются путём пересчёта через кинематические функции, поэтому вычисляются с большими погрешностями.

4. Трёхслойная композитная пластина

Точное решение задачи изгиба шарнирно-опертой многослойной плиты под специальной нагрузкой в рамках пространственной теории упругости получено методом Фурье в [18]. Используя результаты этой работы, проведём анализ полученных численных решений с использованием модифицированного метода КНН в рамках двумерных теорий пластин.

Рассмотрим прямоугольную трёхслойную композитную пластину с ортотропными слоями равной толщины. Размеры пластины в плоскости: а = 1 м, Ь = 3 м. Материал пластины является анизотропным и имеет следующие осредненные характеристики: Е = 172, = 6.8 ГПа, = 0.25, = 3.4, Ск3 = 3.4, Ск3 = 1.3 ГПа (к - номер слоя). Волокна внешних слоев ориентированы вдоль оси х (01 = 03 = 0), внутреннего слоя — вдоль оси у (02 = п/2). Пластина находится под действием поперечной нагрузки вида д(х,у) = д0 вт(пх/а) вт(пу/Ь). В расчётах принималось д0 = 0.1 МПа.

Исследуем зависимость максимального прогиба срединной поверхности от толщины пластины при использовании различных теорий пластин. В табл. 2 приведены величины прогибов в центре пластины, где они для данной задачи принимают максимальные значения. Здесь представлены отклонения величин прогибов, полученные в рамках пространственной теории упругости и при расчётах по разным двумерным теориям. Индексы 3Д, КЬ, Т, ВЬ соответствуют пространственной теории упругости, классической теории Кирхгофа — Лява, теориям Тимошенко и ломаной линии соответственно.

Таблица 2. Сравнение данных расчёта функции прогиба -ш(а/2, Ь/2) трёхслойной прямоугольной пластины в рамках пространственной теории упругости и различных теорий пластин

Б = а/Н "зд "кь "т "вь Егкь Егт Егвь

4 2.618-5 4.680-6 2.200-5 2.544-5 82.12 15.95 2.81

10 1.333-4 7.310-5 1.139-4 1.306-4 45.16 14.55 2.02

20 7.078-4 5.850-4 6.710-4 7.015-4 17.35 5.20 0.89

50 9.427-3 9.140-3 9.230-3 9.421-3 3.05 2.09 0.07

100 7.368-2 7.310-2 7.341-2 7.360-2 0.79 0.37 0.11

Таблица 3. Сравнение данных расчёта напряжений ахх(а/2,Ъ/2, Н/2), ауу (а/2,Ъ/2,Н/6) и аху(0, 0, -Н/2) трёхслойной прямоугольной пластины в рамках пространственной теории упругости и различных теорий пластин

Б = а/Н ахх 3Д аххКЬ аххТ ахх БЬ ЕГКЬ ЕгТ ЕГбь

ахх(а/2,Ъ/2,Н/2)

4 1.824+6 9.968+5 9.808+5 1.576+6 45.35 46.23 13.60

10 7.260+6 6.230+6 6.210+6 7.000+6 14.19 14.46 3.58

20 2.600+7 2.492+7 2.492+7 2.572+7 4.15 4.15 1.08

50 1.570+8 1.558+8 1.558+8 1.568+8 0.80 0.80 0.16

100 6.240+8 6.230+8 6.220+8 6.240+8 0.16 0.32 0.00

ауу (а/2, Ъ/2, Н/6)

4 1.744+5 4.032+4 1.494+5 1.712+5 76.88 14.31 1.83

10 4.180+5 2.520+5 3.740+5 4.100+5 39.71 10.53 1.91

20 1.176+6 1.008+6 1.128+6 1.168+6 14.29 4.08 0.68

50 6.475+6 6.300+6 6.375+6 6.425+6 2.70 1.54 0.77

100 2.530+7 2.520+7 2.510+7 2.530+7 0.40 0.79 0.00

аху (0, 0, -Н/2)

4 4.496+4 1.328+4 3.280+4 4.176+4 70.46 27.05 7.12

10 1.230+5 8.300+4 1.050+5 1.180+5 32.52 14.63 4.07

20 3.720+5 3.320+5 3.520+5 3.640+5 10.75 5.38 2.15

50 2.100+6 2.075+6 2.075+6 2.100+6 1.19 1.19 0.00

100 8.300+6 8.300+6 8.300+6 8.300+6 0.00 0.00 0.00

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

н>, м 1Е-04

5Е-05

0

й =0.1 м

-

** "V чч

¿г ** * * ч \ \\ ч\

0

0.2

0.4 0.6 X, м

0.8

М

6Е-04 4Е-04 2Е-04 0

И =0.05 м

- ГГ- — -.-;

. » ■ » N

* ■чу

0

0.2

Классическая теория Теория

Кирхгофа—Лява Тимошенко

0.4 0.6

х, м

Теория ломаной линии

0.5

Рис. 3. Прогибы срединной поверхности трёхслойной пластины w(x,Ъ/2), вычисленные при использовании различных теорий

уточняет теорию Тимошенко. Следует заметить, что математические постановки задач пространственной теории упругости и двумерные постановки задач теории пластин не являются эквивалентными, т.е. здесь рассматриваются изначально разные задачи. Это связано, в том числе, с постановкой граничных условий, выбором отсчётной поверхности и рядом других факторов. Однако эти различия тем меньше, чем больше значение параметра $.

м 0.04

0.02

А = 0.1 М

-0.02

-0.04

V' (|

1» II

ъ

■ с^Па |1

- |Н Я , ,

г, м 0.02

к =0.05 м

-0.02

-500000

0 500000

(7>у, Па

-1е+06

Классическая теория ______Теория

Кирхгофа—Лява Тимошенко

0 1е+06

Оуу. Па

Теория ломаной линии

Рис. 4. Распределение напряжений ахх(а/2, Ь/2, г) (а) и ауу (а/2,Ь/2,г) (б) по толщине трёхслойной пластины, рассчитанных при использовании различных теорий

а

В табл. З представлены результаты расчётов для компонент тензора напряжений axx(a/2, b/2, h/2), ayy(a/2, b/2, h/6) и axy(0, 0, -h/2), на рис. З, 4 — соответственно прогибы срединной поверхности и распределение напряжений axx(a/2, b/2, z) и ayy(a/2, b/2, z). Как следует из приведённых данных, для тонких пластин все три теории пластин дают близкий результат.

Заключение

В работе предложен и реализован модифицированный метод коллокаций и наименьших невязок для расчёта напряжённо-деформированного состояния прямоугольных анизотропных слоистых пластин, находящихся под действием поперечных нагрузок произвольного вида. Особенность модифицированного метода КНН заключается в специальном выборе точек коллокации в корнях полинома Чебышёва и использовании базиса в виде прямого произведения базисных полиномов Лагранжа. Это позволяет получить решение с высокой точностью, используя лишь одну ячейку для прямоугольной области. Сопоставление найденных численных решений с известными решениями других авторов показало высокую точность численных решений, полученных с помощью модифицированного метода КНН.

Проведено сравнение характеристик НДС в задачах изгиба прямоугольных и анизотропных слоистых пластин, рассчитанных в рамках пространственной теории упругости, классической теории пластин Кирхгофа — Лява, теории Тимошенко и теории ломаной линии Григолюка — Куликова. Показано, что для тонких и очень тонких пластин полученные значения характеристик НДС хорошо согласуются для всех теорий, что позволяет использовать упрощённые теории при решении практических задач. Однако для толстых и средней толщины пластин теории Кирхгофа — Лява и Тимошенко искажают характер распределения напряжений и смещений в слоистых конструкциях и в этом случае необходимо производить расчёт НДС в рамках либо теории ломаной линии, либо пространственной теории упругости.

Список литературы

[1] Голушко С.К., Идимешев С.В., Шапеев В.П. Метод коллокаций и наименьших невязок в приложении к задачам механики изотропных пластин // Вычисл. технологии. 2G13. Т. 18, № б. С. 31-43.

[2] Голушко С.К., Немировский Ю.В. Прямые и обратные задачи механики композитных пластин и оболочек вращения. М.: Физматлит, 2GG8.

[3] Новожилов В.В. Теория тонких оболочек. Л.: Судпромгиз, 1951.

[4] Reissner E., Stavsky Y. Bending and stretching of certain types oh heterogeneous aeolotropic elastic plates // J. of Appl. Mech. 19б1. Vol. 28. P 4G2-4G8.

[5] Лехницкий С.Г. Анизотропные пластинки. 2-е изд. М.: Гостехиздат, 1957. 4б3 c.

[6] Григоренко Я.М., Василенко А.Т. Задачи статики анизотропных неоднородных оболочек. М.: Наука, 1992.

[7] Mindlin R.D. Influence of rotatory inertia and shear on flexural motions of isotropic, elastic plates // J. of Appl. Mech. 1951. Vol. 18. P. 31-38.

[8] Reissner E., Stein M. Torsion and Transverse Bending of cantilever Plates. Technical Note 23б9. National Advisory Committee for Aeronautics. Washington, 1951.

[9] Григолюк Э.И., Куликов Г.М. Многослойные армированные оболочки. М.: Машиностроение, 1988.

[10] Григолюк Э.И., Чулков П.Л. Нелинейные уравнения тонких упругих слоистых анизотропных пологих оболочек с жестким заполнителем // Изв. АН СССР. Механика. 1965. № 5. С. 65-80.

[11] РикАрдс Р.Б. Метод конечных элементов в теории оболочек и пластин. Рига: Зинатне, 1988.

[12] Dey P., Sheikh A.H. A new element for the analysis of composite plates // Finite Elements in Analysis and Design. 2014. Vol. 82. P. 62-71.

[13] Zhang Y.X., Yang C.H. Recent developments in finite element analysis for laminated composite plates // Composite Structure. 2009. Vol. 88. P. 147-157.

[14] Бреббия К., Теллес Ж., Вроубел Л. Методы граничных элементов. М.: Мир, 1987.

[15] Albuquerque E.L., Aliabadi M.H. A boundary element analysis of symmetric laminated composite shallow shells // Comput. Methods in Appl. Mech. and Eng. 2010. Vol. 199, is. 41-44. P. 2663-2668.

[16] Dirgantara T., Aliabadi M.H. Elastoplastic boundary element method for shear deformable shells // Eng. Structures. 2012. Vol 45. P. 62-67.

[17] Reddy J.N. Mechanics of Laminated Composite Plates and Shells. CRC Press, 2003.

[18] Pagano N.J. Exact solution for rectangular bidirectional composites and sandwich plates // J. Composite Materials. 1970. Vol 4. P. 20-34.

[19] Ascher U., Christiansen J., Russel R.D. A collocation solver for mixed order systems of boundary value problems // Math. Comput. 1979. Vol. 33. P. 659-679.

[20] Ascher U., Christiansen J., Russel R.D. Collocation software for boundary value ODEs // ACM. Trans. Math. Software, 1981. Vol. 7, No. 2. P. 45-53.

[21] Russel R.D., Christiansen J. Adaptive mesh selection strategies for solving boundary value problems // SIAM J. Numer. Anal., 1978. Vol. 15. No. 1. P 59-80.

[22] ПлясуновА А.В., Слепцов А.Г. Коллокационно-сеточный метод решения нелинейных параболических уравнений на подвижных сетках // Моделирование в механике. Новосибирск: ВЦ, ИТПМ СО АН СССР. 1987. Т. 1(18), № 4. С. 116-137.

[23] Сёмин Л.Г., Слепцов А.Г., Шапеев В.П. Метод коллокаций-наименьших квадратов для уравнений Стокса // Вычисл. технологии. 1996. Т. 1, № 2. С. 90-98.

[24] Слепцов AT., Шокин Ю.И. Адаптивный проекционно-сеточный метод для решения эллиптических задач // Докл. АН. 1996. Т. 347, № 2. С. 164-167.

[25] Слепцов А.Г., Шокин Ю.И. Адаптивный проекционно-сеточный метод для эллиптических задач // Журн. вычисл. математики и матем. физики. 1997. Т. 37, № 5. С. 572-586.

[26] Беляев В.В., Шапеев В.П. Метод коллокаций и наименьших квадратов на адаптивных сетках в области с криволинейной границей // Вычисл. технологии. 2000. Т. 5, № 4. С. 12-21.

[27] Исаев В.И., Шапеев В.П. Варианты метода коллокаций и наименьших квадратов повышенной точности для численного решения уравнений Навье — Стокса // Журн. вычисл. математики и матем. физики. 2010. T 50, № 10. C. 1758-1770.

[28] Бабенко К.И. Основы численного анализа. 2-е изд., испр и доп. М.; Ижевск: РХД, 2002. С. 847.

[29] Деммель Дж. Вычислительная линейная алгебра. Теория и приложения. М.: Мир, 2001.

Поступила в 'редакцию 1 сентября 2014 г., с доработки —17 сентября 2014 г.

36

C. К. ronymKO, C. B. HgHMemeB, B. n. WaneeB

Development and application of collocations and least residuals method to the solution of problems in mechanics of anisotropic laminated plates

Golushko Sergey K.1'2, Idimeshev Semyon V.1, Shapeev Vasiliy P.3'4

Purpose. This work is aimed at the development of a modified method of collocation and least residuals and its application to solution of problems in mechanics of multilayered composite plates.

Methodology. The basic idea of the modified method of collocations and least residuals lies in reduction of the grid cells number due to the higher accuracy of approximation of the solution in a single cell. For example, in the case of a rectangular initial area it is sufficient to use only one cell which coincides with the whole area.

Higher accuracy of the approximation is achieved by a special choice of the collocation points in the nodes of the Chebyshev polynomial along with the basis functions in the form of a multiplication of Lagrange polynomials.

Findings. The problems of bending of rectangular orthotropic and anisotropic laminated plates under transverse loads of any kind are both considered and solved. Solutions are obtained within the classical Kirchhoff—Love's laminated plate theory, Tymoshenko's theory and the Grigolyuk —Chulkov's theory of a broken line.

A comparison of the numerical solutions with known analytical solutions is carried out for the case of orthotropic and anisotropic plates within 3D theory of elasticity under the special fixing and load.

It is shown that the calculated characteristics of the stress-strain state are in good agreement for all the four theories considered above for the case of thin and very thin laminated plates, so the simplified theory is allowed to solve practical problems. Originality/value. For plates of medium thickness and thick plates, Kirchhoff—Love theory and the theory of Timoshenko distort the nature of the distribution of stresses and displacements in layered structures, so application of a more accurate theory becomes necessary. The presented modified method of collocation and least residuals showed high efficiency for solution of bending of multilayered anisotropic plates problems applied to the canonical squared and rectangular forms.

There is reason to believe that it would be no less effective in dealing with a wider class of problems in the mechanics of composite structures, including multi-layer shells. Keywords: method of collocations and least residuals, Chebyshev polynomials and Lagrange, boundary value problem broken line theory, theory of elasticity, plates theory.

Received 1 September 2014 Received in revised form 17 September 2014

1 Design Technological Institute of Digital Techniques SB RAS, 630090, Novosibirsk, Russia

2Institute of Computational Technologies SB RAS, 630090, Novosibirsk, Russia

3The Institute of Theoretical and Applied Mechanics SB RAS, 630090, Novosibirsk, Russia

4 Novosibirsk State University, 630090, Novosibirsk, Russia Corresponding author: Idimeshev Semyon V., e-mail: idimeshev@gmail.com Aknowlegements: This work were supported by RFBR (project No. 13-01-12032-ofi_m and 1301-00227), grant of Russian Federation President for support of leading scientific schools (NSh-5006.2014.9).