CC BY
13
DOI: 10.25702/^С.2307-5252.2018.9.5.117-126 УДК 537.877+519.6
О.И. Ахметов, И.В. Мингалев, О.В. Мингалев, З.В. Суворова, С.М. Черняков, В.Б. Белаховский
ВЛИЯНИЕ D-ОБЛАСТИ ИОНОСФЕРЫ ВЫСОКИХ ШИРОТ НА РАСПРОСТРАНЕНИЕ ЭЛЕКТРОМАГНИТНЫХ СИГНАЛОВ С ЧАСТОТОЙ 82 ГЦ ВО ВРЕМЯ СОЛНЕЧНОЙ ВСПЫШКИ 15.03.2013
Аннотация
Представлены результаты численных экспериментов прохождения электромагнитных сигналов в волноводе Земля-ионосфера при различных состояниях ионосферы в области высоких широт. Профили концентрации, используемые в численных экспериментах, базируются на данных установки частичных отражений Полярного геофизического института, расположенной на радиофизическом полигоне «Туманный» Мурманской области (69.0 с.ш., 35.7 в.д.), во время и после солнечной вспышки М-класса 15.03.2013.
Ключевые слова:
СНЧ, ИНЧ, численное моделирование, ионосфера.
O. I. Akhmetov, I.V. Mingalev, O. V. Mingalev, Z.V. Suvorova, S.M. Chernyakov, V.B. Belakhovsky
D REGION INFLUENCE OF HIGH LATITUDE IONOSPHERE ON THE SPREAD OF ELECTROMAGNETIC SIGNALS OF FREQUENCY 82 AND 820 Hz DURING A SOLAR FLASH 15.03.2013
Abstract
The results of numerical experiments of the electromagnetic signals propagation in the Earth-ionosphere waveguide with various states of the ionosphere at the high latitudes atmosphere are presented. The concentration profiles used in numerical experiments are based on data from the partial reflections facility at the observatory "Tumanny" (69.0N, 35.7E) of the Polar Geophysical Institute during the M-class solar flash 15.03.2013
Keywords:
SLF, ULF, numerical modeling, ionosphere.
Введение
Наблюдения ионосферных эффектов солнечных вспышек имеют огромное значение для всех видов деятельности, связанной с распространением радиоволн в атмосфере Земли. Состоянием ионосферы определяются возможности связи в диапазоне СНЧ в большей мере, чем иными факторами. Исследования влияния ионосферы методами прямых измерений проводятся десятки лет в большинстве стран мира, в том числе и в области высоких широт. Очевидные выводы об особенностях ионосферного влияния на распространение радиоволн были уже сделаны авторами, первыми взявшимися за эту задачу, дальнейшее же продвижение требует применения нетривиальных методов как для регистрации, так и для обработки и интерпретации результатов. Все это относится
к области D, одной из наиболее сложных и менее всего изученных частей ионосферы. Одним из таких нетривиальных методов изучения эффектов влияния Б-области ионосферы на распространение радиоволн является численное моделирование. Данный метод позволяет получать чистые от постороннего шума эффекты, распределённые во всем интересующем исследователя пространстве, в отличие от локализированных прямых измерений каких-либо величин.
Рис. 2. Концентрация электронов по данным средневолнового радиолокатора вертикального излучения для исследования нижней ионосферы, расположенного на радиофизическом полигоне «Туманный» Полярного геофизического института в Мурманской области (69.0 с.ш., 35.7 в.д.), во время и после солнечной вспышки М-класса 15.03.2013.
В работе представлены результаты численных экспериментов прохождения электромагнитных сигналов в волноводе Земля-ионосфера при различных состояниях ионосферы в области высоких широт. Профили концентрации, используемые в численных экспериментах, получены при помощи доработки профилей, сгенерированных моделью Ш!2016. Доработка профилей базируется на данных средневолнового радиолокатора вертикального излучения для исследования нижней ионосферы, расположенного на радиофизическом полигоне «Туманный» Полярного геофизического института в Мурманской области (69.0 с.ш., 35.7 в.д.), во время и после солнечной вспышки М -класса 15.03.2013 [1], а именно в 08:00 иТ (возмущенные условия) и в 16:00 иТ (спокойные условия). Данные по концентрации электронов средневолнового радиолокатора за 15.03.2013 представлены на рис. 1.
Описание численных экспериментов
Авторами были проведены расчёты распространения радиоволн в волноводе Земля-ионосфера для четырех разных конфигураций ионосферных условий на частоте 82 Гц. Для возмущенной ионосферы использовался профиль, представленный на рис. 2а, соответствующий 8-часовым условиям 15.03.2013, профиль для невозмущенных условий, соответствующий 16-часовым условиям, представлен на рис. 2б.
Концентрация электронов м"3 Рис. 2. Аппроксимированные высотные профили концентрации электронов в 08:00(ИТ) (а) и 16:00(ИТ) (б); аппроксимированный профиль проводимости
литосферы (в).
Две конфигурации представляли собой однородную по горизонту в атмосфере область моделирования с вертикальным профилем концентрации, соответствующим возмущенным и невозмущенным условиям. Еще две конфигурации представляли собой неоднородную по горизонту в атмосфере область моделирования с переходом по функции Гаусса от невозмущенных условий над источником к возмущенным и наоборот. Источник во всех экспериментах представлял собой горизонтальный гармонический диполь
длиной 100 км. Литосфера во всех экспериментах выбрана горизонтально однородной с профилем проводимости, представленным на рис 2в; данный профиль построен на основе работ по исследованию проводимости грунта на Кольском полуострове [2]. Диэлектрическая проницаемость литосферы была принята равной 9.
Описание модели, граничных условий и функции источника
Используемая в представленной работе модель распространения электромагнитных сигналов в различных средах построена на основе схемы с противопотоковой аппроксимацией пространственных производных (метод Годунова с коррекцией потоков). Также используется расщепление по пространственным направлениям и по физическим процессам, причем затухание поля сигнала за счет проводимости и его вращение при наличии холловской проводимости среды учитываются на отдельных шагах расщепления по аналитическим формулам. Схема является монотонной, имеет 2-й порядок точности по времени и 3-й по пространственным переменным, а также является консервативной.
Рассмотрим кратко основные идеи, лежащие в основе численной схемы (более подробно см. в работах [3,4]).
(1)
где
и {вх,ву,в2,ех,еу,е2), Р=(Р1,Р2,Р3,Р4,Р5,РбУ, В = с0В, Е = фцхЕ, =]
Р\ = Рг = Рз = 0,
(РА,Р5,Р6)Т = -[Мх В]- ЧЕ-],
М = Ус +
С&ц)
0=
£0 и - электрическая и магнитная постоянные, £ и ^ - безразмерные относительные диэлектрическая и магнитная проницаемости среды, зависящие от пространственных координат.
Симметричные матрицы Ах, Ау, А2 определяются формулами:
^х =
' 0 ^) А =( 0
^-сПх 0 / у \-cRy
ИуС 0
00, =
0
-с0у
=
'0 0 0 0 ^0 1
-11 00у = 0
0 0 1 0 -1 0 0 0 0 ), 0х = (1 0 0 -1 0 0 0 0 0
Векторное уравнение (1) задает 6-мерную линейную гиперболическую систему уравнений 1-го порядка, записанную в консервативной форме. Ее правая часть ¥ линейно зависит от компонент вектора и. Для численного интегрирования
0
таких систем разработано достаточно много разностных схем, в том числе схемы повышенного порядка точности, которые применяются для уравнений газовой динамики. Наиболее полное описание этих схем содержится в монографиях [5, 6]. Используя метод расщепления по пространственным направлениям и физическим процессам [6, 7], можно построить явные монотонные схемы численного интегрирования системы (1), которые сводятся к последовательному интегрированию 1 -мерных по пространству гиперболических систем уравнений.
Для изотропной среды на каждом временном шаге нужно последовательно проинтегрировать 3 системы уравнений, например в таком порядке:
дг дхк х ■> ди д
~дГ ду
д и д
= Г,, (2)
д, +тЛА*и1 = г*
При этом правые части систем (2) выбираются так, чтобы они не изменялись на своем шаге расщепления и удовлетворяли равенству Р = Рх + Ру + Рг- На каждом из шагов расщепления для двух компонент магнитного поля и двух компонент электрического поля, ортогональных направлению шага, рассчитывается распространение сигнала конечно-разностным способом, а также рассчитывается по аналитическим формулам затухание за счет проводимости третьей компоненты электрического поля. В качестве начальных условий для каждой системы уравнений в (2) берутся значения, рассчитанные в результате предыдущего шага расщепления. Сохранить второй порядок аппроксимации по времени в схеме расщепления можно, если циклически изменять порядок выполнения шагов расщепления. Например, выполняя сначала в следующем порядке шаги по пространственным направлениям: xyz, yxz, zxy, xzy, yzx, zyx. Обоснование этого утверждения содержится, например, в монографиях [6,7].
В случае анизотропной среды с холловской проводимостью (такой средой является плазма в ионосфере и магнитосфере) тензор проводимости представляется в виде суммы его симметричной и антисимметричной частей. В этом случае к трем шагам расщепления схемы (2) добавляется четвертый шаг расщепления. Параметры среды и переменные, зависящие от проводимости и диэлектрической проницаемости, становятся тензорными. Компоненты вектора Б задаются формулами:
(РА, Р5, Р6)т = - (д) КхВ - (ду) %в - (£) ЪгВ -с[МхВ]-Г1Е-1
ч = е^Ч-1'2™-1'2, а = Щ12.
При этом на трех шагах расщепления (2) учитываются только симметричная часть тензора проводимости, а на четвертом шаге учитывается вращение электрического поля за счет антисимметричной части тензора проводимости, которое описывается аналитическими формулами. При этом магнитное поле не изменяется. На этом шаге в каждой точке расчетной сетки аналитически решается система уравнений:
§=[ПХЕ] (3)
здесь П = (/¡уг, 7]Хг, 7}Ху) - угловая скорость, в которой г)у2, Т]Х2, 7)ху являются компонентами антисимметричной части тензора ]. Эта система задает вращение поля Е с вектором угловой скорости П . Соответствующее циклическое изменение последовательности выполнения шагов расщепления обеспечивает второй порядок аппроксимации по времени и в этом случае.
Во всех поставленных в работе численных экспериментах область моделирования представляла собой параллелепипед с основанием 2048x6144 км, высотой в атмосфере 300 км и глубиной в литосфере 100 км. Шаги сетки по горизонтали составляли 8 км, по вертикали в атмосфере 2 км и 1 км в литосфере. Шаг по времени составлял 4-10"6 сек.
'шШНМШнШ
0.3 0.4 0.5
Время сек Время сек
Рис. 3. Пример временных вариаций электрического поля на расстоянии 520 км от источника на высоте 1 км (а); на высоте 30 км (б).
В модельных экспериментах на всех внешних границах действовало граничное условие свободного ухода волны, что достигалось обнулением исходящего за пределы области потока. Применяемая авторами схема в таких условиях имеет достаточно низкие коэффициенты отражения плоской волны от границ области моделирования: для волн, падающих под углом от 80 до 90 градусов, отношение амплитуды отраженной волны к амплитуде падающей плоской волны не превышает 0,01. При угле падения 60 градусов это отношение уже составляет примерно 0,05, при угле падения 45 градусов - примерно 0,16, при угле падения 27 градусов - примерно 0,33, а при угле падения 18,4 градуса -примерно 0,43 [3]. Метод FDTD (finite-differences time-domain method) [8] при применении простых граничных условий, таких как условия Мура (Mur) [9] и Лиао (Liao) [10], дает отражения порядка 0,1..1 %, но только при падении волны на границу под прямым углом. При падении под острым углом коэффициент отражения растет вплоть до 100 % при падении по касательной. Однако при использовании непрерывно действующего источника даже столь малых отражений, которые порождает применяемая схема, достаточно для накопления ошибок в области моделирования, и возникает необходимость в использовании методов подавления, подобных PML (perfectly matched layer), использующихся в FDTD-моделях [11]. Именно такой тип источника применялся в представленных
авторами экспериментах, что привело к необходимости адаптации и применения метода РМЬ. Разделение схемы по пространственным переменным и физическим процессам позволяет применять профиль электрических и магнитных потерь, предложенный Беренгером, непосредственно к потокам противопотоковой схемы на границе области моделирования. Геометрический профиль потерь внутри отдельного слоя имеет вид:
(4)
где д - коэффициент геометрической прогрессии,Дх - шаг по пространству, с0 -скорость света, N - номер PML-слоя, считая от интерфейса счетного региона и границы, г - расстояние от границы, Яо - коэффициент отражения от первого слоя. В представленных численных экспериментах авторы используют профиль потерь, рассчитанный по формуле 4 со следующими параметрами: Яо=0.01 (1 %), коэффициент прогрессии д = 2,15, количество слоев 14. Несмотря на то, что коэффициент отражения от первого слоя не лучше, чем характерный для данной схемы при обнулении исходящих потоков на углах падения 80-90 градусов, а на практике даже хуже вследствие отражений от последующих слоев, основным преимуществом метода РМЬ является его крайне слабая зависимость от угла прихода электромагнитной волны. Данную особенность демонстрирует и адаптированный для противопотоковой схемы вариант.
В качестве источника сигнала во всех представленных экспериментах используется длинная линия (100 км), по которой протекает изменяющийся по закону синуса ток с максимальным значением 100 А. Линия направлена вдоль оси х, а по оси у смещена таким образом, чтобы расстояние до трех боковых границ счетной области от геометрического центра антенны было одинаково. Высота линии над поверхностью Земли 15 м. Распределение поля в 2 ближайших к линии узлах в каждую сторону рассчитывалось аналитически для условий £ = ^ = 1, а = 0.
Рис. 4. Пример записи вариаций магнитной напряженности поля на расстоянии 520 км в спокойных условиях на высоте 30 км.
Результаты и обсуждение
На рис. 3 представлены примеры записи электрического поля на расстоянии 520 км от источника на высотах 1 и 30 км в спокойных условиях. Графики демонстрируют двукратное падение амплитуды поля с изменением высоты на 30 км при распространении в волноводе Земля-ионосфера. На рис. 4 представлен пример записи магнитного поля на той же дистанции на высоте 30 км. Не совсем стабильная синусоида в самом начале записей связана с началом счета и с тем, что аналитическое задание генератора сигнала не вполне соответствует условиям среды. Видно, что основными переносчиками поля являются компоненты Ех, Ну и в гораздо меньшей степени Нх.
Основной интерес в работе представляют результаты по преодолению сигналом областей с горизонтальным градиентом электронной концентрации, представленные на рис. 5. Для наглядности и удобства на рис. 5 представлены огибающие модуля вектора Пойнтинга и его У-компоненты (направлена в сторону распространения). Цветом показаны ионосферные условия: синий -спокойные, красный - возмущенные, зеленый - переход от спокойных условий в возмущенные, а черный - переход от возмущенных к спокойным. Видно, что происходит заметное изменение амплитуды при переходе через область горизонтальных градиентов концентрации электронов. Во всех случаях энергия сигнала после перехода меньше, чем в случае с невозмущенной ионосферой, что может свидетельствовать о возможном прерывании СДВ радиосвязи в более сложных, чем представленные в работе, ионосферных условиях в области Б-слоя.
Рис. 5. Зависимость огибающей модуля вектора Пойнтинга от дистанции до
источника (а); Зависимость огибающей модуля У-компоненты вектора Пойнтинга от дистанции до источника (б). Цветом показаны ионосферные условия: синий - спокойные, красный - возмущенные, зеленый - переход от спокойных условий в возмущенные, черный- переход от возмущенных к
спокойным.
Заключение
В работе показано влияние условий ионосферы на распространение радиоволн частотой 82 Гц в волноводе Земля-ионосфера. Авторы выявили, что при проходе волной участка замагниченного ионосферного волновода с более высокой D-областью в более низкую или наоборот никакого сжатия и последующего расширения волны не происходит, а происходит безвозвратная ускоренная диссипация энергии волны в областях с низкой D-областью, которая при последующем выходе в участок нормальной ионосферы не восстанавливается.
Благодарности. Работа выполнена при финансовой поддержке РНФ (проект 1877-10018).
Литература
1. Терещенко В.Д., Васильев Е.Б., Овчинников Н.А., Попов А.А., Средневолновый радиолокатор Полярного геофизического института для исследования нижней ионосферы, В кн.: Техника и методика геофизического эксперимента, Апатиты: Изд. КНЦ РАН. 2003. c. 37-46.
2. Korja T., Engels M., Zhamaletdinov A.A., Kovtun A. A., Palshin N.A., Smirnov M.Yu., Tokarev A.D., Asming V.E., Vanyan L.L., Vardaniants I.L., and the BEAR Working Group. Crustal conductivity in Fennoscandia—a compilation of a database on crustal conductance in the Fennoscandian shield // Earth Planets Space. 54. 2002. 535-558.
3. Мингалев И.В., Мингалев О.В., Ахметов О.И., Суворова З.В. Явная схема расщепления для уравнений Максвелла // Математическое моделирование. 2018. том 30. № 12. стр. 17-38.
4. Мингалев О.В., Мингалев И.В., Мельник М.Н., Ахметов О.И., Суворова З.В. Новый метод численного интегрирования системы Власова-Максвелла // Математическое моделирование. 2018. том 30. № 10. стр. 21-43.
5. Куликовский А.Г., Погорелов Н.В., Семенов А.Ю. Математические вопросы численного решения гиперболических систем уравнений. 2-е изд., испр. и доп. - М.: ФИЗМАТЛИТ. 2012. 656 с.
6. Бисикало Д.В., Жилкин А.Г., Боярчук А.А.Газодинамика тесных двойных звезд. Москва: ФИЗМАТЛИТ. 2013. 632 с.
7. Ковеня В. М., Яненко Н. Н. Метод расщепления в задачах газовой динамики. Новосибирск: Наука, 1981.
8. Yee Kane. Numerical solution of initial boundary value problems involving Maxwell's equations in isotropic media // IEEE Transactions on Antennas and Propagation. 1966. Vol. 14. PP. 302-307.
9. Mur G. Absorbing boundary conditions for the finite-difference approximation of the time domain electromagnetic field equations // IEEE Trans. Electromagnetic Compatibility. vol. 23. no. 4. pp. 277-382, Nov. 1981.
10. Liao Z. P., Wong H. L., Yang B. P., Yuan Y. F. A transmitting boundary for transient wave analyses. // ScientaSinica (series A), 17. 1984. S. 1063 - 1076.
11. Berenger J-P. A perfectly matched layer for the absorption of electromagnetic waves. Journal of Computational Physics Volume. 114. Issue 2. October 1994, P. 185-200.
Сведения об авторах Ахметов Олег Иршатович
к.ф.-м.н., н. с., Полярный геофизический институт, г. Апатиты E-mail: [email protected]
Мингалев Игорь Викторович
д.ф.-м.н., в. н. с., Полярный геофизический институт, г. Апатиты; E-mail: [email protected]
Мингалев Олег Викторович
к. ф.-м. н., зав. сектора, Полярный геофизический институт, г. Апатиты E-mail: [email protected]
Суворова Зоя Викторовна
программист, Полярный геофизический институт, г. Апатиты E-mail: [email protected]
Черняков Сергей Ыихайлович
н. с., Полярный геофизический институт, г. Мурманск Белаховский Владимир Борисович
к. ф.-м. н., н. с., Полярный геофизический институт, г. Апатиты E-mail: [email protected]
