CC BY
13УДК 519.6:550.34.013.06
АЛГОРИТМЫ ОПРЕДЕЛЕНИЯ УПРУГИХ ПАРАМЕТРОВ ПО ПЛОЩАДНЫМ СИСТЕМАМ НАБЛЮДЕНИЙ (ПРЯМАЯ ЛИНЕЙНАЯ ОБРАБОТКА ДАННЫХ СЕЙСМИЧЕСКИХ
НАБЛЮДЕНИЙ)
С. И. Кабанихин1'2,3, М. А. Шишленин1,2,3, Н. С. Новиков1,3
1 Институт вычислительной математики и математической геофизики СО РАН 2 Институт математики им. С. Л. Соболева 3 Новосибирский государственный университет [email protected], [email protected], [email protected]
В статье предложены алгоритмы определения упругих параметров по площадным системам наблюдений. Рассматривается прямой метод решения коэффициентных обратных задач на основе многомерных аналогов уравнений И. М. Гельфанда, Б. М. Левитана и М. Г. Крейна. Метод заключается в сведении многомерной нелинейной обратной задачи к однопараметрическому семейству линейных интегральных уравнений Фредгольма.
Ключевые слова: уравнение акустики, многомерный аналог уравнения М. Г. Крейна, коэффициентная обратная задача, численные методы.
ELASTIC PARAMETERS ESTIMATION ALGORITHMS FOR AREA OBSERVATION SYSTEMS (DIRECT LINEAR DATA PROCESSING OF SEISMIC OBSERVATIONS)
S. I. Kabanikhin1,2 3, M. A. Shishlenin12 3, N. S. Novikov13
1 Institute of Computational Mathematics and Mathematical Geophysics 2 Sobolev Institute of Mathematics 3 Novosibirsk State University [email protected], [email protected], [email protected]
The authors propose algorithms for estimating elastic parameters in area observation systems. The direct method for solving inverse coefficient problems with the multidimensional analogs of I. M. Gelfand, B. M. Levitan and M. G. Krein equations is considered. The method includes the reduction of a multidimensional nonlinear inverse problem to a one-parameter family of linear integral Fredholm equations.
Keywords: acoustic equation, multidimensional analog of Krein equation, inverse coefficient problem, numerical methods.
Введение
В настоящее время, благодаря площадным системам наблюдений, удалось создать принципиально новый метод решения трёхмерных обратных задач, в котором используются: трёхмерные аналоги уравнений Гельфанда—Левитана—Крейна, параллельные вычисления на высокопроизводительных кластерах, методы Монте-Карло, супербыстрые алгоритмы обращения блочно-теплицевых матриц больших размерностей.
Основной проблемой моделирования трёхмерных упругих сред является большой объём вычислений. Даже для сравнительно небольшой области 8 км3 решение прямой задачи сейсморазведки является сложной проблемой, а если учесть, что большинство современных методов решения обратных задач являются итерационными, то даже количество операций, требуемых для проведения нескольких итераций, может привести к большим ошибкам.
Сейсмические волны отражаются и рассеиваются на неоднородностях и приносят на поверхность Земли информацию об искомых объектах.
В данной статье мы изложим алгоритм, который позволяет избежать многократного решения прямых задач. Вторым достоинством алгоритма является сведение нелинейной обратной задачи к системе линейных интегральных уравнений (многомерный аналог уравнений Гельфанда—Левитана— Крейна). Предлагаемый подход использует многократное возбуждение среды и площадную систему записи сейсмограмм.
Поверхность Земли, под которой находятся искомые неоднородности, покрывается дискретной сеткой размера N х N .В каждой точке сетки устанавливается сейсмоприёмник. Далее в каждой из точек сетки производится возбуждение, а во всех остальных точках сетки записывается отклик среды на это возбуждение. Таким образом производится N возбуждений среды и накапливается N — 1)2 сейсмограмм.
У
В каждой точке сетки с координатами k = (k\,k2) устанавливается приёмник акустических колебаний. Далее в каждой из точек сетки производится возбуждение g(k)(t). Во всех остальных точках сетки записывается отклик среды f (k^(t) на возбуждение g(k)(t). Отметим, что после применения деконволюции, функции g(k\t) можно положить равными дельта-функции Дирака S(t).
Метод Гельфанда—Левитана является одним из наиболее распространённых в теории обратных задач и заключается в сведении нелинейной обратной задачи к однопараметрическому семейству линейных интегральных уравнений Фредгольма.
В 1951 г. И. М. Гельфанд и Б. М. Левитан [1] предложили метод восстановления оператора Штурма-Лиувилля по спектральной функции, а также сформулировали достаточные условия для того, что заданная монотонная функция являлась спектральной функцией оператора. М. Г. Крейн [2] рассмотрел обратную задачу о натянутой струне и сформулировал теоремы о разрешимости обратной задачи. А. С. Благовещенский [3] привёл новое доказательство результатов М. Г. Крейна по теории обратных задач для уравнения струны.
Идеи метода Гельфанда—Левитана активно использовались в обратных динамических задачах сейсмики, начиная с работ А. С. Алексеева [4], G. Kunetz [5], Б. С. Парийского [6, 7]. А. С. Благовещенский [3, 8] разработал динамический (во временной области) вариант метода Гельфанда—Левитана для обратной задачи акустики. А. С. Алексеев и В. И. Добринский [9] использовали дискретный аналог метода Гельфанда—Левитана при исследовании численных алгоритмов решения одномерной обратной динамической задачи сейсмики. Подробный обзор численных методов решения уравнений типа Гельфанда—Левитана изложен в работе Б. С. Парийского [7]. В работе С. И. Кабанихина [10] предложен новый алгоритм решения уравнения Гельфанда—Левитана, основанный на использовании достаточного условия разрешимости обратной задачи. В монографии В. Г. Романова и С. И. Кабанихина [11] динамический вариант метода Гельфанда—Левитана применяется к одномерной обратной задаче геоэлектрики для квазистационарного приближения системы уравнений Максвелла. В работах М. И. Белишева [12, 13], С. И. Кабанихина [10, 14], М. И. Белишева и А. С. Благовещенского [15] получены многомерный аналог уравнения Гельфанда—Левитана.
Первые многомерные аналоги уравнения Гельфанда—Левитана были получены Л. Д. Фаддее-вым [16] и R. G. Newton [17]
Оптимальное по времени восстановление плотности с помощью дополнительной информации в спектральной и временной области измеренной на части границы, предложено М. И. Бели-шевым [18]. Предложенный метод восстановления использует сингулярные гармонические функции.
(£) — иХХ + иЦЦ-Ч\пр(х,у) Vи(£), х > 0, у — (у ,У2) е Ж2, * > 0; (1)
Спектральный вариант был реализован численно [18]. В работах [19, 20] динамический вариант метода в задачах геофизики был модифицирован на основе геометрической оптики. Двумерный аналог уравнения М. Г. Крейна Рассмотрим последовательность прямых задач (£ = (&1,&2) £ ^2):
и++ — ии х х + и
и(£) |*<о = 0, 4£)(+0,у,*) —е1(£у)5(*); (2)
Здесь (£,у) — £ 1 у 1 + £2У2.
Обратная задача: требуется определить функцию р(х,у) из соотношений (1), (2) по дополнительной информации:
и(£)(+0,у,*) — /(£)(у,*), £е X2. (3)
Предположим, что \п р(х,у) является 2ж-периодической функцией.
Рассмотрим вспомогательную последовательность прямых задач (т € X2):
т(т) — ) + и^ - V \п р(х,у) Vт(т\ х > 0, у е Ж2, * е Ж, (4)
Лт)
^1х—0 — е1(т,у)5(*), тхт) |х—0 — 0. (5)
Решение вспомогательной прямой задачи (4), (5) имеет следующий вид [10, 14, 21]: т(т)(х,у,*) — 1е1(т,у^у^Юу)^(х + *) + 5(х - *)] + т(т)(х,у,*)
Здесь т(т)(х,у,*) гладкая функция при х < |*|.
Функция и(£)(х,у,*) после нечётного продолжения по переменной * удовлетворяет условиям:
и(£) — + - V 1пр(х,у) Vи(£), х > 0, у е к2, * е К; (6)
и(£) |х—0 — ! (£)(у,*), их£)| х—0 — 0. (7)
Решения последовательностей прямых задач (4), (5) и (6), (7) связаны соотношением
и(£)(х,у,*) — [ £/т£)(^ 5)т(т)(х,у,5)б5. (8)
т
Здесь х > 0, у е Ж2, * е Ж, £ е X2 и
/ (£)(у,*) —Е /т£)(*)е;(т,у).
т
Применим оператор
х
о*} Рау)* 0
к соотношению (8) и обозначим
х
0
Можно показать [14, 22], что функция Ф(т)(х,*) удовлетворяет многомерному аналогу уравнения М. Г. Крейна:
х
г г е1(£,у)
2Ф\х,*)^] №)'(*- s)Ф(m)(x,s)ds — ^у^у, ^ (9)
Для каждого фиксированного значения х система уравнений (9) является системой линейных интегральных уравнений второго рода.
Решение обратной задачи р(х,у) можно определить по формуле:
р(х,у) = -¿[Е фт(х,х^ 0)е-^у)]-2. (10)
т
Для того, чтобы найти решение р(х,у) обратной задачи в точке хо > 0, необходимо решить систему (9), полагая х = хо, и вычислить решение р(хо,у) по формуле (10).
Дискретные аналоги уравнения Гельфанда—Левитана исследованы в работах [22-26]. Восстановление скорости с(х,у)
Рассмотрим обратную задачу определения скорости с(х,у) из следующей последовательности обратных задач (£ = 0, ± 1, ± 2,...):
^2(х,у)ы{?) = + хе Я, уе Я, I > 0;
а(к) и=о = 0, ^=о = ейу5 (х).
п(и)(0,у,1)= ! (и)(у,1), ах^(+0,у,0=0.
Пусть т(х,у) является решением задачи Коши для уравнения эйконала:
Тх2 + Ту2 = с"2(х,у), х > 0, уе Я; (11)
т|х=о = о, Гх|х=о = с"1 (0,у), уе Я. (12)
Введём новые переменные г = т(х,у), у = у и новые функции:
v(k)(г,y,t) = а№)(х,у,0, Ь(г,у) = с(х,у). (13)
Рассмотрим вспомогательную последовательность задач (т = 0, ± 1, ± 2,...) [22, 27]:
4т) = шт + Ъ^ш^) + qwi¡mm) + рш<т), г > 0, у е Я, ^ Я; (14)
т ( m)(0,y,t) = ешуЬ (t), wгm)(0,y,t) = 0. (15)
Здесь
ф,у) = 2Ъ2Ту, р(г,у) = Ъ2(г,у)(Тхх + Тгг). (16)
Предположим, что с(0,у) = Ъ(0,у) известно и для простоты изложения положим Ъ(0,у) = 1 при у € Я.
В окрестности плоскости t = г решение прямой задачи (14), (15) имеет вид:
ш
(т)
(г,у,0 = 3(т)^,у)8(г^ ^ + д(т)^,у)в(г ^ ^ + Ш(т)(г,уД (17)
Здесь ш(т) является непрерывной функцией, а функции 5(т) и являются решением следующих обратных задач:
25(т) + дЗут^ + р5(т) =0, t > 0, у е Я; (18)
5(т) ^=о = 1е1ту. (19)
2п(т) = с(т) =
^т) + Ъ^ут + рф(т) , t > 0, у е Я; (20)
Я(т)\ t=о = 0. (21)
Тем самым из (17) получим двумерный аналог уравнения М. Г. Крейна (т = 0, ± 1, ± 2,.. .):
г
^(т)(гМк)(^ г) + т^(г,у,0 + ^ [¡Ы'«- з)Ш(т\г,у,з)йз = 0, Щ < г. (22)
«-и «-и ™
Построение ^-приближения двумерного аналога уравнения М. Г. Крейна
Как и предыдущем пункте, рассмотрим метод регуляризации для численного решения двумерной обратной задачи (1)—(3), основанный на проекции изначальной постановки задачи на N-мерное подпространство.
Предположим, что решение задачи (1), (2) можно представить в виде ряда:
и(к)(х,у,£) = ^ и«1\х,1УтУ.
т
Также предположим, что функция 1пр(х,у) представима в виде ряда:
1пр(х,у) = £ ат(х)ету.
т
Тогда обратную задачу (1)—(3) можно записать в виде бесконечной системы одномерных обратных задач на коэффициенты Фурье функций:
гР„(к) гР„(к)
*ип (Ги" а2а^)(х,!У
дР дх2 _ ^ ^ _
" X ~хат(х)~хи{(')- т(х^) + X т(т - п)ат(х)и{П1т(х^), (23)
тт
х € К, £ > 0, к,п е иП] I £=0 = 0, ^ |£=0 = 8кп8 (х); (24)
иП} |х=0 = ¡Пк) (£). (25)
Здесь 8ип — символ Кронеккера.
Обратной задаче (23)—(25) сопоставим ее N-приближение:
я2у (к) я2у (к) Яу (к)
^г- = -ЩТ ~ КУ(к) + А(х)УВ(х)-у^-, х > 0, £ > 0; (26)
гШ (к)
V(к) |£=0 = 0, ^ЩГ£=0=/(к^(х); (27)
V (к)| х=0 = Р (к)(£). (28)
Здесь
V(к)(х,£) = (¿%(х,1)/%+у(х,1),..., 4к)(х,£),..., ¿%(х,1)),
Р (к](£) = (!^(О,...,!%т
I(к) — вектор, к -ый элемент которого, равен 1, а все остальные элементы — нули.
Здесь квадратные матрицы К, А, В имеют размерность (2N + 1) х (2N +1) и состоят из следующих элементов:
Ккт = тЧит, \к\ < N, |т| < N; Акт(х) = 0^- |к - т|)т( к ^т)рк-т(х), Щ,\т\ < N, |к - т\ < N;
Вщт(х) = \к - т\)-^Рк-т(х), |к|,|т| < N, |к - т\ < N.
Выпишем N-приближение для двумерного аналога уравнения М. Г. Крейна (9):
х
2Ф(к) (х,£) - Е [ ¡^'(^ х) Ф(т) (х,^ = в(к), (29)
Н^ ¿х
где < х, Щ < N. Уравнение (29) можно записать в векторно-матричной форме:
2Ф(хД) — з)Ф(х,з)? = в
< х.
(30)
,т
Здесь неизвестная вектор-функция Ф(х^) = (Ф( N)(x,t),. .. ,Ф(0)(х,^),... ,Ф(^(хД)) , вектор правой части в = (в(~N),... С(0),... ^^ и матрица
/ . (- N)'
/—N
Г (-N +1)'
Г (0)'
\ Г (N)' \ /—N
Г (—N)^ /-N+1 . ;(-N +1/ -N +1 . Г (-N) ' . . /—N+1 Г (-N +1)' . . /—N +1
Г (0) ' Г (0) ' . . /—N+1
Г N) ' ^N+1 . Г N) ' . . /—N+1
Г (- N)' \ /N
Г (-N +1)' /N
Г (0) /N
Г N)' /N
(0.
/
Проведём дискретизацию системы уравнений (30). Для этого разобьём отрезок (—х,х) на 2Nx + 1 равных частей. Обозначим Н = х/^, Ф£ = Ф(х,£Н), и пусть t = ¡Н и ¡| < Nx. Тогда получим следующую систему линейных алгебраических уравнений:
- н
Н
N-1
£
к=-N+1
А = с,
п< N
Запишем полученную систему в виде:
АФ = в.
Квадратная матрица А имеет размерность (2Nx + 1) х (2N +1). Рассмотрим восстановление функции плотности (см. рис. 2):
(31)
р(х,у) =
1 + + 1) exp{^160(x^ 0,3)2},
1 + 1.5 cos(7гу + 1) exp{^200(x^ 0,8)2}, 1,
х е [0.1,0.5],у е [-2,2], х е [0.7,0.9],у е [-1,1],
в остальных случаях.
Зафиксируем следующие параметры задачи: х = 1, Nx = 50, Ny = 50.
Приближённые данные брались в виде:
Г (у, 0 = I (у,0 + ea(y,t)(fmax- /min).
Здесь е уровень шума в данных, а(у^) равномерно распределённое случайное число на отрезке [—1,1] для фиксированного значения переменных у и t, fmax и /шш максимум и минимум значений точных данных.
На рисунках 3, 4 восстанавливается гладкая функция плотности р(х,у) (см. рисунок 2) в случае зашумлённых данных е = 0.05.
В работах [28, 29] рассматривается обратная задача для определения коэффициента д(х) гиперболического уравнения ии = Ах,уи — ц(х,у)и. В результате использования многомерного аналога метода Гельфанда—Левитана задача была сведена к семейству линейных интегральных уравнений. В работе [28] полученные уравнения решались методом Монте-Карло — искомое решение представлялось математическим ожиданием некоторой случайной величины. При этом использовался тот факт, что для восстановления коэффициента достаточно вычислить решение уравнения Гельфанда—Леви-тана только в одной точке. В работе [29] был использован стохастический проекционный метод, основанный на сведении интегрального уравнения к линейной системе уравнений и проектировании решения на гиперплоскость, задаваемую случайной строкой уравнения. В отличие от большинства стохастических алгоритмов решения СЛАУ, данный подход не ограничен сходимостью ряда Неймана.
х
х
Г
Рис. 2. Точное решение обратной задачи N = 30
Рис. 3. Приближённое решение обратной задачи в случае 5 источников и 5 приёмников N = 5 и уровнем шума в данных 5 %
Рис. 4. Приближённое решение обратной задачи в случае 10 источников и 10 приёмников N =10 и уровнем шума в данных 5 %
В работе [30] была рассмотрена двумерная обратная задача акустики. Используя метод Гель-фанда—Левитана—Крейна в сочетании с проектированием задачи на конечное подпространство, задаваемое базисом Фурье, задача также сводилась к однопараметрическому семейству линейных интегральных уравнений типа Фредгольма второго рода. Для решения этих уравнений был использован алгоритм обращения матрицы, опирающийся на блочно-тёплицевую структуру системы. Дальнейшее использование структуры задачи позволило получить решение семейства уравнений Крейна за счёт решения только одного уравнения. Это позволило значительно улучшить эффективность метода.
Работа выполнена при финансовой поддержке Министерства образования и науки Российской Федерации, Министерства образования и науки Республики Казахстан (название проекта «Выявление и изучение геофизическими методами проницаемых структур верхней части разреза Семипалатинского испытательного полигона в целях повышения геоэкологической безопасности») и Российским фондом фундаментальных исследований (проекты 14-01-31182, 15-01-09230, 16-01-00755).
ЛИТЕРАТУРА
1. Гельфанд И. М., Левитан Б. М. Об определении дифференциального уравнения по его спектральной функции // Известия АН СССР. 1951. Т. 15, № 4. С. 309-360.
2. Крейн М. Г. Об одном методе эффективного решения обратной краевой задачи // Доклады АН СССР. 1954. Т. 94, № 6. С. 987-990.
3. Благовещенский А. С. Об обратной задаче теории распространения сейсмических волн // Проблемы мат. физики. 1966. Т. 1. С. 68-81.
4. Алексеев А. С. Обратные динамические задачи сейсмики // Некоторые методы и алгоритмы интерпретации геофизических данных. М. : Наука, 1967. С. 9-84.
5. Kunetz G. Essai d'analyse de traces sismiques // Geophysical Prospecting. 1961. Vol. 9. P. 317-341.
6. Парийский Б. С. Обратная задача для волнового уравнения с воздействием на глубине // Некоторые прямые и обратные задачи сейсмологии. 1968. С. 139-169.
7. Парийский Б. С. Экономичные методы численного решения уравнений в свертках и систем алгебраических уравнений с теплицевыми матрицами. М. : ВЦ АН СССР, 1977.
8. Благовещенский А. С. О локальном методе решения нестационарной обратной задачи для неоднородной струны // Тр. МИАН СССР. 1971. Т. 115. С. 28-38.
9. Алексеев А. С., Добринский В. И. Некоторые вопросы практического использования обратных динамических задач сейсмики. Математические проблемы геофизики // АН СССР. Сиб. отд-ние. Вычислительный центр. Новосибирск. 1975. Т. 6, № 2. С. 7-53.
10. Кабанихин С. И. Линейная регуляризация многомерных обратных задач для гиперболических уравнений. Препринт № 27 института математики СО АН СССР. Новосибирск, 1988.
11. Романов В. Г., Кабанихин С. И. Обратные задачи геоэлектрики. М. : Наука, 1991. С. 304.
12. Белишев М. И. Уравнения типа Гельфанда—Левитана в многомерной обратной задаче для волнового уравнения // Зап. научн. сем. ЛОМИ. 1987. Т. 165. С. 15-20.
13. Белишев М. И. Об одном подходе к многомерным обратным задачам для волнового уравнения // ДАН СССР. 1987. Т. 297, № 3. С. 524-527.
14. Кабанихин С. И. О линейной регуляризации многомерных обратных задач для гиперболических уравнений // Доклады РАН. 1989. Т. 309, № 4. С. 791-795.
15. Белишев М. И., Благовещенский А. С. Многомерные аналоги уравнений типа Гельфанда—Левита-на—Крейна в обратной задаче для волнового уравнения // Условно-корректные задачи математической физики и анализа. Новосибирск : Наука, 1992. С. 50-63.
16. Фаддеев Л. Д. Обратная задача квантовой теории рассеяния. II // Итоги науки и техн. Сер. Соврем. пробл. мат. 1974. Т. 3. С. 93-180.
17. Newton R. G. Inverse Schrodinger scattering in three dimensions // Texts and Monographs in Physics. Berlin : Springer-Verlag. 1989.
18. Belishev M. I. Boundary control in reconstruction of manifolds and metrics (the BC method) // Inverse Problems. 1997. Vol. 13, no. 5. P. R1-R45.
19. Belishev M. I. How to see waves under the Earth surface (the BC-method for geophysicists) // Ill-Posed and Inverse Problems. Boston : VSP, Utrecht, 2002. P. 67-84.
20. Belishev M. I. Dynamical Inverse Problem for the Equation utt — Au ^ VpVu = 0 (the BC-method) // CUBO A Mathematical Journal. 2008. Vol. 10, no. 2. P. 17-33.
21. Kabanikhin S. I., Shishlenin M. A. Numerical algorithm for two-dimensional inverse acoustic problem based on Gelfand—Levitan—Krein equation // J. Inverse and Ill-Posed Problems. 2011. Vol. 18, no. 9. P. 979-996.
22. Kabanikhin S. I., Satybaev A. D., Shishlenin M. A. Direct Methods of Solving Inverse Hyperbolic Problems. The Netherlands : VSP, 2004.
23. Gladwell G. M. L., Willms N. B. A discrete Gelfand—Levitan method for band-matrix inverse eigenvalue problems // Inverse Problems. 1989. Vol. 5. P. 165-179.
24. A discrete Gelfand—Levitan theory : Technical report / Institut flir Numerische und instrumentelle Mathematik Universitat Miinster Germany ; Executor: F. Natterer : 1994.
25. Kabanikhin S. I., Bakanov G. B. Discrete analogy of Gelfand—Levitan method // Doklady Akademii Nauk. 1997. Vol. 356, no. 2. P. 157-160.
26. Кабанихин С. И., Баканов Г. Б. Дискретный аналог метода Гельфанда—Левитана в двумерной обратной задаче для гиперболического уравнения // Сиб. мат. журн. 1999. Т. 40, № 2. С. 307-324.
27. Kabanikhin S. I. On linear regularization of multidimensional inverse problems for hyperbolic equations // Sov. Math. Dokl. 1990. Т. 40, № 3. С. 579-583.
28. Kabanikhin S. I., Sabelfeld K. K., Novikov N. S., Shishlenin M. A. Numerical solution of an inverse problem of coefficient recovering for a wave equation by a stochastic projection methods // Monte Carlo Application. 2015. Vol. 21, no. 3. P. 189-203.
29. Kabanikhin S. I., Sabelfeld K. K., Novikov N. S., Shishlenin M. A. Numerical solution of the multidimensional Gelfand—Levitan equation // Journal of Inverse and Ill-Posed Problems. 2015. Vol. 23, no. 5. P. 439-450.
30. Kabanikhin S. I., Novikov N. S., Oseledets I. V., Shishlenin M. A. Fast Toeplitz linear system inversion for solving two-dimensional acoustic inverse problem // Journal of Inverse and Ill-Posed Problems. 2015. Vol. 23, no. 6. P. 687-700.
