Метод наибыстрейшего спуска

Федеральное агентство по образованию

Государственное образовательное учреждение высшего профессионального образования

Тульский государственный университет

Кафедра прикладной математики и информатики

ПОЯСНИТЕЛЬНАЯ ЗАПИСКА

к курсовому проекту по курсу

«ЧИСЛЕННЫЕ МЕТОДЫ»

на тему

«Метод наибыстрейшего спуска»

Автор работы

студент гр. 530432

Руководитель работы

Гаев А.В.

Тула 2006



Реферат

В курсовой работе рассмотрено решение интегральных уравнений методом наибыстрейшего спуска.

На примерах показано использования метода в вычислительной практике.

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

Ключевые слова:

Объем работы страниц. Список использованной литературы содержит источника (ов). Работа содержит приложения.



Введение

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

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

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



1. Основные понятия и сущность метода

Метод наибыстрейшего спуска

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

1.1 Теорема о минимуме квадратичного функционала и ее следствие

В начале гл. 8 мы уже отмечали, что наибольший интерес для нас будут представлять уравнения вида

(1)

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

в G,(2)

,(3)

где Г--граница области G. Оператор -- рассматривался на линеале М₁функций, принадлежащих С⁽²⁾ (G) и удовлетворяющих условию (3) на Г. Мы показали, что М₁ плотен в гильбертовом пространстве L₂(G) и что оператор -- симметричен и положителен на М₁. Было упомянуто, что позднее будет доказана положительная определенность этого оператора.

В дальнейшем будет очень часто встречаться следующая постановка задачи: задано некоторое гильбертово пространство (не обязательно пространство L₂(G), как в случае рассмотренного примера; при решении некоторых задач целесообразно рассматривать более общие пространства или пространства с существенно отличными свойствами). В этом пространстве заданы линеал D_A, а также плотный в H и положительный (или положительно определенный) оператор А, который отображает D_A в Н. Далее, задан некоторый элемент . Мы ищем такой элемент, который

удовлетворяет в H уравнению

Au = f.(4)

Записывая (4), мы имеем в виду, что уравнение Au = f удовлетворяется в рассматриваемом пространстве H, т. е. что Аи--f-- нулевой элемент H. Например, если H=L₂(G), то (4) означает выполнение равенства Au = f почти везде в области G, т. е. везде, кроме, может быть, точек, образующих множество меры нуль.

Сначала мы покажем, что если А положителен в D_A, то уравнение (4) может иметь не более одного решения в H. Прежде чем перейти к доказательству, повторим вкратце некоторые определения.

Пусть D_A-- линеал, плотный в гильбертовом пространстве H. Оператор А, отображающий D_A в H, называется симметричным в D_A, если он линеен в D_A и если равенство

(Аи, v) = (u, Av)(5)

справедливо для любой пары элементов . Оператор А называется положительным в D_A, если он симметричен в D_A и если для всех выполняется неравенство

(Аи,u) 0,(6)

причем

(Аи, u) = 0=>u = 0 в D_A.(7)

Теорема 1. Если А--положительный оператор в D_A, то уравнение (4), в котором , имеет не более одного решения в H.

Доказательство. Допустим, что существуют два решения , т. е. пусть в Я выполнены равенства

(8)

(9)

Вычитая (8) из (9) и принимая во внимание, что А положителен и поэтому заведомо является линейным оператором, так что Аи_г--Аи_х = А*(н₂--щ), получаем А (и₂--uj = 0 в Я. Умножая^х) последнее уравнение на элемент и₂--u₁ЈD_A, получаем (А (и₂~и₁), и₂--M_t) = 0, откуда, согласно (9.7), и₂--и_х = 0 в D_A, или, что то же самое, и_г -- и_г в D_A, что и требовалось доказать.

Следующая теорема чрезвычайно важна для наших целей.

Теорема 2. (Теорема о минимуме квадратичного функционала).

Пусть А-- положительный оператор в D_A, fH. Пусть уравнение (4) имеет решение u_oЈD_A, т. е. выполняется

Au_o = f в Н, u_o D_A.(10)

Тогда квадратичный функционал

Fu = (Au, u) -- 2(f, и)(11)

принимает на и₀ минимальное значение в D_A, т. е. для всех uЈD_A имеет место Fu^Fu₀, причем Fu = Fu₀ только при и = и₀. С другой стороны, пусть среди всех элементов uЈD_A функционал (9.11) принимает свое минимальное значение на элементе и₀. Тогда и_й является решением уравнения (9.4) в Н, т. е. справедливо (9.10).

Доказательство. Во-первых, очевидно, что функционал Fu определен при всех u^D_A.

1. Пусть уравнение (9.10) удовлетворяется в Н при и = и₀, так что = Аи₀. Подставляя это в (9.11), получаем при uЈD_A Fu = (Au, и) -- 2(Ащ, и). Однако легко видеть, что

Fu = (Au, и)--2(Аи₀, и) = (Аи, и) -- (Аи_и, и)--(и, Аи₀) -- = (Au, u)--(Аи₀, и) --(Аи, и_й) = {Аи, и)--(Аи₀, и)--(Аи, и_о) + (Аи_о, и_о)--(Аи_о, и_о) = = (A(u--u₀), и -- и_о) -- (Аи_о, и₀).

(Мы использовали, с одной стороны, симметрию скалярного произведения, т. е. равенство (Au₀, и) = (и, Аи₀), и, с другой стороны, симметрию оператора Л, согласно которой (и, Л«_о) = (Аи, «„).) Следовательно, если выполнено Ли_о = , то

Fu = (A(u -- u₀), u -- u₀) -- (Au₀,u₀).(12)

Член (Лм₀, и₀) в (12) не зависит от и и, таким образом, остается постоянным; затем, оператор Л положителен по предположению, так что для первого члена в правой части (12) мы имеем (Л(и -- и₀), и -- «₀)^0 при любом uЈD_A, причем (А (и -- и_а), и -- м_о) = О в том и только том случае, если и -- и_о = 0 в D_A. Таким образом, из (12) следует, что Fu ^ Fu₀ при любом и Ј D_А и при этом Fu -- Fu_a, только если и = и_и в D_A. Следовательно, если удовлетворяется уравнение Лы_о = , то функционал Fu принимает свое минимальное значение в D_A в точности на элементе и = и_а. Это завершает доказательство первого утверждения теоремы 9.2.

2. Пусть функционал Fu принимает свое минимальное значение в D_A на элементе и₀. Это значит, что если мы выберем произвольный элемент v g D_A и произвольное вещественное число t (так что u_o-{-tv, разумеется, также принадлежит D_A), то

F(u + tv .(13)

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

Следовательно, имеем

(14)

Поскольку ы₀ g Ј) и f€#--фиксированные элементы, то из (14) очевидно, что для фиксированного и g D_A функция F (и₀ - tv) является квадратичной по переменной. Из (13) следует существование у этой функции локального минимума при =0. Это означает, что ее первая производная равна нулю при t = 0,

(15)

или, согласно (14), что

(16)

т. е. что

(Au_o -- f,v) = 0.(17)

Элемент vЈD_A можно выбрать произвольным образом (использовалось лишь то, что он фиксирован), так что от условия (13) можно прийти к уравнению (17) при любом vЈD_A. Линеал D_A по предположению плотен в Н. Следовательно, из (17) следует, что элемент Лм₀ -- ортогонален в Н всем элементам множества, являющегося плотным в Н. Таким образом, в соответствии с теоремой 6.18

Аи₀ --f = 0 в H,,(18)

т. е. и₀ является в Н решением уравнения Au_o -- f, что и требовалось доказать.

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

(19)

с граничными условиями

и (0) = и (I) =0, u'(0) = u'(l) = 0-(20)

Предположим, что функции Е (х), I (х) и их первые и вторые производные, равно как и функция q(x), непрерывны в интервале [0, /] и что



Е(х)>0, 1(х)>0 в [0, /].(9.21)

Уравнение (19) можно интерпретировать как дифференциальное уравнение изгиба балки длины / с модулем упругости Е (х), с моментом инерции / (х) поперечного сечения относительно оси изгиба (эта ось носит также название нейтральной) и с распре-еленной нагрузкой q(x). Условия (20) означают, что балка жестко закреплена («заделана») на обоих концах.

Если Е и / постоянны в [0, /], то (19) можно записать в виде Е1и⁽⁴⁾ = q

Выберем H = L₂ (0, /) и обозначим через D_A линеал функций, которые непрерывны со своими производными до четвертого порядка включительно в интервале [0, I] и удовлетворяют условиям (20). Определим оператор А на D_A выражением

(22)

Задача (9.19), (9.20) может быть записана тогда в виде

Au = q.(23)

Легко доказать, что А положителен на D_д. Сначала докажем, что он симметричен в D_A. Для любых uЈD_A, vЈD_A получаем

(Выражения [(Ј/«")' v]i и [EIu"v']l равны нулю, поскольку вследствие (20) мы имеем v(0) -- v(l) = 0 и v' (0) = v' (I) = 0 соответственно.) Аналогично получаем i

.(25)

Из (24) и (25) немедленно следует симметрия А на D_A. Далее, из (24) получаем

(26)

то при всех .

Ввиду (21) из (26) следует, что (Аи, u)>0 при любом . Кроме того, если (Аи, u) = 0, то из (26) и (21) следует, что и"(х) = 0 на [0, l], откуда и(х) = ах + b на [0, /] и вследствие (20) и(х) = 0 на [0, /].

Это доказывает, что А является положительным на D_A.

Функционал Fu = (Au, и)--2(q, и) имеет в данном случае вид

(27)

и при заданном прогибе uЈD_A равен удвоенной общей потенци-альной энергии Lu рассматриваемой балки:

(28)

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

(29)

Пусть и_о(х)--решение (23), т. е. пусть и

(E]ul)" = q.(30)

Подставляя q из (30) в (27), получаем

оо

Интегрируя по частям аналогично (24), имеем (для заданного и при любом )

так что

(31)

Поскольку (u"--uo)²0 и в то же время выполняется (21), из (31) видно, что Fu минимально в D_A, только если и"--и'^=0 в D_A, или, вследствие (20), только если и = и₀ в D_A. Таким образом, мы можем сформулировать наше первое утверждение, которое соответствует первому утверждению теоремы 2: если является решением уравнения (23), то функционал Fu принимает минимальное значение в D_A в точности на и₀.

Пусть, обратно, -- элемент, минимизирующий в D_A функционал (27). Пусть -- произвольный элемент из D_A и пусть t -- произвольное вещественное число. В этом случае условие минимума



Имеет вид

о

или

.(32)

Деля (32) на 2 и интегрируя первый член по частям (см. (24)), можно свести рассматриваемое условие к виду

Это равенство должно выполняться при всех . Однако линеал D_A плотен в L₂ (0, /), откуда (так же, как и в (18))

(33)

в L₂ (О, I). Но функции, входящие в (33), так же как и их производные, непрерывны в [0, /] по предположению, так что (33) выполняется и в обычном смысле. Таким образом, мы приходим ко второму утверждению, снова соответствующему утверждению теоремы 2: если функция и_о(х), принадлежащая D_A (и удовлетворяющая тем самым условиям (20)), минимизирует в D_A функционал (27), то она удовлетворяет условию (33), т. е. является решением дифференциального уравнения для прогиба рассматриваемой балки.

Оба только что приведенных утверждения представляют собой не что иное, как принцип минимума потенциальной энергии, известный из теории упругости и примененный к рассматриваемой балке. Они показывают тесную связь между теоремой 2 и вариационными принципами механики. Конечно, теорема 2 рассматривает уравнение Au = f безотносительно к его физической интерпретации или к интерпретации соответствующего функционала Fu.

1.2 Понятие энергетического пространства

Пространство называется энергетическим пространством оператора А элементов - энергетической нормой элемента u.

1.3Минимизирующая последовательность и ее сходимость

Пусть F (и) -- некоторый функционал, значения которого ограничены снизу. Тогда существует конечная нижняя грань

d = inf F (u).

Определение. Говорят, что последовательность {и_п}, элементов из Нд является минимизирующей последовательностью для функционала F (и), если выполняется соотношение:

Пусть А -- положительно-определенный оператор и щ -- обобщенное решение уравнения Au = f.

Функционал энергии приводится к виду:

(18)

тогда минимизирующая последовательность для этого функционала характеризуется равенством:

(19)

Из этого замечания вытекает следующая теорема.

Теорема 3. Если оператор А -- положительно-определенный, то всякая последовательность {и_п}, минимизирующая для функционала энергии F (и), сходится по энергии к обобщенному решению соответствующего уравнения Аи = /, т.е. выполняется соотношение:

(20)

И наоборот, если для последовательности {и_п} выполняется соотношение (19), то {и_п} является минимизирующей последовательностью функционала энергии F (и).

Доказательство.

1.Если {и_п} -- минимизирующая последовательность, то выполняется (3.19) и из (3.18) следует, что

,

т.е. выполнено соотношение (20)

2. Если выполняется (20), то из (18) следует

т.е. в На (по энергии). Одновременно Н.

Теорема 3. доказана

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

Существует ряд способов построения -- последовательности:

метод ортонормированных рядов,

метод Ритца,

метод Галёркина,

метод наименьших квадратов,

метод Куранта,

метод быстрейшего спуска.

Метод наименьших квадратов применим и к таким уравнениям вида (1), в которых оператор А не положителен и даже не симметричен. Но в случае положительного А метод наименьших квадратов приводит к -- последовательности.

Отметим, что приведенная теорема о сходимости минимизирующей последовательности для квадратичного функционала F (и) верна и для только положительного оператора А, если при этом уравнение Аи = / имеет решение , принадлежащее энергетическому пространству



2. Разработка алгоритма приближенного решения обыкновенного интегрального уравнения

Пусть А-- положительно определенный оператор в гильбертовом пространстве Н. Пусть , и пусть является в Н решением уравнения

Как мы знаем , элемент и_а минимизирует в H функционал

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

Функционал на каждомпринимает определенное значение. Геометрически можно. представить этот функционал в виде некоторой «поверхности» над пространством , которая имеет «наименьшую координату» в точности при. Отправляясь из начальной точки, мы теперь хотим как можно быстрее достигнуть точки, для которой эта «координата» минимальна. Мы движемся по поверхности именно в «направлении наибыстрейшего спуска».

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

Однако мы имеем

Тогда требования (1) и (2) означают (как легко видеть), что

(поскольку если норма v_x устанавливается из (1), то скалярное произведение правой части (5) будет максимально в точности при . В то же время из (4) и (3) видно, что для такого выбора функционал будет принимать минимальное значение вдоль «прямой линии»при

В качестве второго приближения решения берется элемент

Если , мы продолжаем дальше в том же духе, т. е. мы конструктируем, подобно (6) и (8), элементы

где

и т. д. Так получается последовательность элементов

Пусть т, М -- положительные постоянные, такие что

Тогда можно видеть, что последовательность сходится к решению u₀ данного уравнения как в , так и в пространствесо скалярным произведением , и что справедлива оценка

2.1 Постановка задачи. Разработка алгоритма ее решения

В пространстве рассмотрим интегральное уравнение

где Заметим прежде всего, что оператор

симметричен в, поскольку

Далее,

Но, согласно неравенству Шварца, для любых двух функцийвыполняется

Следовательно, полагая

получим

Еще раз применив неравенство Шварца, можно получить

или, другими словами,



Подстановка в (11) дает

откуда прежде всего вытекает положительная определенность оператора (10) и в, поскольку ясно, что

Если метод наибыстрейшего спуска применяется к решению уравнения (9), то процесс получается сходящимся.

В качестве численного примера возьмем случай на промежутке. Мы получаем уравнение

Поскольку коэффициент 0.1 при втором члене «мал», очевидно, что удобно будет в качестве первого приближения решения взять правую часть данного уравнения, т. е. функцию

В соответствии с (6) мы определим сначала

Далее, согласно (7), имеем

Однако

так как

Подстановка дает

Следовательно, согласно (8), второе приближение решения имеет вид

Согласно (12), это второе приближение в по меньшей мере в раз (т.е. по меньшей мере в четыре раза)

лучше, чем первое приближение

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

Интегрируя и сравнивая коэффициенты при линейно независимых функциях и, получаем

откуда для неопределенных коэффициентов. находим значения

Следовательно, согласно (13), решение уравнения (14) имеет вид

Так как числа и малы по сравнению с числами 1 и 0.2 соответственно, (15) дает действительно очень хорошее приближение решения (16).

3. Численная реализация задачи на ЭВМ

Листинг программы.

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

#include <math.h>

#include <conio.h>

#include "matrix6.cpp"

double norma(TMatrix<double> &x)

{

double res=0;

int m=x.getSizeRow(), n=x.getSizeCol();

for(int i=1; i<=m; i++)

{

for(int j=1; j<=n; j++)

{

res+=x(i,j)*x(i,j);

}

}

return sqrt(res);

}

//решение методом наибыстрейшего спуска

TMatrix<double> linearSolveFastDescent(TMatrix<double> &a,//матрица коэфф.

TMatrix<double> &b,//правая часть

double sigma, //погрешность

long max_step) //максю кол-во шагов

{

int n=a.getSizeRow();

TMatrix<double> x(b), a_tr(a.getTranspose()), r(n,1), s(n,1);

for(long k=0; k<max_step; k++)

{

r=b-a*x;

s=a_tr*r;

x+=(r%r)/(s%s)*s;

if(norma(r)<sigma)

{

break;

}

}

cout<<"число шагов: "<<k<<endl;

if(k==max_step)

{

cout<<"заданная точность не достигнута "<<endl;

}

return x;

}

void main(void)

{

char ch;

TMatrix<double> a_b;

double sigma=.0000001;

long max_step=10000;

cout.setf(ios::showpoint);

for(;;)

{

clrscr();

cout<<"Решение интегрального ур-ния"<<endl

<<"Метод наибыстрейшего спуска "<<endl<<endl

<<"1 - ввести данные из файла"<<endl

<<"2 - ввести погрешность<<endl

<<"3 - ввести максимальное кол-во шагов<<endl

<<"4 - решить интегральное ур-ние"<<endl<<endl

<<"Esc - выход"<<endl<<endl<<endl;

do

{

ch=getch();

}

while(ch!='1' && ch!='2' && ch!='3' && ch!='4' && ch!=27);

switch(ch)

{

case '1':

{

clrscr();

char str[20];

cout<<"введите имя файла:"<<endl;

cin>>str;

ifstream in_file(str);

if(in_file)

{

int n;

in_file>>n;

a_b.setSize(n,n+1);

a_b.readArray(in_file);

a_b.writeArray(cout);

}

else

{

cout<< не могу открыть файл\"”<<str<<"\""<<endl;

}

getch();

break;

}

case '2':

{

clrscr();

cout<<"введите погрешность"<<sigma<<"): ";

cin>>sigma;

break;

}

case '3':

{

clrscr();

cout<<"введите максимальное количество шагов"<<max_step<<"): ";

cin>>max_step;

break;

}

case '4':

{

TMatrix<double> x;

clrscr();

cout<<"Решение интегрального ур-ния“"<<endl;

cout<<"Расширенная матрица системы"<<endl;

a_b.writeArray(cout);

cout<<"Погрешность: "<<sigma<<endl;

cout<<"Максимальное кол-во шагов: "<<max_step<<endl;

int n=a_b.getSizeRow();

TMatrix<double> a(a_b.getPart(1,1,n,n)), b(a_b.getCol(n+1));

x=linearSolveFastDescent(a,b,sigma,max_step);

cout<<"решение: "<<endl;

x.writeArray(cout);

ofstream out_file("result.txt");

if(out_file)

{

out_file<<"Решение интегрального уравнения“"<<endl;

out_file<<"Расширенная матрица системы:"<<endl;

a_b.writeArray(out_file);

out_file<<"Погрешность: "<<sigma<<endl;

out_file<<"Максимальное количество шагов: "<<max_step<<endl;

out_file<<"Решение: "<<endl;

x.writeArray(out_file);

}

else

{

cout<<"не могу создать файл\"result.txt\""<<endl;

}

getch();

break;

}

case 27:

{

return;

}

}}

}