Расчёт параметров электромагнитного поля, рассеянного морской поверхностью при малых углах скольжения

УДК 537.874.31

Ф.Н. Захаров, Ю.П. Акулиничев

Расчёт параметров электромагнитного поля, рассеянного морской поверхностью при малых углах скольжения

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

Ключевые слова: параболическое уравнение, численные методы, метод конформного отображения, коэффициент преломления, коэффициент отражения.

Реальная морская поверхность всегда в той или иной степени неровная, в силу чего отражение от такой поверхности сопровождается явлениями, которые отсутствуют в случае идеально гладких границ раздела. Степень «гладкости» определяется в первую очередь соотношением между длиной волны и геометрическими параметрами неровностей.

Характер рассеяния определяется многими факторами. Кроме размеров неровностей и длины волны падающего излучения, играют роль размеры рассеивающей площади, способ её облучения, а также поляризация первичной волны, отражающие и преломляющие свойства воды и т.д. В зависимости от соотношения между различными параметрами применяют те или иные приближённые методы расчёта рассеянного поля [1, 2]. Два наиболее простых и часто применяемых метода - это метод малых возмущений [3, 4] и метод Кирхгофа [1, 2, 5]. Метод малых возмущений пригоден для описания поверхностей с малыми неровностями и углами наклона в единицы градусов. Взволнованная морская поверхность не всегда удовлетворяет данному требованию. Приближение Кирхгофа применимо для поверхности с радиусами кривизны, значительно превышающими длину волны [2], и даёт хорошие результаты при достаточно больших углах скольжения радиоволны. При малых углах скольжения необходимо учитывать затенения отдельных участков морской поверхности. Существует также метод поверхностных интегральных уравнений [6, 7]. Однако данный метод применим больше для решения задачи дифракции на небольших диэлектрических телах простой формы, а при расчёте поля над импедансными телами он даёт значительные ошибки [8].

Перспективным способом учёта рельефа подстилающей поверхности, свободным от перечисленных недостатков, является метод конформного отображения [9] криволинейной системы координат, согласованной с неровной морской поверхностью, в декартову систему координат над плоской поверхностью. Учёт якобиана такого преобразования в уравнении Гельмгольца приводит к введению искусственного поля коэффициента преломления тропосферы [10]. В этом случае характеристики поля можно вычислить, используя хорошо освоенные сеточные методы численного решения параболического уравнения (ПУ), в частности, эффективный метод, основанный на быстром преобразовании Фурье [11]. Недостаток предлагаемого подхода в том, что конформное отображение существует лишь для двумерного ПУ.

Параболическое уравнение выводится из уравнения Гельмгольца в предположении, что энергия радиоволны, рассеянная в обратном направлении, мала. Оно имеет вид

где и(х,г) = Е(х,^)ехр(-/га.г) - комплексная огибающая монохроматического поля Е(х, г) частоты ю, распространяющегося преимущественно в направлении оси Ох, к = 2лД - волновое число, X - длина радиоволны, п(х, г) - коэффициент преломления среды распространения.

Основным методом численного решения ПУ является сеточный метод, который предполагает, что область расчёта покрывается прямоугольной сеткой с ячейками размером Ах • Дг, а значения напряжённости поля вычисляются в каждом узле этой сетки. Процедура численного решения ПУ заключается в следующем. На первом этапе задаются значения отсчётов начального поля в узлах сетки при х = 0 и затем, шаг за шагом удаляясь от источника, находятся значения поля ип т = и(пДх,тДг) во всех узлах сетки [11]. Неоднородности среды задаются в виде пространственного распределения коэффициента преломления на каждом шаге сетки по дальности и учитываются при помощи метода расщепления [11].

Эквивалентный коэффициент преломления тропосферы. Нерегулярность и хаотичность морской поверхности дают основание рассматривать такое волнение как случайный вероятностный процесс и использовать для его описания спектральный метод. В основе спектрального метода лежит представление о волновом процессе как суперпозиции достаточно большого числа волн малой амплитуды [12]

<х>

£

/ =1

где = у[к^ - частота /-й гармоники с амплитудой а,; к/ = 2%к- ^ - волновое число или пространственная частота; - длина волны /-й гармоники; g = 9,8 м/с2 - ускорение силы тяжести; ф; - начальная фаза /-й гармоники, распределённая равномерно на интервале [0, 2п]. Следует отметить, что в связи с особенностью получения экспериментальных данных модель (2) хорошо описывает характер морской поверхности лишь в направлении, соответствующем направлению движения морской волны, поэтому и направление трассы распространения радиоволн (РРВ) должно быть таким же. Это существенное ограничение, но, к сожалению, более полной модели пока не существует. Основываясь на теории случайных процессов, выражение (2) можно представить в виде

<Х> гС Б а/ соь

£( х) = £ /=1

С Б

где а/ и а/ - взаимно некоррелированные нормальные случайные величины с нулевыми математическими ожиданиями и с равными при одинаковых индексах / дисперсиями. В выражении (3) опущена зависимость амплитуды от времени, т.е. рассматривается зависимость только от пространственной координаты.

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

<Х> г -|

О

/=1

Используя якобиан (4), были получены выражения для эквивалентного коэффициента преломления тропосферы

-6 М Г С Б 1

п

/=1 1 -1 и его пространственной корреляционной функции [13]

М

Я

/=1

где Дх = х1 - х2, М - количество учитываемых гармоник морского волнения, 2\ и - координаты по высоте, к/ = 2п/X/ - волновое число, соответствующее морскому волнению с длиной волны X/, ст/ - дисперсия /-й гармоники морского волнения, N(х,г) - индекс преломления. Выражения (4) и (5) справедливы при выполнении условия [9] t(1 - t) 1 ^ 1, где t = £а/к/ .

Для того чтобы использовать выражение (6) при расчёте среднего поля, необходимо определить дисперсию каждой гармонической составляющей ст2 . Для этого можно воспользоваться энергетическим спектром морского волнения. В качестве типичного частотного спектра реального полно-

стью развитого морского волнения наибольшее распространение получил спектр Пирсона-Московица [12, 14]:

S (со) = î

-ш 5 • exp

-1,25

/ \4 с

V Ш J

где ар = 8,1-10-3 - параметр, определяющий форму высокочастотной ниспадающей ветви спектра; юр = 0,9§ / W10 - частота, соответствующая максимуму спектральной плотности морского волнения при заданной скорости ветра; Wl0 - скорость ветра на высоте 10 м над уровнем моря. На рис. 1 приведена форма спектра морского волнения, рассчитанная по формуле (7) для двух значений скорости ветра.

Таким образом, дисперсия ий гармонической составляющей равна ст2 = 5(ш,) • Дш , где Дш - шаг дискретизации спектра морского волнения:

Дш=.Шпих mn

M

где Ютт, Ютах - минимальная и максимальная частоты, в пределах которых сосредоточено 95% энергии морского волнения. В работе [14] указано, что ютт = 0,2юр, Ютах = 2,5Юр.

Так как в формуле (6) используется значение пространственной частоты к, а не циклической частоты ю, то удобнее перейти к представлению спектра в области пространственных частот:

о

H

с

3

к л

! U С О

0,1

10 м/с

P

\

I

7 м/с -

0 1 2 3

Частота, рад/с Рис. 1. Форма энергетического частотного спектра морского волнения для значений скорости ветра 7 и 10 м/с

Sk (k ) =

p • k"5/2 • exp

-1,25

= Sk (ki У

ki ,

2 2 где ki = i /g - пространственная частота i-й гармоники, kp = ссp /g - пространственная частота,

соответствующая максимуму спектральной плотности морского волнения, yjg

ki -Vk-Î) = Аш •

Тогда корреляционная функция эквивалентного коэффициента преломления будет иметь вид:

M 2 ___

Req (Дг, Z1, Z2) = 4^ ki • S ((g) -Jg (( -Jk~i) • exp [-ki ( + Z2 ) • cos[[ Ax]. i=1

H u

S

m £

1

0,5

-0,5

На рис. 2 приведен пример зависимости коэффициента корреляции от величины разнесения точек по горизонтали Дх на разных высотах над уровнем моря для скорости ветра 10 м/с.

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

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

-1

V» 1 • î

\ 1 î \ t i

> у*.— » V -v- г*"

0

20

60

100

40

Расстояние, м

--У поверхности моря ----На высоте 10 м .....- На высоте 20 м

Рис. 2. Нормированная корреляционная функция эквивалентного коэффициента преломления при скорости ветра 10 м/с

будет уменьшаться при увеличении скорости ветра и увеличении высоты над уровнем моря. Данный вывод подтверждается рис. 3 и 4, на первом изображены нормированные спектры эквивалентного коэффициента преломления на высоте 10 и 20 м над уровнем моря при скорости ветра 10 м/с, а на втором - зависимость ширины спектральной плотности от скорости ветра (сплошные линии).

й 2

а 5

а &

и С о Л

§ т

§ о

а *

О т

К §

т

1

0,8 0,6 0,4 0,2 0

г. 3 1 |/ 1 \ \ \ 1 \ — - Высота 10 м ----Высота 20 м

31 л ч 1 1

1 1 1 \ ч >

ъ _

—.

0

0,02

0,04

0,16

0,18

0,2

0,06 0,08 0,1 0,12 0,14 Пространственная частота, 1/м Рис. 3. Нормированная спектральная плотность мощности эквивалентного коэффициента преломления

0Л2й

£ Е 0,08

а 5

С о

° к § 0,06,

Высота 10 м

Высо та 20 м

1 3 5 7 9 11 13

Скорость ветра, м/с

Рис. 4. Зависимость ширины спектральной плотности от скорости ветра

15

Зависимости на рис. 4 были аппроксимированы экспоненциальной функцией, вида ДО(^10) = Aexp(—bWlo) (штриховые линии), где параметры А и Ь определялись при помощи метода наименьших квадратов. Их значения представлены в таблице.

Значения параметров экспоненциальной функции аппроксимации

Параметр Высота над уровнем моря

10 м 20 м

А, м-1 0,065 0,030

Ь, с/м 0,1 0,07

Для сопоставления значений эквивалентного коэффициента преломления и коэффициента преломления тропосферы на рис. 5 представлены высотные профили СКО эквивалентного индекса преломления для трёх значений скорости ветра.

Используя данные результаты, можно оценить высоту, на которой ещё нужно учитывать этот эквивалентный индекс преломления в сравнении, например, с турбулентными флуктуациями индекса преломления тропосферы. Из рис. 5 видно, что до высоты 200-300 м влияние морского волнения превышает влияние неоднородностей тропосферы (СКО эквивалентного индекса преломления больше нескольких #-ед.). Таким образом, минимальная высота области численного расчёта над морем при использовании сеточного метода решения ПУ должна быть не менее 200-300 м, вследствие влияния морского волнения на процесс РРВ.

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

Ф.Н. Захаров, Ю.П. Акулиничев. Расчёт параметров электромагнитного поля 1000

800 600

й н о

§ 400 т

200

0

100 101 102 103 104 105

СКО экв. инд. преломления, Ы-ед. Рис. 5. Высотные профили СКО эквивалентного индекса преломления

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

Так как рассеяние радиоволн зависит от соотношения высоты морских волн и длины радиоволны X, то вместо скорости ветра Р^ю будем рассматривать нормированную скорость ветра:

^10

и=-

Я, = Я  ехр

ехр

и 2 • 2 -Ьи 81П

Пример оценки модуля коэффициента отражения радиоволны от взволнованной морской поверхности при значениях нормированной скорости ветра и 19,8; 28,3 и 42,5 с-1 (что при частоте 850 МГц соответствует скорости ветра 7, 10 и 15 м/с) представлен на рис. 6. Для сравнения приведены значения коэффициента отражения, рассчитанные по известной формуле [15]

\2 '

и = 19,8 с-1

Метод стат. испытаний Теоретическая форм. Аппроксимация

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

На рис. 6 представлены также результаты аппроксимации экспериментальных кривых функцией вида

1

0,8

а у 0,6 На 0,4

0 0

4 6 8 10 Угол скольжения, град

и = 28,3 с-1

1 Метод стат. испытаний Теоретическая форм. Аппроксимация

12

14

Угол скольжения, град и = 42,5 с-1

Параметр Ь зависит

в 1

5 I 0,8

II 0,6

0,4

т К о о

0,2 0

\\ :  Лу • ! •! Метод стат. испытаний Теоретическая форм. Аппроксимация

а*—п*ч а и а —<

0

2

12

от скорости ветра и равняется 0,60; 0,74 и 2,18 для и, равного 19,8; 28,3 и 42,5 с-1, соответственно.

4 6 8 10 Угол скольжения, град

Рис. 6. Модуль коэффициента отражения от взволнованной морской поверхности

14

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

Угол скольжения, град Рис. 7. Фаза коэффициента отражения от взволнованной морской поверхности для трёх значений нормированной скорости ветра

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

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

Литература

1. Рытов С.М. Введение в статистическую радиофизику. - Ч. II. Случайные поля / С.М. Рытов, Ю.А. Кравцов, В.И. Татарский. - М.: Наука, 1978. - 464 с.

2. Исимару А. Распространение и рассеяние волн в случайно-неоднородных средах: в 2 т. -Т. 2. - М.: Мир, 1981. - 320 с.

3. Fuks I.M. Wave diffraction by a rough boundary of an arbitrary plane-layered medium // IEEE Transactions on Antennas and Propagation. - 2001. - Vol. 49, № 4. - P. 630-639.

4. Брюховецкий А. С. О методе малых возмущений в теории рассеяния волн статистически неровной поверхностью // Радиофизика и радиоастрономия. - 2006. - Т. 11, № 3. - С. 254-263.

5. Басс Ф.Г. Рассеяние волн на статистически неровной поверхности / Ф.Г. Басс, И.М. Фукс. -М.: Наука, 1972. - 424 с.

6. Vasil'ev E.N. The integral equation method in the problem of electromagnetic waves diffraction by complex bodies / E.N. Vasil'ev, V.V. Solodukhov, A.I. Fedorenko // Electromagnetics. - 1991. - Vol. 11, № 2. - P. 161-182.

7. Комаров А.А. Дифракция плоской электромагнитной волны на диэлектрической ступеньке / А.А. Комаров, В.А. Пермяков // Журнал радиоэлектроники [электронный журнал]. - 2013. - № 8. -Режим доступа: http://jre.cplire.ru/jre/aug13/7/text.pdf, свободный (дата обращения: 23.08.2015).

8. Комаров А.А. Применение методов поверхностных интегральных уравнений и параболического уравнения к анализу распространения радиоволн вдоль земной поверхности сложного профиля / А.А. Комаров, М.С. Михайлов, В.А. Пермяков // IV Всероссийские Армандовские чтения: Радиофизические методы в дистанционном зондировании сред [Электронный ресурс]. - Муром: Изд.-полиграф. центр МИ ВлГУ, 2014. - Режим доступа: http://www.mivlgu.ru/conf/armand2014/ rmdzs-2014/pdf/S2_13.pdf, свободный (дата обращения: 23.08.2015).

9. Новиков А.В. Конформное отображение как метод учёта неровностей подстилающей поверхности с малой крутизной при численном решении параболического уравнения / А.В. Новиков, Ю.П. Акулиничев // Доклады ТУСУРа. - 2011. - № 1(23). - С. 55-59.

10. Новиков А.В. Сравнительная оценка влияния неоднородностей тропосферы и неровностей земной поверхности на характеристики электромагнитного поля / А.В. Новиков, Ю.П. Акулиничев // Доклады ТУСУРа. - 2011. - № 1(23). - С. 60-64.

11. Levy M. Parabolic equation methods for electromagnetic wave propagation. - London: The Institution of Electrical Engineers, 2000. - 336 р.

12. Бородай И.К. Мореходность судов / И.К. Бородай, Ю.А. Нецветаев. - Л.: Судостроение, 1982. - 288 с.

13. Захаров Ф.Н. Корреляционная функция эквивалентного индекса преломления над морской поверхностью // Матер. 51-й Междунар. науч. студ. конф. «Студент и научно-технический прогресс»: Математика. - Новосибирск: Новосиб. гос. ун-т, 2013. - 250 с.

14. Федосов А.Н. Моделирование волнения морской поверхности для решения прикладных задач / А.Н. Федосов, И.Л. Филимонов // Збiрник наукових праць Академп вшськово-морських сил iменi П.С. Нахiмова. - 2010. - Вып. 3. - С. 46-55.

15. Dockery G.D. An improved impedance-boundary algorithm for Fourier split-step solutions of the parabolic wave equation / G.D. Dockery, J.R. Kuttler // IEEE Transactions on Antennas and Propagation. -1996. - Vol. 44, № 12. - P. 1592-1599.

Захаров Фёдор Николаевич

Мл. науч. сотрудник НИИ радиотехнических систем ТУСУРа

Тел.: 8 (382-2) 41-38-89

Эл. почта: fzakharov89@gmail.com

Акулиничев Юрий Павлович

Д-р техн. наук, профессор каф. радиотехнических систем ТУСУРа

Тел.: 8 (382-2) 41-36-70

Эл. почта: ayp@rts.rk.tusur.ru

Zakharov F.N., Akulinichev Yu.P.

Calculation of parameters of electromagnetic field scattered by sea surface at small grazing angle

The goal of the paper is calculation the parameters of electromagnetic field, propagating over rough sea surface. Statistical characteristics of the field of artificial refractive index, arising from the conformal mapping of curvilinear coordinate system over an uneven sea surface in a Cartesian coordinate system over a flat surface, were calculated for the Pierson-Moskowitz sea wave spectrum. The average value of complex coefficient of mirror reflection from the sea surface depending on the grazing angle was estimated by repeated simulation of realizations of such a field and the numerical solution of the parabolic equation.

Keywords: parabolic equation, numerical methods, conformal mapping method, refractive index, reflectance.