Доклады Российской академии наук. Науки о Земле, 2020, T. 495, № 2, стр. 51-55

ВВЕДЕНИЕ

Исходными данными для построения плотностных моделей земной коры и верхней мантии являются измеренные значения аномалий силы тяжести. Одним из этапов предварительной обработки данных гравитационных измерений является расчет поправки за топографию [1, 2]. Он состоит в вычислении гравитационного эффекта масс (в точке наблюдения), лежащих выше некоторой референц-поверхности, например, земного эллипсоида. “Классической” аппроксимацией этого этапа (для “плоской” модели) можно назвать вычисление поправки за влияние промежуточного плоскопараллельного слоя, входящей в поправку Буге [1]. Для небольших площадей со слабо выраженными особенностями рельефа данный подход хорошо приближает “поле рельефа”, но в настоящее время область его применимости ограничена в связи с возросшими требованиями к результирующей погрешности в вычисляемом гравитационном поле. При рассмотрении региональных моделей земной коры большой протяженности (порядка 1000 × 1000 км) [3] требуются более точные методы введения поправок: недостаточно точное приближение рельефа и неучет сферической формы планеты могут внести ошибку [4‒6], которая на два десятичных порядка превышает чувствительность современных гравиметров [7]. Отметим, что предварительное вычисление поправок необходимо для получения значений гравитационного поля на регулярной сетке, удобной для реализации численных алгоритмов интерпретации. Цель данной работы: оценить погрешность вычисления поправки за топографию (по формуле плоскопараллельного слоя и с помощью решения прямой задачи гравиметрии от модели рельефа с учетом и без учета “сферичности”) и предложить универсальный метод интерпретации гравитационных данных, не требующий вычисления поправок.

ПОПРАВКА ЗА ВЛИЯНИЕ ПРОМЕЖУТОЧНОГО ПЛОСКОПАРАЛЛЕЛЬНОГО СЛОЯ

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

Расчет ведется по формуле:

где G – гравитационная постоянная, ρ – плотность слоя, h – превышение рельефа над рефренц-поверхностью. Погрешность, которая возникает при ее использовании, складывается из двух частей: плохая аппроксимация реального поля в горных областях и ошибка за неучет “сферичности”. Вторая составляющая начинает проявляться только при расчете поправки для относительно больших по площади территорий [6].

РАСЧЕТ ПОПРАВКИ ЗА ТОПОГРАФИЮ ПУТЕМ РЕШЕНИЯ ПРЯМОЙ ЗАДАЧИ ГРАВИМЕТРИИ С УЧЕТОМ И БЕЗ УЧЕТА “СФЕРИЧНОСТИ”

В случае существенных возмущений рельефа на целевой площади и необходимости учета “сферичности” требуется перейти к численному решению прямой задачи гравиметрии для модели земного рельефа [8]. Ранее П.С. Мартышко, A.I. Chernoskutov с соавторами [9, 10] разработали высокоэффективный алгоритм численного решения прямой задачи для моделей произвольной формы. Алгоритм основан на приближении элементов разбиения модели многогранниками и последующем суммировании их гравитационных эффектов в точке наблюдения. Поле от многогранника вычисляется по аналитической формуле [6]. За счет такого подхода достигается высокая точность решения, которая зависит только от точности приближения моделируемого объекта многогранниками. Данный метод можно применить для вычисления поля от модели, ограниченной поверхностью рельефа. Для этих же целей в [11] предложено использовать аппроксимацию сферическими параллелепипедами.

Определим “эллипсоидальную” трехмерную плотностную модель, поле которой и будем считать “влиянием рельефа”. Пусть “нижняя” граница модели – поверхность ${{H}_{b}}\left( {L,~B} \right)$, “верхняя” граница (со стороны раздела земля‒воздух) – поверхность рельефа ${{H}_{t}}\left( {L,B} \right)$. Функции ${{H}_{b}}\left( {L,~B} \right)$ и ${{H}_{t}}\left( {L,~B} \right)$ определены как превышения над поверхностью некоторого референц-эллипсоида вращения, $B \in \left[ { - \frac{\pi }{2};\frac{\pi }{2}} \right]$ – геодезическая широта, $L \in \left( {\left. { - \pi ;\pi } \right]} \right.$ – долгота, связанные с эллипсоидом. Все точки, расположенные между H_b и H_t, включены в модель. В указанной области $D \subset {{\mathbb{R}}^{3}}$ задано распределение плотности ρ(p), $p \in D$.

Выберем некоторое разбиение $D = \bigcup\nolimits_{i = 1}^N {{{D}_{i}}} $. Аппроксимируя элементы разбиения ${{D}_{i}}$ многогранниками ${{\hat {D}}_{i}}$, получаем модель, к которой можно применить ранее разработанный алгоритм решения прямой задачи гравиметрии [6].

ПРИМЕР ВЫЧИСЛЕНИЯ ПОПРАВКИ ДЛЯ РЕГИОНАЛЬНОЙ МОДЕЛИ

Для проведения численных экспериментов были взяты цифровые модели рельефа и геоида с онлайн-ресурса ICGEM [12] для исследованной ранее территории 60°–68° с.ш. 48°–72° в.д. [9] в разрешении 200 × 600 точек. В этих файлах высоты рельефа в узловых точках заданы превышением над геоидом, а высоты геоида, в свою очередь, заданы как превышение над эллипсоидом WGS84. Этот эллипсоид был использован в качестве референц-эллипсоида и в качестве нижней границы H_b. Отметим, что в качестве H_b может быть выбрана поверхность геоида (если того требует конкретная задача), это не приведет к изменениям в предлагаемом алгоритме или дополнительным вычислительным расходам. Для вычисления поправки необходимо рассчитать гравитационное поле от слоя, ограниченного H_b “снизу” и рельефом “сверху”, для чего понадобится цифровая модель рельефа, заданная как превышения над эллипсоидом. На рис. 1 показана получившаяся поверхность рельефа. Высоты результирующей модели рельефа лежат в диапазоне [0; 1,44] км.

Для целевой территории были вычислены поправка за промежуточный плоскопараллельный слой (1) и поправка за топографию двумя способами: через решение прямой задачи в “плоском” и “сферическом” случаях. Расчеты проводились с постоянной плотностью ρ = 2.67 г/см³. Отметим, что методы поддерживают задание произвольной кусочно-постоянной функции плотности для “модели рельефа”. Минимумы и максимумы результирующих поправок, мГал.: по формуле (1) – [0.11; 161.96], решение прямой задачи (“плоский” случай) ‒ [‒0.3; 154.9], решение прямой задачи (“сферический” случай) ‒ [0.34; 155.92]. На рис. 2 и 3 приведены графики их разности. Уменьшение разности поправок с учетом и без учета “сферичности” при приближении к краям области, которое видно на рис. 3, вызвано понижением уровня рельефа на краях. Во избежание краевого эффекта мы выполняли расчеты для расширенной, по сравнению с целевой территорией, области. Сравнение результатов расчета представлено в табл. 1.

При расчете поправок в приведенном примере использовалось программное обеспечение, основанное на разработанной авторами программе GRAFEN [10]. Время счета составило ~5 мин с использованием двух GPUNVidia Quadro M6000. Программа реализует самый быстрый (из применяемых в настоящее время) алгоритм вычисления гравитационного поля на нерегулярной сетке от объектов произвольной формы [6].

ЗАКЛЮЧЕНИЕ

Проведенные сравнения показывают, что для рассматриваемой территории ошибка, возникающая при учете влияния топографии с использованием формулы (1) и ошибка “за сферичность” примерно одного порядка. В приведенном примере среднеквадратичная ошибка “простого метода” (1) составила 3.1%, что может существенно повлиять на качество интерпретации гравитационного поля, особенно при восстановлении физических параметров “верхнего” гравиактивного слоя. Таким образом, при необходимости высокоточных расчетов гравитационного влияния неоднородного рельефа для больших по площади территорий необходимо учитывать оба фактора. Разработанный П.С. Мартышко с соавторами алгоритм [6, 9] позволяет вычислять значения гравитационного поля на реальном рельефе без дополнительных затрат вычислительных ресурсов: эффективность метода не зависит от регулярности сеток плотностной модели и вычисляемого поля. Следовательно, интерпретацию наблюденных гравитационных данных с учетом “сферичности” можно проводить без предварительного вычисления поправок за топографию, используя для построения плотностных моделей аномалии в свободном воздухе. При этом будут определяться параметры плотностной модели, ограниченной сверху поверхностью рельефа. Такой подход позволит существенно повысить точность определения параметров моделей.