Министерство образования Республики Беларусь
Учреждение образования
БЕЛОРУССКИЙ ГОСУДАРСТВЕННЫЙ УНИВЕРСИТЕТ
ИНФОРМАТИКИ И РАДИОЭЛЕКТРОНИКИ
Кафедра вычислительной математики и программирования
Пояснительная записка
к курсовому проекту
на тему:
«
Решение дифференциальных уравнений по методу Эйлера
»
Минск 2004
СОДЕРЖАНИЕ:
1. Введение
2. Математическое объяснение метода
2.1 Метод Эйлера
2.2 Исправленный метод Эйлера
2.3 Модифицированный метод Эйлера
3. Блок-схема алгоритма программы
4. Описание программы
Список использованной литературы
Приложение 1 (Текст программы)
Приложение 2 (Результаты работы программы)
1. Введение
Уравнение, содержащие производную функции одной переменной, возникают во многих областях прикладной математики. Вообще говоря, любая физическая ситуация, где рассматривается степень изменения одной переменной по отношению к другой переменной, описывается дифференциальным уравнением, а такие ситуации встречаются довольно часто.
Решение обыкновенных дифференциальных уравнений (нелинейных) первого порядка с начальными данными (задача Коши) – классическая область применения численных методов. Имеется много разностных методов, часть из которых возникла в домашинную эпоху и оказалось пригодным для современных ЭВМ.
В этой программе использовался метод Эйлера, один из самых старых и широко известных методов численного интегрирования дифференциальных уравнений. Этот метод имеет довольно большую ошибку; кроме того, он очень часто оказывается неустойчивым – малая начальная ошибка быстро увеличивается с ростом Х. Поэтому чаще используют более точные методы, такие как: исправленный метод Эйлера и модифицированный метод Эйлера. Нужно, однако, заметить, что метод Эйлера является методом Рунге – Кутта первого порядка.
2. Математическое объяснение метода
2.1 Метод Эйлера
Решить дифференциальное уравнение у/
=f(x,y) численным методом - это значит для заданной последовательности аргументов х0
, х1
…, хn
и числа у0
, не определяя функцию у=F(x), найти такие значения у1
, у2
,…, уn
, что уi
=F(xi
)(i=1,2,…, n) и F(x0
)=y0
.
Таким образом, численные методы позволяют вместо нахождения функции
У=F(x) получить таблицу значений этой функции для заданной последовательности аргументов. Величина h=xk
-xk-1
называется шагом интегрирования.
Метод Эйлера относиться к численным методам, дающим решение в виде таблицы приближенных значений искомой функции у(х). Он является сравнительно грубым и применяется в основном для ориентировочных расчетов. Однако идеи, положенные в основу метода Эйлера, являются исходными для ряда других методов.
Рассмотрим дифференциальное уравнение первого порядка
y/
=f(x,y) (1)
с начальным условием
x=x0
, y(x0
)=y0
(2)
Требуется найти решение уравнения (1) на отрезке [а,b].
Разобьем отрезок [a, b] на n равных частей и получим последовательность х0
, х1
, х2
,…, хn
, где xi
=x0
+ih (i=0,1,…, n), а h=(b-a)/n-шаг интегрирования.
В методе Эйлера приближенные значения у(хi
)»yi
вычисляются последовательно по формулам уi
+hf(xi
, yi
) (i=0,1,2…).
При этом искомая интегральная кривая у=у(х), проходящая через точку М0
(х0
, у0
), заменяется ломаной М0
М1
М2
… с вершинами Мi
(xi
, yi
) (i=0,1,2,…); каждое звено Мi
Mi+1
этой ломаной, называемой ломаной Эйлера,
имеет направление, совпадающее с направлением той интегральной кривой уравнения (1), которая проходит через точку Мi
, смотри рисунок 1.
Рисунок 1.
Если правая часть уравнения (1) в некотором прямоугольнике R{|x-x0
|£a, |y-y0
|£b}удовлетворяет условиям:
|f(x, y1
)- f(x, y2
)| £ N|y1
-y2
| (N=const),
|df/dx|=|df/dx+f(df/dy)| £ M (M=const),
то имеет место следующая оценка погрешности:
|y(xn
)-yn
| £ hM/2N[(1+hN)n
-1], (3)
где у(хn
)-значение точного решения уравнения(1) при х=хn
, а уn
- приближенное значение, полученное на n-ом шаге.
Формула (3) имеет в основном теоретическое применение. На практике иногда оказывается более удобным двойной просчет
: сначала расчет ведется с шагом h, затем шаг дробят и повторный расчет ведется с шагомh/2. Погрешность более точного значения уn
*
оценивается формулой
|yn
-y(xn
)|»|yn
*
-yn
|. (4)
Метод Эйлера легко распространяется на системы дифференциальных уравнений и на дифференциальные уравнения высших порядков. Последние должны быть предварительно приведены к системе дифференциальных уравнений первого порядка.
2.2 Исправленный метод Эйлера
В исправленном методе Эйлера мы находим средний тангенс наклона касательнй для двух точек: xm
, ym
и xm
+h, ym
+hy’m
. Последняя точка есть та самая, которая в простом методе обозначалась xm+1
, ym+1
. Геометрический процесс нахождения точки xm+1
, ym+1
можно проследить по рисунку 2. С помощью метода Эйлера находится точка xm
+h, ym
+hy’m,
лежащая на прямой L1.
В этой точке снова вычисляется тангенс угла наклона касательной, на рисунке этому значению соответствует прямая L2
. Усреднение двух тангенсов дает прямую L’3
. Наконец, через точку xm
, ym
мы проводим прямую L3
параллельную L’3
. Точка, в которой прямая L3
пересечется с ординатой, восстановленной из x= xm+1
=xm
+h, и будет искомой точкой y= ym+1
= ym
+hy’m
. Тангенс угла наклона L3
равен:
F(xm
, ym
)=1/2[f(xm
, ym
)+f(xm
+h, ym
+hy’m
)], (5)
где ym
= f(xm
, ym
) (6)
Уравнение линии L3
при этом записывается в виде:
y = ym
+ (x - xm
)*F(xm
) (7)
так что:
ym+1
= ym
+ h*F(xm
) (8)
Соотношения 5, 6, 7 и 8 описывают исправленный метод Эйлера. (рис. 2)
Рисунок 2.
2.3 Модифицированный метод Эйлера
Этот метод более точен. Рассмотрим дифференциальное уравнение (1) с начальным условием y(x0
)=y0
. Разобьем наш участок интегрирования на n равных частей. На малом участке [x0
,x0
+h] интегральную кривую заменимпрямой линией. Получаем точкуМк
(хк
,ук
). (рис. 3)
Y
Nk
/
y=y(x)
Мк
Мк
/
Yk+1
Yk
хк
хк1/2
xk+h
=xk1
X
Рисунок 3.
Через Мк
проводим касательную: y=yк
=f(xk
,yk
)(x-xk
). Делим отрезок (xк
,xк1
) пополам:
xh+k
/
=xk
+h/2=xk+1/2
(9)
yh+k
/
=yk
+f(xk
,yk
)h/2=yk
+yk+1/2
(10)
Получаем точку Nk
/
. В этой точке строим следующую касательную:
y(xk+1/2
)=f(xk+1/2
, yk+1/2
)=αk
(11)
Из точки Мк
проводим прямую с угловым коэффициентом αк
и определяем точку пересечения этой прямой с прямой xк1
. Получаем точку Мк
/
. В качестве ук+1
принимаем ординату точки Мк
/
. Тогда:
yк+1
=yк
+αк
h
xk+1
=xk
+h
αk
=f(xk+h/2
, xk
+f(xk
,yk
)h/2) (12)
yk
=yk-1
+f(xk-1
,yk-1
)h
Эти формулы называются рекуррентными формулами метода Эйлера.
Сначала вычисляют вспомогательные значения искомой функции yк+1/2
в точках xк+1/2
, затем находят значение правой части уравнения (1) в средней точке
y/
k+1/2
=f(xk+1/2
, yk+1/2
) и определяют yк+1
.
3. Блок-схема алгоритма
где A — начальное значение x, B — конечное значение x, F(x) — значение функции в точке xn
, N — количество промежутков, st – выбор операции, C1,C2,C3 – константы для формул, nom - сохраняет номер используемой функции.
На рисунке представлена блок-схема процесса решения дифференциального уравнения методом Эйлера
Подсчитывая каждый раз новое значение уравнения F(x), получаем последовательность значений xn
yn
, n=1,2,…
По этим значениям строим график.
4. Описание программы
Программа весьма проста. В ней много предусмотрено моментов неправильного ввода данных, о которых программа предупреждает пользователя и сразу же просит повторно ввести данные.
С самого начала программа предоставляет пользователю меню выполняемых функций, которые выделяются при помощи стрелок ↑ и ↓ выбор клавишей Enter:
Formyla -> Enter
-> Open in fails
Reshenie
Graphic
Exit
После запуска программы нужно выбрать Formyla -> Enter, эта опция позволит из предложенного списка формул выбрать одну, по которой компьютер будет производить расчет и строить график. Все предложенные формулы имеют номерацию; чтобы выбрать интересующий вас пример нажмите на цифру равную номеру примера, и сразу же появится новое окно, в котором сверху будет записан ваш пример. Также в окне будет этот же пример но с нулями на месте констант. Под примером будет высвечена большая буква С, это используется для ввода констант. Для этого вам нужно нажать номер константы, он появится, и после знака равно запишите чему она равна (вводятся целые и вещественные значения). По окончании набора нажать Enter. Операцию повторять пока не будут введены все числа. По окончании нажать Esc. После появится строчка «уточните границы изменения Х, от A= до B= » здесь нужно занести данные на каком промежутке абсциссы будет рассматриваться функция. Следующая строчка попросит ввести начальные данные y(A)=. Следующей строчкой будет вопрос: «сохранить данные в файле? Да/Нет» ответить на этот вопрос с помощью клавиш Д и Н (рус), после чего программа вернется в первоначальное меню. Если данные были сохранены (в папке с программой появляется файл form.txt), то в следующий раз чтобы не набирать снова выберите в меню опцию Formyla -> Open in fails и на экране появятся введенные данные с пометкой снизу, сообщая что данные были прочитаны из файла.
Следующая опция Reshenie. После нажатия в окне просят ввести N(целое число) – число промежутков, на которые разделится рассматриваемый участок (ось ОХ). После появится таблица рассчитанных данных (номер точки, значение абсциссы, значение ординаты). При нажатии любой клавиши произойдет переход в меню.
Graphic эта опция позволяет визуально видеть решение, а так же на этом графике прописываются все данные: начальная формула, шаг и промежуток построения графика, масштаб, данные об его изменении(клавишами +(увеличить) и -(уменьшить), а также возможность определить точное значение функции в любой точке.
Опция Exit применяется для выхода из программы.
Заключение
Результатом выполнения курсового проекта является готовый программный продукт, позволяющий решать дифференциальные уравнения по методу Эйлера, демонстрирующий возможности численного решения поставленной задачи с заданной степенью точности.
Данная программа решает заданную пользователем дифференциальное уравнение за минимальный промежуток времени. При этом пользователю предоставляется возможность визуально оценить решение, рассматривая график полученного решения.
К достоинствам программы можно отнести также удобный пользовательский интерфейс, возможность ввода пользовательских дифференциальных уравнений, а также давольно высокая стабильность работы. Однако имеются и некоторые недостатки. К недостаткам программы можно отнести: критичность к вводимым пользователем урававней, отсутствие обработки исключительных событий. Это, естественно, ограничивает возможности программы.
Список использованной литературы
1. Д. Мак-Кракен, У. Дорн. Численные методы и программирование на фортране. –М.: Мир,1977.-389,396-408 с.
2. А.А. Самарский. Введение в численные методы. – М.:Наука,1987.-176 с.
3. Алгоритмы вычислительной математики: Лабораторный практикум по курсу «Программирование» для студентов 1 - 2-го курсов всех специальностей БГУИР/А.К. Синицын, А.А. Навроцкий.- Мн.: БГУИР, 2002.- 65-69 с.
4. ГОСТ 2.105-95. Общие требования к текстовым документам.
5. ГОСТ 7.32-91. Система стандартов по информации, библиотечному и издательскому делу. Отчет о НИР. Структура и правила оформления.
Приложение 1. Текст программы.
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <conio.h>
#include <dos.h>
#include <graphics.h>
#include <math.h>
#include <bios.h>
//---------------------------------------------------------------------------
void formyl(int p)
{
if(p==1) printf("n 1. C1*y' = C2*y + C3*x + C4*x*y");
else if(p==2) printf("n 2. y'/(C1-100) = C2*y + C3*x + (C4+x)*y");
else if(p==3) printf("n 3. pow(e,C1)*y' = C2*y + C3*cos(x) + (C4+x+y)");
else if(p==4) printf("n 4. C1*sin(x)*y' = e*C2*y + C3*arcsin(x) + C4*y/x");
else if(p==5) printf("n 5. C1*y' = sin(C2)*y + tg(C3*x) + C4*ln(x)*y");
else if(p==6) printf("n 6. C1*y' = y*C2 + C3*sin(x) + C4*cos(x)*y");
else if(p==7) printf("n 7. (C1+C2+C3+C4)*y' = C2*y + (C3-x) + lg(C4*x)*y");
else if(p==8) printf("n 8. y'/C1 = y/C2 + C3*sin(x) + C4*x*y");
else if(p==9) printf("n 9. sin(C1)*y' = C2*y + |C3|*x + x*y/C4");
}
//---------------------------------------------------------------------------
void formyl2(int p,double C1,double C2,double C3,double C4)
{
if(p==1) printf("%.2f*y'=%.2f*y+%.2f*x+%.2f*x*y",C1,C2,C3,C4);
else if(p==2) printf("y'/(%.2f-100)=%.2f*y+%.2f*x+(%.2f+x)*y",C1,C2,C3,C4);else if(p==3) printf("pow(e,%.2f)*y'=%.2f*y+%.2f*cos(x)+(%.2f+x+y)",C1,C2,C3,C4);
else if(p==4) printf("%.2f*sin(x)*y'=e*%.2f*y+%.2f*arcsin(x)+%.2f*y/x",C1,C2,C3,C4);
else if(p==5) printf("%.2f*y'=sin(%.2f)*y+tg(%.2f*x)+%.2f*ln(x)*y",C1,C2,C3,C4);
else if(p==6) printf("%.2f*y'=y*%.2f+%.2f*sin(x)+%.2f*cos(x)*y",C1,C2,C3,C4);
else if(p==7) printf("(%.2f+%f+%.2f+%.2f)*y'=%.2f*y+(%.2f-x)+lg(%.2f*x)*y",C1,C2,C3,C4,C2,C3,C4);
else if(p==8) printf("y'/%.2f=y/%.2f+%.2f*sin(x)+%.2f*x*y",C1,C2,C3,C4);
else if(p==9) printf("sin(%.2f)*y'=%.2f*y+|%.2f|*x+x*y/%.2f",C1,C2,C3,C4);
}
//---------------------------------------------------------------------------
double formyl3(int p,double h,double x,double y,double C1,double C2,double C3,double C4)
{
double Y;
if(p==1) Y=h*(C2*y+C3*x+C4*x*y)/C1+y;
else if(p==2) Y=h*(C1-100)*(y*C2+C3*x+(C4+x)*y)+y;
else if(p==3) Y=h*(C2*y+C3*cos(x)+C4+x+y)/exp(C1)+y;
else if(p==4) Y=h*(exp(1)*C2*y+C3*asin(x)+C4*y/x)/(C1*sin(x))+y;
else if(p==5) Y=h*(sin(C2)*y+tan(C3*x)+C4*log10(x)*y)/C1+y;
else if(p==6) Y=h*(y*C2+C3*sin(x)+C4*cos(x)*y)/C1+y;
else if(p==7) Y=h*(C2*y+(C3-x)+log10(C4*x)*y)/(C1+C2+C3+C4)+y;
else if(p==8) Y=h*(y/C2+C3*sin(x)+C4*x*y)*C1+y;
else if(p==9) Y=h*(C2*y+abs(C3)*x+x*y/C4)/sin(C1)+y;
return(Y);
}
//---------------------------------------------------------------------------
void main(void)
{
int vv=0,vv1=0; // руководит операциями
int N=0,W; // кол промежутков
int i,j,k; // используются во всех "for"
int nom; // номер примера
int st=4,vst=0; // строчка в меню
double C1,C2,C3,C4; // константы
double M; // масштаб
double xtoch,ytoch; // считает y(x) по графику
double h,H; // шаг
double A=0,B=0,ii,jj,kk; // пределы интегрирования
double x[102],y[102]; // главные переменные x,y
//----Подключение графики и проверка-----------------------------------------
int g_driver=9,g_mode=2, g_error;
initgraph(&g_driver,&g_mode,"c:BORLANDbgi");
g_error=graphresult();
if(g_error!=grOk)
{
puts("error");
printf("n error=%d, reason=%sn", g_error, grapherrormsg(g_error));
getch();
exit(1);
}
closegraph();
//----Проверка или создание файла--------------------------------------------
FILE *fail;
if((fail = fopen("form.txt", "r")) == NULL)
if((fail = fopen("form.txt", "w")) == NULL)
{
clrscr();
printf("Ошибка при открытии файла");
getch();
exit;
}
fclose(fail);
//----Графическое меню-------------------------------------------------------
menu:
initgraph(&g_driver,&g_mode,"c:BORLANDbgi");
cleardevice();
setbkcolor(5);
setfillstyle(7,8); bar(0,0,getmaxx(),getmaxy());
setfillstyle(1,1); bar(200,45,465,145);
setfillstyle(1,15); bar(203,48,462,142);
setfillstyle(1,1); bar(206,51,459,139);
setcolor(11);
settextstyle(7,0,8);
outtextxy(220,40,"Menu");
setcolor(7);
settextstyle(1,0,4);
outtextxy(100,200,"1.Formyla");
outtextxy(100,240,"2.Reshenie");
outt
outtextxy(100,320," Exit");
settextstyle(2,0,5);
outtextxy(500,440,"beta ver. 101 c");
settextstyle(1,0,4);
izm: // Выборпунктаменю
if(st==1) goto d1;
else if(st==11) goto d11;
else if(st==12) goto d12;
else if(st==2) goto d2;
else if(st==3) goto d3;
else if(st==4) goto d4;
d1:
setcolor(2);
outtextxy(100,200,"1.Formyla");
setcolor(7);
outtextxy(100,240,"2.Reshenie");
outtextxy(100,320," Exit");
if(W==0)
{
W=1;
setfillstyle(7,1);
for(i=0;i<35;i++)
{
bar(280+i*7,190,287+i*7,260);
delay(5);
}
setfillstyle(7,8);
}
outtextxy(300,185,"Enter");
outtextxy(300,215,"Open in fails");
vst=getch();
if(vst==72)
{
st=4;
for(i=0;i<35;i++)
{
bar(530-i*7,190,523-i*7,260);
delay(5);
}
goto izm;
}
else if(vst==80)
{
st=2;
for(i=0;i<35;i++)
{
bar(530-i*7,190,523-i*7,260);
delay(5);
}
goto izm;
}
else if(vst==77) { st=11; goto izm; }
else goto d1;
d11:
setfillstyle(7,1);
bar(280,190,530,260);
setfillstyle(7,8);
setcolor(2);
outtextxy(300,185,"Enter");
setcolor(7);
outtextxy(100,200,"1.Formyla");
outtextxy(100,240,"2.Reshenie");
outtextxy(100,320," Exit");
outtextxy(300,215,"Open in fails");
vst=getch();
if(vst==72) { st=12; goto izm; }
else if(vst==80) { st=12; goto izm; }
else if(vst==75) { st=1; goto izm; }
else if(vst==13) { vv=1; vv1=1; closegraph(); goto resh; }
else goto d11;
d12:
setfillstyle(7,1);
bar(280,190,530,260);
setfillstyle(7,8);
setcolor(2);
outtextxy(300,215,"Open in fails");
setcolor(7);
outtextxy(100,200,"1.Formyla");
outtextxy(100,240,"2.Reshenie");
outtextxy(100,320," Exit");
outtextxy(300,185,"Enter");
vst=getch();
if(vst==72) { st=11; goto izm; }
else if(vst==80) { st=11; goto izm; }
else if(vst==75) { st=1; goto izm; }
else if(vst==13) { vv=1; vv1=2; closegraph(); goto resh; }
else goto d12;
d2:
setfillstyle(1,5);
setfillstyle(7,8);
bar(280,160,600,400);
setcolor(2);
outtextxy(100,240,"2.Reshenie");
setcolor(7);
outtextxy(100,200,"1.Formyla");
outtextxy(100,280,"3.Graphic");
vst=getch();
if(vst==13) { vv=2; closegraph(); goto resh; }
else if(vst==72) { st=1; W=0; goto izm; }
else if(vst==80) { st=3; goto izm; }
else goto d2;
d3:
setcolor(2);
outtextxy(100,280,"3.Graphic");
setcolor(7);
outtextxy(100,240,"2.Reshenie");
outtextxy(100,320," Exit");
vst=getch();
if(vst==13) { vv=3; closegraph(); goto resh; }
else if(vst==72) { st=2; goto izm; }
else if(vst==80) { st=4; goto izm; }
else goto d3;
d4:
setfillstyle(1,5);
setfillstyle(7,8);
bar(280,160,600,400);
setcolor(2);
outtextxy(100,320," Exit");
setcolor(7);
outtextxy(100,200,"1.Formyla");
outtextxy(100,240,"2.Reshenie");
outtextxy(100,280,"3.Graphic");
vst=getch();
if(vst==13) { vv=4; closegraph(); goto resh; }
else if(vst==72) { st=3; goto izm; }
else if(vst==80) { st=1; W=0; goto izm; }
else goto d4;
//----Вводданных------------------------------------------------------------
resh:
if(vv==1)
{
if(vv1==1)
{
C1=0;C2=0;C3=0;C4=0;
clrscr();
printf("nn Выбор формулы решенияn");
for(i=1;i<10;i++)
formyl(i);
printf("nn нажмите любой номер: ");
while(2)
{
nom=getch()-48;
if((nom>0)&(nom<10))
break;
}
while(3)
{
clrscr();
printf("n поокончаниивводаконстантнажать - Esc");
printf("nn выбронная формула имеет видn ");
formyl(nom);
printf("nnnn ввод констант в формулу дифференциального уравнения:nn ");
formyl2(nom,C1,C2,C3,C4);
printf("nnn C");
i=getch();
if(i==49)
{
printf("1 = ");
scanf("%lf",&C1);
}
else if(i==50)
{
printf("2 = ");
scanf("%lf",&C2);
}
else if(i==51)
{
printf("3 = ");
scanf("%lf",&C3);
}
else if(i==52)
{
printf("4 = ");
scanf("%lf",&C4);
}
else if(i==27)
{
break;
}
}
AB:
clrscr();
printf("nnn выброннаяформулаимеетвидn ");
formyl(nom);
printf("nnnn ввод констант в формулу дифференциального уравнения:nn ");
formyl2(nom,C1,C2,C3,C4);
printf("nnn уточните границы изменения Х, от A= ");
scanf("%lf",&A);
clrscr();
printf("nnn выбронная формула имеет видn ");
formyl(nom);
printf("nnnn ввод констант в формулу дифференциального уравнения:nn ");
formyl2(nom,C1,C2,C3,C4);
printf("nnn уточните границы изменения Х, от A= %.2f до B= ",A);
scanf("%lf",&B);
if(B<=A)
{
printf("nn неправильно заданны границы!!! ");
delay(1000);
goto AB;
}
clrscr();
printf("nnn выбронная формула имеет видn ");
formyl(nom);
printf("nnnn ввод констант в формулу дифференциального уравнения:nn ");
formyl2(nom,C1,C2,C3,C4);
printf("nnn уточните границы изменения Х, от A= %.2f до B= %.2f",A,B);
printf("n задайте начальное значение Y(%.2f) = ",A);
scanf("%lf",&y[0]);
printf("nn сохранить данные в файле (Да/Нет) ");
while(11)
{
i=getch();
if((i==76)||(i==108)||(i==132)||(i==164))
{
printf("Да");
delay(300);
fopen("form.txt", "w");
fwrite(&nom,1,2,fail);
fwrite(&C1,1,8,fail);
fwrite(&C2,1,8,fail);
fwrite(&C3,1,8,fail);
fwrite(&C4,1,8,fail);
fwrite(&A,1,8,fail);
fwrite(&B,1,8,fail);
fwrite(&y[0],1,8,fail);
fclose(fail);
break;
}
else if((i==89)||(i==121)||(i==141)||(i==173))
{
printf("Нет");
delay(300);
break;
}
}
goto menu;
}
if(vv1==2)
{
fopen("form.txt", "r");
fread(&nom,1,2,fail);
fread(&C1,1,8,fail);
fread(&C2,1,8,fail);
fread(&C3,1,8,fail);
fread(&C4,1,8,fail);
fread(&A,1,8,fail);
fread(&B,1,8,fail);
fread(&y[0],1,8,fail);
fclose(fail);
clrscr();
printf("nnn выбраннаяформулаимеетвидn");
formyl(nom);
printf("nnnn ввод констант в формулу дифференциального уравнения:nn ");
formyl2(nom,C1,C2,C3,C4);
printf("nnn уточните границы изменения Х, от A= %.2f до B= %.2f",A,B);
printf("n задайте начальное значение Y(%.2f) = %.2f",A,y[0]);
printf("nnnn данные были прочитаны из ранее сохраненного файла");
getch();
goto menu;
}
}
//----Решение примера--------------------------------------------------------
else if(vv==2)
{
if((A==0)&(B==0))
{
clrscr();
printf("nnnn Сперва необходимо записать пример!");
delay(2000);
goto menu;
}
clrscr();
printf("nnn для решения этой задачи нужно:nn");
printf(" задать колличество отрезков, на которое расделитсяn промежуток AB=[%.2f,%.2f] N= ",A,B);
scanf("%d",&N);
ii=B; jj=A;
h=(B-A)/N;
N++;
//-------------------------------------
for(i=0;i<N;i++)
{
x[i]=A+h*i;
y[i+1]=formyl3(nom,h,x[i],y[i],C1,C2,C3,C4);
}
clrscr();
printf("nnn Решениезадачи:nn");
printf("
╔════════╦═══════════════════╦════════════════════════╗n");
printf(" ║ i ║ x=h*i ║ y(i) ║n");
printf(" ╠════════╬═══════════════════╬════════════════════════╣n");
for(i=0;i<N;i++)
printf(" ║ %3d ║ %11lf ║ %14f ║ n",i+1,x[i],y[i]);
printf(" ╚════════╩═══════════════════╩════════════════════════╝n");
getch();
}
//----График-----------------------------------------------------------------
else if(vv==3)
{
if(N==0)
{
clrscr();
printf("nnnn Сперва необходимо решить пример!");
delay(2000);
goto menu;
}
initgraph(&g_driver,&g_mode,"c:BORLANDbgi");
cleardevice();
setbkcolor(0);
setcolor(11);
settextstyle(4,0,5);
outtextxy(100,280,"Open graphic...");
setfillstyle(1,7);
bar(218,238,422,252);
setfillstyle(1,5);
bar(220,240,420,250);
setcolor(9);
for(i=0;i<200;i=i+2)
{
line(221+i,241,221+i,249);
delay(10);
}
M=5;
while(31)
{
// Фон
clrscr();
cleardevice();
setbkcolor(0);
setcolor(10);
settextstyle(5,0,4);
// Выводданныхпографику
printf(" выход в меню - Esc");
printf("n график дифференциального уравнения - ");
formyl2(nom,C1,C2,C3,C4);
printf("n построенный с шагом h = %f",h);
printf(" на промежутке от %.2f до %.2f",A,B);
printf("n выполненный в масштабе %.4f: 1",M);
printf("nn найти решение в точке - Enter");
printf("n изменение масштаба увеличить - +");
printf("n уменьшить - -");
// обводка слов
setcolor(3);
line(455,0,630,0);
line(455,15,630,15);
line(455,0,455,15);
line(630,0,630,15);
setcolor(9);
// Оси
setcolor(10);
line(50,250,600,250);
line(580,245,600,250);
line(580,255,600,250);
line(325,50,325,440);
line(325,50,330,70);
line(325,50,320,70);
setcolor(15);
settextstyle(1,0,1);
outtextxy(600,260,"X");
outtextxy(345,50,"Y");
outtextxy(300,260,"O");
// единичныенасечки
if(M>4)
for(i=1;i<70;i++)
{
if(i*M<=255)
{
line(325+M*i,248,325+M*i,252);
line(325-M*i,248,325-M*i,252);
}
if(i*M<=170)
{
line(323,250-M*i,327,250-M*i);
line(323,250+M*i,327,250+M*i);
}
}
// Прорисовка графика
setcolor(9);
for(i=1;i<N;i++)
{
if((M*y[i]>-1000)&(M*y[i]<1000))
{
line(325+M*((i-1)*h+A),250-M*y[i-1],325+M*(i*h+A),250-M*y[i]);
delay(10);
}
}
k=getch();
if(k==43) // масштаб +
{
if (M>=1) M++;
else if((M<1) &&(M>=0.1)) M=M+0.1;
else if((M<0.1) &&(M>=0.01)) M=M+0.01;
else if((M<0.01) &&(M>=0.001)) M=M+0.001;
else if((M<0.001) &&(M>=0.0001)) M=M+0.0001;
else if((M<0.0001) &&(M>=0)) M=M+0.00001;
}
else if(k==45) // масштаб -
{
if (M>1.000001) M--;
else if((M<=1) &&(M>0.100001)) M=M-0.1;
else if((M<=0.100001)&&(M>0.010001)) M=M-0.01;
else if((M<=0.010001)&&(M>0.001001)) M=M-0.001;
else if((M<=0.001001)&&(M>0.000101)) M=M-0.0001;
else if((M<=0.000100)&&(M>0.000010)) M=M-0.00001;
}
else if(k==13) // Нахождение y по значению x
{
printf("nn введите значение абсциссы x= ");
scanf("%lf",&xtoch);
if((A<=xtoch)&&(xtoch<=B))
{
setcolor(11);
for(i=1;i<N;i++)
{
if((x[i-1]<=xtoch)&&(xtoch<=x[i]))
{
H=xtoch-x[i-1];
ytoch=formyl3(nom,H,x[i-1],y[i-1],C1,C2,C3,C4);
printf("ордината точки равна y= %lf",ytoch);
line(325+M*xtoch,250,325+M*xtoch,250-M*ytoch);
break;
}
}
setcolor(9);
}
else
printf("введенаабсцисса, непренадлежащаяпромежуткурешения");
getch();
}
else if(k==27)
break;
}
closegraph();
}
//----Выходизпрограммы-----------------------------------------------------
else if(vv==4)
{
exit(2);
}
//----Ввод программы в бесконечный цикл--------------------------------------
goto menu;
}
Приложение 2. Результаты работы программы
Приведем расчет дифференциального уравнения первого, второго и третьего порядка методом Эйлера
1. Пусть дано дифференциальное уравнение первого порядка: y/
=2x-y
Требуется найти решение на отрезке [0,1] c шагом h=(1-0)/5=0,2
Начальные условия: у0
=1;
Пользуясь рекурентными формулами (4), находим:
1). x1
=0,2; х1/2
=0,1; y(x1
)=y(x0
)+α0
h; y(x1/2
)=y(x0
)+f(x0
,y0
)h/2;
f(x0
,y0
)=2*0-1=-1
y(x1/2
)=1-1*0,1=0,9
α0
=2*0,1-0,9=-0,7
y1
=1-0,1*0,2=0,86
2). y(x2
)=y(x1
)+α1
h; x2
=0,2+0,2=0,4; x1+1/2
=x1
+h/2=0,2+0,1=0,3
y(x1+1/2
)=y(x1
)+f(x1
,y(x1
))h/2
f(x1
,y1
)=2*0,2-0,86=-0,46
y(x1+1/2
)=0,86-0,46*0,1=0,814
α1
=2*0,3-0,814=-0,214
y2
=0,86-0,214*0,2=0,8172
3). x3
=0,4+0,2=0,6; x2+1/2
=x2
+h/2=0,4+0,1=0,5
f(x2
,y2
)=2*0,4-0,8172=-0,0172
y2+1/2
=0,8172-0,0172*0,1=0,81548
α2
=2*0,5-0,81548=0,18452
y3
=0,8172+0,18452*0,2=0,854104
4).x4
=0,8; x3+1/2
=x3
+h/2=0,6+0,1=0,7
f(x3
,y3
)=2*0,6-0,854104=0,345896
y3+1/2
=0,854104+0,345896*0,1=0,8886936
α3
=2*0,7-0,89=0,5113064
y4
=0,854104+0,5113064*0,2=0,95636528
5).x5
=1; x4+1/2
=0,8+0,1=0,9
f(x4
,y4
)=2*0,8-0,956=0,64363472
y4+1/2
=0,956+0,643*0,1=1,020728752;
α4
=2*0,9-1,02=0,779271248
y5
=0,956+0,7792*0,2=1,11221953