Наука и Образование
МГТУ им. Н.Э. Баумана
Сетевое научное издание
ISSN 1994-0448
Наука и Образование. МГТУ им. Н.Э. Баумана. Электрон. журн. 2014. № 12. С. 858-877.
Б01: 10.7463/1214.0738649
Представлена в редакцию: 10.11.2014
© МГТУ им. Н.Э. Баумана
УДК 539.3
Количественная оценка погрешности математической модели Власова для пологой сферической оболочки
Константинов М. В.
1,*
Jkon&t'gbmatejn£kju :МГТУ им. Н.Э. Баумана, Москва, Россия
Сферическая оболочка, как совершенный в весовом отношении сосуд, используется в космических аппаратах, где тонкостенные элементы и оболочки объединяются рамами. Очевидно, что избежать локального воздействия на оболочку и, следовательно, концентрации напряжений в ней невозможно. Стремление к весовому совершенству космического аппарата ведет к уменьшению коэффициента запаса прочности составляющих его частей. Уменьшение запаса прочности возможно только при определении напряженно-деформированного состояния элементов с контролируемой погрешностью. Известные универсальные численные методы конечного элемента и конечных разностей не могут справиться с такой задачей. Задачу необходимо решать аналитически. Одно из решений основано на математической модели Власова для пологой оболочки. Недостаточная обоснованность критерия пологости в части количественной оценки его влияния на погрешность решения создает трудности применения модели Власова. В работе этот недостаток устранен в результате вычислительного сопоставительного анализа решений модели Власова и математической модели механики деформирования оболочки, линейные дифференциальные уравнения которой получены с погрешностью, не превышающей погрешности допущений Кирхгофа теории оболочек. Впервые определены величины погрешности, обусловленной принятыми в теории Власова упрощениями, в зависимости от характера нагрузки и геометрических параметров оболочки. Показано, что критерию пологости Власова соответствует погрешность не более 35%. С увеличением изменяемости напряженного состояния погрешность снижается. В случае локального нагружения оболочки в месте максимальной концентрации напряжений величина погрешности 3-10% (в зависимости от степени локализации нагрузки). Варьирование значений R/h от 100 до 1000 влияния на погрешность модели Власова не выявило. Результаты работы позволяют использовать математическую модель Власова для определения напряженно-деформированного состояния тонкостенных элементов с априори заданной погрешностью, например, при научно-исследовательских и опытно -конструкторских разработках аэрокосмических систем.
Ключевые слова: сферическая оболочка, математическая модель Власова механики деформирования, критерий пологости, погрешность
Введение
Одной из упрощенных математических моделей механики деформирования сферической оболочки является математическая модель, полученная Власовым В.З. [1] из уравнений общей теории при принятии дополнительных гипотез.
Степень пологости оболочки принято оценивать в соответствии с известным определением Власова: «оболочка, очерченная по части сферической поверхности и имеющая подъем, не превышающий 15 диаметра основания оболочки, как показывают исследования, может быть отнесена к числу пологих оболочек».
Указанное определение сформулировано на основе имеющегося опыта и не вполне обосновано в части количественной оценки погрешности, привносимой принятыми дополнительными гипотезами.
В данной работе на примере задачи деформирования сферической оболочки при локальном воздействии аналитически определены значения внутренних силовых факторов с использованием различных математических моделей: модели Власова для пологой оболочки и эталонной модели. В качестве эталонной принята модель [2,3], построенная с погрешностью допущений Кирхгофа теории оболочек, в развитие результатов исследований Института механики АН УССР [4]. Выполнен вычислительный сопоставительный анализ решений математических моделей. Впервые определены величины погрешности, обусловленной принятыми в теории пологих оболочек упрощениями, в зависимости от характера нагрузки и геометрических параметров оболочки.
1. Математическая модель Власова
Рассмотрим изотропную сферическую оболочку постоянной толщины к с радиусом кривизны срединной поверхности Я в географической системе координат (3, ф), рис. 1.
Примем обозначения в соответствии с книгой [5]: е1, е2, ®, Х1, Х2, 1 - параметры деформации оболочки, связанные с перемещениями и, и2, из точек ее срединной по-
Рис. 1. Географическая система координат
верхности по направлениям ортов е - меридиональному, е2 - окружному, е3 - радиальному; ^, 02 - углы поворота нормали к поверхности оболочки в проекциях на плоскости ее и ее , соответственно; N, N, £, й, & - меридиональное и окружное усилия, усилие сдвига и поперечные силы; Мх, М2, Н - изгибающие и крутящий моменты.
Воспользуемся методом Фурье разделения переменных - в предположении симметричности внешней нагрузки на оболочку относительно начала координаты ф, разложим параметры напряженно-деформированного состояния оболочки в тригонометрические ряды вида
ад ад
/(3,ф) = ^/(3)ес8(кф), /(3,ф) = ^/(3)мп(Кф), (1)
К=0 К=1
соответственно для симметричных и антисимметричных по ф величин. Это позволит свести рассматриваемую математическую модель к системе обыкновенных дифференциальных уравнений с постоянными коэффициентами.
Кроме этого выполним переход к безразмерным величинам по схеме
г ^ г Я, и ^ Ц Я, г = 1,2,3, £ ^ ~Ек, Н ^ НЕк2, N ^ N Ек, & ^ & Ек,
Мг ^МгЕк2, хг ^ Хг/Я, г = 1,2, т^ т/Я, где г = Я в1п(3) - радиус параллели.
Тогда, из [1] при отсутствии поверхностной нагрузки на оболочку, получим:
- разрешающие уравнения теории пологих оболочек
V2^2^к - Vф» = о, бУ2У 2кихк + ЕЯ2 V2у« = о, (2)
где ук - вспомогательная функция, введенная Власовым В.З.; Б = Ек3/12(1 -д2) - цилиндрическая жесткость; V = й2/йг2 + й/гйт - к2/~2 - оператор Лапласа;
- связь Ц д. и ^ с внутренними силовыми факторами
~ Г1 й к2 ^~ ~ йг й (к ~ ^
#1,к = — - , ^2,7 = ~Т^Т , = ~[г 1|/Д
йг1 к йг ^ г к ^ /о = Б кй = Б к
<2,к П2 ~ ¥ ки 3,к , <2,к п2 ~ ¥ ки 3,к ;
ЕкЯ ~ ЕкЯ ~
геометрические уравнения
(3)
= ^г«+иък, в2,к = 1(ии + ки~2,к+~из,к), « =-1 (*" йи йг г г
- + и2,к
йг
= %, Х2,к = 4(еи + ке2,к), ~к = -^к+е2,к), еи =-^, е2,к = киъм;
(4)
7--' ? /V/,« 5 \ 1, К 2,к У у К ~ V 1,К 2,к У ? 1 К 7—
аг г г , йг - физические уравнения
= (в1,к +МН2,к), #2,к = Г^Т(в2,к 1,к), М1,к +^Х2,к), (5)
1 -д 1 -д Ек Я
или
M2k + ДГ1.k), Sk = 1 , Hk = (1
Eh R 2(1 + д) Eh R
~ ~ ~ ~ R ~ ~
£1,k = N1.k - Д~2. k , S2, k = N2. k - N , Xl.k = 12 R (M~i,k - M k ) ,
' . . . . . h
(6)
R _ _ R ~
X2, k = 12- (M~2. k - ДЦ.k ), Qk = 2(1 + д)~ , ~k = 12(1 + д) - k. h h
В результате ряда преобразований, которые ввиду громоздкости здесь приводить не
будем, получим решение уравнений (2)
2 ^2
Vk = УCjJk(А/) + C3.k~k, f~3.k =-ЁCjk(V) + C4.k~k, (7)
1=1 1=1
где А2 = ±ik»; k»2 = EhR2¡D ; i - множитель «мнимая единица»; C ■ ^ - произвольные постоянные интегрирования; Jk - функция Бесселя первого рода, определяемая [6] как сумма степенного ряда
ад / 1 \m
Jk (z) = У--(z/2)2m+k , (8)
k У Г(т + k + 1)Г(т +1) ' W
где Г(z) - гамма-функция Эйлера, имеющая аналитическое представление [6,7], например, в виде предела Эйлера-Гаусса
Г( z) = lim-:--(9)
z(z +1)...(z + n)
и, в случае действительного z , являющаяся естественным распространением факториала на вещественные значения аргумента.
Решение (7) позволяет определить остальные величины задачи. Так из (3), например, получим
2 (л А (1 d k2 ^
~ 1 d k
Nu =У 1 ~d-Jkj + СзЛ j=1 у r dr r J
У r dr r j
rk =
2 d2 2
-j, k
1=1
уd~2 1 j
= -УCjk -rrr + А2 Jkj -C3,kk(k- 1)~k-2. (10)
~ 2 d2J (А,r) .
N~2.k =У Cj. k ^1 ) + C^k (k - 1)rk
При этом использовано уравнение Бесселя [6-8] в форме
d2 Jk (А j ~) 1 dJk (А /) ^ 2 k 21
------+
-
у r J
Jk (А/) = 0. (11)
dr ~
Аналогично, из (4) получаются выражения для компонент изгибной деформации ^к, у^2к,
~, подстановка которых в (5) позволяет определить моменты Мх к, М2к, Нк . Из (6) с учетом (10)
=-Z С
j ,k
1=1
d2
(1 + |д) — + X2,
d~
Jk (X /) - Сз,к (1 + (k -1)~
'J2,k У C1,k 1=1
d~2
(1 + +
Jk (X j ~) + C3,k (1 + ^)k(k - 1)~
k-2
k-2
С другой стороны, в1к и в2к определяются геометрическими уравнениями (4). Это позволяет найти перемещения Цк и и 1к. После подстановки (7) и (12) в первые два уравнения
(4) получим
dUu dr
= -(1+ Ц)
2 d2J, (X ,Г)
У q k kr21)+С3 kk (k -1)~
k-2
1=1
Ux k = (1+ Д)
2 f d2
- У с k У сk
k 1=1
-+X2.
vd~2 1y
d~2
Jk (X1-) + СзЛ (k -1)~
- С r
с 4,k
k-1
- c4t 1 r k+1 -1 .
4,k k k u
Откуда, с учетом (11),
Uhk =-(1+ Д)
UXk = (1 +
2 dJ. (X -)
^ k ( 1 )+c3kkrk-1
У с
1=1
1 ,k
d-
k 2
r У qj (X /)+C3,kk~k-1
1=1
~ 1 r ~ 1 ~
k+1
k+1
Здесь постоянная интегрирования, соответствующая перемещению сферы, как жесткого целого, опущена.
Для разделения мнимой и действительной частей решения воспользуемся тригонометрической формой представления комплексного числа и формулой Муавра. Выразим
Л, следующим образом: Л, = д/КУ[со8(утс/4) ± г в1п(утс/4)]. Тогда ряд, определяющий
функцию Бесселя (8), можно представить в виде
Jk (X /) = Re[Jk (X /)] + i Im[ Jk (X 1~)] = JR (~) ± J (~):
I- / \2m+k
~Vk*/2f
1дс j (f ) = y-—--—'--cos
(13)
=0 Г(т + k + 1)Г(т +1)
(2m + k)
Jk (~) =у-
I- / \2m+k
-^—-—'--sin
(2m + k)
т=0 Г(т + к + 1)Г(т +1) Заменим комплексные произвольные постоянные интегрирования С к по схеме
СК ^ (С1,к + гС2,к V2 для . = 1 , СУ.к ^ (С1,к - гС2,к V2 для . = 2 . Здесь С1,к и С2,к - действительные постоянные интегрирования. С учетом этого из (13) получим
у сл J k (X1 ~) = с JR (~) - с2, kj[(r), ЕЕ q k X2 Jk (X 1~)=-k* С J (~)+q k JRk(~)]
1=1
1=1
у см X41 Jk (X 1~) = -k*2 [cj (~) - С2Л J k(-)].
(14)
1=1
Подстановкой (14) выражения для искомых величин приводятся к окончательному виду. Например (10)
—-СиЛ)"-с2дл7)"-—*(СЛ + с/)] - сикк - 1)~к-2 = — -С-Л)" - Л ] +С 2.—[Л)" +—Л ] - сищ - 1)~к-2,
где (*)" — а2(*)/с~2.
В итоге, в матричной форме,
^.к=[ВДСк}, а, к=а к]{Ск}, к=[а~, к]{Ск}, i—1 , 2, ~к=[~]{Ск},
Нк — [йк]{Ск}, и,— — [¿~,к]{Ск}, i — 1.2.3, еи = [еи]{Ск}, О5)
е,к = [еа]{Ск}, X,— = [ЪкКСк}, i — 1. 2, «к = [®—]{Ск}, т = [т—]{Ск},
где
[ к] =
[ ^~2.к]
[~к] = [<2и]
М —] =
[М~2. к ] [Н~к] =
и к] --
[Ц~2. к] -Й. к]
:[- Л ) + кЛ Л ) + —Л - к(к - 1)~к-2 0]; = [/ )" - (Л )" к(к - 1)~к-2 0];
_ . Л ^
~ Л/ I к(к - 1)~к-2 0
к.
— [Л )' - Л )' 0 0]; [0,к] - ~[- л; Л 0 0];
Я к
Я к
1 - д
1 - д
1 - д
дЛ -V1(Л)" -д/к -V1(Л)" 0 -^к(к-1)~
Г—-2
к*
+л)" -/ + лг 0 1 -~к-2
1 -д Я
к* к -(1+д) (1 + д)
к*
Гк-л
к*
1-тдк (к -1)~
К*
г-Л
1
0 — к(к -1)~
к*
к-2
Л)' - (/)' —
--Л ~/ -—~к-1
11
-к+1
1 + д к +1
1 1 1 + д к +1
-к+1
— —
/к7 Л 0 1 ~к
К;.
1
, [ех. -]— -—*[(/)' (Л)' 0 -кгк-1];
К*
[е2. -]—-
-Л -/Я 0 -к~к-1
к ~ и к Г Г
к*
[еи] — [- (1 + д)Л)" +к */— (1 + д)Л )'+—/ - (1 + д)к(к - 1)~к-2 0]; [в,—] — [(1 + д)Л )"-дк Л— - (1 + д)(Л) " -дк Л (1 + д)к(к-1)~к-2 0];
V
V
*
к ] = 2(1 + Д)
Л
- 3 I-
-л
к (к - 1)~к -2 0
[хи] = к*
\-~2,к ] = к* [ х- ] = к*
1
-(Л) -ЛТ о --к(к-1)?
-к-2
к*
(л-)"+к*Л (Лу--Л о 1 к(к- 1)~к-2
/С*
л
_ . 3 г-
к
1
_лк | о -—(—-1)?
К*
?к-2
{С- } = [Си С 2— Сз- с4,-]т; (*) ' = ^ (*)/ .
Поскольку при построении математической модели [1] для пологой оболочки фактически использована гипотеза о преобладании в сечениях оболочки напряжений от усилий над напряжениями от моментов, для сдвигающей силы нет необходимости выполнять приведение по Кирхгофу, т.е.
~ = ~, о*- = 0.-+-
нк.
(16)
Я Бт(3)
Алгоритм решения краевой задачи следующий:
1) вычисляют значения — из (7), функций Л , Лк из (13) и их производных (произ-
водные получают из рекуррентного соотношения алк (г)
к (г )/бг = Л-х(г) - Л-+1(г)]/2
г)]/ 2 и урав-
нения (11));
2) вычисляют значения величин (15)-(16) до постоянных интегрирования;
3) используют граничные условия и формируют систему линейных алгебраических уравнений (СЛАУ);
4) определяют постоянные интегрирования из СЛАУ;
5) вычисляют значения искомых величин задачи;
6) вычисляют к -ю частичную сумму рядов Фурье искомых величин;
7) оценивают погрешность результатов;
8) в случае достижения заданного значения погрешности заканчивают вычисления, иначе - повторяют вычисления для определения следующей частичной суммы.
V
V
V
V
2. Эталонная математическая модель
В качестве эталонной примем математическую модель Виноградова-Менькова [2], опубликованную в [3].
После разделения переменных методом Фурье (1) в [2,3] выполняют комплексное преобразование уравнений равновесия общей теории оболочек и их упрощение путем отбрасывания малых членов. Под малыми понимают члены порядка отношения к!Я по сравнению с единицей, что позволяет оставаться в пределах погрешности гипотез Кирх-
гофа. В результате указанных действий получают упрощенные уравнения равновесия в комплексной форме
а2 Nк/а32+^(3) а~к/а3 - [ Я ~+—7в1П2 (3)]т~— = о, /й3 + 2^(3) ыи +~ к/ 81П(3) - ^(3) N = 0, й§к /й3 + 2^(3)^ + Л~и к/81п(3) - N к/8т(3) = 0,
или
а2 N /аз2+^(3) аык /аз - [ Я ~+к 2Дт2 (3)]Т^ = о,
+ [2с1в(3) + к/мп(3)1~- - [с1в(3) + к/^п(3)]^ = 0, (17)
#2,к/аз+[2с1в(3) - к/81п(3)]/2,к - [с1в(3) - к/мп(3)]^ = 0,
где £ = - ¡к л/12(1 -Д2); N = N к + N2-; / — = N — + ~; -Л- = N — - . Здесь знак «тильда» обозначает комплексную форму величин. Решением уравнений (17) [2,3] являются
ад
N = Си 1Ек(3/2)Б[~,1 -~;1 + к;81П2(32)], /и = СиХвк(3/2)£Ъ,.,=+— (32),
1=0
к-2.
/2,к = С1Л 1вк(32)]Г /=_к ^п2/ (32) + С2, !1Пк+22(32)
(18)
,, /=-к V I , А к совк+2(32)'
где ~ = 1 + 1 - 4 Я/£ У2 ; С14. , С2к - произвольные постоянные интегрирования; Б - гипергеометрическая функция Гаусса [7], которая определяется как сумма степенного ряда
Р( а ,Ъ, с; г) = У Г(а + 1)Г(Ъ + 1>Г(с> * ' ' и Г(а)Г(Ъ)Г(с + 1 )Г(/ +1)
что, в данном случае, равносильно
_ я1п(ла) Г(~ + / )Г(1 - ~ + /)Г(1 + к)
(19)
Р[~,1 - ~;1 + к;в1п2(32)] = У (3/2), ~ = ^П^ Г(~ + 1)Г(1 - ~ + 1)Г(1 + к)
^ 1 ^г I 1т- I 1 »Г/ 1 I 1\
1=0 к Г(1 + к +1 )Г(/ +1)
где Г(г) - гамма-функция Эйлера (9); учтено, что Г(~)Г(1 - а) = л/slп(лa) .
Для определения функциональных коэффициентов степенных рядов, соответствующих /1,к и /2 к , в [2,3] получено рекуррентное соотношение
~ = ~, Ъ~ = (1 +1 )а/ - - (1 +1)~-1]. (20)
0 2 +1 + к 0' 1 2(1 +1) +1 + к
Внутренние силовые факторы и деформации срединной поверхности оболочки определяются [2,3] с помощью разрешающих функций (18) следующим образом:
N1- = ^Яе/ + /2,к), N2— = 1ке(2Л~к - Лк - /и), Я- = 1ке(Лк - Л-),
1 Б
яК
63
к Б
(21)
01— =^-т^ -¡к , 02- = ^М^)
Я slп(3) —
мхк= D Im
Л*
N -1(1 /и)
м 2,k = D Im
К*
Hk = -(1 - DМЛ k -Л k), ^k = -^Re
2Л л Л*
+1(1 -K>C/, k + .Л, k)
— (1+ц)(/, /2, ,)
2
D ^
4k = 72Re k*
N -1(1 + ^)C/,k + /k)
1
> Xu = Im(2Nk - Â,k - hk) = 2k*
X2, k=^M.~, k+./2, k)> ®k=(1k- J~2, k)> ^^^мл k- J~2, k)=
2Л* Л* 2kjls
k
^k " ^k , ôi,k = ôi,k + n . . Пч Hk ■
R sin(3)
где k» = 4DËh
dNk= Cu tgk (3/ 2)J —F[/,1 - /;1 + k ; sin2 (3/ 2)] +
d3
sin(3)
1 ~(1 - ~) 2 1 + k
+ -
sin(3)F[1 + /,2 - /;2 + k ; sin2 (3 2)]}
Значения перемещений определяются [2,3] формулами
Uu =
R sin2(3)
k2 -1
a1,k - 1R sin(3)a2,k k
U2,k =
R sin2(3) I k2 - cos2(3) k (k2 -1) { sin(3)
R sin2 (3)
U3,k = l2
R sin2 (3)
" k2 -1
k2 -1 1
S2,k
- R sin(3)X2,k - COs(3)«1,k 1
Rx 2 k ^2 k + , R COs(3)Œ2 k k
^ S2,k X2,k
ctg(3)+R °1k- sin3)a"
где k =-ч k ctg(3)+k/d3- k®k/ sin(3) ; «2,k = ^k -®k/R •
(22)
(22)
Представим постоянные интегрирования в виде C-к = £ .+ /Ç.^ и введем обозначе-
ния
Ли = tgk (3 2) F[ / , 1 - /;1 + k ; sin2 (3/ 2)],
Ä2,k = tgk (32)£ j+k sin2j (32), Kk = tgk (32)S j=-k sin2j (32),
j=0 j=0
тогда можно записать
Nk = [1 i][Nk]{Ck}, /тЛ = [1 /][/m,k]{Ck}, m = 1,2,
где
[ Nk ] =
Re(\k ) - Im(\k ) 0 0 bn(Au ) Re(\k ) 0 0
, / ] =
Re(Л2,k ) - Im(Л2,k ) 0 0 МЛ2, ) Re(Л2,k ) 0 0
(23)
[ f2,k ] =
sink-2(Ь/ 2)
Яе(Лз,к) - Im(A 3,)
1т(Лз,) Яе(Лз,)
cosk+2(Ь/ 2) 0
sink-2(3/ 2)
, {Ck} 4
Cu
^2, k C 2, k
cosk+2(Ь/ 2) _
Выбранная форма записи позволяет выполнить разделение мнимой и действительной частей решения и для (21)-(22) записать
Nhk = [Nik]{Ck}, Qhk = Q]{Ck}, M,k = [Mk]{Ck}, i = 1, 2, Sk = [Sk]{Ck}, Hk = [Hk]{Ck}, 8i,k = [Si-k]{Ck}, Xi,k = [Xi,k]{Ck}, i = 1,2,
где
% = [®k]{Ck}, Tk=[Tk]{Ck}, uhk = [Uhk]{Ck}, i = 1,2,3, eu = [euщ},
[S*] = [Sk]{Ck}, [Qu]=[Qu]{Ck},
[Nu ] = ^[1 0][[fhk ] + [f2,k ]]; [ N2,k ] = ^[1 0][2[Nk ] - f ] - [fi,k ]];
(24)
[Sk ] = ^[1 0]f ] - [fi,k ]]; [Quk ] =1D [0 1]
2
Rh
d
т[ Nk]
[&,<-] = [0 Nk]; [M,.k]=D[0 1]
R sin(b) k* k*
[ Nk ] - i(1 -v)lfu ] + fk ]]
2
[M 2,k ] = D [0 1] k*
[Si,k ] = - D [1 0]
[s2,k ] = 0]
m[Nk] + ^(1 -^)t/1,k] + [fxk]] ; [Hk] = -(1 -m)D[0 1]f ]-[f2,k]];
2k D
M[Nk]-^(1 + M)f ]+[f2,k]] ; [®k] = (1 + m)£[1 0]fk]-[f2k]]; 2 _ k*
1 1
[Nk]--a+M)fc/i,k]+f ]] ; [Tk]=-^[0 1]f ]- fk]];
[Xu ] = ^[0 1][2[Nk ] - f ] - [f2k ]]; [X2,k ] = ^[0 1]f ] + [fi,k ]];
2k 2k
—[s2k ] = D [1 0]
db 2 k k2
Т [ Nk ] -1(1 + M) db k 2
.„ [f1,k ] + Jr. [f2,k ] db db
d k 1
[a1,k] = -[s1,k] ctg(b)+—[s2,k] —г-т- К]; [a2,k] = [Tk] - -К];
db sin(b) R
n_ Rsin2(3) ] = k2 -1
R sin2(b)
[a1 k ] - 1 R sin(b)[a2 k ] k
U k ] k(k2 -1)
n_ Rsin2(3)
k ] = k2 -1
k2 - cos2(3) sin(b)
[S2,k ] - R sin(b)[X2,k ] - c0s(b)[a1, k]
1
R[X2k ] - [S2k ] ^ R cos(b)[a2k ] k
0
R sin2(£)
[yu J -
к2 -1
J - [^2,к J
1 к ctg(3) + - К, кJ - —— [«2, кJ R sin(^)
[Q*> J - [Qu J[н J;
R srn($)
а[*]/- аналог [*] из (23), составленный из значений производных соответствующих функций.
Алгоритм решения краевой задачи следующий:
1) вычисляют значения X из (17), ~ из (18) и — из (21);
2) с помощью соотношений (19)-(20) в форме (23) вычисляют до постоянных интегрирования значения разрешающих функций (18) и их производных (производные определяют из (17) и (21));
3) вычисляют значения величин (24);
4) используют граничные условия и формируют систему линейных алгебраических уравнений (СЛАУ);
5) определяют постоянные интегрирования из СЛАУ;
6) вычисляют значения искомых величин задачи;
7) вычисляют к -ю частичную сумму рядов Фурье искомых величин;
8) оценивают погрешность результатов;
9) в случае достижения заданного значения погрешности заканчивают вычисления, иначе - повторяют вычисления для определения следующей частичной суммы.
3. Вопросы программной реализации решений математических моделей
Полученные решения рассматриваемых математических моделей могут быть программно реализованы средствами известных математических пакетов, например, Mathcad, MATHLAB или Mathematica.
Указанные математические пакеты являются коммерческими. Кроме этого, в силу своей архитектуры, имеют невысокое быстродействие. В этой связи оправдано создание программного модуля с использованием языка программирования высокого уровня (различные диалекты C, Pascal и т.п.), позволяющего работать со структурированными типами данных, что дает возможность относительно просто реализовать математические операции теории функций комплексного переменного. В данной работе использован язык программирования Free Pascal.
В ходе программной реализации полученных решений, выполнено исследование сходимости степенных рядов, соответствующих разрешающим функциям.
Так для J (X. исследование показало, что для упрощенной оценки количества N учитываемых членов ряда могут быть использованы зависимости N - lg(R/h)^0,2Rh + 70 при 0 < г < 0,7 и N - lg(Rh^0,3 Rh +140 при 0,7 < Г < 1, обеспечивающие вычисление J (X. г) с погрешностью не более 10-6.
Аналогичный уровень погрешности для рядов, определяющих Ык, /1к и /2к, эталонного решения дает учет 1§(Я/к)^/038 Я/к + 600 членов ряда. В общем, число членов ряда, дающих приемлемую погрешность вычисления, уменьшается с увеличением номера гармоники к и/или уменьшением меридиональной координаты 3 . Так для 0<3<45° количество учитываемых членов может быть уменьшено до ^(Я/к)л]0.2 Я/к +100 без роста указанного значения погрешности.
Оценки N количества учитываемых членов ряда позволяют упростить алгоритм вычисления разрешающих функций, исключив проверки сходимости соответствующих рядов.
Полученные решения математических моделей механики деформирования сферической оболочки предполагают вычисление значений гамма-функции Эйлера. Учитывая быстро возрастающий характер гамма-функции, в вычислениях целесообразно использовать логарифм гамма-функции определяемый асимптотической формулой Стирлинга [7]
1п[Г(г)] — (г-1/2)1п(г)-z + 1пл/2тс+1/12г . Можно показать, что при Яе(г) > 15 данная
формула обеспечивает погрешность менее 10-6. Для соблюдения указанного условия всегда есть возможность увеличить действительную часть числа г до 15 и воспользоваться рекуррентной зависимостью 1п[Г(г)] — 1п[Г(г +1)]- 1п(г) .
4. Вычислительный эксперимент
Погрешность аналитического решения математической модели механики деформирования пологой сферической оболочки определяется принятыми при ее построении дополнительными, помимо допущений Кирхгофа теории оболочек, предположениями:
- о преобладании деформаций относительного удлинения и сдвига срединной поверхности над деформациями искривления и кручения;
- о гладкости и малом отличии срединной поверхности оболочки от плоскости (пологость).
К аналогичным дифференциальным уравнениям приводит теория напряженных состояний с большой изменяемостью, где основные разрешающие уравнения получаются исходя из предположения о характере напряженного состояния оболочки - оно должно быть быстро изменяющимся хотя бы по одной из координатных осей [10,11]. Это подтверждается работой [12], где показано, что в пределе (при увеличении локализации и переходе к сосредоточенной нагрузке) напряженно-деформированное состояние оболочки практически перестает зависеть от кривизны ее срединной поверхности.
Следовательно, максимальные значения погрешности математической модели Власова ожидаемы в условиях малой изменяемости напряженного состояния оболочки.
Рассмотрим сферическую оболочку с параметром Щк = 100 и пологостью, характеризуемой углом , под действием сдвигающей нагрузки, распределенной по закону %(ф) = % sin(2ф), рис. 2. Предположим, что осевые перемещения оболочки ограничены.
Рис. 2. Схема нагружения
Граничные условия при 3 = 30: и к ) - и к СО8(&0) = 0; Б* = % при к = 2, Б* = 0 при к * 2; ЫХк = 0 ; ^ = 0.
В данных условиях реализуется мало изменяющееся практически безмоментное напряженное состояние.
На рис. 3 представлены распределения усилий в оболочке при различных углах подъема . Видно, что при = 14° результаты эталонного решения («ЭР») и решения
математической модели с учетом пологости («ПР») практически совпадают. При = 44° наблюдаются отличия.
Выполним серию расчетов при различных значениях и сопоставим значения усилий, полученных на основе математических моделей с учетом ( N ) и без учета ( И*) пологости оболочки. Для этого вычислим величину несоответствия решений
8М =
1 - И/И 1 -100%, I = 1,2 . (25)
На рис. 4 графически представлены максимальные значения взаимного несоответствия результатов решений математических моделей с учетом и без учета пологости.
Верхняя шкала на рис. 4 соответствует отношению величины стрелы подъема Ь к диаметру 2г0 основания оболочки.
2,0
N2
о
2,0
-2,0
N2
------
Л7!
5 10
.......«ПР» -«ЭР»
15 30
.......«ПР» -«ЭР»
(а) (б)
Рис. 3. Распределение усилий вдоль нулевого меридиана при: (а) = 14°, (б) = 44°
15 25 35 45
Рис. 4. Максимальные значения несоответствия решений 8Ni
Согласно определению пологости по Власову [1] пологой может считаться оболочка с соотношением Ц2г0 < 0,2 или < 43,6° .
Поскольку используемая математическая модель без учета пологости оболочки имеет погрешность допущений Кирхгофа, представленная на рис. 4 кривая фактически отражает значения погрешности, привносимой дополнительными гипотезами математической модели для пологой сферической оболочки.
Таким образом, из рис. 4 следует, что критерию пологости соответствует погрешность математической модели для пологой сферической оболочки в пределах 35%.
В условиях напряженного состояния с большой изменяемостью ожидаемо снижение указанной погрешности. Покажем это на следующей задаче.
Рассмотрим задачу деформирования конструкции в виде двух сферических оболочек (чечевицы) с параметром Щк = 100 и пологостью, характеризуемой углом , под действием самоуравновешенных сил Р , распределенных по дугам окружности, рис. 5.
Разложение внешней нагрузки в ряд Фурье: д(ф) = ^qk еов(кф), где ^ = Р/ш0,
к=0,2,4,...
qk = ^вт^/ку.
Граничные условия: и к §1п(30 ) ~изк ео8(30) = 0, Б* = 0, ^к= 0,
^ соб^) + би эт^) = qk^2 .
Выполним серию расчетов при различных значениях и у .
Значения параметров напряженно-деформированного состояния будем определять с точностью до 3 значащей цифры суммированием старших ( к > 0) гармоник рядов Фурье искомых величин.
Результаты расчетов представлены на рис. 6-8 в виде распределений меридиональных стт = (^/к + 6М^к2) и окружных ст, = (N2/И + 6М2/к2) напряжений вдоль нулевого ( ф = 0 ) меридиана.
Как и в предыдущем случае, при = 14° результаты эталонного решения («ЭР») и решения математической модели с учетом пологости («ПР») практически совпадают. Для оболочки на границе пологости в целом наблюдаются несоответствия в решениях. Однако вблизи места приложения нагрузки решения сближаются тем сильнее, чем выше степень локализации внешней нагрузки.
аН/д0 х10 2
От/ /
У
О 5 10 15 .......«ПР» -«ЭР»
Рис. 6. Распределение напряжений вдоль нулевого меридиана оболочки (= 14°, у = 45° )
ай/^ц х10
0,4
0,0
,9/
34 39 .......«ПР» -«ЭР»
Рис. 7. Распределение напряжений вдоль нулевого меридиана оболочки (= 43,6° , у = 45° )
Рис. 8. Распределение напряжений вдоль нулевого меридиана оболочки (= 43,6° , у = 5°)
С точки зрения выполнения условий прочности оболочки, наиболее значима величина погрешности результата, получаемого в месте с максимальным уровнем напряжений. В данном случае это точка с координатами ,0). Вычислим для этой точки, аналогично (25), величину несоответствия значений напряжений, рис. 9.
Рис. 9. Максимальные значения несоответствия напряжений в зоне их концентрации
Из рис. 9 следует, что в условиях концентрации напряжений погрешность математической модели Власова для пологой сферической оболочки снижается, в данном случае до 3-10% (в зависимости от степени локализации внешней нагрузки).
Аналогичные расчеты, выполненные при Rh = 400, 700, 1000, не выявили зависимость величины несоответствий эталонного решения и решения математической модели с учетом пологости от значения параметра R/h .
Заключение
Выполненный вычислительный сопоставительный анализ показал, что критерию Власова соответствует погрешность не более 35%. С увеличением изменяемости напряженного состояния погрешность математической модели Власова для пологой сферической оболочки снижается. Показано, что в случае локального нагружения оболочки в месте максимальной концентрации напряжений погрешность решения не превышает 3-10% (в зависимости от степени локализации внешней нагрузки).
Варьирование относительной толщины оболочки Rjh в диапазоне от 100 до 1000 влияния на погрешность математической модели Власова не выявило.
Список литературы
1. Власов В.З. Избранные труды. В 3 т. Т. 1. Очерк научной деятельности «Общая теория оболочек». Статьи. М.: АН СССР, 1962. 528 с.
2. Меньков Г.Б. Решение задач механики деформирования оболочек методом функционального нормирования: дис. ... канд. физ.-мат. наук. Казань, 1999. 197 с.
3. Виноградов Ю.И., Меньков Г.Б. Метод функционального нормирования для краевых задач теории оболочек. М.: Эдиториал УРСС, 2001. 160 с.
4. Григоренко Я.М., Ильин Л.А., Коваленко А.Д. Теория тонких конических оболочек и ее приложение в машиностроении. Киев: АН УССР, 1963. 287 с.
5. Филин А.П. Элементы теории оболочек. Л.: Стройиздат. Ленинградское отделение, 1975. 256 с.
6. Ватсон Г.Н. Теория бесселевых функций: пер. с англ. М.: Иностранная литература, 1949. 799 с. [Watson G.N. A treatise on the theory of Bessel functions. 2nd ed. Cambridge University Press, 1945. 804 p.]
7. Люк Ю. Специальные математические функции и их аппроксимации: пер. с англ. М.: Мир, 1980. 608 с. [Luke Y.L. Mathematical functions and their approximations. Academic Press, 1975. 568 p.]
8. Polyanin A.D., Zaitsev V.F. Handbook of exact solutions for ordinary differential equations. 2nd ed. CRC Press, 2002. 816 p.
9. Маричев О.И. Метод вычисления интегралов от специальных функций. Минск: Наука и техника, 1978. 310 с.
10. Михайловский Е.И., Новожилов В.В., Черных К.Ф. Линейная теория тонких оболочек. Л.: Политехника, 1991. 656 с.
11. Новожилов В.В. Теория тонких оболочек. СПб.: Санкт-Петербургский ун-т, 2010. 378 с.
12. Новожилов В.В., Черных К.Ф. К расчету оболочек на сосредоточенные воздействия // Исследования по упругости и пластичности: сб. Вып. 2. Л.: ЛГУ, 1963. С. 48-58.
Science and Education of the Bauman MSTU, 2014, no. 12, pp. 858-877.
DOI: 10.7463/1214.0738649
Received:
10.11.2014
Science ^Education
of the Bauman MSTU
ISSN 1994-0448 © Bauman Moscow State Technical Unversity
Vlasov Mathematical Model Inaccuracy Quantitative Assessment for a Shallow Spherical Shell
M. V. Konstantinov1'* 'konstigsmatejiiskju
1Bauman Moscow State Technical University, Moscow, Russia
Keywords: spherical shell, Vlasov mathematical model for deformation mechanics, flatness criterion, inaccuracy
A spherical is used shell as a perfect vessel in a weight ratio in the spacecraft, where frames assemble thin-walled elements and shells. Obviously, it is impossible to avoid local action on the shell and, consequently, presence of stress concentration. Attempt to reach perfect weight characteristics of spacecraft leads to lower safety factor of its parts. Reducing safety margin is possible provided that controlled inaccuracy of element stress-strain state assessment is guaranteed. Known universal numerical finite element method and finite differences method cannot cope with this task. So, analytical methods should be used to solve the problem. One solution is based on Vlasov mathematical model for a shallow shell.
Difficulties in application of Vlasov model are caused by insufficient flatness criterion validity in terms of quantifying its effect on inaccuracy of solution. This disadvantage is eliminated in this paper as a result of a computational comparative analysis of the Vlasov model solution with that of the mathematical model of shell deformation, linear differential equations of which are obtained with an error not exceeding the error of shell theory Kirchhoff assumptions.
For the first time is defined inaccuracy value due to accepted Vlasov theory simplifications depending on the nature of the loading and the geometric shell parameters. It is shown that inaccuracy not more than 35% corresponds to Vlasov flatness criterion. Increasing stress state variability leads to reducing inaccuracy. In the case of local loading the inaccuracy is 3-10% in the area of maximum stress concentration (depending on degree of loading localization). Variation of R/h from 100 to 1000 does not reveal any influence on the inaccuracy of Vlasov model.
The results of work performed allow using Vlasov mathematical model to determine the stress-strain state of thin-walled elements with a prior given error, for example, in research and development activities of aerospace systems.
References
1. Vlasov V.Z. Izbrannye trudy. V 3 t. T. 1. Ocherk nauchnoy deyatel'nosti "Obshchaya teoriya obolochek". Stat'i [Selected works. In 3 vols. Vol. 1. Essay of scientific activity "General theory of shells". Articles]. Moscow, USSR Academy of Sciences Publ., 1962. 528 p. (in Russian).
2. Men'kov G.B. Reshenie zadach mekhaniki deformirovaniya obolochek metodom funktsional'nogo normirovaniya. Kand. dis. [Solving tasks of shells deformation mechanics by functional normalization method. Dr. dis.]. Kazan', 1999. 197 p. (in Russian).
3. Vinogradov Yu.I., Men'kov G.B. Metod funktsional'nogo normirovaniya dlya kraevykh zadach teorii obolochek [Functional normalization method for shells theory boundary value problems]. Moscow, Editorial URSS Publ., 2001. 160 p. (in Russian).
4. Grigorenko Ya.M., Il'in L.A., Kovalenko A.D. Teoriya tonkikh konicheskikh obolochek i ee prilozhenie v mashinostroenii [The theory of thin conical shells and its application for engineering]. Kiev, Ukrainian SSR Academy of Sciences Publ., 1963. 287 p. (in Russian).
5. Filin A.P. Elementy teorii obolochek [Elements of the shells theory]. Leningrad, Stroyizdat Publ., 1975. 256 p. (in Russian).
6. Watson G.N. A treatise on the theory of Besselfunctions. 2nd ed. Cambridge University Press, 1945. 804 p. (Russ. ed.: Watson G.N. Teoriya besselevykh funktsiy. Moscow, Inostrannaya literatura Publ., 1949. 799 p.).
7. Luke Y.L. Mathematical functions and their approximations. Academic Press, 1975. 568 p. (Russ. ed.: Luke Y.L. Spetsial'nye matematicheskie funktsii i ikh approksimatsii. Moscow, Mir Publ., 1980. 608 p.).
8. Polyanin A.D., Zaitsev V.F. Handbook of exact solutions for ordinary differential equations. 2nd ed. CRC Press, 2002. 816 p.
9. Marichev O.I. Metod vychisleniya integralov ot spetsial'nykh funktsiy [The method of calculating the integrals of special functions]. Minsk, Nauka i tekhnika Publ., 1978. 310 p. (in Russian).
10. Mikhaylovskiy E.I., Novozhilov V.V., Chernykh K.F. Lineynaya teoriya tonkikh obolochek [The linear theory of thin shells]. Leningrad, Politekhnika Publ., 1991. 656 p. (in Russian).
11. Novozhilov V.V. Teoriya tonkikh obolochek [The theory of thin shells]. St. Petersburg, SPbSU Publ., 2010. 378 p. (in Russian).
12. Novozhilov V.V., Chernykh K.F. The issue of shells calculation in case of concentrated loading. Issledovaniya po uprugosti i plastichnosti: sb. Vyp. 2 [Investigations on elasticity and plasticity: collected papers. Is. 2]. Leningrad, LSU Publ., 1963, pp. 48-58. (in Russian).