Texto integral
Введение
Геоакустическая эмиссия это процесс генерации упругих волн горными породами в результате динамической перестройки их структуры. На ее динамику влияют различные механические процессы, протекающие в горных породах, в том числе и механические процессы в очаге готовящегося землетрясения. В ряде исследований установлено [1–5], что предсейсмический отклик наблюдается в сигналах геоакустической эмиссии, зарегистрированной на расстоянии первых сотен километров от источника землетрясения. С целью обоснования связи между вариациями геоакустической эмиссии и характером деформационного процесса ранее было проведено моделирование зон геоакустической эмиссии в упругом однородном приближении среды [6, 7]. Зона геоакустической эмиссии это область поверхности земной коры с деформациями порядка от 10−8 до 10−5, вызванными процессами в очаге готовящегося землетрясения. Такие значения деформаций выше приливных, но ниже порога прочности горных пород.
Результаты ранее проведенного моделирования показывают, что уровень расчетных деформаций в пунктах наблюдений превышает приливные, но на порядок ниже зарегистрированных [6]. Это может быть результатом приближения земной коры в виде однородного пространства. В действительности земная кора состоит из слоев горных пород, часть из которых находятся в закритическом состоянии и проявляют пластические и квазипластические свойства [8, 9]. Поэтому цель настоящей работы заключается в моделировании влияния неоднородных включений на пространственное распределение зон геоакустической эмиссии.
Модель очага землетрясения. Согласно концепции, впервые предложенной Б. В. Костровым, тектоническое землетрясение представляет собой разрыв сплошности материала Земли, который возникает под действием упругих сдвиговых напряжений, накопленных в процессе тектонической деформации. Такой разрыв является разрывом скольжения.
В момент землетрясения происходит полное или частичное снятие накопленных напряжений в его очаге. Деформации, возникающие при подготовке землетрясения, обусловлены приращением потенциальной энергии упругих деформаций , вызванным процессом подготовки землетрясения. Эта энергия больше, чем высвободившаяся сейсмическая энергия . Величина , равная отношению этих энергий, определяет эффективность снятия потенциальной энергии упругих деформаций
(1)
Очаг землетрясения можно описать через некоторую систему сил, распределенную по поверхности разрыва. Для описания произвольно ориентированного разрыва скольжения в изотропной среде используют систему, состоящую из девяти пар двойных сил [10, 11].
Математическая постановка задачи
Уравнения и граничные условия. Рассмотрим земную кору в виде упругого изотропного полупространства. Поведение такой среды можно описать при помощи системы, состоящей из уравнений равновесия, определяющих соотношений обобщенного закона Гука, а также выражений, определяющих тензор малых деформаций :
, (2)
где тензор напряжений, вектор массовых сил, коэффициенты Ламе, вектор смещений. Уравнения записаны в тензорной форме, индексами после запятой обозначено дифференцирование по соответствующим пространственным координатам.
Пусть полупространство занимает область . Тогда поверхность Земли задается уравнением . Эта поверхность свободна от напряжений в направлении оси , следовательно на заданы граничные условия вида:
(3)
Напряжения, создаваемые очагом готовящегося землетрясения, стремятся к нулю на бесконечности:
(4)
где расстояние от точки до начала координат.
Компоненты вектора , соответствующие системе из комбинации двойных пар сил, выражаются следующим образом:
(5)
где интенсивность соответствующей пары сил; дельта-функция; точка приложения системы сил.
Аналитическое решение. Для задачи (2) с граничными условиями (3) и (4) известны решения в терминах вектора смещений через функции Грина. Для единичной силы, приложенной к точке упругого полупространства и направленной вдоль оси , функция Грина имеет вид:
(6)
где коэффициент Пуассона, а и :
(7)
Функция Грина для единичной силы, направленной вдоль оси , выражается в виде:
(8)
Ввиду симметричности задачи, функция Грина , соответствующая действию единичной силы вдоль оси , может быть получена из функции заменой осей и . Функции Грина, отвечающие действию двойных сил, могут быть получены дифференцированием функций по пространственным координатам, то есть в виде [12]. При получим решение для пары двойных сил, направленных вдоль соответствующей оси, при для пары двойных сил, направленных вдоль оси с моментом относительно оси с номером, отличным от .
В общем случае, решения для смещений в упругом полупространстве можно получить при помощи формулы Вольтерра [13]:
(9)
где смещение на поверхности разрыва , единичный вектор нормали к поверхности .
Учитывая, что напряжения могут быть выражены через деформации в соответствии с законом Гука, формула Вольтерра для случая однородной и изотропной среды может быть записана в виде:
(10)
где тензор плотности сейсмического момента [10], который отражает механику очага землетрясения.
Таким образом, в случае точечного источника, решение поставленной задачи может быть найдено в следующем виде:
(11)
где площадь поверхности силового воздействия .
Поскольку связь между компонентами тензора деформации и потенциальной энергией упругих деформаций квадратична, то повышающий коэффициент, позволяющий рассчитать напряженно-деформированное состояние земной коры при подготовке землетрясения, будем полагать равным . Удобный, с точки зрения вычислений, вариант оценки этого коэффициента был дан И. П. Добровольским [14]:
(12)
где моментная магнитуда землетрясения.
Неоднородное включение в земной коре. Область земной коры, обладающую неупругими свойствами можно описать через систему сил, распределенную по ее границе. Подобное силовое воздействие характерно для слоев горных пород, находящихся в закритическом состоянии и проявляющих пластические свойства.
Будем полагать, что неоднородное включение имеет сферическую форму радиуса с границей (рис. 1), которая задана функцией :
(13)
где центр сферического включения.
Рис. 1. Сферическое включение Ω, ограниченное поверхностью Γ с центром в точке (c1,c2,c3), на которой заданы простые силы интенсивности P (🡪) по направлению внешней нормали.
[Figure 1. Spherical inclusion Ω, limited by the surface Γ with the center at the point (c1,c2,c3), on which simple forces of intensity P (🡪) are specified in the direction of the outer normal. ]
Пусть на поверхности в каждой точке заданы простые силы интенсивности , действующие по направлению внешней нормали . Решение для поля перемещений может быть найдено в виде поверхностного интеграла от функций Грина:
(14)
Преобразуя интеграл (14) при помощи формулы Остроградского-Гаусса [15], получим
(15)
Для упрощения интегрирования перейдем в сферическую систему координат :
(16)
Производя замену переменных интегрирования (16) в объемном интеграле (15), получим
(17)
где якобиан.
Результаты моделирования
Вычислительный эксперимент №1. Вычисления во всех экспериментах производятся на поверхности . Параметры гипотетического землетрясения следующие: координаты м, м, м, тензор плотности сейсмического момента
(18)
Параметры неоднородного включения: координаты центра м, м, м, радиус м, величина силы Н. Результаты моделирования в виде компонентов вектора смещений на поверхности представлены на рис. 2 4.
Рис. 2. Результат вычислительного эксперимента №1. Компонент u1 вектора смещений на поверхности x3 = 0 в случае однородной среды (a) и при наличии сферической неоднородности (b); 🟠 проекция сферической неоднородности на поверхность x3 = 0.
[Figure 2. Result of the computational experiment №1. Component u1 of the shift vector on the surface x3 = 0 in case of homogeneous environment (a) and with a spherical inhomogeneity (b); 🟠 projection of the spherical inhomogeneity on the surface x3 = 0. ]
Рис. 3. Результат вычислительного эксперимента №1. Компонент u2 вектора смещений на поверхности x3 = 0 в случае однородной среды (a) и при наличии сферической неоднородности (b); 🟠 проекция сферической неоднородности на поверхность x3 = 0.
[Figure 3. Result of the computational experiment №1. Component u2 of the shift vector on the surface x3 = 0 in case of homogeneous environment (a) and with a spherical inhomogeneity (b); 🟠 projection of the spherical inhomogeneity on the surface x3 = 0.]
Рис. 4. Результат вычислительного эксперимента №1. Компонент u3 вектора смещений на поверхности x3 = 0 в случае однородной среды (a) и при наличии сферической неоднородности (b); 🟠 проекция сферической неоднородности на поверхность x3 = 0.
[Figure 4. Result of the computational experiment №1. Component u3 of the shift vector on the surface x3 = 0 in case of homogeneous environment (a) and with a spherical inhomogeneity (b); 🟠 projection of the spherical inhomogeneity on the surface x3 = 0.]
Вычислительный эксперимент №2. Параметры гипотетического землетрясения следующие: координаты м, м, м, тензор плотности сейсмического момента
(19)
Параметры неоднородного включения: координаты центра м, м, м, радиус м, величина силы Н. Результаты моделирования в виде компонентов вектора смещений на поверхности представлены на рис. 5 7.
Рис. 5. Результат вычислительного эксперимента №2. Компонент u1 вектора смещений на поверхности x3 = 0 в случае однородной среды (a) и при наличии сферической неоднородности (b); 🟠 проекция сферической неоднородности на поверхность x3 = 0.
[Figure 5. Result of the computational experiment 2. Component u1 of the shift vector on the surface x3 = 0 in case of homogeneous environment (a) and with a spherical inhomogeneity (b); 🟠 projection of the spherical inhomogeneity on the surface x3 = 0.]
Рис. 6. Результат вычислительного эксперимента №2. Компонент u2 вектора смещений на поверхности x3 = 0 в случае однородной среды (a) и при наличии сферической неоднородности (b); 🟠 проекция сферической неоднородности на поверхность x3 = 0.
[Figure 6. Result of the computational experiment №2. Component u2 of the shift vector on the surface x3 = 0 in case of homogeneous environment (a) and with a spherical inhomogeneity (b); 🟠 projection of the spherical inhomogeneity on the surface x3 = 0.]
Рис. 7. Результат вычислительного эксперимента №2. Компонент u3 вектора смещений на поверхности x3 = 0 в случае однородной среды (a) и при наличии сферической неоднородности (b); 🟠 проекция сферической неоднородности на поверхность x3 = 0.
[Figure 7. Result of the computational experiment №2. Component u3 of the shift vector on the surface x3 = 0 in case of homogeneous environment (a) and with a spherical inhomogeneity (b); 🟠 projection of the spherical inhomogeneity on the surface x3 = 0.]
Вычислительный эксперимент №3. Все параметры эксперимента идентичны предыдущему. Добавлено второе неоднородное включение с центром в точке м, м, м и величиной силы Н. Результаты моделирования в виде компонентов вектора смещений на поверхности представлены на рис. 8 10.
Рис. 8. Результат вычислительного эксперимента №3. Компонент u1 вектора смещений на поверхности x3 = 0 в случае однородной среды (a) и при наличии сферической неоднородности (b); 🟠 проекция сферической неоднородности на поверхность x3 = 0.
[Figure 8. Result of the computational experiment №3. Component u1 of the shift vector on the surface x3 = 0 in case of homogeneous environment (a) and with a spherical inhomogeneity (b); 🟠 projection of the spherical inhomogeneity on the surface x3 = 0.]
Рис. 9. Результат вычислительного эксперимента №3. Компонент u2 вектора смещений на поверхности x3 = 0 в случае однородной среды (a) и при наличии сферической неоднородности (b); 🟠 проекция сферической неоднородности на поверхность x3 = 0.
[Figure 9. Result of the computational experiment №3. Component u2 of the shift vector on the surface x3 = 0 in case of homogeneous environment (a) and with a spherical inhomogeneity (b); 🟠 projection of the spherical inhomogeneity on the surface x3 = 0.]
Рис. 10. Результат вычислительного эксперимента №3. Компонент u3 вектора смещений на поверхности x3 = 0 в случае однородной среды (a) и при наличии сферической неоднородности (b); 🟠 проекция сферической неоднородности на поверхность x3 = 0.
[Figure 10. Result of the computational experiment №3. Component u3 of the shift vector on the surface x3 = 0 in case of homogeneous environment (a) and with a spherical inhomogeneity (b); 🟠 projection of the spherical inhomogeneity on the surface x3 = 0.]
Обсуждение результатов. Как показывают результаты моделирования, наличие неоднородных включения в среде оказывает влияние на смещения поверхности земной коры. Наблюдается, например, «экранирование» части области поверхности земной коры (рис. 8b), которое выражается в уменьшении абсолютного значения компонентов вектора смещений. Это может быть одной из причиной несогласованности результатов наблюдений геоакустической эмиссии и результатов моделирования в случае однородной среды.
Заключение
Произведено моделирование влияния неоднородного строения земной коры на зоны геоакустической эмиссии. Неоднородности описываются системой сил, распределенной по поверхности сферического включения. Решения получены в виде свертки функций Грина для упругого полупространства. Показано, что сферические включения оказывают влияние на поле вектора смещений поверхности земной коры. Может наблюдаться как увеличение, так и уменьшение абсолютного значения компонентов вектора смещений в зависимости от количества неоднородных включений и их расположения относительно очага готовящегося землетрясения.
В настоящей работе величина сосредоточенных сил, распределенных по поверхности неоднородного включения, полагалась постоянной. В реальности же будет наблюдаться зависимость от пространственных координат, т. е. . Построение такой модели может являться дальнейшим направлением развития работы.