Geomechanical Markers of the Stress and Strain State and Interaction of Structures in Inhomogeneous Geoenvironments

Full Text

1. Введение

Геодинамические явления в неоднородной 3D-геосреде рассматриваем как деформационные процессы, которые характеризуются полями перемещений, поворотов и деформаций. Количественные и пространственные характеристики распределения этих полей обеспечиваются соответствующими полями напряжений. Напряженность и соответствующее гравитационное поле, обусловленное плотностной неоднородностью Земли, является одним из основных внутренних источников возникновения тектонических напряжений [1]. Гравитационное поле обеспечивает формирование соответствующего усложненного деформирования геосреды, состоящего из перемещений, сдвигов и поворотов текущих точек по всем координатным направлениям.

Механико-математическое моделирование напряженно-деформированного состояния неоднородных геосред основано на теории подобия реологических процессов в неоднородных многослойных средах, подверженных конечным деформациям [2–8]. В трехмерной постановке теории конечных деформаций и линеаризованной теории устойчивости предложена единая 3D-модель для изучения возможных дислокаций в неоднородных геосредах. Напряженно-деформированное состояние (НДС) разделяется на основное и возмущенное [4, 8]. Возможность существования основного (невозмущенного) состояния с пренебрежимо малыми напряжениями и деформациями, которое существует до начала действия внутренних сил, действующих в земной коре, и обусловлено пренебрежимо малым гравитационным полем, является гипотетическим предположением. Кинематика возмущенного состояния определяется компонентами тензора конечных деформаций Коши–Грина, поступательными перемещениями и перемещениями линейного поворота, напряженное состояние – компонентами несимметричного тензора напряжений Пиола–Кирхгофа. Аналитический алгоритм основан на получении решений относительно параметров, имеющих физический и геометрический смысл. В данной модели – это радиальное перемещение, перемещение поворота и результирующая по главным направлениям деформация. Предложенная методика позволяет на этапе численно-графического анализа встроить неоднородное распределение физико-механических свойств, конкретизировать граничные условия, уравнение состояния, количественно и качественно оценить порядок нелинейности параметров, вклад каждой компоненты в распределение результирующих полей напряженно-деформированного состояния [9–11].

В статье [11] поставлена и решена задача исследования устойчивого равновесия для основного и возмущенного состояний полого шара с граничными условиями на внутреннем основании и на внешней (дневной) поверхности. Внешняя поверхность (ρ = ρ₂) полого шара свободна от нагружения. На сферической поверхности (ρ = ρ₁) в основании полого шара задается граничное условие в напряжениях с учетом следящего давления $\tilde{p}$. Под внутренним следящим давлением понимается нагрузка, которая направлена по нормали к поверхности и не изменяет направление и величину в процессе деформирования [8]. Его значение соответствует давлению на заданной глубине и тем самым моделирует влияние отсутствующей внутренней части шара [12]. Решение задачи выполнено в системе координат O_ρφλ (ρ – радиус, ρ₁ ÷ ρ₂ км; φ – широта, –90° ÷ 90°; λ – долгота, –180° ÷ 180°) в физических составляющих компонент тензора деформаций Грина, несимметричного тензора напряжений Пиола–Кирхгофа, компонент векторов перемещений и поворотов. Для основного и возмущенного состояний определяющие уравнения и их решения получены в аналитической форме для закона состояния в форме нелинейного вязкоупругого потенциала

$$\begin{gathered} {\Phi }^{ij}\left({\epsilon }^{ij},{\xi }^{ij}\right)={\vartheta }_1{\epsilon }^{ij}+{\vartheta }_2{\epsilon }^{ik}{\epsilon }_k^j+{\vartheta }_3{\xi }^{ij}+{\vartheta }_4{\xi }^{ik}{\xi }_k^j+{\vartheta }_5\left({\epsilon }^{ik}{\xi }_k^j+{\xi }^{ik}{\epsilon }_k^j\right)+ \\ +{\vartheta }_6{\epsilon }^{ik}{\xi }_k^n{\epsilon }_n^j+{\vartheta }_7{\xi }^{ik}{\epsilon }_k^n{\xi }_n^j, \end{gathered}$$

где ε^ij – контравариантные компоненты тензора деформаций Коши–Грина, ξ^ij – контравариантные компоненты тензора скоростей деформаций, ϑ_i – функции физических констант (упругости, жесткости, вязкости и других свойств геосред).

Здесь используются результаты построенного решения задачи и комплексы программ для расчета всех компонент и результирующих параметров НДС для закона состояния в форме потенциала Мурнагана [7], физико-механических свойств верхних слоев (мощностью 35 км) океанической версии физической модели Земли РЕМ-О (РЕМ – Parametric Earth Model, О – океаническая модель) [12]. Выбор уравнения состояния и входных данных параметров плотности, упругости и жесткости согласно модели Земли РЕМ–О определяется геофизическим обоснованием данных [12].

Рассматриваем любой вертикальный (горизонтальный) геопрофиль как 2D-фрагмент 3D-шара, каждая точка которого имеет текущие координаты: глубина, широта, долгота. Поэтому на боковых границах 2D-фрагментов дополнительных граничных условий не требуется. Расчет может быть продолжен в любом координатном направлении. Алгоритм предназначен для расчета НДС вертикальных и горизонтальных сечений 3D-шара по заданным координатам. 2D-фрагменты считаются в рамках 3D-модели, поэтому в расчетных значениях компонент параметров в текущих точках учитывается влияние компонент в смежных точках.

Полученное общее решение задачи устойчивого равновесия в аналитической форме [11], предполагает возможность встраивания параметра плотности $\wp \pm \Delta \wp$ модели [13], для построения которой использованы данные глубинных сейсмических зондирований [14–16], геолого-геофизических экспедиций, выполненных в 2005–2010 гг. на НИС «Академик Лаврентьев» [13, 17]. Плотностная модель задает численное значение плотности по определенной сетке глубины и длины профиля, по текущим координатам (глубина, широта, долгота). Используя аналитическое решение задачи при среднем значении плотности $\wp$ по данным модели РЕМ–О, можем рассчитать значения параметров НДС при значениях плотности $\wp \pm \Delta \wp$ по такой же сетке, где $\pm \Delta \wp$ – приращение плотности. Таким образом, в общее решение, построенное для средних значений плотности, упругости и жесткости среды для модели РЕМ–О «встраивается» распределение плотностей с привязкой к текущим координатам. Вклад осадков и водного слоя, стратификация распределения плотностей по глубине профилей в расчетах компонент и интенсивности напряжений учитывается фактическими значениями заданного поля плотностей.

В статье [18] по этому алгоритму представлены некоторые результаты расчета полей интенсивности гравитационных напряжений в земной коре по двум сейсмическим профилям района Центральных Курил: о. Уруп – о. Симушир – о-ва Ушишир и Охотское море – о. Симушир – Курило-Камчатский желоб (это основная часть профиля 2.1–2.10) по обобщенным данным, которые получены Институтом морской геологии и геофизики ДВО РАН [14–16]. Плотностные модели профилей, построенные на базе этих данных, морской и спутниковой гравиметрии соответствуют «усложненной слоисто-блоковой структуре» земной коры, установленной в результате глубинных сейсмических исследований [13]. Гравитационные силы и соответствующие распределения полей интенсивности напряжений в этих условиях, являются следствием плотностной неоднородности тектоносферы. Концентрация изолиний поля интенсивности напряжений соответствует положению глубинных разломов. Полученные расчетные данные полей интенсивности гравитационных напряжений по указанным профилям дают полное совпадение пространственной структуры плотностной неоднородности и полей интенсивности напряжений с учетом локальных аномалий распределения плотностей [18].

Пространственные изменения и локальные возмущения полей напряжений определяются плотностными границами слоев и блоков, выделенных по сейсмическим и гравитационным данным. В рамках алгоритма [9–11] по входным данным плотностных моделей выполнен детальный расчет, характеризующий напряженно-деформированное состояние двух профилей [13, 17], расположенных в квадрате: 45.0° ÷ 50.0° с.ш., 149.0° ÷ 156.0° в.д. Обобщение компонентного анализа расчетных значений перемещений, поворотов, деформаций и соответствующих значений компонент тензора напряжений по профилям 1.1–1.8 (длиной 700 км, глубиной 30 км) [17] и 2.1–2.10 (длиной 500 км, глубиной 35 км) [13] позволило определить общие свойства НДС всего квадрата, которое характеризуется растяжением на фоне преобладающего сжатия. Применение модельных результирующих параметров интенсивности напряжений и деформаций, как геомеханических маркеров, позволило выявить особенности распределения полей напряжений и возможного распределения свойства пластичности геосреды по профилям.

2. Расчетная модель изучения напряженно-деформированного состояния

Динамика неоднородных упругих сред рассматривается в лагранжевых переменных. Напряженно-деформированное состояние устойчивого равновесия определяется основным (невозмущенным) и возмущенным, каждое из которых описываются соответствующими линеаризованными соотношениями [4, 8]. Используем соотношения теории конечных деформаций и их линеаризованные выражения в произвольной ортогональной криволинейной системе координат. В приведенных ниже выражениях и далее все величины, отнесенные к основному состоянию, отмечены индексом “b”, величины без индексов относятся к возмущенному состоянию.

Основные соотношения для возмущенного состояния линеаризуются в окрестности значений соответствующих величин невозмущенного состояния. Ковариантные составляющие тензора конечных деформаций Коши–Грина основного ${\epsilon }_{ij}^b$ и возмущенного ε_ij состояний имеют вид

$$\begin{gathered} 2{\epsilon }_{ij}^b={\nabla }_i{u}_j^b+{\nabla }_j{u}_i^b+{\nabla }_i{u}_m^b\cdot {\nabla }_j{u}_b^m \\ 2{\epsilon }_{ij}=\left(g_i^m+{\nabla }_i{u}_b^m\right)\cdot {\nabla }_j{u}_m+\left(g_j^m+{\nabla }_j{u}_b^m\right)\cdot {\nabla }_i{u}_m, \end{gathered}$$

где ${\nabla }_i{u}_m^b, {\nabla }_j{u}_b^m$ – ковариантные производные от ковариантных и контравариантных составляющих вектора перемещений в основном состоянии, ∇_iu_m – ковариантная производная от ковариантных составляющих вектора перемещений в возмущенном состоянии, $g_i^m$ – метрический тензор в недеформированном состоянии.

В общем виде компоненты симметричного тензора обобщенных напряжений s^ij для сжимаемого нелинейно-упругого тела в невозмущенном $s_b^{ij}$ и возмущенном sⁱⁿ состоянии определены [8] выражениями соответственно

$$\begin{gathered} s_b^{ij}=\frac{1}{2}\left(\frac{\partial }{\partial {\epsilon }_{ij}^b}+\frac{\partial }{\partial {\epsilon }_{ji}^b}\right){\Psi }^b\left(A_1^b,A_2^b,A_3^b\right), s^{in}={\aleph }^{in\alpha \beta }{\nabla }_{\beta }{u}_{\alpha } \\ {\aleph }^{in\alpha \beta }=\frac{1}{4}\left[\left(g_m^{\alpha }+{\nabla }_m{u}_b^{\alpha }\right)\left(\frac{\partial }{\partial {\epsilon }_{m\beta }^b}+\frac{\partial }{\partial {\epsilon }_{\beta m}^b}\right)\left(\frac{\partial }{\partial {\epsilon }_{in}^b}+\frac{\partial }{\partial {\epsilon }_{ni}^b}\right)\right]{\Psi }^b\left(A_1^b, A_2^b, A_3^b\right), \end{gathered}$$ (2.1)

где ℵ^inαβ – тензор четвертого ранга.

Закон состояния ${\Psi }^b\left(A_1^b, A_2^b, A_3^b\right)$ представляет собой упругий потенциал и является дважды непрерывно-дифференцируемой функцией базисных алгебраических инвариантов $A_i^b$ тензора деформаций Грина в невозмущенном состоянии [11]

$A_1^b={\epsilon }_n^{nb}, A_2^b={\epsilon }_n^{mb}{\epsilon }_m^{nb}, A_3^b={\epsilon }_n^{mb}{\epsilon }_m^{kb}{\epsilon }_k^{nb}$

Напряженное состояние описывается несимметричным тензором Пиола–Кирхгофа, который определяется дифференцированием от энергии системы по градиенту вектора перемещений и характеризует элементарную силу, отнесенную к элементарной площадке в отсчетной конфигурации [6, 8]

$t^{im}=\left(g_n^m+{\nabla }_n{u}^{mb}\right)\cdot s^{in}+s_b^{in}\cdot {\nabla }_n{u}^m$, (2.2)

где тензоры обобщенных напряжений в возмущенном sⁱⁿ и невозмущенном $s_b^{in}$ состояниях определяются соответственно формулой (2.1) и выражением

$s_b^{in}=\left(g^{in}\frac{\partial }{\partial A_1^b}+2{\epsilon }_b^{in}\frac{\partial }{\partial A_2^b}+3{\epsilon }_b^{ik}{\epsilon }_k^{nb}\frac{\partial }{\partial A_3^b}\right){\Psi }^b\left(A_1^b, A_2^b, A_3^b\right)$

Линеаризованное уравнение равновесия объемного элемента в контравариантных компонентах несимметричного тензора напряжений Пиола–Кирхгофа имеет вид

${\nabla }_i{t}^{ij}+F^j=0$, (2.3)

где F^j – массовые силы.

Решение 3D-задачи получено в аналитической форме в системе координат O_ρφλ (ρ – радиус, φ – широта, λ – долгота) в физических составляющих компонент тензора деформаций Грина, несимметричного тензора напряжений Пиола–Кирхгофа, физических составляющих вектора линейных перемещений и поворотов. Физические компоненты параметров имеют индексы, определяющие их ориентацию в системе координат O_ρφλ (ρ, φ, λ).

3. Детализация напряженно-деформированного состояния профилей и обобщение этих результатов на заданный регион

Любая геодинамическая модель в приложении должна следовать из анализа структур, развивающихся в заданном геологическом регионе одновременно.

На рис. 1 приведены два профиля 1.1–1.8 и 2.1–2.10, расположенные в квадрате 149.0° ÷ 156.0° в.д., 45.0° ÷ 50.0° с.ш.

Рис. 1. Положение профилей 1.1–1.8 и 2.1–2.10 в квадрате 45.0° ÷ 50.0° с.ш., 149.0° ÷ 156.0° в.д. Цифровая модель рельефа построена по данным http://topex.ucsd.edu [19]

В рамках единой 3D-модели выполнен покомпонентный расчет всех полей, характеризующих напряженно-деформированное состояния заданных профилей. Вычислительный эксперимент выполнен для модели PEM-О, в которой выделен слой, моделирующий земную кору мощностью 35 км (ρ₁ = 6336 км, ρ₂ = 6371 км – верхняя поверхность). В соответствие с геофизическими свойствами PEM-О [12] используем в расчетах средне-взвешенные по мощности слоев значения параметров: упругой постоянной η = 54.777 ГПа, коэффициента жесткости μ = 53.837 ГПа, плотности $\wp$ = 3097.839 кг/м³, ускорения силы тяжести g = 9.839 м/с², внутреннего следящего давления $\tilde{p}$ = 1.82 ГПа.

Физический закон состояния, который определяет способность к деформированию в зависимости от свойств упругости и жесткости среды, принят в форме потенциала Мурнагана [7]

$\Psi =\frac{\eta }{2}{A}_1^2+\mu A_2$ (3.1)

Основное невозмущенное состояние является радиально-симметричным и устойчивым, характеризуется физическими компонентами радиальной деформации 0.0338 ≤ ${\epsilon }_{\left(\rho \rho \right)}^b$ ≤ 0.0370, радиальным напряжением 0.00565 ГПа ≤ $t_{\left(\rho \rho \right)}^b$ ≤ 0.0062 ГПа, удлинением 0.99981 ≤ δ_b ≤ 1.0 и радиальным перемещением –1.2 км ≤ $u_{\rho }^b$ ≤ 0 км [11].

На стадии возмущения в общее решение, построенное для модели РЕМ–О, «встраивается» неоднородное распределение плотностей с привязкой к текущим координатам [13, 17]. Кинематика деформированной геосреды описывается поступательными перемещениями, вектором линейного поворота, тензором конечных деформаций Коши–Грина. Компоненты вектора перемещений определяют вектор линейного поворота и компоненты тензора деформаций по известным формулам [5, 8]. В статье [20] опубликованы подробные численные результаты расчета полей перемещений, поворотов, деформаций и гравитационных напряжений по глубине 35 км для профиля 2.1–2.10 длиной 500 км, которые далее используются.

В табл. 1–3 для сравнения представлены интервальные численные значения компонент перемещений, поворотов, деформаций по двум профилям.

Таблица 1. Значения физических компонент u_ρ, u_φ, u_λ и модуля $\mathbf{\left|}\overrightarrow{\mathbf{u}}\mathbf{\right|}$ векторов перемещений по профилям 1.1–1.8 и 2.1–2.10

Таблица 2. Значения физических компонент w_ρ, w_φ, w_λ вектора поворота по профилям 1.1–1.8 и 2.1–2.10

Таблица 3. Значения интенсивности  и физических компонент тензора деформаций Коши-Грина по профилям 1.1–1.8 и 2.1–2.10

В табл. 1 приведены интервальные значения физических компонент u_ρ, u_φ, u_λ и модуля $\left|\overrightarrow{u}\right|$ векторов перемещений, определяющих поступательные перемещения текущих точек возмущенного напряженно-деформированного состояния по заданным профилям

$\left|\overrightarrow{u}\right|=\sqrt{u_{\rho }^2+u_{\phi }^2+u_{\lambda }^2}$

В текущих точках меридиональные компоненты u_λ > 0 положительные; радиальные и широтные компоненты u_ρ < 0 и u_φ < 0 имеют отрицательные направления. Между абсолютными значениями компонент векторов перемещений выявлена зависимость |u_ρ| << u_λ < |u_φ| по всей глубине профилей 1.1–1.8 и 2.1–2.10. Перемещения в радиальном направлении являются бесконечно малыми величинами |u_ρ| << 1, значит, в результате деформирования мощность слоев профилей остается практически без изменения. Поступательные перемещения текущих точек преобладают в латеральных направлениях и обусловлены компонентами u_φ и u_λ. При этом, абсолютные значения широтных перемещений u_φ больше, чем меридиональные u_λ.

Компоненты w_ρ, w_φ, w_λ определяют линейные повороты (табл. 2) в окрестности текущих точек в результате перемещений u_ρ, u_φ, u_λ на стадии возмущенного состояния. Расчеты показывают, что поворотные перемещения есть во всех координатных направлениях, отсутствуют нулевые значения, их порядок соизмерим со значениями компонент деформаций.

Компоненты тензора деформаций (табл. 3) полностью характеризуют деформацию среды в окрестности текущей точки в результате перемещений u_ρ, u_φ, u_λ, линейных поворотов w_ρ, w_φ, w_λ, косвенно определяя значения деформационных удлинений и сдвигов [5] и как результат соответствующие изменения размеров и формы.

Положительным линейным деформациям ε_(ρρ) соответствуют удлинения, а отрицательным ε_(φφ) и ε_(λλ) – укорочения. Положительные угловые деформации ε_(ρφ) и ε_(φλ) определяют уменьшение углов между положительными направлениями координатных осей (меньше 90°), а отрицательные значения ε_(ρλ) – увеличение тех же углов (больше 90°).

Учитывая то обстоятельство, что при очень медленных процессах в условиях всестороннего сжатия предполагается возможность пластического деформирования массивов пород [21], введем в рассмотрение положительный параметр интенсивности деформаций:

${\epsilon }_{\mathrm{int}}=\frac{\sqrt{2}}{3}\sqrt{{\left({\epsilon }_{\left(\rho \rho \right)}-{\epsilon }_{\left(\phi \phi \right)}\right)}^2+{\left({\epsilon }_{\left(\phi \phi \right)}-{\epsilon }_{\left(\lambda \lambda \right)}\right)}^2+{\left({\epsilon }_{\left(\lambda \lambda \right)}-{\epsilon }_{\left(\rho \rho \right)}\right)}^2+\frac{3}{2}\left({\epsilon }_{\left(\rho \phi \right)}^2+{\epsilon }_{\left(\phi \lambda \right)}^2+{\epsilon }_{\left(\lambda \rho \right)}^2\right)}$ (3.2)

Распределение интервальных значений компонент перемещений, поворотов и деформаций по профилям 1.1–1.8 и 2.1–2.2 имеет одинаковые качественные, но разные количественные характеристики. Максимальные перемещения и максимальные углы поворота в окрестности текущих точек характеризуют и определяют жесткость среды, способность сопротивляться деформации. Компоненты деформации, характеризуют удлинения и сдвиги бесконечно малого объемного элемента и определяют прочность среды, способность выдерживать без разрушения заданную нагрузку [5]. В данных условиях параметр перемещения |u_ρ| << 1 определяет бόльшую жесткость в радиальном направлении, угол поворота w_φ определяет бόльшую сопротивляемость деформации при повороте в плоскости O_ρλ, компонента деформации ε_(ρρ) << 1 определяет бόльшую прочность и способность выдерживать нагрузки в радиальном направлении. Весь комплекс перемещений, поворотов и деформаций обеспечивается соответствующим полем напряжений t^im (2.2).

Для изучения напряженного состояния заданного квадрата разложим симметричную часть sⁱⁿ (2.1) несимметричного тензора напряжений Пиола–Кирхгофа t^im (2.2). Для поверхности второго порядка в каждой текущей точке всегда можно выбрать такие направления координатных осей, для которых на трех взаимно перпендикулярных плоскостях касательные компоненты симметричного тензора напряжений sⁱⁿ обращаются в нуль. Эти плоскости называются главными плоскостями, а напряжения (σ_I, σ_II, σ_III), действующие на них – главными напряжениями, т.е. в каждой текущей точке вместо девяти компонент тензора напряжений получаем три компоненты с важным физическим свойством. Главные напряжения определяют экстремальные направления всего комплекса напряжений. На главных площадках главные напряжения принимают свои экстремальные значения – максимум σ_I, минимум σ_III и минимакс σ_II (σ_I > σ_II > σ_III) [22]. Если в текущей точке известны все компоненты тензора напряжений, то соответствующие главные направления и главные напряжения вычисляются через решение кубического характеристического уравнения, составленного на основе определителя компонент тензора sⁱⁿ [5]. По определению, растягивающее напряжение считается положительным, сжимающее – отрицательным.

На рис. 2 представлено двумерное графическое представление трехмерного напряженного состояния в текущих точках по профилям 1.1–1.8 и 2.1–2.10 (круги Мора).

Рис. 2. Круги Мора в текущих точках по профилям 1.1–1.8 и 2.1–2.10

Предполагаем, что возможно обобщение численного анализа результатов решения характеристического уравнения по определению главных напряжений по двум заданным профилям на весь заданный квадрат. По результатам расчетов напряженное состояние в заданном квадрате характеризуется в одном направлении растяжением на фоне двухосного сжатия: максимальное растяжение (σ_I > 0), сжатие (σ_II < 0), максимальное сжатие (σ_III < 0).

При этом абсолютные значения главных напряжений относительно пропорциональны: σ_I ≈ ½ |σ_III|, |σ_II| ≈ ⅛ |σ_III|. Главное напряжение в зоне сжатия почти в 2 раза превышает такой же параметр в зоне растяжения. Таким образом, можно предположить, что общим свойством напряженно-деформированного состояния всего квадрата 149.0° ÷ 156.0° в.д., 45.0° ÷ 50.0° с.ш. является растяжение на фоне преобладающего сжатия.

4. Детализация и особенности напряженно-деформированного состояния профилей

Кинематика деформированной геосреды по сечению профилей 1.1–1.8 и 2.1–2.10 (табл. 1–3) обеспечивается соответствующим полем напряжений, при этом должна быть выражена согласованность с данными о глубинном строении земной коры, положением плотностных границ и разломов.

Для анализа особенностей распределения полей деформаций и гравитационных напряжений по заданным профилям, рассмотрим распределение модельных положительных параметров интенсивности деформаций (3.2) и интенсивности напряжений (4.1):

$T_{\mathrm{int}}=\frac{1}{\sqrt{2}}\sqrt{{\left(t_{\left(\rho \rho \right)}-t_{\left(\phi \phi \right)}\right)}^2+{\left(t_{\left(\phi \phi \right)}-t_{\left(\lambda \lambda \right)}\right)}^2+{\left(t_{\left(\lambda \lambda \right)}-t_{\left(\rho \rho \right)}\right)}^2+3\left(t_{\left(\rho \phi \right)}^2+t_{\left(\phi \lambda \right)}^2+t_{\left(\rho \lambda \right)}^2+t_{\left(\phi \rho \right)}^2+t_{\left(\lambda \phi \right)}^2+t_{\left(\lambda \rho \right)}^2\right)}$ (4.1)

На рис. 3 показаны графики градиентного распределения интенсивности деформаций ε_int по плоскости профилей 1.1–1.8 (рис. 3, а) и 2.1–2.10 (рис. 3, б). Распределение интенсивности деформаций моделирует степень пластического состояния геосреды профилей согласно заданной плотностной дифференциации.

Рис. 3. Реконструкция поля интенсивности деформаций ε_int по профилям: а: 1.1–1.8 значения интенсивности деформаций 2.07 ⋅ 10⁻⁵ ≤ ε_int ≤ 2.35 ⋅ 10⁻⁵, горизонтальные границы раздела плотностей обозначены штриховой и сплошной линиями, в верхней части разломы соответствуют модели [17]; б: 2.1–2.10 значения интенсивности деформаций 2.04 ⋅ 10⁻⁵ ≤ ε_int ≤ 2.44 ⋅ 10⁻⁵, значения плотностей, границы плотностных блоков и разломы соответствуют модели [13]

Максимальные значения интенсивности деформаций ε_int концентрируются в зонах (окрашены в серый градиентный цвет), в которых возможность перехода в пластическое состояние обусловлена напряженно-деформированным состоянием геосреды.

В плоскости профиля 1.1–1.8 (рис. 3, а) степень пластичности преобладает в верхних слоях срединной части земной коры в окрестности разломов и горизонтальных границ раздела: мантия нижняя кора – верхняя кора [17]. В плоскости профиля 2.1–2.10 (рис. 3, б) степень пластичности преобладает в верхних слоях земной коры в окрестности разломов. Распределение максимальных значений поля ε_int оконтуривает окрестности разломов, плотностных границ и определяет их зону влияния. Изменения значений интенсивности деформаций по всей плоскости профиля 1.1–1.8 определяются интервалом 2.07 ⋅ 10⁻⁵ ≤ ε_int ≤ 2.35 ⋅ 10⁻⁵; профиля 2.1–2.10 – интервалом 2.04 ⋅ 10⁻⁵ ≤ ε_int ≤ 2.44 ⋅ 10⁻⁵. Максимальное значение интенсивности деформаций по сечению профиля 2.1–2.10 равно 2.44 ⋅ 10⁻⁵. Оно больше максимального значения интенсивности деформаций по сечению профиля 1.1–1.8, но порядок значения тот же, 2.44 ⋅ 10⁻⁵ > 2.35 ⋅ 10⁻⁵. Учитывая характеристики главных напряжений, это сравнение, возможно, косвенно определяет положение профиля 2.1–2.10 в зоне сжатия, а профиля 1.1–1.8 в зоне растяжения.

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

На рис. 4 показаны графики изолиний и градиентного распределения интенсивности напряжений T_int по плоскости профилей 1.1–1.8 (рис. 4, а) и 2.1–2.10 (рис. 4, б), которые согласуются с положением глубинных разломов и плотностных границ слоев [13, 18].

Рис. 4. Распределение изолиний и градиентных полей интенсивности напряжений T_int по профилям: а: 1.1–1.8 (длиной 700 км, глубиной 30 км), горизонтальные границы раздела плотностей обозначены штриховой и сплошной линиями, в верхней части разломы отмечены сплошной линией, соответствуют модели [17]; б: 2.1–2.10 (длиной 500 км, глубиной 35 км), границы плотностных блоков и разломы соответствуют модели [13]. Пересечение профилей 1.1–1.8 и 2.1–2.10 обозначено белой сплошной линией

Изолинии концентрируются и преломляются там, где локализованы зоны раздела плотностей и разломы. Максимальное значение поля интенсивности напряжений по профилю 2.1–2.10 (рис. 4, б) достигает 1790 ГПа на юго-востоке на уровне мантии. Максимальное значение поля интенсивности напряжений по профилю 1.1–1.8 (рис. 4, а) достигает 1690 ГПа на уровне мантии в средней части разреза. Пересечение профилей обозначено белой сплошной линией. Интервальные значения интенсивности напряжений для профиля 1.1–1.8 (рис. 4, а) составляют 1516 ГПа ≤ T_int ≤ 1690 ГПа, они меньше значений этого же параметра для профиля 2.1–2.10 (рис. 4, б) 1530 ГПа ≤ T_int ≤ 1790 ГПа. Характеристики главных напряжений σ_i, значения параметров интенсивности напряжений и градиентное распределение интенсивности напряжений (рис. 4) определяют положение профиля 2.1–2.10 в зоне сжатия, а профиля 1.1–1.8 в зоне растяжения. Возможно, что на рис. 4, а на глубине ниже отметки –20 км на протяжении – 40 км < L < 200 км моделируется часть океанической плиты профиля 2.1–2.10, НДС которой характеризуется сжатием с юго-востока на северо-запад (рис. 4, б).

Предполагаем, что динамика возможных структурных изменений обусловлена не только значениями интенсивности напряжений, но и распределением градиентов ΔT_int этой функции относительно смежных точек. На рис. 5 приведены графики полей градиентов интенсивности напряжений ΔT_int по плоскостям профиля 1.1–1.8 (рис. 5, а) и фрагменту профиля 2.1–2.10 длиной 100 км ≤ L ≤ 350 км (рис. 5, б).

Рис. 5. Распределение изолиний, градиентов ΔT_int поля интенсивности напряжений ΔT_int и оконтуренных линией I модельных «блоков», внутри которых ΔT_int ≤ 2.5 ГПа/км, по профилям: а: 1.1–1.8, в верхней части разломы отмечены сплошной линией и соответствуют модели [17]; б: фрагмент профиля 2.1–2.10 длиной 100 км ≤ L ≤ 350 км, границы плотностных блоков и разломы соответствуют модели [13]. Звездочками обозначены гипоцентры Симуширских землетрясений 2006, 2007 гг.

Векторные поля градиентов интенсивности напряжений определяют направление и величину наискорейшего изменения интенсивности $\Delta T_{{\mathrm{int}}_{\max }}$ в расчетных точках профилей.

Предполагаем, что малые значения градиентов интенсивности напряжения в смежных точках определяют более плотную и менее подвижную часть геосреды. На рис. 5 контуром I выделены «пустоты» в векторной сетке – модельные «блоки», внутри которых ΔT_int ≤ 2.5 ГПа/км. И наоборот, менее плотная и более подвижная часть геосреды, например, разломные зоны, характеризуются градиентами интенсивности напряжений ΔT_int ≥ 2.5 ГПа/км.

Межблоковое пространство является неустойчивой составляющей геосреды, характеризуется бóльшими значениями градиентов интенсивности напряжений, которые изменяются по величине и направлению. При этом моделируются зоны гравитационной неустойчивости с повышенной способностью к структурообразованию, которые могут изменять границы глубинных разломов и раздела плотностей [23].

На рис. 5, а распределение векторного поля градиентов интенсивности напряжений моделирует систему из двух горизонтальных массивных модельных «блоков» по длине 700 км, разделенных горизонтальной искривленной границей. Возможно, что это граница раздела плотностей геоструктуры профиля 1.1–1.8. При дополнительном тектоническом воздействии на эту модельную систему: «блок» – межблоковое пространство – «блок» (рис. 5, а) в условиях растяжения, возмущения геосреды и развитие гравитационной неустойчивости затруднительно.

На рис. 5, б представлен наиболее информативный, с точки зрения зафиксированных геособытий, фрагмент 2.1–2.10 длиной 100 км ≤ L ≤ 350 км. Распределение векторного поля градиентов интенсивности напряжений моделирует систему из восьми соразмерных модельных блоков, расположенных на глубине −10 ± 3 км вдоль протяженных разломных зон. Они разделены вертикальными межблоковыми пространствами. При дополнительном тектоническом воздействии на эту модельную систему: «блок» – межблоковое пространство – «блок» (рис. 5, б), которая находится в зоне сжатия, возможно ожидать возмущений геосреды и развития гравитационной неустойчивости. Эти возмущения может служить механизмом, определяющих ход различных геодинамических процессов, в том числе и триггером землетрясений. На рис. 5, а гипоцентры землетрясений (обозначены звездочками) моделируются между соразмерными блоками, в зоне влияния разломов и плотностных границ, которые характеризуются восходящими и нисходящими пиками изолиний интенсивности напряжений, что способствует развитию здесь изгибов слоев с нарушением сплошности или формированием систем трещин и могут быть зоной повышенной сейсмичности. По опубликованным данным гипоцентры Симуширских землетрясений 2006, 2007 гг. располагались на глубине порядка 10 км в зоне разломов основания коры [15, 16].

Параметр интенсивности напряжений предложен как геомеханический маркер возможного взаимодействия неоднородных геоструктур, так как распределение относительно малых по величине градиентов интенсивности напряжений позволяют смоделировать геосреду как систему элементов взаимодействия: «блок» – межблоковое пространство – «блок», где «блок» относительно (локально) устойчивая часть, «межблоковое пространство» – неустойчивая часть геоструктуры. При этом распределение градиентов интенсивности напряжений, которые характеризуются большими значениями, моделируют направление и величину наискорейшего изменения поля интенсивности и последующие возможные геодинамические изменения.

5. Обсуждение результатов и выводы

В едином 3D-пространстве, структурированного соотношениями (2.1)–(2.3), представлено исследование НДС усложненных сред с заданной способностью к деформированию [7], неоднородным распределением физико-механических свойств по данным параметрической модели РЕМ–О [12], распределением плотностей, полученных на основании обобщения данных мониторинга геологической обстановки Центральных Курил [13, 17]. Введенные в рассмотрение положительные параметры интенсивности деформаций ε_int и напряжений T_int являются модельными результирующими характеристиками компонент тензора деформаций и напряжений в текущей точке. Очевидно, что в естественных условиях работает иной механизм определения результирующих характеристик деформирования и напряжений. Но уместность применения в моделировании результирующих параметров деформирования и напряжений обусловлена возможностью получения пространственного распределения параметров, характеризующих особенности НДС в профилях.

Количественное и качественное распределение возможных полей напряжений и деформаций обусловлено расчетной моделью (2.1)–(2.3), граничными условиями, уравнением состояния (4), входными данными полей физико-механических свойств геосреды [12, 13, 17]. Получены согласованные поля распределения значений и градиентов компонент тензоров, результирующих параметров, характеризующих напряженно-деформированное состояние.

Здесь представлены результаты расчетов по двум профилям в квадрате Центральных Курил, в условиях плотностной неоднородности. Использованы результаты структурно-плотностного моделирования земной коры [13, 17]. Плотностная модель определяет важные свойства типов геовещества в текущих точках, распределенных по всей плоскости профилей. Входными данными являются значения плотности в зависимости от распределения давления, температуры и текущих координат. Встраивая эти данные плотностей в расчетную механико-математическую модель, в которой учтены распределения упругости и жесткости получаем реконструкцию возможного динамического 3D-взаимодействия из системы элементов неоднородной геосреды: «блок» – межблоковое пространство – «блок». Вклад осадков и водного слоя, стратификация распределения плотностей по глубине в расчетах компонент НДС учитывается фактическими значениями заданного поля плотностей. На основании обобщения численного анализа результатов НДС по двум профилям 1.1–1.8 (длиной 700 км, глубиной 30 км) и 2.1–2.10 (длиной 500 км, глубиной 35 км) на весь квадрат 149.0° ÷ 156.0° в.д., 45.0° ÷ 50.0° с.ш.,  установлены следующие общие характеристики НДС всего квадрата:

1. Перемещения в радиальном направлении являются бесконечно малыми величинами |u_ρ| << 1, значит, в результате деформирования мощность слоев заданного региона остается практически без изменения. Поступательные перемещения по глубине обусловлены компонентами u_φ и u_λ в латеральных направлениях, при этом преобладают широтные компоненты перемещений;
2. Компоненты линейных поворотов в результате перемещений фиксируются во всех координатных направлениях (отсутствуют нулевые значения), их порядок соизмерим со значениями компонент деформаций;
3. Компоненты тензора деформаций Коши–Грина полностью характеризуют деформацию среды в окрестности текущей точки в результате перемещений, линейных поворотов, косвенно определяя значения деформационных удлинений и сдвигов и как результат соответствующие изменения размеров и формы. Положительным линейным деформациям ε_(ρρ) соответствуют удлинения, а отрицательным ε_(φφ) и ε_(λλ) – укорочения. Положительные угловые деформации ε_(ρφ) и ε_(φλ) определяют угол между положительными направлениями координатных осей меньше 90°, а отрицательные значения ε_(ρλ) – тот же угол больше 90°;
4. Компоненты несимметричного тензора напряжений Пиола–Кирхгофа полностью обеспечивают кинематику, описанную в пп. 1–3 этого списка. Предложенный параметр интенсивности напряжений является результирующей положительной функцией всех девяти компонент тензора напряжений Пиола–Кирхгофа в текущей точке. Чем больше значение параметра интенсивности напряжений в определенном направлении, тем более вероятно в этом месте относительное усиление движения и сопутствующего деформирования. Направление усиления/ослабления движения определяется вкладом компонент тензора t_(φλ) и t_(λφ). Вклад касательных t_(ρλ) и t_(λρ) значительно меньше, чем значения компонент t_(φλ), t_(λφ), t_(ρρ), t_(φφ), t_(λλ), t_(ρφ) и t_(φρ) [20, 23];
5. На основании алгоритма построения главных напряжений по данным расчетов текущих значений компонент тензора напряжений Пиола–Кирхгофа по двум профилям, получена общая характеристика напряженного состояния всего квадрата 149.0° ÷ 156.0° в.д., 45.0° ÷ 50.0° с.ш.: растяжение на фоне преобладающего сжатия (максимальное растяжение–сжатие–максимальное сжатие).

Для выявления особенностей распределения НДС по двум профилям 1.1–1.8 и 2.1–2.10 используем те же расчетные данные компонент перемещений, поворотов, деформаций и напряжений. Эти данные информативны в плане общих характеристик НДС заданного квадрата, при этом особенности распределения полей НДС по геоструктуре профилей согласуются и обеспечиваются общими характеристиками, приведенными в списке пп. 1–5.

Установлены следующие особенности распределения НДС по профилю 1.1–1.8:

1. Рельеф поверхности профиля определяется напряженно-деформированным состоянием глубинных структур.
2. На рис. 4а пространственное распределение изолиний интенсивности напряжений моделирует положение профиля в зоне растяжения, на глубине ниже –20 км на протяжении –40 км < L < 200 км моделируется часть океанической плиты профиля 2.1–2.10, НДС которой характеризуется сжатием с юго-востока на северо-запад.
3. Распределение маркера интенсивности напряжений моделирует устойчивую слоисто-блоковую структуру по плоскости всего профиля.
4. Распределения кинематических параметров определяют возможные качественные изменения вещественных элементов геосреды профиля (список пп. 1–5). Параметр интенсивности деформаций, как геомеханический маркер степени пластического деформирования неоднородной геосреды по глубине профиля, моделирует концентрацию максимальных значений в зонах локализации разломов и плотностных границ.

Установлены следующие особенности распределения НДС по профилю 2.1–2.10:

1. На рис. 4, б пространственное распределение изолиний интенсивности напряжений и более насыщенная окраска на юго-восточной границе профиля моделирует потенциальную возможность бокового давления «поддвигающейся Тихоокеанской плиты, создающей зону сжатия, повышение плотности геологической среды в пределах уклона Курильского желоба и, как следствие, возникновение зоны повышенных напряжений. По мере удаления от указанной зоны в сторону Охотского моря, величина интенсивности напряжений, естественно, должна уменьшаться…» [18], что подтверждается более светлой окраской верхней части вдоль всего профиля на рис. 4, б.
2. Распределение градиентов маркера интенсивности напряжений моделирует устойчивые/неустойчивые структурные элементы в плоскости профиля относительно возможного динамического взаимодействия: устойчивые структурные элементы – «блоки», которые разделены зонами с повышенной активностью – «межблоковыми пространствами». Динамическое воздействие на межблоковые пространства в системе из восьми соразмерных модельных блоков, расположенных на глубине −10 ± 3 км (рис. 5, б), моделирует возможность сейсмических деформаций. Распределение, детализация слоисто-блоковой структуры (расположение, соизмеримость блоков взаимодействия) в гипоцентрах Симуширских землетрясений, а также концентрация полей деформаций и напряжений согласуются с плотностной дифференциацией земной коры и свидетельствуют о прошлой и возможной в будущем сейсмической активности в этом регионе.
3. Рельеф поверхности профиля определяется напряженно-деформированным состоянием глубинных структур.
4. Распределения кинематических параметров определяют возможные качественные изменения вещественных элементов геосреды профиля (список пп. 1–5). Применен параметр интенсивности деформаций, как геомеханический маркер степени пластического деформирования неоднородной геосреды по глубине профиля. Максимальные значения маркера концентрируются в зонах локализации разломов по всей плоскости профиля, при этом распределение максимальных значений оконтуривают окрестности разломов, плотностных границ и определяют их зону влияния.

Полученные результаты позволяют реконструировать и объяснить особенности структурной эволюции неоднородных геосред в поле собственной гравитации.

В рамках расчетной 3D-модели неоднородность свойств (упругость, жесткость) геосреды может быть учтена через приращение ±Δ. Закон состояния, определяющий способность к усложненному деформированию, учитывающий взаимодействие компонент деформаций и скоростей деформаций, может быть принят в форме Φ^ij (ε^ij, ξ^ij). Но автор не располагает данными неоднородного распределения свойств упругости, жесткости, вязкости и других свойств геосред.

Работа выполнена в рамках госбюджетной темы «Изучение структуры, физических и вещественных характеристик и геодинамики литосферы, сейсмической активности и закономерностей размещения полезных ископаемых в регионе дальневосточных морей и северо-западном секторе Тихого океана», регистрационный номер: 124022100082-4.

About the authors

E. B. Osipova

Author for correspondence.
Email: osipov@poi.dvo.ru
Russian Federation, Vladivostok

References

Supplementary files

Supplementary Files

Action

1. JATS XML

Download

2. Fig. 1. Position of profiles 1.1-1.8 and 2.1-2.10 in the square 45.0° ÷ 50.0°N, 149.0° ÷ 156.0°E. The digital elevation model was built using data from http://topex.ucsd.edu [19]

Download (51KB)

Indexing metadata

3. Fig. 2. Mora circles at current points along profiles 1.1-1.8 and 2.1-2.10

Download (13KB)

Indexing metadata

4. Fig. 3. Reconstruction of the εint strain intensity field along the profiles: a: 1.1-1.8 strain intensity values 2.07 ⋅ 10-5 ≤ εint ≤ 2.35 ⋅ 10-5, horizontal density boundaries are indicated by dashed and solid lines, faults in the upper part correspond to the model [17]; b: 2. 1-2.10 deformation intensity values 2.04 ⋅ 10-5 ≤ εint ≤ 2.44 ⋅ 10-5, density values, density block boundaries and faults correspond to the model [13]

Download (53KB)

Indexing metadata

5. Fig. 4. Distribution of isolines and gradient fields of Tint stress intensity along the profiles: a: 1.1-1.8 (length 700 km, depth 30 km), horizontal boundaries of density interfaces are marked by dashed and solid lines, in the upper part faults are marked by solid line, correspond to the model [17]; b: 2.1-2.10 (length 500 km, depth 35 km), boundaries of density blocks and faults correspond to the model [13]. The intersection of profiles 1.1-1.8 and 2.1-2.10 is marked by the white solid line

Download (83KB)

Indexing metadata

6. Fig. 5. Distribution of isolines, gradients ΔTint of the stress intensity field ΔTint and model "blocks" delineated by line I, within which ΔTint ≤ 2.5 GPa/km, along the profiles: a: 1.1-1.8, in the upper part the faults are marked by a solid line and correspond to the model [17]; b: fragment of profile 2.1-2.10 with the length of 100 km ≤ L ≤ 350 km, the boundaries of density blocks and faults correspond to the model [13]. Asterisks indicate hypocentres of the Simushir earthquakes of 2006, 2007

Download (101KB)

Indexing metadata