Результат интеллектуальной деятельности: ГИБРИДНЫЙ СПОСОБ ДЛЯ ПОЛНОВОЛНОВОЙ ИНВЕРСИИ С ИСПОЛЬЗОВАНИЕМ СПОСОБА ОДНОВРЕМЕННЫХ И ПОСЛЕДОВАТЕЛЬНЫХ ИСТОЧНИКОВ

Перекрестные ссылки на родственные заявки

Данная заявка притязает на приоритет предварительной заявки на патент США 61/386828, поданной 27 сентября 2010 года, которая содержится в данном документе по ссылке во всех разрешенных юрисдикциях. Данная заявка связана с конкретной заявкой, озаглавленной "Simultaneous Source encoding and Source Separation as the Practical Solution for Full Wavefield Inversion", и притязает на приоритет предварительной заявки на патент США 61/386831, поданной 27 сентября 2010 года, и заявки на патент США 12/903749, поданной 13 октября 2010 года. Эти родственные заявки также содержатся по ссылке в данном документе во всех разрешенных юрисдикциях.

Область техники, к которой относится изобретение

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

Уровень техники

Геофизическая инверсия [1, 2] нацелена на нахождение модели свойств геологической среды, которая оптимально объясняет наблюдаемые данные и удовлетворяет геологическим и геофизическим ограничениям. Существует большое число известных способов геофизической инверсии. Эти известные способы разделяются на одну из двух категорий, итеративную инверсию и неитеративную инверсию. Ниже приведены определения того, что обычно подразумевается под каждой из двух категорий:

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

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

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

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

Способы оптимизации функции затрат являются или локальными или глобальными [3]. Глобальные способы просто заключают в себе вычисление функции затрат S(M) для совокупности моделей {M₁, M₂, M₃,...} и выбор набора из одной или более моделей из этой совокупности, которые приблизительно минимизируют S(M). Если требуется дополнительное улучшение, этот новый выбранный набор моделей затем может быть использован в качестве основы, чтобы формировать новую совокупность моделей, которые могут быть снова протестированы относительно функции затрат S(M). Для глобальных способов каждая модель в тестовой совокупности может рассматриваться в качестве итерации, или на верхнем уровне каждый тестируемый набор совокупностей может считаться итерацией. Известные глобальные способы инверсии включают в себя метод Монте-Карло, метод моделирования отжига, генетические и эволюционные алгоритмы.

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

Алгоритм 1. Алгоритм для выполнения локальной оптимизации функции затрат.

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

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

где:

S является функцией затрат,

M является вектором N параметров (m₁, m₂,..., m_N), описывающим модель геологической среды,

g является индексом сейсмограммы,

w_g является функцией источника для сейсмограммы g, которая является функцией от пространственных координат и времени, для точечного источника она является дельта-функцией от пространственных координат,

N_g является числом сейсмограмм,

r является индексом приемника в сейсмограмме,

N_r является числом приемников в сейсмограмме,

t является индексом временной выборки в трассе,

N_t является числом временных выборок,

W является функцией критериев минимизации (обычно выбирается W(x)=x², которая представляет собой критерии по методу наименьших квадратов (L2)),

ψ_calc является расчетными данными сейсмического давления из модели M,

ψ_obs является измеренными данными сейсмического давления.

Сейсмограммы могут быть любым типом сейсмограммы, которая может быть моделирована в одном проходе программы сейсмического прямого моделирования. Обычно сейсмограммы соответствуют сейсмическому взрыву, хотя взрывы могут быть более общими, чем точечные источники. Для точечных источников индекс g сейсмограммы соответствует местоположению отдельных точечных источников. Для плоских сейсмоволн источники g должны соответствовать различным направлениям распространения плоских сейсмоволн. Эти данные обобщенных источников, ψ_obs, могут или быть обнаружены в полевых условиях или могут быть синтезированы из данных, обнаруженных с использованием точечных источников. Расчетные данные ψ_calc, с другой стороны, обычно могут вычисляться непосредственно посредством использования функции обобщенных источников при прямом моделировании. Для многих типов прямого моделирования, включающих в себя моделирование на основе конечных разностей, время вычисления, необходимое для обобщенного источника, примерно равно времени вычисления, необходимому для точечного источника.

Уравнение (1) может быть упрощено до:

причем теперь подразумевается сумма по приемникам и временным выборкам, и

Инверсия пытается обновлять модель M таким образом, что S(M) является минимумом. Это может быть выполнено посредством локальной оптимизации функции затрат, которая обновляет данную модель M^(k) следующим образом:

где k является номером итерации, α является скалярным размером обновления модели, и ∇_MS(M) является градиентом функции несоответствия, принятой относительно параметров модели. Возмущения модели или значения, посредством которых обновляется модель, вычисляются посредством умножения градиента целевой функции на длину α шага, которая должна вычисляться многократно.

Из уравнения (2) следующее уравнение может извлекаться для градиента функции затрат:

Таким образом, чтобы вычислять градиент функции затрат, следует отдельно вычислять градиент доли каждой сейсмограммы в функции затрат, затем суммировать эти доли. Следовательно, трудоемкость вычислений, требуемая для вычисления ∇_MS(M), в N_g раз превышает трудоемкость вычислений, требуемую для того, чтобы определять долю одной сейсмограммы в градиенте. Для геофизических задач, N_g обычно соответствует числу геофизических источников и составляет порядка 10000-100000, значительно повышая затраты на вычисление ∇_MS(M).

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

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

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

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

Уменьшение вычислительных затрат, полученное посредством инвертирования комбинированных источников, по меньшей мере, частично смещается вследствие того факта, что инверсия комбинированных данных обычно формирует менее точную инвертированную модель. Эти потери в точности обусловлены тем фактом, что информация теряется, когда суммируются отдельные источники, и, следовательно, суммированные данные не ограничивают инвертированную модель так строго, как несуммированные данные. Эти потери информации во время суммирования могут быть минимизированы посредством кодирования каждой записи по взрыву перед суммированием. Кодирование перед комбинированием сохраняет значительно больше информации в данных одновременных источников, и, следовательно, лучше ограничивает инверсию [4]. Кодирование также дает возможность комбинирования близко расположенных источников, тем самым давая возможность комбинирования большего числа источников для данной вычислительной области. Для этой технологии могут быть использованы различные схемы кодирования, включающие в себя кодирование на основе сдвига по времени и случайное фазовое кодирование. Оставшаяся часть этого раздела "Уровень техники" кратко анализирует различные опубликованные технологии на основе геофизических одновременных источников, как кодированных, так и некодированных.

Van Manen [6] предлагает использование способа сейсмической интерферометрии для того, чтобы ускорять прямое моделирование. Сейсмическая интерферометрия работает посредством размещения источников по всей границе интересующей области. Эти источники моделируются по отдельности, и записывается волновое поле во всех местоположениях, для которых требуется функция Грина. Функция Грина между любыми двумя записанными местоположениями затем может быть вычислена посредством взаимной корреляции трасс, обнаруженных в двух записанных местоположениях, и суммирования по всем граничным источникам. Если данные, которые должны быть инвертированы, имеют большое число источников и приемников, которые находятся в интересующей области (в противоположность наличию одних или других на границе), это представляет собой очень эффективный метод для вычисления требуемых функций Грина. Тем не менее для случая сейсмических данных редко, когда как источник, так и приемник для данных, которые должны быть инвертированы, находятся в интересующей области. Следовательно, это улучшение имеет очень ограниченную применимость для задачи сейсмической инверсии.

Berkhout [7] и Zhang [8] предлагают то, что инверсия, в общем, может быть улучшена посредством инвертирования некодированных одновременных источников, которые суммируются когерентно, чтобы формировать некоторый требуемый фронт волны в некоторой области геологической среды. Например, данные точечных источников могут быть суммированы со сдвигами по времени, которые являются линейной функцией от местоположения источников, чтобы формировать нисходящую плоскую сейсмоволну под некоторым конкретным углом относительно поверхности. Эта технология может применяться ко всем категориям инверсии. Проблема этого способа заключается в том, что когерентное суммирование сейсмограмм источников обязательно уменьшает объем информации в данных. Таким образом, например, суммирование для того, чтобы формировать плоскую сейсмоволну, удаляет всю информацию в сейсмических данных, связанную со временем пробега в зависимости от выноса источников от приемников. Эта информация является критически важной для обновления медленно варьирующейся фоновой модели скорости, и, следовательно, способ Berkhout не подчиняется оптимальным ограничениям. Чтобы преодолевать эту проблему, могут инвертироваться множество различных когерентных сумм данных (например, множество плоских сейсмоволн с различными направлениями распространения), но в таком случае теряется эффективность, поскольку затраты на инверсию являются пропорциональными числу различных инвертированных сумм. В данном документе, такие когерентно суммированные источники называются обобщенными источниками. Следовательно, обобщенный источник может быть либо точечным источником, либо суммой точечных источников, которая формирует фронт волны некоторой требуемой формы.

Van Riel [9] предлагает инверсию посредством некодированного накопления или частичного накопления (относительно выноса источников от приемников) входных сейсмических данных, а затем задания функции затрат относительно этих накопленных данных, которые должны быть оптимизированы. Таким образом, эта публикация предлагает улучшение инверсии на основе функции затрат с использованием некодированных одновременных источников. Как и для способа инверсии одновременных источников по Berkhout [6], накопление, предлагаемое посредством этого способа, уменьшает объем информации в данных, которые должны быть инвертированы, и, следовательно, инверсия подчиняется субоптимальным ограничениям по сравнению с исходными данными.

Mora [10] предлагает инвертирование данных, которые являются суммой широко разнесенных источников. Таким образом, эта публикация предлагает повышение эффективности инверсии с использованием моделирования на основе некодированных одновременных источников. Суммирование широко разнесенных источников обеспечивает преимущество сохранения гораздо большего объема информации, чем когерентная сумма, предложенная Berkhout. Тем не менее суммирование широко разнесенных источников подразумевает, что апертура (инвертируемая область модели), которая должна быть использована в инверсии, должна быть увеличена, чтобы приспосабливать все широко разнесенные источники. Поскольку время вычисления является пропорциональным площади этой апертуры, способ Mora не формирует такой прирост эффективности, который может достигаться, если суммированные источники находятся рядом друг с другом.

Ober [11] предлагает ускорение сейсмической миграции, частного случая неитеративной инверсии, посредством использования одновременных кодированных источников. После тестирования различных способов кодирования Ober выяснил, что результирующие мигрированные изображения значительно уменьшают отношение "сигнал-шум" вследствие того факта, что функции широкополосного кодирования обязательно являются только приблизительно ортогональными. Таким образом, при суммировании более 16 взрывов, качество инверсии не является удовлетворительным. Поскольку неитеративная инверсия не является слишком затратной для начального уровня, и поскольку требуется инверсия высокого отношения "сигнал-шум", эта технология не используется широко на практике в геофизической отрасли.

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

Кодирование источников, предложенное Krebs и др. в публикации PCT заявки на патент номер WO 2008/042081, которая содержится в данном документе по ссылке во всех разрешенных юрисдикциях, представляет собой очень эффективный с точки зрения затрат способ для того, чтобы инвертировать данные полного волнового поля. (Идентичный подход одновременной инверсии кодированной сейсмограммы должен работать для приемников либо через обратимость источников-приемников, либо посредством кодирования фактических местоположений приемников в сейсмограммах общих источников данных). Для стационарных приемников прямые и сопряженные вычисления должны выполняться только для одного эффективного источника; см. публикацию PCT заявки на патент номер WO 2009/117174, которая содержится в данном документе по ссылке во всех разрешенных юрисдикциях. Учитывая тот факт, что сотни взрывов записываются для типичных геометрий двумерного обнаружения и тысячи в случае трехмерных съемок, экономия вычислительных ресурсов от этого способа является довольно значительной. На практике предположение о стационарных приемниках не является строго достоверным для наиболее распространенных геометрий обнаружения полевых данных. В случае данных морского сейсморазведочного кабеля, источники и приемники перемещаются для каждого нового взрыва. Даже при съемках, в которых местоположения приемников являются стационарными, практика зачастую такова, что не все приемники "прослушивают" каждый взрыв, и приемники, которые прослушивают, могут варьироваться в зависимости от взрыва. Это также нарушает "предположение о стационарных приемниках". Помимо этого вследствие логистических проблем, трудно записывать данные близко к источнику, и это означает то, что данные минимальных выносов типично отсутствуют. Это является истинным и для морских, и для наземных съемок. Оба из этих факторов означают, что для сейсмограммы одновременных источников, каждое местоположение приемника представляет собой отсутствующие данные для некоторых исходных взрывов. В общих словах, в инверсии одновременных кодированных источников, для данной одновременной кодированной сейсмограммы, данные требуются во всех местоположениях приемников для каждого взрыва, и это может упоминаться как предположение о стационарных приемниках инверсии одновременных кодированных источников. В WO 08/042081, некоторые раскрытые варианты осуществления могут работать лучше, чем другие, когда не удовлетворяется предположение о стационарных приемниках. Следовательно, должно быть преимущественным иметь приспособление или регулирование для прямого применения инверсии одновременных кодированных источников (и/или приемников), что должно повышать производительность, когда скомпрометировано предположение о стационарных приемниках. Настоящее изобретение предоставляет способы для осуществления этого. Haber и др. [15] также описывает подход для задачи перемещения приемников в инверсии одновременных кодированных источников с использованием способа стохастической оптимизации и применяет его к задаче удельного сопротивления по постоянному току.

Сущность изобретения

В одном варианте осуществления изобретение представляет собой машинореализованный способ для инверсии полного волнового поля измеренных геофизических данных, чтобы определять модель физических свойств для области геологической среды, содержащий: (a) использование компьютера для того, чтобы инвертировать выбранное околоповерхностное временное окно поступлений из измеренных геофизических данных посредством инверсии одновременных кодированных источников и/или приемников, чтобы получать первую модель физических свойств для области геологической среды; (b) использование компьютера для того, чтобы инвертировать измеренные геофизические данные или выбранное глубокое временное окно поступлений из измеренных геофизических данных посредством итеративной инверсии последовательных источников, которая может использовать только разреженную дискретизацию измеренных данных, чтобы получать вторую модель физических свойств для области геологической среды, при этом первая модель физических свойств используется в качестве начальной модели, и набор местоположений источников используется для того, чтобы обновлять вторую модель физических свойств в итеративной инверсии последовательных источников; и (c) вывод, отображение или сохранение в устройстве хранения данных второй модели физических свойств области геологической среды.

В другом варианте осуществления изобретение представляет собой машинореализованный способ для инверсии полного волнового поля измеренных геофизических данных, чтобы определять модель физических свойств для области геологической среды, содержащий: (a) использование компьютера для того, чтобы инвертировать выбранное околоповерхностное временное окно поступлений из измеренных геофизических данных посредством инверсии одновременных кодированных источников и/или приемников, чтобы получать первую модель физических свойств для области геологической среды; (b) использование первой модели физических свойств для того, чтобы моделировать, с использованием компьютера, синтетические данные для более длинных выносов, соответствующих поступлениям из места глубже упомянутого околоповерхностного временного окна; (c) использование компьютера для того, чтобы инвертировать измеренные геофизические данные, при этом данные с более длинными выносами дополняются, причем упомянутая инверсия является инверсией одновременных кодированных источников и/или кодированных приемников, чтобы получать вторую модель физических свойств области геологической среды, при этом упомянутые дополненные данные с более длинными выносами являются суммой синтетических данных для более длинных выносов и измеренных данных при более длинных выносах; и (d) вывод, отображение или сохранение в устройстве хранения данных второй модели физических свойств области геологической среды.

Краткое описание чертежей

Вследствие ограничений на основе правил выдачи патентов на использование цвета, прилагаемые чертежи являются черно-белыми репродукциями цветных оригиналов. Копии этой заявки на патент или публикации с цветными чертежами могут быть получены из Патентного ведомства США (заявка на патент США № 12/903744) после запроса и выплаты необходимых платежей.

Настоящее изобретение и его преимущества должны лучше пониматься посредством рассмотрения нижеприведенного подробного описания и прилагаемых чертежей, на которых:

Фиг.1 является принципиальной схемой, показывающей окно данных, которое может использоваться для инверсии одновременных источников;

Фиг.2 является блок-схемой последовательности операций способа, показывающей базовые этапы в одном варианте осуществления настоящего изобретаемого способа, в котором кодирование одновременных источников используется в околоповерхностном временном окне, и разреженная инверсия последовательных источников используется для более глубоких окон в подходе I и в подходе II, околоповерхностная модель используется для того, чтобы вычислять отсутствующие трассы данных, которые сначала кодируются, чтобы выполнять инверсию одновременных источников;

Фиг.3 показывает "истинную" модель скорости, используемую в примере, чтобы формировать измеренные данные;

Фиг.4 показывает две сейсмограммы по взрывам данных последовательных источников, вычисленные с использованием модели скорости на фиг.3 и кодированные с использованием бинарного кодирования, описанного в WO 2008/042081, указывающего околоповерхностное временное окно 0-3 секунд;

Фиг.5 показывает начальную модель скорости для инверсии полного волнового поля в тестовом примере;

Фиг.6 показывает модель, полученную из инвертирования одновременных источников посредством инверсии низкочастотных данных (пиковая частота 7 Гц) от околоповерхностного окна (0-3 секунды); как показано на фиг.4, предположение о стационарных приемниках является достоверным для временного окна в 0-3 секунды;

Фиг.7 показывает модель, полученную с использованием инверсии одновременных источников низкочастотных данных (пиковая частота 7 Гц) с окном данных от 0-4 секунд; поскольку предположение о стационарных приемниках более не является достоверным, артефакты могут быть четко видны в инвертированной модели;

Фиг.8 показывает модель, полученную с использованием инверсии одновременных источников низкочастотных данных (пиковая частота 7 Гц) с окном данных от 0-5 секунд; поскольку предположение о стационарных приемниках более не является достоверным, артефакты могут быть четко видны в инвертированной модели;

Фиг.9 показывает модель, полученную в тестовом примере посредством использования инверсии последовательных источников низкочастотных данных (пиковая частота 7 Гц) с окном данных от 0-4 секунд с использованием разреженного набора источников с разделением между взрывами в 1,2 км; поскольку используются последовательные источники, модель не имеет артефактов по сравнению с фиг.7;

Фиг.10 показывает модель, полученную в тестовом примере с использованием инверсии последовательных источников низкочастотных данных (пиковая частота 7 Гц) с окном данных от 0-5 секунд с использованием разреженного набора источников с разделением между взрывами в 1,2 км; поскольку используются последовательные источники, модель не имеет артефактов по сравнению с фиг.8;

Фиг.11 показывает модель, полученную в тестовом примере с использованием инверсии последовательных источников вседиапазонных частотных данных (пиковая частота 40 Гц) с окном данных от 0-5 секунд с использованием разреженного набора источников с разделением между взрывами в 1,2 км; поскольку используются последовательные источники, модель не имеет артефактов по сравнению с фиг.8;

Фиг.12 показывает ограниченную по минимальному и максимальному выносу сейсмограмм по взрывам, вычисленную в тестовом примере с использованием истинной модели на фиг.3; самый минимальный доступный вынос составляет 200 м, и самый максимальный доступный вынос составляет 5000 м;

Фиг.13 показывает ограниченную по минимальному выносу сейсмограмм по взрывам на фиг.12, заполненную с трассой из самого минимального доступного выноса в x=200;

Фиг.14 иллюстрирует инверсию с использованием кодированных данных с ограниченными по минимальному выносу сейсмограммами по взрывам, показанными на фиг.12, при этом модель, показанная на итерации 5, получена с помощью инверсии одновременных источников с низкочастотными данными (пиковая частота 7 Гц);

Фиг.15 иллюстрирует инверсию с использованием кодированных данных с ограниченными по минимальному выносу сейсмограммами по взрывам, показанными на фиг.12, при этом модель, показанная на итерации 10, получена с помощью инверсии одновременных источников с низкочастотными данными (пиковая частота 7 Гц); он четко показывает то, что модель расходится со значительными артефактами;

Фиг.16 иллюстрирует инверсию с использованием кодированных данных с сейсмограммами по взрывам с заполнением минимальных выносов, показанными на фиг.13 (поскольку заполнение минимальных выносов аппроксимировано, градиент в водном слое подавляется, чтобы уменьшать воздействие аппроксимации на инверсию одновременных источников), при этом модель, показанная на итерации 5, получена с помощью инверсии одновременных источников с низкочастотными данными (пиковая частота 7 Гц); и

Фиг.17 иллюстрирует инверсию с использованием кодированных данных с сейсмограммами по взрывам с заполнением минимальных выносов, показанными на фиг.13 (поскольку заполнение минимальных выносов аппроксимировано, градиент в водном слое подавляется, чтобы уменьшать воздействие аппроксимации на инверсию одновременных источников), при этом модель, показанная на итерации 200, получена с помощью инверсии одновременных источников с низкочастотными данными (пиковая частота 7 Гц).

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

Подробное описание примерных вариантов осуществления

Один вариант осуществления настоящего изобретения представляет собой гибридное комбинирование инверсии одновременных кодированных источников с традиционной инверсией последовательных источников. Этот вариант осуществления использует кодирование одновременных источников для околоповерхностного временного окна данных и использует разреженные последовательные источники, чтобы инвертировать более глубокую часть данных. Krebs и др. [5, 16] показывают, что функция затрат кодированных одновременных источников может быть вычислена более эффективно, чем традиционные функции затрат, при одновременном предоставлении точных инверсий. Функция затрат одновременных источников задается здесь следующим образом (следует сравнить с вышеприведенным уравнением (2)):

причем подразумевается суммирование по приемникам и временным выборкам, аналогично уравнению (2), и:

задает сумму по сейсмограммам посредством подгрупп сейсмограмм,

S_sim является функцией затрат для данных одновременных источников,

G является группами одновременных обобщенных источников, и

N_G является числом групп,

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

Внешнее суммирование в уравнении (6) выполняется по группам одновременных обобщенных источников, соответствующих типу сейсмограммы (например, точечных источников для сейсмограмм по общим взрывам). Внутреннее суммирование, g, выполняется по сейсмограммам, которые группируются для одновременного вычисления. Для некоторых способов прямого моделирования, таких как моделирование на основе конечных разностей, вычисление прямой модели для суммированных обобщенных источников (внутренняя сумма g∈G) может быть выполнено за идентичное количество времени, что и вычисление для одного источника. Следовательно, как показано в Krebs и др. [5], δ(M,Σc_g⊗w_g) может быть вычислено очень эффективно с использованием алгоритма 3.

Алгоритм 3. Алгоритм для вычисления функции затрат кодированных одновременных источников.

Как также показано в Krebs и др. [5], этот алгоритм может вычислять S_sim(M) на коэффициент N_g/N_G раз быстрее, чем S(M) из уравнения (2).

В геометрии стационарных приемников сейсмограммы по взрывам имеют полное покрытие по приемникам. Для еще более околоповерхностного окна это может достигаться, если могут быть заполнены минимальные выносы. Еще более околоповерхностные окна задаются как соответствующие временам, меньшим, чем времена поступления для более быстрых режимов при выносах наибольшей длины. Размер околоповерхностного окна, в котором является применимым подход на основе кодированных одновременных источников, зависит от приповерхностной скорости и максимальных доступных выносов. Фиг.1 является принципиальной схемой, показывающей окно данных, которое может использоваться для инверсии одновременных источников. Линии с меньшими наклонами задают нижнюю часть окна для случая (100) большей околоповерхностной скорости. Линии с большими наклонами задают нижнюю часть окна для случая (101) меньшей околоповерхностной скорости. Вертикальная ось представляет собой время, а горизонтальная ось представляет собой вынос источников от приемников, причем нулевой вынос для выбранного приемника находится в середине схемы, при этом положительные выносы задаются вправо, а отрицательные выносы - влево. Схематический вид, показанный на фиг.1, предлагает, что для небольших приповерхностных скоростей применимое окно (101) является более длинным по сравнению со случаем, когда приповерхностные скорости являются большими (100).

Отрицательные выносы в сейсмограммах по взрывам могут быть заполнены с использованием обратимости (201). Отсутствующие минимальные выносы могут быть заполнены оцененными значениями, например, из соседних трасс (201). Заполнение минимальных выносов в околоповерхностном окне предоставляет полное покрытие приемников, чтобы использовать кодирование одновременных источников. Другая альтернатива заключается в том, чтобы просто исключать все приемники, которые представляют минимальный вынос для некоторого взрыва. Существуют различные подходы предшествующего уровня техники для того, чтобы заполнять данные минимальных выносов, к примеру, способы восстановления, типично используемые в связанном с поверхностью множественном ослаблении. Это обеспечивает соответствие сейсмограмм (202) по взрывам с околоповерхностным окном предположению о стационарных приемниках, предпочтительному для кодирования (203) одновременных источников. Следовательно, для околоповерхностного окна можно формировать кодированные данные из одновременных источников (204) и инвертировать (205) для модели (206) геологической среды, которая соответствует этим данным. Поскольку время вычисления значительно меньше на несколько порядков величины, более сложный алгоритм моделирования (например, эластичная инверсия полного волнового поля) потенциально может быть использован при необходимости. После инверсии одновременных источников с помощью данных околоповерхностного окна следующий этап представляет собой использование околоповерхностной модели в качестве начальной модели (207), чтобы инвертировать данные (208) более глубоких окон, поскольку предположение о стационарных приемниках нарушается для более глубоких окон.

Настоящий изобретаемый способ включает в себя два варианта осуществления для инвертирования данных более глубоких окон. Блок-схема последовательности операций способа этих гибридных подходов представляется на фиг.2. В первом подходе, показанном на фиг.2, традиционные технологии последовательных взрывов используются для того, чтобы инвертировать (209) данные (208) более глубоких окон. Тем не менее более разреженный набор источников является достаточным, чтобы обновлять более глубокую часть модели (210) и, следовательно, преимущественно может быть использован. То, сколько источников требуется, зависит от длины пространственной волны модели и максимальных доступных выносов. С точки зрения физики разрешение геологической среды ухудшается с глубиной. Например, в ограниченных по выносу данных доля в более глубоком окне преобладает из отражений, а не из передачи. Таким образом, в зависимости от требуемой длины шкалы модели, которая должна быть обновлена, выбор того, сколько последовательных источников необходимо, может быть зависимым от задачи. Например, если цель состоит в том, чтобы только получать сглаженную модель скорости, которая лучше объясняет кинематику (времена пробега событий), то пользователь изобретения может выбирать разреженный набор источников для инверсии. Если цель состоит в том, чтобы получать целеориентированное обновление модели высокого разрешения, то пользователь может включать больше источников выше целевой области и разреженных источников в любом другом месте.

Во втором подходе (подходе II на фиг.2), преимущество обусловлено тем фактом, что околоповерхностная модель (206), являющаяся результатом инверсии околоповерхностных окон, должна быть удовлетворительной для прогнозирования данных (211) длинных выносов во времена, превышающие околоповерхностное окно. Следовательно, околоповерхностная модель, полученная из одновременной инверсии, используется для того, чтобы вычислять данные при больших отсутствующих выносах (211). Преимущество этого подхода состоит в том, что с доступными трассами с максимальным выносом, полученными с использованием прямого моделирования, может применяться предположение о стационарных приемниках. Таким образом, с доступными данными вместе с расчетными максимальными выносами, сейсмограммы по взрывам могут быть кодированы (212) для инверсии (214) более глубоких окон (213), как указано на фиг.2.

Гибридный подход настоящего изобретения является применимым не только к данным сейсморазведочного кабеля, но также и к множеству других геометрий обнаружения, в которых нарушается предположение о стационарных приемниках. Например, обнаружения на основе морского донного кабеля (OBC) типично являются коммутируемыми (стационарные приемники для поднабора взрывов), и они не соответствуют идеализированному случаю, когда все приемники являются стационарными и ведут запись для всех взрывов. Аналогично для наземных обнаружений, вследствие логистических проблем трудно достигать геометрии стационарных приемников.

В уравнении (6) может быть использовано множество типов функций c_g кодирования, включающее в себя, но не только:

- линейное, случайное, на основе линейной частотной модуляции и на основе модифицированной линейной частотной модуляции частотно-зависимое фазовое кодирование, как представлено в Romero и др. [13];

- частотно-независимое фазовое кодирование, как представлено в Jing и др. [14];

- кодирование на основе случайного сдвига по времени;

- мультиплексирование с частотным разделением каналов (FDMA), мультиплексирование с временным разделением каналов (TDMA) и мультиплексирование с кодовым разделением каналов (CDMA), используемые в связи.

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

Следует отметить, что технология на основе одновременных кодированных источников может использоваться для многих типов функций затрат при инверсии. В частности, она может использоваться для функций затрат на основе норм, отличных от L2, поясненной выше. Она также может быть использована в более сложных функциях затрат по сравнению с функцией затрат, представленной в уравнении 2, включающих в себя регуляризованные функции затрат. В завершение, способ одновременных кодированных источников может быть использован с любым типом глобального или локального способа инверсии функции затрат, включающим в себя метод Монте-Карло, метод моделирования отжига, генетический алгоритм, эволюционный алгоритм, градиентный линейный поиск, сопряженные градиенты и метод Ньютона.

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

Некоторые вариации варианта осуществления, описанного выше, включают в себя:

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

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

- Это изобретение не ограничено однокомпонентными точечными приемниками. Например, приемники могут быть матрицами приемников, или они могут быть многокомпонентными приемниками.

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

- Обнаружение данных одновременных кодированных источников может приводить к значительному снижению затрат на обнаружение геофизических данных.

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

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

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

Пример

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

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

Фиг.4 показывает пример двух ограниченных по максимальному выносу характерных сейсмограмм по взрывам, кодированных с использованием бинарного кодирования, описанного в публикации PCT заявки на патент номер WO 2008/042081. Ограниченная по выносу природа сейсмограмм по взрывам, которая является очевидной из чертежа, приводит к недостоверности предположения о стационарных приемниках. Тем не менее также можно видеть то, что предположение о стационарных приемниках является достоверным для околоповерхностного временного окна от 0 до приблизительно 3 секунд. Фиг.4 показывает то, что приемники для левого взрыва между x=800 и x=2000 не имеют доли от правого взрыва, т.е. энергия из правого взрыва не достигает этих приемников ранее 3 секунд. Аналогично, приемники для правого взрыва между x=2800 и x=4200 не имеют доли от левого взрыва в околоповерхностном временном окне между 0 и 3 секундами. Самый максимальный вынос, используемый в этом примере, составляет 5 км. Данные формируются с использованием разнесения между взрывами 80 м, и разделение приемников составляет 10 м. Этот пример использует данные минимальных выносов, вычисленные с использованием модели. Чтобы формировать начальную модель скорости для инверсии, сейсмограммы на основе общей глубинной точки (CDP) формируются с ограниченными по выносу в 5 км последовательными взрывами, и глубинная миграция Кирхгофа используется для того, чтобы выравнивать сейсмограммы по общим изображениям. Начальная модель скорости, показанная на фиг.5, получается с использованием глубинной миграции Кирхгофа. Поскольку глубинная миграция соответствует кинематической части данных, начальная модель является сглаженной моделью.

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

Фиг.6 показывает модель, полученную посредством инвертирования низкочастотных данных (с пиковой частотой 7 Гц). Поскольку предположение о стационарных приемниках не нарушается до 3 секунд данных (показано на фиг.4), восстановленная модель достигает стабильной сходимости с помощью инверсии одновременных источников. Для данных с временным окном до 3 сек максимальный вынос составляет приблизительно 5 км. Таким образом, инверсия восстанавливает еще более околоповерхностную часть модели хорошо по сравнению с более глубокой частью. Фиг.7 и 8 показывают результат инвертирования данных с временным окном в 4 сек и 5 сек, соответственно. Из результатов очевидно, что инвертированная модель имеет артефакты, поскольку кодированные данные в этих временных окнах не имеют данных больших выносов, и как результат, вычисляется некорректный градиент для обновления модели. Поскольку предположение о стационарных приемниках нарушается в этих временных окнах, инверсия одновременных источников формирует артефакты в модели.

Модель, полученная из инверсии одновременных источников околоповерхностного окна (фиг.5), используется в качестве начальной модели для инверсии последовательных источников для более глубокого временного окна в 4 сек. Разреженный набор источников (20 источников из 383 источников, используемых для кодирования данных) используется для инверсии окон в 4 сек. Результирующая модель показывается на фиг.9. С использованием модели из фиг.9 в качестве начальной модели, инвертируются данные с окном в 5 сек. Для инверсии данных в 5 сек используются только 10 последовательных источников. Результирующая модель показывается на фиг.10. Для ограниченных по выносу данных апертура становится ограниченной по глубине, следовательно, разреженный набор источников является достаточным, чтобы получать обоснованную модель. Использование разреженных источников предоставляет значительное снижение трудоемкости вычислений.

Фиг.11 показывает результат инвертирования окна данных в 5 сек с вседиапазонными данными (пиковая частота 40 Гц) для 10 последовательных источников.

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

Фиг.12 показывает характерную сейсмограмму по взрыву с отсутствующими выносами, которые кодируются для одновременной инверсии. Поскольку существует отсутствующая информация минимальных выносов в кодированных данных, не только одновременная инверсия формирует артефакты в восстановленной модели, но также трудно обеспечивать соответствие данных. Восстановленная модель на итерации 5 и 50 показывается на фиг.14 и 15, соответственно. Инвертированная модель имеет значительные артефакты, и решение полностью расходится для инверсии низкочастотных данных. Фиг.15 показывает то, что вследствие отсутствующей информации минимальных выносов в кодированных данных существуют большие артефакты скорости в водном слое. Чтобы улучшать решение, интерполяция или регуляризация данных может быть использована для того, чтобы заполнять трассы с минимальным выносом. Это является общей процедурой для подготовки данных для задач множественного исключения. Другой фактор, который приводит к артефактам в инверсии, представляет собой отсутствующие данные по большой энергии вследствие поступлений прямой волны. Типично данные минимальных выносов с поступлениями прямой волны имеют значительную энергию, следовательно, с точки зрения обеспечения соответствия данных они имеют большую долю в несоответствии данных. Отсутствующая информация в измеренных кодированных данных обостряет проблему и приводит к тому, что инверсия расходится, что вызывает значительные артефакты. При признании этого аспекта обеспечения соответствия данных, проблема может в определенной степени обходиться посредством маскирования или подавления градиента в околоповерхностной части модели, которая обычно ассоциирована со скоростью потока воды. Это демонстрируется в следующем примере.

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

Фиг.16 и 17 показывают результат инверсии для итерации 5 и 200, соответственно, с использованием инверсии одновременных источников низкочастотных данных с подавлением градиента в водном слое. Фиг.17 четко показывает то, что инверсия является стабильной, и восстановленная модель является аналогичной инверсии на фиг.5. Результаты должны дополнительно улучшаться, если минимальные выносы интерполируются с использованием способов регуляризации данных, таких как интерполяция Радона. Все представленные примеры демонстрируют подход "I", описанный на фиг.2.

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

Библиографический список

1. Tarantola, A., "Inversion of seismic reflection data in the acoustic approximation", Geophysics 49, 1259-1266 (1984 год).

2. Sirgue, L. и Pratt G. "Efficient waveform inversion and imaging: A strategy for selecting temporal frequencies", Geophysics 69, 231-248 (2004 год).

3. Fallat, M. R., Dosso, S. E., "Geoacoustic inversion via local, global and hybrid algorithms", Journal of the Acoustical Society of America 105, 3219-3230 (1999 год).

4. Hinkley, D. и Krebs, J., "Gradient computation for simultaneous source inversion", публикация PCT заявки на патент номер WO 2009/117174.

5. Krebs, J. R., Anderson, J. A., Neelamani, R., Hinkley, D., Jing, C., Dickens, T., Krohn, C., Traynin, P., "Iterative inversion of data from simultaneous geophysical sources", публикация PCT заявки на патент номер WO 2008/042081.

6. Van Manen, D. J., Robertsson, J.O.A., Curtis, A., "Making wave by time reversal", SEG International Exposition and 75^th Annual Meeting Expanded Abstracts, 1763-1766 (2005 год).

7. Berkhout, A. J., "Areal shot record technology", Journal of Seismic Exploration 1, 251-264 (1992 год).

8. Zhang, Y., Sun, J., Notfors, C., Gray, S. H., Cherris, L., Young, J., "Delayed-shot 3D depth migration", Geophysics 70, E21-E28 (2005 год).

9. Van Riel, P. и Hendrik, W. J. D., "Method of estimating elastic and compositional parameters from seismic and echo-acoustic data", Патент США номер 6876928 (2005 год).

10. Mora, P., "Nonlinear two-dimensional elastic inversion of multi-offset seismic data", Geophysics 52, 1211-1228 (1987 год).

11. Ober, C. C., Romero, L. A., Ghiglia, D. C., "Method of Migrating Seismic Records", Патент США номер 6021094 (2000 год).

12. Ikelle, L. T., "Multi-shooting approach to seismic modeling and acquisition", Патент США номер 6327537 (2001 год).

13. Romero, L. A., Ghiglia, D. C., Ober, C. C., Morton, S. A., "Phase encoding of shot records in prestack migration", Geophysics 65, 426-436 (2000 год).

14. Jing X., Finn, C. J., Dickens, T. A., Willen, D. E., "Encoding multiple shot gathers in prestack migration", SEG International Exposition and 70^th Annual Meeting Expanded Abstracts, 786-789 (2000 год).

15. Haber, E., Chung M. и Herrmann, "An effective method for parameter estimation with PDE constraints with multiple right hand sides", Preprint - UBC http://www.math.ubc.ca/~haber/pubs/PdeOptStochV5.pdf (2010 год).

16. Jerome R. Krebs, John E. Anderson, David Hinkley, Ramesh Neelamani, Sunwoong Lee, Anatoly Baumstein и Martin-Daniel Lacasse, "Full-wavefield seismic inversion using encoded sources", Geophysics 74-6, WCC177-WCC188 (2009 год).