УДК 531.383:532.516: 681.03.06
И.А. Ковалева
МОДЕЛИРОВАНИЕ ДИНАМИКИ НЕЛИНЕЙНЫХ ВОЛН В СООСНЫХ ФИЗИЧЕСКИ НЕЛИНЕЙНЫХ ОБОЛОЧКАХ, СОДЕРЖАЩИХ ВЯЗКУЮ НЕСЖИМАЕМУЮ ЖИДКОСТЬ МЕЖДУ НИМИ
Исследование посвящено анализу распространения нелинейных волн деформаций в физически нелинейных соосных упругих цилиндрических оболочках,
содержащих вязкую несжимаемую жидкость между ними. Волновые процессы в упругой цилиндрической оболочке без взаимодействия с жидкостью ранее исследованы с позиций теории солитонов. Наличие жидкости потребовало разработки новой математической модели и компьютерного моделирования процессов, происходящих в рассматриваемой системе.
Цилиндрические соосные оболочки, колебания, волны деформации, гидроупругость, вязкая несжимаемая жидкость, солитон, базис Грёбнера
I.A. Kovaleva NONLINEAR WAVES DYNAMICS MODELLING IN COAXIAL PHYSICALLY NONLINEAR SHELL CONTAINING VISCOUS INCOMPRESSIBLE FLUID BETWEEN
The present investigation is devoted to the analyses of non-linear deformation waves propagation in physically non-linear coaxial elastic cylinder covers containing viscous incompressible liquid in between. The wave processes in elastic cylinder cover without interaction with liquid were investigated earlier on the basis of soliton theory.
The presence of liquid demanded working out a new mathematical model and computer modelling of the processes taking place in the system.
Cylinder coaxial shells, oscillations, deformation waves, hydroelasticity, viscous incompressible liquid, solitary wave, Groebner basis
1. Введение
Приведение систем алгебраических, дифференциальных и разностных уравнений к канонической форме, называемой базисом Грёбнера [1], представляет собой качественный аналитический метод исследования соответствующих математических моделей.
В частности, при поиске частных решений дифференциальных уравнений методом неопределённых коэффициентов возникают переопределённые системы алгебраических уравнений. Построение базиса Грёбнера позволяет проверить совместность системы, определить, обладает ли система конечным или бесконечным числом решений, а в ряде случаев построить решения в явном виде.
Не для всех моделей, описываемых уравнениями в частных производных, удаётся построить аналитические решения и в этом случае для их исследования можно применять численные эксперименты на соответствующих разностных схемах. Так, для построения разностных схем из первоначально заданных базовых разностных соотношений, аппроксимирующих исходную систему дифференциальных уравнений, строится базис Грёбнера разностного идеала. Из этого базиса иногда в нелинейном и всегда в линейном случае можно извлечь разностную схему, которую невозможно построить традиционными методами генерации разностных схем. Зачастую такие разностные схемы обладают уникальными свойствами, хорошо передающими физику процессов, описываемых исходными дифференциальными уравнениями [2].
Кроме того, знание базиса Грёбнера даёт возможность проверить совместность исходных разностных соотношений, определить произвол в решении, посчитав полином Гильберта и, применяя специальный вид допустимого упорядочения при его построении, получить другое представление первоначальных разностных соотношений.
В представленной работе данная техника будет использована в качестве примера для анализа распространения нелинейных волн деформаций в упругих физически нелинейных соосных цилиндрических оболочках, содержащих вязкую несжимаемую жидкость между ними.
2. Постановка задачи гидроупругости и построение математической модели методом возмущений. Волновые процессы в упругой цилиндрической оболочке без взаимодействия с жидкостью ранее исследованы в [3, 4] с позиций теории солитонов. Рассмотрим бесконечно длинные соосные упругие цилиндрические оболочки, между которыми находится вязкая несжимаемая жидкость. Уравнения движения вязкой несжимаемой жидкости и уравнения неразрывности в цилиндрической системе координат r, З, x записываются в случае осесимметричного течения в виде [5]
= V
дУ. + +у эvL +!_ др
дt г дг х дх р дг
дУх .. дУх дУх 1 др
х-+Уг—х+Ух^т^+—^~ =v
(■л2
д2У 1 дУ д2У У
дt
дУ
дг
дх р дх
2
дУ
у дг
^ 2
- +------------------г--------------------г.
г дг дх 2 г
1 дУ дУ л
+ —
V
дг г дг дх
(1)
г Уг дУх
г + —+ - х
= 0.
дг г дх
На границе с оболочкой выполняются условия прилипания жидкости
У
дW
(0
дt
-,Ух
ди
(‘)
дt
при г = Я{ — W
(0
(2)
Здесь t - время; Уг,Ух - проекции вектора скорости жидкости на оси цилиндрической системы координат; р - давление; р - плотность; V - кинематический коэффициент вязкости; и(,) -
продольное упругое перемещение оболочек по оси х; W(,) - прогиб, положительный к центру кривизны оболочки; Я1 - внутренний радиус внешней оболочки; Я2 - внешний радиус внутренней оболочки (Я1 = Я2 + 3); 3 - толщина слоя жидкости в кольцевом сечении трубы; г = 1 относится к
внешней, а г = 2 - к внутренней оболочке.
Записывая уравнения движения элемента цилиндрических оболочек в перемещениях для модели Кирхгофа-Лява, считаем материал нелинейно-упругим с кубической зависимостью интенсивности напряжений <Г1 от интенсивности деформаций е1 [6,7]
01 = Ее1 — ше1. (3)
Здесь Е - модуль Юнга; ш - константа материала, определяемая из опытов на растяжение или сжатие.
Уравнения динамики физически нелинейных оболочек с учетом (3) записываются в виде [7]
Е(г) д ,,ди (г) п) W (,\Г1 4 ш(°г ди (г\2 ди (г) W (г) № (°
х1+3 ^[(— )3 -
<(-
1 — А дх дх
д 2и(г)
р) иО) д и __Л) р0 'Ч
дt2
Чх •
о ) 0
Е(1) к
3
4Т*/(0
я
Е({)кр 1
дх
дх Н'~ + ^)2]1> —
(
0
12(1 — а0°2) дх4
1 — У02 Я(0
А
ди(г)
(I) ди 0
W
(г) Л
ди(г) W(0
(г)
дх Я
Я
1Г+(ЪгУ]})+р0г)к0г)
\\
д 2W(0
дt
2
Здесь р0г) - плотность материала оболочки; А) - коэффициент Пуассона; Я(п - радиусы
дх
(г)
Я
(г)
{1 —
4 ш{1^ ди (г'\2
3 Е
(г) «[(:
дх
)2 —
(4)
(— 1)м
?(0
срединной поверхности оболочек; к0г) - толщины оболочек (к|(1)/2 = Я(1) — Я1,к0(2)/2 = Я2 — Я(2));
(г)
сУ -
скорость звука в материале оболочки; Чх , Чп - напряжения со стороны жидкости, находящейся внутри кольцевого сечения.
Если снести напряжения на невозмущенную поверхность оболочек (!) << Я{), то можно
считать, что поверхностные напряжения со стороны жидкости определяются формулами
чхх) =
, дУх дУг
рv| —х + —г 1 дг дх
Чп =
г = Я;
р + 2рv
ЭУ,
дг
(5)
г = Я;
Принимая за характерную длину волны I и считая, что соосные оболочки изготовлены из одного материала, то есть опуская индекс г у Е,ш,р0, А0, ^, перейдем к безразмерным переменным для исследования уравнений (4)
w(! ) = ),
и
(г)
: ищи
(г)
г * = с° г, I
* _ х
х = 7 ’
Е
Р0
(1 — А02)
(6)
п
с
0
Положим
= є = 0(1),
І
Н
(і)
Я
(і)
= О (є),
Ж
т
Я (І)
Е_
т
= О(є), = О(є),
Я
(і)
= О(є1/2),
(7)
где £ << 1 - малый параметр в задаче (4).
Применим метод двухмасштабных асимптотических разложений, вводя независимые переменные в виде
с* & & &
£ = х — сг , г = я , (8)
где с - безразмерная неизвестная скорость волны, 1 - внутренняя переменная, а зависимые
переменные представлены в виде разложения по малому параметру £ :
= и^ + £иЦ +___, и3) = и 30 + £и31 +____ (9)
Подставляя (6), (8), (9) в уравнения (4) с учетом оценок (7), получим в нулевом приближении
по £ линейную систему уравнений, из которой следует связь
итЯ
и определяется безразмерная скорость волны
(і)
(і) и30
А
диш
д#
1 — А)-
(10)
(11)
Из следующего приближения по Є, учитывая (10) и (11), находится система уравнений,
являющихся составными для и
(і) . 10 :
Э 2и« 1 Г Я(і) ^
- +-
дфт є
2т ( ит Еє { І
А
V1 — Ао д4и:
4,/і ) 10
д#4
(1
А + А
У1
А
Г ди !•') ^
д
2„«
V
д# і д#2
(12)
І2
ЄитР0Н0 —0 2д/1 — Аі2
і
> (і)
(і) Я(і) д^
^ — А0----------------------
І д#
(—1)
і—1
В случае отсутствия жидкости правая часть уравнений (12) равна нулю и система распадается на два одинаковых уравнения, представляющих модифицированные уравнения Кортевега - де Вриза для
ди
(і)
10
ЖтІ
и
(і)
(і) 30
д£ АИш^
Для определения правой части уравнения (12) ведем безразмерные переменные и параметры
V
І
ж
Я1
г — Я2 * —0 * х
--------г = — г, х = -, р
8 І І
РУ—0ІЖт
83
Р;
¥
8 = 0(1), 1 = Жт = жтъ...
Я
8
8 = _8Яі= Я _ 1
І ~ Я2 І ~¥ І ~ 1
Я2 8
О
(13)
у¥ у
іу<< 1, Л<< 1.
Подставляя (13) в уравнения гидродинамики (1) и граничные условия (2), представим безразмерные скорость и давление в виде разложения по малому параметру Л :
V0 + IV1 + ..., vr = V0 + Лv1 + ..., Р = Р° + ЛР1 + к.
(14)
В нулевом приближении по (31 ~ 0 - гидравлическая теория смазки), считая
(31 )(8с^')<< 1 (ползущие течения), и в нулевом приближении по Л получаем уравнения
гидродинамики (классические уравнения гидродинамической теории смазки)
дР0
дг *
= 0,
дР0
дv0 + сМ
дх дг*2 дг * дх
=0
(15)
І
2
—
2
І
2
2
1
V
х
и граничные условия
ди<31)
0
vr =
дґ * ди®
дґ *
vx =0 при г = 1, V0 = 0 при г* = 0.
(16)
Из решения задачи (15), (16) следует, что
Р0
= 12/ Л'^
ди
дґ дґ
йх
йх*,
Г дVx0 ] =6/ г ди3(2) ди31) Л
[ дг" ] г*=1 { дґ * дґ * I
дvx0 =6/ ^ ди3(1) ди32) Л
дг * г*=0 ' дґ * дґ * )
йх\
йх*.
(17)
Учитывая, что были введены переменные (8), (9), и имея соотношения (10) и (11), из (17) получим
,(1) (1) п(2) (2)
и10 - Я и10 !
С принятой точностью по Є,у,А. из (5) найдем
(І)
чх)
и, следовательно, << дп и в правой части уравнения (12) остается выражение
ру
™тС0
81 дг*
ру
*=1,0
™тС0 1Р 82 8 ’
£-=СІ8'
Чп
і
6А
рі У Я(і)Я(2) ди® Я(і)Я(1) ди® ,
р0к0° 8с0є 82 д£ 82 д£ (
(-1)
і—1
С принятой точностью по /, £ положим
Я(1) » Я(2) = я, н(1]
■к
(2)
00
Подставляя (19) в уравнение (12), окончательно получим
V
д 2"(0
и
10 +1Г Я12 А
л/ї-А)" д 4и(і)
и
10
д^дт є{ і
д£4
2т Г и.
(1 -А +А02 У1-
Г ди(і) Л д2 и(і)
- ЕЄ {?1 (1 -Д> +А>
а
10
д2и
10
д$ I д£2
2Г (1)
6АІ
ди^ ди10
рі у ' Я
р0к0 8:0є{ 81{ д% д%
{ ^ у
(2) ^
(-1У =0.
у
Легко видеть, что замена
ди (1) ди (2)
2!^- сф{1\ ^ = ер(2\
позволяет записать систему уравнений (20) в виде
(р^1)г +ф1\щ — бф~1) \ — (^(1) — ^(2))(-1)( = 0.
Постоянные с, с1, с2 определяются при подстановке (21) в (20) и имеют вид
С2 = 6 А
рі У
р0к0 8с0Є { 8\
Я1, с
С2Є —
2
і
1/3
3
А
2т и1 1 -А +А02
1/2
(18)
(19)
(20)
(21)
3
Ч
п
г
2
2
2
с
с
В случае отсутствия жидкости последние два слагаемых (р(1) - р(2)) в уравнениях (22) исчезают и система распадается на два независимых уравнения МКдВ (модифицированные уравнения Кортевега - де Вриза) и имеет точное частное решение в виде кинк - антикинк
р = ±к ґапк(кц + 2к ґ).
(23)
Эти решения при г = 0 можно взять в качестве начальных условий при решении задачи Коши для системы уравнений (22).
3. Компьтерное моделирование. Метод конечных объёмов сводится к дискретизации исходных уравнений, представленных в интегральной форме, в противоположность методу конечных разностей, который обычно применяется к исходным уравнениям в их дифференциальной форме. При этом если исходная система обладала законами сохранения, построенная разностная схема будет обладать хорошими консервативными свойствами просто по построению. Кроме того, при этом подходе упрощается вывод разностных соотношений на границах вычислительной области.
Если исходные уравнения содержат производные выше первого порядка, то метод конечных объёмов нуждается в модификации. Это модификация получила название интегро-интерполяционного метода, недостатками которого являются отход от работы только с интегральными соотношениями и прямая замена производных их конечными разностями. Если на этом этапе добавить интегральные соотношения, связывающие искомые функции с их производными, а затем использовать алгоритм Бухбергера построения базисов Грёбнера или инволютивный алгоритм, то можно получить соотношения, связывающие только искомые функции [2, 8].
Запишем уравнение (22) в интегральной форме
/д^(2р°)3 -р)пп)йґ + р°)йц-/£(р(1) -р(2))(-1)‘йґйп
0
(24)
для любой области ^ . Для перехода к дискретной формулировке сопоставим и(!). = р(1 \гп ) и
выберем в качестве базового контур, показанный на рис. 1.
Рис. 1. Базовый контур для уравнения (24) Добавим интегральные соотношения
и = и(1 \г ,п;.+1) — и(1 \г ),
п+1и
1п1
п+1и{ї)тйц = и(і)п(ґ,П;.+1) - и(!)п(ґ,П;').
(25)
Используя для интегрирования по времени и по четным производным по П формулу трапеций, а по нечетным производным по 1) - формулу среднего значения и полагая гп+1 — гп = Т, П^+1 —Ц. = к , перепишем соотношения (24), (25) в виде
3П+1
2и(!)Л] + 2и(г)Л] - 2и(г)Л]+2 - 2и(г)3]+2
- (и(г
(і) п (і) п+1 (і) п
(і) п+1
и‘щ^ и’пп] - и’пп]+2 - и’пп]+2 ))' 2 +
+ (и(і)п+1 - и(і)п+1) • 2к - (и(і)п+1 + и(і)п+1)(-1)г • кт = 0,
(и (і) п + и (і) п ) к = и (і)п, и (і)п
(и п }+1 + и п } ) • 2 = и 1+1 - и 1,
(і) п
(і) п
(і) п
и 'пп 1+1 • 2к = и 'п 1+2 - и п 1 •
г
Вводя сеточные операторы сдвига 6{ ,6^ по переменным соотвественно, запишем
уравнения в операторной форме:
2 Т
— (1 + 6 — 62 —6г6п') ° (—3и (!) + и (!)ПП) ' ^ +
+ (6п6г — 6п) ° и(0 ' 2Ь — (6п6г + 6п)° и(0(—1)!'' Нт = 0,
Ь
(6ц +1) ° и(!)п ' 2 = (6п — 1) ° и(0,
6П °и' 2Ь = 6 — 1) °и%.
Выбирая допустимое лексикографическое упорядочение сначала по функциям
и (1)пп ^ и (2)пп ^ и % ^ и (2)п ^ и(1) ^ и(2), затем по переменным 6( ^ 6П , можно построить базис
Грёбнера или инволютивный базис [8]. В результате получим в качестве отдельных элементов авторедуцированного базиса Грёбнера следующие разностные схемы для уравнений (22), аналогичные схеме Кранка-Николсона для уравнения теплопроводности
и+1 п зП+1 зП+1 зП зП
и(г) и(г)п (и(г )■ 1 и(г )■ 1) + (и(г)■ 1 и() • 1)
и 1 — и 1 — 2(и . +1 — и .—1) +(и .+1 — и .—1) +
Т 4h
(и(!) П+2 - 2и(!
+
( (i)n+1 О (i)П+1 . О (i)"+1 (i)^ \ / (i)n О (i)n , О (i)n (i)n \
(u j+2 — 2u j+1 + 2u j-1 — и j-2) + (u j+2 — 2u j+1 + 2u j-1 — и j-2)
4h3
(i)n+1 (i)
- u j + u j (-1)i =0.
2
Полученные неявные разностные схемы имеют кубическую нелинейность для следующего временного слоя. При построении решения использована следующая линеаризация:
vl+1 = vl+1 - vl + vl = (vk+1 - vk )(v2+1 + vk+1vk + v2) + vl ~ vk+1 • 3v2 - 2v3.
Шаг по времени t брался равным половине шага по переменной П. Программа расчёта была написана на языке Python с использованием пакета SciPy [9] и занимает около 150 строк, включая сам расчёт и построение графиков.
Результаты проведённого компьютерного моделирования представлены на рис. 2, 3 и позволяют сделать следующие выводы. В первом случае отсутствие жидкости между оболочками приводит к возникновению волны деформации (кинк) только во внешней оболочке, а во внутренней оболочке, как и в начальный момент, деформация равняется нулю (см. рис. 2).
Наличие жидкости между оболочками приводит к возникновению волны деформации и во внутренней оболочке, в которой в начальный момент деформации равнялись нулю (см. рис. 3).
Это процесс происходит за счёт «перекачки» энергии (через слой жидкости) от волны (кинк) во внешней оболочке и сопровождается падением амплитуды волны во внешней оболочке и, как следствие, снижением скорости её распространения.
В результате во внешней и внутренней оболочках устанавливается волна деформации постоянной амплитуды и скорости распространения.
Проведенное моделирование позволяет сделать вывод, что рассматриваемая механическая система начинает вести себя как единый трёхслойный пакет с двумя несущими слоями (внешняя и внутренняя оболочки), по которым распространяются волны деформации, и заполнителя - слоя вязкой несжимаемой жидкости.
4. Заключение
Проведенное моделирование с использованием компьютерной алгебры позволило выявить особенности поведения волн деформаций в физически нелинейных упругих соосных цилиндрических оболочках, содержащих вязкую несжимаемую жидкость между ними.
Использование базиса Грёбнера для генерации разностной схемы при численном решении задачи Коши для системы двух нелинейных уравнений в частных производных третьего порядка по пространственной переменной позволило получить результат расчета без осцилляций, вызываемых численной реализацией. Численная схема была протестирована на точном решении при отсутствии жидкости (см. рис. 2).
n
я
Рис. 2. График численного решения уравнений (22) при отсутствии жидкости и с начальным условием (23) для р(1) с к = 0.2, и для р(2) = 0.0
— *=0.00 - - * =20.04 ' - - * =40.08 Ь =60.12 Ь =80.16 -- « = 100.20
/ / ; ; / / / : / /
/ ' ' '
//• / '/
“э. 00
V
Рис. 3. График численного решения уравнений (22) при наличии жидкости и с начальным условием (23) для р(1) с к = 0.2, и для р(2) = 0.0
Полученный расчет показал влияние вязкой несжимаемой жидкости на поведение нелинейной волны деформации в соосных оболочках (см. рис. 3). В результате возникает нелинейная волна деформации во внутренней оболочке, в которой ее не было в начальный момент времени и амплитуды волн деформации в соосных оболочках со времением начинают совпадать. Эти амплитуды в два раза меньше исходной аплитуды волны деформации внешней оболочки в начальный момент времени.
Эту конструкцию можно толковать как трехслойный пакет, заполнителем которого является жидкость.
В заключение приношу благодарность профессору Л. И. Могилевичу за постановку задачи и внимание к работе.
Работа выполнена при поддержке гранта РФФИ № 13-01-00049а.
ЛИТЕРАТУРА
1. Buchberger B. Groebner bases and applications / B. Buchberger, F. Winkler // London Mathematical Society Lecture Notes Series. Cambridge University Press, 1998. Vol. 251. P. 552.
2. Gerdt V.P. Groebner bases and generation of difference schemes for partial differential equations / V.P. Gerdt, Yu.A. Blinkov, V.V. Mozzhilkin // Symmetry, Integrability and Geometry: Methods and Applications. 2006. Vol. 2. P. 26. http://www.emis.de/journals/SIGMA/2006/Paper051/index.html.
3. Землянухин А.И. Нелинейные волны деформаций в цилиндрических оболочках / А.И. Землянухин, Л.И. Могилевич // Изв. вузов. Прикладная нелинейная динамика. 1995. Т. 3. N° 1. С. 52-58.
4. Землянухин А.И. Нелинейные волны в цилиндрических оболочках: солитоны, симметрии, эволюция / А.И. Землянухин, Л.И. Могилевич. Саратов: СГТУ, 1999. С. 132.
5. Лойцянский Л.Г. Механика жидкости и газа / Л.Г. Лойцянский. М.: Дрофа, 2003. С. 840.
6. Каузерер К. Нелинейная механика / К. Каузерер. М.: Иностр. лит., 1961. С. 240.
7. Вольмир А.С. Оболочки в потоке жидкости и газа: задачи гидроупругости / А.С. Вольмир. М.: Наука, 1979. С. 320.
8. Блинков Ю.А. Генерация разностных схем для уравнения Бюргерса построением базисов Грёбнера / Ю.А. Блинков, В.В. Мозжилкин // Программирование. 2006. Т. 32. № 2. С. 71-74.
9. SciPy. http://www.scipy.org/
Ковалева Ирина Александровна - Irina A. Kovaleva -
аспирант кафедры «Теплогазоснабжение, Postgraduate
вентиляция, водообеспечение и прикладная Department of heat and gas supply,
гидрогазодинамика» Саратовского ventilation, water supply
государственного технического and applied fluid dynamics,
университета имени Гагарина Ю.А. Gagarin Saratov State Techmcal University
Статья поступила в редакцию 12.07.12, принята к опубликованию 06.11.12