ПРИЛОЖЕНИЕ:
ОТВЕТ(задание 3)
Вспомогательные функции:
function fun1(t,y1,y2:real):real; //дано
begin
fun1:=(y2-y1)*t;
end;
function fun2(t,y1,y2:real):real; //дано
begin
fun2:=(y2+y1)*t;
end;
//функция реализующая формулы (12), все обозначения совпадают
function tform1.yh1(h,t0,y10,y20:real):real;
var
k: array [1..4,1..2] of real;
begin
k[1,1]:=h*fun1(t0,y10,y20);
k[1,2]:=h*fun2(t0,y10,y20);
k[2,1]:=h*fun1(t0+(1/2)*h,y10+(1/2)*k[1,1],y20+(1/2)*k[1,2]);
k[2,2]:=h*fun2(t0+(1/2)*h,y10+(1/2)*k[1,1],y20+(1/2)*k[1,2]);
k[3,1]:=h*fun1(t0+(1/2)*h,y10+(1/2)*k[2,1],y20+(1/2)*k[2,2]);
k[3,2]:=h*fun2(t0+(1/2)*h,y10+(1/2)*k[2,1],y20+(1/2)*k[2,2]);
k[4,1]:=h*fun1(t0+h,y10+k[3,1],y20+k[3,2]);
k[4,2]:=h*fun2(t0+h,y10+k[3,1],y20+k[3,2]);
yh1:=y10+(1/6)*(k[1,1]+2*k[2,1]+2*k[3,1]+k[4,1]);
//yh2:=y20+(1/6)*(k[1,2]+2*k[2,2]+2*k[3,2]+k[4,2]);
end;
//построение полинома Эрмита по 2 точкам и 4-ём значениям в этих точках; см. формулу (13)
function tform1.ermit(t0,t1,ft0,ft1,dft0,dft1,trm:real):real;
begin
ermit:=(sqr(trm-t1)/sqr(t0-t1))*(ft0*(1-2*(trm-t0)/(t0-t1))+dft0*(trm-t0))+
(sqr(trm-t0)/sqr(t1-t0))*(ft1*(1-2*(trm-t1)/(t1-t0))+dft1*(trm-t1));
end;
//процедура решения ЗК по МТРК с автоматическим выбором переменного шага интегрирования по заданной величине Е ЛП, поиск значений в точках выдачи результата
procedure TForm1.Button5Click(Sender: TObject);
var
kol,i,j,k,m:integer;
levgran, pravgran,l,t0:real;
y,y2:array [0..1000, 1..2] of real;
sum,h,x,tr:array [0..1000] of real;
begin
kol:=strtoint(edit3.text); //kol-количество точек разбиения
levgran:=strtofloat(edit1.text); //levgran - t0
pravgran:=strtofloat(edit2.text); //pravgran - t0+T
l:=(pravgran-levgran); //длина промежутка = Т
m:=strtoint(edit13.text); //m-кол-во точек выдачи результата
for i:=1 to 1000 do h[i]:=l/kol; //h[i] - начальный шаг интегрирования
stringgrid1.RowCount:=kol+1; //размер таблицы выдачи результата
stringgrid1.Cols[0][0]:=edit4.text; //заполняем первую строку таблицы (это Н.У.З.К.)
stringgrid1.Cols[1][0]:=edit5.text; //заполняем первую строку таблицы (это Н.У.З.К.)
stringgrid1.Cols[2][0]:=edit7.text; //заполняем первую строку таблицы (это Н.У.З.К.)
t0:=strtofloat(edit4.text); //t0
y[0,1]:=strtofloat(edit5.text); //Н.У.
y[0,2]:=strtofloat(edit7.text); //Н.У.
y2[0,1]:=y[0,1]; //Н.У.
y2[0,2]:=y[0,2]; //Н.У.
for i:=1 to m do tr[i]:=t0+(l/m)*i; //фиксируем точки выдачи результатов(tr[i])
i:=0;
x[0]:=t0;
repeat //начало цикла реализующего поиск приближенного решения
i:=i+1;//переход к следующей точке в которой будет происходить поиск приб. реш.
//далее цикл из двух шагов: 1-ый нахождение y[i] ч/з y[i-1] c заданным наперёд начальным шагом инт-я h; 2-ой нахождение y[i] ч/з y[i-1], но с шагом, который обеспечивает точность Е ЛП
for j:=0 to 1 do begin
y[i,1]:=yh1(h[i],x[i-1],y[i-1,1],y[i-1,2]); //находим по МТРК y[i] (П5); см. что такое yh1
y[i,2]:=yh2(h[i],x[i-1],y[i-1,1],y[i-1,2]); //находим по МТРК y[i] (П5); см. что такое yh2
if j<>1 then begin
y2[i,1]:=yh1(h[i]/2,x[i-1],y[i-1,1],y[i-1,2]); //находим по МТРК y2[i] в точке x[i-1]+h/2 (П6)
y2[i,2]:=yh2(h[i]/2,x[i-1],y[i-1,1],y[i-1,2]);
y2[i,1]:=yh1(h[i]/2,x[i-1]+h[i]/2,y2[i,1],y2[i,2]); //находим по МТРК y2[i] в точке x[i-1]+h/2+h/2 (П6)
y2[i,2]:=yh2(h[i]/2,x[i-1]+h[i]/2,y2[i,1],y2[i,2]);
sum[i]:=sqrt(sqr(y2[i,1]-y[i,1])+sqr(y2[i,2]-y[i,2])); //модуль разности для формулы (11)
h[i]:=(h[i]/2)*exp(ln((15*exp(ln(10)*strtofloat(edit12.text))/sum[i]))*(1/4)); //формула (11)
end;
end;
x[i]:=x[i-1]+h[i]; //новая точка в которой нам предстоит искать приближ. значение реш-я
for k:=1 to m do begin //цикл для нахождения значения полинома Эрмита в точке выдачи результата
if (x[i-1]<tr[k])and(tr[k]<x[i]) then begin //проверка, если точка выдачи оказывается м/ду точками разбиения
stringgrid2.Cols[0][k-1]:=floattostr(tr[k]); //заполнение результата в таблицу
stringgrid2.Cols[1][k-1]:=floattostr(ermit(x[i-1],x[i],y[i-1,1], y[i,1],fun1(x[i-1],y[i-1,1],y[i-1,2]),fun1(x[i],y[i,1],y[i,2]),tr[k])); //заполнение результата в таблицу; см. function ermit
stringgrid2.Cols[2][k-1]:=floattostr(ermit(x[i-1],x[i],y[i-1,2], y[i,2],fun2(x[i-1],y[i-1,1],y[i-1,2]),fun2(x[i],y[i,1],y[i,2]),tr[k]));
end;
if tr[k]=x[i] then begin //проверка, если точка выдачи совпадлает с точкой разбиения
stringgrid2.Cols[1][k-1]:=floattostr(y[i,1]);
stringgrid2.Cols[2][k-1]:=floattostr(y[i,2]);
end
end;
until x[i]>t0+l; //завершения основного цикла, происходит в случае выхода за пределы заданного промежутка интегрирования
//Далее: на всякий случай заполним таблицу со значениями во всех точках дробления
stringgrid1.RowCount:=i+1;
for j:=1 to i do
begin
stringgrid1.Cols[0][j]:=floattostr(x[j]);
stringgrid1.Cols[1][j]:=floattostr(y[j,1]);
stringgrid1.Cols[2][j]:=floattostr(y[j,2]);
end;
end;
Рассмотрим процедуру, реализующую нахождение глобальной погрешности.
procedure TForm1.Button1Click(Sender: TObject);
var
kol,i,j:integer;
levgran, pravgran,h,l,t0,max:real;
y,y2,gp,gp2:array [0..10000, 1..2] of real;
sum:array [0..10000] of real;
begin
kol:=strtoint(edit3.text);
levgran:=strtofloat(edit1.text);
pravgran:=strtofloat(edit2.text);
l:=(pravgran-levgran);
h:=l/kol;
stringgrid1.RowCount:=kol+1;
stringgrid1.Cols[0][0]:=edit4.text;
stringgrid1.Cols[1][0]:=edit5.text;
stringgrid1.Cols[2][0]:=edit7.text;
t0:=strtofloat(edit4.text);
y[0,1]:=strtofloat(edit5.text);
y[0,2]:=strtofloat(edit7.text);
y2[0,1]:=y[0,1];
y2[0,2]:=y[0,2];
//всё что выше повторяет предыдущую процедуру
max:=1;
//далее, если убрать комментарий, программа будет находить приближенное решение с наперёд заданной точностью, но выбор шага интегрирования будет происходить не автоматически, а по итогам оценки ГП
//while max>exp(ln(10)*strtofloat(edit11.text)) do begin
//алгоритм нахождения ГП полностью соответствует (П4)
for i:=1 to kol do //1-ое: ищем значения y в точках xn с шагом h
begin
y[i,1]:=yh1(h,t0+h*(i-1),y[i-1,1],y[i-1,2]);
y[i,2]:=yh2(h,t0+h*(i-1),y[i-1,1],y[i-1,2]);
end;
for i:=1 to kol*2 do //2-ое: ищем значения y в точках xn с шагом h/2
begin
y2[i,1]:=yh1((h/2),t0+(h/2)*(i-1),y2[i-1,1],y2[i-1,2]);
y2[i,2]:=yh2((h/2),t0+(h/2)*(i-1),y2[i-1,1],y2[i-1,2]);
end;
for i:=1 to kol do //подсчитываем ГП в каждой точке дробления
begin
gp[i,1]:=(y2[i*2,1]-y[i,1])/(1-1/16);
gp[i,2]:=(y2[i*2,2]-y[i,2])/(1-1/16);
gp2[i,1]:=(y2[i*2,1]-y[i,1])/(16-1);
gp2[i,2]:=(y2[i*2,2]-y[i,2])/(16-1);
sum[i]:=sqrt(sqr(gp2[i,1])+sqr(gp2[i,2])); //фиксируем модуль ГП в каждой точке
end;
max:=sum[1];
for i:=1 to kol-1 do //выбираем максимальную ГП
begin
if max<sum[i+1] then begin j:=i; max:=sum[i+1]; end;
end;
stringgrid1.RowCount:=kol+1;
h:=h/2;
kol:=kol*2;
//end;
//далее аналогично предыдущей процедуре
edit6.Text:=floattostr(h*2);
edit8.Text:=floattostr(max);
edit9.Text:=floattostr(t0+h*2*j);
for i:=1 to trunc(kol/2) do
begin
stringgrid1.Cols[0][i]:=floattostr(t0+h*2*i);
stringgrid1.Cols[1][i]:=floattostr(y2[i*2,1]+gp2[i,1]);
stringgrid1.Cols[2][i]:=floattostr(y2[i*2,2]+gp2[i,2]);
end;
end;
Уважаемый посетитель!
Чтобы распечатать файл, скачайте его (в формате Word).
Ссылка на скачивание - внизу страницы.