Физика Земли, 2021, № 5, стр. 133-145
Эффекты ударной волны, генерируемой на поверхности Земли при падении космических тел размером от 20 м до 3 км
Д. О. Глазачев 1, *, О. П. Попова 1, **, Е. Д. Подобная 1, ***, Н. А. Артемьева 1, ****, В. В. Шувалов 1, *****, В. В. Светцов 1, ******
1 Институт динамики геосфер имени академика М.А. Садовского РАН
г. Москва, Россия
* E-mail: GlazachevD@gmail.com
** E-mail: olga@idg.chph.ras.ru
*** E-mail: epodobnaya@gmail.com
**** E-mail: natalia.artemeva@gmail.com
***** E-mail: valeryvshuvalov@gmail.com
****** E-mail: svettsov07@rambler.ru
Поступила в редакцию 19.02.2021
После доработки 28.03.2021
Принята к публикации 01.04.2021
Аннотация
Разрушения на поверхности Земли, вызываемые ударной волной, являются одним из наиболее важных и опасных последствий падений астероидов и комет. Избыточное давление и скорость ветра за фронтом ударной волны, приводящие к тем или иным опасным последствиям, могут быть оценены на основе специально разработанных моделей взаимодействия космических объектов с атмосферой и поверхностью Земли. Ударная волна является также причиной сейсмических эффектов, но в данной работе рассматривается только прямое действие ударной волны. Систематическое численное моделирование взаимодействия космических объектов с атмосферой было ранее проведено для большого количества различных сценариев в рамках гидродинамической модели. Анализ результатов этого моделирования позволяет предложить аппроксимационные соотношения, которые дают возможность оценить избыточное давление, скорость ветра за ударной волной и их распределение на поверхности, если известны параметры импактора, его скорость и угол наклона траектории. Эти соотношения учитывают пространственную неоднородность распределения избыточного давления на поверхности Земли. Полученные соотношения подобия были протестированы на данных Челябинского и Тунгусского событий, и было показано хорошее соответствие с данными разрушений. Полученные аппроксимационные соотношения используются в интернет-калькуляторе (http://www.AsteroidHazard.pro), позволяющем быстро и достаточно точно оценивать различные эффекты падений.
ВВЕДЕНИЕ
Земля постоянно сталкивается с космическими телами, однако размеры большинства тел слишком малы для возникновения опасных эффектов. Небольшие объекты – метеороиды с размером от микроскопической пылинки до 1 м, как правило, разрушаются и сгорают в атмосфере. Фрагменты наиболее крупных из них могут достигать поверхности как метеориты, но опасности для людей они не несут. При падении большей части тел с размерами 1–10 м атмосфера Земли все еще является хорошей защитой, хотя железные объекты такого (и даже меньшего) размера долетают до поверхности и вызывают образование кратеров и кратерных полей таких, как, например, Сихоте-Алинский железный дождь [Krinov, 1971; Nemtchinov, Popova, 1997] или метеорит Мораско [Bronikowska et al., 2017].
С увеличением размера космические объекты становятся гораздо более опасными. Они способны причинить значительный ущерб, что и показало падение 19-метрового Челябинского астероида [Popova et al., 2013; Brown et al., 2013]. Челябинское событие вызвало массовые повреждения в жилых домах и общественных зданиях, за медицинской помощью только официально обратилось более 1600 человек, финансовые затраты на восстановление разрушенной инфраструктуры были оценены в 1 млрд рублей [Акимов и др., 2015; Kartashova et al., 2018]. Удары тел типа Челябинского (и крупнее) случаются не каждый год, однако такие падения невозможно предсказать заранее. Кроме этого, для надежной оценки последствий возможного удара необходимо проводить довольно сложные и времязатратные расчеты, причем нужно как можно точнее определить параметры космического тела, траекторию и скорость входа. Такие параметры, как правило, наиболее точно определяются лишь в последний момент (перед входом в атмосферу), когда на проведение детальных расчетов уже нет времени. Учитывая сложность и трудоемкость оценки последствий через непосредственное моделирование взаимодействия метеороида с атмосферой и поверхностью Земли, представляется актуальной и важной возможность простой быстрой оценки поражающих факторов, на основе предварительно подготовленной (рассчитанной) сетки моделей последствий. Быструю оценку можно получить по аппроксимационным формулам, построенным по данным моделей. В данной работе представлены результаты расчетов такой сетки и приведены аппроксимационные формулы.
Одним из наиболее важных и опасных последствий падений астероидов и комет являются разрушения на поверхности Земли, вызываемые ударной волной в нижних слоях атмосферы. Например, в случае Челябинского события 15.02.2013 г. именно избыточное давление в ударной волне было главным поражающим фактором. Кроме того, ударная волна, распространяясь в атмосфере, доходит до Земли и порождает сейсмические волны аналогично сосредоточенному воздушному взрыву [Svetsov et al., 2017; Хазинс и др., 2018]. Но в данной работе рассматривается только прямое действие ударной волны. Избыточное давление и скорость ветра за фронтом ударной волны, приводящие к тем или иным опасным последствиям, могут быть оценены на основе реальных данных о ядерных взрывах [Collins et al., 2005; 2017] или на основе специально разработанных моделей взаимодействия космических объектов с атмосферой и поверхностью Земли [Shuvalov et al., 2013].
МОДЕЛЬ ВЗАИМОДЕЙСТВИЯ КОСМИЧЕСКИХ ТЕЛ С АТМОСФЕРОЙ И ПОВЕРХНОСТЬЮ ЗЕМЛИ
Моделирование процессов, сопровождающих падение космического тела на Землю, проводилось в два этапа. На первом этапе моделировалось движение метеороида в атмосфере с учетом его деформации, торможения, разрушения и испарения [Shuvalov, Artemieva, 2002; Shuvalov, Trubetskaya, 2007; Shuvalov et al., 2013; 2017a; 2017b]. Падающее тело предполагается квазижидким (с нулевой прочностью), а его деформация под действием аэродинамической нагрузки описывается уравнениями гидродинамики. Заметная деформация метеороидов размером больше нескольких десятков метров начинается на высотах (около 30–40 км), где аэродинамические нагрузки значительно превышают прочность падающих тел, и поэтому можно считать, что деформируемое тело в начале процесса деформации уже полностью разрушено [Svetsov et al., 1995]. Кроме того, многие крупные астероиды изначально представляют собой сильно разрушенные объекты типа кучи щебня (rubble-pile asteroids) с низкой прочностью [Sánchez, Scheeres, 2014].
Полностью разрушенное тело не вполне эквивалентно жидкому из-за присутствия внутреннего трения. Наличие внутреннего трения (касательных напряжений, пропорциональных давлению) может повлиять на характер деформации разрушенного метеороида, в частности, замедлить фрагментацию и торможение метеороида [Ivanov et al., 1997]. Влияние внутреннего трения и акустической флюидизации было рассмотрено в работе [Шувалов, Трубецкая, 2010]. При использовании жидкостной модели нужно иметь в виду, что реальная высота торможения метеороида может оказаться на несколько километров меньше, чем полученная в расчетах по жидкостной модели, что приведет к некоторому увеличению максимального избыточного давления.
Моделирование на первом этапе проводилось в двумерном приближении в системе координат, связанной с падающим телом, которое обдувается воздухом. Плотность воздуха в набегающем потоке изменяется в соответствии со стратификацией атмосферы, а скорость потока равна скорости тела. При торможении тела скорость набегающего потока соответственно уменьшается. Для моделирования переноса излучения в метеорных парах и воздухе в процессе полета и разрушения метеороида применялось приближение лучистой теплопроводности. Росселандовы пробеги излучения для воздуха при типичных температурах за фронтом головной ударной волны (1–5 эВ) не превышают 1 м на высотах ниже 30–40 км [Авилова и др., 1970], поэтому это приближение вполне применимо для крупных метеороидов. Излучение, падающее на тело, приводит к его испарению. Давление поступающих паров считается равным давлению газа вблизи поверхности тела, а температура и плотность паров определяются из кривой фазового равновесия. Расчеты первого этапа заканчивались, когда метеороид достигал поверхности Земли или терял свою начальную кинетическую энергию.
Полученные распределения газодинамических и термодинамических параметров использовались в качестве начальных данных для второго этапа расчетов, на котором в системе координат, связанной с поверхностью Земли, моделировалось образование кратера, взаимодействие продуктов выброса с атмосферой, распространение ударной волны в атмосфере Земли на большие расстояния [Artemieva, Shuvalov, 2019]. Оба этапа расчетов были реализованы с помощью численного комплекса SOVA [Shuvalov, 1999]. Независимо от угла наклона траектории на втором этапе решалась трехмерная задача. В трехмерной задаче перенос излучения учитывался в приближении лучистой теплопроводности на высотах до 80 км и в приближении объемного высвечивания при малых оптических толщинах на больших высотах.
В рамках этой модели были проведены серийные расчеты, рассматривались тела с размерами от 20 до 3000 м, входящие в атмосферу под углами от 15° до 90° со скоростями входа от 20 до 70 км/с, частично описанные в работе [Artemieva, Shuvalov, 2019]. Скорость и угол входа задавались в точке входа на высоте 100 км, угол отсчитывается от горизонтали. Всего было рассмотрено 39 вариантов для астероидов (с плотностью 3320 кг/м3), 39 вариантов для комет (с плотностью 1000 кг/м3) и 3 варианта для астероидов с меньшей плотностью (2650 кг/м3). Кинетическая энергия импакторов лежит в диапазоне от 0.5 Мт ТНТ до 4.2 Тт ТНТ. В качестве примера на рис. 1 приведено распределение максимального избыточного давления ΔP на поверхности Земли (ΔP = P – P0, где P0 – нормальное атмосферное давление в 1 атм).
Результатом моделирования является распределение максимального избыточного давления на поверхности Земли, учитывающее как падающую, так и отраженную волны. Дополнительных препятствий, которые могли бы быть источниками дополнительных отраженных волн, не вводилось, свойства подстилающей поверхности, которые могли бы повлиять на амплитуду отраженной волны, не учитывались. Результаты моделирования можно рассматривать как максимальную оценку вблизи эпицентра и минимальную на больших расстояниях (см. обсуждение в работе [Artemieva, Shuvalov, 2019]). В данной работе длительность воздействия избыточного давления не оценивалась.
Результаты серийного численного моделирования дали возможность построить аппроксимационные соотношения для избыточного давления и скорости ветра, зависящие только от параметров падающего тела и позволяющие быстро оценивать последствия падений.
ЭФФЕКТИВНАЯ ВЫСОТА ЭНЕРГОВЫДЕЛЕНИЯ
В зависимости от размера объекта, материала импактора и угла наклона траектории к горизонту могут реализовываться различные сценарии входа в атмосферу.
Вхождение в атмосферу высокоскоростного космического тела можно классифицировать как “воздушный взрыв”, если основная часть энергии этого тела выделяется в атмосфере и основная масса тела не достигает поверхности (или достигает ее со скоростью намного меньшей, чем скорость входа), но возникающие ударные волны достаточно сильны, чтобы достичь поверхности Земли и вызвать те или иные повреждения [Shuvalov et al., 2013]. Воздушные взрывы характерны для вертикального входа астероидов размером менее 100 м и комет менее 150 м. Удары космических тел размером более 300 м, как правило, приводят к возникновению кратера. Уменьшение угла наклона траектории сдвигает эти границы к большим размерам. Например, комета диаметром 300 м при угле входа <30° теряет большую часть своей энергии в атмосфере, приводя к мощному воздушному взрыву, и не образует кратер. Возможны и промежуточные случаи, когда часть энергии выделяется в атмосфере, но тем не менее энергии хватает и для образования кратера.
Эффективная высота Heff – высота, на которой импактор теряет большую часть энергии и заметно тормозится, дает возможность определить границу между различными сценариям [Shuvalov et al., 2016]. В этой работе было показано, что в первом приближении избыточное давление в ударной волне от падающего космического тела можно оценить при помощи рассмотрения взрыва на высоте Heff, энергия которого определяется энергией падающего тела. Конечно, такое приближение не дает возможность учесть пространственную неоднородность в распределении избыточного давления, которая оказывается существенной при входе космических тел в атмосферу под углами менее 60°.
Эффективная высота определяется как
(1)
${{H}_{{eff}}} = \,\,~{{\left( {H - 1.3H{\text{ln}}\left( {D\frac{{\sin \alpha }}{H}{{{\left( {\frac{{{\rho }}}{{{{{{\rho }}}_{0}}}}} \right)}}^{{{2 \mathord{\left/ {\vphantom {2 3}} \right. \kern-0em} 3}}}}} \right)} \right)} \mathord{\left/ {\vphantom {{\left( {H - 1.3H{\text{ln}}\left( {D\frac{{\sin \alpha }}{H}{{{\left( {\frac{{{\rho }}}{{{{{{\rho }}}_{0}}}}} \right)}}^{{{2 \mathord{\left/ {\vphantom {2 3}} \right. \kern-0em} 3}}}}} \right)} \right)} {1000}}} \right. \kern-0em} {1000}},$Если для выбранных характеристик импактора значение Heff равно нулю, то можно считать, что это кратерообразующее событие; если Heff > 5–7 км, то это событие типа метеорного взрыва; импакты с высотой энерговыделения от 0 до 5−7 км являются переходными случаями (как правило, при размерах D ≈ 100–300 м). Стоит отметить, что выражение (1) для Heff имеет точность порядка нескольких км [Shuvalov et al., 2016].
Зависимость имеет сложный характер, например, удар астероида с энергией 83 Мт ТНТ с D = 100 м, скоростью V = 20 км/с и углом входа α = 45° приводит к образованию кратера, а падение кометы с большей энергией в 307 Мт ТНТ с D = 100 м, V = 70 км/с и α = 30° является метеорным взрывом. Чем меньше угол входа, тем больший диаметр импактора требуется для образования кратера, эффективная высота для кометных тел выше, чем для астероидов с аналогичными параметрами (рис. 2а).
В широко известном интернет-калькуляторе для быстрой оценки последствий импактов Earth Impact Effects Program [Collins et al., 2005; 2017] также оценивается эффективная высота энерговыделения Zb. В этом случае высота Zb определяется на основе полуэмпирической модели разрушения объекта (pancake модели), описывающей поперечное расширение фрагментировавшего метеороида. Моделирование методом Монте-Карло при заданных распределениях по скоростям, углу входа, прочности и плотности космических тел дало возможность получить зависимость Zb от энергии импактора [Collins et al., 2017]. На рис. 2б приведено сравнение зависимостей Zb и Heff от размера тела при близких значениях других параметров. Для астероидов величина Heff оказывается ниже, чем Zb, для комет зависимости близки друг к другу.
СТЕПЕННОЙ ЗАКОН ДЛЯ АППРОКСИМАЦИИ ИЗБЫТОЧНОГО ДАВЛЕНИЯ
Распределения максимального избыточного давления ΔP (аналогичные представленному на рис. 1) для 81 варианта дают возможность для поиска простых аппроксимационных зависимостей. Все имеющиеся распределения максимального давления за фронтом ударной волны ΔP были переинтерполированы на равномерную сетку в полярных координатах, был определен центр и усредненные профили распределений. Эти профили использовались для поиска аппроксимаций.
Исследуемые распределения максимального избыточного давления имеют большой разброс значений. Например, для небольшого импактора (20–30 м в диаметре) максимальное избыточное давление не превышает 2.5% от нормального давления атмосферы и составляет всего 0.025 атм, в то время как для километровых тел значения ΔP могут превышать несколько тысяч атмосфер [Artemieva, Shuvalov, 2019]. Действительно, интересующая нас область, где требуется наиболее точное описание избыточного давления, лежит в области от 0.005 (примерный порог разрушения остекления [Gi et al., 2018]) до примерно 6–10 атм, при которых происходит разрушение мостов, груженых железнодорожных вагонов [Glasstone, Dolan, 1977]. Разрушение стен из бетонных блоков и кирпича (толщиной 24–36 см) происходит при избыточном давлении ΔP ~ 0.2–0.35 атм (ΔP = 20–35 кПа). Детально описывать область избыточных давлений >100 атм не имеет смысла, это будут области полного разрушения. Поэтому при поиске аппроксимаций использовалась весовая функция, зависящая от максимума распределения и позволяющая обращать внимание прежде всего на диапазон избыточных давлений 0.005–100 атм. Отметим, что в большинстве расчетов избыточные давления ниже 0.1–0.05 атм определялись с невысокой точностью (до нескольких раз).
Следует отметить, что для метеорных взрывов максимальное значение избыточного давления на поверхности Земли наблюдается на некотором расстоянии Y0 от точки пересечения траектории с поверхностью Земли (без учета торможения) (рис. 1). Для сдвига Y0 ниже предложено аппроксимационное соотношение, но поиск зависимости для максимального избыточного давления проводился по центрированным распределениям с началом координат в точке максимального значения ΔP.
Аппроксимации для ΔP искались в виде степенного закона с учетом эффективной высоты энерговыделения Heff :
где: r – расстояние от центра распределения; Heff – высота источника; p (power) – степень, характеризующая уменьшение величины ΔP с расстоянием; a – максимальная амплитуда.Одновременный поиск параметров p и a может оказаться затруднителен в ряде случаев в силу разного характера переменных (степени и линейного множителя). При логарифмическом преобразовании зависимость приобретает линейный вид по параметрам $p$ и $\lg a$:
(3)
${\text{lg}}\left( {\Delta P} \right) = \lg a - 0.5p\lg \left( {{{r}^{2}} + H_{{eff}}^{2}} \right).$При определении оптимальных значений параметров $\left( {p,a} \right)$ и $(p,\lg a)$ использовались две функции невязки, относительная невязка RD для решения в виде (2) для ΔP и абсолютная невязка AD для решения в виде (3) для${\text{\;lg}}\left( {{{\Delta }}P} \right)$:
(4)
$RD = \sqrt {\frac{1}{n}\mathop \sum \limits_{i = 1}^n {{{\left( {\frac{{\Delta {{P}_{i}} - \Delta P{{M}_{i}}}}{{\Delta {{P}_{i}}}}} \right)}}^{2}}} ,$(5)
$AD = \sqrt {\frac{1}{n}\mathop \sum \limits_{i = 1}^n {{{\left( {\lg (\Delta {{P}_{i}}) - \lg \left( {\Delta P{{M}_{i}}} \right)} \right)}}^{2}}} ,$На рис. 3а приведено поле невязки в координатах (p, a). Область минимума имеет вытянутый характер, т.е. можно получить очень хорошее (хотя и не идеальное) решение для несколько различающихся наборов параметров (p, a). При этом для некоторого зафиксированного значения степени p можно подобрать такой коэффициент a, при котором будет получаться хорошая аппроксимация усредненного профиля. Это дает возможность найти универсальный показатель степени p, одинаковый для всех вариантов, и позволяет значительно снизить сложность разрабатываемой аппроксимационной формулы.
Для каждого значения параметра p (power) от 1 до 2 с шагом 0.01 оценивалось среднее значение невязки по всем вариантам распределений усредненного профиля избыточного давления (рис. 3б) и для оптимального для каждого варианта параметра a. Минимум невязки наблюдается при значении параметра p = 1.56 ± 0.02. Это и есть оптимальная степень для описания избыточного давления как метеорных взрывов, так и кратерообразующих ударов. Такое изменение величины с расстоянием занимает промежуточное положение между падением давления с расстоянием от линейного (p = 1) или точечного источников (p = 2), что кажется достаточно разумным, поскольку сложный характер энерговыделения при импакте представляет собой нечто среднее между точечным и линейным взрывами.
Комбинация параметров ударника в различных степенях рассматривалась как вероятная зависимость для каждого исследуемого параметра, например, для значения коэффициента a, позволяющего описывать все рассчитанные варианты, для сдвига Y0, для параметра эллиптичности (см. далее) и других необходимых величин. Угол входа был включен через тригонометрическую функцию. В некоторых случаях использовалась логарифмическая функция кинетической энергии. Коэффициенты и показатели степени для этих комбинаций были оценены путем минимизации функции невязки для всех имеющихся расчетных вариантов или для отдельных групп вариантов, например, метеорные взрывы и кратерообразующие варианты, астероидные тела и кометные могли рассматриваться как отдельные однородные группы. Окончательный выбор между зависимостями был сделан с учетом ошибки подбора (минимальное значение функции неопределенности) и простоты.
АППРОКСИМАЦИОННЫЕ СООТНОШЕНИЯ
Построенные аппроксимационные соотношения для распределения максимального избыточного давления зависят только от параметров ударника, которые включают его диаметр D (м), скорость V (км/с), угол входа α к горизонту (градусы) и плотность ρ (кг/м3). Аппроксимации основаны на численных результатах для комет и астероидов с заданными плотностями, поэтому для ударников с другими плотностями используется линейная интерполяция по плотности.
Общая формула для всех тел как выделяющих свою энергию в атмосфере, так и кратерообразующих ударов можно записать в виде:
(6)
$\begin{gathered} {{\Delta }}P = El \times {{10}^{5}} \times 0.14E_{{kt}}^{{0.43}} \times \\ \times \,\,{{\left( {\frac{1}{{1000 \times \sqrt {H_{{eff}}^{2} + {{{\left( {y - {{Y}_{0}}} \right)}}^{2}} + {{x}^{2}}} }}} \right)}^{{1.56}}}, \\ \end{gathered} $Соотношение для Y0 зависит от плотности импактора и растет с увеличением эффективной высоты и уменьшением угла входа:
(7)
${{Y}_{0}} = \left\{ {\begin{array}{*{20}{c}} {0.57H_{{eff}}^{{1.3}}{{{({\text{ctg}}{\kern 1pt} {\kern 1pt} \alpha )}}^{{0.71}}},\,\,\,\,~~\rho = 1000\,\,~{{{\text{кг}}} \mathord{\left/ {\vphantom {{{\text{кг}}} {{{{\text{м}}}^{{\text{3}}}}}}} \right. \kern-0em} {{{{\text{м}}}^{{\text{3}}}}}}} \\ {0.69H_{{eff}}^{{1.3}}{{{({\text{ctg}}{\kern 1pt} {\kern 1pt} \alpha )}}^{{0.71}}},\,\,\,\,~~\rho = 3320\,\,~{{{\text{кг}}} \mathord{\left/ {\vphantom {{{\text{кг}}} {{{{\text{м}}}^{{\text{3}}}}}}} \right. \kern-0em} {{{{\text{м}}}^{{\text{3}}}}}}~} \end{array}} \right..$Разница между приближением (7) и результатами численного моделирования в среднем составляет около 2 км, большие значения могут наблюдаться в случае малых углов входа α.
Эллиптичность $El$ описывает гетерогенность распределения избыточного давления. Она является поправкой, превращающей круговую симметрию в эллиптическую для двух полупространств. Поля распределений обладают симметрией в областях до и после точки максимального энерговыделения. Введение${\text{\;}}El$ позволяет увеличить или уменьшить оценку избыточного давления в точке пространства, тем самым пододвигая или отдаляя изолинию постоянного уровня к или от точки центра распределения.
Область под траекторией (до точки максимального энерговыделения) контролируется параметром Elb (0 <φ < π), а область вперед по направлению траектории после точки максимального энерговыделения контролируется параметром Elf (–π < φ <0), где φ полярный угол ($\phi = {\text{arctg}}\left( {{y \mathord{\left/ {\vphantom {y x}} \right. \kern-0em} x}} \right)$):
(8)
$\begin{gathered} E{{l}_{b}} = \frac{{nb}}{{\sqrt {{{n}^{2}}{\text{sin}}{{{\left[ \phi \right]}}^{2}} + {{{\text{b}}}^{2}}{\text{cos}}{{{\left[ \phi \right]}}^{2}}} }},~ \\ E{{l}_{f}} = \frac{{nf}}{{\sqrt {{{n}^{2}}{\text{sin}}{{{\left[ \phi \right]}}^{2}} + {{f}^{2}}{\text{cos}}{{{\left[ \phi \right]}}^{2}}} }}. \\ \end{gathered} $Для кратерообразующих ударов и метеорных взрывов соотношения получались отдельно. Коэффициенты n, b и f, входящие в уравнение (8), для случаев метеорных взрывов обозначаются индексом “ab” и имеют следующий вид:
(9)
${{n}_{{ab}}} = 1 - 0.002H_{{eff}}^{{1.16}}{{\left( {\cos \alpha } \right)}^{{2.7}}}\ln \left( {~{{E}_{{kt}}}} \right),$(10)
${{b}_{{ab}}} = \left\{ {\begin{array}{*{20}{c}} {0.19{{D}^{{0.4}}}{{{\left( {\sin \alpha } \right)}}^{{0.32}}},~\,\,\,\,~\rho = 1000\,\,~{{{\text{кг}}} \mathord{\left/ {\vphantom {{{\text{кг}}} {{{{\text{м}}}^{{\text{3}}}}}}} \right. \kern-0em} {{{{\text{м}}}^{{\text{3}}}}}}} \\ {0.39{{V}^{{0.26}}}{{{\left( {\sin \alpha } \right)}}^{{0.26}}},\,\,\,\,~\rho = 3320\,\,~~{{{\text{кг}}} \mathord{\left/ {\vphantom {{{\text{кг}}} {{{{\text{м}}}^{{\text{3}}}}}}} \right. \kern-0em} {{{{\text{м}}}^{{\text{3}}}}}}} \end{array}} \right.,$(11)
${{f}_{{ab}}}\, = \,\left\{ \begin{gathered} 0.5{{D}^{{0.1}}},\,\,\,\,~\rho = 1000\,\,~{{{\text{кг}}} \mathord{\left/ {\vphantom {{{\text{кг}}} {{{{\text{м}}}^{{\text{3}}}}}}} \right. \kern-0em} {{{{\text{м}}}^{{\text{3}}}}}} \hfill \\ 0.23\,\, + \,\,0.25H_{{{\text{eff}}}}^{{{\text{0}}{\text{.37}}}}{\kern 1pt} {{\left( {\sin \alpha } \right)}^{{0.52}}},\,\rho \, = \,3320~\,\,{{{\text{кг}}} \mathord{\left/ {\vphantom {{{\text{кг}}} {{{{\text{м}}}^{{\text{3}}}}}}} \right. \kern-0em} {{{{\text{м}}}^{{\text{3}}}}}} \hfill \\ \end{gathered} \right..$Для кратерообразующих ударов коэффициенты n, b и f обозначаются индексом “с” и имеют следующий вид:
(12)
${{n}_{c}} = \left\{ {\begin{array}{*{20}{c}} {1.1{{{\left( {\sin {{\alpha }}} \right)}}^{{ - 0.34}}},~\,\,\,\,~\rho = 1000~\,\,{{{\text{кг}}} \mathord{\left/ {\vphantom {{{\text{кг}}} {{{{\text{м}}}^{3}}}}} \right. \kern-0em} {{{{\text{м}}}^{3}}}}} \\ {0.56{{{\left( {\sin {{\alpha }}} \right)}}^{{ - 0.27}}},\,\,\,\,~~\rho = 3320~\,\,{{{\text{кг}}} \mathord{\left/ {\vphantom {{{\text{кг}}} {{{{\text{м}}}^{3}}}}} \right. \kern-0em} {{{{\text{м}}}^{3}}}}} \end{array}} \right.,$(13)
${{b}_{c}} = \left\{ {\begin{array}{*{20}{c}} {1.7{{{\left( {\sin {{\alpha }}} \right)}}^{{ - 0.54}}},\,\,\,\,~\rho = 1000\,\,~{{{\text{кг}}} \mathord{\left/ {\vphantom {{{\text{кг}}} {{{{\text{м}}}^{{\text{3}}}}}}} \right. \kern-0em} {{{{\text{м}}}^{{\text{3}}}}}}} \\ {0.84{{{\left( {\sin {{\alpha }}} \right)}}^{{ - 1.2}}},\,\,\,\,~\rho = 3320~\,\,{{{\text{кг}}} \mathord{\left/ {\vphantom {{{\text{кг}}} {{{{\text{м}}}^{{\text{3}}}}}}} \right. \kern-0em} {{{{\text{м}}}^{{\text{3}}}}}}} \end{array}} \right.,$(14)
${{f}_{c}} = \left\{ {\begin{array}{*{20}{c}} {0.65 + 0.064{{{\left( {\sin {{\alpha }}} \right)}}^{{ - 2.8}}},~\,\,~\rho = 1000~\,\,{{{\text{кг}}} \mathord{\left/ {\vphantom {{{\text{кг}}} {{{{\text{м}}}^{{\text{3}}}}}}} \right. \kern-0em} {{{{\text{м}}}^{{\text{3}}}}}}} \\ {0.44 + 0.0083{{{\left( {\sin {{\alpha }}} \right)}}^{{ - 3.8}}},\,\,~\rho = 3320\,\,{{{\text{кг}}} \mathord{\left/ {\vphantom {{{\text{кг}}} {{{{\text{м}}}^{{\text{3}}}}}}} \right. \kern-0em} {{{{\text{м}}}^{{\text{3}}}}}}} \end{array}} \right..$Стоит отметить наличие некоторых недостатков в соотношениях. Например, для вертикального удара не получится круговой симметрии, и описываемые параметры выглядят не одинаково. Процесс их получения был следующим. Составлялся список из параметров основных ударов (плотность ρ, диаметр D, скорость V, кроме угла входа α) и их производных характеристик (кинетическая энергия $~{{E}_{{kt}}}$ и ее логарифм, эффективная высота Heff). Угол входа α был представлен набором тригонометрических функций. Далее методом перебора рассматривались все возможные комбинации этих параметров в некоторых степенях, с прибавлением констант. Эти параметры подбирались так, чтобы минимизировать относительную невязку на всем поле распределения избыточного давления. Таким образом, соотношения $El$ также несут в себе общую поправку на аппроксимацию среднего профиля распределения для каждого конкретного варианта, а не только его пространственную неоднородность. Эта процедура применялась как в полной выборке всех вариантов метеорных взрывов и кратерообразующих ударов, но также и только для комет, только астероидов и специфических выборок по диапазонам Heff (равные 0 км, от 0–5 км, больше 5 км, больше 10 км). Из отсортированных по значению ошибок всех этих тестов выбирались наиболее простые и лучше всего описывающие поля комбинации. Можно видеть, что разделение на группы по плотности космического тела оказалось самым лучшим.
На рис. 4 показано сравнение распределений максимального избыточного давления ΔP для 30-метрового астероида со скоростью 20 км/с, углом входа 45°, полученных по аппроксимации (6) и в результате численного моделирования. Численное моделирование и аппроксимационные соотношения оценивают избыточное давление как 0.01–0.1 во всей области. Соответствующие относительные погрешности в большей части рассматриваемой области не превышают 33%, в небольшой части области достигают 100% и лишь на небольших участках вблизи траектории, в области низких значений, оказывается немного хуже.
Для кратерообразующего удара 1-км кометы соответствующие распределения приведены в нижнем ряду рис. 4. Как и для случая метеорного взрыва в большей части области различие в давлениях не превышает 100%, хуже всего описывается область относительно невысоких давлений под траекторией, где ошибка достигает 800%. В целом аппроксимация (6) удовлетворительно описывает результаты численного моделирования. С учетом приближений, принятых в расчетной модели, и неопределенности информации об ударнике, такая точность кажется приемлемой.
В упоминавшемся ранее интернет-калькуляторе [Collins et al., 2005; 2017] используются аппроксимации избыточного давления, построенные по результатам ядерных испытаний. Предложенное в работе [Collins et al., 2017] соотношение для метеорных взрывов выглядит следующим образом:
(15)
$\begin{gathered} p\left( r \right) = 3.14 \times {{10}^{{11}}} \times {{\left( {{{r}^{2}} + Z_{b}^{2}} \right)}^{{ - \frac{{2.6}}{2}}}} + {{\left( {{{r}^{2}} + Z_{b}^{2}} \right)}^{{ - \frac{{1.13}}{2}}}}, \\ {{r}^{2}} = {{x}^{2}} + {{y}^{2}}, \\ \end{gathered} $Сравнение предсказываемых значений по уравнению (15) и результатов численного моделирования приведено на рис. 5. В качестве эффективной высоты Zb использовалось модальное значение (Zb 50%), приведенное в работе [Collins et al., 2017]. Для воздушного взрыва (верхний ряд рис. 5) предложенное соотношение (15) дает удовлетворительное описание, на большей части распределения относительная ошибка не превышает 100%, область с ошибкой менее 30% оказывается меньше, чем на рис. 4, вероятно из-за отсутствия пространственной неоднородности в выражении (15), предполагающем сферическую симметрию. Следует отметить, что для метеорных взрывов калькулятор Impact Earth выдает диапазон значений давления (различающихся в 2 раза) в центральной области с учетом неопределенности в величине Zb. Для данного набора параметров ударника (D, V, α, ρ) нижняя граница диапазона давлений близка к приведенной нами оценке и удовлетворительно описывает результаты численного моделирования, верхняя граница диапазона приводит к заметному росту погрешностей в центральной области.
Для кратерообразующих ударов соотношения, предложенные в работе [Collins et al., 2005], не модифицировались и по-прежнему используются в калькуляторе. Сравнение результатов для удара 1-км кометы с результатами численного моделирования приведено на рис. 5 в нижнем ряду. Ошибки в оценках выше и составляют 33–100% для области вперед по траектории, а в области под траекторией значительно выше – до 5 раз и выше.
В работе [Collins et al., 2017] и в (15) используется двухкомпонентная зависимость для избыточного давления от расстояния. В ближней зоне ${{\Delta }}P\sim {{~E_{{kt}}^{{0.86}}} \mathord{\left/ {\vphantom {{~E_{{kt}}^{{0.86}}} {{{r}^{{2.6}}}}}} \right. \kern-0em} {{{r}^{{2.6}}}}}$, а на больших расстояниях . Для точечного источника формула Садовского предлагает трехкомпонентную зависимость, которую можно свести к виду:
(16)
${{\Delta }}P\sim \frac{{A{{E}_{{kt}}}}}{{{{r}^{3}}}} + \frac{{BE_{{kt}}^{{0.66}}}}{{{{r}^{2}}}} + \frac{{CE_{{kt}}^{{0.33}}}}{r},$Для цилиндрического источника выражение для избыточного давления предложено в работе [Цикулин, 1969]. Оно имеет двухкомпонентный вид, и его можно свести к следующему виду:
(17)
${{\Delta }}P\sim \frac{{G{\text{\;}}{{E}_{{kt}}}}}{{{{r}^{2}}}} + \frac{{FE_{{kt}}^{{0.375}}}}{{{{r}^{{0.75}}}}},$К опасным факторам ударной волны относится не только избыточное давление, но и скорость ветра за фронтом. Уровень избыточного давления и скорость ветра за фронтом ударной волны связаны уравнениями Гюгонио. Следовательно, распределение максимальной массовой скорости ветра может быть получено из распределения избыточного давления:
(18)
${{V}_{{{\text{max}}}}} = \frac{{330}}{{{\gamma }}}\frac{{{{\Delta }}P}}{{{{P}_{0}}}}{{\left( {1 + \frac{{{{\gamma }} + 1}}{{2{{\gamma }}}}\frac{{{{\Delta }}P}}{{{{P}_{0}}}}} \right)}^{{ - 1/2}}},$Распределения ветра, полученные по этой формуле с использованием выражения (6), сравнивались с численными расчетами, в которых определялась максимальная горизонтальная составляющая скорости. Значительными можно считать скорости ветра в диапазоне от 20 м/с (начинают падать сухие ветки с деревьев), в этом диапазоне средняя относительная ошибка нашей модели (18) составляет около 2–4.
Строго говоря, формула (18) хорошо работает на больших расстояниях, где фронт ударной волны становится вертикальным. Особенно большую погрешность формула (18) дает в зоне “телеграфного леса” для воздушных взрывов, где почти все деревья остаются стоять вертикально с оторванными кронами, в этой области скорость газа за фронтом ударной волны направлена почти вертикально. Аппроксимационные соотношения были применены к двум наиболее известным событиям – Челябинскому астероиду 2013 г. (рис. 6) и Тунгусскому событию 1908 г. Максимальное избыточное давление составляет около 0.03 атм (3 кПа) для Челябинска и 0.3 атм (30 кПа) для Тунгуски.
Сравнение полученных распределений избыточного давления для астероида Челябинск с данными о разрушениях и травмах показывает неплохое согласие [Kartashova et al., 2018]. Сравнение результатов численного моделирования [Shuvalov et al., 2017b] с оценками масштабных соотношений приведено на рис. 6 и демонстрирует хорошее согласие. Следует отметить, что моделируемые области, полученные с помощью аппроксимационного соотношения (6), уже в направлении, перпендикулярном траектории, чем полученные с помощью численного моделирования. Это показатель того, что параметр El недооценен для косых ударов относительно небольших тел. Кроме того, значения относительного избыточного давления в этом случае низкие (на нижней границе значений относительного давления), и их точность ниже, чем в других случаях.
В случае Тунгусского события можно рассматривать широкий диапазон параметров импактора (т.е. кинетическая энергия от 10 до 15 Мт, плотность космического объекта от 1000 до 3320 кг/м3; угол входа в диапазоне 20°–45° и диапазон скоростей 15–40 км/с). Сравнение результатов аппроксимации (6) с показаниями очевидцев демонстрирует удовлетворительное согласие: избыточные давления на расстояниях 300–400 км оказываются достаточными для повреждения остекления, что согласуется с недавно опубликованными данными о разбитых стеклах [Jenniskens et al., 2019].
КАЛЬКУЛЯТОР ЭФФЕКТОВ ПАДЕНИЙ КОСМИЧЕСКИХ ТЕЛ
В интернете доступны несколько веб-калькуляторов для оценки последствий импактов космических тел. Самым популярным из них является ImpactEarth (https://impact.ese.ic.ac.uk/ImpactEarth/), который подробно описан в работах [Collins et al., 2005; 2017]. Как уже отмечалось выше, этот калькулятор основан на упрощенных моделях и данных о ядерных взрывах. Можно использовать оценки некоторых эффектов; другие оцениваются с существенными погрешностями, некоторые эффекты вообще не рассчитываются. Например, этот калькулятор не оценивает излучения, возникающего при пролете и разрушении космического тела в атмосфере. Другой калькулятор, доступный в интернете (http://simulator. down2earth.eu), фактически представляет собой тот же калькулятор.
Наш онлайн-калькулятор доступен на http://AsteroidHazard.pro и реализует представленные аппроксимационные соотношения для эффектов ударной волны, которые основаны на серийном численном моделировании падения космических тел. Помимо воздействия ударных волн, калькулятор также дает оценку радиационных эффектов [Svetsov, Shuvalov, 2018; 2019; Podobnaya et al., 2018], параметров образования кратеров и выбросов [Artemieva et al., 2019], сейсмических эффектов [Svetsov et al., 2017] и атмосферных возмущений в верхних слоях атмосферы [Shuvalov, Khazins, 2018].
ЗАКЛЮЧЕНИЕ
Анализ результатов численного моделирования падений космических тел в широком диапазоне их свойств позволил предложить простые соотношения, которые дают возможность оценивать избыточное давление в любой точке пространства на поверхности Земли для падающего тела с заданными свойствами и могут использоваться для прогнозирования последствий удара. Предложенные соотношения зависят только от параметров импактора (размер, плотность, скорость и угол входа в атмосферу), не требуют длительных вычислений и представляют собой относительно простые алгебраические выражения.
Все аппроксимации, описанные выше, реализованы в онлайн-калькуляторе (http://AsteroidHazard.pro), который можно использовать для быстрой оценки опасных последствий. Кроме того, в калькуляторе возможен расчет других эффектов (тепловые потоки, размер кратера, выбросы из кратера, сейсмическое воздействие и ионосферные возмущения).
Список литературы
Авилова И.В., Биберман Л.М., Воробьев В.С., Замалин В.М., Кобзев Г.А., Лагарьков А.Н., Мнацаканян А.Х., Норман Г.Э. Оптические свойства горячего воздуха. М.: Наука. 1970. 320 с.
Акимов В.А., Глазачев Д.О., Емельяненко В.В., Краминцев А.П., Нароенков С.А., Попова О.П., Пучков В.А., Рыхлова Л.В., Савельев М.И., Светцов В.В., Трубецкая И.А., Шувалов В.В., Шугаров А.С., Шустов Б.М. Астероидно-кометная опасность: стратегия противодействия / Под общ. ред. Пучкова В.А. М.: ФГБУ ВНИИ ГОЧС (ФЦ). 2015. 272 с.
Садовский М.А. Механическое действие воздушных ударных волн по данным экспериментальных исследований // Физика взрыва. М. : АН СССР. 1952. № 1. С. 20–110.
Шувалов В.В., Трубецкая И.А. Влияние внутреннего трения на деформацию разрушенного метеороида // Астрономический вестник. 2010. Т. 44. № 2. С. 117–122.
Хазинс В.М., Шувалов В.В., Светцов В.В. Сейсмическая эффективность ударов космических тел // Астрономический вестник: Исследования Солнечной системы. 2018. Т. 52. № 6. С. 566–576.
Цикулин, М.А. Ударные волны при движении в атмосфере крупных метеоритных тел. М.: 1969. Наука. 86 с.
Artemieva N., Shuvalov V. Atmospheric shock waves after impacts of cosmic bodies up to 1000 m in diameter // Meteorit. Planet. Sci. 2019. V. 54. № 3. P. 592–608.
Artemieva N., Shuvalov V.V., Khazins V.M. Upper atmosphere effects after the entry of small cosmic bodies: Dust trains, plumes, and atmospheric disturbances // Icarus. 2019. V. 327. P. 60–71.
Bronikowska M., Artemieva N.A., Wünnemann K. Reconstruction of the Morasko meteoroid impact ‑ Insight from numerical modeling // Meteoritics & Planetary Science. 2017. V. 52. № 8. P. 1704–1721.
Brown P.G., Assink J.D., Astiz L., Blaauw R., Boslough M.B., Borovička J., Brachet N., Brown D., Campbell-Brown M., Ceranna L., Cooke W., De Groot-Hedlin C., Drob D.P., Edwards W., Evers L.G., Garces M., Gill J., Hedlin M., Kingery A., Laske G., Le Pichon A., Mialle P., Moser D.E., Saffer A., Silber E., Smets P., Spalding R.E., Spurný P., Tagliaferri E., Uren D., Weryk R.J., Whitaker R., Krzeminski Z. A 500-kiloton airburst over Chelyabinsk and an enhanced hazard from small impactors // Nature. 2013. V. 503. № 7475. P. 238–241.
Collins G.S., Melosh H.J., Marcus R.A. Earth Impact Effects Program: A Web-based computer program for calculating the regional environmental consequences of a meteoroid impact on Earth // Meteorit. Planet. Sci. 2005. V. 40. № 6. P. 817–840.
Collins G.S., Lynch E., McAdam R., Davison T.M. A numerical assessment of simple airblast models of impact airbursts // Meteorit. Planet. Sci. 2017. V. 52. № 8. P. 1542–1560.
Gi N., Brown P., Aftosmis M. The frequency of window damage caused by bolide airbursts: A quarter century case study // Meteorit. Planet. Sci. 2018. V. 47. P. 1413–1431.
Glasstone S., Dolan P.J. The Effects of Nuclear Weapons. 3rd edition. Washington, D.C.: US Department of Defense. 1977. 653 p.
Ivanov B.A., Deniem D., Neukum G. Implementation of dynamic strength models into 2D hydrocodes: Applications for atmospheric break up and impact cratering // Int. J. Impact Eng. 1997. V. 20. P. 411–430.
Jenniskens P., Popova O.P., Glazachev D.O., Podobnaya E.D., Kartashova A.P. Tunguska eyewitness accounts, injuries, and casualties // Icarus. 2019. V. 327. P. 4–18.
Kartashova A.P., Popova O.P., Glazachev D.O., Jenniskens P., Emel’yanenko V.V., Podobnaya E.D., Skripnik A.Y. Study of injuries from the Chelyabinsk airburst event // Planet. Space Sci. 2018. V. 160. P. 107–114.
Krinov E.L. New Studies of the Sikhote-Alin Iron Meteorite Shower // Meteoritics. 1971. V. 6. № 3. P. 127–138.
Nemtchinov I.V., Popova O.P. An Analysis of the 1947 Sikhote-Alin Event and a Comparison with the Phenomenon of February 1, 1994 // Sol. Syst. Res. 1997. V. 31. № 5. P. 12.
Podobnaya E.D., Glazachev D.O., Popova O.P., Svetsov V.V., Shuvalov V.V. Scaling Relations for Radiation Effects due to Impacts of Large Cosmic Objects // Meteorit. Planet. Sci. 2018. V. 53. P. 6030.
Popova O.P., Jenniskens P., Emel’yanenko V., Kartashova A., Biryukov E., Khaibrakhmanov S., Shuvalov V., Rybnov Y., Dudorov A., Grokhovsky V.I., Badyukov D.D., Yin Q.Z., Gural P.S., Albers J., Granvik M., Evers L.G., Kuiper J., Kharlamov V., Solovyov A., Rusakov Y.S., Korotkiy S., Serdyuk I., Korochantsev A.V., Larionov M.Y., Glazachev D., Mayer A.E., Gisler G., Gladkovsky S.V., Wimpenny J., Sanborn M.E., Yamakawa A., Verosub K.L., Rowland D.J., Roeske S., Botto N.W., Friedrich J.M., Zolensky M.E., Le L., Ross D., Ziegler K., Nakamura T., Ahn I., Lee J.I., Zhou Q., Li X.H., Li Q.L., Liu Y., Tang G.Q., Hiroi T., Sears D., Weinstein I.A, Vokhmintsev A.S., Ishchenko A.V., Schmitt-Kopplin P., Hertkorn N., Nagao K., Haba M.K., Komatsu M., Mikouchi T. Chelyabinsk Airburst, Damage Assessment, Meteorite Recovery, and Characterization // Science. 2013. V. 342. № 6162. P. 1069–1073.
Sánchez P., Scheeres D.J. The strength of regolith and rubble pile asteroids // Meteorit. Planet. Sci. 2014. V. 49. № 5. P. 788–811.
Shuvalov V.V. Multi-dimensional hydrodynamic code SOVA for interfacial flows: Application to the thermal layer effect // Shock Waves. 1999. V. 9. № 6. P. 381–390.
Shuvalov V.V., Artemieva N.A. Numerical modeling of Tunguska-like impacts // Planet. Space Sci. 2002. V. 50. № 2. P. 81–192.
Shuvalov V.V., Trubetskaya I.A. Aerial bursts in the terrestrial atmosphere // Sol. Syst. Res. 2007. V. 41. № 3. P. 220–230.
Shuvalov V.V., Khazins V.M. Numerical Simulation of Ionospheric Disturbances Generated by the Chelyabinsk and Tunguska Space Body Impacts // Sol. Syst. Res. 2018. V. 52. № 2. P. 129–138.
Shuvalov V.V., Svettsov V.V., Trubetskaya I. A. An estimate for the size of the area of damage on the Earth’s surface after impacts of 10–300-m asteroids // Solar System Research. 2013. V. 47. № 4. P. 260–267.
Shuvalov V.V., Popova O.P., Svettsov V.V., Trubetskaya I.A., Glazachev D.O. Determination of the Height of the “Meteoric Explosion” // Sol. Syst. Res. 2016. V. 50. № 1. P. 1–12.
Shuvalov V.V., Svettsov V.V., Artem’eva N.A., Trubetskaya I.A., Popova O.P., Glazachev D.O. Asteroid Apophis: Evaluating the impact hazards of such bodies // Sol. Syst. Res. 2017a. V. 51. № 1. P. 44–58.
Shuvalov V.V., Svetsov V.V., Popova O.P., Glazachev D.O. Numerical model of the Chelyabinsk meteoroid as a strengthless object // Planet. Space Sci. 2017b. V. 147. P. 38–47.
Svetsov V.V., Shuvalov V.V. Thermal radiation and luminous efficiency of superbolides // Earth Planet. Sci. Lett. 2018. V. 503. P. 10–16.
Svetsov V.V., Shuvalov V.V. Thermal radiation from impact plumes // Meteorit. Planet. Sci. 2019. V. 54. P. 126–141.
Svetsov V.V., Nemtchinov I.V., Teterev A.V. Disintegration of large meteoroids in Earth’s atmosphere: Theoretical models // Icarus. 1995. V. 116. P. 131–153.
Svetsov V.V., Artemieva N.A., Shuvalov V.V. Seismic Efficiency of Meteor Airbursts // Dokl. Earth Sci. 2017. V. 475. № 2. P. 935–938.
Дополнительные материалы отсутствуют.