Главная Коллекция "Otherreferats" Программирование, компьютеры и кибернетика Одношаговые методы решения обыкновенного дифференциального уравнения

Одношаговые методы решения обыкновенного дифференциального уравнения

Характеристика практического применения одношаговых методов Эйлера, Рунге-Кутта 2-го и 4-го порядков. Результаты расчетов и анализ их применения. Специфика кода программы решения перечисленных методов на языке программирования Borland C++ Builder 6.

Рубрика Программирование, компьютеры и кибернетика
Вид курсовая работа
Язык русский
Дата добавления 01.12.2009
Размер файла 40,7 K

Отправить свою хорошую работу в базу знаний просто. Используйте форму, расположенную ниже

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

13

Введение

Данная курсовая работа включает в себя три одношаговых метода решения обыкновенного дифференциального уравнения.

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

1. Метод Эйлера.

2. Метод Рунге-Кутта 2-го порядка (метод Гюна).

3. Метод Рунге-Кутта 4-го порядка.

Также данная курсовая работа включает в себя: описание метода, применение метода к конкретной задаче (анализ), код программы решения вышеперечисленных методов на языке программирования Borland C++ Builder 6.

Описание метода

Рассмотрим некоторые численные методы решения задачи Коши (начальной задачи) обыкновенных дифференциальных уравнений первого порядка. Запишем данное уравнение в общем виде, разрешенном относительно производной (правая часть уравнения не зависит от первой производной):

(1)

Необходимо найти значения функции y в заданных точках сетки x1,x2,…,xn, если известны начальные значения (x0,y0), где y0=y(x0) есть значение функции y(x) в начальной точке x0.

Преобразуем уравнение (1) умножением на : . И проинтегрируем левую и правую части между k-ым и k+1-ым узлами сетки:

(2)

Мы получили выражение для построения решения в k+1 узле интегрирования через значения x и y в k-ом узле сетки. Сложность заключается в том, что интеграл в правой части есть интеграл от неявно заданной функции, нахождение которого в аналитическом виде в общем случае невозможно. Численные методы решения ОДУ различным способом аппроксимируют (приближают) значение этого интеграла для построения формул численного интегрирования ОДУ.

Наша задача состоит в том, чтобы рассмотреть и реализовать в программном виде методы Эйлера, Рунге-Кутта 2-го порядка (метод Гюна) и Рунге-Кутта 4-го порядка. Они дают начальное представление о подходах к решению данной задачи в рамках численного решения задачи Коши.

Рассмотрим и опишем отдельно каждый метод:

1. Метод Эйлера.

Первым и наиболее простым способом численного решения задачи Коши для ОДУ первого порядка является метод Эйлера. В его основе лежит аппроксимация производной отношением конечных приращений зависимой y и независимой x переменных между узлами равномерной сетки:

,

где yk+1 это искомое значение функции в точке xk+1.

Если теперь преобразовать это уравнение, и учесть равномерность сетки интегрирования, то получится итерационная формула, по которой можно вычислить yk+1 , если известно yk в точке хк:

(3)

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

Также метод Эйлера можно представить в графической интерпретации. Действительно, исходя из вида решаемого уравнения (1) следует, что значение F(xk,yk) есть значение производной функции y(x) в точке x=xk - , и, таким образом, равно тангенсу угла наклона касательной, проведенной к графику функции y(x) в точке x=xk.

Из прямоугольного треугольника изображенного на рисунке можно найти => - формула Эйлера. Таким образом, суть метода Эйлера заключается в замене функции y(x) на отрезке интегрирования прямой линией, касательной к графику в точке x=xk. Если искомая функция сильно отличается от линейной на отрезке интегрирования, то погрешность вычисления будет значительной.

Ошибка метода Эйлера прямо пропорциональна шагу интегрирования (ошибка ~ h).

Процесс вычислений строится следующим образом. При заданных начальных условиях x0 и y0 можно вычислить:

y1 = y0 + F(x0,y0)hx1 = x0 + h

y2 = y1 + F(x1,y1)hx2 = x1 + h

y3 = y2 + F(x2,y2)hx3 = x2 + h

…………………………………

yn = yn-1 + F(xn-1,yn-1)hxn = xn-1 + h

Таким образом, строится таблица значений функции y(x) с определенным шагом h по x на отрезке [x0, xn]. Ошибка в определении значения y(xk) при этом будет тем меньше, чем меньше выбрана длина шага h (что определяется точностью формулы интегрирования).

При больших h метод Эйлера весьма неточен. Он дает все более точное приближение при уменьшении шага интегрирования. Если отрезок [xk, xk+1] слишком велик, то каждый участок [xk, xk+1] разбивается на n отрезков интегрирования и к каждому из них применяется формула Эйлера с шагом , то есть шаг интегрирования h берется меньше шага сетки, на которой определяется решение.

2. Метод Рунге-Кутта 2-го порядка (метод Гюна).

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

(4)

Данная формула оказывается неявной относительно yk+1 (это значение есть и в левой и в правой части выражения), то есть является уравнением относительно yk+1. Однако, можно поступить иначе и приблизительно вычислить значение функции в узле k+1 с помощью обычной формулы Эйлера: , которое можно использовать при вычислении (4).

Таким образом, метод Рунге-Кутта 2-го порядка (метод Гюна) получается с пересчетом, т.е. как бы усовершенствованием метода Эйлера. Для каждого узла интегрирования производится следующая цепочка вычислений:

(5)

Благодаря более точной формуле интегрирования, погрешность метода Рунге-Кутта 2-го порядка (метод Гюна) пропорциональна уже квадрату шага интегрирования.

3. Метод Рунге-Кутта 4-го порядка.

Повышение точности модифицированных методов Эйлера было достигнуто за счёт дополнительных по сравнению с обычным методом Эйлера вычислений функции f(x,y) из правой части дифференциального уравнения. При этом вычислять частные производные от функции f не требовалось. На этой идее дополнительных вычислений правой части основаны методы Рунге-Кутты высокой точности. В этих методах правая часть дифференциального уравнения вычисляется в нескольких точках, составляется линейная комбинация вычисленных значений, которая используется при определении значения yk+1.

Дальнейшее улучшение точности решения ОДУ первого порядка возможно за счет увеличения точности приближенного вычисления интеграла в выражении.

Выше было описано какое преимущество дает переход от метода Эйлера к методу Гюна при аппроксимации этого интеграла.

При использовании различных методов приближенного вычисления, получаются выражения для методов Рунге-Кутта различного порядка точности.

Алгоритм Рунге-Кутта четвертого порядка (погрешность порядка h4):

(6), где k0 =hF(xk, yk)

k 1=hF(xk+h/2, yk+k0/2)

k 2=hF(xk+h/2, yk+k1/2)

k 3=hF(xk+h, yk+k2)

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

Применение метода к конкретной задаче (анализ).

Дано дифференциальное уравнение вида 4xy'' + 2y' + y = 0, x[1,2], где y(1)=1,38177, y'(1)=-0,15058.

Точное решение этого уравнения имеет вид .

Задача состоит в том, чтобы решить это дифференциальное уравнение, применив одношаговые методы решения ОДУ: метод Эйлера, метод Гюна, метод Рунге-Кутта 4-го порядка.

Изучив методы и применив их к данному дифференциальному уравнению, приходим к такому выводу: при решении данного уравнения тремя известными способами при шаге h=0,1 все три метода дают приблизительно одинаковый результат, отличается лишь погрешность (но отличается значительно!).

Сравнительная Таблица 1

x

Метод Эйлера

Метод Гюна

Метод Рунге-Кутта 4-го порядка

y

отн. погр.

y

отн. погр.

y

отн. погр.

1

1,38177

3,29067604005725E-6

1,38177

3,29067604005725E-6

1,38177

3,29067604005725E-6

1,1

1,366712

0,00127814487613993

1,365453334

1,94798193216625E-5

1,365430962

2,89225942433997E-6

1,2

1,348952475

0,00216968182269993

1,346797805

1,50127182925036E-5

1,346780296

2,49669502851158E-6

1,3

1,328894037

0,00275219839584001

1,326127541

1,42980303934863E-5

1,326139732

2,10696659125027E-6

1,4

1,306861050

0,00308541704017657

1,303709146

6,64873545153573E-5

1,303773908

1,72461604735239E-6

1,5

1,283119920

0,00321459803483112

1,279765447

0,000139875317013123

1,279903972

1,35042681162576E-6

1,6

1,257893007

0,00317480492688739

1,254485219

0,00023298398264003

1,254717218

9,84768342856233E-7

1,7

1,231368459

0,00299371987609308

1,228030252

0,000344486721386327

1,228374111

6,27778015693109E-7

1,8

1,203707344

0,00269355862283386

1,200540610

0,000473175694179785

1,201013506

2,79460818010651E-7

1,9

1,175048956

0,0022924073958248

1,172138608

0,000617939973481596

1,172756609

6,0254226042205E-8

2

1,145514818

0,0018051777630581

1,142931890

0,000777750138060593

1,143710032

3,91481384411801E-7

Причем, если в случае применения метода Эйлера при h=0,1 погрешность весьма велика (результат можно посмотреть в Таблице 1), то в случае применения метода Рунге-Кутта 2-го порядка погрешность значительно уменьшается. В случае применения метода Гюна результаты намного точнее, чем результаты, полученные методом Эйлера, причем при том же шаге интегрирования h=0,1(сравнить результаты можно в Таблице 1).

Отметим существенное увеличение точности вычислений метода Гюна по сравнению с методом Эйлера. Так, для узла x=0,1 относительное отклонение значения функции, определенного методом Гюна, оказывается во много раз меньше. Следовательно, при использовании метода Гюна при одинаковой точности вычислений понадобится на много меньше времени ЭВМ, чем при использовании метода Эйлера.

Теперь рассмотрим и сравним метод Рунге-Кутта 4-го порядка с двумя другими методами на интервале x[1,2] с шагом h=0,1 и его относительное отклонение от точного решения. Как видно из Таблицы 1, точность решения, полученного методом Рунге-Кутта 4-го порядка, намного превышает точность решения, полученного методом Гюна, не говоря уже о методе Эйлера. При шаге h = 0,1 он позволил точно определить четыре значащие цифры решения, тогда как для достижения такой точности с помощью метода Эйлера необходимо взять h = 0,0001, что требует неимоверное количество вычислений функции F(x,y).

Высокая точность, вместе с достаточной простотой реализации делает метод Рунге-Кутта 4-го порядка одним из весьма распространенных численных методов решения задачи Коши ОДУ.

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

//---------------------------------------------------------------------------

#include <vcl.h>

#pragma hdrstop

#include "Unit1.h"

#include <math.h>

//---------------------------------------------------------------------------

#pragma package(smart_init)

#pragma resource "*.dfm"

TForm1 *Form1;

int n=10;

double x1=1,x2=2;

double x[100],y[100],p_y[100];

double h=(x2-x1)/n;

//---------------------------------------------------------------------------

__fastcall TForm1::TForm1(TComponent* Owner)

: TForm(Owner)

{

}

//---------------------------------------------------------------------------

double dif_ur(double x,double y,double p_y)

{

return ((-2*p_y-y)/(4*x));

}

void __fastcall TForm1::Button_EylerClick(TObject *Sender)

{

Label_TypeMethod->Caption="Yeea?";

y[0] = 1.38177;

x[0] = x1;

p_y[0] = -0.15058;

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

{

x[i] = x[i - 1] + h;

p_y[i] = p_y[i - 1] + h * (-2*p_y[i-1]-y[i-1])/(4*x[i-1]);

y[i]=y[i-1]+h*p_y[i-1];

}

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

{

StringGrid1->Cells[0][j]=FloatToStr(x[j]);

StringGrid1->Cells[1][j]=FloatToStr(y[j]);

StringGrid1->Cells[2][j]=FloatToStr(sin(sqrt(x[j]))+cos(sqrt(x[j])));

StringGrid1->Cells[3][j]=FloatToStr(fabs(StrToFloat(StringGrid1->Cells[2][j])-y[j]));

}

}

//---------------------------------------------------------------------------

void __fastcall TForm1::Button_ExitClick(TObject *Sender)

{

Close();

}

//---------------------------------------------------------------------------

void __fastcall TForm1::Button_GunClick(TObject *Sender)

{

Label_TypeMethod->Caption="A?i";

y[0] = 1.38177;

x[0] = x1;

p_y[0] = -0.15058;

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

{

x[i]=x[i-1]+h;

p_y[i]=p_y[i-1]+(h/2.0)*(dif_ur(x[i-1],y[i-1],p_y[i-1])+dif_ur(x[i-1]+h,y[i-1],p_y[i-1]+h*dif_ur(x[i-1],y[i-1],p_y[i-1])));

y[i]=y[i-1]+(h/2.0)*(p_y[i-1]+p_y[i]);

}

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

{

StringGrid1->Cells[0][j]=FloatToStr(x[j]);

StringGrid1->Cells[1][j]=FloatToStr(y[j]);

StringGrid1->Cells[2][j]=FloatToStr(sin(sqrt(x[j]))+cos(sqrt(x[j])));

StringGrid1->Cells[3][j]=FloatToStr(fabs(StrToFloat(StringGrid1->Cells[2][j])-y[j]));

}

}

//---------------------------------------------------------------------------

void __fastcall TForm1::RadioButton1Click(TObject *Sender)

{

Button_Eyler->Visible=true;

Button_Gun->Visible=false;

Button_Runge->Visible=false;

}

//---------------------------------------------------------------------------

void __fastcall TForm1::RadioButton2Click(TObject *Sender)

{

Button_Eyler->Visible=false;

Button_Gun->Visible=true;

Button_Runge->Visible=false;

}

//---------------------------------------------------------------------------

void __fastcall TForm1::RadioButton3Click(TObject *Sender)

{

Button_Eyler->Visible=false;

Button_Gun->Visible=false;

Button_Runge->Visible=true;

}

//---------------------------------------------------------------------------

void __fastcall TForm1::Button_RungeClick(TObject *Sender)

{

Label_TypeMethod->Caption="?oiaa-Eoooa";

y[0] = 1.38177;

x[0] = x1;

p_y[0] = -0.15058;

double F1,F2,F3,F4;

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

{

x[i]=x[i-1]+h;

F1=h*dif_ur(x[i-1],y[i-1],p_y[i-1]);

F2=h*dif_ur(x[i-1]+h/2,y[i-1]+p_y[i-1]*h/2+F1*h/8,p_y[i-1]+F1/2);

F3=h*dif_ur(x[i-1]+h/2,y[i-1]+p_y[i-1]*h/2+F1*h/8,p_y[i-1]+F2/2);

F4=h*dif_ur(x[i-1]+h,y[i-1]+p_y[i-1]*h+h*F3/2,p_y[i-1]+F3);

p_y[i]=p_y[i-1]+(F1+2*F2+2*F3+F4)/6;

y[i]=y[i-1]+h*(p_y[i-1]+(F1+F2+F3)/6);

y[i]=y[i-1]+h*p_y[i];

}

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

{

StringGrid1->Cells[0][j]=FloatToStr(x[j]);

StringGrid1->Cells[1][j]=FloatToStr(y[j]);

StringGrid1->Cells[2][j]=FloatToStr(sin(sqrt(x[j]))+cos(sqrt(x[j])));

StringGrid1->Cells[3][j]=FloatToStr(fabs(StrToFloat(StringGrid1->Cells[2][j])-y[j]));

}

}

//---------------------------------------------------------------------------

Результаты расчета

Метод Эйлера

x

exact solution

y

error

1

1,38177329067604

1,38177

3,29067604005725E-6

1,1

1,36543385512386

1,366712

0,00127814487613993

1,2

1,3467827931773

1,348952475

0,00216968182269993

1,3

1,32614183910416

1,3288940375

0,00275219839584001

1,4

1,30377563353274

1,30686105057292

0,00308541704017657

1,5

1,27990532272839

1,28311992076322

0,00321459803483112

1,6

1,25471820307239

1,25789300799928

0,00317480492688739

1,7

1,2283747392501

1,23136845912619

0,00299371987609308

1,8

1,20101378595756

1,20370734458039

0,00269355862283386

1,9

1,17275654886199

1,17504895625781

0,0022924073958248

2

1,14370964075811

1,14551481852117

0,0018051777630581

Метод Рунге-Кутта 2-го порядка (метод Гюна)

x

exact solution

y

error

1

1,38177329067604

1,38177

3,29067604005725E-6

1,1

1,36543385512386

1,36545333494318

1,94798193216625E-5

1,2

1,3467827931773

1,34679780589559

1,50127182925036E-5

1,3

1,32614183910416

1,32612754107377

1,42980303934863E-5

1,4

1,30377563353274

1,30370914617822

6,64873545153573E-5

1,5

1,27990532272839

1,27976544741138

0,000139875317013123

1,6

1,25471820307239

1,25448521908975

0,00023298398264003

1,7

1,2283747392501

1,22803025252871

0,000344486721386327

1,8

1,20101378595756

1,20054061026338

0,000473175694179785

1,9

1,17275654886199

1,17213860888851

0,000617939973481596

2

1,14370964075811

1,14293189062005

0,000777750138060593

Метод Рунге-Кутта 4-го порядка.

x

exact solution

y

error

1

1,38177329067604

1,38177

3,29067604005725E-6

1,1

1,36543385512386

1,36543096286444

2,89225942433997E-6

1,2

1,3467827931773

1,34678029648227

2,49669502851158E-6

1,3

1,32614183910416

1,32613973213757

2,10696659125027E-6

1,4

1,30377563353274

1,30377390891669

1,72461604735239E-6

1,5

1,27990532272839

1,27990397230158

1,35042681162576E-6

1,6

1,25471820307239

1,25471721830405

9,84768342856233E-7

1,7

1,2283747392501

1,22837411147208

6,27778015693109E-7

1,8

1,20101378595756

1,20101350649674

2,79460818010651E-7

1,9

1,17275654886199

1,17275660911622

6,0254226042205E-8

2

1,14370964075811

1,14371003223949

3,91481384411801E-7

Подобные документы

Работы в архивах красиво оформлены согласно требованиям ВУЗов и содержат рисунки, диаграммы, формулы и т.д.
PPT, PPTX и PDF-файлы представлены только в архивах.
Рекомендуем скачать работу.

© 2000 — 2023, ООО «Олбест» Все права защищены