Содержание к диссертации
Введение
1 Постановка задачи геологического моделирования. Современные подходы к построению полей геологических пара метров. 13
1.1 Виды и основные способы представления трехмерных геоло гических моделей 13
1.1.1 Цифровая трехмерная адресная геологическая модель (ЦТАГМ) 13
1.1.2 Цифровая трехмерная адресная фильтрационная модель (ЦТАФМ) 16
1.2 Типы исходных данных и априорной информации для трех мерного моделирования 16
1.2.1 Типы эмпирических исходных данных о геологических объектах 17
1.2.2 Примеры априорных сведений о геологических объектах 19
1.3 Современные подходы к построению трехмерных геологиче ских моделей 20
1.3.1 Построение трехмерных геологических моделей методами послойной интерполяции 24
1.3.2 Построение трехмерных геологических моделей с использованием методов крайгипга 29
1.3.3 Последовательные методы моделирования полей геологических параметров 47
1.3.4 Использование метода имитации отжига для построения и корректировки геологических моделей. 5G
2 Построение полей геолого-геофизических параметров с использованием В-сплайнов ( авторский подход применения априорной информации в виде 2D карт). 62
2.1 Интерполяционная и аппроксимационная постановки задачи моделирования поля геологического параметра G3
2.2 Поиск нормального решения задачи моделирования 3D поля геологического параметра 72
2.3 Основные результаты численных экспериментов 80
3 Новый метод построения 3D геологических моделей. Использование крайгинга с возможностью учета 2D геологи ческих карт . 88
3.1 Постановка задачи. Теоретическое обоснование подхода к решению 88
3.2 Некоторые практически важные обобщения метода
3.3 Результаты численных экспериментов , 99
4 Построение 3D литолого-фациальных моделей методом имитации отжига с использованием разнородной априор ной информации. 103
4.1 Гиббсовские распределения, потенциальные поля и марковские случайные функции 103
4.2 Динамические методы Монте-Карло: моделирование гибб-совских распределений и метод имитации отжига
4.2.1 Метод Гиббса ... 10G
4.2.2 Метод имитации отжига 109
4.3 Пример построения целевой функции для моделирования по ля литологии методом имитации отжига 111
4,3.1 Результаты численных экспериментов 11G
4.4 Обоснование авторского алгоритма корректировки литолого-фациальных моделей методом имитации отжига 124
4.4.1 Численные эксперименты 131
Заключение. 139
Литература
- Цифровая трехмерная адресная фильтрационная модель (ЦТАФМ)
- Поиск нормального решения задачи моделирования 3D поля геологического параметра
- Некоторые практически важные обобщения метода
- Динамические методы Монте-Карло: моделирование гибб-совских распределений и метод имитации отжига
Введение к работе
Высокая эффективность эксплуатации нефтегазового месторождения возможна только при условии правильного планирования и своевременного проведения мероприятий по его разработке. Б настоящее время решения о проведении таких мероприятий принимаются на основе анализа особого рода информационных моделей - трехмерных геологических моделей разрабатываемых нефтегазоносных объектов. В связи с этим, очень большое практическое значение имеет качество таких моделей, их адекватность реальности.
Трехмерная геологическая модель, как частный случай информационной модели, состоит из позиционной и атрибутивной составляющей. Позиционная составляющая представляет собой геометрический каркас -трехмерную сетку, определяющую положение исследуемого объекта в пространстве. Атрибутивная составляющая - поля геологических параметров - описывает внутреннее строение объекта, оценивает свойства геологических пород в различных его точках. Задачи построения этих составляющих решаются, как правило, отдельно. При этом, если при построении геометрического каркаса сейсмические данные обеспечивают довольно точное решение (в пределах известной погрешности), то данные для построения полей геологических параметров гораздо более; скудны. Результатов непосредственного исследования нефтегазоносного объекта, как правило, недостаточно для получения полного представления о его строении. Более того, на практике имеющимся данным об объекте может соответствовать множество качественно различных моделей (см., например, [А.Д. Бекман, 2005], [А.Д. Бекмап, 2005(2)]). Такие проблемы характерны для многих задач информационного моделирования (см. [В.Г. Гитис, Б.В. Ермаков, 2004]), поэтому при решении таких задач считается необходимым, чтобы модель удовлетворяла не только результатам непосредственного наблюдения моделируемого объекта, но и разного рода дополнительным данным о нем: априорным знаниям, экспертным представленням и гипотезам. Часто данные такого рода являются трудноформализуемыми и представляют собой многолетний опыт исследователей, использование которого позволяет
закладывать в модели те или иные качественные особенности и тем самым улучшать их адекватность действительности.
По причинам актуальности задачи трехмерного геологического моделирования и практической необходимости совершенствования программных средств для такого моделирования, предметом настоящего исследования были выбраны алгоритмы построения полей геологических параметров. В то же время, объектом исследования стали различные формы задания априорных данных и способы их учета в процессе построения этих полей. Это обусловлено, с одной стороны, подчеркнутой выше важностью такого рода данных, а с другой стороны, тем, что проблеме учета таких данных в геологическом моделировании, по мнению автора, в современных исследованиях уделяется недостаточное внимание.
Моделирование полей параметров геологических объектов имеет множество особенностей, не характерных для задач математического моделирования, В частности, одной из таких особенностей, и пожалуй самой важной, является то, что в настоящее время полностью отсутствует детерминистская математическая модель, описывающая взаимосвязь между свойствами геологических пород в различных точках объекта. Более того, закономерности такой взаимосвязи не выделены достаточно точным образом, в связи с чем не возможно привлечь для моделирования геологической среды аппарат дифференциального исчисления. Это приводит к невозможности обоснованного использования для решения рассматриваемой задачи математических моделей в виде дифференциальных уравнений, с успехом применяющихся для нужд механики, теории упругости, гидродинамики и т.д. Привлекательность численных методов, основанных на таких моделях, с одной стороны в том, что они представляют собой хорошо изученный математический аппарат. С другой стороны, они позволяют точно (или в пределах известной погрешности) восстановить свойства моделируемого объекта во всей изучаемой области по начальным и граничным условиям. Отказ (пусть даже и вынужденный) от такого подхода приводит к полной неопределенности относительно свойств объекта вне скважин. Другими словами, множество возможных вариантов модели представляет собой множество всех возможных функций, принимающих известные нам значения в скважинах и любые значения вне скважин. Тем не менее, нужды практики требуют от исследователя выбрать наиболее правдоподобный (или один из наиболее правдоподобных) вариант модели объекта в качестве основы для принятия решений. Именно этой цели и служат все известные
численные методы моделирования полей геологических параметров.
По глубокому убеждению автора, единственным возможным критерием оценки качества модели, т.е. степени ее соответствия реальному объекту, может служить ее соответствие имеющейся информации об объекте. Исходя из этого, представляется разумным для решения каждой конкретной задачи использовать такой численный метод, который бы позволил учесть всю имеющуюся информацию об объекте. С точки зрения достоверности, априорную информацию целесообразно классифицировать на строгую, т.е. требующую точного соответствия, и нестрогую, требующую соответствия приближенного. Каждый из классов информации будет использован в рамках численного метода по-своему. Так, строгая информация будет использована для сужения множества возможных вариантов модели, из которых численный метод будет выбирать единственное решение. Примером строгих априорных знаний может служить интервал возможных значений исследуемого свойства, полученный в результате статистических исследований породы, составляющей моделируемый объект. Благодаря такому знанию, при моделировании указанного свойства, множество возможных решений может быть сужено за счет исключения из пего функций, реализующих значения вне заданного интервала.
Если строгая информация не определяет решение однозначно и в за-уженом множестве возможных вариантов остается некоторый произвол, то выделить единственное решение можно с помощью нестрогой информации. При этом, большинство известных методов реализуют один из двух подходов:
Детерминистский подход. Вся имеющаяся нестрогая информация объединяется в один функционал, заданный на множестве возможных вариантов модели, ограниченном строгими сведениями. Функционал позволяет сравнивать различные варианты модели по степени соответствия заложенной в функционал информации: обычно он задается таким образом, чтобы его наименьшие значения соответствовали вариантам, наилучшим с точки зрения имеющейся нестрогой информации. (Примеры такого рода функционалов приводятся в главах 1 - 3,) Если при этом функционал имеет не единственный минимум, т.е. не определяет решение задачи однозначно, то для обеспечения единственности решения могут потребоваться дополнительные предположения о решении, например, предположения о том, что модель
является непрерывно дифференцируемой функцией, удовлетворяет некоторому дифференциальному уравнению и т.п. Такие предположения служат для обеспечения корректности математической постановки задачи, на основании которой выбирается численный метод ее решения. Разумеется, желательно, чтобы эти предположения несли в себе какую-то информацию об объекте, но если такой информации нет, то они могут выбираться произвольно и выбирать любой вариант модели из множества, ограниченного имеющимися строгими и нестрогими сведениями.
Стохастический подход, предполагает задание па множестве возможных решений, ограниченном строгой информацией, вероятностного распределения. Распределение задается таким образом, чтобы серо ятность вариантов решения была "пропорциональна" их правдоподобности. Алгоритм выбора в этом случае представляет собой некоторый специализированный метод Монте-Карло, генерирующий небольшую выборку заданного распределения. Полученная выборка анализируется экспертами и один из вариантов объявляется актуальной моделью.
Любой известный численный метод может быть рассмотрен с позиций приведенной выше схемы. Тем не менее, зачастую в рамках распространенных численных методов используется довольно мало информации о моделируемом объекте. Приведем несколько примеров:
При использовании для нужд геологического моделирования В-сплайпов (см., например, [A.M. Волков, 1988]) кроме использования скважишюй информации предлагается в качестве регуляризатора использовать функционал, соответствующий уравнению Лапласа, Привлечение такой информации служит для выделения единственного решения из множества возможных, а также выражает мнение исследователя о том, что исследуемое поле ведет себя подобно мыльной пленке. В качестве строгой информации, ограничивающей множество возможных моделей, используется информация о том, что модель принадлежит пространству В-сплайнов, а, значит, представляет собой дважды непрерывно дифференцируемую функцию. Последнее также может являться экспертным мнением относительно моделируемого поля.
Один из традиционных геостатистических методов, простой край-гинг, помимо скважшшых данных использует предположение, что исследуемое поле представляет собой стационарную в широком смысле случайную функцию с заданной функцией ковариации. Исходя из этого предположения, множество возможных вариантов модели сужается до множества функций-оценок, заданных на некоторой конечной сетке, и представляющих в каждом узле этой сетки линейную комбинацию известных скважшшых значений. Критерием выбора единственного решения из этого множества является требование наименьшего среднеквадратического отклонения оценки от неизвестной функции в каждой исследуемой точке (ем. [М. Давид, 1980], [Ж. Матероп, 1968]).
Численный метод Sequential Gaussian Simulation (SGSIM, см., например, [P. Goovacrts, 1997], глава 8), часто используемый в зарубежных программных продуктах для геологического моделирования, относится к разряду стохастических и базируется на тех же предположениях относительно исследуемого поля, что и крайгинг. Как и в случае краіігшіга, множество возможных моделей сужается до множества функций, заданных на некоторой конечной сетке и удовлетворяющих скважииным данным. Вероятностное распределение формулируется на основании требования наилучшего удовлетворения заданной функции ковариации (или вариограммы). Выборка вероятных вариантов получается с помощью специализированного метода типа Монте-Карло. Заметим, что и гипотеза стационарности случайной функции и вариограмма вносятся в задачу извне для возможности осуществления численного метода. На практике проверка гипотезы стационарности является практически невозможной, а вариограмма выбирается из нескольких шаблонов, часто не имеющих отношения к исследуемому объекту.
Любой из численных методов, приведенных в примерах может быть дополнен новой информацией. В частности, в главах 2 и 3 даны модификации этих численных методов, позволяющие учитывать (строго или пестрого) априорную информацию в виде двумерных геологических карт. Относительно двумерных геологических карт требуется заметить, что они являются традиционным для отечественной науки типом геологических моделей. Методы геологического картирования совершенствовались деся-
тилетиями и несут в сеое опыт многих поколении геологов, который может оказаться весьма полезным при построении и трехмерных геологических моделей. Кроме того, двумерный способ представления информации является более удобным для восприятия человеком нежели трехмерный, В связи с этим, методики экспертного контроля трехмерной геологической модели зачастую базируются на анализе ее сечений и построенных по ней двумерных карт. Еще одним доводом за использование априорных двумерных карт при построении трехмерных геологических моделей является то. что, согласно существующим регламентам (см,, например, [Регламент по созданию постоянно действующих геолого-технологических моделей...]), подсчет запасов месторождений углеводородов проводится по двумерным моделям, в то время как для целей гидродинамического моделирования требуются модели трехмерные, причем соответствующие двумерным по основным интегральным характеристикам, например, объемам коллскторских пород, объемам нефте- и газонасыщенных пород, запасам углеводородов в пластовых условиях и т.д. Таким образом, создание численных методов, позволяющих строить трехмерные геологические модели в соответствии с имеющимися двумерными, является практически важной задачей.
Представленные в настоящей работе алгоритмы моделирования трехмерных полей геологических параметров были реализованы автором в модуле ''Gektra" программного комплекса "BASPRO Оптима" и используются для решения практических задач моделирования нефтегазоносных пластов в научно-исследовательском институте "ТюменНИИГипроГаз".
Цифровая трехмерная адресная фильтрационная модель (ЦТАФМ)
Все исходные данные о геологических объектах, которые доступны исследователю на практике и могут быть использованы в процессе геологического моделирования, целесообразно разделить на две группы: 1. Данные, полученные эмпирическим путем. 2. Априорные сведения и гипотезы.
Данные первой группы на практике более ценны, так как они являются результатами экспериментального исследования объекта. Данные, полученные таким путем, являются наиболее точными сведениями об объекте. Даже если эти данные являются приближенными, то, как правило, их погрешность известна. Данные такого рода должны иметь больший приоритет, нежели априорные сведения, так как неудовлетворение экспериментальным данным делает модель заведомо неадекватной. Ниже рассмотрим названные трупы данных подробнее.
Типы эмпирических исходных данных о геологических объектах.
Современные технические средства и методы исследования геологических объектов подразумевают два основных источника прямых эмпирических данных: сейсмические исследования и исследования в скважинах.
Сейсмические методы исследований" территории способны с некоторой точностью определить глубины залегания различных геологических объектов, а также форму границ между ними. Результатами таких исследований являются поверхности границ крупных геологических объектов - пластов, пачек. Данные такого рода в основном необходимы па этапе построения геометрического каркаса модели. Они определяют форму поверхностей, ограничивающих каркас будущей модели сверху и снизу. Кроме того, они позволяют выявить такие особенности строения как тектонические разломы, выклинивание геологических объектов, границы несогласного залегания и т.д. Такие особенности также должны учитываться при построении ГК.
Исследования скважин представлены двумя основными группами методов исследований: петрографические исследования образцов породы, полученных при бурешши скважины, и геофизические исследования свойств пород вдоль ствола скважин. На этапе моделирования полой геологических параметров исследователь, как правило, располагает уже результатами интерпретации сважинных исследований (РИ-ГИС), которые имеют следующий вид: траектория каждой скважины разделена на интервалы глубин, для каждого из которых определена порода и ее геолого-геофизические свойства, такие как пористость, проницаемость, характер насыщения и т.п. (рис. 1.2) Подробный обзор методов исследований в скважинах можно найти, например, в [В.М. Бондаренко, Г.В. Демура, A.M. Ларионов, 1986].
Благодаря результатам интерпретации исследований скважин исследователь обладает полной информацией о свойствах пород геологических объектов вдоль траектории скважин. Из этого следует, что для ячеек ГК модели, пересекаемых скважинами, известны представляющие их породы и все интересующие свойства этих пород. Таким образом, задача восстановления полей параметров геологической модели сводится к экстраполяции этих данных на всю область модели наиболее правдоподобным образом. Методов такой экстраполяции в настоящее время разработано достаточно много, подробнее о них будет сказано в пунктах 1.3.1 - 1.3.4. Однако, можно выделить общую закономерность для всех методов: каждый метод основывается на привлечении того или иного вида априорной информации. Это связано с тем, что задача экстраполяции скважинных данных на все ячейки ГК в общем случае с очевидностью является некорректной по
Адамару, так как имеет иеединствешюс решение. Обладая знанием о значениях исследуемого свойства в одних точках пространства мы не можем в общем случае ничего определенного сказать о значениях этого свойства в других точках пространства. Для того, чтобы это стало возможным требуются знания или предположения о том, что значения свойства в различных точках пространства связаны между собой, подчиняются какому-либо общему правилу пли закону. Выбирая различным образом это организующее правило мы можем приходить к различным корректным постановкам задачи экстраполяции. Например, можно предположить, что решение является многочленом трех переменных, или, что оно описывается рядом Фурье с некоторым числом членов. Каждая из этих постановок будет предполагать единственное решение задачи, причем для различных постановок, очевидно, различными будут и решения. В связи с этим, желая получить наиболее адекватную действительности модель, исследователь должен выбирать дополнительные сведения о поведении свойств объекта в пространстве исходя не столько из математических соображений, сколько из геологических. Другими словами, важно учесть при построении модели как можно больше имеющихся априорных сведений о модели. Некоторые из возможных видов априорных сведений, наиболее часто используемых на. практике, приведены ниже.
Поиск нормального решения задачи моделирования 3D поля геологического параметра
В случае, когда исследуемая задача является некорректной в смысле неединственности ее решения, добиться ее корректности можно, определив дополнительный критерий выбора из множества возможных решений единственного. Академик А.Н. Тихонов (см,, например, [А.Н. Тихонов, В,Я, Арсенин, 1979]) предлагал в качестве такого критерия в рамках многих практических использовать требование минимальности нормы решения в выбранном функциональном пространстве. Такие решения принято называть нормальными. Выбор нормы функционального пространства определяет в этом случае известные a-priori желаемые свойства решения. В частности, нормой может случить расстояние до некоторого эталонного" элемента пространства.
В рассматриваемой задаче в качестве дополнительной информации о решении можно использовать следующие соображения геологов: так как каждый слой модели несет в себе информацию о породах, формировавшихся в определенную геологическую эпоху, в схожих условиях, следовательно, внутри слои значения параметров пород должны быть близки - мало отличаться от константы. Такие требования к модели выдвигают, например в [А.Н. Сидоров, Н.Г. Хорошев, 1987]. В роли значения такой константы в некоторых методах оценки качества моделей (см,, например, [Пьяиков В.Н.. Сыртланов В,Р., Майсюк Д.М., 2003]) может выступать выборочное среднее по скважшшым измерениям в слое. Таким образом, представляется разумным к уже введенным в функционал требованиям добавить еще одно - наименьшее отклонение среднего значения каждого слоя от соответ ствуюіцеґо выборочного среднего по скважинам: I J м (jjYIYlYl аткфтіхі], Vij) - йк)1 = min, !=1 j=l m=\
Такая формулировка дополнительного условия направлена на наилучшее удовлетворение статистическому критерию оценки качества модели, описанному в [Пьяпков В.Н., Сыртлаиов В.Р.. Майсюк Д.М-, 2003]. С другоіі стороны, дополнительное условие может быть сформулировано и иначе: искомая функция в каждой исследуемой точке должна как можно меньше отличаться от выборочного среднего в слое: / J м { a mixij-yij) -ukf = min. (2.11)
В такой формулировке условие более соответствует приведенным выше геологическим соображениям, чем в предыдущем случае, когда к константе стремится лишь среднее значение функции в рассматриваемой области. В дальнейшем будем рассматривать именно такой подход к постановке задачи.
Замечание 2.2.1. Искомая модель в каждом из подпространств О С D интересует пас. не, для всех пар точек (х.у), а только для I х J точек, которые, соответствуют центрам столбцов модели и образуют сетку. Можно ввести в рассмотрение функциональное пространство D, являющееся проекцией D на эту сетку. В полученном конечномерном пространстве D выражение (2.11) при (фиксированном наборе значений ilk будет являться квадратом евклидовой нормы. Таким образом, использование слагаемого вида (2.11) приводит к нахождению нормальных решений исходных задал в интерполяционной и аппроксилищиопиой постановках. Корректность задач при использовании новых постановок будет доказана ниже.
Запишем функционал, соответствующий интерполяционной постановке задачи, дополненный новыми требованиями:
Ниже покажем, что полученная задача корректна. Для доказательства по- требуется следующая лемма.
Лемма 2.2.1. Пусть Р - матрица размерностью W х Л/, причем М W, a Q симметричная невырожденная матрица положительно определенной квадратичной формы, размерностью М х AL Тогда, в случае линейной независимости столбцов матрицы Р, матрица PQP также является невырожденной.
Доказательство: Представим матрицу Q в виде произведения Q — S х S1, где S - также имеет размерность М х М. Очевидно, что S - невырожденная, так как в противном слу.чае согласно теореме об определителе произведения вырожденной была бы и матрица Q. Следовательно, имеет место равенство:
PQPT = PSSTPT - PS{PS)T.
Будем смотреть на матрицу PQP1, как на матрицу некоторой, положительно определенной квадратичной форыы. Тогда существует невырожденная матрица N, приводящая квадратичную форму к каноническому виду Н: Н = NPQPrNT = NPS{PS)Ti\rT = NPS(NPSf. Обозначим В = NPS и предположим, что R = В В является вырожденной. Тогда, так как R является матрицей канонического вида, в пей должна быть по меньшей мере одна нулевая строка. Обозначим номер нулевой строки через п. Выразим элементы этой строки;l Заметим, что при j п равенство пулю элемента r„j возможно только в том случае, если равны нулю все элементы п-той строки матрицы В. Отсюда следует, что ранг матрицы В меньше W. Однако это не возможно, так как ранг матрицы Р равен W, а матрицы Лг и S являются невырожденными. Полученное противоречие доказывает невырожденность матрицы Я, а следовательно, и матрицы PQP . Лемма доказана.
Некоторые практически важные обобщения метода
Будем рассматривать задачу построения поля некоторого свойства геологического объекта. Предположим, что имеются следующие данные об объекте моделирования:
Геометрический каркас с размерностью сетки (/+1) х (J+l) х (К+1). Для простоты изложения ниже будем считать, что геометрический каркас разбивает исследуемую область на пропорциональные слои, то есть для каждой фиксированной точки (&\т/) мощность (протяженность по вертикали) всех К слоев предполагается одинаковой (см, рис. 3.1). В общем случае при решении рассматриваемой задачи от этого условия можно отказаться. Изменения алгоритма решения, связанного с отказом от такого условия будут указаны отдельно в пункте 3.2.
Скважипные данные: uklv,w — l...\V,k = 1.,, К, Величина ukw представляет собой среднее значение искомого параметра наблюда емое в W- ЮЇІ скважине на глубине k-того слоя модели. Траектории скважин предполагаются вертикальными, для u -той скважины она может быть задана в параметрическом виде
Пару (xw,yw) будем ниже называть координатами u- -той скважины. Общее число скважин будем обозначать через W.
Двумерная модель (геологическая карта) исследуемого свойства. На практике двумерная модель представляет собой массив значений V — {vjj},i — 1.../, j — 1... J, величина Vij задает вертикальное среднее значение моделируемого свойства для (i,j)-ioro столбца модели. Однако, для удобства изложения будем считать, что двумерная модель представляет собой непрерывную дифференцируемую функцию V(o y). Предполагается, что при построении карты использовались достаточно точные данные, мощные математические методы, а также представления опытного эксперта-геолога. Будем считать построенную карту объективной двумерной макромодельго свойства трехмерного тела, включающую в себя все имеющиеся эмпирические и неэмпирические сведения.
Статистические закономерности пространственного распространения значении исследуемого свойства: средние значения моделируемого свойства для каждого слоя модели тк,к — 1...К, функция кова-риации K(h).
Модель поля будем искать в виде семейства функций двух переменных {й (х. г/)} =1, подразумевая при этом, что при фиксированном к функция й (х,у) описывает поле в k-том слое геометрического каркаса. Функции йк(х.у) должны быть построены так. чтобы наилучшим образом удовлетворять имеющимся скважинным данным, априорной двумерной модели и данным статистического характера.
Возможность использования априорных данных статистического характера при решении рассматриваемой задачи связана с методом восстановления неизвестного поля. Для этого предполагается, что в каждом слое модели неизвестная функция исследуемого свойства U (х.у) является стационарной в широком смысле гауссовой случайной функцией, с математическим ожиданием т и ковариационной функцией K{h). На практике ковариационная функция может быть каким-либо образом оценена по имеющейся двумерной модели или скважинным данным, либо получена в ходе более ранних исследований сходных по геологической истории объектов. Исходя из указанных предположений функции й (х,у) в каждой фиксированной точке {х.у) можно искать как оценку неизвестной случайной функции U , т.е. как условное математическое ожидание:
Такой подход является широко распространенным в картопостроении, именно он приводит к численному методу "простой крайгинг" (Simple Kriging). Однако, для того чтобы использовать априорную информацию в виде двумерной карты, данный подход необходимо изменить. Действительно, для соответствия двумерной карте требуется, чтобы среднее значение свойства в каждом столбце модели равнялось a-priori заданному числу, но выполнение такого условия не может обеспечить метод, в рамках которого значения свойства в каждом слое восстанавливаются независимо от других слоев. Изложенный ниже авторский подход использует модифицированный метод простого крайпшга для достижения соответствия двумерной модели (см, также, [А.Д. Бекмап, 2005J).
Рассмотрим отдельно задачу восстановления значений свойства для (/,,/)-того столбца модели. Пусть центр проекции этого столбца на плоскость хОу имеет координаты (х, у). Таким образом, необходимо найти значения йк(х,у). к = 1...К, как наилучшие оценки неизвестных случайных величин /%г\y),k — l...K, при условии 1 к - Y,{ikd,y) = Vi3 = V(x,y). (3.1) k=l Предположим для начала, что математическое ожидание функций (Jk(x,y) равно нулю. Оценки {йк(х,у),к — 1... К] будем искать в виде линейной комбинации известных значений в точках, соответствующих скважинам: йк(х,у) = 2akwj{x,y)ukw и-=1 или. в матричной записи, uk(xty) = ak{uk)T. (3.2)
Определим неизвестные коэффициенты аки,(х, у) из условий минимизации среднеквадратических отклонений Е(йк(х,у) — (Jk(x,y))2 для всех к и выполнения условия соответствия двумерной и трехмерной моделей (3.1). Последнее условие может быть удовлетворено либо точно, либо приближенно. В первом случае подход будем называть интерполяционным, во втором - аппроксимационным. Рассмотрим интерполицоппый подход.
Динамические методы Монте-Карло: моделирование гибб-совских распределений и метод имитации отжига
По формуле Байеса хиіх) іх--У) Ріу) также является распределением вероятности на X и задаст вероятность попасть в состояние у, если на предыдущем шаге распределение было v{x). Заметим, что формула представляет собой матричное умножение строки на матрицу.
Аналогично можно показать, что матричное произведение двух марковских ядер также будет марковским ядром, а значит формула vP...P{\j), где множитель Р входит п раз, задает распределение вероятности появления состояния у, если известно распределение вероятностей и шагов назад было и{х).
Определение 10. Однородная цепь Маркова па конечном пространстве X определяется начальным распределением i/(x) и переходными вероятностями Р(х.у). Вероятность того, что цепь в моменты времени 0....,п будет последовательно принимать состояния a.Jo,a. i,a. 2,. .хп равняется и{х{))Р{х{).Х\)Р{х\,Х )).. ,Р{хп-\,хп). Бесконечная последова-телъпоежь (XQ, ... ,хп,...) называется траекторией марковской цепи.
В случае, если переходные вероятности зависят не только от состояния системы, но и от номера шага (момента времени), цепь называется неоднородной. Вероятность того, что неоднородная марковская цепь в моменты времени 0,..., гс будет последовательно принимать состояния хц, хих2.... х„ равняется U(XQ)P\{Х0,ХЛ)Р {Х\,Х2) ... Р„(ж„_і.хп).
Для неоднородной марковской цепи распределение вероятности па X па п-том шаге будет иР{ (х), где произведение понимается как матричное. В этом случае для любого состояния х мы знаем, с какой вероятностью система придет в это состояние, но нам не важно по какой именно траектории.
Рассмотрим еще несколько определений.
Определение 11. Распределение вероятностей /.t(x) называют инвариантным (стационарным) относительно мартовского ядра Р. если выполнено /хР" — //.
Определение 12. Марковское ядро, некоторая степень которого явля-Є7ПСЯ строго полооїсителгшой называют примитивным..
Рассмотрев пространство всех распределений вероятности на X и введя в нем метрику JJ\, можно доказать следующую теорему (см.,например, теорему 4.1 в [Г. Винклер, 2002]):
Теорема 4.2.1. Для примитивного марковского ядра Р па конечном пространстве существует единственное инвариантное распределение /.(, и При П —5- ОС иРп - /і. Сходимость равномерная по всем v. Покажем, как смоделировать поведение системы с распределением Т вероятностей вида Щх) — —%—. При этом будем придерживаться обозна чений предыдущего пункта. Переформулируем определение 3 для одноточечных множеств и гибб-совских случайных полой;
Определение 13. Локальными характеристиками случайного поля Щх) называются марковские ядра (переходные вероятности), определяемые для ячейки s Є S формулой: 0, иначе. = E(s:J54, .,e-"(-).
Другими словами, локальную характеристику Us(x, у) можно воспринимать, как вероятность перехода из состояния х в состояние у, отличное от х в одной только ячейке s. Она же определяет и вероятность обратного перехода. Используя локальные характеристики для всех ячеек s Є 5, можно определить марковское ядро при котором может быть достигнуто любое состояние системы: Р(х,у) = Щф,у)-...-Щм}{х1у) (4.4) Доказано (см., например, [S. Geman, D. Gemaii, 1984]), что для марковской цепи с ядром (4.4) справедливо следующее соотношение: lim иРп{х) = Щх) (4.5) п сс
Другими словами, относительная частота появления каждой конфигурации х Є X в марковской цепи с ядром (4.4) будет стремиться к Щх) при стремлении количества итераций к бесконечности. Доказательство справедливости этого соотношения следует из теоремы 4.2.1, если учесть известный факт, что гиббеовские поля инвариантны относительно своих локальных характеристик. Численный метод, использующий формулу 4.5 для численного моделирования гиббсовскнх полей, известен как метод Гиббсл или стохастическая релакелцгія.
Замечание 4.2.1. Порядок обхода ячеек (может также употребляться термин схема сканирования) в формуле (4.4) мооїсет быть произвольным. Как показано в [Г. Вииклер, 2002J схема сканирования мооїсет даже носить случайный характер. В последнем случае соотношение (4.5) останется справедливым, если вероятность оказаться, в любой ячейке $ Є S при, обходе является ненулевой.
Изложенный выше метод моделирования гиббсовскнх полей базируется на использовании однородных цепей Маркова, Ниже приведем подход, известный как метод имитации отоісига (Simulated Annealing), использующий неоднородные марковские цени для нахождения глобальных минимумов функций.
Рассмотрим следующее параметрическое семейство гиббсовскнх полей: П{Х) = Е е-""" (46) Параметр fl принято называть обратной температурой системы. Достаточно легко показать (см., например [Г, Винклер, 2002]), что lim П (яО = { № ХЄМ] ,з ос I 0, иначе.
При этом, для х Є М функция Р — П3{х) монотонно возрастает, для остальных - с некоторого момента становится убывающей. Такое свойство гиббсовскнх полей позволяют находить минимумы функции Н(х), Для этого можно построить специального вида неоднородную цепь Маркова, предельным распределением которой будет равномерное на множестве , Применяя ее к произвольному состоянию системы можно получать ее состояния, принадлежащие множеству М. Ниже кратко покажем, как это можно сделать.
