Министерство образования и науки РФ
Глазовский инженерно-экономический институт(филиал) федерального государственного учреждения высшего образования
“Ижевский государственный технический университет имени М.Т. Калашникова”
(ГИЭИ (филиал) ФГБОУ ВО “ ИжГТУ имен М.Т. Калашникова)
Кафедра основы архитектуры, устройство и функционирование ВС
КУРСОВАЯ РАБОТА
Численное интегрирование дифференциальных уравнений методом Эйлера
Выполнила
студентка группы Д19-790-31, обучающаяся по специальности «Информационные системы» Ельцова С.Г.
Научный руководитель: Савельева Т.А.
г. Глазов
2021
ВВЕДЕНИЕ
Уравнение, содержащие производную функции одной переменной, возникают во многих областях прикладной математики. Вообще говоря, любая физическая ситуация, где рассматривается степень изменения одной переменной по отношению к другой переменной, описывается дифференциальным уравнением, а такие ситуации встречаются довольно часто. численное дифференцирование метод эйлера
Решение обыкновенных дифференциальных уравнений (нелинейных) первого порядка с начальными данными (задача Коши) - классическая область применения численных методов. Имеется много разностных методов, часть из которых возникла в домашинную эпоху и оказалось пригодным для современных ЭВМ.
В этой программе использовался метод Эйлера, один из самых старых и широко известных методов численного интегрирования дифференциальных уравнений. Этот метод имеет довольно большую ошибку; кроме того, он очень часто оказывается неустойчивым - малая начальная ошибка быстро увеличивается с ростом Х. Поэтому чаще используют более точные методы, такие как: исправленный метод Эйлера и модифицированный метод Эйлера. Нужно, однако, заметить, что метод Эйлера является методом Рунге - Кутта первого порядка.
1. Численное дифференцирование
Численное дифференцирование применяется в тех случаях, кода функция f(x)задана таблично и методы дифференциального исчисления неприменимы.
В основе численного дифференцирования лежит следующий прием: исходная функция f(x) заменятся на рассматриваемом отрезке [a, b] Интерполяционным полиномом Pn (x) И считается, что f'(x). И P'n(x) Примерно равны, т.е. f'(x)=P'n(x), (a ? x ? b)
Для численного дифференцирования используется интерполяционный многочлен с равноотстоящими узлами, так как это значительно упрощает формулы численного дифференцирования.
1.1 Постановка задачи численного дифференцирования
Множество прикладных задач, решаемых средствами классического математического анализа, сводится нахождению производной функции, первообразной или определённого интеграла от заданной функции.
При аналитическом задании указанной функции ее дифференцирования и(или) интегрирование стремится выполнить аналитически.
При табличном здании функции, для которой требуется осуществить дифференцирование или интегрирование, возможности аналитических методов вообще неясны без уточнения постановки задачи.
Как известно, производная аналитически заданной функции f(x) в точке x=a опредеяется следующим образом
f'(a)=
1.2 Метод Эйлера
Простейшим численным методом решения задачи Коши является метод Эйлера, называемый иногда методом ломанных Эйлера.
Угловой коэффициент касательной и интегральной кривой в точке P0(x0,y0) есть
y0'= f (x0, y0). (1)
Найдем ординату y1 касательной, соответствующей абсциссе x1 = x0 + h. Так как уравнение касательной к кривой в точке P0 имеет вид y-y0=y0'(x-x0), то
y1 = y0 + hf (xo, y0) (2)
Угловой коэффициент в точке P1(x1, y1) также находится из данного дифференциального уравнения y1'=f (x1, y1). На следующем шаге получаем новую точку P2(x2, y2), причем
x2 = x1+h (3)
y2 = y1 +hf (x1, y1). (4)
Продолжая вычисления в соответствии с незамеченной схемой, получим формулы Эйлера для m приближенных значений решения задачи Коши начальными данными (x0, y0) на сетке отрезка [a,b] с шагом h:
x1 = xi-1+h, y1 = yi-1 + hf(xi-1, yi-1) (i = 1,2,…, m) (5)
Графической иллюстрацией приближенного решения является ломаная, соединяющая последовательно точки P0, P1, P2,…, Pm, которую называют ломанной Эйлера (Рис.1)
Оценим погрешность метода Эйлера на одном шаге. Для этого запишем разложение точного решения задачи Коши в точке x1 по формуле Тейлора с остаточным членом в форме Лагранжа:
g (x1) = g (x0 + h) = g (x0) + g'(x0) *h + Ѕ g” (z)* h2 = (6)
= y0 + hf (x0, y0) + 1/2g”(z)* h2 = y1 + Ѕ g”(z) *h2, z € (x0, x1).
Погрешность метода на одном шаге имеет порядок h2, так как
d0 = |g(x1)-y1| = Ѕ |g” (z)|* h2<= Ѕ * h2. (7)
После m шагов погрешность вычисления значения ym в конечной точке отрезка возрастает не более чем в m раз. Погрешность метода Эйлера можно оценить неравенством
d Ѕ h2m = Ѕ * (b-a) h (8)
или представить в виде d = Ch, где C € [0, Ѕ
Это означает, что метод Эйлера имеет первый порядок точности. В частности, при уменьшении шага h в 10 раз погрешность уменьшится примерно в 10 раз.
Практическую оценку погрешности решения, найденного на сетке с шагом h/2, в точке x1 € [a, b] производят с помощью приближенного равенства - правила Рунге :
|g (zi) - yi (h/2) (| yi (h) -yi (h/2)|) / (2p - 1)(9)
Где p - порядок точности численного метода. Таким образом, оценка полученного результата по формуле (9) вынуждает проводить вычисления дважды: один раз с шагом h, другой- с шагом h/2.
В основе метода ломанных Эйлера лежит идея графического построения решения дифференциального уравнения, однако этот метод дает одновременно и способ нахождения искомой функции в численной (табличной) форме.
Пусть дано уравнение
y' = f (x, y) (10)
с начальным условием
y(x0) = y0 (11)
(т .е. подавлена задача Коши). Вначале найдем простейшим способом приближенное значение решения в некоторой точке x1 = x0 + h, где h достаточно малый шаг. Заметим, что уравнение (10) совместно с начальным условием (11) задают направление касательной к искомой интегральной кривой в точке M0 (x0, y0). Двигаясь вдоль этой касательной (рисунок 2), получим приближенное значение решения в точке x1:
y1 = y0 + hf(x0, y0) (12)
Приближенным значением в точке M1(x1, y1), можно повторить описанную выше процедуру: построить прямую, проходящую через эту точку под углом, определяемым условием tgв = f (x1, y1), и по ней найти приближенное значение решения в точке x2 = x1 + h. Заметим, что, в отличие от ситуации, изображенной на рисунке 2, эта прямая не есть касательная к реальной интегральной кривой, поскольку точка ~M1 нам недоступна. Однако представляется интуитивно ясным, если h достаточно мало, то получаемые приближения будут близки к точным значениям решения.
Продолжая эту идею, построим систему равноотстоящих точек xi = x0 + ih (i = 0, 1, 2, …, n). Получение таблицы значений искомой функции y(x) по методу Эйлера заключается в циклическом применении пары формул:
yi = hf (xi, yi); yi+1 = yi + yi (i= 0, 1, 2, …, n). (13)
Геометрическая иллюстрация метода Эйлера приведена на рис 3. Вместо интегральной кривой в реальности получается совокупность прямых (так называется ломаная Эйлера).
Методы численного интегрирования дифференциальных уравнений, в которых решение получается от оного узла к другому, называются пошаговыми. Методы Эйлера - простейший представитель семейства пошаговых методов.
Особенностью любого пошагового метода является то, что, начиная со второго шага, исходное значение yi в формуле (13) само является приближенным, т.е., погрешность на каждом следующем шаге систематически возрастает.
Наиболее используемым эмпирическим методом оценки точности, как метода Эйлера, так и других пошаговых методов приближенного численного интегрирования обыкновенных дифференциальных уравнений является способ двойного прохождения заданного отрезка - с шагом h и с шагом h/2.
Метод Эйлера легко распространяется на системы дифференциальных уравнений и на дифференциальные уравнения высших порядков. Последние должны быть предварительно приведены к системе дифференциальных уравнений первого порядка.
1.3 Исправленный метод Эйлера
В исправленном методе Эйлера мы находим средний тангенс наклона касательный для двух точек: xm , ym и xm +h, ym+hy'm . Последняя точка есть та самая, которая в простом методе обозначалась xm+1 , ym+1 . Геометрический процесс нахождения точки xm+1 , ym+1. С помощью метода Эйлера находится точка 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) = Ѕ [f(xm , ym)+f(xm +h, ym +hy'm)], (14)
где ym = f(xm , ym) (15)
Уравнение линии L3 при этом записывается в виде:
y = ym + (x - xm)*F(xm) (16)
так что:
ym+1 = ym+ h*F(xm) (17)
Соотношения 14, 15, 16 и 17 описывают исправленный метод Эйлера. (рис 4)
Рис 4.
1.4 Блок-схема алгоритма
Рис 5.
где A -- начальное значение x, B -- конечное значение x, F(x) -- значение функции в точке xn , N -- количество промежутков, st - выбор операции, C1,C2,C3 - константы для формул, nom - сохраняет номер используемой функции.
На рисунке представлена блок-схема процесса решения дифференциального уравнения методом Эйлера
Подсчитывая каждый раз новое значение уравнения F(x), получаем последовательность значений xn, yn , n=1,2,…
По этим значениям строим график.
2. Численное дифференцирование
Численное дифференцирование применяется в тех случаях, кода функция f(x)задана таблично и методы дифференциального исчисления неприменимы.
В основе численного дифференцирования лежит следующий прием: исходная функция f(x) заменятся на рассматриваемом отрезке [a, b] Интерполяционным полиномом Pn (x) И считается, что f'(x). И P'n(x) Примерно равны, т.е. f'(x)=P'n(x), (a ? x ? b)
Для численного дифференцирования используется интерполяционный многочлен с равноотстоящими узлами, так как это значительно упрощает формулы численного дифференцирования.
Практическая часть
Приложение 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