Г.М. Кобельков - Курс лекций по численным методам (1160467), страница 15

Просмтор этого файла доступен только зарегистрированным пользователям. Но у нас супер быстрая регистрация: достаточно только электронной почты!

Текст из файла (страница 15)

Метод балансаПро него всё равно ничего доказано не было. Он есть в книжке Бахвалова – Жидкова – Кобелькова на страницах 464–466.Точнее говоря, на странице 464 всё начинается с последнего абзаца и заканчивается формулой (8) на странице 466. Приведённаятам схема и есть «ответ».Ну ладно, настало время написать этот бред.Рассмотрим уравнение′′ k(x)y (x) − p(x)y(x) = f (x),y(0) = a,y(l) = b.(126)Здесь p(x) > 0, а k(x) > k0 > 0, причём k(x) ∈ C3 , а p, f ∈ C2 за исключением конечного числа точек.Казалось бы, можно раскрыть производную произведения, но этого как раз делать не нужно, а то сходитьсябудет совсем плохо.Обозначим w(x) := k(x)y ′ (x).

Будем считать, что узлы сетки xi не попадают на точки разрыва нашихфункций. Имеем(ky ′ )′(ky ′ )′xnxn+ 1− (ky ′ )′2≈xn− 12h=k(xn+ 12 ) ·yn+1 −ynh− k(xn− 12 ) ·hyn −yn−1h.(127)Но это нам ничего не даст. А именно, данная схема будет иметь порядок O(h), и аппроксимация будет порядкаO(1). Это полная лажа, надо придумывать что-то ещё (кому нужна аппроксимация с точностью до константы?).Таким образом, это уравнение доставляет пример того, когда схема «аппроксимация + устойчивость» не даётсходимости.Предлагаемый ниже метод называется методом баланса и позволяет построить схему порядка h2 . Правда,доказать, что она имеет именно такой порядок, весьма непросто, и мы этого делать не будем.Ну-с, приступим: имеемxn+ 1Z 2w(xn+ 21 ) − w(xn− 12 ) =p(x)y(x) + f (x) dx.(128)xn− 12Далее, y(x) = y(xn ) + O(h), если x − xn = O(h). Тогдаxn+ 1Z2xn− 12p(x)y(x) + f (x) + O(h) dx = h pn yn + f n + O(h2 ),55(129)гдеxn+ 1pn :=1hZxn+ 12p(x) dx,f n :=Z1hxn− 12f (x) dx.(130)xn− 122Делить на h в этих формулах не так уж страшно, так как в пределе при h → 0 получим pn и fn соответственно.Рассмотрим схемуw(xn+ 12 ) − w(xn− 12 )− pn y(xn ) = f n + O(h).(131)hА вот теперь начинаем хитрить.

Введём новую переменнуюt(x) :=Zxdx,k(x)tn :=Zxndx.k(x)(132)00Тогдаw(xn− 12 ) = k(x)dydx=xn− 1dydt2Введём ещё две функции1k + :=h=tn− 12xZn+1xny(xn ) − y(xn−1 )+ O(tn − tn−1 ).tn − tn−11k − :=hdx,k(x)Zxnxn−1dx.k(x)(133)(134)Окончательно, схема имеет вид1hyn+1 − ynyn − yn−1−hk +hk −− pn yn = f n .(135)Почему-то в конспекте чёрточки над fn не было.

Что-то мне подсказывает, что это больше похоже на опечатку Арушаняна.Иначе какого чёрта мы эту функцию вводили?!Далее читателю предлагается поверить в тот факт, что эта схема имеет порядок O(h2 ). Можно про этопрочесть в книжке БЖК, но на лекциях этого не было, да и в программу не вошло.Кстати, как вычислять те интегралы, которые тут фигурируют в большом количестве, не было сказано ни слова. Видимо, поквадратурным формулам.

С другой стороны, если использовать только узлы xi , то будет как-то не очень точно, а если использоватьбольше узлов, то это фактически означает уменьшение шага h.5.9. Метод конечных элементов (проекционный метод)Тот метод, который мы рассмотрим сейчас, является очень частным случаем метода конечных элементов.Точнее говоря, даже сам метод конечных элементов не есть чёткий алгоритм решения той или иной дифференциальной задачи, а скорее представляет собой некоторую идею (или технологию). Ничего принципиальнонового он в себе не несёт.Ключевая идея метода (т. е.

её вольная интерпретация автором конспекта) излагается ниже мелким текстом.Философское отступление о разных базисах в функциональных пространствахРазложение функций по какому-нибудь базису было придумано очень давно. Но вот какой базис брать — этот вопрос частобывает не праздным. Например, в некоторых случаях хорош тригонометрический базис, однако в задачах сжатия информацииуже давно поняли, что он не очень удачен. А вот зато базис Хаара (или как там по-русски пишут его фамилию?), о которомречь пойдет ниже, себя оправдывает. Вообще, идея использования в качестве базиса функций-«всплесков» («wavelets»), то естьфункций, которые имеют ярко выраженный экстремум, оказалась очень плодотворной, например, в задачах сжатия.

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

Этот метод можно использовать и длямногомерных задач, и даже для областей произвольной формы, но вид получаемой системы будет зависеть от способа нумерациинаших «конечных элементов», то есть кусочков области, где сосредоточены носители базисных функций.Не претендуя на оригинальность и новизну, рассмотрим до боли знакомую одномерную задачу′′−u + pu = f,u′ + αu = 0,0u(l) = 0.56(136)Но мы будем искать не классическое решение, а обобщённое. А именно, перейдём к интегральной задаче,для начала спроецировав всё на сетку на отрезке [0, l] с шагом h и N узлами:a(u, v) :=Zl(u′ v ′ + puv) dx + αu(0)v(0) = (f, v),(137)0где(f, v) :=N−1Xhui vi .(138)i=1Будем говорить, что мы нашли приближённое решение, если данное уравнение выполняется при любой функции v. От «настоящего» обобщённого решения оно отличается лишь тем, что мы рассматриваем функции изподпространства сеточных (На самом деле не совсем сеточных, а кусочно-линейных, продолженных линейнымобразом между узлами) функций.Введём в нашем пространстве вот такой прикольный базис:(1, x = xi ;vi =(139)0, x 6= xi .RНоситель k-й функции есть отрезок [xi−1 , xi+1 ].

Очевидно, что если |i − j| > 2, то vi (x)vj (x) dx = 0. Отсюдаследует, что при разложении по наших функций по базису {vj } мы получим ленточную (точнее даже трёхдиагональную матрицу). Как она получается, наверное, читатель догадается сам.Для формалистов: на странице примерно 479 книжки БЖК (и далее ещё на 10 страницах) излагается «варианционно-разностныйметод Ритца», в котором делается всё примерно то же самое.

Там же (на стр. 480) явно вычислены коэффициенты разложения.Пункты 4, 5 и 6 данного параграфа представляют собой обобщение рассмотренного нами случая (для сведения к частному случаюнадо в формуле (13) на странице 481 взять k(x) ≡ 1. В пункте 6 рассматривается в точности данная задача.Тут была ещё теорема о том, что этот метод вообще сходится, и что матрица получится положительно определённой. Нодоказательство было неправильное, а в книжке написано хреново.

Но кому интересно, почитайте пункт 5 на странице 484.5.10. Интегральные уравнения второго родаБудем рассматривать интегральные уравнения следующего видаZ1K(x, y)u(y) dy = f (x).(140)0λZ1K(x, y)u(y) dy + u(x) = f (x).(141)0Они называются интегральнымиуравнениями первого и второго рода соответственно.RRПредполагаем, чтоk dx dy < ∞.Эта задача, вообще говоря, является некорректно поставленной. Нужно определить, на каком классе функций мы ищем решения.Первый метод решения состоит в приближении конечномерными интегральными операторами. Можно разложить ядро в ряд по системе функций ϕj :Kn =nXaij ϕi (x)ϕj (y).(142)i,j=1Если система полная, то kKn − Kk → 0.

Подынтегральное выражение и правую часть тоже разложим в ряды:un =nXci ϕi (x).(143)nXdi ϕi (x).(144)i=1fn =i=157Можно доказать, что un сходятся к точному решению u. Тогда уравнение запишется в видеλnXi,j=1aijZ1ϕi (x)ϕj (x)0nXcl ϕl (y) dy +nXi=1l=1ci ϕi (x) =nXdi ϕi (x).Рассмотрим это уравнение при каждом фиксированном i. Сумму по l можно вытащить:nn Z1XXλaij ϕj (x)ϕl (y) dy  ui + ui = fi .i,j=1(145)i=1(146)l=1 0Получилась система линейных уравнений.

Одно плохо: даже если система ортогональна, то матрица получаетсяполная. То есть придётся решать методом Гаусса, по-честному. А вот когда уравнение имеет вид свёртки, тостановится легче: матрица становится «полосатой». При этом можно применить быстрое преобразование Фурье,и тогда свёртка перейдёт в произведение, и всё будет хорошо.Второй способ — замена интеграла квадратурной формулой.

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

Характеристики

Список файлов лекций