Рандомизированные методы решения краевых задач математической физики тема диссертации и автореферата по ВАК РФ 01.01.07, кандидат наук Моцартова, Надежда Сергеевна

• Моцартова, Надежда Сергеевна
• кандидат наук
• 2013

Оглавление диссертации кандидат наук Моцартова, Надежда Сергеевна

Оглавление

Введение

1 Рандомизированные алгоритмы для решения больших систем линейных алгебраических уравнений

1.1 Рандомизация с помощью случайных разреженных матриц

1.2 Итерационные методы решения систем линейных алгебраических уравнений

1.2.1 Метод, основанный на преобразовании спектра

1.2.2 Нестационарный итерационный процесс со случайными параметрами

1.2.3 Рандомизация метода Гаусса-Зейделя

1.3 Дискретные варианты алгоритма блуждания по границе

1.3.1 Алгоритм изотропного случайного блуждания по границе

1.3.2 Дискретная версия метода случайного блуждания по границе

1.4 Численные результаты

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

2.1 Сингулярное разложение матриц

2.2 Варианты рандомизации сингулярного разложения матриц

2.2.1 Метод разреженного сингулярного разложения матриц

2.2.2 Рандомизированный метод главных компонент

2.3 Моделирование случайных полей

2.3.1 Разложение Кархунена-Лоева

2.3.2 Дискретизация разложения Кархунена-Лоева

2.3.3 Численные эксперименты по моделированию неоднородных случайных полей

2.3.4 Численные результаты сравнения методов рандомизированного сингулярного разложения

2.4 Алгоритмы блуждания по границе

3 Стохастические граничные методы фундаментальных решений и их

приложения

3.1 Формулировка метода фундаментальных решений

3.2 Уравнение Лапласа

3.3 Уравнения Ламе

3.4 Аппроксимация функции Грина

3.5 Дискретизация интеграла Пуассона

3.6 Метод, основанный на обращении интегральной формулы Пуассона

3.7 Сравнительные эксперименты

3.8 Задача о нахождении электроемкости для цепочки сфер

Заключение

Список литературы

Введение диссертации (часть автореферата) на тему «Рандомизированные методы решения краевых задач математической физики»

Введение

Статистическое моделирование (или метод Монте-Карло) - это численный метод решения математических задач при помощи моделирования случайных величин [1]. Датой рождения метода Монте-Карло принято считать 1949 год, когда была опубликована статья Метрополиса и Улама [2]. Однако до появления вычислительных машин этот метод не имел широкого применения из-за трудоемкости моделирования большого числа случайных величин. Таким образом, возникновение и интенсивное развитие методов Монте-Карло стало возможным благодаря появлению вычислительной техники.

Метод Монте-Карло естественно применяется при моделировании любого процесса, на протекание которого влияют случайные факторы, например при решении задач статистической физики, теории турбулентности, теории переноса. Статистическое моделирование используется и для решения многих математических задач, не связанных с какими-либо случайностями, но для которых можно искусственно придумать вероятностную модель, таких как вычисление числа 7Г (предложенная Бюффоном в 1777 году), решение линейных интегральных и алгебраических уравнений, обращение матриц [3-5]. В первой работе на эту тему [2] была представлена схема Неймана-Улама, область применения которой ограничена условием сходимости ряда Неймана. Основная идея метода Неймана-Улама состояла в построении случайных оценок членов ряда Неймана на основе цепи Маркова в соответствии с ядром данного линейного уравнения. Применения методов Монте-Карло для решения интегральных уравнений были традиционно ограничены задачами переноса излучения и уравнениями диффузии [6,7]. Для решения применялся метод случайного блуждания по решетке (одни из первых работ по случайному блужданию для решения задач в частных производных в конечномерном пространстве [8,9]). В работе [10] представлены варианты различных модификаций первоначальной схемы Неймана-Улама, направленные на улучшение свойств дисперсии. Следует отметить большое количество других методов [11,12], представляющих собой рандомизированные оценки итераций на основе цепи Маркова. В работе [13] был разработан метод случайного блуждания по сферам, основанный на методе релаксации и методе Гаусса-Зейделя с черно-красной нумерацией.

В 1982 году К. Сабельфельд предложил в своей работе [14,19] первый метод Монте-Карло для решения интегральных уравнений с расходящимся рядом Неймана. Основная идея заключалась в перегруппировке членов ряда Неймана интегрального уравнения и подборе весов, таких, чтобы преобразованный ряд сходился и было возможным его устойчивое вычисление. Это было достигнуто с помощью конформного преобразования спектрального параметра интегрального уравнения, хорошо известного из линейной алгебры [15,16]. Фактически каждое конформное преобразование порождает соответствующий итерационный процесс [17].

Мотивацией при разработке метода, представленного в работе [14], было сконструировать метод Монте-Карло, который позволяет решать граничные интегральные уравнения, представимые в виде потенциалов простого и двойного слоя, плот-

ность которых удовлетворяет интегральным уравнениям со спектральным радиусом равным единице. Метод преобразования спектрального параметра был применен для решения таких задач. В итоге был сконструирован метод блуждания по границе, который наиболее эффективен для выпуклых областей. В случае невыпуклых областей данный метод также применим, но дисперсия существенно возрастает с увеличением "невыпуклости области". Метод блуждания по границе применим для решения всех классических внутренних и внешних граничных задач теплопроводности, упругости и электростатики [18-20]. Эргодическая версия алгоритма случайного блуждания по границе была предложена в работе [21] для решения задачи Робена. Множество других применений алгоритма блуждания по границе было разработано, например, в работе [22] рассматривались задачи теории упругости, а в [23] данный метод используется для вычисления электроемкости. Все эти алгоритмы хорошо применимы для выпуклых областей и сталкиваются с проблемой сильного увеличения дисперсии в случае невыпуклой области. Некоторые преобразования были предложены для того, чтобы как-то преодолеть быстрое увеличение дисперсии, например метод стабилизации в [19], раздел 6.2, или ветвление изотропного случайного блуждания по границе [24].

Новая дискретная версия алгоритма блуждания по границе предложена в нашей работе [25]. Применение данного алгоритма для решения граничных интегральных уравнений можно разделить на две основные части: 1) построение сетки на границе области и дискретного аналога интегрального уравнения, 2) применение метода рандомизированной оценки большой матрицы набором случайных матриц меньшей размерности. Основная особенность данного алгоритма состоит в том, что дисперсия незначительно возрастает при переходе от выпуклой области к невыпуклой. Идея рандомизации больших матриц набором случайных матриц меньшей размерности не является новой. Этот метод применяют для построения аппроксимации меньшего ранга [26-30], для рандомизации итерационного метода Ланцоша [31], при вычислении произведения матриц [32,33]. Во всех этих работах было подтверждено, что использование рандомизации является эффективным при работе с очень большими системами.

С увеличением мощности вычислительной техники размерность задач и требуемый объем вычислений постоянно увеличиваются. Рандомизированные методы все чаще появляются в задачах, где требуется работать с данными большой размерности, например, при решении многомерных задач для уравнений в частных производных, трехмерных интегральных уравнений теории потенциала, при моделировании случайных полей, решении обратных задач томографии или кристаллографии, при решении уравнения Шредингера, при моделировании турбулентных потоков и т.д. Область применения методов статистического моделирования для решения задач большой размерности ограничена следующими условиями: (1) дисперсия оценки мала, (2) не требуется большой точности, (3) вычислительная сложность случайной оценки растет медленно с возрастанием размерности т. Условие (3) достаточно часто выполнено, а вот условия (1) и (2) являются основной проблемой, так как скорость сходимости методов Монте-Карло мала, примерно как а/ у/П, где а - среднеквадра-

тическое отклонение и^- размер выборки. Поэтому точность методов Монте-Карло на практике обычно не превосходит 1—2%. Уменьшение дисперсии часто достигается с помощью некоторого преобразования, результатом применения которого является случайная оценка с меньшей дисперсией. Уменьшение размерности задачи является одним из возможных приемов для уменьшения дисперсии. Например, дисперсия будет уменьшена, если при вычислении многократного интеграла возможно интегрирование только по части переменных.

За последние пару десятков лет, рандомизированные методы, основанные на явлении концентрации меры [34], интенсивно развивались. Первые результаты для концентрации были получены П. Леви в его книге по функциональному анализу [35]. Явление концентрации мер приобрело большую популярность в математическом сообществе и нашло многочисленные приложения в функциональном анализе, геометрии, теории вероятности, комбинаторике и технических науках. Разработано много статистических методов для работы с большими матрицами, которые основаны именно на этой концепции. В линейной алгебре особенно стоит отметить следующие методы: рандомизированный проекционный метод, который был признан мощным инструментом для уменьшения размерности [36-39]. Также следует отметить фундаментальный результат Джонсона и Линденштрауса [40]: для любого 0 < е < 1, множества X, состоящего из т точек в пространстве М^, и для любого целого числа к > О(^р), существует отображение /: —» такое что Ми, V 6 X выполнено

(1 - е)||и - < ||/(и) - /(«)|| < (1 + е)||и - г>||2-

В последнее время появилось множество работ по рандомизации различных матричных операций: произведение матриц, сингулярное разложение и др., [25,26,28-32,4153]. В данной работе для нахождения решения применялись несколько рандомизированных алгоритмов для работы с матрицами большой размерности: метод разрежения, рандомизированный метод сингулярного разложения и рандомизированный метод главных компонент.

В последние десятилетия значительные усилия исследователей были направлены на развитие так называемых бессеточных численных методов. Среди этих методов большую популярность завоевал метод фундаментальных решений - граничный бессеточный метод аппроксимации решений однородных эллиптических уравнений. Для аппроксимации решения однородного эллиптического уравнения Ь[и] = 0 метод фундаментальных решений использует линейную комбинацию сингулярных решений £(х, у) этого же уравнения с координатами сингулярностей у (координатами источников) за пределами области решения. Таким образом, применение метода фундаментальных решений предполагает, что фундаментальные решения рассматриваемого уравнения известны в явном аналитическом виде. Хотя это обстоятельство сужает сферу применения, тем не менее, в последние десятилетия этот метод нашел обширные приложения в инженерных расчетах. При использовании даного метода исследователи столкнулись со следующими особенностями: (1) полученные системы уравнений являются плохо обусловленными, (2) в случае сложных невыпуклых областей его применение приводит к решению прямоугольных систем большой

размерности. Метод фундаментальных решений был впервые предложен В. Д. Купрадзе и М. А. Алексидзе в работах [54], [55] как метод обобщенных рядов Фурье. Сходимость была изучена в работе [56]. Более общий подход, основанный на аппроксимации граничного условия линейной комбинацией полного набора функций, использовался в методе Треффца [57], [58], [59]. Применение сингулярного разложения матриц в методе фундаментальных решений описано в работах [60], [61], [62]. В диссертации предложен новый стохастический граничный алгоритм решения многомерных краевых задач на основе рандомизированных вариантов метода фундаментальных решений, при этом координаты источников выбираются случайно распределенными, и задача сводится к решению большой прямоугольной системы с помощью рандомизированного сингулярного разложения матриц.

Перейдем к изложению содержания диссертации по главам.

Первая глава диссертации посвящена применению алгоритмов аппроксимации исходной матрицы последовательностью разреженных матриц к решению систем линейных алгебраических уравнений. В разделе (1.1) приводится стандартная схема рандомизации с помощью случайных разреженных матриц и оценки несмещенности. В (1.2) кратко описываются итерационные методы, к которым будет применен данный алгоритм: (1) метод, основанный на смещении спектра, (2) нестационарный итерационный процесс со случайными параметрами, (3) метод Гаусса-Зейделя. В этой главе также предложен новый алгоритм двойной рандомизации для метода Гаусса-Зейделя. В разделе (1.3) на примере решения двумерной задачи Дирихле для уравнения Лапласа сравниваются стандартный алгоритм случайного блуждания по границе и разработанная дискретная версия метода случайного блуждания по границе. Раздел (1.4) содержит численные результаты применения метода аппроксимации матрицы последовательностью разреженных матриц к разным итерационным методам, поведение дискретного алгоритма случайного блуждания по границе для различных областей и влияние выбора строк на сходимость методов. Результаты главы опубликованы в работе [25].

Вторая глава посвящена рандомизированному сингулярному разложению и его применению к моделированию случайных процессов и полей. Первые два раздела (2.1) и (2.2) являются вводными и содержат хорошо известную в литературе информацию. В разделе (2.1) дается общая теория сингулярного разложения матрицы и приводится формулировка теоремы Экарта-Янга. В следующем разделе (2.2) приводятся два алгоритма рандомизированного сингулярного разложения: метод разреженного сингулярного разложения матриц и рандомизированный метод главных компонент, также приводятся оценки погрешности и трудоемкости обоих методов. Моделирование случайных процессов и полей рассмотрено в разделе (2.3). При моделировании случайных полей использовалась следующая схема: рассматривались случайные процессы и поля с заданной корреляционной матрицей, с помощью алгоритма рандомизированного сингулярного разложения находились собственные числа и вектора корреляционного оператора. Далее, используя разложение Кархунена-Лоэва из раздела (2.3.1) (в разделе (2.3.2) приведен его дискретный аналог), полу-

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

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

í=l

где г-ранг системы, а, > а2 > ... > о,, сингулярные значения, и и^ е М'п, г/') £ ]1{а, г = 1 ,...,г левые и правые сингулярные вектора соответственно. Аппроксимируя первых 5, в <С г сингулярных числа и левых/правых сингулярных вектора, перейти к

новой системе с матрицей 5 = А — ^ 0{ Преобразованная система уравнений

»=1

с матрицей системы 5 уже будет пригодна для применения метода простой итерации. Результаты главы опубликованы в работе [63].

Третья глава посвящена методу фундаментальных решений. В разделе (3.1) содержится общая информация о классическом методе фундаментальных решений, который был разработан В. Д. Купрадзе и М. А. Алексидзе. Применение метода для решения задачи Дирихле для уравнения Лапласа содержится в разделе (3.2) и для уравнений Ламе в разделе (3.3). В разделе (3.4) приведен новый метод, который позволяет аппроксимировать функцию Грина и находить решение неоднородной задачи. Используя интеграл Пуассона, было реализовано два новых подхода: дискретная аппроксимация интеграла, описанная в разделе (3.5), и его представление в виде ряда, представленная в разделе (3.6). Численные эксперименты, представленные в разделе (3.7), показали, что второй метод приводит к уменьшению числа обусловленности задачи. В разделе (3.8) приведены численные результаты решения задачи о нахождении электроемкости для цепочки сфер. На примере решения задачи о нахождении электроемкости была построенна аппроксимация нормальных производных для трехмерной внешней задачи с использованием подхода, предложенного в работе В. Д. Купрадзе и М. А. Алексидзе [54]. Результаты главы опубликованы в работе [64].

Глава 1

Рандомизированные алгоритмы для решения больших систем линейных алгебраических уравнений.

Рассмотрим методы построения решения системы линейных алгебраических уравнений х = Ах 4- Ь размерности п с матрицей {Л, р(Л) < 1}, где р{А) - спектральный радиус матрицы А. Для данной системы применим метод простой итерации хк+1 = Ахк + ь

Сначала рассмотрим стандартный метод Монте-Карло, основанный на цепи Маркова [80]. Моделируется цепь Маркова (г'о, г'ь ..., г'дт} соответствующая конечному фазовому пространству (1,2,..., т). Она определяется вероятностями перехода гг0 из состояния г в у, начальными вероятностями 7Г, = Р{го = г) и вероятностями обрыва рьу при переходе из состояния г в Оценка по столкновениям имеет вид

При условии, что 7г( ф 0, если Ь1 ф 0, и ф 0, если а1и Ф 0, выполнено = Если требуется найти оценку только для компоненты Хк, то = 1 и = 0 для V? ф к. На каждом шаге стандартного метода Монте-Карло используется только одна компонента матрицы а,^.

Рассмотрим модификацию стандартного метода, для которой на каждом шаге будет производиться умножение столбца или строки матрицы А. На первом шаге итераций строится случайная оценка £1 вектора А Ь (£о = Ь)

где Л« - г'-й столбец матрицы А, {р,,г = 1 ,...,п} - вероятностное распределение столбцов матрицы А и и - случайный индекс, имеющий распределение р, и £1 - случайный вектор, равный Предполагается выполнение условия р„ ф 0, если Ьу Ф 0. Здесь и далее математическое ожидание обозначается Е.

На втором шаге аппроксимируется вектор А2Ь — А(АЬ) и рассматривается в качестве несмещенной случайной оценки вектора А Ь

а2ь = е 'а^А^Ь; = Е

. Р-у .

где 1/, 7 - независимые случайные индексы, распределенные в соответствии с р. Для оценки А:-ой итерации ^ используется оценка с предыдущего шага и т.д.

Этот подход требует минимум памяти и трудоемкость ведет себя как 0(п2), погрешность убывает как сг/у/М, где <т2 - дисперсия оценки, N - статистика. Сходимость данной схемы можно улучшить, используя оптимальное распределение столбцов Рг, которое будет приведено в разделе 1.1.

В этой главе предлагается более общий тип случайных оценок, которые используют несмещенные случайные матричные оценки для каждой итерации: А Ь = Е 5 V, где в - случайная разреженная матрица и V - случайная оценка вектора Ъ. Итерационный процесс с построенными оценками для решения системы линейных алгебраических уравнений организуется таким образом, чтобы матрично-векторные оценки на разных шагах итераций были независимы. Отметим, что в работе [33] был предложен следующий способ построения случайной оценки матрицы А: матрица 5 состоит из I случайно выбранных столбцов исходной матрицы А с некоторыми весами, остальные компоненты матрицы равны нулю. Дисперсия такой оценки может быть уменьшена с помощью увеличения числа I, размера подматрицы 5.

1.1 Рандомизация с помощью случайных разреженных матриц.

Рассматривается система линейных алгебраических уравнений следующего вида:

х = Ах + Ь, (1.1)

х = (х\,..., хп)т, Ь = (6Ь ■ • ■, Ьп)т е Яп, А = {Агу}"у=1 и определитель матрицы А не равен нулю.

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

итерационного процесса

х(т+1) = Лх(т) + Ь; х(0) = Ь; т = о, 1,2..........(1.2)

Запишем аппроксимацию решения конечным отрезком ряда Неймана

м

= ^ЛЧэ . (1.3)

1=0

Начнем с рассмотрения подхода представленного в работах [33,52]. Пусть С? = - несмещенная оценка матрицы А, определенная как случай-

ная матрица, такая что Ед^ = € {1,... ,п}, и пусть ...,

последовательность независимых реализаций случайной матрицы (7. Определим случайный процесс следующим образом

£(т+1) = С(т)£(т) + Ь -щ = 0, 1, . . . , М — 1, (1.4)

где = Ь. В силу того, что т = 0,1,... независимы, из равенства (1.4)

получим Е£{м) = х<м\

Рассмотрим частный случай построения случайной матрицы С в виде разреженной матрицы. Мы будем строить матрицу О по столбцам: выберем случайное множество 3, состоящее из I целых чисел равномерно распределенных от 1 до п без повторения и I п. Это означает, что выбирается как целое число, равномерно распределенное среди 1,2,..., п, затем ]2 - целое число, равномерно распределенное среди оставшихся п — 1 чисел, и т.д. Определим матрицу (7

Gík = {у А%к для к е 3; 0 иначе|

для г = 1,2,..., п. Таким образом, (7 - случайная матрица размерности п, которая имеет ровно I ненулевых столбцов, и для любых г, к Е {1,2,..., п} справедливо

ЕС* = в* р{к ез} = -хА1кр{к е 3} = А1к.

Здесь и далее вероятность случайного события обозначается Р. Следует отметить, что для вычисления вектора необходимы только I компонент вектора В

свою очередь, для получения необходимо знать I компонент и т.д. Со-

ответственно, на каждом шаге итераций количество операций будет равно I2. Для аппроксимации решения необходимо NМ12 операций, где N - объем статистики и М - длина отрезка ряда Неймана.

Описанный выше процесс можно рассматривать как рандомизацию вычисления произведения матрицы на вектор. Идея алгоритма рандомизации произведения матриц может быть объяснена через рассмотрение следующего вычислительного метода Монте-Карло для набора интегралов

Л = J Мх)д(х)ёх , (1.5)

с

где ¡1, г = 1,2..., п семейство функций, определенных на области интегрирования С. Пусть плотность вероятности р(х) определена на С и удовлетворяет условию согласования: р(х) ф 0 если /г(х)д(х) ф 0 для г = 1,2, ...п. Выберем I независимых

случайных точек к = 1,2,...в области С в соответствии с плотностью р(х). Несмещенная оценка для интегралов 1г определяется следующим образом

Х^уЕШ)^, г = 1,2... ,п. (1.6)

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

Е ^ = /Е ~ Е ' (1-т)

.=1 £ г=1 ,=1

где сг - заданные вещественные константы и дисперсия Хг обозначается Ихг- Эта проблема была сформулирована и решена в книге Ермакова С.М. и Михайлова Г.А. [65]. Приведенный выше метод может быть применен для случайной оценки произведения матрицы на вектор, где используется одно и тоже фиксированное случайное распределение столбцов матрицы. Рассматривая дискретный аналог случайной оценки интегралов (1.5), получаем оценку для произведения матрицы на вектор, где ¿ = 1,2,... есть номера строк и интегрирование по х соответствует сумме столбцов.

Минимум выражения (1.7) достигается, когда достигается минимум выражения

в

Рассмотрим случайную оценку

'/=(Ес»Л2Ю)1/2^Ц интеграла 1 = | ( £ с2/2) ^М* .

Здесь £ - случайная точка из области С, распределенная с вероятностью р(х). Из того, что Б = Ет]2 и От] = 5 — I2 > 0, следует, что

(ЕсЦ^яШх)2 . а '=1

Прямой подстановкой можно проверить, что равенство достигается для

Р(*) = у( Е^Л2)1/2Ф). (1.8)

<=1

_ п

отсюда следует, что минимальное значение дисперсии равно = /2 — •

1=1

Таким образом, если необходимо получить случайную оценку для всех интегралов

/,, г = 1,2..., п одновременно, то можно использовать один и тот же набор случайных точек £2 • • • * £п с плотностью распределения Значение оптимальной плотности дается формулой (1.8).

Рассмотрим метод, представленный выше, в более общем случае. Пусть необходимо построить несмещенную оценку произведения двух матриц, тп х п матрицы А и пхр матрицы В. Схема, представленная ниже, достаточно часто используется в анализе данных и при построении случайной низкоранговой аппроксимации [32,50,51].

Произведение матрицы А и В может быть представлено в виде суммы АВ = ]С"=1 гДе Л(т) - столбец матрицы А с номером т, и В^ - строка матрицы

В с номером т. Пусть также {рг,р2, ■ ■ ■ ,Рп} - некоторый набор вероятностней, определяющий случайное распределение индексов в диапазоне от 1 до п. Матрицы Б и Л - рандомизированные несмещенные оценки матриц А и В соответственно. Введем обозначение Пи для элемента г'-й строки и ^'-го столбца матрицы О.

Рандомизированная оценка произведения матриц А В строится следующим образом:

1. Для г от 1 до / выбираем независимо случайные числа гт из интервала от 1 до п

в соответствии с вероятностным распределением {рь, к = 1,..., п}. Столбец матрицы 5 выбирается как ^4(1т)/\Да7> строка матрицы Я выбирается равной ВМ/уДр^.

2. Несмещенная оценка для произведения матриц АВ - матрица 57?. Приведем доказательство несмещенности оценки Б К

Е[(ЯН),;]=(АВ),„ ¡,1 = 1,...,п . (1.9)

Действительно, случайные величины

т_х.....,

\ 1Р,Т ) %д 1Р<т

являются независимыми по построению. Вычисляя математическое ожидание (5Я),, = ]Г)т=1Х(т) получим

Е[Хт} = р^Рк = ](АВ)г],

I

отсюда

E[{SR)l3) = £ E[Xr] = (АВ)Ч .

Т = 1

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

Выбор вероятностного распределения столбцов оказывает большое влияние на поведение дисперсии случайной оценки. Дисперсия (SR)l;l представима в следующем

виде

A\Bl

D[(SRU = ) £ - ){AB)l . (1.10)

' k=i Pk

Действительно, компоненты (SR)V представимы в виде суммы I независимых случайных величин, следовательно, дисперсия (Sfí)tJ вычисляется как сумма дисперсий. Из определения дисперсии получаем следующее выражение

" А2 R2 1

откуда вытекает равенство (1.10).

Используя схему, приведенную в (1.5), получим оптимальное вероятностное распределение {рк}- Значение весов {p¿} определяются из задачи минимизации величины £?[||Л.В — S-ñH^], что соответствует выбору с, = 1 в формуле (1.7). Здесь || * ||f - норма Фробениуса, которая для пхтг матрицы X определяется как ЦХЦ^ = Т" г2

Выбор весов pi в виде

_ \A(k)\\BW\

минимизирует дисперсию погрешности, которая в этом случае запишется как

Е[\\АВ - SR\\2F] = у (¿\AW\- \\\АВ\\1 . (1.12)

Действительно, несмещенность оценки SR получается с использованием формулы (1.10)

E[\\AB-SR\\2f]

1=1 J=1

= 7 ¿ ~\Aw\2\B{k)\2 - J\\ab\\F ■

í k=l Pk l

Заметим, что последнее слагаемое }||Л5|||. не зависит от выбора рк, и задача минимизации сводится к минимизацией функции

f(p i.....Pn) = Í2^M\B^\2

К=1 К

при условии, что YlkPk = 1- Решение может быть получено с использованием стандартного метода множителей Лагранжа, который в данном случае приводит к необходимости минимизации функционала

п

F(pi, ...,рп) = f{Pi, ■ ■ ■ ,Рп) + - l) •

¡=1

Здесь ц - множитель Лагранжа. Решая систему уравнений §^- = 0, г = 1,2,... ,п п учитывая, что рп = 1 — Рк получим оптимальные значения вероятностей pk в виде (1.11).

1.2 Итерационные методы решения систем линейных алгебраических уравнений

Вначале рассмотрим кратко итерационные методы, к которым будет применен алгоритм построения разреженной оценки матрицы.

1.2.1 Метод, основанный на преобразовании спектра

Список литературы диссертационного исследования кандидат наук Моцартова, Надежда Сергеевна, 2013 год

Литература

[1] И. M. Соболь Метод Монте-Карло М., «Наука», 1968

[2] N. Metropolis and S. Ularri. The Monte Carlo method, J. Amer. Statist. Association, vol. 44 (1949), pp. 335-341.

[3] G. E. Forsythe, R. A. Leibler, "Matrix inversion by a Monte Carlo method MTAC, v.4, 1950, p. 127-129.

[4] J. H. Curtiss, "Sampling methods applied to differential and differetce equations Seminar on Scientific Computation, Nov. 1949, p. 87-109. IBM, New York.

[5] J. H. Curtiss, "Monte Carlo mathods for the iteration of linear operators Journalof Mathematics and Physics, vol. 32 (1954), p. 209-232.

[6] G. I. Marchuk, G. A. Mikhailov, M. A. Nazaraliev, et al, Monte Carlo Methods in Atmospheric Optics, Springer, 1992.

[7] Jerome Spanier and Ely M. Gelbard. Monte Carlo Principles and Neutron Transport Problems. New York, NY : Dover, 2008. - 256 p.

[8] R. Courant, К. Friedrichs, and II. Levy. Ueber die partiellen differenzenglei- chungen der mathemischen physik, Math. Ann. 100 (1928), pp. 32-74

[9] I. Petrovsky. Ueber das Irrfahrtproblem. Mathem. Annalen, 1934 (109), 1, 425-444.

[10] Михайлов Г.А. "Оптимизация весовых методов Монте-Карло "Наука, 1987.

[11] Mikhailov G.A. New Monte Carlo methods with estimating derivatives. Utrecht, VSP, 1995.

[12] Ермаков С.M. Методы Монте-Карло в вычислительной математике. Невский диалект, Бином. Лаборатория знаний. Санкт-Петербург 2009

[13] К. Sabelfeld, I. Shalimova., and A. Levykin. Random Walk on Fixed Spheres for Laplace and Lamé equations. Monte Carlo Methods and Applications, v. 12 (2006), N1, 55-93.

[14] Сабельфельд К. К. Векторные алгоритмы метода Монте-Карло для решения эллиптических уравнений второго порядка и уравнений Ламе Доклады АН СССР. 1982. - Т. 262, № 5. - С. 1076-1080.

[15] Kantorovih L.W. and Krylov V.I. Approximate methods of higher analysis. Interscience, New York, 1964.

[16] Кублановская В. H. Применение аналитического продолжения посредством замены переменных в численном анализе. Тр. Мат. ин-та им. В.А.Стеклова. - 1959.-Т.53, стр. 145-185.

[17] Sabelfeld, К.К. Monte Carlo Methods in Boundary Value Problems. Springer-Verlag, Berlin - Heidelberg - New York, 1991.

[18] Курбанмурадов 0., Сабельфельд К.К. Решения многомерных задач теории потенциала методом блуждания по границе. В сб.: Численные методы механики сплошной среды, Новосибирск, 1984, т. 15, с. 77-102.

[19] Sabelfeld, К.К. and Siinonov, N.A. Random Walks on Boundary for Solving PDEs. VSP, The Netherlands, Utrecht, 1994.

[20] О. А. Курбанмурадов, К. К. Сабельфельд, Н. А. Симонов Алгоритмы случайного блуждания по границе Под ред. Г. А. Михайлова; АН СССР, Сиб. отд-ние, Новосибирск ВЦ СО АН СССР 1989

[21] Сабельфельд К. К. Эргодический процесс блуждания по границе для решения задачи робена Теория и приложения статистического моделирования. Новосибирск: ВЦ СО АН СССР, 1989.- С. 118-123.

[22] D. Shia and Y. Hui. A Monte Carlo solution method for linear elasticity. International Journal of Solids and Structures, 37 (2000), issue 42, pages 6085-6105.

[23] M. Mascagni and N. Simonov.The Random Walk on the Boundary Method for Calculating Capacitance. Journal of Computational Physics, 2004, vol. 195, N 2, 465-473.

[24] Копылов Ю. H. Ветвящийся процесс случайного блуждания по границе для невыпуклых областей Теория и применения статистического моделирования. — Новосибирск: ВЦ СО АН СССР, 1991. - С. 52-57.

[25] Sabelfeld К. and Mozartova N. Sparsified Randomization Algorithms for large systems of linear equations and a new version of the Random Walk on Boundary method. Monte Carlo Methods and Applications, 2009, vol 15, issue 3, 257-284.

[26] A. Frieze, R. Kannan, and S. Vempala. Fast Monte Carlo algorithms for finding low-rank approximations. J. ACM, 51 (2004), no. 6, 1025-1041.

[27] Dimitris Achlioptas and Frank McSherry. Fast computation of low rank matrix approximations. Proceedings of the 33rd Annual Symposium on Theory of Computing, 2001, and Journal of the ACM (JACM), vol. 54, issue 2 (2007).

[28] Per-Gunnar Martinsson, Vladimir Rokhlin, Mark Tygert. A randomized algorithm for the approximation of matrices. Technical Report YALE U/DCS/TR-1352June 7, 2006

[29] Franco Woolfe, Edo Libei ty, Vladimir Rokhlin, Mark Tygert. A fast randomized algorithm for the approximation of matrices. Applied and Computational Harmonic Analysis, 25 (2008), 335-366.

[30] M. Kobayashi, G. Duplet, O. King, and H. Samukawa. Estimation of singular values of very large matrices using random sampling. Computers and mathematics with applications. 42 (2001), 1331-1352.

[31] W Eberly, E Kaltofen. On Randomized Lanczos Algorithms. International Conference on Symbolic and Algebraic Computation archive Proceedings of the 1997 international symposium on Symbolic and algebraic computation, pp. 176 - 183, 1997.

[32] Petors Drineas, Ravi Kannan, and Michael W. Mahoney. Fast Monte Carlo algorithms for matrices I: approximating matrix multiplication. SIAM J. Comput., vol.36 (2006), N 1, 132-157.

[33] Yu.V. Bulavsky and S.A. Terrmikov. Randomized method of successive approximations. In: Mathematical Methods in Stochastic Simulation and Experimental Desig. - St. Petersburg University Publishing House. 1996, pp. 64-68.

[34] Ledoux M. The concentration of measure phenomeon. Mathematical surveys and monographs, N89, AMS. 2000

[35] Levy P. Probl'eme concretes d'analyse fonctionelle. Paris: GauthierVillars, 1951.

[36] Halko N, Martinsson P , and Tropp J. Finding structure with randomness: stochastic algorithms for constructing approximate matrix decompositions. 2009. arXiv:0909.4061vl

[37] Vempala S. The Random projection method. DIMACS Series in Discrete Math, and Theoretical Comp. Science Vol. 65, American Math. Soc., Providence, R.I., 2004.

[38] Butnariu D. The expected-projection method: its behavior and applications to linear operator equations and convex optimization. Journal of Applied Analysis 1995; 1: 93-108.

[39] Diaconis P, Freedman D. Klifue K, Saloff-Coste L. Stochastic alternating projections. http://www-stat.stanfoul.edu kdkhare/altprojnew.pdf.

[40] W. B. Johnson and J. Lindenstrauss, Extensions of Lipschitz maps into a Hilbert space, Contemp. Math. 26 (1984), 189-206.

[41] Gene H. Golub and Charles F. Van Loan, Matrix Computations, Johns Hopkins University Press, Baltimore, Maryland, third ed., 1996.

[42] Cornelius Lanczos. An iteration method for the solution of the eigenvalue problem of linear differential and integral operators. Journal of Research of the National Bureau of Standards, vol. 45 (1950), N4, 255 - 282.

[43] Neil Muller, Lourenceo Magaia, B. M. Herbst. Singular Value Decomposition, Eigenfaces, and 3D Reconstructions. SIAM Review, 46 (2004), No. 3, pp. 518-545.

[44] Christos H. Papadimitriou, Prabhakar Raghavan, Hisao Tamaki, and Santosh Vempala. Latent semantic indexing: a probabilistic analysis. Journal of Computer and System Sciences, Volume 61, Number 2, October 2000 , pp. 217-235(19).

[45] Vladimir Rokhlin, Arthur Szlam, and Mark Tygert. A randomized algorithm for principal component analysis. SIAM J. Matrix Anal. Appl, 2009 - arxiv.org.

[46] G. W. Stewart. On the Early History of the Singular Value Decomposition. SIAM Review, Vol. 35, No. 4, 1993.

[47] Petros Drineas, Alan Frieze, Ravi Kannan, Santosh Vempala, and V. Vinay. Clustering Large Graphs via the Singular Value Decomposition, Machine Learning, 56 (2004), no. 1-3, 9-33.

[48] Petros Drineas, Ravi Kannan, Pass Efficient Algorithms for Approximating Large Matrices, Proceedings of the 14th Annual Symposium on Discrete Algorithms (Baltimore, MD), 2003, pp. 223-232.

[49] Gregory Beylkin and Martin J. Mohlenkam. Algorithms for numerical analysis in high dimension. SIAM .Journal on Scientific Computing, Volume 26 , Issue 6 (2005), 2133-2159.

[50] Drineas P, Drinea E, and Huggins P. An experimental evaluation of a Monte Carlo algorithm for singular value decomposition. Lecture notes in computer science 2003; 2563: 279-296.

[51] Petros Drineas, Ravi Kannan. Fast Monte Carlo algorithms for approximate matrix multiplication. Proceedings of the 42nd IEEE symposium on Foundations of Computer Science, Page: 452, 2001, ISBN:0-7695-1390-5.

[52] E. Cohen and D. Lewis. Approximating matrix multiplication for pattern recognition tasks. J. Algorithms, 30 (1999), N 2, 211-252.

[53] Edo Liberty, Franco Woolfe, Per-Gunnar Martinsson, Vladimir Rokhlin, and Mark Tygert. Randomized algorithms for the low-rank approximation of matrices. Yale Dept. of Computer Science Technical Report 1388.

[54] В. Д. Купрадзе, M. А. Алексидзе Метод функциональных уравнений для приближенного решения некоторых граничных задач, Ж. вычисл. матем. и матем. физ., 4:4 (1964), 683-715

[55] М.А. Алексидзе, К вопросу о практическом применении одного нового приближенного метода. Дифференциальные уравнения 2 (1966), 1625-1629. MR0203975(34:3822)

[56] Barnett А.Н. and Betcke Т. Stability and Convergence of the Method of Fundamental Solutions for Heimholt z problems on analytic domains. Journal of Computational Physics 2008; 227, issue 1 1, 7003-7026.

[57] Trefftz E. Ein Gegenstück zum Ritzschen Verfahren, 2er Intern. Kongr. für Techn. Mech., Zürich, 1926, 131—137.

[58] Betcke T. and Trefethon L.N. Reviving the method of particular solutions, SIAM Rev 2005; 47: 469U-491

[59] Kita E and Kamiya N. Trefftz method: an overview. Advances in Engineering software 1995; 24: 3-12.

[60] Ramachandran P. Method of fundamental solutions: singular value decomposition analysis. Communications in Numerical Methods in Engineering 2002; 18: 789-801.

[61] Smyrlis Y. The method of fundamental solutions: a weighted least-squares approach. BIT 2006; 46: 1, 163U194.

[62] Cao Y, Schultz W and Beck R. Three-dimensional desingularized boundary integral methods for potential problems. International Journal for Numerical Methods in Fluids 1991; 12: 785-803.

[63] K. Sabclfcld, N. Mozart ova Sparsified Randomization algorithms for low rank approximations and applications to integral equations and inhomogeneous random field simulation. Mathematics and Computers in Simulation, Volume 82, Issue 2, October 2011, Pages 295-317.

[64] K. Sabelfeld, N. Mozartova Stochastic boundary collocation and spectral methods for solving PDEs. Monte Carlo Methods and Applications.-202.-Vol. 18 ДОЗ.-Р.217-264.

[65] Ермаков С. M., Михаилов Г. А. Статистическое моделирование.— Москва: Наука, 1982.

[66] Владимиров B.C. Уравнения математической физики. М.: Наука,1981. 512 с.

[67] Ю. В. Воробьев, "Случайный итерационный процесс", Ж. вычисл. матем. и ма-тем. физ., 4:6 (1964), 1088-1093

[68] Ю. В. Воробьев, "Случайный итерационный процесс в методе переменных направлений", Ж. вычисл. матем. и матем. физ., 8:3 (1968), 663-670

[69] А. Н. Коновалов, "Введение в вычислительные методы линейной алгебры "Наука, Новосибирск,1993.

[70] Ermakov S.M., Wagner W. Monte Carlo difference schemes for the wave equations. Monte Carlo Methods and Applications, 8 (2002), N 1, 1-29.

[71] Rosea V, Leitao V. Quasi-Mote Carlo mesh-free integration for meshless weak formulations. Engineering Analysis with Boundary Elements 2008; 32: 471-479.

[72] Cisilino A.P., Sensale B. Application of a simulated annealing algorithm in the optimal placement of the source points in the method of the fundamental solutions. Computational Mechanics 2002; 28: 129-136.

[73] Nir Ailon and Bernard Chazelle. Approximate nearest neighbors and the fast JohnsonLindenstrauss Transform. Proceedings of the thirty-eighth annual ACM symposium on Theory of computing, 2006, 557-563.

[74] Dimitris Achlioptas. Database-friendly random projections: Johnson-Lindenstrauss with binary coins. Journal of' Computer and System Sciences 66 (2003) 671-687.

[75] Eckhart С and Young G. A principal axis transformation for non-Hermitian matrices, Bulletin of the American Mathematical Siciety 1939; 45: 118-121.

[76] Vladimir Rokhlin, Artlms Szlam, and Mark Tygert A randomized algorithm for Principal Component, Analysis UCLA Computational and Applied Math. Technical Report 08-60 September 5, 2008

[77] Pearson K., On lines and planes of closest fit to systems of points in spacc, Philosophical Magazine, (1901) 2, 559-572

[78] P. Kramer, O. Kurbanniuradov, and K. Sabelfeld. Comparative Analysis of Multiscale Gaussian Random Field Simulation Algorithms. Journal of Computational Physics, v. 226 (2007), 897-924.

[79] Прохоров Ю.В., Розанов IO.А. Теория вероятностей. Основные понятия, предельные теоремы, случайные процессы. М.: Наука, 1973

[80] Г. А. Михайлов, А. В. Войтишек Численное статистическое моделирование. Методы Монте-Карло. Академия, 2006

[81] Rokhlin, V., Rapid solution of integral equations of classical potential theory, J. Сотр. Phys., 60 (1985). 187-207.

[82] Ермаков С.М. Сипии A.C. Новая схема метода Монте-Карло для решения задач мат. физики, ДАН СССР, 285, N3, 1985.

[83] Jorge Buescu and A.C. Paixao. Eigenvalues of positive definite integral operators of unbounded intervals. Positivity, 10 (2006), 627-646.

[84] J. Buescu, F. Garcia, I Lourtie, A. Paixao. Positive definiteness, integral equations and Fourier transfoims. Journal of Integral equations and Applications, 16 (2004), 33-52

[85] I.M. Novitsky. Representation of kernels of integral operators by bilinear series. Siberian Math. J., 25 (1981), N3, 774-778. Translated from the Russian: Sibirsk. Mat. Zh., 25 (1984), N5. 114-118.

[86] A.M. Yaglom. Con elation theoiy of stationary and related random functions I. Basic results. Springer-Veilag. New York - Heidelberg - Berlin, 1987.

[87] K.K Phoon, H.W Huang, S.T. Quek. Simulation of strongly non-Gaussian processes using Karhunen-Loeve expansion. Probabilistic engineering Mechanics, 20 (2005), 188-198.

[88] B.B. Воеводин, Ю.А. Кузнецов Матрицы и вычисления. Москва "Наука 1984.

[89] Walker A.J. New fast method for generating discrete random numbers with arbitrary friquency distribution. Elect ionic Letters, 1974; 10; 127-128.

[90] Пригарин С. M. Методы численного моделирования случайных процессов и полей. - Новосибирск: 1Ьд. ИВМиМГ СО РАН, 2005.

[91] Кроновер Р. М. Фрак™ ды п хаос в динамических системах. М., 2000

[92] Strang G. The fundamental Theorem of linear algebra. The American Mathematical Monthly 1993; 100: 848-855.

[93] Chein-Shan Liu Impioving the Ill-conditioning of the Method of Fundamental Solutions for 2D Laplai e Equation. CMES,vol.851, no.l,pp.l-17,2009

[94] Katsurada, M. (1989). A mathematical study of the charge simulation method II. J. Fac. Sei., Univ. of Tokvo. Sect. 1A, Math. 36, 135-162.

[95] Kitagawa, T. (1988). On the numerical stability of the method of fundamental solution applied to the Diiiihlet pioblem. Japan J. Appl. Math. 5, 123-133.

[96] E. Trefftz, Ein Gegenstück zum Ritzschen Verfahren, 2 Intern. Kongr. für Techn. Mech., Zürich, 1926. pp 131-137.

[97] M.A. Алексидзе. К вопросу о практическом применении одного нового приближенного метода. Дифференциальные уравнения 2 (1966), 1625-1629. MR0203975(34:3822)

[98] Sabelfeld K.K. Expansion of random boundary excitations for some elliptic PDEs. Monte Carlo Methods and Applications 2007; 13: 403-451.

[99] F.E. Browder, Approximation by solutions of partial differential equations, Amer. J. Math.84(1962), 134-160. MR0178247(31:2505)

[100] A. H. Cho, M. A. Golberg, A. S. Muleshkov, and X. Li, Trefftz methods for time dependent partial diff'oiinitial equations, Comput. Mat. Cont. 1(2004), no. 1, 1-37.

[101] A. Doicu, Y. Eremin. and T. Wriedt, Acoustic and Electromagnetic Scattering Analysis using Discrete Sources, Academic Press, New York, 2000.

[102] G. Fairweather, A. Karageorghis, and P. A. Martin, The methos of undamental solutions for scattering and radiation problems, Eng. Anal. Bound. Elem. 27 (2003), 759-769.

[103] О. В. Бесов Тригонометрический ряд Фурье. Учебно-методическое пособие (для студентов 2-го курса). МФТИ. М., 2004.

[104] Jack F. Douglas, Husan-Xiang Zhou, and Joseph B. Hubbard Hydrodinamic friction and the capacitance of arbitrary shaped objects. Physical Review E, volume 49, number 6, June 1994.

[105] P.C. Hansen, The L-curve and its use in the numerical treatment of inverse problems. Advances in Computational Bioengineering, 2000, WIT Press

[106] Ландау Л. Д., Лифшпц Е. М. Теоретическая физика. В 10-ти т. Т. УШ.Электродипа.мпка < и.юшных сред: Учеб. пособие. — 4-е изд., исир. и доп. — М.; Наука. Гл. род. фи j.-мат. лит., 1987.