Разработка и исследование методов численного решения задачи распространения электромагнитных волн в неоднородной тропосфере тема диссертации и автореферата по ВАК РФ 05.13.18, кандидат наук Лытаев Михаил Сергеевич

• Лытаев Михаил Сергеевич
• кандидат наук
• 2019

Оглавление диссертации кандидат наук Лытаев Михаил Сергеевич

Введение

Глава 1. Анализ методов решения задачи распространения

электромагнитных волн вблизи поверхности Земли

1.1 Исследование исходной математической модели

1.2 Анализ методов моделирования тропосферного распространения электормагнитных волн

1.3 Исследование методов численного решения параболического уравнения

1.4 Поиск аналогичных задач, моделей и методов в других научных областях

1.5 Изучение опыта проведения натурных экспериментов

1.6 Классификация программных комплексов расчета тропосферного распространения

1.7 Постановка научно-технической задачи исследования

1.8 Выводы

Глава 2. Разработка численного метода решения уравнения

Гельмгольца в неоднородной тропосфере

2.1 Формализация постановки задачи

2.2 Аппроксимация оператора пошагового распространения

2.3 Построение разностной схемы по поперечной координате

2.4 Дискретизация условия третьего рода на нижней границе

2.5 Анализ возможности сквозного интегрирования через неоднородности

2.6 Построение дискретных прозрачных граничных условий

2.7 Анализ размера расчетной сетки и порядка аппроксимации

2.8 Обобщение метода на трехмерное уравнение Гельмгольца

2.9 Проведение и анализ результатов численных экспериментов

2.10 Выводы

Стр.

Глава 3. Разработка метода моделирования направленных

источников излучения

3.1 Постановка задачи

3.2 Редукция к объемному источнику

3.3 Учет переменного индекса преломления

3.4 Особенности построения численного решения

3.5 Выводы

Глава 4. Разработка метода моделирования рассеяния

электромагнитных волн на массиве из тонких ребер

4.1 Постановка задачи

4.2 Вывод модельного уравнения

4.3 Исследование свойств системы модельных уравнений

4.4 Усеченное модельное уравнение

4.5 Разрешимость модельного уравнения

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

4.7 Проведение и анализ результатов численных экспериментов

4.8 Выводы

Глава 5. Синтез архитектуры программного комплекса

5.1 Определение состава комплекса и используемых программных средств

5.2 Разработка объектно-ориентированной модели предметной области

5.3 Разработка вспомогательных модулей

5.4 Пример использования программного интерфейса

5.5 Выводы

Заключение

Список сокращений и условных обозначений

Список литературы

Приложение А. Акт о внедрении результатов работы

Введение диссертации (часть автореферата) на тему «Разработка и исследование методов численного решения задачи распространения электромагнитных волн в неоднородной тропосфере»

Введение

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

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

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

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

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

Степень разработанности темы. Классические фундаментальные результаты теории распространения электромагнитных волн вблизи поверхности Земли получены в исследованиях А. Зоммерфельда, В.А. Фока [1], М.А. Леонтовича, Е.Л. Фейнберга [2], Г.И. Макарова [3] и других ученых. К сожалению, далеко не все классические результаты могут быть доведены до практического использования. Методы, основанные на интегральных уравнениях и разложении по собственным функциям, могут дать решения некоторых простых модельных задач, однако они не применимы в сложных условиях распространения. Для расчетов в сложных областях была разработана менее строгая геометрическая теория дифракции [4; 5]. Широко используется аппроксимация неровной поверхности Земли массивом их тонких ребер [6] или клинов, хотя даже в такой упрощенной постановке строгое решение представляется труднореализуемым.

Необходимо отдельно выделить разработанный советскими учеными М.А. Леонтовичем и В.А. Фоком [1; 7] метод параболического уравнения (ПУ). С одной стороны, метод ПУ основан на строгих уравнениях Максвелла, с другой — этот метод показал себя достаточно удобным для численной реализации. Основные результаты западных ученых по данной тематике посвящены особенностям построения эффективных численных методов. R.H. Hardin и F.D. Tappert разработали метод расщепления Фурье для численного решения ПУ с переменными коэффициентами [8]. В работах J.R. Kuttler, G.D. Dockery [9] и A.E. Barrios [10] данный метод адаптирован для задачи тропосферного распространения. Метод ПУ получил широкое распространение для решения аналогичных задач в гидро-

акустике [11; 12], оптике и геофизике. Это способствовало более интенсивному теоретическому исследованию и разработке численных методов для решения ПУ Теоретический анализ метода ПУ, его обобщение на случай произвольного угла распространения и произвольного коэффициента преломления рассмотрены в работах A. Bamberger [13; 14] и L. Fishman [15]. Строгое математическое обоснование на основе теории операторов изложено в работах J.P. Corones [16; 17]. Особенности построения решения в неограниченных областях отражены в работах А.В. Попова [18; 19], M.F. Levy [20], M. Ehrhardt [21]. Следует отметить монографию M.F. Levy [22], описывающую важнейшие достижения в области прогнозирования радиовидимости методом ПУ к 2000 году. На базе метода ПУ реализован ряд программных комплексов: AREPS, TEMPER [23], CARPET, PETOOL [24].

В настоящее время ведутся разработки моделей тропосферного распространения в сложных условиях среды. В работах M.D. Collins [25] и O. Ozgun [26] на основе метода ПУ разработан численный метод учета обратного рассеяния от препятствий. Особенности прохождения волн через лесные массивы рассмотрены в работах В.А. Пермякова, М.С. Михайлова [27] и P. Holm [28]. R. Janaswamy предложил векторное трехмерное ПУ для решения задачи в условиях городской застройки [29]. Метод трехмерного ПУ получил развитие в работах L. Zhou [30], М.С. Михайлова [31] и других ученых. Особенности применения математического моделирования для решения конкретных практических задач изложены в монографиях R. Janaswamy [32], J.S. Seybold [33], N. Blaunstein [34] и других исследователей. Последние достижения отечественных ученых изложены в работах Ф.Н. Захарова [35], П.Н. Дагурова [36], В.А. Пермякова [27], Ю.П. Акулиничева [37].

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

Для достижения указанной цели были поставлены и решены следующие задачи:

1. Исследовать существующие математические модели, численные методы и комплексы программ, предназначенные для прогнозирования зон радиовидимости в различных условиях.

2. Разработать пошаговый численный метод решения уравнения Гельм-гольца в неограниченной области.

3. Разработать метод моделирования направленных источников излучения.

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

5. Разработать комплекс программ для проведения вычислительных экспериментов.

Объект исследования: Математическая модель распространения волн.

Предмет исследования: Методы решения задачи распространения электромагнитных волн в неоднородной тропосфере.

Научная новизна:

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

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

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

4. В отличие от аналогов, предназначенных для решения данного класса задач, разработанный комплекс программ основан на детерминированных моделях теории распространения волн, объектно-ориентированной модели предметной области и организован в виде программной библиотеки.

Теоретическая и практическая значимость. Теоретическая значимость изложенных в работе положений состоит в следующем:

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

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

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

Практическая значимость полученных результатов состоит в следующем:

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

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

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

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

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

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

Основные положения, выносимые на защиту:

1. Численный метод решения уравнения Гельмгольца в неоднородной тропосфере.

2. Метод моделирования направленных источников излучения.

3. Метод моделирования рассеяния электромагнитных волн на массиве из тонких ребер.

4. Комплекс программ для расчета характеристик тропосферного распространения электромагнитных волн.

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

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

- 9-й конференция «Информационные технологии в управлении» (ИТУ-2016), АО «Концерн «ЦНИИ «Электроприбор», Санкт-Петербург, 4-6 октября 2016 г.

- XV Санкт-Петербургская международная конференция «Региональная информатика (РИ-2016)», Санкт-Петербург, 26-28 октября 2016 г

- VI, VII Международная конференция «Актуальные проблемы инфоком-муникаций в науке и образовании» (АПИНО), СПбГУТ им. проф. М.А. Бонч-Бруевича, Санкт-Петербург, 1-2 марта 2017 г., 28 февраля - 1 марта 2018 г.

- The 8th International Symposium «Information Fusion and Intelligent Geographie Information Systems» (IF&IGIS'17), Шанхайский морской университет, Шанхай, КНР, 10-12 мая 2017 г

- Annual international conference «Days on Diffraction 2017», ПОМИ РАН, Санкт-Петербург, 19-23 июня 2017 г

- Городской семинар по вопросам теории распространения волн под руководством В.М. Бабича, ПОМИ РАН, Санкт-Петербург, 24 октября 2017 г, 18 сентября 2018 г.

- The 10th International congress «Ultra Modern Telecommunications and Control Systems (ICUMT2018)», Москва, 5-9 ноября 2018 г

- Advances in Wireless and Optical Communications (RTUWO'18), Рижский Технический Университет, Рига, Латвия, 15-16 ноября 2018 r.

- IEEE Conference of Russian Young Researchers in Electrical and Electronic Engineering (ElConRus), ЛЭТИ, Санкт-Петербург, 28-31 января 2019 r.

Публикации. Основные результаты по теме диссертации изложены в 15 работах [38—52], из них 4 в рецензируемых научных изданиях из перечня ВАК, 4 в научных изданиях, индексируемых международными базами цитирования Web of Science и Scopus, 5 в других изданиях и материалах конференций, 2 свидетельства о государственной регистрации программы для ЭВМ.

Соответствие паспорту научной специальности. Содержание диссертации отвечает пунктам 1 (положения 2 и 3), 3 (положение 1) и 4 (положение 4) паспорта специальности 05.13.18 - математическое моделирование, численные методы и комплексы программ.

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

Глава 1. Анализ методов решения задачи распространения электромагнитных

волн вблизи поверхности Земли

Любые электромагнитные взаимодействия могут быть полностью описаны системой уравнений Максвелла с соответствущими начальными и граничными условиями (ГУ). Решение такой задачи будет полностью соответствовать требованиям и определяется единственным образом. К сожалению, в подавляющем большинстве реальных задач точное замкнутое решение получить не удается. Реальные задачи распространения электромагнитных волн решаются при помощи численного моделирования или проведения натурных экспериментов. На первом этапе строится упрощенная математическая модель, наиболее полно описывающая поставленную задачу Такая модель может основываться на строгих уравнениях Максвелла (детерминированные модели) или статистическом анализе измерений (эмпирические модели). Возможны и комбинированные модели. Возникает задача валидации математической модели на предмет соответствия реальной проблеме [53]. На втором этапе разрабатывается численный алгоритм для реализации полученной модели на компьютере. Для численных методов необходимо определить границы их применимости и точность в зависимости от входных параметров. На последнем этапе при помощи построенных численных методов создается программная реализация модели. Здесь необходимо учитывать конкретные особенности эксплуатации программ и их целевую аудиторию. Проверка полученной триады «модель—алгоритм—программа» [54] осуществляется при помощи сравнения с имеющими замкнутое решение простыми модельными задачами, другими аналогичными моделями или результатами натурных экспериментов. При этом одновременно проверяется как адекватность модели, так и корректность использованных численных схем и их программной реализации.

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

основе проведенного анализа выбраны перспективные задачи, решаемые в данной диссертации.

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

Далее во всей работе предполагается гармоническая зависимость электромагнитной волны от времени вида вхр(—гш1), где ш — угловая частота. Тогда уравнения Максвелла для векторов электрического поля Е и магнитного поля Н в декартовых координатах с учетом материальных уравнений запишутся следующим образом

1.1 Исследование исходной математической модели

1.1.1 Уравнения Максвелла для гармонической волны

Ух Е =

(1.1)

V х Н

4п

j — ¿кеЕ,

(1.2)

с

У • (еЕ) = 0,

У • (^Н) = 0,

где

- к = ш/с = — волновое число;

- Л — длина волны в атмосфере;

2п Л

- j(x, у, г) — вектор плотности электрического тока;

- £ = £(ш,х,у, г) — относительная диэлектрическая проницаемость среды;

- ^ — относительная магнитная проницаемость.

В задачах тропосферного распространения обычно полагается, что ^ = 1 [2]. Запишем ток j в виде суммы токов источника и индуцированных токов

j = jsrc + jшd. (1.3)

Для индуцированных токов применим закон Ома

jwd = стЕ, (1.4)

где ст = ст(ш, х, у, г) — проводимость. С учетом соотношений (1.3) и (1.4), уравнение (1.2) может быть перезаписано следующим образом

4п

Ух Н = —'^с - 1Ы2Е, (1.5)

с

где коэффициент преломления п = п(ш, х, у, г) определяется следующим образом

2 4Ыст

П = £ +--. (1.6)

Ш

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

Применяя оператор Ух к обоим частям уравнений (1.1) и (1.5), приходим к векторным уравнениям

ЛЕ - У (У • Е) + гк2п2Е = - ^^гс) ^

АН + к^п2 ^Н - п У ^ ^г) х (Ух Н) = -п2—У х

п

1

"2 j

Здесь использовано известное тождество

Ух Ух = У(У^) - Л.

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

ЛЕ + гк2п2 Е = ^^ (1.7)

с

4п

ЛН + к2п2Н =--V х jsrc. (1.8)

с

Уравнения (1.7) и (1.8) могут быть связаны через ГУ Волновое поле в данном случае порождается известным током источника в правой части. Тут следует заметить, что знание функции jsгc подразумевает моделирование внутренней структуры излучающей антенны и протекающих в ней токов, что представляет собой отдельную трудоемкую задачу При этом излучающие антенны обычно задаются своей диаграммой направленности в дальней зоне. В связи с этим актуальной является задача моделирования направленных источников излучения в уравнении Гельмгольца по известной диаграмме направленности антенны (ДНА).

1.1.2 Вывод двумерного уравнения Гельмгольца

Важным упрощением уравнений Максвелла (1.1) и (1.2) является двумерная задача, когда геометрия не зависит от какой либо пространственной координаты, в нашем случае координаты у. Приравняв в уравнениях (1.1) и (1.5) все частные производные по у к нулю, можно получить две независимые системы скалярных уравнений

дЕу дх

—

< дЕу — -гкЯ7

дх

дНг дх

(1.9)

дх

— 4-П]у - гкп2ву

дНу

— - -Пзх + гкПЕх

= - -?3х + гкп2Ех

дх

дНу дх

дЕг дЕх =

дх дх = АкИУ

(1.10)

х

Система (1.9) соответствует горизонтальной поляризации. Подставляя соответствующие компоненты, данная система сводится к двумерному скалярному уравнению Гельмгольца относительно Еу

д2Еу д 2Е

+ —2 + к2п2Еу = дг2 у

4пк

~3у.

дх2 дг* " с

Несколько более сложным является случай вертикальной поляризации, которому соответствует система (1.10). Для простоты опустим внешние токи (] = 0). Тогда система (1.10) может быть преобразована к следующему уравнению

д_ дх

1_ дИу п2 дх

_д_ дг

1_ дИу п2 дг

+ к2 Ну = 0.

Чтобы свести данное уравнение к более простому виду, введем в рассмотрение следующую замену переменных

ф = Н.

п

При этом предполагается, что функция п дифференцируема необходимое число раз, а случаи резкого скачкообразного изменения п рассматриваются отдельно. Тогда приходим к следующему уравнению Гельмгольца относительно функции ф

д2ф д2ф /, 2 2 д

ТГ2 + ТТ2 + к2п2 + п— дх2 дг2 V дх

1 дп п2 дх

д

+ п— дг

1 дп

п2 дг

ф = 0.

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

?Ф + Й + ^ = 0,

дх2 дг2

(1.11)

где

ф(х,г) =

Еу (х,г) для горизонтальной поляризации,

И(х,г)/п(х,г) для вертикальной поляризации.

1.1.3 Сферичность земной поверхности и коэффициент преломления

тропосферы

Все рассмотренные ранее уравнения выписаны в декартовой системе координат. Вместе с тем, при моделировании загоризонтного распространения волн необходимо учитывать кривизну поверхности Земли. Замена переменных на сферические координаты сопряжена с рядом проблем. Гораздо более удобной является замена переменных, предложенная Пекерисом [55]. Записывая уравнения Максвелла в сферических координатах (г,0,ф) и применяя замену

ж = RQ

z = R log (1 + R)

уравнение на сферической поверхности Земли может с хорошей точностью быть представлено уравнением Гельмгольца (1.11), где функция n2(n,z) заменяется модифицированным коэффициентом преломления

2 / \ 2/ \ 2z m (x,z) = n (x,z) + —,

R

где R — радиус Земли.

В радиофизике коэффициент преломления тропосферы чаще всего представляется при помощи так называемых N и M профилей, которые определяются выражениями

N = (n - 1)106,

М = (т - 1)106.

Отметим, что в определении М-профиля по умолчанию учтена кривизна Земли.

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

М

Рисунок 1.1 — Схематичное изображение типовых М-профилей коэффициента

преломления тропосферы.

Профиль коэффициента преломления тропосферы выражается через основные метеорологические параметры следующей формулой [57]

п = 1 + 7.76 х 10-5Р + 0.373е/Т2, (1.12)

где Р — атмосферное давление (мбар), Т — температура (К), в - давление водяного пара (мбар). Из-за различного нагрева земной поверхности, испарения воды из океанов, озер, рек, горизонтальных и вертикальных ветров атмосферное давление, температура и содержание водяного пара варьируются по горизонтали и вертикали по всей тропосфере, что приводит к изменению показателя преломления. Указанные атмосферные параметры претерпевают сильное пространственное изменение по высоте, в то время как горизонтальные изменения обычно выражены достаточно слабо. Основываясь на многочисленных измерениях температуры, влажности и давления, Международный Союз Электросвязи (МСЭ) рекомендует [58] использовать следующее выражение для рефракции в стандартной атмосфере

N (г) = 315 х ехр(—1.36 х 10-5 г).

Для высот, не превышающих нескольких километров, коэффициент преломления в стандартной атмосфере можно считать линейно убывающим со средним градиентом dN/ йг = -39 №ед./км (или йИ/йг = 118 М-ед./км).

Такой простой профиль коэффициента преломления обычно имеет место над сушей. Поверхность воды характеризуется образованием волноводов, приводящих к весьма нетривиальным эффектам распространения. Резкий отрицательный градиент влажности на высоте до нескольких десятков метров над уровнем моря, в соответствии с формулой (1.12), приводит отрицательной рефракции в приводном слое. Далее спад влажности происходит более плавно и начинают доминировать эффекты, связанные с кривизной Земли. Уменьшение и последующее увеличение с высотой коэффициента преломления способствует появлению т.н. волновода испарения [59]. Волновод испарения проявляется почти повсеместно над водоемами и оказывает существенное влияние на работу морских радаров. Еще большее влияние на распространение над морской поверхностью оказывают т.н. поверхностные волноводы. Они могут приводить как к захвату волн и распространению их на значительные расстояния (до нескольких тысяч км), так и к появлению зон тени вблизи излучателя. На рисунке 1.1 схематично изображены М-профили типовых волноводов. Аномальные тропосферные эффекты обычно имеют место на высотах до нескольких сотен метров, а далее не оказывают существенного влияния.

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

Список литературы диссертационного исследования кандидат наук Лытаев Михаил Сергеевич, 2019 год

ИСТОЧНИК

тонкие ребра

Н,

непроницаемая поверхность л/ л

Рисунок 4.1 — Схематичное описание задачи дифракции на тонких ребрах.

4.2 Вывод модельного уравнения

Совершая синус преобразование Фурье по переменной х, уравнение (4.1) запишется следующим образом

_ к2 [

+ (к2 — р2 + ¿к2а)"ф +—п(х,р — р' )"ф(х,р' )(р' = Ь(х)С((р), (4.3)

(х2 у/2п

где

—оо

—

"ф(х,р) = у -П у -ф(х,г) sm(pz)(х, о

—

(((р) = у - У ((х) sm(pх)(х,

......7-1......

-=Л ....../П -

Заметим, что использование синус преобразования Фурье автоматически гарантирует выполнение условия (4.2) на поверхности Земли.

Функция Грина С(х,х',р) для соответствующего однородного уравнения (при п(х,г) = 0) определяется следующим образом

+ (к2 — р2 + гк2а)С = 6(х — х'),

ах2

С

ас

(х

= (3

х—Ух'+

ас

хх

ах

х— х +

= 1.

х— х —

Принимая во внимание принцип предельного поглощения [161], функция Грина запишется в следующей замкнутой форме

С(х,х',р) = — ^ в—у^рР)х—х' \

где

у(р) = а(р) + %((р),

а(р) =

' У(р2 - - к2)2 + а2 — (к2 — -р2)

2 '

Мр2- — к2)2 + а2 + (к2 — р2)

((Р) = — у 2 ■

Применяя построенную функцию Грина к уравнению (4.3), приходим к следующему интегральному уравнению

2 +то +то

"ф(х,р) + ! J С(х,х',р)П(х'р — р)"ф(х',р')(р (х' = <С(х,0,р)<5(р). (4.4)

—то —то

Дальнейшие рассуждения основаны на предположении, что толщина препятствий К много меньше длины волны как внутри препятствия, так и в остальном пространстве. Сначала предположим, что функция *ф(х',р') при фиксированном р' не меняет свой знак на отрезке [хз- — К/2,хз + К/2]. Тогда, пользуясь интегральной

теоремой о среднем, внутренний интеграл в уравнении (4.4) может быть представлен следующим образом

+то

/б(хх',р)и(х'р — р' Жхр'=

—то

+оо

т 1 „

^^Ь, (р — р') С(х,х',р)п (хХ )*ф(х',р' )(хХ =

3=1 —ТО

хо+Н/2

т

£ ь,(р — р) у с(х,х',р)^(х',р'=

3 1 х1 — Н/2

х о+Л./2

т

(р — р' )С(х,хч ,р) ^(х'р' )(х',

3=1 х—Н/2

где точка хС1 принадлежит отрезку [х, — Н/2,х, + Н/2] в соответствии с теоремой о среднем. В общем случае, когда функция *ф(х',р') меняет знак на отрезке [х, — Н/2,х, + Н/2], последний разбивается на соответствующее число отрезков, на которых *ф(х',р') не меняет знак, и указанные выше рассуждения повторяются к каждому из них. Пользуясь малостью Н и заменяя каждую функцию С(х,хС1 ,р) на С(х,х,,р), можно записать приближенное соотношение

+то х о+Н/2

т

I 6(х,х'р)й(х'р—р)Мх'1р)1(х' * ^ь'(р—р)д{х'х';рч Мх!р)(х'

—То 3 1 х—Н/'2

Теперь уравнение (4.4) можно заменить на приближенное

+то X 1+Н/2

2т

к2

"ф(х,р) + / ^ 0(х,х,,р) у Ьз(р — р') у ^(х',р')(х'(р' =

х—Н/2

((р)О(х0,р'). (4.5)

3=1

Введем в рассмотрение новые функции

хз+Н/2 +то

фз (р)= / Г'з(р — р )ар'(х''

хз —Н/2 —то

с помощью которых можно переписать уравнение (4.5) следующим образом

к

^(х,р) = С(х,хз ,р)фз (р) + Я(р)С(х,0,р).

(4.6)

з=1

Выражение (4.6) дает возможность интерпретировать функции фз- (р) как плотности распределения точечных источников.

Подставив выражение (4.6) в уравнение (4.5) и приравняв соответствующие слагаемые, получим следующую систему интегральных уравнений относительно функций фз (р)

т + то

Ф1 (р) + Е I К1 ,з(р,р)Фз(р)(р = ^1(Р)

з=1 —то

т + то

ф2М + Е I К2,з (р,р)Фз (р)(р = 92(р)

з=1 —то

(4.7)

т + то

фт(р) + Е I Кт,з (р,р)фз (р)(р = 9т(р)

з=1 —то

с ядром

к2

Кг.з (рр) = Ьг(р — р)Хгз (р) '

хг + Н/2

Хг,з (р)= J С(х',хз ,р)(х' =

х—Н/2

1

2 (е—У(р')^/2 — 1) ,

гр , - гр .

гз

У2(р') I е—у(р')х—хзI [е—у(р')н/2 — еу(р')н/2] , хг = х

и правой частью

+то

9г(р)= 3г (р — р' )Уг(р )(р

—оо

2т

xi+h/2

1

yi(p) = Q(p') J G(x',0,p')dx' = ^pjQ(P)e-y(p')xi (e-y(p')h/2 - ey(p')h/^ .

Xi-h/2

Принимая во внимание малость h, последние два выражения можно заменить на приближенные

_ / a bhk2 sin ((p - p') Hi) (p^x_x ,

Kij (p,p') «---. vv y e Y(p )X Xj

i J 2n (p - p') y(p') '

o-h Г í Sin ((p - p') H) Q(n')e-y(p)Xi

9i(p) ~ -bh\¡ nj 2 (p - p) Y(p') Q(p )e •

—oo

4.3 Исследование свойств системы модельных уравнений

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

ф(р) = (ф1,ф2, . . . ,фш) , д(р) = {91,92, . . . ,9т) и матричный оператор

\

Кт

Тогда система уравнений (4.7) может быть записана в операторной форме

m р m р

Аф = (S J Kíj(pp(p'w-.-Y^ J

. j 1 — ГУ} j 1 — ПО

вия

(I + А) ф = g. (4.8)

Легко заметить, что для каждой из функций Ьу выполнены следующие усло-

\bj(pi) - Ь3(p2)| < ci\pi -P2I, Vpi,p2 e R

sup \bj (£)\ < C2,

где ci, c2 — некоторые константы.

Потребуем выполнение следующего ограничения на функцию объемного источника Q(z)

sup\Q(p)\ < + то.

p

Данное условие выполнятся для любого разумного выбора источника, в том числе для точечного источника, которому отвечает дельта функция.

Далее будет поставлен вопрос о нахождении решения модельного уравнения (4.8) в банаховом пространстве B, состоящем из ограниченных векторознач-ных функций и снабженном нормой

|Ф\\ = max

i=1..m

SUp \фг(р)\ р^+то

Р_ 2

Лемма 1. Следующие несобственные интегралы сходятся

+то +то

J \Хгj(Р)\dp' <сз, J \Yi(p')\dp < C4,

—то —то

где c3 и c4 — некоторые константы.

Доказательство. Необходимо обратить внимание, что Re(y(p')) = a(p') ~ при \p'\ ^ то. Тогда искомый результат следует непосредственно из явного вида соответствующих подынтегральных функций. □

Лемма 2. Функция g(p) в правой части уравнения (4.8) ограничена.

Доказательство. Для каждой из функций gi(p) имеет место следующая оценка

+то

\gi (p)\ < sup \ bi(^)\ - J \ji(p' )\dp < C2 • C4.

—то

Таким образом, каждая из функций gi (p) ограничена, откуда следует ограниченность g(p). □

Лемма 3. Оператор А отображает равномерно ограниченные функции в равномерно ограниченные.

Доказательство. Справедливость утверждения вытекает из следующих оценок

+

тто

"С f*

\Аф\\ = max sup \ V / KhJ(p,p)ф(p)dp'\ < i=1..m \р\<то p1_J

то

т

max У^C2 - C3 - sup ф(p)\ < C2 - C3 - \ \ф\\.

i=1..mj=1 р<то

□

Лемма 4. Оператор А отображает ограниченные функции в равностепенно непрерывные.

Доказательство. Для каждой из функций yi(p) справедливы соотношения

+

/т

J2(Ki,j (pip ) - Ki j (p2,p ))Ф3 (p )dp \ <

—то j=1

+

/m

\Ki, j (p1'p) - Ki, j (p2p )\dp - 11ф|1 ^

—то j=1

+

1 f m

72= / ~bi(p' - p1) -~bi(p' - p2)\-\Xij(p)\dp' - ||ф|| <

v я_оо j=1

m - Cu , ,, ,,

72n\p1 — p2\- C3

из которых следует справедливость утверждения.

□

Лемма 5.

+то

—oo

sup \Ki j J (pp )\dp' < +то.

\р\<+то

Лемма 6. Функция g(p) непрерывна.

4.4 Усеченное модельное уравнение

Наряду с уравнением (4.8), рассмотрим приближенное к нему уравнение, определенное на конечном интервале интегрирования

(I + Ап) ф = g, (4.9)

где оператор Ап определяется следующим образом

/ +п +п \

/ т т \

Апф = / Ки(р,р')ф3(р'/ Ктл(р,р')фз(р')ф' I . (4.10)

Х^1 —п 3=1 -п )

Оператор Ап определен в банаховом пространстве Сп непрерывных на интервале [—п,п] векторозначных функций с нормой

|ф|| = тах

г=1..т

SUp |фг(р)|

\р\<п

Из Леммы 4 следует, что Ап : Сп ^ Сп.

Установим теперь разрешимость усеченного уравнения (4.9) в пространстве

Сп

Теорема 1. Пусть п — некоторое фиксированное число. Допустим, что однородное уравнение, соответствующее (4.9), имеет только тривиальное решение. Тогда усеченное модельное уравнение (4.9) имеет единственное решение в пространстве Сп для произвольной функции g(p), принадлежащей пространству Сп.

Доказательство. Из указанных выше Лемм 1-4 и теоремы Арцела-Асколи следует компактность оператора Ап в Сп. Таким образом, искомый результат непосредственно следует из теории Рисса-Шаудера [162] и соответствующей теоремы Фредгольма. □

Замечание 1. Условие Теоремы 1 для случая одного препятствия (т = 1) проверяется непосредственно при помощи резольвенты Фредгольма [163].

Замечание 2. Следует напомнить, что физический смысл спектрального оператора р выражается следующем выражением

p = k sin 6,

где 6 - угол между направлением распространения и положительным направлением оси x. При этом значения \p\ > k отвечают эванесцентным волнам, которые обычно игнорируются при тропосферном распространении радиоволн, но могут быть учтены в рамках предложенного метода. Это наблюдение дает возможность обеспечить физически обоснованное усечение системы (4.9), положив n = k. Действительно, игнорируя эванесцентные волны и полагая функцию "ф(х,р) равной нулю при \p\ > k, область интегрирования в (4.3) можно заменить на конечную, что приводит к усеченному уравнению (4.9).

4.5 Разрешимость модельного уравнения

В данном разделе устанавливаются условия, при которых имеет место разрешимость модельного уравнения (4.8). Обозначим полученное в Теореме 1 решение через фп(р). Рассмотрим уравнение (4.9) при фиксированном п и р е (—то, +то). Тогда его решение в пространстве В запишется в виде следующей функциональной последовательности

—Апфп(р)+g(p), р е [—п,п].

Из ограниченности фп(р) и Леммы 5 следует ограниченность Fn(p) (при фиксированном п). Кроме того, как нетрудно заметить, Fn(p) является непрерывной функцией.

Рассмотрим операторы Т и Тп, определяемые следующим образом

Тф = (1 + А) ф, (4.11)

Тпф = (I + Ап) ф, (4.12)

где Т, Тп : В ^ В. Обозначим через ||| • ||| операторную норму.

Теорема 2. Предположим, что существуют сколь угодно большие п, такие, что ||| (Тп)—1 ||| < Ь, где число Ь не зависит от п. Тогда уравнение (4.8) имеет единственное решение.

Доказательство. Из Леммы 5 следует выполнение следующего предельного соотношения

|||Т — Тп||| ^ 0прип ^ ос.

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

|||Т — Тп|||-||| (Тп)—1 ||| < 1,

из которого следует непрерывная обратимость оператора Т [162]. Кроме того, выполняется предельное соотношение |||Т—1 — Т—1||| ^ 0 при п ^ о.

□

Следствие 1. Отметим, что указанное решение непрерывно и является пределом последовательности Fn(p) при п ^ оо.

Замечание 3. Результат, полученный в Теореме 2 носит условный характер. Условие Теоремы 2 можно проверить непосредственно для случая одного препятствия, однако для случая нескольких препятствий конструктивного способа проверки на данный момент не существует.

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

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

2 1.5

1

0.5 0

-2 -1.5 -1 -0.5 0 0.5 1 1.5 2

p/k

Рисунок 4.2 — Функция y(p) при значении диссипативного параметра

а = 1E — 4.

ль - 1/А

Таким образом, увеличивая показатель преломления препятствия и одновременно уменьшая его толщину в той же пропорции, мы оставляем последний параметр также незначительным по сравнению с длиной волны внутри препятствия. С другой стороны, в рассматриваемом случае горизонтальной поляризации, бесконечное увеличение диэлектрической проницаемости препятствия по отношению к его значению во внешней среде приводит, в соответствии с теорией Френеля, к эффекту полного отражения всех волн от препятствия. При этом можно показать [164], что значение поля на границе будет равно нулю. Отметим, что в процессе построения численного решения параметры b и h остаются конечными.

Далее рассматривается усеченное уравнение (4.9) на интервале [—k,k]. Для численного решения системы интегральных уравнений (4.9) используется метод квадратур [159] (Nystrom method). Рассмотрим разбиение интервала интегрирования т : —k = p0 < p\ < ... < ps = k. Размер разбиения можно легко оценить, проанализировав осциллирующие части ядра и правой части уравнения (4.9). Равномерного разбиения обычно достаточно, но при необходимости точно учитывать вертикально распространяющиеся волны, необходимо иметь в виду поведение функции у(р) в окрестностях точекp = ±k. Графики функций y(p) и Y—l(p) изображены на рисунках 4.2 и 4.3 соответственно. Заметим, что абсолютная величина функции y(p) никогда не обращается в ноль, что гарантируется условием а > 0.

I 30

Е

Е 20 со

10

о о

-1.5

-0.5

0 р/к

0.5

1.5

Рисунок 4.3 — Функция у 1 (р) при значении диссипативного параметра

а = 1Е - 4.

Интегралы в (4.9) приближенно записываются с использованием Гауссовой квадратуры

ь

Кг.з (Рр)фп ,3 (Р') ~ ^ П)кКгЛ (Р,Рк)ф3 (Рк) '

ь

к=1

где ,шк - веса квадратуры. Используя данное приближение, запишем следующую систему

Фг(Р) + Е Е WkКг,3 (РРк)Фз (Рк) = 9г(Р)-3=1 к=1

(4.13)

Подставляя точки из разбиения т в систему (4.13), получаем систему линейных алгебраических уравнений размера те х те. Далее подставим полученное численное решение системы в выражение (4.6) и получим приближение функции ■ф(жр) в узлах т. На последнем этапе вычисляется обратное дискретное синус преобразование Фурье [83] от функции ■ф(жр), что дает искомое численное решение для функции -ф(х,г). Наиболее эффективно использовать квадратуру Филона-Кленшоу-Кертиса [146].

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

0

1

2

4.7 Проведение и анализ результатов численных экспериментов

В первом примере проведем моделирование для случая одного препятствия. Используются следующие вычислительные параметры: К = 10-7 м, Ь = 105 м, а = 10-4. Выбранные параметры К и Ь гарантируют непроницаемость препятствия в соответствии с указанными выше соображениями. Высота препятствия равна 50 м, расстояние от источника — 500 м. Точечный источник расположен на высоте 50 м и излучает гармонические радиоволны горизонтальной поляризации на частоте 300 МГц. Двумерное распределение интенсивности электромагнитного поля изображено на Рисунке 4.4. Интерференция падающих и отраженных волн четко наблюдается перед препятствием. Как и следовало ожидать, за препятствием образовалась зона тени.

Графики решения интегрального уравнения ф(р) и правой части д(р) для данного примера изображены соответственно на рисунках 4.5 и 4.6. Хорошо видно, что значения данных функций пренебрежимо мало при значениях спектрального параметра, соответствующего углам распространения |0| > 30°. Действительно, точечный источник излучает волны равномерно по всем направлениям, но как видно на рисунке 4.4, отражение от препятствия происходит на углах |0| < 30°.

4.7.1 Сравнение с методом двунаправленного ПУ

В следующем примере рассматривается дифракция на двух препятствиях одинаковой высоты. Гауссова антенна с шириной диаграммы направленности в = 10° расположена на высоте 50 м над подстилающей поверхностью. Для моделирования направленного источника использован метод, предложенный в главе 3. Препятствия расположены на расстояниях 200 и 800 м от источника излучения. Высота препятствий равна 50 м. В данном примере проводится сравнение предложенного подхода с методом двунаправленного ПУ Использована реализация из работы [24]. Для двунаправленного ПУ использовались следующие расчетные параметры: шаг по продольной координате Ах = 1 м, шаг по высоте Аг = 0.5 м, порог сходимости £ = 0.05. Результаты моделирования, полученные предло-

150

100

со

I-

8 Л CQ

50

0

п

ы

-10 -15 -20 -25 -30 -35 -40

0

200

800

1000

400 600

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

Рисунок 4.4 — Рассеяние электромагнитных волн на одном препятствии. Пространственное распределение электромагнитного поля 10 log |"ф|.

phi

0

-0.6 -0.4 -0.2

0 p/k

0.2 0.4 0.6

Рисунок 4.5 — График функции ф(р) решения интегрального уравнения.

женным методом и методом двунаправленного ПУ изображены на рисунке 4.7. Распределение амплитуды поля на высоте 51м, рассчитанное обоими методами, изображено на рисунке 4.8. Легко заметить, что оба метода дают едва различимые результаты. Хорошо видны эффекты обратного рассеяния от препятствий и

РЭ

р/к

Рисунок 4.6 — График функции д(р) в правой части интегрального уравнения.

множественные переотражения. Отраженное поле от второго препятствия обеспечивает освещение зоны тени за первым препятствием.

Далее рассматривается дифракция на пяти препятствиях с разной высотой и взаимным расположением. Высоты препятствий равны соответственно 25, 50, 75, 50 и 125 м, считая слева направо. Они расположены соответственно на расстояниях 100,200, 300, 600, и 800 м от источника. Распределение поля, рассчитанное предложенным методом, показано на рисунке 4.9. Сравнение с двунаправленным ПУ приведено на рисунке 4.10. Как и в предыдущем примере, результаты, полученные обоими методами, фактически совпадают.

Очевидно, что вычислительное время обоих методов существенно возрастает с увеличением числа препятствий. Размер системы уравнений (4.13) зависит от числа препятствий т как те х те. Между тем, метод двунаправленного ПУ требует значительно большего числа итераций с увеличением числа препятствий, особенно когда препятствия расположены в непосредственной близости друг от друга. В таблицах 13 и 14 указано вычислительное время, затраченное на расчет рассмотренных выше примеров.

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

100

со i-

8 Л CQ

50

0

0

200

150

100

со i-о о л m

50

0

0

200

400 600

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

800

I

0 -5 -10 -15 -20 -25 J -30 -35 -40

1000

п

U

0 -5 -10 -15 -20 -25 J -30 -35 -40

800

1000

400 600

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

Рисунок 4.7 — Рассеяние электромагнитных волн на двух препятствиях, рассчитанное предложенным методом (сверху) и методом двунаправленного ПУ (снизу). Пространственное распределение электромагнитного поля 10 log |"ф|.

f=300 МГц.

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

z = 51 m

Двунаправленное ПУ Предложенный метод

LQ

ч:

со"

5

С <

0

100 200

300

400 500 600 700 800 900 100С Расстояние, м

Рисунок 4.8 — Рассеяние электромагнитных волн на двух препятствиях. Распределение электромагнитного поля 10 log \"ф\ на высоте 51м. f=300 МГц.

150

100

со

I-

о о л m

50

0

I

о

-5

-10

-15

-20

-25

-30

-35

-40

0

200

800

1000

400 600

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

Рисунок 4.9 — Рассеяние электромагнитных волн на пяти препятствиях, рассчитанное предложенным методом. Пространственное распределение электромагнитного поля 10 log \-ф\. f=300 МГц.

х=175 м

х=700 м

140

120

100

2

го 1- 80

о

о

.0 ш 60

40

20

0

Двунаправленное ПУ Предложенный метод

-30 -20 -10 Амплитуда, дБ

140 120 100 го 80

I-

о

Ъ 60 Ш

40 20 0

х

Двунаправленное ПУ Предложенный метод

-40 -30 -20

Амплитуда, дБ