Научная статья на тему 'Оценка точности конечно-разностного метода Рунге-Кутты 4-го порядка при решении задач динамики тонкостенных оболочек методом конечных элементов'

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

CC BY
229
16
i Надоели баннеры? Вы всегда можете отключить рекламу.

Аннотация научной статьи по физике, автор научной работы — Муравьев В. В.

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

i Надоели баннеры? Вы всегда можете отключить рекламу.

Похожие темы научных работ по физике , автор научной работы — Муравьев В. В.

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.
i Надоели баннеры? Вы всегда можете отключить рекламу.

Estimation of Precision of 4-th Order Finite Difference Method by Runge-Kutta in Solving Problems of Thinwalled Shells Dynamics by Finite Element Method

A way of analysis of the numerical solution of problems of thin-walled shells dynamics is suggested, which implies the application of finite element method of and 4-th order finite difference method by Runge-Kutta. Functions predicting the numerical solution and having the form of analytical solution functions are constructed. Recommendations on selecting the step of integration for time are developed. For free and forced oscillations the error of the obtained solution is estimated. The numerical solution predictions are compared to the analytical solutions and results of numerical experiments.

Текст научной работы на тему «Оценка точности конечно-разностного метода Рунге-Кутты 4-го порядка при решении задач динамики тонкостенных оболочек методом конечных элементов»

3.Волков И. К., Загоруйко Е. А. Исследование операций. - М.: Изд-во МГТУ им. Н.Э. Баумана, 2000. - 436 с.

4. П у г а ч е в В. С. Теория вероятностей и математическая статистика. - М.: Наука, 1979. -496 с.

5. Chebotarev A. N., Vinogradova E. I., Volkov I. K. Effect of the state of the geomagnetic field and solar activity on chromosome aberrations frequencies dynamics // Russ. J. Num. Anal. Math. Modelling. - 2003. - Vol. 18, № 6. - P. 467-478.

6. Альберт А. Регрессия, псевдоинверсия и рекуррентное оценивание. - М.: Наука, 1977. -224 с.

7. Зельнер А. Байесовские методы в эконометрии. - М.: Статистика, 1980. -438 с.

8.Т ю р и н Ю. Н. О предельном распределении статистики Колмогорова-Смирнова для сложной гипотезы // Изв. АН СССР. Математика. - 1984. - Т. 48, № 6. - С. 1314-1343.

9. З а к с Ш. Теория статистических выводов. - М.: Мир, 1975. - 776 с.

Статья поступила в редакцию 1.06.2006

Елена Ивановна Виноградова родилась в 1980 г., окончила МГТУ им. Н.Э. Баумана в 2002 г. Старший преподаватель кафедры "Математическое моделирование" МГТУ им. Н.Э. Баумана. Автор 5 научных работ в области идентификации математических моделей по данным эксперимента.

Ye.I. Vinogradova (b. 1979) graduated from the Bauman Moscow State Technical University in 2002. Senior teacher of "Mathematical Simulation" department of Bauman Moscow State Technical University. Author of 5 publications in the field of identification of mathematical models according to experimental data._

УДК 519.62:519.622.2:531.312

В. В. Муравьев

ОЦЕНКА ТОЧНОСТИ КОНЕЧНО-РАЗНОСТНОГО МЕТОДА РУНГЕ-КУТТЫ 4-ГО ПОРЯДКА ПРИ РЕШЕНИИ ЗАДАЧ ДИНАМИКИ ТОНКОСТЕННЫХ ОБОЛОЧЕК МЕТОДОМ КОНЕЧНЫХ ЭЛЕМЕНТОВ

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

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

численными методами, в том числе с использованием метода конечных разностей по времени. Возможно использование условно устойчивых схем (методы центральных разностей, Рунге-Кутты, Вильсона) и безусловно устойчивых (методы Хубольта, Адамса-Мултона) [1], [2]. Вопросам точности и сходимости посвящены работы Б. Парлетта [3], О.Зенкевича [4], Г.Н.Кувыркина [5], Г.И. Марчука [6].

Широко используется явная условно устойчивая схема метода Рунге-Кутты 4-го порядка. Важно выбрать такой шаг по времени, который обеспечит на заранее известном интервале времени требуемую точность получаемого решения.

Рассмотрим сплошное неоднородное тело объемом V, плотностью р, на которое действуют нестационарные объемные силы /, а на части поверхности $1 — нестационарные поверхностные силы /и. На части полной поверхности заданы кинематические граничные условия. Перемещения считаем малыми, деформации упругими. Необходимо определить напряженно-деформированное состояние с известной точностью на интервале времени. Перемещения щ и их первые производные по времени в начальный момент известны. Воспользуемся вариационным принципом возможных перемещений в формулировке

ац5et]dV - J J J (-püi)8uidV-

(V) (V)

игШБ - JJJ ^пЖ = 0. (1)

(V)

Представим тело конечными элементами. Записывая формулу (1) для выбранного конечного элемента и применяя стандартные процедуры метода конечных элементов, получаем следующее уравнение движения:

[К] ш = - [М] Ш + } , (2)

где [К] — матрица жесткости; [М] — матрица масс конечно-элементной модели; (д} — вектор узловых перемещений, (^} — вектор узловых сил. Начальные условия

(д}|4=о = (до}, Ш^о = Ш. Представим матрицу [М] согласно работе [4] в виде

[М] = [Ь] [Ь]т,

где [Ь] — нижняя треугольная матрица, определенная, например, по алгоритму Холецкого [3]. Умножим обе части равенства (2) слева на [Ь]-1, получим

[Ь]-1 [К] ([Ь]т)-1 [Ь]т (д} = - И-1 [Ь] [Ь]т(д} + [Ь]-1 (^} . (3)

Так как ([Ь]т)-1 = ([Ь]-1)т, то запишем [Ь]-1 [К] [Ь]-1Т = [В], где [В] - симметричная матрица с ленточной структурой, соответствующей ленточным структурам [М] и [К]. Обозначим (д} = [Ь]т (д} и = [Ь]-1 (^}. Учитывая, что [Ь]т(д} = |д|, перепишем выражение (3) в виде

[В] (д} = - {$} + {. (4)

Для простоты дальнейшего изложения символ" над векторами (д} и указывать не будем.

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

Для дифференциального уравнения вида (д} = (/ (¿, (д})} соотношения метода Рунге-Кутты 4-го порядка записываются для вектора узловых перемещений и вектора узловых скоростей. Дифференциальное уравнение второго порядка представляется в виде системы двух дифференциальных уравнений первого порядка:

Ш = (V},

(5)

(V} = (/ (*, (д})}.

Для рассматриваемого случая имеем вектор (/ (¿к, (дк})} = = — [В] (дк} + (^к}. В общем виде соотношения метода Рунге-Кутты 4-го порядка для системы (5) представлены в работе [2]. Предположим неизменность вектора } на к-м шаге по времени. Опуская выкладки, выразим значения векторов узловых перемещений и узловых скоростей на следующем шаге по времени:

{qk+i} = ( [B]2 ^ - [B] Ar ) {qk} +

{qk+i} = [E] + [B ]

6

+ [E] + [B ]

AT 4

AT 4 24

- [B]

At 2\

2

) fe} +

24

+ ([E] At - [B] AT3 ){Fk} ,

At 2 \

- [B] — J {qk} +

+ ([E] At - [B] AT3 ) {qk} +

+ ([E ]AT2 - [B]At4 ) {Fk} .

2

24

Входящую в формулы (4) и (6) матрицу [В] можно представить в соответствии со спектральной теоремой в виде [3]

[B] = [G] [A] [G]

i

(7)

где [С] — матрица форм колебаний, причем [С]т = [С]-1; а [А] — матрица собственных значений V квадратной матрицы [В] размерности п х п, ] = 1, 2, ...п — номер собственного значения. Собственные значения матрицы [В] удовлетворяют условию det([K] — V [М]) = 0. Подставим равенство (7) в систему (6) и введем новые векторы (дк},

[С]-1 {&} и

{qk} и {Fk},

такие, что

{Qk} = [G] 1 {qk}, {q4 =

Fk }• = [G] 1 {Fk}. Уравнения (6) примут вид

, At 3 л .

Qk+i, = i- v At )qk, +

> At 4

+ 1 + vj--v,—

1 j 24 j 2

At 2\

Qk, + At - v

At 3

Fk

(8)

2 At

Qk+ij = ( 1+

4

- v,

Ar)

Qkj +

At3\Q /At2 At4\ Q + I At - v, — ) qj + ( —--v,-^- ) fj.

6

2

24

Отсюда следует, что исследование выражения (6) на сходимость есть исследование на сходимость п систем (8) или п одномерных задач.

2

6

Для простоты дальнейшего изложения символ ~ над векторами указываться не будет, как и индекс ] у компонентов вектора узловых перемещений.

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

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

А^ = Да

+ Дфазы + Д

частоты •

(9)

амплитуды I "фазы

Рассмотрим задачу о свободных колебаниях. Погрешность в первую очередь будет определяться расхождением частот получаемого численного и аналитического решений и затуханием амплитуды колебаний численного решения. Будем искать численное решение системы

(8) в виде дк = д (4) = и дк = д(^) = ,

где Д— — отклонение частоты численного решения от аналитического.

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

После подстановки в выражение (8) дк и дк получаем систему

' +Лт) =

= 1 - v

De

Дт2

2

+ v2

'Дт4

24

C + Дт - V.

Дт3

D e

l(^ÜJ+Aw)(tk+Ar) =

= V

,Дт3

6

- V. Дт C + 1 - V

Дт2 ' 2

+ v2

'ДГ4

24

De'

(10)

После сокращения обеих частей второго уравнения на ег^Ьк находим ДегА^к. По формуле Лагранжа для остатка ряда можно

Дт2 . Дт3 о Дт4

записать е'^Ат = 1 + i^Vj"Дт -_Дт 5

V.

2

- iV з^з-^- + V

24

+

+ ZV.

V.-e'^^1, где G [0, Дт]. Опираясь на численные экс-

120

перименты, считаем, что Д— много меньше a/V.. Представляем

6

егДадДт _ i + ¿AwAr. Тогда будет справедливо соотношение

( V2 /V"Ar5 \

" v " _eVj«1 + AwAreVjЛт

1__120

л AT3

v/V" At - VV"

у V " 6 у

В знаменателе дроби пренебрегли слагаемыми высшего порядка ма-

Ат3

лости относительно У'^/У! . Подстановка полученного выражения в первое уравнение системы (10) и сокращение правой и левой частей на ег(^+Аш)^ приводит к уравнению

, ч Ат 2 Ат з

еЧ^+д-)Ат = 1 + г^Ат - у - гу^-Г- +

2 ^ ' 6 Ат 4 Ат5

+ V? — - ¿*? /V"—е^1 - *А—Ате^Дт,

которое с использованием формулы Лагранжа для остатка ряда позволяет получить соотношение для А-—:

Aw = , (11)

2 Ат4 е^«1 + ег^«2

^ 120 2ег^ Дт

где € [0, Ат]. Из равенства (11) видно, что в первом приближении действительная часть А— отрицательна и по модулю не превышает

у?/^ Ат 4 .

значения -, а мнимая часть А— положительна и по моду-

120 ^ уЗАт 5

лю не превышает значения -. Отсюда следуют выражение для

120

численного решения и предельные значения погрешностей:

% = д ) = Се-^^^Ь, (12)

у3Ат5

А амплитуды 1 (13)

120

у?/У! Ат4

А = _£ (14)

"частоты 120 V1 V

Помимо того, ошибки в определении амплитуды и фазы получаемого решения могут появиться при записи конечно-разностной схемы (8) для первого шага интегрирования по времени. Если аналитическое решение имеет вид д = С0 соз(/у1£ + £0), а численное дк = д (£к) = С соэ(^4 + то будет справедливым записать начальные условия (для £ = ¿0 = 0)

С0 соэ(£0) = С1 соэ(^1),

/У1Со 81п(^о) = э1п(^1).

Используем основное тригонометрическое тождество и тот факт, что ^ есть сумма ^^ и действительной части выражения (11). Получим для величины С равенство

/ V 4Дт 8

С? = С? (1 - ^ 81П2(^О)

v4Дт 8

Начальная погрешность амплитуды имеет порядок --, что

2 * 120

значительно меньше погрешности по равенству (13) при £ = Дт. Ошибка определения фазы имеет тот же порядок малости, что и погрешность амплитуды. Пренебрегаем ею, сравнив с погрешностью по соотношению (14).

Предполагая, что на интервале времени [0, ¿пред] погрешность не превышает предельного значения Дпред и, пользуясь ее малостью, исходя из выражений (9), (13) и (14) получаем неравенство

V2^ Ar

.4

120 ДТ + ^ £пред < Дпред. (15)

Если нас интересует только точность нахождения амплитудных значений перемещений, то получаем неравенство

V 3Дт 5

3

120 £пред ^ Дпред. (16)

Из неравенств (15) и (16) следуют ограничения, накладываемые на Дт.

Рассмотрим задачу о вынужденных колебаниях. Решение системы (2), соответствующей системе (8), зависит от вида функции = ^). Пусть функция ^) имеет вид ^(£) = С э1п(ш£), когда

ш = . Частное решение неоднородного дифференциального урав-

С

нения имеет вид д (£) = --зт(ш£). Внешнюю нагрузку можно

С

представить в виде F(£) = — (ег^ — е-г^). Частное решение неодно-

2

родного дифференциального уравнения от F(£) = Сег^ есть функция С

д(£) =--. Соответствующее ему численное решение системы

V —

(8) будем искать в виде дк = д ) = . Подставляя в первое и

второе уравнения, получаем

,Дт 3

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

+Дт) = A ^v2

- v,Ar +

j 6 j

( Ar4 Ar 2\ / Ar 3\

+A i 1 + v2— - Vj ——— j + C (Ar - v,— j

A^Mík+Дт) = a 1 + v2

24

At 2

- V,

+

+A AT - v.

At 3

~6~

+ C

At2 2

- V,

At4

24

Теперь умножим второе уравнение системы на ф и приравняем правые части. После приведения подобны слагаемых, получаем соотношение

A - V, И At - V,

At 3

gí^tfc =

= C -At + V.

At 3

+ г ^

At 2

At 4

Из этого соотношения следует, что ^ = ш и A=

C

24

1+г-

At

если опускать величины высшего порядка малости относительно

т. Полное решение системы (8) имеет вид дк = д(£к) = +

2 _ _

+ + . Начальные условия (£к0 = 0) точно соответ-

ствуют аналитическому решению, поэтому

Ao = - г

ш

At c

4 V2 - ш2

ш \ шАт C

1 + — и Bo = -г—--2-2

W 4 V2 - ш2

1

ш

Полное решение системы (8) представляется в виде

C

qfc =

1 + г-

.шАт

1 - 2(1 + /V1^-

1

A1 /Vle

. (17)

Последнее выражение показывает, что у перемещения относительная погрешность по амплитуде

Аамплитуда ^

шАт

ш

1 + — V

и по фазе

Афаза < шАт.

(18)

(19)

В частном случае при F (t) = C sin (wt) и w = //Vj частное решение неоднородного дифференциального уравнения имеет вид

C п q (t) = --1 sin ( //Vjt — — ). Численное частное решение будем искать в виде qk = q (tk) = и qk = q(tk) = . Подстановка в

и

2

6

6

2

2

2

систему (8) выражений для перемещений, скоростей и внешней силы, приводит к системе

f'

ВегМ^+Ат) = а/+ в/ег^й — —

V '

/' / — 1

АегМ^+Дт) = А/ег^ — — Вег^ — -_,

V? V? '

, , 2 Дт4 Дт2 Г. 2 Дт3 Л „ где / = 1 + V?2—--V? —г—, а /' = V2—--V? Дт. Выражая из первого

24 2 6

уравнения системы В и подставляя его во второе уравнение системы, получаем, что^ = и

А (V? (е — /)2 + /'2) = с / — (/ — 1) (е^т — /^ .

Опустим слагаемые высшего порядка малости относительно V? Дт2. Тогда

А = —60-^- (1 — - V? Дт Л / (1 + - V? Дт2 ^Дт4 V 12 3 ) / V 36 ?

Полное численное решение можно представить с учетом равенства (11)в виде

. 120 ^ -г л/^ 120 ^

дк = д(^) = Аег^ + Асе V ) + Вое

Начальные условия соответствуют аналитическому решению, а пото-

Дт4 + ^ 1

му полное численное решение при —1-< 1 представляется

120

как

120^- (1 — V?Дт2' v3Дт 4 V 12 ?

д- «-?--—о-^-. (20)

1 \ /V2/^Дт4 \ К }

1 + — ^-Дт2 вт -

36 ? У V 240

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

Погрешности численного решения в определении амплитуды

^Дт 8

1440

Другой частный случай, когда ^(¿) = С. Частное решение неС

однородного дифференциального уравнения д (¿) = —. Соответствующее ему численное решение системы (8) будем искать в виде

дк = д (¿к) = А и дк = д(4) = 0. Подстановка в уравнения (8) дает С

величину А = —. Численное решение совпадает с аналитическим.

V

Д амплитуда ^ ,, , „ ^fc. (21)

Рис. 1. Расчетная схема квадратной шарнирно опертой по углам пластины

Для проверки условий (12)-(14), (17)-(21) решалась задача о поперечных колебаниях квадратной шарнирно опертой по углам пластины. Ее расчетная схема приведена на рис. 1.

Пластина выполнена из однородного изотропного материала со следующими характеристиками: модуль упругости Е = 73,1 ГПа, коэффициент Пуассона V = 0,3, плотность р = 2821 кг/м3. Размер пластины в плане Ь = 0,305 м, толщина к = 0,00328 м.

Оценивалось влияние плотности конечно-элементной разбивки на точность определения собственных частот и форм. Результаты сравнивались с решением, приведенным в работе [7]. Удовлетворительная точность была достигнута на конечно-элементной модели, содержащей 81 узел и 128 элементов Зенкевича [8]. Узлы распределены по пластине равномерно. Численно методом Якоби на 50 итерациях было получено значение первой частоты собственных колебаний пластины 61,05 Гц, которому соответствует собственное значение VI = 147140,1 с-2.

Для проверки соотношения (12) решалась задача о свободных колебаниях в главных координатах. Начальные условия: вектор узловых перемещений совпадает с первой формой собственных колебаний, вектор узловых скоростей равен нулю. Поперечные перемещения центральной точки пластины были найдены несколькими способами и приведены на рис. 2. Они были получены в численном эксперименте с шагом интегрирования по времени Дт = 0,002 с. Этим результатам соответствует на рис. 2. график "численное решение". Перемещениям, найденным из аналитического решения и по соотношению (12), соответствуют графики "аналитическое решение" и "прогноз численного решения", соответственно. Результаты говорят об удовлетворительном соответствии наших предположений о рассогласовании частоты, фазы и амплитуды в случае свободных колебаний результатам численных

Рис. 2. Зависимость перемещений от времени в случае свободных колебаний:

о — численное решение; • — прогноз численного решения; х — аналитическое решение

экспериментов даже при относительно больших шагах интегрирования по времени (« 10 % периода). Последнее объясняет значительное расхождение численного и аналитического решений.

Соотношение (17) проверялось при частоте возбуждения, близкой к резонансной. Решалась задача о вынужденных колебаниях в главных координатах. Внешняя нагрузка, приведенная к первой главной координате, имеет вид = 18^(1,01/^). Внешняя нагрузка, приведенная к остальным главным координатам, имеет вид = 0, г = 2, 3,... п. Поперечные перемещения центральной точки пластины были получены несколькими способами. Во-первых, они были найдены из аналитического решения — график "аналитическое решение" на рис.3. Во-вторых, перемещения были получены в численном эксперименте с шагом интегрирования по времени Дт = 0,001 с — график "численное решение" на рис. 3. В-третьих, они были найдены по соотношению (17) при Дт = 0,001 с — график "прогноз численного решения" на рис. 3. Совпадение результатов, полученных по соотношению (17) и в численном эксперименте, удовлетворительное.

Соотношение (17) проверялось при частоте возбуждения, значительно отличающейся от частоты резонанса. Решалась задача о вынужденных колебаниях в главных координатах. Внешняя нагрузка, приведенная к первой главной координате, имеет вид = 1 8^(2/^). Внешняя нагрузка, приведенная к остальным главным координатам, имеет вид = 0, г = 2, 3,... п. На рис. 4 приведены поперечные

1,50Е-05

tD С

а bf в и

&0.00Е+00

tD '

-1.00Е-05

-5,00Е-06

5,00Е-06

1.00Е-05

-1,50Е-05

время, с

Рис. 3. Зависимость перемещений от времени в случае синусоидальной внешней нагрузки с частотой, близкой к резонансной:

□ — численное решение;--аналитическлое решение; х — прогноз численного

решения

перемещения центральной точки пластины, полученные несколькими способами. График "аналитическое решение" показывает перемещения, которые были найдены из аналитического решения. Графики "численное решение 1" и "численное решение 2" соответствуют перемещениям из численных экспериментов с шагами интегрирования по времени Дт = 0,001 с и Дт = 0,0001 с соответственно. По соотношению (17) для Дт = 0,001 с и Дт = 0,0001 с были построены графики "прогноз численного решения 1" и "прогноз численного решения 2" соответственно. Из представленных на рис. 4 результатов видно удовлетворительное соответствие соотношения (17) численному эксперименту в широком диапазоне шагов интегрирования по времени.

Соотношение (20) проверялось при частоте возбуждения, равной резонансной частоте. Решалась задача о вынужденных колебаниях в главных координатах. Внешняя нагрузка, приведенная к первой главной координате, имеет вид Fl(t) = 1sin(^Vii). Внешняя нагрузка, приведенная к остальным главным координатам, имеет вид F(t) = 0, i = 2, 3,... n. Поперечные перемещения центральной точки пластины были получены, во-первых, из аналитического решения — график "аналитическое решение" на рис. 5. Во-вторых, в численных экспериментах с шагами интегрирования по времени Дт = 0,001 с и

Рис. 4. Зависимость перемещений от времени в случае синусоидальной внешней нагрузки с частотой, в два раза большей резонансной:

□ — численное решение 1; х — прогноз численного решения 1;----прогноз

численного решения 2;--аналитическое решение; о — численное решение 2

Рис. 5. Зависимость перемещения от времени при синусоидальной внешней нагрузке с резонансной частотой:

□ — численное решение 1; х — прогноз численного решения 1;----прогноз

численного решения 2;--аналитическое решение; А — численное решение 2

Дт = 0,002 с — графики "численное решение 1" и "численное решение 2". Помимо того, по соотношению (20) для Дт = 0,001 с и Дт = 0,002 с были построены графики "прогноз численного решения 1" и "прогноз численного решения 2". Видно удовлетворительное соответствие спрогнозированных значений результатам численных экспериментов.

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

СПИСОК ЛИТЕРАТУРЫ

1. Постнов В. А. Численные методы расчета судовых конструкций. - Л.: Судостроение, 1977. - 280 с.

2. Б у с л о в В. А., Яковлев С. Л. Численные методы решения уравнений: Курс лекций. - СПб: Кафедра вычислительной физики физического факультета СПбГУ, 2001.-44 с.

3. Парлетт Б. Симметричная проблема собственных значений. - М.: Мир, 1983.-384 с.

4. Зенкевич О., Морган К. Конечные элементы и аппроксимации. - М.: Мир, 1986.- 318 с.

5. К у в ы р к и н Г. Н. Термомеханика деформируемого твердого тела при высокотемпературном нагружении. - М.: Изд-во МГТУ им. Н.Э. Баумана, 1993. -354 с.

6. М а р ч у к Г. И. Методы вычислительной математики. - Новосибирск: Наука, 1989.- 608 с.

7. R e i d R. Comparison of methods in calculating frequencies of corner supported rectangular plate // NASA technical Node D - 3030.

8. Зенкевич О. Метод конечных элементов в механике. - М.: Мир, 1975. -541 с.

Статья поступила в редакцию 29.06.2006

Василий Викторович Муравьев родился в 1979 г., окончил в 2003 г. МГТУ им. Н.Э. Баумана. Аспирант кафедры "Ракеты-носители и космические аппараты" МГТУ им. Н.Э. Баумана. Специализируется в области численных методов расчета многослойных тонкостенных конструкций.

V.V. Muravyev (b. 1979) graduated from Bauman Moscow State Technical University in 2003, Post-graduate of "Rocket and spacecrafts" department of Bauman Moscow State Technical University. Specializes in the field of numeric methods of calculation thin multilayer skins.

i Надоели баннеры? Вы всегда можете отключить рекламу.