Физика Земли, 2022, № 2, стр. 3-26

СТАНОВЛЕНИЕ МЕТОДА

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

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

Если рассматривать истокообразные аппроксимации как средство моделирования и исходить из принципов классификации математических моделей гравитационного поля Земли, приведенных в работе [Гравиметрия и геодезия, 2010], то в данной статье речь пойдет преимущественно о детерминированных, дискретных, 2D- (площадных) компьютерных моделях наблюденных аномальных полей. Эти модели предназначены для проведения вычислительных экспериментов, сжатия и отображения информации. Они реализуются в виде комплекса программных средств и связанных с ним баз данных [Непоклонов, 1998]. Область действия рассматриваемых моделей может быть различной: для глобальных моделей – 10⁹ км² и более; для региональных моделей – от 10⁸ до 10⁶ км²; для локальных моделей – от 10⁵ до 10² км² и менее.

Исторически возникновение метода истокообразной аппроксимации геопотенциальных полей относится ко времени проведения активных исследований для задачи определения фигуры Земли по наблюденным значениям силы тяжести, с 30-х по 50-е годы XX столетия. Концептуальная основа буквально всех существующих методов истокообразной аппроксимации геопотенциальных полей во многом заимствована из этой задачи геодезической гравиметрии. Однако, по мнению В.М. Гордина и С.А. Тихоцкого [Гордин, Тихоцкий, 1999], сами идеи, связанные с истокообразной аппроксимацией полей, появились в физике гораздо раньше, еще в XVII–XVIII веках.

М.С. Молоденский, подводя итоги исследований по определению фигуры Земли, рассматривал возмущающий потенциал $W$ на поверхности реальной (нерегуляризованной, не редуцированной к уровню геоида) Земли в виде суммы

в которой первое слагаемое $U$ представляет собой потенциал теоретической Земли (сфероида или референц-эллипсоида) и определяется обобщенной формулой Стокса, а второе – возмущающий потенциал $T$ – приближается потенциалом простого слоя [Молоденский, 1945]. Продифференцировав обе части выражения (1.1), его можно переписать в виде: или где: $n$ – нормаль к поверхности референц-эллипсоида из заданной точки ${{P}_{i}}$ на поверхности Земли; $\Delta g = g - {{\gamma }}$ – аномальное значение силы тяжести; $g$ – наблюденное значение силы тяжести; ${{\gamma }}$ – нормальное значение силы тяжести [Молоденский, 1948; 1949].

Перефразируя М.А. Алексидзе [Алексидзе, 1987, с. 58], можно сказать, что М.С. Молоденский концептуально сформулировал подход к решению задачи истокообразной аппроксимации возмущающего потенциала T в виде вопроса: “Нельзя ли подобрать такие значения аномалий возмущающего потенциала на сфере внутри геоида, чтобы, будучи продолженными на поверхность реальной Земли, они полностью совпали с измеренными значениями потенциала?”. Положительный ответ на этот вопрос дал сам М.С. Молоденский в работе [Молоденский, 1945] в виде сложного линейного интегрального уравнения типа Фредгольма, при этом оставив его решение другим исследователям. Впрочем, несколько позже, вводя понятие “квазигеоида”, М.С. Молоденский приводит приближенный способ решения полученного уравнения [Молоденский, 1948; 1949].

В обозначенный период задачу М.С. Молоденского в основном пытались решить путем прямого продолжения (проецирования) гравитационного поля с поверхности Земли на внутреннюю сферу, в результате получив некую гармоническую функцию ${{\varphi }}$. Несмотря на то, что ${{\varphi }}$ иногда называли поверхностной плотностью, ее рассматривали лишь как проекцию на сферу гармонической функции ${{V}_{z}}$ или $\Delta g$, заданной измеренными значениями на физической поверхности Земли. Большой цикл работ на данную тему опубликовал выдающийся норвежский ученый А. Bjerhammar, ссылки на некоторые из них приводятся в списке литературы [Bjerhammar, 1962; 1963; 1965; 1969a; 1969b; 1975; 1985]. Однако выполненные им исследования, как и у М.С. Молоденского, лежали в области геодезической гравиметрии, а не прикладной (разведочной) геофизики, и были также посвящены задаче определения фигуры Земли.

Бесспорно, А. Bjerhammar является одним из основоположников метода истокообразной аппроксимации. М. Hotine в своей вышедшей посмертно монографии [Hotine, 1969, с. 323] так описывает суть предложенного подхода: “Реализующий существенную простоту классического подхода, А. Bjerhammar в одном из самых современных методов использует сферическую поверхность относимости, полностью расположенную внутри Земли. Затем он использует интеграл продолжения поля вверх… Во-первых, чтобы определить значения гравитационной аномалии на поверхности относимости по ее значениям, измеренным на топографической поверхности, и во-вторых, чтобы определить значения внешнего поля по значениям на поверхности относимости. Метод применим без изменений для редуцирования силы тяжести или любой другой производной интеграла Пуассона…”.

Существенный вклад в решение внутренней краевой задачи Дирихле – определение значений потенциальной функции на некоторой внутренней поверхности в теле Земли по ее значениям на поверхности Земли – внесли работы В.Д. Купрадзе и М.А. Алексидзе [Купрадзе, 1967; Алексидзе, 1966; 1987; Купрадзе, Алексидзе, 1964]. В них приводится доказательство возможности существования приближенного решения задачи М.С. Молоденского со сколь угодно малой погрешностью, которое послужило впоследствии отправной точкой для теории метода истокообразной аппроксимации. В 1968 г. М.А. Алексидзе в соавторстве с Г.И. Санадзе представляют версию метода истокообразной аппроксимации с использованием точечных источников, распределенных в некоторой замкнутой области внутри Земли, но не предлагают при этом конкретной вычислительной схемы [Алексидзе, Санадзе, 1968]. Также еще в 1963 г. И.Д. Савинским были получены теоретические результаты, которые также можно считать прообразом современных методов истокообразной аппроксимации [Савинский, 1963].

В своей работе [Молоденский, 1967, с. 287] М.С. Молоденский писал: “Итак, основная геодезическая задача (изучение фигуры и внешнего гравитационного поля Земли) четко отделяется от основной геофизической задачи (изучение внутреннего строения Земли)”. Однако следует заметить, что еще в 1956 г. В.А. Кузивановым было высказано предложение об использование результатов М.С. Молоденского в области разведочной геофизики [Кузиванов, 1956]. Подход, который развивали в своих работах М.С. Молоденский и А. Bjerhammar, а также в общем виде обозначаемый в работе других авторов, был прообразом метода аппроксимации геопотенциального поля потенциальными гармоническими функциями, названный позднее истокообразной аппроксимацией.

В отличие от предшественников, В.И. Аронов разглядел в этом подходе предпосылки для создания принципиально нового метода трансформации геопотенциальных полей, но применительно к решению задач прикладной геофизики, и в частности – гравиразведки. Изучая особенности интерпретации данных гравиметрической съемки в горной местности, связанные с влиянием резко расчлененного рельефа земной поверхности²², В.И. Аронов [Аронов, 1963], опираясь на ранние работы М.С. Молоденского, В.А. Кузиванова, Б.А. Андреева, В.Н. Страхова, предлагает новый способ редуцирования гравитационных аномалий с земной поверхности на горизонтальную плоскость относимости, проходящую через максимальную высотную отметку рельефа. В этой [Аронов, 1963] и последующих работах [Аронов, 1965a; 1965б] приводится не только сам способ, но и реализующая его вычислительная схема, где в качестве базового эквивалентного источника используется участок простого слоя в виде квадратной площадки (трехмерный случай) либо пластинки (двумерный случай). Как пишет сам В.И. Аронов [Аронов, 1970, с. 79]: “Идея метода заключается в том, что на внутренней плоскости $z{\kern 1pt} '$ подбирается такая функция ${{U}_{z}}\left( {\xi ,z{\kern 1pt} '} \right)$, которая при аналитическом продолжении в узлы расчетной сети будет совпадать с заданными значениями ${{U}_{z}}\left( {\xi ,z} \right)$ сколь угодно точно. Дальнейшие расчеты ведутся по обычным формулам, основанным на использовании интеграла Пуассона для решения внешней задачи Дирихле”. И далее [Аронов, 1970, с. 80]: “Фактически решение сводится к аппроксимации исходной функции системой гармонических функций. … С физической точки зрения это решение может быть интерпретировано как подбор на плоскости $z{\kern 1pt} '$ эквивалентного простого слоя с кусочно-непрерывной плотностью (постоянной на элементарном отрезке интегрирования $\Delta x$), который создает на физической поверхности наблюденные аномалии³³”.

В работе В.И. Аронова совместно с В.М. Гординым [Аронов, Гордин, 1971, с. 1] прослеживается эта же логика: “Основная идея метода заключается в замене реального распределения возмущающих масс простым слоем с кусочно-постоянной плотностью, расположенным на внутренней плоскости $z{\kern 1pt} '$ и создающим в гравиметрических пунктах наблюдаемые аномалии $\Delta g$. Исходное поле при этом фактически интерполируется системой гармонических функций – аномалий от элементарных квадратов с площадями $h \times h$, залегающих на глубине $z{\kern 1pt} '$”.

Казалось бы, метод создан, имеется достаточно простая вычислительная схема с использованием сеточной модели, то есть состоялся переход от интегрального представления к сумме и решению системы линейных алгебраических уравнений (СЛАУ). Что еще нужно? Чего не хватает? Как оказалось, не хватает простоты. Аномальные гравитационные эффекты квадратных пластинок – тех самых “простых” потенциальных гармонических функций, суммой которых предполагалось аппроксимировать аномальное поле – имеют довольно сложное аналитическое выражение, содержащее трансцендентную функцию arctg. Легко представить, насколько трудоемким было построен аппроксимации с многократным использованием подпрограммы расчета функции arctg даже для небольшого числа точек задания поля, учитывая, что происходило это в эпоху ранних вычислительных машин, не обладающих большой производительностью. Возможно, что работа [Алексидзе, Санадзе, 1968] прошла для остальных исследователей незамеченной или В.И. Аронов с коллегами считали оправданным применение в качестве эквивалентного источника именно пластинки для решения своих задач.

В книге [Аронов, 1977] В.И. Аронов уже рекомендует использовать при аппроксимации более простые гармонические функции. Сложно отследить тот момент, когда он обратил внимание на аномальный эффект точечной массы как на самый практичный для использования в роли истокообразной функции. Еще в работе [Гордин и др., 1980] ее авторы в качестве эквивалентного источника использовали притяжение горизонтальной вещественной линии. Однако в наиболее полной монографии В.И. Аронова [Аронов, 1990], посвященной, в том числе, истокообразной аппроксимации потенциальных полей, в качестве аппроксимирующей функции уже фигурирует поле точечного источника.

Точечные источники также начали использоваться в геодезической гравиметрии. Дж. Бальмино представил аномальную часть гравитационного потенциала шарообразной Земли потенциалом притяжения 126 точечных масс, расположенных на глубинах от 1000 до 2000 км [Бальмино, 1975]. Пространственное расположение этих масс согласовано с наиболее интенсивными аномалиями силы тяжести, зафиксированными на земной поверхности.

По всей видимости, из-за “железного занавеса”, препятствующего регулярному обмену знаниями и опытом между зарубежными и советскими учеными, незамеченными для российских геофизиков того времени прошли работы [Dampney, 1966a; 1966b; 1969], австралийского геофизика C.N.G. Dampney, и в особенности его магистерская диссертация [Dampney, 1966b], представленная к защите в ноябре 1966 г. и защищенная им в 1967 г. в Университете Тасмании. В ней он впервые предлагает использовать сеточное распределение точечных источников в качестве эквивалентной модели для аппроксимации потенциального поля и приводит полное теоретическое обоснование метода, а кроме того – рабочий алгоритм, программу для ЭВМ “Эллиот-503” и пример ее практического применения для продолжения гравитационного поля в верхнее полупространство. Позднее, в 1969 г., выходит его наиболее известная работа по этой теме – “The equivalent source technique” [Dampney, 1969]. Именно на нее, как на первоисточник, ссылаются все зарубежные геофизики, имеющие отношение к развитию или применению метода истокообразной аппроксимации.

Нет оснований считать, что зарубежные геофизики не были знакомы с трудами советских коллег. Работы C.N.G. Dampney изобилуют ссылками на ранние работы В.А. Романюка [Romanyuk, 1959], В.Н. Страхова [Strakhov, 1961], в том числе в соавторстве с М.И. Лапиной [Strakhov, Lapina, 1963], А.В. Тархова и А.А. Сидорова [Tarkhov, Sidorov, 1960], но, к сожалению, не содержат ни одного упоминания работ В.И. Аронова. C.N.G. Dampney занимался этой тематикой в то же время, что А. Bjerhammar и В.И. Аронов, но развивал свою линию исследования, опираясь в основном на работы англоязычных авторов, среди которых, в частности – E.C. Bullard, R.I.B. Cooper [Bullard, Cooper, 1948], R.A. Smith [Smith, 1960] и др. В ранней работе [Dampney, 1966a] он, как и В.И. Аронов, еще использует в качестве базисной функции при аппроксимации поле горизонтальной пластинки. Однако именно ему принадлежит первенство в переходе к более простой в вычислительном плане сеточной модели точечных эквивалентных источников, которая по сей день является самой популярной у геофизиков при решении прикладных задач [Dampney, 1966b].

И все же нельзя сказать, что идея истокообразной аппроксимации применительно к трансформациям аномалий геопотенциальных полей возникла только в 60-х годах ХХ века. Еще в 1956 г., опередив В.И. Аронова и C.N.G. Dampney на целое десятилетие, А.К. Маловичко в своей книге [Маловичко, 1956, с. 67] отмечал: “Любое аномальное поле можно интерпретировать одним или несколькими аномальными телами с центрами на соответствующих глубинах. Если такая интерпретация лишена геологического смысла, то результаты ее с большим эффектом могут быть использованы в качестве подсобного решения. При помощи найденных выпуклых тел можно вычислить значения вертикальной составляющей притяжения в точках вспомогательной плоскости и тем самым найти аналитически продолженные значения аномалий. Точность такого решения зависит всецело от точности аппроксимации наблюденного аномального поля действием подобранных выпуклых тел. … Из обширного списка тел, для которых решена обратная задача и которые могут быть использованы в качестве подсобных решений, в практическом отношении оправданным будет применение только самых простых тел, вроде шара, горизонтальной полосы и пр. … Разумеется, наибольшего внимания заслуживает шар”.

А.К. Маловичко не представил простой вычислительной схемы для подбора параметров подсобных тел при построении аппроксимации, но его идея была своевременной. Ее возникновение было предопределено развитием аппаратурной базы и появлением первых гравиметров, заменивших вариометры и градиентометры и существенно увеличивших производительность и объемы выполнения гравиметрических съемок. Гравиметрия как метод перестала ассоциироваться только с решением геодезической задачи об определении фигуры Земли и постепенно эволюционировала в гравиметрический метод, применяющийся для геокартирования, поисков и разведки полезных ископаемых – в гравиразведку [Сорокин, 1953]. Очевидно, что идея использования подсобных тел для продолжения аномалий силы тяжести в верхнее полупространство в то время “витала в воздухе” и материализовалась не только в работах А.К. Маловичко, но и других отечественных и зарубежных геофизиков. К сожалению, многие из этих работ были вскоре забыты.

Что касается C.N.G. Dampney, то отвечая на вопросы, которые задавал A. Roy в дисскуссионном разделе журнала “Geophysics” [Roy, 1970], он на первый план он поставил работу [Zidarov, 1965] в качестве основной предпосылки для разработки его “The equivalent source technique”. Действительно, именно в этой работе D. Zidarov показал принципиальную возможность подбора физических параметров совокупности фиксированных в пространстве точечных масс, суммарное поле которых отвечает наблюденному полю реальных геологических тел. Сама процедура подбора, представляющая собой по сути решение обратной задачи гравиметрии, сводится к решению СЛАУ. Причем, D. Zidarov распространил задачу на случай магнитного и электрического полей, предложив для этой цели соответствующие типы элементарных источников.

Если попытаться разобраться в том, кому принадлежит авторство отдельных аспектов теории, то, на наш взгляд, в некоторых вопросах C.N.G. Dampney опередил В.И. Аронова, например, раньше исследовав и опубликовав теоретические результаты, касающиеся оптимальной глубины расположения эквивалентных источников, связав эту глубину с обусловленностью решаемой СЛАУ. Однако позднее В.И. Аронов не только дал развернутое изложение теоретической составляющей метода истокообразной аппроксимации, но и подробно осветил многие прикладные аспекты его применения [Аронов, 1990].

В работе А.К. Маловичко [Маловичко, 1956, с. 139] выдвигается другая, не менее революционная идея, относящаяся к методу подсобных тел: “При решении задачи практически используется аномальное поле только в радиусе пяти–шести значений ${{z}_{0}}$ вокруг данной точки⁴⁴. Это, очевидно, и определяет размеры области, в точках которой измеренные значения поля должны быть аппроксимированы полем подобранного тела”. Идея А.К. Маловичко впоследствии будет отмечена другими геофизиками и с развитием компьютеров приведет к существенному сокращению вычислительных затрат и ускорению решения СЛАУ путем разрежения матрицы коэффициентов (значений истокообразной функции) за счет выделения ленточной части. Такой подход использован в диссертационной работе [Пугин, 2007] при решении задачи истокообразной аппроксимации большой размерности.

РАЗВИТИЕ МЕТОДА ИСТОКООБРАЗНЫХ АППРОКСИМАЦИЙ ГЕОПОТЕНЦИАЛЬНЫХ ПОЛЕЙ. АНАЛИЗ РАБОТ СОВЕТСКИХ (РОССИЙСКИХ) И ЗАРУБЕЖНЫХ ИССЛЕДОВАТЕЛЕЙ

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

Применительно к истории развития истокообразных аппроксимаций вторую половину ХХ века условно можно разделить на следующие временные периоды:

1) ориентировочно с 1960 по 1980 гг. – докомпьютерная эпоха или эпоха ручного счета, которая характеризуется отсутствием персональных компьютеров и необходимостью обрабатывать геофизические данные вручную. Это процесс лишь частично был автоматизирован с применением электронных калькуляторов или с использованием первых серийных ЭВМ, имеющих низкую производительность. В этот период в разумные сроки (от нескольких часов до суток) могли быть решены лишь задачи небольшой размерности и вычислительной сложности (массивы данных, содержащие от нескольких десятков до нескольких сотен значений поля);

2) с 1980 по 2000 гг. – ранняя компьютерная эпоха, время широкого внедрения первых персональных компьютеров, когда появляются возможности для автоматизированной обработки в те же сроки сравнительно больших массивов данных и задач средней размерности (массивы до нескольких тысяч или даже первых десятков тысяч значений поля);

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

ПЕРВЫЙ ПЕРИОД – ДОКОМПЬЮТЕРНАЯ ЭПОХА ИЛИ ЭПОХА “РУЧНОГО СЧЕТА”

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

Одним из первых практических приложений метода истокообразной аппроксимации, помимо аналитического продолжения потенциальных полей в верхнее полупространство, стало решение задачи интерполяции данных, без которой не удавалось решить другую не менее важную задачу гравиметрии – вычисление гравиметрических уклонений отвеса. Именно этим задачам посвящена одна из первых совместных публикаций В.И. Аронова и С.М. Гордина [Аронов, Гордин, 1971]. Нужно отметить, что метод обеспечил погрешность вычислений составляющих уклонения отвесной линии ∼1" в центральной и ближней зоне [Гравиметрия и геодезия, 2010]. Позднее это приведет к целому циклу статей, посвященных трехмерной (учитывающей рельеф поверхности наблюдений) интерполяции методом истокообразной аппроксимации, в основу которой положено решение СЛАУ. Результат интерполяции представляет собой гармоническую функцию, а сама интерполяция осуществляет одновременно фильтрацию негармонической помехи.

Предложенное В.А. Кузивановым теоретическое решение задачи редуцирования гравиметрических измерений на горизонтальную плоскость было программно реализовано в 70-х годах ХХ века группой исследователей под руководством О.К. Литвиненко [Литвиненко и др., 1973]. В процессе вычислений итерационно определялась плотность простого слоя, расположенного на поверхности Земли. По приведенным оценкам погрешность метода обычно составляет 1–1.5%, но при углах наклона поверхности измерений более 36° возрастает до 5%.

В 1976 г. В.И. Ароновым и В.О. Михайловым “при реализации метода интерполяции системой гармонических функций” окончательное решение рассматривалось “как сумма решений на ряде последовательных шагов, при каждом из которых решение является устойчивым” [Аронов, 1990, с. 88]. Это достигается путем размещения аппроксимирующих масс на ряде поверхностей, глубины до которых от поверхности измерений поля соизмеримы со средними расстояниями между точками в соответствующих выборках значений поля. Максимальное значение глубины в этом случае соответствует радиусу наибольшего “белого пятна” (лакуны) в области задания поля. Это действительно был революционный шаг, поскольку впервые было предложено адаптировать сеточное распределение эквивалентных источников к структуре исходных данных, что делало метод еще более гибким для решения прикладных задач.

Под обеспечением устойчивости интерполяции понималось, прежде всего, решение проблемы минимизации количества ложных высокочастотных аномалий (артефактов), возникающих между точками наблюдений. Данное явление также было описано C.N.G. Dampney и в зарубежной литературе называется “алиасингом” (aliasing – в ориг.). В отечественной геофизике термин “алиасинг” употребляется преимущественно сейсморазведчиками.

За рубежом в эти же годы A. Bjerhammar продолжает развивать свой метод продолжения поля на внутреннюю сферу применительно к задаче определения фигуры Земли и вводит дискретное распределение потенциальной гармонической функции $\Delta g$ на ней с использованием функции Дирака, по сути, переходя от континуального к сеточному распределению [Bjerhammar, 1975]. Этот переход был оправдан тем, что существенно упрощал практические расчеты.

В 1972 г. M.H.P. Bott, A. Ingles представляют метод аппроксимации аномального магнитного поля, названный ими матричным, использующий горизонтальный эквивалентный слой из плотной упаковки прямоугольных призм [Bott, Ingles, 1972]. Авторы показывают возможность псевдогравитационного преобразования с помощью введенной ими конструкции из аппроксимирующих призм.

В 1973 г., уже со ссылкой на работу C.N.G. Dampney [Dampney, 1969], появляется работа D. Emilia [Emilia, 1973], знаковая по той причине, что впервые подход с элементарными эквивалентными источниками был распространен на случай магнитного поля. В качестве эквивалентных источников применялись магнитные диполи.

На теоретическом примере D. Emilia показал, что хотя полный вектор аномального магнитного поля $\Delta {\mathbf{T}}$, строго говоря, не является потенциальной гармонической функцией, но учет магнитного наклонения и склонения при аппроксимации эквивалентными источниками позволяет добиться вполне приемлемого совпадения модельного поля и его производных с соответствующими функциями, полученными путем решения прямой задачи магниторазведки. Причем сопоставление проводилось в том числе и при наличии (D. Emilia пишет – при доминировании) остаточной намагниченности у аномалиеобразующих тел, отличающейся по направлению вектора от современного нормального магнитного поля Земли. Этот пример показал, что нет никаких принципиальных препятствий для использования техники эквивалентных источников в прикладных задачах магниторазведки, если погрешности вычислений невелики, а результат может быть корректно интерпретирован с геологических позиций.

Используя подходы C.N.G. Dampney и D. Emilia в 1977 г. B.K. Bhattacharyya, K.C. Chan представили свою модификацию аппроксимационного метода в попытках подавить влияние расчлененного рельефа, сложенного неоднородно намагниченными породами [Bhattacharyya, Chan, 1977]. Эквивалентные источники располагались непосредственно на поверхности наблюдений. Рассмотрен также случай, когда поверхность наблюдений – горизонтальная плоскость, находящаяся выше максимальной высоты рельефа местности, но источники расположены на земной поверхности. Авторы показывают, что редуцирование поля на горизонтальную плоскость выше рельефа местности с использованием их аппроксимации существенно подавляет эффект разновысотности, а если редуцирование производить на поверхность в верхнем полупространстве, конкордантную рельефу (полученную увеличением отметок высот на некоторую константу), то частично подавляется и влияние самого магнитного рельефа. К сожалению, в статье не раскрывается тема алиасинга, который неизбежно возникает при расположении эквивалентных источников, независимо от их типа, в непосредственной близости от точек задания поля.

ВТОРОЙ ПЕРИОД – РАННЯЯ КОМПЬЮТЕРНАЯ ЭПОХА

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

В работе В.И. Старостенко рассматривается “построение числовой модели гравитационных полей”, которое подразумевает “отнесение поля к узлам равномерной сети, приведение его к плоскости, расположенной на определенной высоте, выполнение некоторых трансформаций поля” [Старостенко, 1978]. Предлагается устойчивый алгоритм построения аппроксимационной конструкции, представленной точечными массами, использующий метод регуляризации А.Н. Тихонова. В процессе аппроксимации используется скользящее окно, параметры которого определяются экспериментально, затем для каждого местоположения окна решается нормальная СЛАУ с числом неизвестных, не превышающем 121. В работе [Старостенко, 1978] приводятся 5 примеров, свидетельствующих о высокой эффективности предложенного алгоритма, в т.ч. в плане подавления интенсивных аномалий-помех, осложняющих исходные данные.

Неоценимый вклад в развитие истокообразных аппроксимаций внес В.М. Гордин, начинавший научную работу под руководством В.И. Аронова и являвшийся соавтором многих его работ. Как пишет сам В.М. Гордин в своей книге “Избранные труды и воспоминания”, изданной лишь спустя год после его смерти: “Наивысшим достижением этих исследований⁵⁵ явился гибрид идеи истокообразных аппроксимаций с оптимальной фильтрацией по критерию Колмогорова–Винера. Полученные результаты были “доведены до числа” и легли в основу статьи “Физические аспекты фильтрации и разделения геофизических аномалий”, написанной в соавторстве с Б.О. Михайловым и В.О. Михайловым [Гордин и др., 1980]”.

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

Еще в 70-х годах ХХ века была обоснована возможность решения задачи разделения потенциальных полей, обусловленных разными источниками, на основе выметания масс с применением сеточных аппроксимирующих конструкций. Задача сводится к решению интегрального уравнения Фредгольма 1-го рода, реализованному в виде итерационного алгоритма с регуляризацией по А.Н. Тихонову [Непомнящих и др., 1978]. Эквивалентное перераспределение физического параметра осуществлялось по шеститочечной схеме Д. Зидарова. Было создано программное обеспечение, апробированное на модельных и практических примерах, разработаны методические рекомендаций по его применению. В диссертации Л.В. Ли показано, что построенные таким образом результативные сеточные модели также вполне удовлетворительно характеризуют морфологию, условия залегания и размеры реальных источников [Ли, 1984].

Тем временем в параллельной “зарубежной реальности” в 1981 г. T. Nakatsuka [Nakatsuka, 1981] исправляет существенную ошибку, допущенную авторами в работе [Bhattacharyya, Chan, 1977]. При выводе финального уравнения для аппроксимации магнитного поля B.K. Bhattacharyya и K.C. Chan использовали допущение о постоянстве направления вектора магнитных моментов аппроксимирующих эквивалентных диполей, располагаемых на произвольной поверхности S, полученной параллельным переносом вниз поверхности наблюдений, заявляя в то же время, что диполи ориентированы по нормали к поверхности S. Очевидно, что с изменением угла наклона этой поверхности изменяются направления магнитных моментов диполей. T. Nakatsuka справедливо отметил, что в таком случае направление вектора магнитных моментов эквивалентных источников будет переменным, и будет весьма сильно изменяться в случае значительного перепада высот поверхности S. С его слов, введя допущение о постоянстве направления вектора магнитных моментов, авторы эффективно ввели некую аппроксимацию “слукавив” (не говоря об этом прямо), что поверхность S, оказывается, не так уж сильно криволинейна и слабо отличается от плоскости. Вводя другое допущение о постоянном вертикальном направлении данного вектора, T. Nakatsuka исправляет их оплошность и дает более точное и быстрое в плане сходимости решение задачи истокообразной аппроксимации, однако тема алиасинга здесь также не раскрывается.

В 1981 г. von R.R.B. Frese с коллегами [von Frese et al., 1981], опираясь на работы C.N.G. Dampney и D. Emilia, адаптируют истокообразные аппроксимации на случай задания результатов измерения гравитационного и магнитного полей в сферических координатах и затем применяют метод для аппроксимации региональных аномалий гравитационного поля на земной поверхности и спутниковых данных магнитного поля. Аппроксимация здесь выполняется с целью последующих аналитических вычислений производных гравитационного и магнитного потенциала в верхнем полупространстве и приведения их к одной поверхности. Эта работа была написана на основе научного отчета [von Frese et al., 1980], выполненного для агентства NASA и послужившего основой докторской диссертации von R.R.B. Frese. Обсуждая проблемы вычислительного характера, авторы, со ссылкой на свою более раннюю работу, рассматривают возможность размещения эквивалентных источников в модели в соответствии с особенностями аппроксимируемых аномалий – конструктивная идея, создающая предпосылки для развития алгоритмов с предварительным анализом структуры аномальных полей и последующей адаптацией к ней пространственного распределения эквивалентных источников. Впоследствии уже в 1998 г. von R.R.B. Frese выпускает статью с корректировками [von Frese, 1998] к работе [von Frese et al., 1981], исправляя ошибки, допущенные при выводе некоторых уравнений.

В 1984 г. R.O. Hansen и Y. Miazaki [Hansen, Miazaki, 1984] анонсируют возможности истокообразной аппроксимации в задаче пересчета потенциальных полей на произвольную поверхность в верхнем полупространстве. Они отмечают слабое место подхода предшественников [Bhattacharyya, Chan, 1977; Nakatsuka, 1981] – расположение эквивалентных источников непосредственно на поверхности наблюдений. R.O. Hansen, Y. Miazaki утверждают, что такая аппроксимация весьма ограничена в плане пересчета поля вверх, и предлагают располагать источники на поверхности, конкордантной дневному рельефу, но лежащей всюду ниже него. Вместо “дискретных” диполей в качестве эквивалентных источников предложено использовать поле однородно намагниченных горизонтальных пластин. Авторы приводят интересный практический пример редукции данных аэромагнитной съемки, выполненной с обтеканием рельефа местности, на горизонтальную плоскость, расположенную на средней высоте съемочных полетов.

J.B.C. Silva в 1986 г. предлагает использовать сетку эквивалентных источников (магнитных диполей) для приведения аномального магнитного поля к полюсу, показывая, что такой подход существенно упрощает решение этой задачи и делает ее устойчивым [Silva, 1986].

В 1990 г. выходит наиболее полная монография В.И. Аронова [Аронов, 1990], знаковая для отечественных геофизиков, обобщающая все теоретические, методические и практические аспекты метода истокообразной аппроксимации, разработанные к этому времени, и предвосхищающая многие вопросы, которые будут позднее поставлены и изучены. Несмотря на то, что технология эквивалентных источников уже более 20 лет применялась на практике, в зарубежной и отечественной геофизике настолько всеобъемлющих публикаций по методу истокообразных аппроксимаций не было ни до 1990 г., ни после. Следует отметить, что в этой монографии рассматривается весьма широкая область практического применения метода – кроме работы с потенциальными геофизическими полями, рассмотрены задачи интерполяции структурных поверхностей, геометризации залежей углеводородов и подсчета запасов, оптимизация размещения скважин при разведочных работах на нефть и газ.

Анализ содержания двух следующих работ показывает, насколько внушительным по времени было отставание исследований зарубежных геофизиков от группы В.И. Аронова. В 1991 г. J. Xia и D.R. Sprowl [Xia, Sprowl, 1991] вновь обращаются к проблеме пересчета гравитационного поля с земной поверхности на горизонтальную плоскость с целью подавления искажающего влияния разновысотности точек измерений. При этом они опираются на технологию C.N.G. Dampney [Dampney, 1969] и, все также, располагают эквивалентные источники в одном слое, обосновывая выбор оптимальной глубины их расположения с позиций минимизации эффекта алиасинга.

Годом позже L. Cordell [Cordell, 1992], разрабатывая алгоритм интерполяции измеренных потенциальных полей в узлы регулярной сети и опираясь при этом на работу J. Xia и D.R. Sprowl, предлагает распределять источники, адаптируя их положение к амплитуде поля и его морфологии. Алгоритм является итерационным и весьма оригинален: предварительно задается некий коэффициент пропорциональности, связывающий глубину погружения источника с расстоянием до ближайшей точки измерения поля. Определяется максимальная амплитуда поля в пределах площади исследований, затем – глубина до эквивалентного источника под этой точкой, исходя из этого коэффициента. Физический параметр (например, масса) источника вычисляется, следуя предположению, что аномалия в точке над ним обусловлена исключительно данным элементарным источником. Далее гравитационный эффект источника вычитается из измеренного поля во всех точках и определяется новое максимальное значение поля.

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

В 1992 г. В.Н. Страхов публикует развернутую статью с исследованием алгоритмов трансформаций геопотенциальных полей, основанных на истокообразной аппроксимации. В статье рассмотрена проблема оптимизации параметров аппроксимационной конструкции для более точного приближения исходного поля [Страхов, 1992].

В этом же году, работая с данными морской магнитной съемки, J.B.C. Silva в соавторстве с С.A. Mendonça реализуют весьма экономичный в плане вычислительных затрат алгоритм истокообразной аппроксимации с однослойным расположением эквивалентных источников, используя, подобно тому, как это было в подходе L. Cordell, технологию последовательных приближений. На первой итерации находятся минимальное и максимальное значения поля и вносятся в выборку “определяющих” значений, тогда как остальные значения поля принимаются избыточными. Эквивалентные источники располагаются только под точками экстремумов поля и по этим “определяющим” значениям находятся их физические параметры. Далее во всех остальных точках вычисляются модули отклонений наблюденного и модельного полей, и под точкой измерения, где имеет место наибольшее значение модуля вновь располагается эквивалентный источник. Значение измеренного поля в этой точке добавляется в выборку “определяющих” и участвует теперь в нахождении физических параметров всех источников. В практическом примере, содержащем 3137 наблюдений морской гравиметрической съемки на шельфе северной Бразилии, размер “определяющей” выборки составил 294 значения. Алгоритм продемонстрировал высокую экономичность (по оценкам авторов – ускорение до двух порядков) и в плане вычислений, и с точки зрения малого числа эквивалентных источников, необходимых для аппроксимации поля. Авторы назвали его “determination of equivalent data (DED)” – определение эквивалентных данных [Mendonça, Silva, 1992; 1994] и использовали для решения задачи интерполяции потенциальных полей, представив годом позднее исчерпывающий сравнительный анализ своего алгоритма и метода минимальной кривизны [Mendonça, Silva, 1995].

В этот период времени под руководством И.Д. Савинского была создана компьютерная технология “Альфа-4”, все алгоритмы которой основаны на едином математическом аппарате решения интегрального уравнения Фредгольма 1-го рода, сводящимся, как и у многих других исследователей, к решению СЛАУ. В некоторых случаях применяется регуляризация. С физической точки зрения, это аппроксимация поля совокупностью элементарных источников, суммарный эффект которых автор называет представлением поля в виде полинома. Система предусматривает выполнение псевдогравитационного преобразования магнитного поля, пересчет полей на нижележащую наклонную плоскость, вычисление распределения плотности гравитирующих тел и градиента этого распределения, интерпретацию трехкомпонентной скважинной магниторазведки, оценку ресурсов полезных ископаемых путем определения интегральных характеристик источников аномалий [Савинский, 1995]. Приведены примеры успешного использования технологии “Альфа-4” на Березовском железорудном месторождении, месторождении хромитов на Полярном Урале, медно-никелевом месторождении Норильск-2, при обработке практических материалов среднемасштабных гравиметрических и аэромагнитных съемок [Савинский, 1994]. Однако примеров использования этой технологии другими исследователями авторам найти не удалось.

В 1995 г. Т. Nakatsuka предлагает интересную технику для подавления краевых эффектов при трансформациях магнитного поля на основе истокообразной аппроксимации [Nakatsuka, 1995]. Применяя в качестве эквивалентных источников прямоугольные призмы, он предлагает разместить часть призм за пределами области задания поля. Параметры источников определяются при этом путем решения недоопределенной СЛАУ, которая возникает за счет превышения числа источников над числом точек.

В 1996 г. выходит обширная монография R. Blakely [Blakely, 1996], в которой в числе прочего, со ссылкой на многие перечисленные выше работы, кратко излагаются основы аппроксимации с применением эквивалентных источников. Эта работа носит обзорный характер и не привнесла ничего нового в метод истокообразной аппроксимации, но часто используется вместо ссылки на исследования C.N.G. Dampney в публикациях зарубежных геофизиков, как обобщающий источник.

В 1997 г. выходит обстоятельная и разносторонняя по содержанию статья Г.Я. Голиздры [Голиздра, 1997], в которой приводится наиболее полный обзор российских и зарубежных работ, касающихся аналитического продолжения потенциальных полей, опубликованных до этого времени. В ней приводится большой список работ В.Н. Страхова по теории методов аналитического продолжения потенциальных полей, переносить который целиком в данную статью мы не сочли возможным. Интересно, что в одной из работ Владимир Николаевич называет продолжение поля в верхнее полупространство с использованием аналитических аппроксимаций не аналитическим, а аппроксимационным, имея в виду его приближенный характер.

Следует отметить, что с 1997 г. В.Н. Страховым и его учениками развивается метод линейных интегральных представлений [Страхов и др, 1997; 1999], являющийся обобщением и развитием классического метода линейных интегральных уравнений. Этот метод обеспечивает решение многих задач гравиметрии и магнитометрии на общей аппроксимационной основе, в т.ч. нахождение линейных аналитических аппроксимаций элементов аномальных геопотенциальных полей: S- и F-аппроксимации. Теоретические работы В.Н. Страхова определили стратегию развития теории интерпретации геофизических полей на многие десятилетия вперед, их перечень можно найти в автобиографической монографии [Страхов, 2008]. К сожалению, объем статьи вынуждает авторов избирательно подойти к выбору ссылок на работы этого выдающегося ученого.

В 1999 г. в сборнике материалов первой всероссийской конференции “Геофизика и математика” выходит статья В.М. Гордина и С.А. Тихоцкого по истории развития метода истокообразной аппроксимации [Гордин, Тихоцкий, 1999].

В конце 90-х годов с появлением в России достаточно мощных персональных компьютеров метод истокообразной аппроксимации начинает широко использоваться для решения прикладных задач в области гравиразведки, магниторазведки [Долгаль, Христенко, 1997; Долгаль, 1999a], а затем и электроразведки методом естественного электрического поля (ЕП) [Долгаль, 1999б]. Причем в статье [Долгаль, 1999a] впервые предлагается оригинальный подход к решению проблемы краевых эффектов и боковых источников при трансформациях геопотенциальных полей с использованием истокообразных аппроксимаций. Данные эффекты возникают в силу ограниченных размеров площади задания поля. Обычно их пытаются подавить “гипотетическими” процедурами, основанными не на фактических измерениях поля за пределами площади, а на каких-либо допущениях о поведении поля и его источниках. Здесь же предлагалось расширить пределы площади за счет привлечения фактических данных мелкомасштабных гравиметрических и магнитных съемок, выполненных ранее и сохраненных в банках данных специализированных предприятий.

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

ТРЕТИЙ ПЕРИОД – ЭПОХА КОМПЬЮТЕРНЫХ ТЕХНОЛОГИЙ

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

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

В начале 2000-х годов иностранные исследователи все еще пытаются решить задачу интерполяции геопотенциального поля в узлы регулярной сети, используя однослойную конструкцию эквивалентных источников, но изменяя глубину их размещения в зависимости от расстояния между точками задания данных [Cooper, 2000], а также пытаясь уйти от процедур обращения матриц при определении параметров источников к более экономичным методам [Cooper, 2002].

В английском языке регулярная, т.е. квадратная или прямоугольная сеть точек задания параметра называется “грид” (grid – в ориг.), а процедура перевода данных из нерегулярной сети в грид – “гриддинг” (gridding – в ориг.). Этот термин уже прочно заимствован и достаточно часто используется в российском геофизическом диалекте научного языка.

В 2001 г. в работе [Li, 2001] используется подход C.N.G. Dampney для аппроксимации данных гравитационной градиентометрической съемки, но подобно Т. Nakatsuka с применением в качестве эквивалентного слоя плотной упаковки кубических тел.

В 2000 г. В.И. Аронов в статье [Аронов, 2000] рассматривает проблему трехмерной аппроксимации с целью интерполяции с использованием координат x, y, z исходных данных (не только значений потенциальных полей). Представлены примеры параметрической интерполяции – восстановления гладкой границы произвольного трехмерного тела по координатам вершин аппроксимирующего многогранника, заданного на его поверхности. Приведена структурная карта по подошве меловых отложений юго-западной части Прикаспийской впадины и карта концентрации азота на данном стратиграфическом уровне, построенная путем 3D-интерполяции. Выявлены достаточно заметные различия этих карт. Вполне можно согласиться с выводами автора о том, что “переход к решению трехмерных задач обработки геолого-геофизических данных открывает новые перспективы при подсчете запасов и изучении внутренней структуры месторождений, определении эксплуатационных параметров”. Не менее важно замечание автора о том, что отсутствие эффекта алиасинга в данном случае достигается не только за счет регуляризации, сколько за счет специфических особенностей самого алгоритма, позволяющих ему адаптироваться к структуре сети исходных данных и отношению сигнал/помеха.

В 2000 г. выходит большая обзорная работа выдающихся ученых [Страхов и др., 2000] о достижениях в области интерпретации геопотенциальных полей в ХХ веке, в которой в качестве инструмента для устойчивого и точного нахождения трансформаций предлагается метод, истории развития которого посвящена настоящая работа. Отмечается что истокообразная аппроксимация является универсальным подходом, пригодным для численной реализации как корректных, так и некорректных (например – вычисление производных) трансформаций.

В это же время А.С. Долгаль выпускает серию работ прикладного характера, в которых истокообразные аппроксимации успешно применяются для обработки потенциальных полей при решении рудных задач в горной местности [Долгаль, 2000; 2001а; 2002б; Долгаль и др., 2002]. В них приводятся примеры комплексной интерпретации наземной гравиметрической и аэромагнитной съемок, построения аналитической модели рельефа земной поверхности, а также развивается подход В.И. Аронова с многоуровневым размещением эквивалентных источников при интерполяции измеренных значений геопотенциальных полей. Местоположение источников в каждом слое в проекции на горизонтальную плоскость определяется путем анализа поля в пределах скользящего окна заданного размера, а их физические параметры – исходя из решения соответствующей СЛАУ, в которой правая часть определена из скользящего среднего (среднеарифметического значения поля в окне).

В 2002 г. выходит монография А.С. Долгаля [Долгаль, 2002б], отражающая основное содержание докторской диссертации [Долгаль, 2002c]. В ней описан оригинальный подход, который более чем в 1.5 раза ускоряет итерационный процесс вычисления физических параметров эквивалентных источников в слое при квадратной сети задания поля. Предлагается модифицированный метод Зейделя с предварительной профильной фильтрацией поля по строкам грида на основе той же истокообразной аппроксимации и закреплением некоторого числа “каркасных точек” при решении СЛАУ. По сути, осуществляется декомпозиция задачи на ряд подзадач меньшей размерности. Также предложен подход к комплексной интерпретации данных наземной гравиразведки и аэромагнитной съемки, приведенных к единой горизонтальной плоскости при поисках платино-медно-никелевого оруденения в Норильском районе.

В период с 2000 по 2011 гг. выходит серия работ Е.Г. Булаха и его учеников [Булах, Шиншин, 2000; 2002; Булах и др., 2011], в которых предлагается использовать для истокообразной аппроксимации источник, собранный из трех попарно ортогональных стержней, ориентированных вдоль осей декартовой системы координат и образующих симметричную структуру относительно точки их пересечения. В этом случае для построения аналитической аппроксимации требуется также определить вектор физических параметров стержней. Также рассматриваются аппроксимации в других классах источников, например, для тел Л.Н. Сретенского, обладающих средней плоскостью [Bulakh, Markova, 2008].

В этот же период выходят две статьи: в одной из них, обзорной, аналитические аппроксимации рассматриваются для различных классов тел с точки зрения построения вспомогательных конструкций для последующих трансформаций и с позиций построения геологически содержательных решений [Bulakh, 2007]; в другой предлагается использовать эквивалентные источники для аппроксимации рельефа земной поверхности с целью интерполяции высот [Булах, Маркова, 2009]. Все перечисленные выше задачи были включены в последнюю монографию Е.Г. Булаха, вышедшую спустя год после смерти автора [Булах, 2010].

В 2001 г. на 28-й сессии Международного семинара им. Д.Г. Успенского С.А. Тихоцким и Д.Ю. Шур, в развитие подхода В.И Аронова, представлен алгоритм построения многоуровневых истокообразных аппроксимаций для решения задач магнитной картографии и трансформаций магнитного поля [Тихоцкий, Шур, 2001].

В 2002 г. В.Н. Страхов и И.Э. Степанова публикуют две ключевые работы [Страхов, Степанова, 2002a; 2002б] по методу S-аппроксимаций, представляющему собой частный случай метода линейных интегральных представлений, в котором потенциальное поле предлагается приближать суммой полей простого и двойного слоя. Метод был предложен в локальном и региональном вариантах и впоследствии вошел в диссертационную работу И.Э. Степановой [Степанова, 2003]. Чуть позднее выходит ее работа [Степанова, 2004], в которой приводятся результаты опробования S-аппроксимаций на материалах гравиметрической и магнитной съемок.

В 2004 г. на очередной сессии Международного семинара им. Д.Г. Успенского В.М. Гордин, С.А. Тихоцкий и Д.Ю. Шур делают доклад о возможности и адекватности восстановления гармонической компоненты скалярных аномалий магнитного поля по аномалии модуля напряженности магнитного поля с применением истокообразных аппроксимаций [Гордин и др., 2004], а в 2006 г. выходит расширенная статья и ее англоязычная версия [Gordin et al., 2006]. Условия сходимости итерационного процесса и его устойчивость исследованы в широком диапазоне амплитуд аномального поля на численных примерах. Установлена эффективность предложенного алгоритма при решении часто встречающихся практических задач.

В 2005 г. коллектив авторов, в составе которого в том числе С.А. Тихоцкий и Д.Ю. Шур, публикуют интересную работу [Соловьев и др., 2005] с алгоритмом определения вектора магнитного момента аномалообразующих тел на основе кластерного анализа множества эквивалентных источников, адаптированного к распределению особых точек магнитного поля.

В это же время А.С. Долгалем разработан алгоритм, реализующий интегральный метод решения обратной задачи гравиметрии на основе сеточных распределений эквивалентных источников [Долгаль, 2005].

В 2006 г. японские геофизики T. Nakatsuka и Sh. Okuma используют технику эквивалентных источников для аппроксимации результатов вертолетной аэромагнитной съемки с целью гриддинга и последующих трансформаций поля [Nakatsuka, Okuma, 2006]. В ней со ссылкой на работу 1995 г. [Nakatsuka, 1995] повторно уделяется внимание методике расширения эквивалентного слоя путем выноса дополнительных источников за пределы площади задания поля с последующим решением недоопределенной СЛАУ для нахождения их параметров. По нашему мнению, данный метод способен лишь несколько сгладить краевые эффекты, хотя сами авторы утверждают о полном их подавлении.

В том же 2006 г. J. Li и I. Morozov на ежегодной конференции SEG выступают с докладом, где критикуют дискретное распределение точечных источников, и предлагают перейти к непрерывному эквивалентному слою [Li, Morozov, 2006]. При этом в численной реализации метода непрерывный слой так же дискретизируется и представляется в виде плотной упаковки прямоугольных призм. На наш взгляд, преимущества такой конструкции в плане повышения точности вычислений сомнительны, а сам процесс аппроксимации усложняется.

R. Jia и R.W. J. Groom в 2007 г. на этой же конференции делают доклад о применении техники эквивалентных источников для вычисления градиентов магнитного поля, где сравнивают этот метод с аппроксимацией дискретных значений поля рядом Фурье [Jia, Groom, 2007], еще раз говоря о необходимости корректного выбора глубины размещения источников. Интересно, что при создании эквивалентного слоя они также расширяют его на 20% за пределы площади задания поля, по всей видимости, для подавления краевых эффектов. K. Davis и Ya Li в свою очередь представляют способ построения эквивалентного слоя с одновременным использованием полного вектора аномального магнитного поля и его градиентов [Davis, Li, 2007].

В России в 2006–2008 гг. появляется серия статей, посвященных быстрым алгоритмам построения многоуровневых аппроксимаций, ключевой процедурой в которых является декомпозиция сложной вычислительной задачи большой или сверхбольшой размерности на ряд подзадач меньшей размерности с применением бинарных деревьев – метода квадродерева и быстрого вейвлет-преобразования на основе базисных функций Хаара [Пугин, 2006; 2007; Dolgal, Pugin, 2006; Долгаль, Пугин, 2006a; 2006б; Долгаль, 2007; Dolgal, Simanov, 2008]. В алгоритмах применяются процедуры сжатия информации с отбрасыванием избыточных данных о поле, что приводит к существенному уменьшению размерности СЛАУ на каждой стадии декомпозиции задачи. Однако отбрасывание проходит по “мягкому” варианту, т.е. участки поля, отброшенные из аппроксимации на начальной стадии декомпозиции, в обязательном порядке анализируются вновь и поступают в аппроксимацию на последующих стадиях во избежание неконтролируемого увеличения погрешности приближения. Разрастание аппроксимационной конструкции идет от нижних уровней глубин эквивалентных источников к верхним. В совокупности с истокообразной аппроксимацией, данные методы позволили разработать точные и весьма экономичные решения в плане вычислительных затрат, которые впоследствии вошли в диссертацию одного из авторов [Пугин, 2007].

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

В 2007 г. В.Н. Страхов публикует работу [Страхов, 2007], в которой рассматривает широкий спектр вопросов, касающихся построения линейных аналитических аппроксимаций при решении задач геодезии, геоинформатики и гравиметрии, а в 2009 г. выходит обширная и весьма интересная монография [Страхов и др., 2009], в которой подробно излагаются теоретические основы аппроксимационного подхода на основе интегральных представлений, а также методы F- и S-аппроксимации и результаты их опробования на модельных и практических примерах.

Спустя 11 лет после выхода работы [Долгаль, 1999a], предложенный в ней подход был реализован в виде компьютерной технологии для учета влияния боковых источников при трансформациях геопотенциальных полей [Долгаль, Пугин, 2008; Веселкова, Пугин, 2009; 2010]. Предлагалось расширить площадь задания поля за счет вовлечения материалов мелкомасштабных съемок из банков цифровых данных, а аппроксимационную конструкцию строить в виде нескольких слоев эквивалентных источников, нижние из которых аппроксимируют мелкомасштабные данные на большей площади, а верхние – разность между исходными данными (имеющими более крупный масштаб и точность) и результатами аппроксимации мелкомасштабной съемки. Позднее способ трансформаций с учетом влияния боковых источников был оформлен в виде патента [Бычков и др., 2010] и применялся для интерпретации результатов аэромагнитных съемок 1 : 100 000 – 1 : 250 000 масштабов на северо-западе Сибирской платформы [Долгаль и др., 2020].

В период с 2008 по 2013 гг. выходит большая серия статей группы пермских геофизиков, работающей под руководством А.С. Долгаля, с примерами практического применения различных алгоритмов истокообразной аппроксимации для решения прикладных задач гравиметрии и магнитометрии [Долгаль, Христенко, 2008; Шархимуллин, Долгаль, 2010; Долгаль и др., 2011a; 2011б], в том числе для трехмерной интерполяции и фильтрации геопотенциальных полей [Батырева, 2008; Батырева, Долгаль, 2009; Новикова и др., 2013].

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

В 2010 г. Ya Li и D.W. Oldenburg [Li, Oldenburg, 2010] предлагают быстрый алгоритм построения слоя эквивалентных источников для работы с данными магнитной съемки, имеющими большую размерность. Источники по-прежнему предлагается располагать в один конкордантный поверхности наблюдений слой либо размещать их на горизонтальной плоскости. Перед решением СЛАУ нужно выполнять разрежение матрицы ее коэффициентов на основе вейвлет-анализа и последующий синтез с отбрасыванием малых вейвлет-коэффициентов. На наш взгляд, данный подход оправдан лишь в равнинных территориях, поскольку любая предварительная фильтрация потенциального поля, измеренного в горной местности, не учитывающая рельеф земной поверхности наряду со значениями поля, может приводить к существенным искажениям результатов вследствие неучтенного эффекта разновысотности наблюдений. В этом кроется существенный недостаток алгоритма в сравнении с более ранними способами, учитывающими фрактальную структуру данных, предложенными в работах А.С. Долгаля и его учеников.

В том же 2010 г. выходит работа [Salem et al., 2010], в которой авторы предлагают перед построением эквивалентной модели источников выполнять фильтрацию высокочастотной компоненты магнитного поля с применением фильтра нижних частот, что близко к подходу, предложенному А.С. Долгалем в работе [Долгаль, 2002], однако предлагаемый иностранными авторами подход менее автоматизирован и более сложен в плане вычислений. В качестве эквивалентных источников использованы горизонтальные цилиндры бесконечного простирания, что обусловлено профильным вариантом алгоритма.

В 2011 г. выходит статья [Степанова, 2011], в которой рассматриваются аппроксимации рельефа земной поверхности и расчет топопоправок в рамках метода линейных интегральных представлений. Позднее выходит серия статей на данную тему [Stepanova et al., 2019; 2020a; 2020b]. Определенный интерес представляет статья [Gudkova et al., 2020], посвященная истокообразной аппроксимации рельефа поверхности и гравитационного поля Марса.

В 2012 г. П.С. Мартышко, Н.В. Федорова, Д.В. Геймайдинов провели анализ магнитного поля Урала (координаты 52°–64° с.ш., 54°–66° в.д.) на основе дискретизации интегрального уравнения Фредгольма 1-го рода и решения СЛАУ (размер сетки 512 × 512 точек) с регуляризацией по Лаврентьеву на суперкомпьютере “Уран”. С помощью пересчета поля на высоты от 1 до 200 км были выделены его разночастотные составляющие, отражающие распределение намагниченности в земной коре. Следует отметить, что карты магнитных аномалий на высотах 100–200 км хорошо согласуются с данными спутниковых магнитометрических наблюдений [Мартышко и др., 2012]. В дальнейшем был выполнен подбор магнитного поля сингулярными источниками для обширной территории Уральского региона с последующим вычислением ряда трансформант, в т.ч. для приведения к полюсу региональных магнитных аномалий [Мартышко и др., 2016].

В 2013 г. Jr. Oliveira с соавторами [Oliveira et al., 2013] предлагают оригинальный способ определения физических параметров эквивалентных источников в слое на основе метода локальных полиномов. В скользящем окне при каждом его положении определяется распределение физических параметров элементарных источников, попадающих в пределы окна, на основании подобранного полинома с некоторыми коэффициентами. В последующем полиномы “сшиваются” в единое распределение физических параметров источников для всего аппроксимирующего слоя.

В 2014 г. отечественными исследователями был предложен еще один оригинальный алгоритм построения истокообразной аппроксимации геопотенциального поля, реализующий концепцию “случайного вброса” или случайного поиска положения источника в сеточной модели [Пугин и др., 2014]. Впоследствии вышла статья с расширенным составом авторов, обобщающая все ранее предложенные ими алгоритмы построения многоуровневых аппроксимационных конструкций [Balk et al., 2016].

В том же 2014 г. коллектив авторов статьи [Kara et al., 2014] демонстрируют на практическом примере применение истокообразных аппроксимаций для трансформаций аномалий силы тяжести в условиях горного рельефа на полуострове Biga Peninsula в Турции.

Технология эквивалентных источников была модифицирована и использована для аппроксимации глобального магнитного поля Земли, полученного в результате спутниковых измерений CHAMP в период 2009–2010 гг. [Kother et al., 2015]. Аппроксимация производилась в сферических координатах с расположением эквивалентных источников на глубине 100 км ниже земной поверхности в узлах грида с икосаэдрическими ячейками. Предварительно из геомагнитного поля было исключено влияние земного ядра и магнитосферных источников с применением модели CHAOS-4.

В 2015 г. И.Э. Степановой и Д.Н. Раевским предложен модифицированный метод S-аппроксимаций геопотенциальных полей [Raevskiy, Stepanova, 2015a; 2015b], вошедший впоследствии в диссертационную работу [Раевский, 2016]. Эффективность метода подтверждена на модельных примерах, в работе [Stepanova et al., 2017] представлен практический пример трансформации гравитационного поля северной части Тихого океана (13341 пункт измерений, квадратная сеть точек с шагом 0.5°). Позднее был опубликован цикл статей, посвященных повышению эффективности данного метода, в том числе применительно к интерпретации больших объемов данных гравиметрической и магнитной съемок с использованием предложенной модификации метода S-аппроксимаций [Stepanova et al., 2016; 2017; Koneshov et al., 2017; Stepanova et al., 2020], а также для вычисления высот относительно геоида [Koneshov, Stepanova, 2019].

С 2018 по 2020 гг. выходит серия статей, посвященных различным модификациям метода линейных интегральных представлений и их совместному использованию [Stepanova et al., 2018a; 2018b; 2019].

В 2016 г. C. Martinez и Ya Li [Martinez, Li, 2016] приводят пример использования технологии эквивалентных источников для аппроксимации данных гравитационной градиентометрической съемки. Располагая эквивалентные источники в один слой и решая линейную обратную задачу с регуляризацией достаточно сложным путем, авторы вынуждены искать дополнительные критерии, например, использовать L-кривую Хансена, для обеспечения приемлемой точности аппроксимации аномальных градиентов тяжести, не допуская при этом увеличения вычислительных затрат.

В опубликованной в 2017 г. статье M. Pilkington и O. Boulanger сравнивают различные алгоритмы продолжения потенциального поля с поверхности измерений в верхнее полупространство на модельном и практическом примере, в том числе и с применением истокообразных аппроксимаций [Pilkington, Boulanger, 2017]. К сожалению, из российских разработок упоминается только один метод, предложенный В.Н. Страховым в 1965 г., который не использует истокообразную аппроксимацию.

В 2017 г. выходят две статьи [Siqueira et al., 2017a; 2017b], в которых предлагается для ускорения сходимости итерационного процесса нахождения масс эквивалентных источников в однослойной конструкции использовать некоторую систему параметров – весовых коэффициентов, зависящих от амплитуды аномалии поля в данной точке.

В 2017 г. С.А. Тихоцкий и Д.Ю. Шур [Тихоцкий, Шур, 2017] расширили доказательство эквивалентности метода истокообразной аппроксимации оптимальному фильтру Колмогорова–Винера при определенном выборе автокорреляционной функции на более общий случай нерегулярной сети и криволинейной поверхности измерений [Гордин и др., 1980].

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

В том же 2018 г. выходит статья [Ali et al., 2018], в которой на практическом примере продемонстрировано использование техники эквивалентных источников с однослойным расположением для пересчета поля на горизонтальную плоскость, расположенную всюду ниже земной поверхности, но выше эквивалентных источников. Истокообразная аппроксимация использовалась для вычисления производных гравитационного поля в нижнем полупространстве и определения по ним границ различных структурных элементов земной коры на территории межостровного моря Сулу в Тихом океане.

В этом же году M.H. Dransfield и T. Chen [Dransfield, Chen, 2019] приводят пример успешного практического применения истокообразной аппроксимации для трансформаций данных гравитационной градиентометрической съемки, выполненной с борта вертолета над вулканом ASO в Японии.

В 2019 г. Е.И. Есин с соавторами [Есин и др., 2019] на примере метода S-аппроксимаций показывают, что в условиях ограниченных площадей однослойные эквивалентные модели приводят к существенным искажениям трансформант при пересчетах магнитного поля в верхнее полупространство, особенно при наличии сильномагнитных боковых источников, на что также указывалось в работе, приведенной абзацем выше. Авторы предлагают вместо расширения площади использовать для аппроксимации данные аэросъемок на нескольких высотных уровнях, справедливо утверждая, что это повышает корректность построения эквивалентной модели. Однако при этом возникает необходимость получения данных на верхних уровнях измерений с меньшей погрешностью, чем на нижних уровнях, чего сложно добиться на практике.

В этом же году А.С. Долгаль для подавления влияния практической эквивалентности при решении обратной задачи гравиразведки монтажным методом предложил выполнять предварительное разделение гравитационного поля на составляющие, обусловленные объектами с разнознаковыми эффективными плотностями. Используется система точечных эквивалентных источников, последовательно погружающихся на глубину, аппроксимирующих интерпретируемое поле. Экспериментально установлено, что при некоторой глубине расположения источников возможно достаточно точное восстановление “положительной” и “отрицательной” компонент аномального поля [Долгаль, 2019].

На ежегодной конференции SEG 95 J. Capriotti и Ya Li представляют доклад [Capriotti, Li, 2019], в котором предлагают свой способ аппроксимации полного вектора гравитационного поля (всех трех компонент в декартовой системе координат, а не только вертикальной составляющей) с использованием регуляризации методом L – кривой Хансена при решении плохо обусловленной СЛАУ.

На этой же конференции другой коллектив авторов [Reis et al., 2019] предлагает способ определения направления вектора намагничения пород, создающих измеренное магнитное поле, посредством истокообразной аппроксимации с построением однослойной конструкции эквивалентных источников. В работе [Piauilino et al., 2019] повторно излагается способ быстрого построения аппроксимирующего эквивалентного слоя с применением регуляризатора, ускоряющего сходимость итерационного процесса нахождения физических параметров источников, с демонстрацией эффективности алгоритма на модельном и практическом примерах. Следует заметить, что успешное решение подобной задачи было реализовано в 2006 г. российскими геофизиками в программе REIST [Бабаянц и др., 2006].

В 2020 г. А.Ю. Давыденко представил устойчивый метод определения магнитной восприимчивости и вектора остаточной намагниченности трехмерных объектов на основе линейной регрессионной модели по методу эластичной сети [Давыденко, 2020].

На другой конференции представлен еще один доклад группы авторов [Piauilino et al., 2019], в котором они сравнивают два метода построения аналитической аппроксимации гравитационного поля: метод истокообразной аппроксимации в подходе C.N.G. Dampney и аппроксимацию рядом Фурье. В последующем обе аппроксимационные конструкции используются для вычисления тензора градиентов силы тяжести с целью сопоставления результатов. По мнению авторов, для модели локального источника первый метод является предпочтительным. При преобразовании суммы локальной и региональной аномалий оба метода дают неприемлемые результаты.

В одной из работ К. Mendonça предлагает способ решения обратной задачи для построения аппроксимирующего слоя эквивалентных источников, пытаясь разделить аномалиеобразующие объекты, имеющие различную геологическую природу [Mendonça, 2020]. В последующем подпространства, содержащие эквивалентные источники с положительными и отрицательными значениями эффективной плотности, используются раздельно для вычисления суммарной массы тел, создающих положительные и отрицательные аномалии гравитационного поля. Приведен пример успешного практического применения метода при интерпретации материалов аэрогравиметрической съемки в провинции Карахас (Бразилия).

В 2020 г. D. Takahashi с соавторами приводят подробное описание нового алгоритма, реализующего быстрое построение эквивалентного слоя [Takahashi et al., 2020]. Авторы применяют к матрице значений истокообразной функции двумерное дискретное преобразование Фурье, перегруппировывая ее элементы и получая матрицы иного типа, новую систему уравнений, которая решается непосредственно в частотной области. При этом матрица коэффициентов истокообразной функции рассматривается как блочная матрица Теплица, что отвечает лишь идеализированному случаю задания потенциального поля на горизонтальной плоскости. Однако именно такой случай – аппроксимация данных аэрогравиметрической съемки, выполненной на постоянной барометрической высоте, приводится в качестве примера, демонстрирующего очень высокое быстродействие алгоритма.

Наконец в 2020 г. появляется статья [Zuo et al., 2020], в которой авторы предлагают использовать многослойную конструкцию эквивалентных источников для аппроксимации данных аэромагнитной съемки, выполненной с обтеканием рельефа, с последующим продолжением поля вниз вплоть до земной поверхности. Детали алгоритма не раскрываются, приводится лишь функционал сложного вида, который необходимо минимизировать. Примечательно, что возврат к разноглубинному расположению эквивалентных источников при истокообразной аппроксимации поля в работах зарубежных авторов происходит лишь спустя много лет после выхода работы [Cordell, 1992], в то время как российские геофизики, опираясь на результаты В.И. Аронова стали широко эксплуатировать данный подход, начиная с конца 90-х годов ХХ века.

Метод трехмерной интерполяции гравиметрических и магнитных данных был предложен S.R. Soler и L. Uieda в работе [Soler, Uieda, 2020]. Авторы используют разбиение области задания поля на ряд блоков с вычислением средних для каждого из них значения поля. Эти средние значения аппроксимируются набором точечных источников, находящихся под центрами блоков на разных глубинах. Метод обладает высокой точностью и быстродействием. Этот подход был развит на случай больших объемов данных (миллионы точек наблюдений) и показал хорошие результаты [Soler, Uieda, 2021]. Рекомендовано выбирать размер блоков, близкий к шагу результативной сети гридирования. Наилучшим компромиссом между скоростью расчетов и точностью результатов является взаимное перекрытие блоков на 50% площади.

В 2015–2020 гг. группа пермских геофизиков расширила область практического применения истокообразной аппроксимации геофизических полей [Долгаль и др., 2017; Долгаль, 2020]. Новую геологическую информацию удалось получить при пересчете геопотенциальных полей Норильского района и ультраосновного массива Кондер на сглаженные синтезированные поверхности, расположенные выше земной поверхности [Долгаль, Дудин, 2015]. В качестве альтернативы сплайн-аппроксимации при построении огибающих сигнала в методе эмпирической модовой декомпозиции (EMD) рассматривается аномальный эффект эквивалентных источников [Долгаль, Христенко, 2017]. Истокообразные аппроксимации выбраны в качестве инструмента оценки различий между радиальной и вертикальной составляющими силы тяжести для больших территорий (с размерами от 100 км в поперечнике), характеризующей влияние сферичности Земли при выполнении интерпретационных построений в прямоугольной системе координат [Дол-галь и др., 2018]. Аппроксимационный пересчет полей в верхнее полупространство с одновременным их дифференциальным преобразованием начал использоваться при интерпретационной томографии с целью приближенного выделения эффекта от горизонтального слоя горных пород. Эти технологии активно применяются при изучении подводных вулканов Курильской островной дуги по данным гидромагнитной съемки [Блох и др., 2018; 2019; Долгаль и др., 2021].

В 2020 г. вышел справочник [Параметры Земли…, 2020] с параметрами Земли 1990 г. в котором приводятся глобальная, региональная и локальная модели аномального гравитационного поля Земли (АПГЗ), построенные с использованием метода истокообразной аппроксимации. Описание аналитических моделей содержится в приложениях в виде параметров точечных источников (аномальных масс и геоцентрических координат). Вообще же одним из основных приложений метода точечных масс являлся учет АГПЗ для повышения точности баллистических расчетов. Несмотря на то, что уже получены коэффициенты разложения возмущающего потенциала Земли в ряд по сферическим функциям до 360-й и более высокой степени, метод точечных масс продолжает успешно применяться при расчете траектории полета ракеты. Данной теме посвящен отдельный раздел энциклопедического издания – книги “Гравиметрия и геодезия” [Гравиметрия и геодезия, 2010].

Эпоха развития и совершенствования компьютерных технологий истокообразной аппроксимации в прикладной геофизике продолжается. Благодаря новым вычислительных возможностям, можно ожидать совершенно неожиданных подходов к решению проблемы моделирования геопотенциальных полей. Подводя промежуточные итоги, отметим во многом “параллельное” развитие отечественных и зарубежных технологий, обусловленное недостаточным обменом информацией. Причем складывается впечатление, что российские ученые лучше знакомы с англоязычными работами, нежели их иностранные коллеги – с российскими. Можно выделить основные тенденции развития метода истокообразной аппроксимации:

1) расширение области применения многоуровенных аппроксимационных конструкций, адаптированных к структуре сети задания и морфологическим особенностям поля;

2) совершенствование методов решения СЛАУ с числом неизвестных О(10⁶), направленное на снижение времени расчета параметров аппроксимационной конструкции и повышение точности результатов;

3) декомпозиция общей задачи истокообразной аппроксимации на ряд подзадач меньшей размерности;

4) обсуждение проблемы краевых эффектов и борьба с их проявлениями;

5) определение физических параметров эквивалентных источников без решения СЛАУ (одномерная оптимизация, статистические зависимости);

6) применение для 3D-интерполяции данных полевых магнитометрических и гравиметрических измерений и фильтрации осложняющих их помех;

7) использование для расчета трансформант геофизических полей ‒ продолжение в верхнее полупространство, подавление эффекта разновысотности путем пересчета на горизонтальную плоскость, пересчет на различные криволинейные поверхности – драппинг и т.п.;

8) разработка алгоритмов, учитывающих сферообразную форму Земли, предназначенных для построения глобальных и региональных моделей геофизических полей;

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

Следует также отметить неуклонно возрастающее количество публикаций в области практического применения метода истокообразной аппроксимации и его существенную технологическую адаптацию к решению широкого круга задач прикладной геофизики.