Физика Земли, 2021, № 3, стр. 63-81

ВВЕДЕНИЕ

Около 40% всех ощущаемых землетрясений, в той или иной мере обусловленных человеческой деятельностью, связаны с извлечением твердых полезных ископаемых, а еще 15% – с добычей нефти и газа [Foulger et al., 2018]. В результате существенного увеличения плотности и качества инструментальных наблюдений на сегодняшний день in situ регистрируется множество сейсмических событий, вплоть до столь малых магнитуд, как ${{M}_{w}} = - 4.1$ [Kwiatek et al., 2011]. Параметры микро- и наносейсмичности изучаются во многих исследованиях, связанных, прежде всего, с горными работами [Gibowicz, Kijko, 1994; Oye et al., 2005; Yamada et al., 2007; Wesseloo, 2018; Ma et al., 2018; и др.].

Исследования, проведенные в горных массивах с разным напряженным состоянием и тектонической обстановкой, показали, что можно обнаружить два типа сейсмического отклика на проведение горных работ [Кочарян и др., 2018; De Santis et al., 2020]. К первому относятся микрособытия, непосредственно вызванные добычей. Их источники по пространственным и временным параметрам приурочены непосредственно к фронту ведения горнопроходческих работ. Расположение очагов этих событий отмечает либо зоны концентрации напряжений, либо, напротив, разгруженные области массива [Gibowicz, Kijko, 1994; Urbancic, Trifu, 2000; Young et al., 2004; Orlecka-Sikora et al., 2012]. Второй тип – инициированные подвижки по существующим тектоническим нарушениям [Кочарян и др., 2018; Lu et al., 2018; 2019; Sainoki, 2013]. Такие события могут быть удалены от места проведения работ на несколько километров, могут происходить с задержкой по времени по отношению к производству взрывов. Именно крупные динамические подвижки по разломам – события второго типа – приносят наибольший ущерб [Ловчиков, 2013; Lu et al., 2018].

Анализ очаговых параметров землетрясений, инициированных горными работами, обычно проводится на основе подходов, развитых в “большой” сейсмологии и дает довольно много дополнительной информации для изучения масштабного эффекта землетрясений [McGarr, 1999; Kanamori, Heaton, 2000; Boettcher et al., 2009; Kwiatek et al., 2011]. При этом остается открытым вопрос, насколько изучение параметров слабой сейсмичности помогает интерпретировать данные, получаемые при регистрации средних и крупных коровых землетрясений. В рамках предположения о том, что процесс разрыва самоподобен [Aki, 1967], такие параметры очага, как приведенная сейсмическая энергия $e = \frac{{{{E}_{s}}}}{{{{M}_{0}}}}$ (${{E}_{s}}$ – энергия, излученная при землетрясении, ${{M}_{0}}$ – сейсмический момент) и кажущееся напряжение ${{\sigma }_{a}}$, должны быть одинаковы для землетрясений любого размера. Хотя исследованиями закономерностей сейсмического масштабирования занимались и занимаются многие исследовательские группы, геофизическое сообщество пока не пришло к единому мнению [Walter et al., 2006]. Если в одних случаях авторы подтверждают гипотезу самоподобия и постоянства кажущихся напряжений [Kanamori, Anderson, 1975; McGarr, 1999; Ide et al., 2003; Baltay et al., 2010; и др.], то другие исследования, напротив, свидетельствуют в пользу зависимости величины приведенной сейсмической энергии от размера землетрясения [Abercrombie, 1995; Kwiatek et al., 2011; Кочарян и др., 2016; и др.]. В любом случае, наблюдается очень большой разброс в рассчитываемых величинах $e$. Так, например, в работе [Кочарян и др., 2016] был проанализирован большой объем наблюдательных данных – более 1500 событий в диапазоне моментных магнитуд от ${{M}_{w}}$ = –3.5 до ${{M}_{w}}$ = 9.2. Анализ был ограничен землетрясениями, очаг которых располагался на глубине не более 50 км. Были привлечены только надежные оценки сейсмического момента и излученной энергии, а при анализе очаговых параметров землетрясений малых магнитуд либо использовались данные, полученные с использованием высокочастотной регистрирующей аппаратуры, либо проводилась коррекция результатов по методу, предложенному [Ide, Beroza, 2001]. Анализ составленной базы данных показал, что при изменении сейсмического момента на 19 порядков – в диапазоне от 10³ до 10²² Н ⋅ м – подавляющее большинство значений параметра $e$ лежит в диапазоне $e$ ~ 10^–6–10^–3, а среднее значение составляет $e$ = 2.82 × 10⁵. При рассмотрении всего комплекса данных явной зависимости отношения от масштаба землетрясения не обнаруживается: дисперсия очень велика и вариации приведенной энергии при одном и том же значении сейсмического момента превышает 3 порядка. Однако, судя по данным [Кочарян и др., 2016], средняя приведенная энергия неодинакова в разных диапазонах магнитуд. Многими авторами отмечается, что для каждого отдельного подмножества данных, особенно для слабых землетрясений, часто наблюдаются очевидные отклонения от пропорциональной связи сейсмической энергии и сейсмического момента [Walter et al., 2006]. И если авторы ряда исследований последовательностей слабых природных землетрясений склоняются к соотношению ${{E}_{s}}\sim {{M}_{0}}$ [Ide, Beroza, 2001; Matsuzawa et al., 2004], то для рядов наблюдений индуцированной сейсмичности наклон регрессионных зависимостей ${{E}_{s}}({{M}_{0}})$ иногда сильно отличается от единицы [Маловичко, Маловичко, 2010; Gibowicz et al., 1991]. Это может означать, что разломная зона или участок массива, к которым приурочено расположение очагов группы слабых землетрясений, определяют свои закономерности масштабирования, определяемые местными условиями [Walter et al., 2006; Беседина и др., 2015; Кочарян, 2016]. При этом совместное рассмотрение большого количества рядов данных, относящихся к разным участкам, создает видимость постоянства величины средней приведенной энергии [Walter et al.; 2006; Kwiatek et al., 2011; Кочарян и др., 2016].

Физические причины вариации эффективности излучения сейсмических колебаний с ростом размера землетрясения до сих пор обсуждаются. Среди них рассматриваются зависимость модуля деформации массива горных пород и жесткости нарушений сплошности от масштаба [Беседина и др., 2013]. Помимо этого, можно выделить также два фактора, приводящих к возникновению зависимости излучательной эффективности сейсмического источника от масштаба события и не связанных с физикой очага. Во-первых, это эффект влияния ограниченной частотной характеристики регистрирующей аппаратуры [Ide, Beroza, 2001]. При использовании недостаточно высокочастотной аппаратуры этот эффект оказывается весьма значительным для событий малой энергии. Второй эффект связан с поглощением и рассеянием высоких частот при распространении колебаний. Для слабых событий наблюдается недооценка рассчитанной сейсмической энергии, быстро увеличивающаяся с ростом расстояния [Беседина и др., 2013].

Несмотря на большое количество данных, до настоящего момента нет единого мнения и о причинах большого разброса значений отношения ${{{{E}_{s}}} \mathord{\left/ {\vphantom {{{{E}_{s}}} {{{M}_{0}}}}} \right. \kern-0em} {{{M}_{0}}}}$. Хотя, как правило, величина приведенной излученной энергии для индуцированных землетрясений лежит в диапазоне, характерном для “нормальных” природных землетрясений ${{10}^{{ - 6}}} < e < {{10}^{{ - 4}}}$, в некоторых случаях зарегистрированные величины попадают в диапазон ${{{{E}_{s}}} \mathord{\left/ {\vphantom {{{{E}_{s}}} {{{M}_{0}}}}} \right. \kern-0em} {{{M}_{0}}}} < {{10}^{{ - 6}}}$, т.е. эффективность подвижки как источника сейсмических колебаний оказывается крайне низкой [Oye, 2005; Маловичко, Маловичко, 2010; Kwiatek et al., 2011]. Между тем ущерб, наносимый землетрясением, определяется именно величиной излученной энергии, так что важно изучать причины наблюдаемых эффектов.

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

ДАННЫЕ И МЕТОДЫ

Коробковское железорудное месторождение

Коробковское железорудное месторождение сложено из 8 крутопадающих залежей железистых кварцитов длиной 380–3500 м и горизонтальной мощностью 50–1050 м [Адушкин и др., 2017б]. Согласно результатам измерений, горизонтальные напряжения на глубине 300 м от поверхности составляют 13–15 МПа, превышая примерно в 1.5 раза вертикальные [Адушкин и др., 2017б; Кочарян и др., 2018]. Существенное превышение величины горизонтальных напряжений над вертикальными стало основанием для того, чтобы надзорные органы объявили массив удароопасным, хотя до настоящего времени горных ударов здесь не наблюдалось. Поскольку для небольших глубин, на которых разрабатывается Коробковское месторождение, уровень напряжений достаточно высок, было решено организовать временные сейсмические наблюдения на участке, который считается наиболее опасным при проведении горных работ.

Разработка Коробковского железорудного месторождения ведется подземным способом. Средние значения физико-механических характеристик разрабатываемых кварцитов (плотность $\rho $, модуль Юнга $E$, коэффициент Пуассона ν, модуль сдвига $\mu $ и скорость распространения продольных и поперечных волн ${{C}_{\alpha }}$ и ${{C}_{\beta }}$ соответственно), отобранных в горных выработках Коробковского месторождения, приведены в табл. 1 по [Кочарян и др., 2018].

Таблица 1.  

Физико-механические характеристики разрабатываемых кварцитов

$\rho $, кг/м3	$E$, ГПа	ν	$\mu $, ГПа	${{C}_{\alpha }}$, м/с	${{C}_{\beta }}$, м/с
3400–3700	0.65–1.16	0.13–0.33	29–48	4400–6300	2900–3600

На шахте применяется этажно-камерная система разработки. Камеры находятся на глубине 245–285 м от дневной поверхности. При подготовке взрывов участок породы шириной 10–15 м внутри камеры разбуривается рядами вертикальных скважин длиной 50–55 м и диаметром 105 мм. При каждом массовом взрыве с интервалом 30–50 с взрываются 3–4 камеры, расположенные на расстоянии около 1 км друг от друга. Применяется короткозамедленное взрывание. Суммарная масса ВВ в камере составляет, в среднем, 10–20 тонн. Внутри каждой камеры последовательно, через определенный интервал времени (~20 мс) взрываются группы скважин, составляющих одну ступень замедления массой 0.5–1.5 т ВВ.

Пример записи вертикальной скорости смещения грунта на свободной поверхности непосредственно над одной из камер показан на рис. 1. При этом взрыве суммарная масса ВВ в камере была 10.3 т, а максимальная масса ВВ в одной ступени замедления составляла 0.518 т. В начале взрыва интервалы замедления между ступенями составляли 20 мс, а в конце 500 мс, так что в конце записи легко выделяются волны от отдельных ступеней замедления. Длительность колебаний от взрыва одной ступени составляет около 200 мс. Суммарная длительность колебаний не превышает 5–6 с.

Рис. 1.

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

Величина максимальной массовой скорости (PGV) при короткозамедленном взрыве определяется суммарным весом заряда в группе замедления, расстоянием до взрыва, а на малых расстояниях еще и геометрией заряда. На достаточно большом удалении от взрыва для оценки воздействия на массив можно воспользоваться эмпирической зависимостью, построенной нами для короткозамедленных взрывов на шахте им. Губкина комбината КМАруда [Кочарян и др., 2019]:

(1)

${{V}_{m}} = 4110 \times {{\left( {\frac{{{{q}^{{1/3}}}}}{R}} \right)}^{{1.73}}},$

где: ${{V}_{m}}$ – измеряется в мм/с; $R$ – эпицентральное расстояние в метрах, а максимальное количество ВВ в одной ступени $q$ – в кг.

Зависимость величины PGV от приведенного расстояния приведена на рис. 2. Можно видеть, что на малых расстояниях, сопоставимых с длиной скважин, где сказывается геометрия заряда, величина PGV лежит заметно ниже соотношения (1).

Рис. 2.

Зависимость максимальной массовой скорости волны от приведенного расстояния до взрыва. Точки – результаты измерений на глубине 300 м. Линия – эмпирическая зависимость (1), построенная по результатам измерений на свободной поверхности параметров колебаний от взрывов на шахте им. Губкина комбината КМАруда [Кочарян и др., 2019].

Целью наблюдений, описываемых в настоящей статье, являлось изучение сеймоакустической реакции напряженного массива на воздействие промышленного взрыва. Для Коробковского месторождения это первый опыт многоканальной регистрации с помощью сейсмических станций, организованных по принципу малоапертурной группы. Предварительно в течение года проводились измерения с одной и двумя измерительными точками [Besedina et al., 2018; Besedina, Ostapchuk, 2018; Besedina, 2019].

Измерения проводились 06.07.2019 г. в течение 16 ч в период прекращения горных и строительно-монтажных работ, которые останавливаются примерно за 2 ч до проведения взрыва и возобновляются в полной мере через сутки. Параметры взрывов приведены в табл. 2.

Таблица 2.  

Параметры взрывов 06.07.2019 г.

№ взрыва	1		2	3
№ камеры	1.1	1.2	2	3
Масса ВВ в камере Q, кг	10 353	2081	5222	4838
Масса ВВ в ступени замедления q, кг	1200	552	864	800
Число ступеней	17	5	6	6
Продолжительность детонации, мс	450	100	125	125
Расстояние до точки Т3(0,0), м	310	100	4380	2500

Сейсмическая группа с апертурой ~130 м, состоящая из 4-х акселерометров Bruel&Kjaer 8306 (рабочая полоса частот 0.08 Гц–1 кГц), была установлена на глубине около 300 м от поверхности. Три датчика T₁–T₃ были установлены по профилю вдоль одного орта, четвертый датчик T₄ располагался в соседнем орте (рис. 3). Акселерометры устанавливались на вертикальной стенке выработки. Регистрация велась на 14-разрядное АЦП E14-440 L-Card с частотой опроса 10 кГц. Чувствительность измерительного канала составляла 7.8 (мкм/с²)/отсчет. Дополнительно на отдельный канал велась запись временной метки для контроля частоты опроса АЦП.

Рис. 3.

Схема Коробковского месторождения в зоне смыкания Главной и Сретенской залежей рудных полей и более подробная схема части исследуемой области, оборудованной акселерометрами T₁–T₄ (треугольники) в разломной зоне шахты. Сплошные линии соответствуют зонам тектонических нарушений. Серым цветом выделены рудные поля. Цифры – номера разрабатываемых камер, подорванных 06.07.2019 г., а их положение показано соответствующими прямоугольниками. Жирные линии внутри прямоугольников соответствуют зонам взрыва. Цифры соответствуют номерам камер в табл. 2. Серые точки – сейсмические события, зарегистрированные после взрывов. Начало координат (0,0) соответствует измерительной точке T₃. Ось Y направлена на север. Большие круги (сплошная линия) соответствуют динамическим деформациям порядка 10^–6 от соответствующего взрыва. Маленький круг (пунктирная линия) соответствует динамической деформации примерно 3 × 10^–6.

Расположение группы было выбрано в окрестности зоны смыкания Главной и Стретенской залежей [Адушкин и др., 2017б], которая характеризуется как зона тектонических нарушений общей мощностью около 100 м. В строении участка выделяются четыре зоны дробления с видимой мощностью 13, 34, 3 и 4 м.

Использовалась система координат с началом в точке расположения датчика Т₃ и направлением оси ординат на север. Координаты измерительных точек в этой системе координат приведены в табл. А1 Приложения. Основные зоны разломов на участке показаны на рис. 3a, 3б.

Обработка данных и локация событий

За два часа до проведения взрывов было зарегистрировано три сейсмических события, которые представляли собой технологические взрывы малой мощности. Сразу после проведения короткозамедленных взрывов, в течение нескольких часов фиксировались волны от множества сейсмических событий разной магнитуды, вне всякого сомнения индуцированных взрывами. Предварительная обработка исходных рядов данных состояла в удалении низкочастотной компоненты и высокочастотных резонансных пиков с помощью полосового фильтра Баттерворта 2-го порядка в диапазоне частот 10 Гц–1 кГц. Обнаружение сигналов проводилось по каждому каналу отдельно. Стандартная методика детектирования сигналов STA/LTA показала много ложных срабатываний из-за постоянно присутствующих техногенных помех. Поэтому в работе использовалась методика выделения, отработанная ранее при проведении предварительных измерений при мониторинге поствзрывной активности в исследуемой области [Besedina et al., 2018]. С помощью быстрого преобразования Фурье в скользящем окне длительностью 0.05 с и перекрытием 0.025 с рассчитывались спектральные характеристики исследуемого временного интервала. Обнаружение событий производилось по превышению выбранным параметром порогового значения, в качестве которого была выбрана сумма спектров ускорения в скользящем окне в заданном диапазоне частот $\sum\nolimits_{{{f}_{1}}}^{{{f}_{2}}} {\left| {A(f)} \right|} $, где $A(f)$ – ускорение, рассчитанное в скользящем окне 0.05 с, a ${{f}_{1}}$ и ${{f}_{2}}$ – частоты, определяющие границы частотного диапазона при обнаружении сигнала. Для измерительных каналов Т₁–Т₃${{f}_{1}}$ = 10 Гц, ${{f}_{2}}$ = 1 кГц. На канале датчика Т₄ присутствуют стационарные помехи на частотах 580–740 Гц, так что этот диапазон частот был исключен из анализа; обнаружение сигнала проводилось по параметру $\sum\nolimits_{{{f}_{1}}}^{{{f}_{3}}} {\left| {A(f)} \right|} + \sum\nolimits_{{{f}_{4}}}^{{{f}_{2}}} {\left| {A(f)} \right|} $, где ${{f}_{1}}$ = 10 Гц, ${{f}_{2}}$= 1 кГц, ${{f}_{3}}$ = 580 Гц, ${{f}_{4}}$ = 740 Гц.

Для каждого обнаруженного сигнала определялось время прихода волны на станцию. Локация сейсмических событий проводилась с помощью метода засечек по 4 станциям. На рис. 4 в качестве примера показаны волновые формы одного из выделенных сигналов, зарегистрированных измерительными точками T₁–T₄. Пунктирная линия соответствует временам первых вступлений. В финальную выборку вошли 94 из 130 событий (рис. 3б) с эпицентральными расстояниями 25–400 м. Для остальных событий сильная зашумленность одного и/или более каналов не позволила провести локацию.

Рис. 4.

Волновая форма события № 43 (табл. А2), записанная всеми сейсмическими станциями. Пунктирная линия соответствует первому вступлению сейсмической волны. Максимальные амплитуды для каждой станции указаны над волновой формой в единицах мм/с², тогда как на рисунке приведены сейсмограммы, нормированные на максимальную амплитуду.

При локации событий было определено время в очаге, что позволило оценить скорость распространения продольных волн на участке измерений в предположении однородности среды ${{C}_{\alpha }} = 5570 \pm 630$ м/с. С учетом коэффициента Пуассона $\nu = 0.27$ была рассчитана скорость распространения поперечных волн, которая составила ${{C}_{\beta }} = 3130 \pm 630$ м/c. Полученные значения скоростей ${{C}_{\alpha }}$ и ${{C}_{\beta }}$ находятся в хорошем соответствии со значениями, полученными при оценке свойств отобранных образцов (табл. 1).

Определение парамеров источника

Оценка очаговых параметров источника проводилась по результатам анализа параметров низкочастотной области спектра сейсмических колебаний. Для оценки величины скалярного сейсмического момента использовалось соотношение [Кейлис-Борок, 1957] в модификации [Gibowicz, Kijko, 1994]

(2)

${{M}_{0}} = \frac{{4\pi \rho {{C}^{3}}R{{\Omega }_{0}}}}{{{{F}_{c}}{{R}_{c}}{{S}_{c}}}},$

где: ρ – плотность массива (при расчетах $\rho = 3400\,\,{{{\text{кг}}} \mathord{\left/ {\vphantom {{{\text{кг}}} {{{{\text{м}}}^{{\text{3}}}}}}} \right. \kern-0em} {{{{\text{м}}}^{{\text{3}}}}}}$); С – скорость распространения волн; R – расстояние от источника до точки наблюдения; Ω₀ – амплитуда “полки” в низкочастотной части спектра смещения, которая оценивается при помощи аппроксимации зарегистрированного спектра одним из модельных; ${{F}_{c}}$ – коэффициент, учитывающий направленность излучения; ${{R}_{c}}$ – коэффициент, определяемый углом выхода луча на свободную поверхность; ${{S}_{c}}$ – станционная поправка. Для большинства событий на малых эпицентральных расстояниях объемные волны не разделяются, поэтому расчет очаговых параметров проводился по всей волновой записи сейсмического сигнала. Как отмечено в работе [Gibowicz, Kijko, 1994], если механизм источника не определен, оценку сейсмического момента можно производить по среднеквадратичному значению коэффициента, учитывающего направленность излучения, рассчитанному исходя из сферичной формы источника, используя значения ${{F}_{\alpha }}$= 0.52 и ${{F}_{\beta }}$ = 0.63 [Boore, Boatwright, 1984]. Следуя [Gibowicz, Kijko, 1994], можно пренебречь станционной поправкой и коэффициентом, который учитывает угол выхода луча на свободную поверхность (${{R}_{c}}$ = 1), для датчиков, расположенных в подземных выработках.

Для аппроксимации низкочастотной части спектра смещения мы использовали модель Дж. Брюна [Brune, 1970], которая описывает характеристики сигнала, излучаемого при внезапном высвобождении сдвигового напряжения на круговой поверхности в упругом изотропном пространстве. Для этой модели спектр смещения описывается функцией:

(3)

$\Omega (f) = \frac{{{{\Omega }_{0}}}}{{1 + {{{{f}^{2}}} \mathord{\left/ {\vphantom {{{{f}^{2}}} {f_{0}^{2}}}} \right. \kern-0em} {f_{0}^{2}}}}},$

где ${{f}_{0}}$ – угловая частота излучения источника. Построение аппроксимирующего спектра согласно модели производилось при минимуме среднеквадратичных отклонений от исходных данных амплитудного спектра смещения. По значениям сейсмического момента можно оценить моментную магнитуду [Kanamori, 1977; Hanks, Kanamori, 1979]:

(4)

${{M}_{w}} = {2 \mathord{\left/ {\vphantom {2 3}} \right. \kern-0em} 3}(\lg {{M}_{0}} - 9.1).$

Расчет излученной сейсмической энергии проводился путем интегрирования спектра скорости [Gibowicz, Kijko, 1994]:

(5)

${{E}_{S}} = 4\pi \rho C{{\left\langle {{{F}_{c}}} \right\rangle }^{2}}\frac{{\sum\limits_n {{{J}_{{{{c}_{n}}}}}R_{n}^{2}} }}{{\sum {F_{{{{c}_{n}}}}^{2}} }},$

где: $n$ – количество точек измерения; ${{J}_{{{{c}_{n}}}}}$ – поток энергии; ${{R}_{n}}$ – расстояние между источником до n-й измерительной станции. Поток энергии определяется как [Gibowicz, Kijko, 1994]

(6)

${{J}_{n}} = 3 \times 2\int\limits_{{{f}_{1}}}^{{{f}_{2}}} {{{{\left| {V(f)} \right|}}^{2}}} df,$

где ${{f}_{1}}$ = 10 Гц, ${{f}_{2}}$ = 1 кГц, $\left| {V(f)} \right|$ – модуль спектра скорости; поскольку система измерений однокомпонентная, использовался множитель “3”. Если механизм очага неизвестен, то коэффициент, учитывающий направленность излучения $\frac{{\sum {F_{{{{c}_{n}}}}^{2}} }}{{4\pi {{{\left\langle {{{F}_{c}}} \right\rangle }}^{2}}}}$ из (5), становится равным $\frac{n}{{4\pi }}$, и формулу (5) можно переписать в виде:

(7)

${{E}_{S}} = 4\pi \rho C\frac{{\sum\limits_n {{{J}_{{{{c}_{n}}}}}R_{n}^{2}} }}{n}.$

При расчете излученной энергии использовались значения потока $J$ с учетом потерь энергии за счет ограниченного частотного диапазона измерительной аппаратуры. Как показано в работе [Ide, Beroza, 2001], в тех случаях, когда полоса регистрации обрезана сверху частотой ${{f}_{M}}$, отношение зарегистрированной энергии ${{E}_{{reg}}}$ к излученной E_s определяется выражением

(8)

$\begin{gathered} R = \frac{{{{E}_{{reg}}}}}{{{{E}_{s}}}} = \\ = \frac{2}{\pi }\left[ {{\text{arctg}}({{{{f}_{M}}} \mathord{\left/ {\vphantom {{{{f}_{M}}} {{{f}_{0}}}}} \right. \kern-0em} {{{f}_{0}}}}) - {{({{{{f}_{M}}} \mathord{\left/ {\vphantom {{{{f}_{M}}} {{{f}_{0}}}}} \right. \kern-0em} {{{f}_{0}}}})} \mathord{\left/ {\vphantom {{({{{{f}_{M}}} \mathord{\left/ {\vphantom {{{{f}_{M}}} {{{f}_{0}}}}} \right. \kern-0em} {{{f}_{0}}}})} {{{{(1 + {{{{f}_{M}}} \mathord{\left/ {\vphantom {{{{f}_{M}}} {{{f}_{0}}}}} \right. \kern-0em} {{{f}_{0}}}})}}^{2}}}}} \right. \kern-0em} {{{{(1 + {{{{f}_{M}}} \mathord{\left/ {\vphantom {{{{f}_{M}}} {{{f}_{0}}}}} \right. \kern-0em} {{{f}_{0}}}})}}^{2}}}}} \right], \\ \end{gathered} $

где ${{f}_{0}}$ – угловая частота спектра излученного сигнала. Для используемой измерительной сети частота ${{f}_{M}}$ равна 1 кГц для всех датчиков. Данные измерительной точки Т₄, на которой присутствуют стационарные помехи, были исключены из рассмотрения, чтобы не вносить дополнительных инструментальных ошибок.

Для оценки величины сброса напряжения в массиве использовалось соотношение для модели кругового разлома R. Madariaga [Madariaga, 1976]:

(9)

$\Delta \sigma = \frac{7}{{16}}{{M}_{0}}{{\left( {\frac{{2\pi {{f}_{0}}}}{{K{{C}_{\beta }}}}} \right)}^{3}},$

где K принимает значения ${{K}_{\alpha }} = 2.01$ и ${{K}_{\beta }} = 1.32$ для P- и S-волн соответственно. Полученные оценки дают хорошую сходимость с результатами наблюдений в шахтах [Gibowicz et al., 1990; Gibowicz, Kijko, 1994]. Если механизм очага не определен, то могут быть использованы средние значения коэффициентов ${{K}_{\alpha }}$ и ${{K}_{\beta }}$ [Gibowicz, Kijko, 1994].

Кажущееся напряжение оценивалось по обычному соотношению:

(10)

$\sigma {}_{a} = \mu \frac{{{{E}_{s}}}}{{{{M}_{0}}}},$

где модуль сдвига среды $\mu \approx 33$ ГПа.

Исследования по определению добротности исследуемой области ранее не проводились. Опираясь на данные [Oye et al., 2005], где для рудного месторождения наилучшая сходимость с модельным решением была получена для значения Q = 350, нами были рассмотрены две модели среды с добротностью $Q$ равной 200 и 400 путем введения множителя $\exp \left( { - \frac{{\pi Rf}}{{QC}}} \right)$. На рис. 5 представлены спектры смещения для среды с учетом добротности $Q$ = 200 и $Q$ = 400 и среды без затухания. Эффекты введения затухания видны только на высоких частотах, поэтому далее все расчеты проведены для модели среды без учета затухания.

Рис. 5.

Амплитудный спектр смещения, рассчитанный для трассы без учета затухания (1) и с постоянными Q = 200 (2) и Q = 400 (3) для события № 13 (табл. A2 ). Приближение по модели ω² показано пунктирной линией. Эффект Q невелик и наблюдается только на высоких частотах.

РЕЗУЛЬТАТЫ

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

Рис. 6.

(а) – Кумулятивное количество событий, обнаруженных до и после взрыва; нулевая точка соответствует времени взрыва; сплошная линия – логарифмическое приближение (11); (б) – взаимосвязь между количеством событий, магнитудой ${{M}_{w}}$ и плотностью вероятности распределения количества событий по магнитуде (гистограмма, серый цвет). Пунктирная линия соответствует параметру b-value, равному 1.

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

На рис. 7 приведена зависимость магнитуд зарегистрированных событий от расстояния. Линия определяет приблизительный порог чувствительности сейсмической группы. Это означает, что все зарегистрированные события были инициированы, в основном, взрывом в камере 1.2 и, возможно, в камере 1.1. Динамические деформации от взрывов в камерах №№ 2 и 3 были значительно меньше по амплитуде (рис. 8).

Рис. 7.

Магнитуды зарегистрированных событий в зависимости от расстояния до точки T₃. Прямая линия – оценка порога чувствительности сейсмической сети.

Рис. 8.

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

Для модели ω² проведена оценка моментных магнитуд ${{M}_{w}}$ выделенных событий в диапазоне от –2.7 до –1.4 с угловой частотой от 120 до 760 Гц. Полученные значения магнитуд и угловых частот соответствуют, согласно модели в работе [Madariaga, 1976], радиусам очага в первые метры. Величина сброса напряжения варьируется от 0.002 до 0.25 MПa (рис. 9). Величина рассчитанной сейсмической энергии варьируется от 10^–3 до 1 Дж при изменении сейсмического момента примерно на два порядка от 10⁵ до 10⁷ Н ⋅ м, что соответствует вариации кажущегося напряжения от 0.0003 до 0.017 MПa (рис. 10a). Несмотря на значительный разброс, можно отметить близкий к линейному тренд роста средней величины приведенной энергии ${{{{E}_{s}}} \mathord{\left/ {\vphantom {{{{E}_{s}}} {{{M}_{0}}}}} \right. \kern-0em} {{{M}_{0}}}}$ примерно на порядок в диапазоне изменения момента от 1.5 × 10⁵ до 1.5 × 10⁶ Н ⋅ м (рис. 10б). При этом в диапазоне ${{M}_{0}}$ (10⁶–10⁷) такой тенденции не отмечено.

Рис. 9.

Соотношение между скалярным сейсмическим моментом ${{M}_{0}}$ и угловой частотой источника ${{f}_{0}}$ для зарегистрированных слабых сейсмических событий. Значения постоянного сброса напряжения рассчитаны для модели из источника [Madariaga, 1976].

Рис. 10.

(а) – Излученная сейсмическая энергия ${{E}_{s}}$ в зависимости от скалярного сейсмического момента ${{M}_{0}}$ (моментной магнитуды ${{M}_{w}}$); (б) – приведенная сейсмическая энергия ${{{{E}_{s}}} \mathord{\left/ {\vphantom {{{{E}_{s}}} {{{M}_{0}}}}} \right. \kern-0em} {{{M}_{0}}}}$ в зависимости от скалярного сейсмического момента ${{M}_{0}}$. Уровни постоянного кажущегося напряжения ${{\sigma }_{a}}$ показаны пунктирными линиями с учетом модуля сдвига 33 ГПа.

Отношение кажущегося напряжения ${{\sigma }_{a}}$ к величине сброса напряжения $\Delta \sigma $ можно использовать для оценки скорости распространения разрыва ${{V}_{r}}$ через параметр излучательной эффективности [Husseini, 1977].

Рис. 11.

Соотношение между моментной магнитудой и ${{{{\sigma }_{a}}} \mathord{\left/ {\vphantom {{{{\sigma }_{a}}} {\Delta \sigma }}} \right. \kern-0em} {\Delta \sigma }}$ со средним значением, равным 0.11.

Мониторинг динамических перемещений по тектоническим трещинам, проведенный на том же участке, не выявил значимых остаточных смещений, инициируемых взрывами. Методика и результаты этих измерений приведены в работах [Адушкин и др., 2017а; Кочарян и др., 2019]. Величина относительного смещения берегов трещины после прохождения сейсмических колебаний от взрыва невелика и составляет от единиц до 10–20 мкм.

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

Расположение очагов зарегистрированных событий показано на рис. 12. Можно видеть, что практически все точки находятся внутри зоны, в которой уровень динамических деформаций от взрывов достигает величины 10^–6. Судя по тому, что до воздействия активность практически отсутствует, все эти события инициированы взрывом. Большая часть очагов приурочена к зоне разгрузки вокруг разрабатываемой камеры (схематически показана прямоугольником на рис. 12).

Рис. 12.

Расположение очагов зарегистрированных событий. Серые точки – сейсмические события, зарегистрированные после взрывов. Точка (0,0) соответствует измерительной точке T₃. Прямоугольник – взрывная камера; толстая линия внутри прямоугольника – зона взрывных скважин. Большие круги (сплошная линия) соответствуют уровню динамических деформаций 10^–6 от взрыва, маленький круг (пунктирная линия) – динамической деформации 3 × 10^–6.

Результаты мониторинга участка разломной зоны Коробковского месторождения говорят о том, что в массиве наблюдаются только события, связанные с поствзрывной активностью. C учетом предыдущих работ по оценке напряженно-деформированного состояния массива [Кочарян и др., 2018] и по мониторингу сейсмической активности [Besedina et al., 2018], это свидетельствует о спокойной геодинамической обстановке в пределах рудного поля.

Как отмечалось выше, совместное рассмотрение большого количества рядов данных, относящихся к разным участкам, создает видимость постоянства величины средней приведенной энергии [Walter et al., 2006; Kwiatek et al., 2011; Кочарян и др., 2016]. Однако наклон регрессионных зависимостей ${{E}_{s}}({{M}_{0}})$ для отдельных рядов наблюдений индуцированной сейсмичности иногда сильно отличается от единицы. Например, на месторождении в Канаде (события с магнитудами от –2.2 до –0.4) [Urbancic, Young, 1993] и при мониторинге наведенной сейсмичности в Китае (магнитуды 0.1–4.2) [Hua et al., 2013] получены одинаково сильные зависимости для масштабных соотношений ${{{{E}_{s}}} \mathord{\left/ {\vphantom {{{{E}_{s}}} {{{M}_{0}}}}} \right. \kern-0em} {{{M}_{0}}}}{{ \sim }}M_{0}^{{0.63}}$. Обработка результатов наблюдений за микросейсмичностью (155 событий магнитудой от –3.6 до –1.9), проводившихся при проходке участка шахты в граните [Gibowicz et al., 1991], дает соотношение ${{{{E}_{s}}} \mathord{\left/ {\vphantom {{{{E}_{s}}} {{{M}_{0}}}}} \right. \kern-0em} {{{M}_{0}}}}\sim M_{0}^{{0.33}}$. В Германии, при анализе индуцированных сейсмических событий с магнитудами от –2 до 0 получены зависимости приведенной энергии от масштаба ${{{{E}_{s}}} \mathord{\left/ {\vphantom {{{{E}_{s}}} {{{M}_{0}}}}} \right. \kern-0em} {{{M}_{0}}}}{{ \sim }}M_{0}^{{0.90}}$ [Jost et al., 1998]. При этом разброс данных, как правило, весьма велик, так, что регрессии отражают лишь общий тренд для каждого ряда. На рис. 13 сопоставлены результаты, полученные нами на Коробковском месторождении, с некоторыми мировыми данными по регистрации слабой сейсмичности.

Рис. 13.

Зависимость приведенной сейсмической энергии ${{{{E}_{s}}} \mathord{\left/ {\vphantom {{{{E}_{s}}} {{{M}_{0}}}}} \right. \kern-0em} {{{M}_{0}}}}$ от сейсмического момента ${{M}_{0}}$: 1 – [Yamada et al., 2007]; 2 – [Urbancic, Young, 1993]; 3 – [Oye et al., 2005]; 4 – [Gibowicz et al., 1991]; 5 – [Kwiatek et al., 2011]; 6 – [Маловичко, Маловичко, 2010]; 7 – [Jost et al., 1998]; 8 – Коробковскоe месторождение.

Некоторые возможные причины отклонения от закона подобия подробно описаны, например, в работе [Беседина и др., 2013], где, наряду с влиянием ограниченной частотной характеристики регистрирующей аппаратуры [Ide, Beroza, 2001], рассматривается эффект, связанный с поглощением и рассеянием высоких частот при распространении колебаний. Эти причины отклонения от закона подобия не связаны с физикой очага динамического события и не отражают эффективность последнего как источника сейсмических колебаний.

Отклонения от закона подобия могут быть связаны и с зависимостью от масштаба модуля деформации массива горных пород и жесткости нарушений сплошности. Это следует из рассмотрения баланса энергии при динамической подвижке по трещине. Рассматривая баланс энергии (рис. 14а), следуя, например, [Kanamori, Brodsky, 2002], предположим для простоты линейный закон фрикционного разупрочнения (рис. 14), т.е. постоянную скорость снижения сопротивления сдвигу в процессе перемещения $\left| {\frac{{d\sigma }}{{dD}}} \right| = {{k}_{f}} = {\text{const}}$, где ${{k}_{f}}$ – сдвиговая жесткость разгрузки (скорость снижения сопротивления сдвигу) разломной зоны, $D$ – относительное перемещение берегов, $\sigma $ – сдвиговое напряжение.

Рис. 14.

(а) – Баланс энергии при динамической подвижке по трещине: ${{E}_{G}}$ – энергия, затраченная на трещинообразование, ${{E}_{s}}$ – энергия, излученная сейсмическими волнами, $\Delta {{E}_{f}}$ – энергия, диссипированная в процессе скольжения; (б) – трансформация жесткости разломной зоны в окрестности предела прочности.

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

Соответственно, излучательная эффективность подвижки определяется выражением:

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

Одним из условий формирования динамического срыва по существующему разлому является соотношение [Scholz, 1998]:

Оценить величину скорости снижения сопротивления сдвигу без проведения специальных измерений можно лишь опосредованно, используя результаты лабораторных [Pyrak-Nolte, 1996] и in situ [Kocharyan et al., 1997; Кочарян, 2014] экспериментальных определений “динамической” жесткости трещин и разломов по параметрам сейсмических волн, прошедших через нарушение сплошности. В работе [Кочарян, 2014] приводится эмпирическое выражение для оценки сдвиговой жесткости природных трещин и небольших разломов длиной $L$ менее 500–1000 м in situ:

Соотношение (19) демонстрирует тот же закон масштабирования, что и многочисленные данные экспериментального определения статической сдвиговой жесткости, согласно которым жесткость обратно пропорциональна размерам образца [Barton, 2007]. При этом, судя по данным лабораторных опытов, динамическая жесткость по крайней мере в несколько раз выше соответствующей величины сдвиговой жесткости при статическом нагружении [Pyrak-Nolte, 1987; 1996; Barton, 2007]. Так что для оценок можно полагать, что статическая жесткость для нарушения сплошности исследуемого диапазона описывается соотношением:

Очевидно, что если определяемая в эксперименте величина ${{k}_{s}} = \frac{{d\sigma }}{{dD}}$ относится к восходящей ветви зависимости $\sigma (D)$, то скорость снижения сопротивления сдвигу ${{k}_{f}}$ при скольжении представляет собой абсолютную величину производной на запредельном участке зависимости $\sigma (D)$ после достижения предела прочности (рис. 14б). При этом по мере приближения к пределу прочности величина жесткости заметно снижается [Johnson, Xia, 2005; Кочарян и др., 2018]. Опираясь на многочисленные данные лабораторных экспериментов [Marti et al., 2017], можно заключить, что в окрестности экстремума форма зависимости $\sigma (D)$ часто близка к симметричной, так что в грубом приближении ${{k}_{f}} \approx {{k}_{{s0}}} < {{k}_{s}}.$

В соответствии с (20) можно принять:

Соотношение (21) дает разумную оценку величины сброса напряжения. В самом деле:

Соответственно, из (21) и (22) имеем $\Delta \sigma \,\,{{\sim }}\,\,{\text{1}}{\kern 1pt} - {\kern 1pt} {\text{5 МПа}}$.

Преобразуя (17) и (21), имеем:

В тех случаях, когда мелкие события связаны с образованием новых трещин, можно считать справедливым приближение Кейлиса-Борока–Эшелби [Keilis-Borok, 1959], согласно которому скачок напряжений при смещении по круговой трещине пропорционален амплитуде смещения Δu и обратно пропорционален радиусу r трещины, что приводит к независимости величины приведенной сейсмической энергии от масштаба при постоянстве упругих характеристик среды.

Наблюдаемая на шахте им. Губкина комбината КМАруда и на многих других горных предприятиях тенденция к росту среднего значения e (см. рис. 13) с увеличением размера события, возможно, связана с эффектом снижения значения модуля сдвига породы при увеличении характерного размера задачи. Хорошо известно, что на относительно небольших глубинах результаты измерения скорости распространения колебаний in situ неизменно зависят от базы измерений. Так, например, результаты измерений скорости распространения поперечных волн в массиве кварцита на глубине около 1000 м демонстрируют устойчивое снижение этого параметра от значения, характерного для измерений на образцах ${{C}_{\beta }} \approx 3000\,\,{{\text{м}} \mathord{\left/ {\vphantom {{\text{м}} {\text{с}}}} \right. \kern-0em} {\text{с}}}$, до величины ${{C}_{\beta }} \approx 2600\,\,{{\text{м}} \mathord{\left/ {\vphantom {{\text{м}} {\text{с}}}} \right. \kern-0em} {\text{с}}}$, полученной при измерениях на базах 100–500 м [Кочарян, 2013]. В механике горных пород хорошо известен эффект снижения прочностных характеристик и деформационных модулей массива горных пород по сравнению с характеристиками образца, причем эта разница возрастает в зонах влияния выработки и в разломных зонах [Barton, 2007]. Так, по данным [Barton, 1999], на глубине 300 м эффективный модуль деформации ненарушенной породы достигает величины $E$ ~ 80 ГПа, модуль горной породы в окрестности мелких разломов $E$ ~ 20 ГПа, а в зоне влияния крупного разлома $E$ ~ 10 ГПа.

Гипоцентры большинства индуцированных событий, зарегистрированных на шахте Губкина, расположены либо в зоне влияния разрабатываемой камеры, либо в зоне влияния крупного разлома. Чем больше масштаб события (и, соответственно, длина разрыва), тем выше вероятность того, что очаг, либо его значительная часть, расположен на участке с относительно низким эффективным модулем, и, следовательно, эффективность сейсмического источника будет достаточно высокой. Заметим, что при проведении наблюдений в глубоких шахтах ($H$ ~ 3.5 км) эффекта увеличения кажущихся напряжений с ростом масштаба обнаружено не было [Kwiatek et al., 2011]. Возможно, из-за большого литостатического давления величина $\mu $ in situ достигает значений, близких к наблюдаемым в образцах.

Обращают внимание полученные крайне низкие значения величины приведенной сейсмической энергии ${{{{E}_{S}}} \mathord{\left/ {\vphantom {{{{E}_{S}}} {{{M}_{0}}}}} \right. \kern-0em} {{{M}_{0}}}}\,\,{{\sim }}\,\,{\text{1}}{{{\text{0}}}^{{ - {\text{8}}}}}{\kern 1pt} - {\kern 1pt} 5 \times {{10}^{{ - 7}}}$. Это позволяет отнести деформационные события в массиве Коробковского месторождения к категории медленных землетрясений [Peng, Gomberg, 2010], что подтверждается и приведенной выше оценкой скорости распространения разрыва ${{V}_{R}} \cong (0.1{\kern 1pt} - {\kern 1pt} 0.7){{C}_{\beta }}$.

ВЫВОДЫ

Мы рассмотрели результаты мониторинга индуцированной взрывом сейсмичности на Коробковском месторождении Курской магнитной аномалии. Зарегистрированные события с моментными магнитудами от –2.7 до –1.4 оказались сгруппированы вокруг разрабатываемой камеры и в окрестности крупной разломной зоны. Очаги практически всех зарегистрированных событий находятся внутри области с уровнем динамических деформаций от взрывов более 10^–6.

Для зарегистрированных событий были определены очаговые параметры и установлены эмпирические корреляционные зависимости между скалярным сейсмическим моментом и излученной энергией. Отмечается близкий к линейному тренд роста средней величины приведенной энергии ${{{{E}_{s}}} \mathord{\left/ {\vphantom {{{{E}_{s}}} {{{M}_{0}}}}} \right. \kern-0em} {{{M}_{0}}}}$: примерно на порядок в диапазоне изменения момента от 1.5 × 10⁵ до 1.5 × 10⁶ Н ⋅ м. При этом в диапазоне ${{M}_{0}}$ (10⁶–10⁷ Н ⋅ м) такой тенденции не отмечено. Можно предположить, что этот эффект связан с тенденцией снижения модуля сдвига породы при увеличении размера сейсмического события.

Рассчитанные значения величин приведенной сейсмической энергии ${{{{E}_{\beta }}} \mathord{\left/ {\vphantom {{{{E}_{\beta }}} {{{M}_{0}}}}} \right. \kern-0em} {{{M}_{0}}}}{{\sim 1}}{{{\text{0}}}^{{ - {\text{8}}}}}{\kern 1pt} - {\kern 1pt} 5 \times {{10}^{{ - 7}}}$ и скорости распространения разрыва позволяют предполагать, что промышленный взрыв инициировал рой сейсмических событий с низкой излучательной эффективностью – медленных микроземлетрясений.

ПРИЛОЖЕНИЕ A.