CC BY
27УДК 532.59:629.12
С.А. Химич, Ю.Ф. Орлов
ДВИЖЕНИЕ НАД ГРАНИЦЕЙ РАЗДЕЛА ДВУХ СРЕД КРЫЛА КОНЕЧНОГО РАЗМАХА
Нижегородский государственный технический университет им. Р.Е. Алексеева
Рассмотрена задача о генерации поверхностных волн при движении над поверхностью раздела жидкостей с различной плотностью на малой высоте крыла конечного размаха в квадрупольном приближении. При этом двумерное сингулярное интегральное уравнение теории крыла вырождается в краевую задачу для уравнения Пуассона для определения плотности двойного слоя, распределённого по поверхности крыла. Эта задача решалась численно. Рассмотрен пример движения плоского крыла, имеющего размеры реального экраноплана.
Ключевые слова: граница раздела, крыло экраноплана, квадрупольная теория крыла, волны под крылом.
Задачи о движении крыла экраноплана интенсивно исследовались в 60-х годах прошлого века [1]. Первые решения задач методом функциональных параметров получены А.Н. Панченковым [1, 2]. Позже им же создана квадрупольная теория крыла вблизи экрана [3], которая конструктивно использует малость отстояния крыла от экрана, позволяющую перейти от моделирования несущей поверхности слоем диполей к слою квадруполей, распределенных на экране под крылом. Это вырождение дает возможность перейти от сингулярного интегрального уравнения к дифференциальному соотношению в плоских потоках и к уравнению Пуассона для области течения под крылом в пространственном случае. Для задач со свободной поверхностью первые такие алгоритмы построены в [4, 5]. Обширный круг задач по рассматриваемой теме представлен также в монографии М.А. Басина и В.П. Шадрина [6].
Рассмотрим задачу о волнах при движении на малой высоте над поверхностью раздела жидкостей с различной плотностью крыла конечного размаха. Пусть малоискривленное крыло S движется с постоянной скоростью V0 на малой высоте h над границей раздела Г
двух жидкостей различных плотностей, Pj - плотность верхней жидкости, р2 - плотность нижней жидкости. Введем систему координат, связанную с крылом, S - проекция поверхности крыла S на координатную плоскость xy, S0 - проекция поверхности S на границу раздела жидкостей (рис. 1).
Рис. 1. Выбор системы координат
© Химич С.А., Орлов Ю.Ф., 2011.
Основное двумерное сингулярное интегральное уравнение линейной задачи о движении крыла вблизи границы раздела жидкостей имеет вид [1, 2, 3, 4]
4 У = ¥(д), (1)
z
>2 ,
где N0 =-f [] dx; А, у = -^ íy(p)d q) ds; p = е Sp ; q = (х, J, z) eSp; y(p) -
-L 4ж ¿ & p P
плотность двойного слоя на поверхности S ; F(q) - проекция на нормаль к S скорости набегающего потока жидкости; G (P, q) - фундаментальное решение уравнения Лапласа.
Введем параметр относительной плотности р = р: / р2 << 1, который предполагается малой величиной, функции потенциала ускорений (х, у, z) - для верхней жидкости, 02 (х, у, z) - для нижней жидкости. Линеаризованное граничное условие на поверхности раздела Г жидкостей имеет вид
р (01** -|0„ + V 0„)-02** +|Ах - V 02z = 0 , (2)
g 2
где V = —-, g = 9.81 м/c - ускорение свободного падения; | - коэффициент диссипативных
Vo
сил. По условию (2) фундаментальные структуры определяются в виде [3]:
Gl(p,q) =1 л —г+ G11(P,q)+ G12(P,q)' (3)
R(p,q) Rl(p,q)
G2 ^ q) = q) + G22 (p, q), (4)
где G (p, q), G2 (p, q) - потенциалы единичного источника, движущегося в верхней жидкости.
ж/2
G i(p, q) = Re iva • (1 - a) • í exp[--— (- (z + О + 2h) + /ш)]-sec2 0 d0,
¿9 cos20
-ж /2
G12(p,q) = -Re^ -V.p.fíexp[r(-(z + C + 2h) +/ш)]гcos20drd0, 2ж r cos 0- va
-ж 0 ж/2
G21 (p, q) = - Re iva • (1 - a) • f exp[--— (z-£-2h + /ш)] • sec2 0 d0,
J me 2 A
-ж /2 ^^ "
,2 ,
G22(p,q) = Re^-02-F.p. f f (Zf + rCOS 0drd0, (5)
2л JJ, r cos 0- ra
-л 0
R(p,q) = V(--^)2 + (У-Л)2 + (*-C)2 , R(p,q) = V(--^)2 + (У-Л)2 + (* + C + 2h)2 ,
- 1 -P
a =——; ю = (--^)cos 0 + (y ^)sin 0.
1+ P
Поставленная задача будет решена, если будет решено уравнение (1) с ядром (3), (5). В [2] впервые получено решение (1) методом функциональных параметров. В [4] близкое к (1) уравнение о движении тела по поверхности тяжелой жидкости после асимптотического анализа структур типа (5) сведено к двумерному уравнению сингулярному по поперечной координате и уравнению с периодическим ядром по продольной координате. Получено решение первого приближения по шкале va. Существенного упрощения расчетов можно добиться используя квадрупольные вырождения фундаментальных структур (3) [3].
Поскольку (3) - гармонические функции, которые вместе со всеми своими производными равномерно непрерывны всюду в трехмерном евклидовом пространстве, исключая точки (£, л, C) и (£, л, - C - 2"), то их разложения в ряды Тейлора относительно точек
р0 — - к), лежащих на границе раздела, - равномерно сходящиеся ряды. Выполнив необходимые вычисления, получим:
8 °1(р- ч) = 8^(Ро» ч)+(2+к) 8 2 °1(Ро' ч)+(1+И)28 ^^ ч) + (6)
8£ 8С ( ) 8£2 ( ) 8С3 • 2! "''
^(А» Ч) _ 8Оп(Ро, Ч) , 8^12(Ро» Ч) где — -г ,
8С 8С 8С
= О8 2(1/Г0) , 8 2°11(Ро» Ч) | 8 ^(Ро» Ч)
8С2 8С2 8С2 8С2 ' ( ^
&30(Р°2Ч) = о 8 3(1/ Го) | 8 ^Ро» Ч) | 8 3^12(Ро» Ч) 8С3 8С3 8С3 8С3 '
Го — ^(X-02+(у-л)2+(2-к)2 .
Будем использовать асимптотическую квадрупольную теорию крыла, которая основана на аппроксимации фундаментальной структуры двумя членами ряда (6) [3, 4, 5]. С учетом (7) получим:
8С
Это выражение можно переписать в виде
8G1(p1q) = 8Gn(Po, q) | 8Gi2 (Po , q) , 2(r , 7i)82(1/ ^o) + 8C 8C 8C ( ) 8C2
i (r i tq8 2[Gii^Po-q)+Gi2( Po> q)]
(8)
8 Gl(p'q) =82[f Gn(Po, q)dC + j G^fo, q)dC+ 2( z + h) / —o +
8C J J (9)
+ (z + h) (Gu(Po, q) + Gi2( Po, q))]/ 8C2.
Введем безразмерные величины соотношениями: x = xa, y = y b, z = z a, X = b / a,
— ^^ ^
^ = , ^ = ^b, С = Ca, г = —, G = ~, G2 = —, у = 2Xy, h = 2a h , ш0 = va = -
II ---13-" -----J — o /->7-^2'
a a a 2F—
где Fr = ,V - число Фруда, 2a - хорда, 2b - размах, X - относительное удлинение крыла.
V2ag
Выделив реальную часть в выражениях (5), получим (черточки над безразмерными координатами далее опущены):
г/2 ехР [ZWoa • (z + С + 4h)] __
п г \ — /1 —ч f cos 0 . Qo a ш Gii(P,q) = -Qoa •(1 -a)• I -^--sin(—T^)d0,
cos2 0 cos2 0
-л /2
tt , Л 1 - a гг — cos 0-exp[-— • (z + С + 4h )]
Gi2(P>q) = — V.P. I I-^-=--cos(—Q)d—d0.
2л JnJ0 — cos 0-шо a
г/2exp^-Q°a • (z -C- 4h )] __
n f \ - n f cos 0 Qo a ш
G2i(P,q) =Qoa •(1 -a)• I -^--sin(-TK)d0,
cos 2 0 cos 2 0
-л / 2
ч 1 - a rr — cos2 0^ exp[— • (z-C-4h)] _
G22(P,q) = —— V.P. II-^^—=—-•cos(—ш)d—d0,
2г JnJ0 — cos 0-шо a
- 1 -Р — где a =—з, ю = (х-^)cos0 + X(y-^)sin0 . 1+ Р
В силу гармоничности функций —, G i(p, q), Gi (Р, q) и использования (9) в (1), по-
r0
лучим краевую задачу:
д2
(^ + X-2 ^Г) ф(х, y) = ^^, (10)
д2 .2 д2 Л F(х,y)
,—T + X 2-5") Ф(х, У) = — /
дх2 д y2 4 H X2
где H = h/2b; х е[0,2]; y е[-1,1]; Ф(х,± 1) = Ф(2,у) = 0; дФ(0,у) = 0;
д х
Ф(х,y) = N,,Ly) + ± j ,).[G-<p.q)tG-2(p,q) + д(G„<p.q> +g.2(p,q))]ds(11)
I 0 h д^
Уравнение Пуассона (10) описывает движение жидкости в ограниченной области между крылом и поверхностью раздела жидкостей. Исходя из аппроксимации (9), предполагается, что возмущения за пределы этой области не распространяются. Поэтому в (11) опера-
2
тор N0 следует заменить на оператор N00 = j"[] dx, определяющий отсутствие возмущений
х
перед поверхностью S 0. Вместо (11) можно записать
у(х, у) = (1 - a )Гу--дФху), (12)
дх
где Lу = -(Л) Гу(£,Л)ds;
Sp
q) = (1- a) • <~и (р, q); G12 С?, q) =(1 - а) ■ <12 ^ q) •
Из соотношения (12) плотность двойного слоя у(х, у) может быть найдена и использована для нахождения потенциала скоростей нижней жидкости:
ф(х,у,г) = -(-1)|у(£, л)• ,q) ds • (13)
Sp
При вычислении формы свободной поверхности жидкости может быть использована известная формула теории волн малой амплитуды [7]:
С,(х,у,г) = й-2аф%М, (14)
дх
где С(х,У,г) = С(х,у,г)/а - относительная высота волны в долях полухорды крыла.
Итак, форма свободной поверхности жидкости под низколетящим крылом может быть найдена в результате последовательного решения краевой задачи (10), интегрального уравнения (12) и вычисления потенциала скоростей нижней жидкости и формы поверхности по формулам (13) и (14) соответственно.
Полученная последовательность задач (10), (12), (13), (14) решалась численно. Для решения краевой задачи (10) строится прямоугольная сетка размером п х т точек
= К,У1 : 0 = х0 < Х1 <... < Хп = 2; -1 = у0 < у <... < Ут = 1}. Шаги сетки по х и у зада-
2 2
ются формулами к = —, к = — соответственно. Во внутренних узлах этой сетк:и значения
п т
функции Ф( х, у) могут быть вычислены по итерационной формуле [8]:
h2 h2 k2 Fl , 9-' =-^-— ■ (15)
2k2 + — X2
значения же функции ф( x, y) в граничных узлах вычисляются по формулам [8]:
ф =ф =ф =o
o фг, m фn, j o ,
h2 h2 k2 Ft j
2Г- Ф,,+ (ф o.j „+Ф o.j-О - -4Hf.
фo,j =-—- ■ (16)
2k + 2T
Значения 8ф(x*y) во внутренних узлах сетки могут быть найдены по центрирован-
8x
ной формуле порядка O(h4) [8]:
ЯсТ) - ф „ + 8ф , - 8ф , + ф „ _ _
^ = фi+2j + 8фi+1j 8фi-1j +фi-2j , i = 2, n - 2, j = o, n . (17) 8x i, j 12h
Значения 8ф(x, y) в приграничных узлах сетки могут быть найдены по центрированной
8x
формуле порядка O(h2) [8]:
8ф ф 2, -ф o, —
^^ = —^-о- , j = o, n, (18)
8x 1, - 2h
ao =ф,,,-ф^, -=o;^. (19)
8x n-1, - 2h
Значения же 8ф(x, y) в граничных узлах сетки могут быть вычислены по формулам правой
8x
2
и левой разностей порядка O(h ) соответственно [8]:
8ф - 3ф0- + 4ф1 ,-ф 2- —
^ =-о--у-г-,- = o, n, (20)
8x o, - 2h
8ф 3ф n,j - 4ф n-1,j +ф n-2, j . Г"
- =--- , j = o, n. (21)
8x n, - 2h
TT " J, 8ф( x, y) o Дискретное множество значений функции -, заданных на сетке So , можно
8x h k
интерполировать с помощью бикубических сплайнов, положив g(x, y) = 8ф(x, y) , по форму-
8x
ле [9, c. 146-149]:
(x - x)3 (x - x )3
g(x, y) = gxx(xi-1, У) ' + gxx(xi, У) '" +
6h 6h
+ (g (x, _1, y) - g"( x-',y) h2) ^ + (g (x,, y) - b&M) ^ 6 h 6 h
где
(22)
(y - y)3 (y - y ,)3
gxx(x,-1, y) = Ki_1, -_1—^ + K_1, . 'y y--1) +
K i ,-ik2 y - y K i , k2 y - y.
+ (M, _1,j _1 - ) ^ + (M, _Lj - ^) ^;
* (х у) - К "У)3 + К (У~Уз-1)3 +
8Лхг,у) - Ки>-1 + К' 6* +
КJ-1 к\У]-У К'>1к\У-У]-
+(м,) —1—+(м, 1 —) —г~; 6 к 6 к
(У -У)3 (У-У- ,)3
19 У) - Ыг 1 . У) + N 1 , У^-1) +
5 4 г-Ъ?) г-1, у-1 бк '-1"/ 6к
+ (8-1, у-1--1--+ ( 8-1,1--)—т^-;
6 к 6 к
(У -у)3 (У - У ,)3
8(х., У) - Ы1Г У) + Ыг , У^-1) +
^ г ^ г, 1 -1 6к г, 1 6к
Ы, 1 -1к \ У7~У , , Ы, 1к \ У~У}-1
+ (8г, 7-1--1-^— + (8г, 7--^--т~;
6 к 6 к
N - 8уу (х,, У} ) > М - 8хх (х,, У} ) > К - ЯххУУ (хг, У; ) •
Так как оператор Ьу сжимающий [5], то, положив в качестве нулевого приближения -/ ч - д. дФ(х, у)
у(х, у) интерполированное множество значений функции--, для уравнения (12)
дх
можно использовать расчетную схему, основанную на вычислении последовательных приближений:
дФ( х, у)
У 0 -
дх дФ( х, у)
У и+1 - (1- а) ЬУ п- я
дх
, 1 чг2-/е ч Л,-2И,х,У,г) + G12(£,л,-2й,х,У,г) где Ьу - -(—) | | у(£, л) • [ '----+
-10
+ д(^11 Л, - , x, У, г) + ^12 ^ - , x, У, г)) ^
дС •
Было установлено, что для сочетаний жидкостей воздух-вода достаточно вычислить одно-два приближения.
После нахождения дискретного множества значений функции у(х, у), его необходимо интерполировать по формуле (22), положив 8(х, у) - у(х, у), и полученный результат использовать в (13) для нахождения множества значений функции ф(х, у, г) потенциала скоростей нижней жидкости:
/ N / 1 ^дG2(^,л,-2й,х,У,г)
ф(х,у,г) - -(—) II у(£, л)—2 ' •
4л -/о дС
( \
После этого необходимо найти множество значений функции ф хУ, 2 по формулам (17)-
дх
(21) и использовать полученный результат для нахождения высоты волны в долях полухорды крыла по формуле (14).
В качестве примера был выполнен расчет волнообразования для крыла с плоским профилем размахом 40 м, хордой 10 м, углом атаки 3° при его движении на высоте 2 м над разделом сред воздух-вода со скоростью 100 м/с. Результаты расчета представлены на рис. 2 и рис. 3. На рис. 2 изображена плотность двойного слоя, характеризующая нагрузку на крыле. На
рис. 3 представлен рельеф свободной поверхности под крылом. Её форма близка к волнам при движении сосредоточенного давления [4].
Рис. 2. Плотность двойного слоя
Рис. 3. Форма свободной поверхности жидкости под крылом
Таким образом, в задачах о движении на малой высоте над поверхностью раздела жидкостей с различной плотностью несущей поверхности оказывается эффективным квад-рупольное вырождение фундаментальных структур. Построенная численная схема решения задачи может быть просто обобщена на случай движения слабо искривленных тонких крыльев. Также просто можно вычислить силы и моменты на крыле.
Библиографический список
1. Панченков, А.Н. Гидродинамика подводного крыла / А.Н. Панченков. - Киев: Наук. думка, 1965. - 552 с.
2. Панченков А.Н. Линейные задачи гидродинамики крыла над поверхностью раздела жидкостей разных плотностей // Гидродинамика больших скоростей. Киев: Наук. думка, 1965. Вып. 1. С. 7-20.
3. Панченков, А.Н. Квадрупольная теория крыла вблизи твердой границы // Асимптотические методы в динамике систем. - Новосибирск: Наука, 1980. С. 5-116.
4. Орлов, Ю.Ф. Потенциал ускорений в гидродинамике корабельных волн / Ю.Ф. Орлов. - Новосибирск: Наука, 1979. - 214 с.
5. Орлов Ю.Ф. Алгоритмы расчета формы свободной поверхности тяжелой жидкости под низколетящим крылом / Ю.Ф. Орлов // Асимптотические методы в механике. - Новосибирск: Наука, 1983.
6. Басин, М.А. Гидроаэродинамика крыла вблизи границы раздела сред / М.А. Басин, В.П. Шадрин. - Л.: Судостроение, 1980. - 304 с.
7. Сретенский, Л.Н. Теория волновых движений жидкости / Л.Н. Сретенский. - М.: Наука, 1977. - 816 с.
8. Мэтьюз Джон, Г. Численные методы. Использование MATLAB: [пер. с англ.] / Г. Мэтьюз Джон, Финк Д. Куртис. - 3-е изд. - М.: Издательский дом «Вильямс», 2001. - 720 с.
9. Марчук, Г.И. Методы вычислительной математики / Г.И. Марчук. - М.: Наука, 1977. - 456 с.
Дата поступления в редакцию 02.08.2011
S.A. Khimich, Y.F. Orlov
THE MOTION OF THE WING OF THE FINITE WINGSPAN ABOVE THE INTERFACE OF TWO MEDIUMS
The problem of surface waves generation at the motion of the wing of the finite wingspan at the small height above the interface of two mediums is considered here. Two-dimensional singular integral equation of the wing theory degenerates to the boundary problem for the Poisson equation for finding the double-layer density distributed on the wing surface. This problem was solved numerically. The case of the motion of plain wing which has the size of real ekranoplan wing is considered.
Key words: interface, ekranoplan wing, quadrupole theory of the wing, waves under the wing.
