Геомагнетизм и аэрономия, 2020, T. 60, № 4, стр. 528-541

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

В. В. Хегай *

Институт земного магнетизма, ионосферы и распространения радиоволн им. Н.В. Пушкова РАН (ИЗМИРАН)
г. Москва, г. Троицк, Россия

* E-mail: hegai@izmiran.ru

Поступила в редакцию 24.05.2019
После доработки 02.09.2019
Принята к публикации 23.01.2020

Полный текст (PDF)

Аннотация

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

1. ВВЕДЕНИЕ

В ионосфере, при обычных условиях в области хорошей погоды, величина поля, связанная с атмосферным электричеством, составляет ~1 мкВ/м, т.е. в 108 раз меньше, чем обычное поле у поверхности Земли в области хорошей погоды (~100 В/м). Однако, локализованные источники в области ниже верхней границы тропосферы (высота z ≤ ≤ 20 км), связанные с полями грозовых облаков (см. работы Park and Dejnakarintra [1973], Kim and Hegai [2015]), или крупномасштабные источники сейсмического происхождения (см. публикации Hegai et al. [2015], Denisenko et al. [2018]), могут, при определенных обстоятельствах, привести к величинам поля в ионосфере от 0.2 до 1 мВ/м. Такие поля незначимы для области тропосферы, но играют роль в перераспределении плазмы в ионосфере.

Тензор проводимости атмосферного воздуха $\hat {\sigma },$ в значительной степени определяющий величину проникшего в ионосферу поля от приземных источников, находится в сложной зависимости от высоты, а именно, от поведения заряженных частиц на разных высотах. В нижних слоях атмосферы за счет столкновения заряженных частиц с преобладающим количеством нейтральных проводимость изотропна, но с ростом высоты изотропность тензора $\hat {\sigma }$ нарушается. В более разряженной атмосфере, а тем более в ионосфере, количество заряженных частиц увеличивается по сравнению с нейтральной составляющей. Магнитное поле Земли приводит к упорядоченному дрейфу ионов и электронов по спирали вдоль силовых линий. Однако, до определенной высоты концентрация нейтральных компонент еще значительна, и влияние последних на дрейф ионов велико. Здесь возникают, помимо продольной проводимости σ0, так называемые проводимость Педерсена (Pedersen) σp ≡ σ1 и проводимость Холла (Hall) σh ≡ σ2. В общем случае, с учетом наклона геомагнитных силовых линий, в правой декартовой системе координат с осью z, направленной вверх, и вектором геомагнитной индукции B, лежащим в плоскости zy, тензор проводимости $\hat {\sigma }$ имеет вид (см., например, статью [Park and Dejnakarintra, 1973])

$\begin{gathered} \hat {\sigma } = \\ = \left[ {\begin{array}{*{20}{c}} {{{\sigma }_{1}}}&{{{\sigma }_{2}}{\text{sin}}I}&{{{\sigma }_{2}}{\text{cos}}I} \\ { - {{\sigma }_{2}}{\text{sin}}I}&{{{\sigma }_{1}}{{{\sin }}^{2}}I + {{\sigma }_{0}}{\text{co}}{{{\text{s}}}^{2}}I}&{({{\sigma }_{1}} - {{\sigma }_{0}}){\text{cos}}I{\text{sin}}I} \\ { - {{\sigma }_{2}}{\text{cos}}I}&{({{\sigma }_{1}} - {{\sigma }_{0}}){\text{cos}}I{\text{sin}}I}&{{{\sigma }_{1}}{{{\cos }}^{2}}I + {{\sigma }_{0}}{\text{si}}{{{\text{n}}}^{2}}I} \end{array}} \right], \\ \end{gathered} $
где I – угол наклона геомагнитных силовых линий, т.е. угол между вектором B и осью y.

На уровнях ионосферы проводимость вдоль геомагнитных силовых линий (σ| |) много больше, чем проводимость перпендикулярно им (σ), и, следовательно, вектор электрического поля в ионосфере перпендикулярен геомагнитным силовым линиям. При расчете проникновения электростатического поля от локализованных источников с поверхности земли в ионосферу (см., например, работы [Park and Dejnakarintra, 1973; Kim and Hegai, 2015; Hegai et al., 2015]) наклоном геомагнитных силовых линий часто пренебрегают (что оправдано на средних и в высоких широтах), т.е. работают в приближении вертикальных силовых линий геомагнитного поля (I = 90°), при этом тензор проводимости приобретает вид

$\hat {\sigma } = {\text{\;}}\left[ {\begin{array}{*{20}{c}} {{{\sigma }_{1}}}&{{{\sigma }_{2}}}&0 \\ { - {{\sigma }_{2}}}&{{{\sigma }_{1}}}&0 \\ 0&0&{{{\sigma }_{0}}} \end{array}} \right].$

В этом случае, в силу анизотропии проводимости, связанной с геомагнитным полем, если на уровне поверхности земли имеется только вертикальная компонента электрического поля (Ez) от аксиально-симметричного локализованного источника, то (в цилиндрической системе координат) на уровне ионосферы значимой будет только компонента поля, перпендикулярная геомагнитным силовым линиям (Er). Из работы [Denisenko et al., 2018] (см., в частности, формулу (34) этой работы) следует, что учет наклона магнитного поля приводит к тому, что величина компоненты Er в ионосфере будет ослабляться на фактор sin I по мере смещения источника поля на поверхности земли от более высоких широт к экватору, т.е. величину Er, рассчитанную в приближении вертикальных силовых линий геомагнитного поля (I = 90°), необходимо умножить на sin I для получения корректных ее значений.

Целью настоящего исследования является моделирование распределений сейсмогенного электрического поля по данным измерений в приземном слое атмосферы средних широт, которые приведены в работе [Hao et al., 2000] и расчет его величины в ионосфере аналогично тому, как это сделано в работе [Hegai et al., 2015]. При этом распределение возмущения вертикальной компоненты электрического поля ΔEz на поверхности земли (в предположении его цилиндрической симметрии относительно вертикальной оси, проходящий через эпицентр землетрясения) модифицируется согласно данным измерений Hao et al. [2000] и учитываются результаты работы Denisenko et al. [2018], указанные выше, т.е. расчетные величины Er в ионосфере корректируются в соответствии с наклоном геомагнитного поля в расчетной точке.

2. МОДЕЛИРОВАНИЕ СЕЙСМОГЕННОГО ЭЛЕКТРИЧЕСКОГО ПОЛЯ В АТМОСФЕРЕ У ПОВЕРХНОСТИ ЗЕМЛИ

В работе [Hegai et al., 2015] распределение возмущения вертикальной компоненты сейсмогенного квазистатического поля (ΔEz), связанного с приближающимся землетрясением, вблизи поверхности земли задавалось в цилиндрической системе координат (r, φ, z), начало которой находилось на горизонтальной плоскости, расположенной на высоте z = 0 км (ось z направлена вертикально вверх). Это распределение имело следующий вид:

(1)
$\Delta {{E}_{z}} = {{E}_{0}}{\text{exp}}{\kern 1pt} [--{{({r \mathord{\left/ {\vphantom {r {{{R}_{0}}}}} \right. \kern-0em} {{{R}_{0}}}})}^{2}}{\text{ln}}10],$
где E0 и R0 – соответственно максимальное значение Ez и радиус области возмущения поля, на котором исходное поле падает в 10 раз. При этом величина R0 = RD, т.е. характерный масштаб изменения поля для различных магнитуд сейсмических событий (землетрясений), полагался равным характерному радиусу зоны подготовки землетрясения на поверхности земли RD, определенному теоретически в работе [Dobrovolsky et al., 1979] следующим образом:
(2)
${{R}_{D}} = {{10}^{{0.43}}}^{M},$
где M – магнитуда землетрясения, а RD выражается в км. С другой стороны, согласно работе [Hao et al., 2000], результаты измерений которой далее используются при моделировании сейсмогенных аномалий электрического поля на поверхности земли, предшествующих землетрясениям, соответствующий характерный размер зоны (RH) эмпирически определен как
(3)
${{R}_{H}} \cong {\text{exp}}(M),$
где M – также величина магнитуды землетрясения, и RH тоже выражается в км.

Нетрудно заметить, что RD = 100.43Мe0.99M ≅ ≅ exp(M), так что в дальнейшем будем полагать, что радиус зоны подготовки землетрясения на поверхности земли задается выражением (3) и далее везде R0 = exp(M).

Из рассмотрения Table 4 работы [Hao et al., 2000] видно, что сейсмогенные аномалии вертикальной компоненты квазистатического электрического поля (ΔEz) в приземном слое атмосферы наблюдались одновременно на двух станциях Baijiatuan и Baodi за 4 (No. 4, M = 4.1 ), 5 (No. 3, М = 4.7) и 15 (No. 7, М = 4.3) дней до землетрясения соответствующей магнитуды. Эпицентральные расстояния до станций наблюдения (r) лежали в диапазоне от 28 до 95 км, а абсолютные величины аномальных значений Ez составляли от 344 В/м до 721 В/м. Отрицательный знак аномалий соответствует направлению вектора напряженности квазистатического сейсмогенного электрического поля “вверх” от поверхности земли.

На рисунке 1 представлены аналитические аппроксимации (сплошные линии) абсолютных величин этих данных (кривые No. 3, No. 4 и No. 7) по измерениям аномальных значений ΔEz (точки) в зависимости от эпицентрального расстояния r согласно соотношению

(4)
$\begin{gathered} \Delta {{E}_{z}} = {{E}_{0}}{\text{exp}}\left\{ {--([1 + } \right. \\ + \,\,{{{{{({r \mathord{\left/ {\vphantom {r {{{R}_{0}}}}} \right. \kern-0em} {{{R}_{0}}}})}}^{2}}{{]}^{{{1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-0em} 2}}}}--1)} \mathord{\left/ {\vphantom {{{{{({r \mathord{\left/ {\vphantom {r {{{R}_{0}}}}} \right. \kern-0em} {{{R}_{0}}}})}}^{2}}{{]}^{{{1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-0em} 2}}}}--1)} {{{({{{{R}_{E}}} \mathord{\left/ {\vphantom {{{{R}_{E}}} {{{R}_{0}}}}} \right. \kern-0em} {{{R}_{0}}}})\} } \mathord{\left/ {\vphantom {{({{{{R}_{E}}} \mathord{\left/ {\vphantom {{{{R}_{E}}} {{{R}_{0}}}}} \right. \kern-0em} {{{R}_{0}}}})\} } {{{{[1 + {{{({r \mathord{\left/ {\vphantom {r {{{R}_{0}}}}} \right. \kern-0em} {{{R}_{0}}}})}}^{2}}]}}^{{{1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-0em} 2}}}}}}} \right. \kern-0em} {{{{[1 + {{{({r \mathord{\left/ {\vphantom {r {{{R}_{0}}}}} \right. \kern-0em} {{{R}_{0}}}})}}^{2}}]}}^{{{1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-0em} 2}}}}}}}}} \right. \kern-0em} {{{({{{{R}_{E}}} \mathord{\left/ {\vphantom {{{{R}_{E}}} {{{R}_{0}}}}} \right. \kern-0em} {{{R}_{0}}}})\} } \mathord{\left/ {\vphantom {{({{{{R}_{E}}} \mathord{\left/ {\vphantom {{{{R}_{E}}} {{{R}_{0}}}}} \right. \kern-0em} {{{R}_{0}}}})\} } {{{{[1 + {{{({r \mathord{\left/ {\vphantom {r {{{R}_{0}}}}} \right. \kern-0em} {{{R}_{0}}}})}}^{2}}]}}^{{{1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-0em} 2}}}}}}} \right. \kern-0em} {{{{[1 + {{{({r \mathord{\left/ {\vphantom {r {{{R}_{0}}}}} \right. \kern-0em} {{{R}_{0}}}})}}^{2}}]}}^{{{1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-0em} 2}}}}}}}}, \\ \end{gathered} $
где E0 и R0 – соответственно максимальное значение ΔEz и радиус зоны подготовки землетрясения на поверхности земли, а RE определяет масштаб уменьшения поля по мере удаления от эпицентра. В случае, когда имеются одновременные данные измерений ΔEz в двух пространственно разнесенных точках и R0 задано (R0 = RH), из (4) можно определить E0 и RE.

Рис. 1.

Аналитические аппроксимации абсолютных величин сейсмогенных аномалий ΔEz в приземном слое атмосферы в зависимости от эпицентрального расстояния r (сплошные линии) экспериментальных данных (точки) работы [Hao et al., 2000] (Table 4 , No. 3, 4, 7) в соответствии с формулой (4). Номера кривых указаны рядом с ними.

Из рисунка видно, что аналитическое выражение (4) позволяет хорошо описывать измерения ΔEz, выполненные одновременно в двух разнесенных точках, при этом однозначно определяются величины E0 и RE. Для представленных кривых описывающие их параметры имеют следующие значения:

$\begin{gathered} {\text{No}}{\text{. }}3{\text{ (}}М = 4.7{\text{)}} \Rightarrow \\ \Rightarrow [{{E}_{0}} = 890{\text{ }}{{\text{В}} \mathord{\left/ {\vphantom {{\text{В}} {\text{м}}}} \right. \kern-0em} {\text{м}}},\,\,{{R}_{0}} = 110\,\,{\text{км,}}\,\,{{R}_{E}} = 66{\text{ км}}], \\ \end{gathered} $
$\begin{gathered} {\text{No}}{\text{. }}4{\text{ (}}М = 4.1{\text{)}} \Rightarrow \\ \Rightarrow [{{E}_{0}} = 1265{\text{ }}{{\text{В}} \mathord{\left/ {\vphantom {{\text{В}} {\text{м}}}} \right. \kern-0em} {\text{м}}},\,\,{{R}_{0}} = 60\,\,{\text{км,}}\,\,{{R}_{E}} = 40{\text{ км}}], \\ \end{gathered} $
$\begin{gathered} {\text{No}}{\text{. }}7{\text{ (}}М = 4.3{\text{)}} \Rightarrow \\ \Rightarrow [{{E}_{0}} = 802\,\,{{\text{В}} \mathord{\left/ {\vphantom {{\text{В}} {\text{м}}}} \right. \kern-0em} {\text{м}}},\,\,{{R}_{0}} = 74\,\,{\text{км,}}\,\,{{R}_{E}} = 130\,\,{\text{км}}]. \\ \end{gathered} $

Эти аппроксимации дают количественное представление о возможных соотношениях между размерами радиуса зоны подготовки землетрясения на поверхности земли R0 и характерным масштабом изменения поля RE, а также о предполагаемой аномальной величине ΔEz в эпицентре назревающего землетрясения на поверхности земли. Ясно, что получить фактические измерения ΔEz точно в эпицентре будущего землетрясения можно только чисто случайно. Из полученных данных следует, что величина отношения RE/R0 лежит в интервале [0.6; 1.75], при этом средняя величина 〈E0〉 ≈ 1 кВ/м уже для диапазона магнитуд землетрясений М ∈ [4.1; 4.7].

На рисунке 2 дано сравнение кривой No. 3 (сплошная линия, для полученного значения E0 = = 890 В/м) с аппроксимацией сейсмогенного поля (1), использованного в статье [Hegai et al., 2015] (штрихпунктирная линия) и, аналогично, данной в работе [Denisenko et al., 2013] (штриховая линия, см. также формулу (12) для E2(y) работы [Denisenko et al., 2013]) в соответствии с выражением

(5)
$\Delta {{E}_{z}} = {{E}_{0}}{{[1 + {\text{cos}}({{\pi r} \mathord{\left/ {\vphantom {{\pi r} {{{R}_{0}}}}} \right. \kern-0em} {{{R}_{0}}}})]} \mathord{\left/ {\vphantom {{[1 + {\text{cos}}({{\pi r} \mathord{\left/ {\vphantom {{\pi r} {{{R}_{0}}}}} \right. \kern-0em} {{{R}_{0}}}})]} {2,}}} \right. \kern-0em} {2,}}$
при одинаковых значениях R0 = 110 км во всех трех случаях.

Рис. 2.

Сравнение аналитических аппроксимаций сейсмогенных аномалий ΔEz в приземном слое атмосферы: формула (4) – сплошная линия, соответствующая кривой No. 3 на рис. 1; формула (1) – штрихпунктирная линия; формула (5) – штриховая линия. Вертикальная стрелка отмечает положение r = R0 = 110 км, а косые крестики на ней указывают величину ΔEz в этой точке на соответствующих кривых. Остальные обозначения те же, что и на рис. 1.

Из рисунка следует, что лучшая аппроксимация экспериментальных данных исследования [Hao et al., 2000] соответствует предложенной здесь формуле (4), при этом на характерном расстоянии r = R0 = 110 км формула (1) занижает значение ΔEz в 3.55 раза по сравнению с формулой (4), а формула (5) дает нулевую величину поля. Как будет показано далее, этот факт имеет существенное значение при расчетах величины электрического поля, проникающего от поверхности земли на высоты ионосферы.

Так как проникновение сейсмогенного поля в ионосферу наиболее эффективно ночью в периоды низкой солнечной активности (см. работы [Hegai et al., 2015; Denisenko et al., 2018]), то в последующих расчетах его проникновения из указанных выше случаев будет рассмотрен только случай сейсмогенной аномалии ΔEz, наблюдавшейся за 5 дней до землетрясения 10.11.1986 г. с М = 4.7 (географические координаты эпицентра φ = 40.05° N, λ = 158.72° E), соответствующий аппроксимирующей кривой No. 3 на рис. 1. В период его подготовки усредненный суточный поток солнечного излучения F10.7 с 4 по 10 ноября лежал в диапазоне от 70 до 82 единиц SFU, тогда как в других случаях (кривые No. 4 и 7) его величина колеблется в диапазоне от 97 до 280 единиц SFU. Кроме того, из рассмотренных выше случаев это землетрясение имеет наибольшую магнитуду.

Далее, из всех рассмотренных в работе [Hao et al., 2000] сейсмогенных аномалий ΔEz, предваряющих землетрясения, всего лишь 3 из них (см. статью [Hao et al., 2000], Table 1 ) No. 17, No. 25 и No. 42 наблюдаются перед землетрясениями с М ≥ 6.0. При этом только для No. 42 (землетрясение 03.05.1996 г., М = 6.0, географические координаты эпицентра φ = 40.8° N, λ = 109.7° E, величина сейсмогенной аномалии ΔEz = – 633 В/м за 5 дней до толчка) солнечная активность находится на низком уровне (F10.7 ∈ [70; 80] SFU). В этом случае имеются наблюдения аномально высоких по абсолютной величине значений ΔEz только на ОДНОЙ станции – Baotou, лежащей на эпицентральном расстоянии Re = 470 км. Отметим здесь, что сейсмогенная аномалия ΔEz наблюдается на большем расстоянии, чем величина радиуса зоны подготовки на поверхности земли R0(М = 6.0) ≅ ≅ 400 км.

Аппроксимации возможного распределения величины сейсмогенного поля Ez в зависимости от расстояния до эпицентра r перед землетрясением 03.05.1996 г. c М = 6.0 согласно (4) даны на рис. 3 для двух возможных значений отношения RE/R0. Там же представлены числовые значения параметров соответствующих кривых для формулы (4), а опорное измеренное значение ΔEz (r = = 470 км) = 633 В/м с учетом знака поля в выбранной цилиндрической системе координат отмечено косым крестиком. В соответствии приведенными ранее возможными величинами отношения RE/R0, взяты два случая: RE/R0 = 0.6 (сплошная линия) и RE/R0 = 1.75 (штриховая линия).

Рис. 3.

Аппроксимации возможного распределения величины возмущения сейсмогенного поля ΔEz в зависимости от расстояния до эпицентра r перед землетрясением 03.05.1996 г. c М = 6.0 согласно формуле (4) для разных величин отношения характерных масштабов RE и R0 по опорному измеренному значению ΔEz (r = 470 км) = 633 В/м (отмечено косым крестиком) с учетом знака поля в выбранной цилиндрической системе координат: RE/R0 = 0.6 (сплошная линия) и RE/R0 = 1.75 (штриховая линия).

Из рассмотрения рис. 3 хорошо видно, что увеличение RE – характерного масштаба изменения величины сейсмогенного поля ΔEz относительно R0 приводит к уменьшению максимального значения ΔEz, достигаемого в эпицентре землетрясения (E0), в случае, если принимаемая модель аппроксимации по формуле (4) справедлива.

3. МОДЕЛИ ПРОВОДИМОСТИ

Для расчета проникновения квазиэлектростатического поля в ионосферу для двух выбранных случаев сейсмогенных аномалий ΔEz у поверхности земли за 5 дней до землетрясения (перед землетрясениями EQ-1{10.11.1986 г. с М = 4.7} и EQ-2{03.05.1996 г. c М = 6.0} соответственно) необходимо задать конкретные модели проводимости в диапазоне высот от 0 до 90 км для σ0(z) и σ1(z) ≡ σp(z) в предположении горизонтальной стратификации среды, а также и величину интегральной проводимости Педерсена (ΣP), максимально близко к географическим координатам эпицентров рассматриваемых землетрясений и времени появления описанных сейсмогенных аномалий электростатического поля. В проведенных далее расчетах полагается, так же как и в статье [Denisenko et al., 2018] (см. формулу (16) этой работы для ΣP при I = 90°),что интегральная проводимость Педерсена определяется формулой

(6)
${{{\Sigma }}_{P}} = \int\limits_{{{z}_{{up}}}}^{{{z}_{\infty }}} {{{\sigma }_{p}}dz} ,$
при этом zup = 90 км, а z = 1000 км. Как указано в работе [Denisenko et al., 2018], величина проводимости ΣP практически не изменяется при выборе z = 400 или 500 км, что связано с очень быстрым падением σp(z) с увеличением высоты z, начиная с уровня 400 км. По другим оценкам (см. статью Hegai et al. [2015]), σp(z) становится пренебрежимо малой величиной и не вносит значимого вклада в интегральную проводимость ΣP выше уровня z = 600 км.

Ниже 70 км, в зоне изотропной проводимости, мы принимаем модельный вертикальный профиль ночной проводимости на средних широтах, приведенный на рис. 4.

Рис. 4.

Модельный кусочно-непрерывный вертикальный профиль ночной проводимости на средних широтах в области изотропной проводимости от 0 до 70 км. Числа у соответствующих отрезков прямых указывают характерный масштаб (в км) изменения экспоненциально растущей проводимости. Точкой отмечено значение проводимости ночного среднеширотного профиля проводимости из работы [Park and Dejnakarintra, 1976] на уровне 70 км, а затемненный прямоугольник под осью абсцисс дает разброс значений проводимости у поверхности земли в области хорошей погоды согласно Baumgaertner et al. [2013].

Этот профиль разделен на три высотных зоны: z ∈ [0; 20] км – зона 1, z ∈ [20; 60] км – зона 2 и z ∈ [60; 70] км – зона 3. При этом σ01 = σp1 = = b1 exp(z/h1), σ02 = σp2 = b2 exp[(z z1)/h2], σ03 = = σp3 = b3 exp[(z z2)/h3], а индексы 1, 2 и 3 относятся к значениям параметров проводимости в нижней, средней и верхней зонах изотропной проводимости. Интервал высот [0, z1 = 20] определяет зону, в которой параметры (b1, h1) таковы, что вертикальный профиль проводимости соответствует значениям средней глобальной модели проводимости работы [Baumgaertner et al., 2013] в области хорошей погоды. Во второй [z1 = 20; z2 = 60] и третьей [z2 = 60; z3 = 70] зонах используемый далее вертикальный профиль проводимости близок к модели, приведенной в работе [Rycroft et al., 2007]. Значения соответствующих параметров локальной проводимости в изотропной зоне следующие: h1 = 5.87 км, h2 = 9.65 км, h3 = 3.36 км; b1 = 1 × 10–13 См/м, b2 = b1 exp(z1/h1) См/м, b3 = = b2[exp(z2 z1)/h2] См/м.

В четвертой [z3 = 70; z4 = 90], анизотропной зоне, σ0(z) и σ1(z) ≡ σp(z) экстраполируются экспоненциально от своих значений на z3 = 70 км до своих полуночных по местному времени (00 LT) значений на уровне z4 = 90 км (в соответствии с географическими координатами эпицентра рассматриваемого землетрясения и датой наблюдения предшествующей ему сейсмогенной аномалии ΔEz у поверхности земли).

В диапазоне высот 90 ≤ z ≤ 1000 км величины локальных проводимостей определяются следующими соотношениями

(7)
${{\sigma }_{0}} = {{e}^{2}}\left[ {\frac{{{{N}_{e}}}}{{{{m}_{{e~}}}{{\nu }_{e}}}} + \sum\limits_i {\frac{{{{N}_{i}}}}{{{{m}_{{i~}}}{{\nu }_{i}}}}} } \right],$
(8)
${{\sigma }_{p}} = {{e}^{2}}\left[ {\frac{{{{N}_{e}}{{\nu }_{e}}}}{{{{m}_{{e~}}}(\omega _{e}^{2} + \nu _{e}^{2})}} + \sum\limits_i {\frac{{{{N}_{i}}{{\nu }_{i}}}}{{{{m}_{{i~}}}(\omega _{i}^{2} + \nu _{i}^{2})}}} } \right],$
где $\sum\nolimits_i {} $ – это знак суммирования по индексу i; индексы e и i обозначают электроны и i-ю ионную составляющую ионосферной плазмы; Ne и Ni – концентрации электронов и ионов соответственно; e – заряд электрона; me и mi – электронная и ионная массы; νe и νi – частоты столкновений с нейтральными частицами плазмы, а ωe и ωi – электронная и ионная гирочастоты. Частоты νe и νi вычисляются в соответствии с работой [Schunk, 1988]. Все необходимые входные параметры вычисляются по эмпирической модели ионосферы IRI-2016 (https://ccmc.gsfc.nasa.gov/modelweb/models/iri2016_vitmo.php) и модели нейтральной атмосферы NRLMSIS-00 (https://modelweb.gsfc.nasa. gov/atmos/nrlmsise00.html). Нужные для дальнейших расчетов значения интегральной проводимости $\sum _{P}^{{{\text{av}}}}$ = (ΣP + $\sum _{P}^{*}$)/2 получаются в соответствии с формулой (6). Здесь ΣP – интегральная проводимость Педерсена для географических координат эпицентра рассматриваемого землетрясения в момент 00 LT соответствующей даты, а $\sum _{P}^{*}$ – интегральная проводимость Педерсена для магнито-сопряженной с ним точки.

На рисунках 5 и 6 приведены модельные высотные профили проводимости от 0 до 90 км, принятые при расчетах проникновения электрических полей сейсмогенных аномалий ΔEz у поверхности земли, связанных с землетрясениями EQ-1{10.11.1986 г. с М = 4.7} и EQ-2{03.05.1996 г. c М = 6.0}, проявившихся 05.11.1986 г. и 28.04.1996 г. соответственно. Числовые значения $\sum _{P}^{{{\text{av}}}}$ также указаны на рисунках.

Рис. 5.

Модельные кусочно-непрерывные вертикальные профили ночных проводимостей σ0 и σ1 от 0 до 90 км в точке географических координат эпицентра землетрясения EQ-1{10.11.1986 г. с М = 4.7} для 00:00 LT 05 ноября1986 г. Числа у соответствующих отрезков прямых указывают характерный масштаб (в км) изменения экспоненциально растущей проводимости, а штриховые горизонтальные линии отмечают границы высотных зон. Среднее значение интегральной проводимости Педерсена для двух полушарий $\sum _{P}^{{{\text{av}}}}$ = 0.156 См также указано на рисунке.

Рис. 6.

То же, что и рис. 5, но для географических координат эпицентра землетрясения EQ-2{03.05.1996 г. c М = 6.0} для 00:00 LT 28 апреля 1996 г. Среднее значение интегральной проводимости Педерсена для двух полушарий $\sum _{P}^{{{\text{av}}}}$ = = 0.0738 См.

Сравнение рисунков 5 и 6 показывает, что в верхней (четвертой зоне) профили проводимости различаются незначительно, зато $\sum _{P}^{{{\text{av}}}}$ в ноябре 1986 г. (рис. 5) более чем в два раза превосходит свои апрельские значения 1996 г. (рис. 6). Это, как будет видно из дальнейшего рассмотрения, приводит к большим различиям в величинах проникающего в ионосферу поля сейсмогенных электрических аномалий.

4. ОСНОВНЫЕ УРАВНЕНИЯ И ПОСТАНОВКА ЗАДАЧИ

Для расчета проникновения сейсмогенного электростатического поля (СЭП) в ионосферу воспользуемся формализмом, развитым в работе [Park and Dejnakarintra, 1973] при исследовании проникновения в ионосферу электростатического поля грозового облака и использованным также в исследовании [Hegai et al., 2015]. В стационарных условиях исходные уравнения, определяющие распределение электрического потенциала Φ на высотах z ≥ 0 км, имеют вид

(9)
$\nabla \,\,\cdot\,\,{\mathbf{J}} = 0,$
(10)
${\mathbf{J}} = \hat {\sigma }{\mathbf{E}},$
(11)
${\mathbf{E}} = --\nabla \Phi ,$
где J – плотность электрического тока; $\hat {\sigma }$ – тензор электрической проводимости среды, а E и Φ – напряженность электростатического поля и его потенциал соответственно.

Пренебрегая кривизной Земли и используя (как и в работе [Hegai et al., 2015]) цилиндрическую систему координат (r, φ, z), начало которой помещено в эпицентр последовавшего землетрясения, а ось z направлена вверх, представим сейсмогенное возмущение вертикальной компоненты приземного атмосферного электрического поля в соответствии с новой аппроксимирующей формулой (4), предложенной в настоящей работе (см. выше),

(12)
$\begin{gathered} \Delta {{E}_{z}} = {{\left. {--({{\partial \Phi } \mathord{\left/ {\vphantom {{\partial \Phi } {\partial z}}} \right. \kern-0em} {\partial z}})} \right|}_{{z = 0}}} = {{E}_{0}}{\text{exp\{ }}--{\kern 1pt} ({{[1 + {{({r \mathord{\left/ {\vphantom {r {{{R}_{0}}}}} \right. \kern-0em} {{{R}_{0}}}})}^{2}}]}^{{{1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-0em} 2}}}} - \\ {{--\,\,1)} \mathord{\left/ {\vphantom {{--\,\,1)} {{{({{{{R}_{E}}} \mathord{\left/ {\vphantom {{{{R}_{E}}} {{{R}_{0}}}}} \right. \kern-0em} {{{R}_{0}}}})\} } \mathord{\left/ {\vphantom {{({{{{R}_{E}}} \mathord{\left/ {\vphantom {{{{R}_{E}}} {{{R}_{0}}}}} \right. \kern-0em} {{{R}_{0}}}})\} } {{{{[1 + {{{({r \mathord{\left/ {\vphantom {r {{{R}_{0}}}}} \right. \kern-0em} {{{R}_{0}}}})}}^{2}}]}}^{{{1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-0em} 2}}}}.}}} \right. \kern-0em} {{{{[1 + {{{({r \mathord{\left/ {\vphantom {r {{{R}_{0}}}}} \right. \kern-0em} {{{R}_{0}}}})}}^{2}}]}}^{{{1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-0em} 2}}}}.}}}}} \right. \kern-0em} {{{({{{{R}_{E}}} \mathord{\left/ {\vphantom {{{{R}_{E}}} {{{R}_{0}}}}} \right. \kern-0em} {{{R}_{0}}}})\} } \mathord{\left/ {\vphantom {{({{{{R}_{E}}} \mathord{\left/ {\vphantom {{{{R}_{E}}} {{{R}_{0}}}}} \right. \kern-0em} {{{R}_{0}}}})\} } {{{{[1 + {{{({r \mathord{\left/ {\vphantom {r {{{R}_{0}}}}} \right. \kern-0em} {{{R}_{0}}}})}}^{2}}]}}^{{{1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-0em} 2}}}}.}}} \right. \kern-0em} {{{{[1 + {{{({r \mathord{\left/ {\vphantom {r {{{R}_{0}}}}} \right. \kern-0em} {{{R}_{0}}}})}}^{2}}]}}^{{{1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-0em} 2}}}}.}}}} \\ \end{gathered} $
Если предположить, что геомагнитное поле B является вертикальным и тензор электрической проводимости $\hat {\sigma }$ зависит только от z, то с учетом азимутальной симметрии задачи из уравнений (9)–(11) можно получить следующее уравнение для потенциала Φ:
(13)
$\begin{gathered} {{{{\partial }^{2}}\Phi } \mathord{\left/ {\vphantom {{{{\partial }^{2}}\Phi } {\partial {{r}^{2}}}}} \right. \kern-0em} {\partial {{r}^{2}}}} + ({1 \mathord{\left/ {\vphantom {1 r}} \right. \kern-0em} r}){{\partial \Phi } \mathord{\left/ {\vphantom {{\partial \Phi } {\partial r}}} \right. \kern-0em} {\partial r}} + \\ + \,\,({1 \mathord{\left/ {\vphantom {1 {{{\sigma }_{p}}}}} \right. \kern-0em} {{{\sigma }_{p}}}}){{\partial ({{{{\sigma }_{0}}\partial \Phi } \mathord{\left/ {\vphantom {{{{\sigma }_{0}}\partial \Phi } {\partial z}}} \right. \kern-0em} {\partial z}})} \mathord{\left/ {\vphantom {{\partial ({{{{\sigma }_{0}}\partial \Phi } \mathord{\left/ {\vphantom {{{{\sigma }_{0}}\partial \Phi } {\partial z}}} \right. \kern-0em} {\partial z}})} {\partial z}}} \right. \kern-0em} {\partial z}} = 0. \\ \end{gathered} $
Напомним здесь, что σp ≡ σ1 – локальная проводимость Педерсена, а σ0 – продольная (вдоль B) проводимость.

Уравнение (13) может быть решено аналитически, если проводимости σp и σ0 зависят от высоты z экспоненциально. Тогда, в области изотропии проводимости, σ0 = σp = b exp(z/h) (где b и h – константы), получим:

(14)
$\begin{gathered} \Phi = \int\limits_0^\infty {{{J}_{0}}(kr)[{{A}_{{\,1}}}(k)\exp ({{c}_{1}}z) + } \\ + \,\,{{B}_{1}}(k)\exp ({{c}_{2}}z)]dk, \\ \end{gathered} $
где J0 является функцией Бесселя первого рода нулевого порядка; c1= – l/(2h) – [l/(4h2) + k2]1/2; c2 = – l/(2h) + [l/(4h2) + k2]1/2, а A1 и B1 – коэффициенты.

Для анизотропной зоны, где σ0 = b0 exp(z/h0), а σp = bp exp(z/hp), решение уравнения (13) дается интегральной формулой

(15)
$\begin{gathered} \Phi = \int\limits_0^\infty {{{J}_{0}}(kr)[{{A}_{2}}(k){{I}_{\nu }}(kf) + } \\ + \,\,{{B}_{2}}(k){{K}_{\nu }}(kf)]{{f}^{\nu }}dk, \\ \end{gathered} $
где Jν и Kν – модифицированные функции Бесселя первого и второго рода порядка ν соответственно, при этом ν = hp/(hph0); f = = 2νh0(bp/b0)1/2 exp[–z/(2νh0)], а A2 и B2 – коэффициенты. Коэффициенты A1, B1, A2, и B2 могут быть определены из граничных условий, а компоненты электрического поля получаются из соотношений
(16)
${{E}_{r}} = {{--\partial \Phi } \mathord{\left/ {\vphantom {{--\partial \Phi } {\partial r}}} \right. \kern-0em} {\partial r}},$
(17)
${{E}_{z}} = {{--\partial \Phi } \mathord{\left/ {\vphantom {{--\partial \Phi } {\partial z}}} \right. \kern-0em} {\partial z}}.$
Как указано выше, мы полагаем, что вектор геомагнитного поля B параллелен оси z, откуда следует, что Ez | | B, в то время как ErB.

Выше 90 км продольная (вдоль геомагнитного поля) проводимость σ0 достаточно высока и намного больше поперечной проводимости σp, поэтому силовые линии геомагнитного поля являются почти эквипотенциальными линиями для случая перпендикулярных электростатических полей с характерными масштабами более нескольких десятков километров. Это позволяет рассматривать ионосферную область от ~90 до ~1000 км как тонкий проводящий слой с интегральной (вдоль геомагнитной силовой линии) проводимостью Педерсена ΣP . Тут следует напомнить, что локальная проводимость σp пренебрежимо мала выше уровня 600 км. Таким образом, уравнение непрерывности электрического тока на уровне z = 90 км можно записать в следующем виде:

(18)
${{\sigma }_{0}}{{E}_{z}} = {{\nabla }_{ \bot }}\left( {2\sum _{P}^{{{\text{av}}}} \cdot \,\,{{{\mathbf{E}}}_{ \bot }}} \right),$
где символ ${{\nabla }_{ \bot }}$ обозначает двумерный оператор градиента (в плоскости, перпендикулярной B), фактор 2 перед $\sum _{P}^{{{\text{av}}}}$ учитывает вклад в интегральную проводимость Педерсена магнито-сопряженной области ионосферы.

Здесь необходимо отметить, что соотношения, аналогичные (18), ранее уже использовались как верхнее граничное условие при решении задачи проникновения СЭП в ионосферу в работах [Denisenko et al., 2008; Ampferer et al., 2010; Hegai et al., 2015; Denisenko et al., 2018]. В выбранной (см. выше) цилиндрической системе координат уравнение (18) может быть записано в явном виде как

(19)
${{{{\sigma }_{0}}\partial \Phi } \mathord{\left/ {\vphantom {{{{\sigma }_{0}}\partial \Phi } {\partial z}}} \right. \kern-0em} {\partial z}} = 2\sum _{P}^{{{\text{av}}}}[{{{{\partial }^{2}}\Phi } \mathord{\left/ {\vphantom {{{{\partial }^{2}}\Phi } {\partial {{r}^{2}}}}} \right. \kern-0em} {\partial {{r}^{2}}}} + ({1 \mathord{\left/ {\vphantom {1 r}} \right. \kern-0em} r}){{\partial \Phi } \mathord{\left/ {\vphantom {{\partial \Phi } {\partial r}}} \right. \kern-0em} {\partial r}}].$
Соотношения (12) и (19) представляют собой нижнее и верхнее граничные условия соответственно для расчета электростатического потенциала Φ в соответствии с (14), (15).

Итак, для того чтобы рассчитать распределение электростатического потенциала Φ, однозначно определяющего величину СЭП соотношениями (16) и (17), накладываются следующие граничные условия:

$\begin{gathered} 1.\,\,{{--\partial \Phi } \mathord{\left/ {\vphantom {{--\partial \Phi } {\partial z}}} \right. \kern-0em} {\partial z}} = {{E}_{0}}{\text{exp}}\{ --({{[1 + {{({r \mathord{\left/ {\vphantom {r {{{R}_{0}}}}} \right. \kern-0em} {{{R}_{0}}}})}^{2}}]}^{{{1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-0em} 2}}}}-- \\ {{{{ - \,\,1)} \mathord{\left/ {\vphantom {{ - \,\,1)} {({{{{R}_{E}}} \mathord{\left/ {\vphantom {{{{R}_{E}}} {{{R}_{0}}}}} \right. \kern-0em} {{{R}_{0}}}})}}} \right. \kern-0em} {({{{{R}_{E}}} \mathord{\left/ {\vphantom {{{{R}_{E}}} {{{R}_{0}}}}} \right. \kern-0em} {{{R}_{0}}}})}}\} } \mathord{\left/ {\vphantom {{{{ - \,\,1)} \mathord{\left/ {\vphantom {{ - \,\,1)} {({{{{R}_{E}}} \mathord{\left/ {\vphantom {{{{R}_{E}}} {{{R}_{0}}}}} \right. \kern-0em} {{{R}_{0}}}})}}} \right. \kern-0em} {({{{{R}_{E}}} \mathord{\left/ {\vphantom {{{{R}_{E}}} {{{R}_{0}}}}} \right. \kern-0em} {{{R}_{0}}}})}}\} } {{{{[1 + {{{({r \mathord{\left/ {\vphantom {r {{{R}_{0}}}}} \right. \kern-0em} {{{R}_{0}}}})}}^{2}}]}}^{{{1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-0em} 2}}}}\,\,{\text{при}}\,\,z = 0}}} \right. \kern-0em} {{{{[1 + {{{({r \mathord{\left/ {\vphantom {r {{{R}_{0}}}}} \right. \kern-0em} {{{R}_{0}}}})}}^{2}}]}}^{{{1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-0em} 2}}}}\,\,{\text{при}}\,\,z = 0}}, \\ \end{gathered} $
$\begin{gathered} {\text{2}}{\text{.}}\,\,{{{{\sigma }_{0}}\partial \Phi } \mathord{\left/ {\vphantom {{{{\sigma }_{0}}\partial \Phi } {\partial z}}} \right. \kern-0em} {\partial z}} = 2\sum _{P}^{{{\text{av}}}}[{{{{\partial }^{2}}\Phi } \mathord{\left/ {\vphantom {{{{\partial }^{2}}\Phi } {\partial {{r}^{2}}}}} \right. \kern-0em} {\partial {{r}^{2}}}} + \left( {{1 \mathord{\left/ {\vphantom {1 r}} \right. \kern-0em} r}} \right){{\partial \Phi } \mathord{\left/ {\vphantom {{\partial \Phi } {\partial r}}} \right. \kern-0em} {\partial r}}] \\ {\text{при}}\,\,z = 90\,\,{\text{км}}, \\ \end{gathered} $
$3.\,\,\Phi \,\,\,\,{\text{непрерывен всюду}}\,\,{\text{при}}\,\,z \geqslant 0.$

5. РЕЗУЛЬТАТЫ

5.1. Оценка максимальных величин ожидаемых значений $E_{r}^{{{\text{(ion)}}}}$max в ионосфере, обусловленных сейсмогенными аномалиями квазистатического электрического поля ΔEz на поверхности земли, наблюдавшимися за пять дней до землетрясений EQ-1 и EQ-2

На рисунке 7 показаны распределения радиальных компонент поля $E_{r}^{{{\text{(ion)}}}}$ на уровне ионосферы в зависимости от эпицентрального расстояния r и исходных характеристик сейсмогенных аномалий ΔEz (указаны на рисунке над кривыми) на поверхности земли, которые могли иметь место за пять дней до землетрясений EQ-1 (штрихпунктирная линия) и EQ-2 (штриховая или сплошная линии) соответственно. Точками отмечена кривая, которая сопоставляется со сплошной линией и моделирует источник сейсмогенного возмущения ΔEz по формуле квази-Гауссоподобного распределения (1), использованного ранее в [Hegai et al., 2015], когда E0 = 2415 В/м и R0 = 400 км.

Рис. 7.

Возможные величины радиальной компоненты поля $E_{r}^{{{\text{(ion)}}}}$ в ионосфере за пять дней до землетрясений EQ-1 (штрихпунктирная линия) и EQ-2 (штриховая или сплошная линии) соответственно, в зависимости от эпицентрального расстояния r и исходных характеристик сейсмогенной аномалии ΔEz на поверхности земли (указаны на рисунке над кривыми). Точками показана кривая, которая сопоставляется со сплошной линией и моделирует источник сейсмогенного возмущения ΔEz по формуле (1), когда E0 = 2415 В/м и R0 = 400 км.

Сопоставление сплошной линии на рис. 7 с линией, данной точками, обнаруживает, что максимальное значение радиальной компоненты $E_{r}^{{{\text{(ion)}}}}$max в ионосфере в первом случае будет в ~1.5 раза выше, чем во втором. Это означает, что представление сейсмогенной аномалии ΔEz на поверхности земли с помощью формулы (4), аппроксимирующей данные измерений работы [Hao et al., 2000], дает при одинаковых моделях проводимости существенно более высокую величину поля в ионосфере по сравнению с использованным в работе [Hegai et al., 2015] выражением (1) для описания ΔEz, полученным из общих соображений.

Далее, сравнение сплошной линии (EQ-2, M = = 6.0, вариант I) со штрихпунктирной (EQ-1, M = = 4.7) показывает, что меньший размер источника R0 на поверхности земли (110 км против 400 км), меньшее значение E0 (890 В/м против 2415 В/м) и высокая величина $\sum _{P}^{{{\text{av}}}}$ (0.156 См против 0.0738 См, что в ~2.11 раза выше) при близких высотных профилях проводимости в диапазоне высот 0–90 км приводит к тому, что максимальное значение в ионосфере $E_{r}^{{{\text{(ion)}}}}$max(EQ-1)|I= 90° ≅ 0.009 мВ/м, тогда как $E_{r}^{{{\text{(ion)}}}}$max(EQ-2)|I= 90° ≅ 0.2 мВ/м, что в ~22 раза больше. С учетом результатов работы [Denisenko et al., 2018], т. к. I(EQ-1) ≅ 57.6° и I(EQ-2) ≅ 59.4°, а sin 57.6° ≅ 0.844 и sin 59.4° ≅ ≅ 0.861, получим:

$\begin{gathered} E_{r}^{{({\text{ion}})}}{\text{max}}{{\left. {({\text{EQ - 1}})} \right|}_{{I = {\text{57}}.{\text{6}}^\circ }}} \cong \\ \cong 0.{\text{844}}E_{r}^{{({\text{ion}})}}{\text{max}}{{\left. {({\text{EQ - 1}})} \right|}_{{I = {\text{9}}0^\circ }}} \cong \\ \cong 0.00{\text{8}}\,\,{{{\text{мВ}}} \mathord{\left/ {\vphantom {{{\text{мВ}}} {\text{м}}}} \right. \kern-0em} {\text{м}}},\,\,\,\,{{r}^{{{\text{max}}}}} \cong {\text{14}}0\,\,{\text{км,}} \\ \end{gathered} $
$\begin{gathered} E_{r}^{{({\text{ion}})}}{\text{max}}{{\left. {({\text{EQ - 2}})} \right|}_{{I = {\text{59}}.{\text{4}}^\circ }}} \cong \\ \cong 0.{\text{861}}E_{r}^{{({\text{ion}})}}{\text{max}}{{\left. {({\text{EQ - 2}})} \right|}_{{I = {\text{9}}0^\circ }}} \cong \\ \cong 0.17\,\,{{{\text{мВ}}} \mathord{\left/ {\vphantom {{{\text{мВ}}} {\text{м}}}} \right. \kern-0em} {\text{м}}},\,\,\,\,{{r}^{{{\text{max}}}}} \cong 480\,\,{\text{км}}\,\,\,\,{\text{(вариант I)}}{\text{.}} \\ \end{gathered} $
Второй возможный вариант для EQ-2 (штриховая кривая, вариант II) дает максимальные возможные значения радиальной компоненты $E_{r}^{{{\text{(ion)}}}}$max(EQ-2)|I= 90° ≅ 0.16 мВ/м и $E_{r}^{{{\text{(ion)}}}}$max(EQ-2)|I= 59.4° ≅ 0.14 мВ/м соответственно, а rmax ≅ 800 км. В этом случае величина максимума электрического поля уменьшается на 0.03 мВ/м, а его положение перемещается на 320 км дальше по сравнению с вариантом I (rmax ≅ 480 км). Из рисунка 7 также видно, что при увеличении RE (RE = 1.75R0) скорость падения величины поля по мере роста эпицентрального расстояния уменьшается по сравнению с вариантом I (RE = 0.6R0).

Таким образом, получается, что для заданных геофизических условий (на момент проявления сейсмогенной аномалии ΔEz на поверхности земли) в ионосфере за 5 дней до землетрясения EQ-2{03.05.1996 г. c М = 6.0; φe = 40.8° N, λe= = 109.7° E} можно было ожидать значений $E_{r}^{{({\text{ion}})}}$max(EQ-2)|I= 59.4° ≅ 0.14–0.17 мВ/м.

5.2. Оценка величины возможных значений $E_{r}^{{{\text{(ion)}}}}$max в ионосфере, перед землетрясением с магнитудой М = 7.2, произошедшим 20.03.2008 г. в приграничной области китайских провинций Синьцзян и Сицзян (Xinjiang-Xizang border region)

Указанное в заголовке подраздела 5.2. землетрясение (здесь и далее – EQ-3) с магнитудой 7.2 (глубина гипоцентра 10 км, т.е. коровое) произошло 20.03.2008 г. в 22:33 UT (или 21.03.2008 г. в 03:33 LT), при этом географические координаты эпицентра оказались следующими: φe = 35.49° N, λe = 81.47° (так же средние широты, при этом I ≅ 54°, а sin I ≅ 0.81). В марте 2008 г. усредненный суточный поток солнечного излучения F10.7 ∈ [70; 80] SFU, т.е солнечная активность была низкой. Измерений вертикальной компоненты электрического поля в области подготовки этого землетрясения (R0 ≅ 1300 км), центрированной координатами его эпицентра, не проводилось, однако предполагая, что 21.03.2008 г. в 00:00 LT (за 3.5 ч до землетрясения) в эпицентральной зоне имела место сейсмогенная электрическая аномалия ΔEz, можно оценить возможные проявления сейсмогенного электрического поля в ионосфере $E_{r}^{{{\text{(ion)}}}}$ на основании имеющихся характеристик самого землетрясения и значений величин геофизических параметров в области эпицентра незадолго до толчка.

На рисунке 8 представлена (аналогично рис. 5 и рис. 6) соответствующая моменту времени 21.03.2008 г. 00:00 LT и географическим координатам эпицентра землетрясения высотная модель проводимости, на которой в нижней части профиля (первая зона) штрихпунктирной и штриховой линиями указаны нижняя (LB) и верхняя (UB) границы разброса проводимости в области хорошей погоды по модели Baumgaertner et al. [2013] соответственно, остальные обозначения те же, что и на рисунках 5, 6.

По входным данным модели проводимости, представленным на рис. 8, был проведен расчет возможного распределения нормированных на E0 значений радиальной компоненты $E_{r}^{{{\text{(ion)}}}}$/E0 в зависимости от радиального расстояния от эпицентра r. Это распределение показано на рис. 9 для отношения RE/R0 = 1, что чуть меньше середины интервала [0.6; 1.75], оцененного по данным измерений в работе [Hao et al., 2000] (см. выше, раздел 2). Штрихпунктирная, сплошная и штриховая линии соответствуют своим моделям проводимости, приведенным в нижней части рис. 8 (первая высотная зона).

Рис. 9.

Распределения возможных нормированных на E0 значений радиальной компоненты $E_{r}^{{{\text{(ion)}}}}$/E0 (предположительно связанных с развитием завершающей фазы процесса подготовки EQ-3 и возникновением сейсмогенного возмущения ΔEz на поверхности земли в области будущего эпицентра за 3.5 ч до подземного толчка), в зависимости от радиального расстояния от эпицентра r для отношения RE/R0 = 1. Штрихпунктирная, сплошная и штриховая линии соответствуют своим моделям проводимости, приведенным в нижней части рис. 8 (LB, AV и UB).

Видно, что при неизменном положении максимума (rmax ≅ 1980 км), наибольшая величина отношения max($E_{r}^{{{\text{(ion)}}}}$/E0) составляет последовательно, в соответствии со своими профилями проводимости в нижней части:

max($E_{r}^{{{\text{(ion)}}}}$/E0)LB ≅ 1.82 × 10–7; max($E_{r}^{{{\text{(ion)}}}}$/E0)AV ≅ ≅ 3.63 × 10–7; max($E_{r}^{{{\text{(ion)}}}}$/E0)UB ≅ 7.27 × 10–7, т.е. каждое последующее значение больше предыдущего примерно в два раза, что характерно и для значений соответствующих им профилям проводимости на поверхности земли σ0(z = 0). Этот результат согласуется с представленным ранее в работе [Ampferer et al., 2010] выводом о пропорциональном уменьшении максимальной величины проникшего в ионосферу сейсмогенного поля при уменьшении проводимости у поверхности земли при прочих равных условиях. Из приведенных выше оценок max($E_{r}^{{{\text{(ion)}}}}$/E0) также следует, что для средних модельных значений σ0(z = 0) = = 1.0 × 10–13 С/м согласно модели Baumgaertner et al. [2013] и величине E0 = 1500 В/м = 1.5 кВ/м у поверхности земли, максимальная величина радиальной компоненты проникшего в ионосферу сейсмогенного электростатического поля (СЭП) могла бы составить $E_{r}^{{{\text{(ion)}}}}$max)|I= 54° = 1500 × 3.63 × × 10–7 sin 54° ≅ 0.44 мВ/м для землетрясения EQ-3 20.03.2008 г. с магнитудой М = 7.2 и учетом наклона геомагнитных силовых линий.

В работе [Choudhury et al., 2013] описана сейсмогенная аномалия ΔEz на поверхности земли, появившаяся за ~14 ч до землетрясения 12.06.2010 г. в районе Никобарских о-ов (Индия) с магнитудой M = 7.8 (по уточненным данным, в каталоге USGS (https://earthquake.usgs.gov/earthquakes/ browse/significant.php?year=2010), магнитуда этого землетрясения M = 7.5), которая существовала в течение ~41 мин. Она имела амплитуду ΔEz = = 1385 В/м на эпицентральном расстоянии Re = = 1763 км. Тогда, если в соответствии с (3), R0 = = RH ≅ exp(7.5) = 1808 км, получим на основе формулы (4) такие оценки: E0 = 3747 В/м ≅ ≅ 3.75 кВ/м (если RE = 0.6R0 = 1085 км) и E0 = = 2876 В/м ≅ 2.88 кВ/м (если RE = R0 = 1808 км). На этом фоне величина E0 = 1.5 кВ/м для рассмотренного нами землетрясения EQ-3 с M = 7.2 не кажется экстремальной. С другой стороны, полученная из дрейфовых измерений скоростей плазмы на спутнике FORMOSAT-5 величина сейсмогенного электрического поля в ионосфере, перпендикулярного магнитным силовым линиям (над эпицентральной зоной Тайваньского землетрясения 31.03.2002 г. с магнитудой M = 6.8 за пять дней до подземного толчка), как указывают авторы работы [Liu and Chao, 2017], оказалась равной 0.91 мВ/м.

Таким образом, полученная выше оценка величины возможного сейсмогенного поля в ионосфере перед землетрясением EQ-3 с магнитудой М = 7.2, произошедшим 20 марта 2008 г. в приграничной области китайских провинций Синьцзян и Сицзян, в соответствии с которой $E_{r}^{{{\text{(ion)}}}}$max|I= 54° ≅ ≅ 0.44 мВ/м, согласуется с результатами, представленными в работе [Liu and Chao, 2017].

5. ОБСУЖДЕНИЕ РЕЗУЛЬТАТОВ

Как видно, использование оптимальной (с точки зрения описания данных работы [Hao et al., 2000] по измерениям электрических полей сейсмогенных аномалий ΔEz у поверхности земли перед землетрясениями) формулы (4) в ~1.5 раза увеличивает проникшее в ионосферу поле $E_{r}^{{{\text{(ion)}}}}$max по сравнению со случаем, когда вместо нее используется выведенная ранее из общих соображений формула (1) при прочих равных условиях. Это происходит потому, что падение амплитуды возмущения поля ΔEz с ростом эпицентрального расстояния r по формуле (1) происходит быстрее, чем по формуле (4). Кроме этого, формула (4) позволяет учесть соотношение между характерным размером зоны подготовки землетрясения (R0) и характерным масштабом изменения сейсмогенного поля на поверхности земли (RE). По оценкам, сделанным на основании данных работы [Hao et al., 2000], отношение RE/R0 ∈ ∈ [0.6; 1.75].

Далее, для реально имевших место сейсмогенных аномалий ΔEz у поверхности земли за 5 дней до землетрясений EQ-1 (с М = 4.7) и EQ-2 (с М = = 6.0) в области средних широт и конкретно заданных геофизических условий (время, географические координаты эпицентров, модели проводимости) получены следующие оценки максимальных величин проникшего в ионосферу поля (с учетом наклона геомагнитных силовых линий в соответствии с работой [Denisenko et al., 2018]): $E_{r}^{{{\text{(ion)}}}}$max(EQ-1M= 4.7)|I= 57.6° ≅ 0.008 мВ/м и $E_{r}^{{{\text{(ion)}}}}$max(EQ-2M= 6.0)|I= 59.4° ≅ 0.14–0.17 мВ/м. В случае землетрясения EQ-3 (c М = 7.2), за ~3.5 ч до него, ГИПОТЕТИЧЕСКИ можно было бы ожидать в ионосфере величину $E_{r}^{{{\text{(ion)}}}}$max(EQ-3M= 7.2)|I= 54.0° ≅ ≅ 0.44 мВ/м. Полученные величины следует удвоить, если проводимость σ0(z = 0) = $\sigma _{0}^{{{\text{UB}}}}$ = = 2.0 × 10–13 См/м, т.е. соответствует верхней (UB) границе разброса проводимости в области хорошей погоды по модели Baumgaertner et al. [2013] (эта же величина приводится также в статье [Cole and Pierce, 1965]). Такая ситуация вполне реалистична, т. к. выход радона в эпицентральных областях будущих землетрясений (см., например, монографию [Pulinets and Boyarchuk, 2004]) увеличивает степень ионизации приземного слоя атмосферы, а следовательно и его проводимость.

Наконец, в работе [Denisenko et al., 2018] для оценки возможной максимальной величины перпендикулярного магнитным силовым линиям электрического поля в ионосфере сейсмогенного источника (см. формулы (34) и (37) этой работы) были использованы “входные” параметры из работы [Hegai et al., 2015]: E0 = 1000 В/м, σ0(z = 0) = = 1.0 × 10–13 С/м, $\sum _{P}^{0}$ = 0.05 См ($\sum _{P}^{{{\text{av}}}}$ в наших обозначениях), x0 = 1100 км (R0 в наших обозначениях в связи с другой геометрией источника сейсмогенного поля, что, согласно формуле (3), соответствует землетрясению с магнитудой M = ln R0 = = ln 1100 ≅ 7.0). По оценкам авторов работы, по их модели такие значения входных параметров дают величину $E_{ \bot }^{{{\text{(ion)}}}}$max|I= 35.0° ≅ 0.9 мВ/м. При этом авторы указывают, что ВСЕ взятые ключевые входные параметры являются экстремальными и ставят вопрос: “Могут ли все они быть достигнуты в реальных условиях одновременно?”. Ответ на этот вопрос состоит в том, что вовсе НЕ ВСЕ выбранные параметры являются экстремальными, хотя такое их сочетание является достаточно редким событием, как и крупные землетрясения (с М ≥ 6.0), и такая их комбинация МОЖЕТ возникнуть в реальных условиях одновременно.

Во-первых, как следует из проведенного выше в этой работе рассмотрения, величина E0 ≅ 1 кВ/м у поверхности земли, согласно измерениям Hao et al. [2000], достигается уже перед землетрясениями с М ∈ [4.1; 4.7], а перед землетрясением с М = 6.0 величина E0 ≅ 1.3–2.4 кВ/м (см. раздел 2).

Во-вторых, величина σ0(z = 0) = 1.0 × 10–13 С/м, как это следует из работы [Baumgaertner et al., 2013] является СРЕДНИМ значением модели проводимости в области хорошей погоды.

В-третьих, в частности, для рассмотренного в этой работе конкретного примера среднеширотного землетрясения EQ-3 с магнитудой М = 7.2 (что соответствует R0 ≅ 1300 км), произошедшего 20.03.2008 г. в приграничной области китайских провинций Синьцзян и Сицзян, за ~3.5 ч до подземного толчка и географических координат эпицентра, согласно формулам (6) и (8) по современным эмпирическим моделям ионосферы IRI-2016 и нейтральной атмосферы NRLMSIS-00, интегральная проводимость Педерсена $\sum _{P}^{{{\text{av}}}}$ ≅ 0.0655 См (ночь, 00 LT, F10.7 ∈ [70; 80] SFU – низкая солнечная активность). Это близко к значению $\sum _{P}^{0}$ = = 0.05 См, принятому в указанной авторами исследования [Denisenko et al., 2018] работе [Hegai et al., 2015], т.к. δ = 100 × ($\sum _{P}^{{{\text{av}}}}$$\sum _{P}^{{\text{0}}}$)/$\sum _{P}^{{{\text{av}}}}$ ≅ 24%. Здесь следует отметить также, что для крупномасштабных плазменных “пузырей” (bubbles), простирающихся до средних широт, $\sum _{P}^{{\text{0}}}$ может иметь значение менее 0.01 См согласно исследованию [Huang et al., 2007].

При этом получается, что если E0 ≅ 1.5 кВ/м, то $E_{r}^{{{\text{(ion)}}}}$max|I= 54.0° ≅ 0.44 мВ/м (и вдвое больше, если σ0(z = 0) = 2.0 × 10–13 С/м, что отвечает верхней (UB) границе разброса проводимости в области хорошей погоды по модели Baumgaertner et al. [2013]).

С учетом того, что число сильных землетрясений с М ≥ 6.0 в году имеет тенденцию увеличиваться в периоды низкой солнечной активности (см. работу [Хегай, 2018]), вероятность благоприятного сочетания нужных геофизических параметров (способствующих проникновению в ионосферу полей электрических сейсмогенных аномалий у поверхности земли ΔEz) и землетрясений, с “подходящими” характеристиками, будет возрастать.

5. ВЫВОДЫ

1. Построена аналитическая модель возмущений вертикальной компоненты сейсмогенного электрического поля ΔEz у поверхности земли перед некоторыми землетрясениями (по данным измерений в приземном слое атмосферы на средних широтах, приведенных в работе [Hao et al., 2000]) и проведен расчет его величины на уровне ионосферы с учетом наклона геомагнитных силовых линий в соответствии с результатами исследования [Denisenko et al., 2018].

2. Полученная аналитическая аппроксимация сейсмогенного возмущения вертикальной компоненты электрического поля ΔEz на поверхности земли адекватно описывает данные измерений Hao et al. [2000]. Показано, что использование этой аппроксимации (формула (4) настоящей работы) приводит к увеличению максимальной амплитуды проникшего в ионосферу поля сейсмогенного источника $E_{r}^{{{\text{(ion)}}}}$max в ~1.5 раза (при прочих равных условиях) по сравнению с ранее использованным в работе [Hegai et al., 2015] представлением ΔEz в виде квази-Гауссоподобного распределения (формула (1) настоящей работы), выведенным из общих соображений.

3. Для реально имевших место сейсмогенных аномалий ΔEz у поверхности земли за 5 дней до землетрясений EQ-1 (с магнитудой М = 4.7) и EQ-2 (с магнитудой М = 6.0) в области средних широт и соответственно определенных геофизических условий (время, географические координаты эпицентров, модели проводимости) получены следующие значения максимальных величин проникшего в ионосферу поля, перпендикулярного геомагнитным силовым линиям):

$\begin{gathered} E_{ \bot }^{{({\text{ion}})}}{\text{max}}{{\left. {({\text{EQ - }}{{{\text{1}}}_{{M = 4.7}}})} \right|}_{{I = {\text{57}}.{\text{6}}^\circ }}} \cong 0.00{\text{8}}\,\,{{{\text{мВ}}} \mathord{\left/ {\vphantom {{{\text{мВ}}} {\text{м}}}} \right. \kern-0em} {\text{м}}} \\ {\text{и}}\,\,\,\,E_{ \bot }^{{({\text{ion}})}}{\text{max}}{{\left. {({\text{EQ - }}{{{\text{2}}}_{{M = 6.0}}})} \right|}_{{I = {\text{59}}.{\text{4}}^\circ }}} \cong 0.14{\kern 1pt} - {\kern 1pt} 0.17\,\,{{{\text{мВ}}} \mathord{\left/ {\vphantom {{{\text{мВ}}} {\text{м}}}} \right. \kern-0em} {\text{м}}}. \\ \end{gathered} $
В случае землетрясения EQ-3 (c магнитудой М = 7.2), за ~3.5 ч до него, ГИПОТЕТИЧЕСКИ можно было бы ожидать в ионосфере значения
$E_{ \bot }^{{({\text{ion}})}}{\text{max}}{{\left. {({\text{EQ - }}{{{\text{3}}}_{{M = 7.2}}})} \right|}_{{I = 54.0^\circ }}} \cong 0.44\,\,{{{\text{мВ}}} \mathord{\left/ {\vphantom {{{\text{мВ}}} {\text{м}}}} \right. \kern-0em} {\text{м}}}.$
Все полученные величины $E_{ \bot }^{{{\text{(ion)}}}}$max должны быть удвоены, если величина проводимости σ0 у поверхности земли соответствует верхней границе разброса проводимости в области хорошей погоды по модели Baumgaertner et al. [2013]. Представленные результаты хорошо согласуются с результатами выполненных ранее теоретических рассмотрений в работах [Hegai et al., 2015; Denisenko et al., 2018] (и уточняют их), а также экспериментального исследования [Liu and Chao, 2017].

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

  1. Хегай В.В. О корреляции между сильными землетрясениями и уровнем солнечной активности в 21, 22 и 23 циклах / Астрономия-2018. Т. 2. Солнечно-земная физика – современное состояние и перспективы. Ред.  Обридко В.Н. М.: Тровант. С. 262–265. 2018. https://doi.org/10.31361/eaas.2018-2.066

  2. Ampferer M., Denisenko V.V., Hausleitner W., Krauss S., Stangl G., Boudjada M.Y., Biernat H. K. Decrease of the electric field penetration into the ionosphere due to low conductivity at the near ground atmospheric layer // Ann. Geophysicae. V. 28. Iss. 3. P. 779–787. 2010.

  3. Baumgaertner A.J.G., Thayer J.P., Neely III R.R., Lucas G. Toward a comprehensive global electric circuit model: Atmospheric conductivity and its variability in CESM1 (WACCM) model simulations // J. Geophys. Res. Atmos. V. 118. P. 9221–9232. 2013. https://doi.org/10.1002/jgrd.50725

  4. – Choudhury A., Guha A., De B.K., Roy R. A statistical study on precursory effects of earthquakes observed through the atmospheric vertical electric field in northeast India // Ann. Geophys. [S.l.] V. 56. № 3. P. R0331. 2013. https://doi.org/10.4401/ag-6235

  5. Cole Jr.R.K., Pierce E.T. Electrification in the Earth’s atmosphere for altitudes between 0 and 100 kilometers // J. Geophys. Res. V. 70. № 12. P. 2735–2749. 1965. https://doi.org/10.1029/JZ070i012p02735

  6. Denisenko V.V., Boudjada M.Y., Horn M., Pomozov E.V., Biernat H.K., Schwingenschuh K., Lammer H., Prattes G., Cristea E. Ionospheric conductivity effects on electrostatic field penetration into the ionosphere // Nat. Hazards Earth Syst. Sci. V. 8. Iss. 5. P. 1009–1017. 2008. https://doi.org/10.5194/nhess-8-1009-2008

  7. Denisenko V. V.,Ampferer M., Pomozov E.V., Kitaev A.V., Hausleitner W., Stangl G., Biernat H. K. On electric field penetration from ground into the ionosphere // J. Atmos. Solar-Terr. Phys. V. 102. P. 341–353. 2013. https://doi.org/10.1016/j.jastp.2013.05.019

  8. Denisenko V.V., Nesterov S. A., Boudjada M.Y., Lammere H. A mathematical model of quasistationary electric field penetration from ground to the ionosphere with inclined magnetic field // J. Atmos. Solar-Terr. Phys. V. 179. P. 527–537. 2018. https://doi.org/10.1016/j.jastp.2018.09.002

  9. Dobrovolsky I.R., Zubkov S.I., Myachkin V.I. Estimation of the size of earthquake preparation zones // Pure Appl. Geophys. V. 117. № 5. P. 1025–1044. 1979.

  10. Hao J., Tang T.M., Li D.R. Progress in the research of atmospheric electric field anomaly as an index for short-impending prediction of earthquakes // J. Earthquake Pred. Res. V. 8. Iss. 3. P. 241–255. 2000.

  11. Hegai V.V., Kim V.P., Liu J.Y. On a possible seismo-magnetic effect in the topside ionosphere // Adv. Space Res. V. 56. Iss. 8. P. 1707–1713. 2015. https://doi.org/10.1016/j.asr.2015.07.034

  12. – https://ccmc.gsfc.nasa.gov/modelweb/models/iri2016_vitmo.php

  13. – https://earthquake.usgs.gov/earthquakes/browse/significant.php?year=2010

  14. – https://modelweb.gsfc.nasa.gov/atmos/nrlmsise00.html

  15. Huang C.-S., Foster J.C., Sahai Y. Significant depletions of the ionospheric plasma density at middle latitudes: A possible signature of equatorial spread F bubbles near the plasmapause // J. Geophys. Res. V. 112. P. A05315. 2007. https://doi.org/10.1029/2007JA012307

  16. Kim V.P., Hegai V.V. On the electric fields produced by dipolar coulomb charges of an individual thundercloud in the ionosphere // J. Astron. Space Sci. V. 32. № 2. P. 141–144. 2015. https://doi.org/10.5140/JASS.2015.32.2.141

  17. Liu J.-Y., Chao C.K. An observing system simulation experiment for FORMOSAT-5/AIP detecting seismo-ionospheric precursors // Terr. Atmos. Ocean Sci. V. 28. № 2. P. 117–127. 2017. https://doi.org/10.3319/TAO.2016.07.18.01(EOF5)

  18. Park C.G., Dejnakarintra M. Penetration of thundercloud electric fields into the ionosphere and magnetosphere: 1. Middle and subauroral latitudes // J. Geophys. Res. V. 78. № 28. P. 6623–6633. 1973. https://doi.org/10.1029/JA078i028p06623

  19. – Park C.G., Dejnakarintra M. Thundercloud Electric Fields in the Ionosphere / Electrical Processes in Atmospheres. Eds. Dolezalek H., Reiter R., Landsberg H.E. Darmstadt: Dr. Dietrich Steinkopff Verlag GmbH & Co. KG., 1976. P. 544–551. https://doi.org/10.1007/978-3-642-85294-7_84

  20. Pulinets S.A., Boyarchuk K.A. Ionospheric precursors of earthquakes. Berlin: Springer. 315 p. 2004.

  21. Rycroft M.J., Odzimek A., Arnold N.F., Füllekrug M., Kułak A., Neubert T. New model simulations of the global atmospheric electric circuit driven by thunderstorms and electrified shower clouds: The roles of lightning and sprites // J. Atmos. Solar-Terr. Phys. V. 69. Iss. 17–18. P. 2485–2509. 2007. https://doi.org/10.1016/j.jastp.2007.09.004

  22. – Schunk R.W. A mathematical model of the middle and high latitude ionosphere // Pure Appl. Geophys. V. 127. № 2/3. P. 255–303. 1988.

Дополнительные материалы отсутствуют.