Большое число задач нефтяной геологии и разведочной геофизики решается с помощью карт. Карты служат примером графических моделей геологических поверхностей, отражают их пространственную изменчивость и являются непрерывным, гладким восполнением дискретных наблюдений. В настоящее время существует большое число алгоритмов и программ картопостроения, использующих разнообразные интерполяционные или аппроксимационные методы приближения. Не вдаваясь в подробный анализ существующих работ по моделированию геологических поверхностей, отметим, что большинство из них не содержит принципиально новых идей, и задача картопостроения в них рассматривается как формальное приближение исходных данных некоторыми функциями. Разнообразие заключается в выборе класса функций и критериев приближения.
Развивается и принципиально иной подход к проблеме построения геологических поверхностей, использующий дополнительную информацию о свойствах восстанавливаемой поверхности. К ним отнесем крайгинги, приближение сплайн-функциями, метод конечных элементов. Дополнительная информация о свойствах поля в этих подходах представлена различными способами. Например, в крайгингах это вариограммы, являющиеся оценкой автокорреляционной функции некоторого случайного поля. В методе сплайн-функций это вариационный функционал. В методе конечных элементов, как и в методе сплайн-функций, в качестве дополнительной информации о свойствах поля задаются линейные операторы, описывающие свойства картируемой поверхности и соответствующие им функционалы. Существует формальная аналогия между крайгингами и сплайнами, методами приближения сплайн-функциями и конечными элементами. В этой работе мы остановимся на методе конечных элементов и покажем, какие содержательные задачи могут быть решены в рамках этого подхода.
Математическая постановка задачи
Пусть H (Ω) – гильбертово пространство, которое ассоциируется с пространством состояний полей геологических параметров, Ω — ограниченная замкнутая область в R2. Пусть А, B – линейные непрерывные операторы, действующие из H (Ω) в L2 (Ω):
(1)
где v∈ H (Ω) — некоторая известная функция, А, B – дифференциальные операторы общего вида
(2)
Коэффициенты операторов могут быть постоянными, а также известными функциями координат точек области. Оператор В также может зависеть от неизвестных числовых коэффициентов λ=(λ1, …,λ r), которые являются модельными параметрами процесса, приведшего к реализации конкретного геологического поля. Равенство (1) содержит априорную информацию о свойствах поля геологических параметров, т.е. является математической моделью геологических знаний. Уравнение (1) – есть некоторая идеализация: в действительности, точное равенство никогда не выполняется, так как поведение реального поля зависит от бесконечного числа факторов. Всегда существует отличная от нуля функция h∈ L2 (Ω), норма которой
показывает степень отклонения пробной функции u∈ H (Ω) от равенства (1). Чем больше это отклонение, тем меньше «вероятность», что функция u может быть выбрана в качестве допустимой модели поля геологического параметра.
Пусть Li, Ni – непрерывные линейные функционалы в пространстве H (Ω), i=1,..n. Запишем
(3)
Эту величину назовем модельным измерением, а функционалы L, N – моделью измерительного прибора, где p∈ H (Ω) – некоторая фиксированная функция. В рассматриваемой постановке эти функционалы имеют вид дифференциальных форм, вычисленных в точках с координатами Qi={x1i, x2i}. Кроме этого функционал N может зависеть от ряда неизвестных заранее числовых параметров μ=(μ1, …,μ m) . Таким образом, величина модельного измерения зависит от пробной функции u и вектора μ. Положим, что реальные измерения поля геологического параметра z содержат случайную ошибку, имеющую нормальный закон распределения ξ∈ N (0,σ2). Оценку математического ожидания наблюденных значений поля будем отыскивать в виде (3).
Задачу (1) и (3) с учетом сделанных допущений можно свести к задаче минимизации функционала:
(4)
где S — стабилизирующий функционал, определяющий гладкость решения задачи, α≥0, β, ρ — положительные числа. В качестве S выберем положительную линейную комбинацию квадратичных функционалов:
(5)
Задача (4) решается методом конечных элементов. В качестве базиса выбраны бикубические В-сплайны, позволяющие свести вариационную задачу (4) к решению системы линейных алгебраических уравнений с разряженной матрицей ленточного типа с окаймлением. В некоторых случаях система уравнений имеет плохо обусловленную матрицу, что ведет к неустойчивости. Для того, чтобы регуляризировать решение, добавим к (4) функционалы вида ||λ-λ0||2l2 и ||μ-μ0||2l2 , где λ0, μ0 – априори известные значения модельных параметров, полученные, например, по другим площадям.
Примеры построения моделей полей геологических параметров
Реализация предложенной постановки осуществлена в программах Medium (Windows) и ГеоФлюид (Unix, Linux). Рассмотрим примеры, в которых используются различные предположения относительно свойств восстанавливаемого поля и способов его измерения.
Простейший случай: построение структурной поверхности по значениям поля в некоторых точках с использованием сглаживающего функционала. В этом случае задача (4) будет иметь вид:
(6)
Построенные карты в зависимости от вида стабилизатора представлены на рис. 1а, б.
Рис. 1. Пример построения карт с разными стабилизаторами: а) минимум площади поверхности, б) минимум кривизны.
При построении первой карты в качестве стабилизатора использовался функционал S1 — минимум кривизны, во втором примере S2 – соответствующий квадрату площади поверхности. Эти функционалы отвечают разным свойствам поверхности. Если минимум кривизны поверхности достигается на всех плоскостях, то минимум площади достигается на функциях, имеющих постоянное значение. Для такого функционала выполняются минимаксные условия, то есть функция u (x, y), минимизирующая функционал (6), может принимать максимальные и минимальные значения только в точках наблюдения. Таким образом, по одним и тем же исходным данным можно построить разные поверхности, используя информацию о вероятных свойствах картируемого поля.
Пример конформного разреза
Термин конформный понимается как сходство форм объектов. Разрез осадочного чехла назовем конформным, если геологические границы подобны. Математически это может быть записано равенством
(7)
где λ — некоторое число. Если 0<λ<1, то границы будут выглядеть, как на рис.2, то есть вышезалегающие границы более пологие, чем нижние.
Рис. 2. Модель конформного разреза
На рис. 3 показан пример построения структурной карты по данным бурения с учетом сейсмической поверхности (изолинии синего цвета), автоматически оценивался коэффициент λ, который оказался равным 0.609. Это значит, что с этим коэффициентом постро енная поверхность лучше согласуется с данными бурения и сейсмическими построениями. Видно, что структурная поверхность повторяет форму сейсмического горизонта.
Рис. 3. Структурная карта, построенная по данным бурения (изолинии черного цвета) с учетом конформности сейсмическому отражающему горизонту (изолинии синего цвета)
Модель смещения свода
Рассмотрим более сложный пример. Известно, что иногда в осадочной толще наблюдается смещение свода структур. Это можно объяснить наличием тренда скорости осадконакопления. Модель смещения свода может быть построена как суперпозиция двух моделей: чистого сдвига и конформности.
Рис. 4. Модель смещения сводов структур
На рис. 4 представлен разрез, соответствующий описываемой модели, которая. математически может быть записана в виде равенства
(8)
где λ — конформный коэффициент, v — вектор смещения структуры. Подставляя в качестве модели знания соотношение (8), мы можем построить поверхность Z1 по данным в скважинах при условии, что за дана поверхность Z2. Одновременно с задачей картопостроения отыскивается вектор смещения и коэффициент конформности.
Рис. 5. Пример построения структурной карты с учетом модели смещения сводов
На рис. 5 черные линии показывают структурную карту по горизонту Z2, построенную по данным скважин. Точки — это исходные данные, по которым строилась структурная карта по горизонту Z1 (синие линии). Абсолютные отметки в скважинах были взяты по горизонту Z2 в точках, сдвинутых от скважин на величину e=(-1,1) км по осям Х, У, соответственно, что равносильно сдвигу поверхности Z1 на величину – e. По этим данным и модели (8) была построена структурная карта, и оценен вектор смещения v, который оказался равным (1.066,-0.993), что хорошо согласуется с его теоретическим значением. Построенная структурная карта по поверхности Z1 сдвинута относительно поверхности Z2 на вектор v. Особенно это хорошо видно в северо-западной части площади, где отсутствуют скважинные данные и карта строится только в соответствии с моделью сдвига.
Использование знаний о свойствах картируемого поля позволяет улучшить прогнозирующие свойства модели. Рассмотрим нетривиальный пример деформации тяжелой, упругой мембраны, закрепленной на прямоугольном контуре.
Интересующая нас задача такова:
(9)
где Ω — область на плоскости x, y: {(x, y): 0≤ x≤1, -1≤ y≤1}. Эта задача возникает при определении положения равновесия тонкой мембраны с закрепленными краями, смещаемой силой 2 единицы на единицу площади в положительном направлении оси z. С помощью стандартного метода разделения переменных можно показать, что решение задачи (9) записывается в виде
(10)
Данная формула будет использована при оценке точности описанной ниже численной процедуры.
Построим в области Ω прямоугольную сетку с шагом по x и y, равным 0.05. На границе области в узлах сетки положим нулевые значения функции u. В соответствии с целевым функционалом (4) задача (9) сводится к минимизации следующего функционала:
Решение этой задачи осуществлялось Картопостроителем программы ГеоФлюид. Для сравнения были взяты вычисленные суммы первых 150 членов ряда Фурье (10). Результаты приведены в табл. 1.
Таблица 1.
Пример использования разнородной информации
Исходные данные: координаты и абсолютные отметки в скважинах (точки черного цвета), координаты и значения тангенсов углов наклона структуры в случайных направлениях (точки синего цвета). Производные по направлениям вычислялись с известной структурной поверхности в точках, полученных с помощью датчика случайных чисел, в случайных направлениях. Всего таких точек 100000, равномерно охватывающих область картирования. Формально задача записывается в следующем виде
(12)
где — абсолютные отметки в скважинах, — единичный вектор, вдоль которого вычислялась производная, — значение тангенса угла наклона структуры в данном направлении.
На рис. 6 представлена структурная карта, построенная по разнородным данным, отвечающая решению задачи (12). Для сравнения синим цветом показаны изолинии исходной структурной поверхности, с которой снимались производные по направлениям. Можно видеть, что соответствие карт в области, где отсутствуют прямые наблюдения (значения абсолютных отметок), вполне удовлетворительное.
Рис. 6. Структурная карта, построенная одновременно по абсолютным отметкам и значениям градиентов
Таким образом, метод построения карт, основанный на конечных элементах, позволяет использовать дополнительную информацию о свойствах картируемого поля. Эта информация представляется в виде линейных дифференциальных уравнений в частных производных, которые служат моделями реальных физических процессов. Одновременно с задачей построения поверхности восстанавливаются неизвестные коэффициенты модели, например, коэффициент конформности, вектор смещения свода структуры. Приведено численное сравнение результатов решения уравнения Пуассона как способ тестирования метода на уравнениях эллиптического типа. Проверено, что метод картопостроения, основанный на конечных элементах, позволяет одновременно использовать разнородные по способам получения данные, например, абсолютные отметки в скважинах и производные по направлениям. Выбор в качестве базиса бикубических В-сплайнов позволяет быстро и эффективно решать все перечисленные задачи.
Литература
1. Аронов В.И. Принципы и методы использования косвенных данных при построении карт //Математические методы исследования в геологии / Экспресс – информация.
— М.: ВИЭМС. -1981.- № 4.
— С. 6–13.
2. Волков А.М. Решение практических задач на ЭВМ. -М.: Недра.-1980.-244 с.
3. Волков А.М., Пятков В.И., Торопов С.В. Сплайн-аппроксимация в задачах картирования // Про блемы нефти и газа Тюмени.
— 1978.- Вып. 39.- С.74–78.
4. Волков В.А. Моделирование геологических поверхностей в связи с задачами размещения скважин и установления достаточности разведки нефтяных и газовых месторождений // Обзорная серия: Математические методы исследований в геологии.
— М.: ВИЭМС.- 1977.- С.71.
5. Волков А.М. Геологическое картирование нефтегазоносных территорий с помощью ЭВМ. -М.: Недра.-1988.- 221 с.
6. Касти Дж., Калаба Р. Методы погружения в прикладной математике.- М.: Мир.- 1976.- 223 с.
7. Сидоров А.Н. Метод построения оптимальных карт // Тр. ЗапСибНИГНИ.- Тюмень.- 1984.- Вып.192.- С.32–39.
8. Сидоров А.Н. Математические методы обработки и интерпретации геолого– геофизической информации на примере построения карт геологических параметров // Проблемы нефти и газа Тюмени. – Тюмень.- 1979.
— Вып. 42.- С. 59–64.
9. Сидоров А.Н., Торопов С.В. Моделирование структурных поверхностей, осложненных дизъюнктивными нарушениями, методом граничных интегральных уравнений //Строение земной коры Западной Сибири/ Тр.ЗапСибНИГНИ. .-Тюмень.
— 1989.- С.90–96.
Автор: Плавник А.Г., Сидоров А.Н., Шутов М.С. (ГУП ХМАО НАЦ РН им. В.И.Шпильмана)
