Реферат: Расчет дифференциального уравнения первого, второго и третьего порядка методом Эйлера
<img src="/cache/referats/3320/image001.gif" v:shapes="_x0000_s1241"> Міністерство освіти України ДАЛПУ Кафедра автоматизаціїтехнологічнихпроцесів і приладобудування
КУРСОВА РОБОТА
з курсу“Математичне моделювання наЕОМ”
на тему“Розв’язок диференціальногорівняння
виду апу(п)+ап-1у(п-1)+…+а1у1+а0у=кх при заданих
початковихумовах з автоматичним вибором кроку
методом Ейлера”
Виконаластудентка групи БА-4-97
Богданова Ольга Олександрівна
Холоденко Вероніка Миколаївна
Перевірила ЗаргунВалентина Василівна
1998
<img src="/cache/referats/3320/image003.jpg" v:shapes="_x0000_s1235">
Блок-схема алгоритма
Блок-схема алгоритма
<img src="/cache/referats/3320/image004.gif" v:shapes="_x0000_s1060">начало
<img src="/cache/referats/3320/image005.gif" v:shapes="_x0000_s1132">
<img src="/cache/referats/3320/image006.gif" v:shapes="_x0000_s1061">у/=f(x,y)
y(x0)=y0
x0,x0+a
<img src="/cache/referats/3320/image007.gif" v:shapes="_x0000_s1126">
<img src="/cache/referats/3320/image008.gif" v:shapes="_x0000_s1064">
h, h/2
<img src="/cache/referats/3320/image007.gif" v:shapes="_x0000_s1127">
<img src="/cache/referats/3320/image009.gif" v:shapes="_x0000_s1067">
k:=0
<img src="/cache/referats/3320/image010.gif" v:shapes="_x0000_s1128 _x0000_s1134 _x0000_s1135">
<img src="/cache/referats/3320/image011.gif" v:shapes="_x0000_s1070"> xk+1/2:=xk+h/2
yk+1/2:=yk+f(xk, yk)h/2
αk:=f(xk+1/2, yk+1/2)
xk+1:=xk+h
yk+1:=yk+αkh
<img src="/cache/referats/3320/image012.gif" v:shapes="_x0000_s1097 _x0000_s1098 _x0000_s1121">
<img src="/cache/referats/3320/image013.gif" v:shapes="_x0000_s1133"><img src="/cache/referats/3320/image014.gif" v:shapes="_x0000_s1099"><img src="/cache/referats/3320/image015.gif" v:shapes="_x0000_s1095"> нет k:=n
<img src="/cache/referats/3320/image016.gif" v:shapes="_x0000_s1129"> да
<img src="/cache/referats/3320/image017.gif" v:shapes="_x0000_s1100"> x0,y0,
x1, y1…
xn,yn
<img src="/cache/referats/3320/image018.gif" v:shapes="_x0000_s1131">
<img src="/cache/referats/3320/image019.gif" v:shapes="_x0000_s1101"> конец
ПОСТАНОВКА ЗАДАЧИИ МЕТОД РЕШЕНИЯ
Решитьдифференциальное уравнение у/=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)<span Times New Roman";mso-hansi-font-family:«Times New Roman»; mso-ansi-language:EN-US;mso-char-type:symbol;mso-symbol-font-family:Symbol">»
yiвычисляются последовательнопо формулам уi+hf(xi, yi) (i=0,1,2…).При этом искомая интегральная кривая у=у(х),проходящая через точку М0(х0, у0), заменяетсяломаной М0М1М2… с вершинами Мi(xi, yi) (i=0,1,2,…);каждое звено МiMi+1 этой ломаной, называемой ломаной Эйлера, имеетнаправление, совпадающее с направлением той интегральной кривой уравнения (1), котораяпроходит через точку Мi.
Если правая часть уравнения (1) внекотором прямоугольнике R{|x-x0|<span Times New Roman";mso-hansi-font-family:«Times New Roman»; mso-ansi-language:EN-US;mso-char-type:symbol;mso-symbol-font-family:Symbol">£
a, |y-y0|<span Times New Roman"; mso-hansi-font-family:«Times New Roman»;mso-ansi-language:EN-US;mso-char-type: symbol;mso-symbol-font-family:Symbol">£b}удовлетворяет условиям:
|f(x, y1)- f(x, y2)| <span Times New Roman";mso-hansi-font-family:«Times New Roman»; mso-ansi-language:EN-US;mso-char-type:symbol;mso-symbol-font-family:Symbol">£
|df/dx|=|df/dx+f(df/dy)|<span Times New Roman";mso-hansi-font-family: «Times New Roman»;mso-ansi-language:EN-US;mso-char-type:symbol;mso-symbol-font-family: Symbol">£
M (M=const),
то имеет место следующая оценка погрешности:
|y(xn)-yn| <span Times New Roman";mso-hansi-font-family:«Times New Roman»; mso-ansi-language:EN-US;mso-char-type:symbol;mso-symbol-font-family:Symbol">£
hM/2N[(1+hN)n-1], (3)где у(хn)-значение точного решенияуравнения(1) при х=хn,а уn — приближенное значение, полученное на n-омшаге.
Формула (3) имеет в основном теоретическоеприменение. На практике иногда оказывается более удобным двойной просчет: сначаларасчет ведется с шагом h, затем шаг дробят и повторный расчет ведется с шагом h/2. Погрешность болееточного значения уn* оценивается формулой
|yn-y(xn)|<span Times New Roman";mso-hansi-font-family:«Times New Roman»; mso-ansi-language:EN-US;mso-char-type:symbol;mso-symbol-font-family:Symbol">»
|yn*-yn|.Метод Эйлералегко распространяется на системы дифференциальных уравнений и надифференциальные уравнения высших порядков. Последние должны бытьпредварительно приведены к системе дифференциальных уравнений первого порядка.
Модифицированный метод Эйлераболее точен.
Рассмотрим дифференциальное уравнение (1) y/=f(x,y)
с начальным условием y(x0)=y0.Разобьем наш участок интегрирования на n
<img src="/cache/referats/3320/image020.gif" v:shapes="_x0000_s1138"> равных частей.На малом участке [x0,x0+h]
у интегральнуюкривую заменим прямой
<img src="/cache/referats/3320/image021.gif" v:shapes="_x0000_s1187"><img src="/cache/referats/3320/image022.gif" v:shapes="_x0000_s1176"><img src="/cache/referats/3320/image023.gif" v:shapes="_x0000_s1175"><img src="/cache/referats/3320/image024.gif" v:shapes="_x0000_s1147"> Nk/ y=y(x) линией. Получаем точку Мк(хк, ук).
<img src="/cache/referats/3320/image025.gif" v:shapes="_x0000_s1219"><img src="/cache/referats/3320/image026.gif" v:shapes="_x0000_s1215"><img src="/cache/referats/3320/image021.gif" v:shapes="_x0000_s1202">
<img src="/cache/referats/3320/image027.gif" v:shapes="_x0000_s1234"><img src="/cache/referats/3320/image028.gif" v:shapes="_x0000_s1231"><img src="/cache/referats/3320/image028.gif" v:shapes="_x0000_s1229"><img src="/cache/referats/3320/image029.gif" v:shapes="_x0000_s1200"><img src="/cache/referats/3320/image030.gif" v:shapes="_x0000_s1188"><img src="/cache/referats/3320/image031.gif" v:shapes="_x0000_s1159"><img src="/cache/referats/3320/image031.gif" v:shapes="_x0000_s1158"><img src="/cache/referats/3320/image031.gif" v:shapes="_x0000_s1157"><img src="/cache/referats/3320/image026.gif" v:shapes="_x0000_s1156"><img src="/cache/referats/3320/image031.gif" v:shapes="_x0000_s1155"> Мк Мк/
<img src="/cache/referats/3320/image030.gif" v:shapes="_x0000_s1199"><img src="/cache/referats/3320/image021.gif" v:shapes="_x0000_s1189"><img src="/cache/referats/3320/image032.gif" v:shapes="_x0000_s1153"> yk+1
<img src="/cache/referats/3320/image029.gif" v:shapes="_x0000_s1213"><img src="/cache/referats/3320/image026.gif" v:shapes="_x0000_s1196"><img src="/cache/referats/3320/image033.gif" v:shapes="_x0000_s1195"><img src="/cache/referats/3320/image026.gif" v:shapes="_x0000_s1174"><img src="/cache/referats/3320/image026.gif" v:shapes="_x0000_s1152"> yk
<img src="/cache/referats/3320/image034.gif" v:shapes="_x0000_s1232"><img src="/cache/referats/3320/image026.gif" v:shapes="_x0000_s1214"><img src="/cache/referats/3320/image026.gif" v:shapes="_x0000_s1212"><img src="/cache/referats/3320/image026.gif" v:shapes="_x0000_s1204">
<img src="/cache/referats/3320/image035.gif" v:shapes="_x0000_s1140"><img src="/cache/referats/3320/image031.gif" v:shapes="_x0000_s1194"> хкхк1/2xk+h=xk1 х
Через Мкпроводим касательную: у=ук=f(xk,yk)(x-xk).
Делимотрезок (хк, хк1) пополам:
xNk/=xk+h/2=xk+1/2
yNk/=yk+f(xk,yk)h/2=yk+yk+1/2
Получаемточку Nk/.В этой точке строим следующую касательную:
y(xk+1/2)=f(xk+1/2,yk+1/2)=αk
Из точки Мкпроводим прямую с угловым коэффициентом αки определяем точкупересечения этой прямой с прямой Хк1. Получаем точку Мк/. В качестве ук+1принимаем ординату точки Мк/.Тогда:
ук+1=ук+αкh
xk+1=xk+h
(4) αk=f(xk+h/2, yk+f(xk,Yk)h/2)
yk=yk-1+f(xk-1,yk-1)h
(4)-рекурентныеформулы метода Эйлера.
Сначалавычисляют вспомогательные значения искомой функции ук+1/2в точках хк+1/2, затем находят значениеправой части уравнения (1) в средней точке y/k+1/2=f(xk+1/2, yk+1/2) и определяют ук+1.
Для оценки погрешности в точке хк проводят вычисления ук с шагомh,затем с шагом 2hи берут 1/3разницы этих значений:
| ук*-у(хк)|=1/3(yk*-yk),
гдеу(х)-точное решение дифференциального уравнения.
Таким образом, методом Эйлера можно решатьуравнения любых порядков. Например, чтобы решить уравнение второго порядка y//=f(y/,y,x) c начальными условиями y/(x0)=y/0,y(x0)=y0, выполняется замена:
y/=z
z/=f(x,y,z)
Темсамым преобразуются начальные условия:y(x0)=y0, z(x0)=z0, z0=y/0.
РЕШЕНИЕ КОНТРОЛЬНОГО ПРИМЕРА
Приведем расчет дифференциального уравнения первого,второго и третьего порядка методомЭйлера
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)+α0h; y(x1/2)=y(x0)+f(x0,y0)h/2;
f(x0,y0)=2<span Times New Roman";mso-hansi-font-family:«Times New Roman»; mso-ansi-language:EN-US;mso-char-type:symbol;mso-symbol-font-family:Symbol">*
0-1=-1y(x1/2)=1-1<span Times New Roman";mso-hansi-font-family:«Times New Roman»; mso-ansi-language:EN-US;mso-char-type:symbol;mso-symbol-font-family:Symbol">*
0,1=0,9α0=2<span Times New Roman";mso-hansi-font-family:«Times New Roman»; mso-ansi-language:EN-US;mso-char-type:symbol;mso-symbol-font-family:Symbol">*
0,1-0,9=-0,7y1=1-0,1<span Times New Roman";mso-hansi-font-family:«Times New Roman»; mso-ansi-language:EN-US;mso-char-type:symbol;mso-symbol-font-family:Symbol">*
0,2=0,862). y(x2)=y(x1)+α1h; 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<span Times New Roman";mso-hansi-font-family:«Times New Roman»; mso-ansi-language:EN-US;mso-char-type:symbol;mso-symbol-font-family:Symbol">*
0,2-0,86=-0,46y(x1+1/2)=0,86-0,46<span Times New Roman";mso-hansi-font-family:«Times New Roman»; mso-ansi-language:EN-US;mso-char-type:symbol;mso-symbol-font-family:Symbol">*
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
2. Даноуравнение второго порядка:
y//=2x-y+y/
Находим решение на том жеотрезке [0,1] c шагом h=0,2;
Замена: y/=z
z/=2x-y+z
Начальные условия: у0=1
z0=1
1).x1=0,2; x1/2=0,1
y(z1)=y(z0)+α0h z(x1,y1)=z(x0,y0)+β0h
y(z1/2)=y(z0)+f(z0,y0)h/2 z(x1/2,y1/2)=z(x0,y0)+f(x0,y0,z0)h/2
f(z0,y0)=f10=1 f(x0,y0,z0)=f20=2*0-1+1=0
y1/2=1+1*0,1=1,1 z1/2=1+0*0,1=1
α0=z0=1 β0=2*0,1-1,1+1=0,1
y1=1+0,2*1=1,2 z1=1+0,2*0,1=1,02
2).x2+0,4; x1+1/2=0,3
f11=z1=1,02 f21=2*0,2-1,2+1,02=0,22
y1+1/2=1,2+1,02*0,1=1,1 z1+1/2=1,02+0,22*0,1=1,042
α1=z1+1/2=1,042 β1=2*0,3-1,302+1,042=0,34
y2=1,2+1,042*0,2=1,4084 z2=1.02+0,34*0,2=1,088
3).x3=0,6; x2+1/2=0,5
f12=z2=1,088 f22=2*0,4-1,4084+1,088=0,4796
y2+1/2=1,4084+1,088*0,1=1,5172 z2+1/2=1,088+0,4796*0,1=1,13596
α2=z2+1/2=1,13596 β2=2*0,5-1,5172+1,13596=0,61876
y3=1,4084+1,136*0,2=1,635592 z3=1,088+0,61876*0,2=1,211752
4).x4=0,8; x3+1/2=0,7
f13=z3=1,211752 f23=2*0,6-1,636+1,212=0,77616
y3+1/2=1,636+1,212*0,1=1,7567672 z3+1/2=1,212+0,776*0,1=1,289368
α3=z3+1/2=1,289368 β3=2*0,7-1,7568+1,289=0,9326008
y4=1,6+1,289*0,2=1,8934656 z4=1,212+0,93*0,2=1,39827216
5).x5=1; y4+1/2=0,9
f14=z4=1,39827216 f24=2*0,8-1,893+1,398=1,10480656
y4+1/2=1,893+1,398*0,1=2,0332928 z4+1/2=1,398+1,105*0,1=1,508752816
α4=z4+1/2=1,508752816 β4=2*0,9-2,03+1,5=1,27546
y5=1,893+1,5*0,2=2,195216163 z5=1,398+1,275*0,2=1,65336416
3. Чтобырешитьуравнение третьего порядка
y///=2x-y-y/+y//
наотрезке [0,1], с шагом h=0,2 и начальными условиями
y0//=1
y0/=1
y0=1
необходимосделать 3 замены: y/=a y0/=a0=1
y//=a/=b y0//=b0=1
b/=2x-y-a+b
1).x1=0,2; x1/2=0,1
y(a1)=y(a0)+a0h y(a1/2)=y(a0)+f10h/2
a(b1)=a(b0)+β0h a(b1/2)=a(b0)+f20h/2
b(x1,y1,a1)=b(x0,y0,a0)+γ0h b(x1/2,y1/2,a1/2)=b(x0,y0,a0)+f30h/2
f10=f(a0,y(a0))=1 y1/2=1+1*0,1=1,1
f20=f(b0,a(b0))=1 a1/2=1+1*0,1=1,1
f30=f(x0,y0,a0,b0)=-1 b1/2=1-1*0,1=0,9
α0=a1/2=1,1 y(a1)=1+1,1*0,2=1,22
β0=b1/2=0,9 a(b1)=1+0,9*0,2=1,18
γ0=2*0,1-1,1-1,1+0,9=-1,1 b(x1,y1,a1)=1-1,1*0,2=0,78
2).x2=0,4; x1+1/2=x1+h/2=0,3
f11=a1=1,18 y1+1/2=1,22+1,18*0,1=1.338
f21=b1=0,78 a1+1/2=1,18+0,78*0,1=1,258
f31=2*0,2-1,22-1,18+0,78=-1,22 b1+1/2=-1,22*0,1+0,78=0,658
α1=a1+1/2=1,258 y2=1,22+1,258*0,2=1,4716
β1=b1+1/2=0,658 a2=1,18+0,658*0,2=1,3116
γ1=2*0,3-1,338-1,258+0,658=-1,338 b2=0,78-1,338*0,2=0,5124
3).x3=0,6; x2+1/2=0,5
f12=a2=1,3116 y2+1/2=1,47+1,3*0,1=1,60276
f22=b2=0,5124 a2+1/2=1,3116+0,5*0,1=1.36284
f32=2*0,4-1,47-1,31+0,512=-1,4708 b2+1/2=0,4-1,4*0,1=0,36542
α2=1,36284 y3=1,4716+1,3116*0,2=1,744168
β2=0,36542 a3=1,3116+0,3654*0,2=1,384664
γ2=2*0,5-1,6-1,36+0,365=-1,60018 b3= 0,51-1,60018*0,2=0,192364
4).x4=0,8; x3+1/2=0,7
f13=1,384664 y3+1/2=1,74+1,38*0,1=1,8826364
f23=0,192364 a3+1/2=1,38+0,19*0,1=1,4039204
f33=2*0,6-1,7-1,38+0,19=-1,736488 b3+1/2=0,19-1,7*0,1=0,0187152
α3=1,4039204 y4=1,74+1,4*0,2=2,0249477
β3=0,0187152 a4=1,38+0,9187*0,2=1,388403
γ3=2*0,7-1,88-1,4+0,0187=-1,8678416 b4=0,192-1,87*0,2=-0,1812235
5).x4=1; x4+1/2=0,9
f14=1,388403 y4+1/2=2,02+1,388*0,1=2,16379478
f24=-0,1812235 a4+1/2=1,4-0.181*0,1=1,370306608
f34=2*0,8-2,02-1,388-0,18=-1,9945834 b4+1/2=-0,18-1,99*0,1=-0,38066266
α4=1,3703 y5=2,02+1,37*0,2=2,2990038
β4=-0,38066 a5=1,388-0,38*0,2=1,3122669
γ4=2*0,9-2,16-1,37-0,38=-2,114764056 b5=-0,181-2,1*0,2=-0,6041734
Программа на Turbo Pascal
usescrt,pram,kurs1_1;
var
yx,xy,l,v,p,ff,ay,by,x:array[0..10] of real;
y,a,b:array[0..10,0..1]of real;
i,n,o:integer;
c,d,h,k:real;
label
lap1;
begin
screen1;
clrscr;
writeln('введите наивысший порядок производной не большетрех ');
readln(n);
ifn=0 thenbegin
writeln('это прямолинейная зависимость и решается безметода Эйлера ');
gotolap1;end;
writeln('введите коэффициенты {a0,a1}');
fori:=0 ton do
readln(l[i]);
if(n=1) and (l[1]=0) or (n=2) and (l[2]=0) or (n=3) and (l[3]=0) then begin
writeln('делениена ноль');
goto lap1;
end;
writeln('введите коэффициент при x');
readln(k);
writeln('введите отрезок ');
readln(c,d);
o:=5;
h:=abs(d-c)/o;
writeln('шаг=',h:1:1);
writeln('задайте начальные условия y(x)= ');
fori:=0 ton-1 do
readln(v[i]);
ifn=3 thenbegin
yx[0]:=v[0];
ay[0]:=v[1];
by[0]:=v[2];
p[0]:=(k*c-l[0]*v[0]-l[1]*v[1]-l[2]*v[2])/l[3];
x[0]:=c;
gotoxy(32,1);
<img src="/cache/referats/3320/image026.gif" v:shapes="_x0000_s1254"><img src="/cache/referats/3320/image026.gif" v:shapes="_x0000_s1252"><img src="/cache/referats/3320/image034.gif" v:shapes="_x0000_s1248"><img src="/cache/referats/3320/image026.gif" v:shapes="_x0000_s1245"><img src="/cache/referats/3320/image034.gif" v:shapes="_x0000_s1244"><img src="/cache/referats/3320/image036.gif" v:shapes="_x0000_s1243"> write(' ');
gotoxy(32,2);
<img src="/cache/referats/3320/image037.gif" v:shapes="_x0000_s1256"><img src="/cache/referats/3320/image032.gif" v:shapes="_x0000_s1253"><img src="/cache/referats/3320/image032.gif" v:shapes="_x0000_s1251"><img src="/cache/referats/3320/image038.gif" v:shapes="_x0000_s1247"><img src="/cache/referats/3320/image039.gif" v:shapes="_x0000_s1246"> write(' x y a b ');
gotoxy(32,3);
<img src="/cache/referats/3320/image026.gif" v:shapes="_x0000_s1340"><img src="/cache/referats/3320/image034.gif" v:shapes="_x0000_s1339"><img src="/cache/referats/3320/image034.gif" v:shapes="_x0000_s1336"><img src="/cache/referats/3320/image034.gif" v:shapes="_x0000_s1335"><img src="/cache/referats/3320/image034.gif" v:shapes="_x0000_s1332"> write('',c:7:7,' ',yx[0]:7:7,' ',ay[0]:7:7,' ',by[0]:7:7,' ');
for i:=0 to o-1 dobegin
x[i]:=x[i]+h/2;
y[i,1]:=yx[i]+(h/2)*ay[i];
a[i,1]:=ay[i]+(h/2)*by[i];
b[i,1]:=by[i]+(h/2)*p[i];
ff[i]:=(k*x[i]-l[0]*y[i,1]-l[1]*a[i,1]-l[2]*b[i,1])/l[3];
xy[i]:=x[i]+h/2;
yx[i+1]:=yx[i]+h*a[i,1];
ay[i+1]:=ay[i]+h*b[i,1];
by[i+1]:=by[i]+h*ff[i];
x[i+1]:=x[i]+h/2;
p[i+1]:=(k*xy[i]-l[0]*yx[i+1]-l[1]*ay[i+1]-l[2]*by[i+1])/l[3];
end;
for i:=0 to o-1 dobegin
gotoxy(32,4+i);
<img src="/cache/referats/3320/image026.gif" v:shapes="_x0000_s1275"><img src="/cache/referats/3320/image026.gif" v:shapes="_x0000_s1274"><img src="/cache/referats/3320/image038.gif" v:shapes="_x0000_s1273"><img src="/cache/referats/3320/image040.gif" v:shapes="_x0000_s1272"><img src="/cache/referats/3320/image040.gif" v:shapes="_x0000_s1270"> write('',xy[i]:7:7,' ',yx[i+1]:7:7,' ',ay[i+1]:7:7,' ',by[i+1]:7:7,' ');
end;
gotoxy(32,4+o);
<img src="/cache/referats/3320/image026.gif" v:shapes="_x0000_s1269"><img src="/cache/referats/3320/image026.gif" v:shapes="_x0000_s1268"><img src="/cache/referats/3320/image026.gif" v:shapes="_x0000_s1267"><img src="/cache/referats/3320/image034.gif" v:shapes="_x0000_s1266"><img src="/cache/referats/3320/image026.gif" v:shapes="_x0000_s1262"><img src="/cache/referats/3320/image041.gif" v:shapes="_x0000_s1261"> write(' ');
end;
ifn=2 thenbegin
x[0]:=c;
yx[0]:=v[0];
ay[0]:=v[1];
p[0]:=(k*c-l[0]*yx[0]-l[1]*v[1])/l[2];
gotoxy(32,1);
<img src="/cache/referats/3320/image026.gif" v:shapes="_x0000_s1280"><img src="/cache/referats/3320/image034.gif" v:shapes="_x0000_s1279"><img src="/cache/referats/3320/image026.gif" v:shapes="_x0000_s1278"><img src="/cache/referats/3320/image034.gif" v:shapes="_x0000_s1277"><img src="/cache/referats/3320/image042.gif" v:shapes="_x0000_s1276"> write(' ');
gotoxy(32,2);
<img src="/cache/referats/3320/image034.gif" v:shapes="_x0000_s1292"><img src="/cache/referats/3320/image026.gif" v:shapes="_x0000_s1290"><img src="/cache/referats/3320/image034.gif" v:shapes="_x0000_s1288"><img src="/cache/referats/3320/image026.gif" v:shapes="_x0000_s1286"> write(' x y a ');
gotoxy(32,3);
<img src="/cache/referats/3320/image033.gif" v:shapes="_x0000_s1296"><img src="/cache/referats/3320/image033.gif" v:shapes="_x0000_s1295"><img src="/cache/referats/3320/image033.gif" v:shapes="_x0000_s1294"><img src="/cache/referats/3320/image040.gif" v:shapes="_x0000_s1293"> write('',c:7:7,' ',yx[0]:7:7,' ',ay[0]:7:7,' ');
for i:=0 to o-1 dobegin
x[i]:=x[i]+h/2;
y[i,1]:=yx[i]+(h/2)*ay[i];
a[i,1]:=ay[i]+(h/2)*p[i];
ff[i]:=(k*x[i]-l[0]*y[i,1]-l[1]*a[i,1])/l[2];
xy[i]:=x[i]+h/2;
yx[i+1]:=yx[i]+h*a[i,1];
ay[i+1]:=ay[i]+h*ff[i];
x[i+1]:=x[i]+h/2;
p[i+1]:=(