Полином Гира (экстраполяция методом Гира)
{****************************************************************************}
{ }
{ Вариант №1.23 Экстраполяция методом Гира гр. ##-### Ф.И.О. }
{ }
{ P^2 P = 50-200 [Ом], }
{ f(w) = ----- * ( 1 + Q^2 * (2*dw/wp))^(-1/2), Rn= 0.1-5 [Ом], }
{ Rn Q = 50-1000, }
{ wp- частота резонанса, }
{ dw- шаг, }
{ }
{****************************************************************************}
{ Аппроксимация
Аппроксимацией называется замена одних математических объектов
(чисел, функций) другими, более простыми в том или ином смысле и
близкими к исходным.
При решении задачи аппроксимации предполагается, что известны
значнения аппроксимируемой функции f(x) (снятые эксперементально или
расчитанные по формуле) в узловых точках (обычно их называют узлами
интерполяции) X1, X2, ..., Xn. Необходимо найти функцию P(x), с помощью
которой можно было бы определить значения исходной функции для других
(отличных от узлов) значений аргумента x с некоторой достаточно малой
погрешностью.
Различают две задачи аппроксимации: интерполяцию и экстраполяцию.
Интерполяция - это отыскание промежуточных значений функции внутри
заданного интервала [X1, Xк] для значений x не совпадающих с узловыми
точками.
Экстраполяция - это распространение выводов, полученных из наблюдения
над одной частью явления на другую ее часть, то есть это определение
значений функции вне заданного интервала [Xн, Xк]. }
PROGRAM GIR;
{$N-}
USES Crt;
CONST
Max_kol = 25; {Максимальное количество узлов}
Tochka = 10; {Количество знаков после запятой}
TYPE
_Real = Real;
_INTEGER = Byte;
_Mas = array [1..Max_kol] of _Real;
VAR
H, F_up, F_down, eps_pol, eps_tek : _Real;
Xpol, Ypol, V : _Mas;
n, m : _INTEGER; {n-количество узлов}
Pr : Char; {переключатель}
{----------------------------------------------------------------------------}
{Расчет табличных данных для заданной функции}
PROCEDURE TABL_FUNC(Var Xpol, Ypol: _Mas; Var n: _INTEGER; Var dw, eps_pol: _Real);
VAR
P, Rn, Q, wp, wn : _Real;
I : _INTEGER; {счетчик}
BEGIN
Write('Сопротивление P ............. = '); ReadLn(P);
Write('Сопротивление нагрузки Rn ... = '); ReadLn(Rn);
Write('Параметр Q .................. = '); ReadLn(Q);
Write('Начальна частота wp ......... = '); ReadLn(wn);
Write('Шаг изменения частота dw .... = '); ReadLn(dw);
Write('Количество узлов (< ', Max_kol:2, ') ..... = '); ReadLn(n);
Write('Точность eps_pol ............ = '); ReadLn(eps_pol);
WriteLn;
{расчет функции в узлах}
wp:= wn;
FOR I:= 1 TO (n+1) DO
begin
Xpol[I]:= wp;
Ypol[I]:= Sqr(P) / (Rn * Sqrt(1 + Sqr(Q) * Sqr(2*dw/wp)));
wp:= wp + dw;
end;
END;
{--------------------------< Конец ПП TABL_FUNC >----------------------------}
{Ввод значений аппроксимируемой функции}
PROCEDURE VVOD_TABL(Var Xpol, Ypol: _Mas; Var n: _INTEGER; Var dw, eps_pol: _Real);
VAR
wp, wn : _Real;
I : _INTEGER; {счетчик}
BEGIN
Write('Начальный узел ............ = '); ReadLn(wn);
Write('Шаг изменения узла ........ = '); ReadLn(dw);
Write('Количество узлов (< ', Max_kol:2, ') ... = '); ReadLn(n);
Write('Погрешность eps .............. = '); ReadLn(eps_pol);
WriteLn;
{ввод значений функции в узлах}
WriteLn(' Введите значение функции в узлах:');
wp:= wn;
FOR I:= 1 TO n DO
begin
Xpol[I]:= wp;
Write('Y(', wp:Tochka, ')= '); ReadLn(Ypol[I]);
wp:= wp + dw;
end;
{задание (n+1)-узла}
Xpol[n+1]:= wp;
END;
{---------------------------< Конец ПП VVOD_TABL >---------------------------}
{Расчет коэффициентов полинома}
PROCEDURE Kff_pol(X: _Mas; H: _Real; n, m: _INTEGER; Var V: _Mas);
TYPE
_Matric = array [1..Max_kol] of _Mas;
_Vector = array [1..Max_kol] of Byte;
VAR
G : _Matric;
P : _Vector;
I, J : _INTEGER; {счетчики}
{----------------------------------------------------------------------------}
{рашение матричной системы уравнений, методом "LU-разложений"}
PROCEDURE GVP(G: _Matric; P: _Vector; m: _INTEGER; Var V: _Mas);
VAR
K, I, J : _INTEGER; {счетчики}
L, U : _Matric; {вспомогательные матрицы}
SOS : _Mas; {вспомогательный вектор}
alf, S : _Real; {вспомогательные переменные}
BEGIN
{обнуление элементов вспомогательных матриц L и U}
FOR I:= 1 TO m DO
For J:= 1 To m Do
begin
L[I, J]:= 0; { L[I, J] = 0, где I=1..m, J=1..m }
U[I, J]:= 0; { U[I, J] = 0, где I=1..m, J=1..m }
end;
{разложение матрицы G на матрицу L и U}
FOR K:= 1 TO m DO
begin
{формирование матрицы L}
For I:= K To m Do
begin
{накопление суммы элементов в S}
S:= 0;
for J:= 1 to (K - 1) do
begin { k-1 }
alf:= L[I, J] * U[J, K]; { S = E (L[I, J] * U[J, K]), }
S:= S + alf; { j=1 }
end;
L[I, K]:= G[I, K] - S; { L[I, K] = G[I, K] - S, где k=1..m, i=k..m }
end;
{формирование матрицы U}
U[K, K]:= 1;
For I:= (K + 1) To m Do
begin
{накопление суммы элементов в S}
S:= 0;
for J:= 1 to (K - 1) do
begin { k-1 }
alf:= L[K, J] * U[J, I]; { S = E (L[K, J] * U[J, I]), }
S:= S + alf; { j=1 }
end;
{ G[I, K] - S }
alf:= G[I, K] - S; { U[K, I] = -------------, где k=1..m, i=(k+1)..m }
U[K, I]:= alf/L[K, K]; { L[K, K] }
end;
end;
{вычисление вспомогательного вектора SOS}
FOR K:= 1 TO m DO
begin
{накопление суммы элементов в S}
S:= 0;
For I:= 1 To (K - 1) Do
begin { k-1 }
alf:= L[K, I] * SOS[I]; { S = E (L[K, I] * SOS[I]), }
S:= S + alf; { i=1 }
end;
{ P[K] - S }
alf:= P[K] - S; { SOS[K] = ----------, где k=1..m }
SOS[K]:= alf/L[K, K]; { L[K, K] }
end;
{вычисление вектора V}
FOR K:= m DOWNTO 1 DO
begin
{накопление суммы элементов в S}
S:= 0;
For I:= (K + 1) To m Do
begin { m }
alf:= U[K, I] * V[I]; { S = E (U[K, I] * V[I]), }
S:= S + alf; { i=k+1 }
end;
V[K]:= SOS[K] - S; { V[K] = SOS[K] - S, где k=1..m }
end;
END;
{------------------------------< конец ПП GVP >------------------------------}
BEGIN
{формирование матрицы G}
FOR I:= 1 TO m DO { | X[n+1] - X[n - (m-1)] | (m-1) }
For J:= 1 To m Do { G[I, J] = |-----------------------| }
G[I, J]:= exp( (I-1)*Ln((X[n+1] - X[n-(J-1)]) / H) ); { | H | }
{формирование вектора P = [1, 0, 0, 0, ..., 0]}
P[1]:= 1;
FOR I:= 2 TO m DO P[I]:= 0;
{обращение к ПП GVP}
GVP(G, P, m, V);
END;
{---------------------------< Конец ПП Kff_pol >-----------------------------}
{Расчет полинома}
FUNCTION F_pol(V, Y: _Mas; n, m: _INTEGER): _Real;
VAR
I: _INTEGER; {счетчик}
S: _Real; {сумматор}
BEGIN
S:= 0;
FOR I:= 1 TO m DO S:= S + V[I]*Y[n+1-I];
F_pol:= S;
END;
{-----------------------------< Конец ПП F_pol >-----------------------------}
{Вывод значения функции в (n+1) узле}
PROCEDURE VIVOD1(Xpol, Ypol, eps_tek: _Real; n, m: _INTEGER);
BEGIN
WriteLn;
WriteLn;
WriteLn('Значение функции в (n+1) узле: Y(', Xpol:Tochka, ')=', Ypol:Tochka);
WriteLn('Количество используемых узлов', m:3, ' из', n:3);
WriteLn('Погрешность расчета =', eps_tek:Tochka);
END;
{-----------------------------< Конец ПП VIVOD1 >----------------------------}
{Вывод значения функции в (n+1) узле для тестовой функции}
PROCEDURE VIVOD2(Xpol, Ypol, eps_tek, Y: _Real; n, m: _INTEGER);
VAR
eps: _Real;
BEGIN { Y_полином - Y_точное }
eps:= Abs((Ypol - Y)/Y); { EPS_точное = ---------------------- }
{ Y_точное }
WriteLn;
WriteLn(' Значение функции в узле X(n+1)=', Xpol:Tochka , ' :');
WriteLn('Y полинома ... =', Ypol:Tochka);
WriteLn('Y точное ..... =', Y:Tochka);
WriteLn(' Погрешность расчета :');
WriteLn('текущая ...... =', eps_tek:Tochka);
WriteLn('точная ....... =', eps:Tochka);
WriteLn;
WriteLn('Количество используемых узлов', m:3, ' из', n:3);
END;
{-----------------------------< Конец ПП VIVOD2 >----------------------------}
{основной блок программы}
BEGIN
ClrScr;
WriteLn(' ExtraPol_GIRA 1.00 Copyright (c) 1999-00 DGur.');
WriteLn;
WriteLn;
{Ввод исходных данных}
Write('Тестовый расчет? [Y/N]'); ReadLn(Pr);
WriteLn;
{Ввод значений функции в узлах}
IF (UpCase(Pr) = 'Y')
{обращение к ПП TABL_FUNC}
THEN TABL_FUNC(Xpol, Ypol, n, H, eps_pol)
{обращение к ПП VVOD_TABL}
ELSE VVOD_TABL(Xpol, Ypol, n, H, eps_pol);
{итерационный процесс}
m:= 2;
F_down:= 0;
REPEAT
{обращение к ПП Kff_pol}
Kff_pol(Xpol, H, n, m, V);
{обращение к ПФ F_pol}
F_up:= F_pol(V, Ypol, n, m);
{точность расчета}
IF Not(F_up = 0) THEN eps_tek:= Abs((F_up - F_down)/F_up);
F_down:= F_up;
Inc(m);
UNTIL (eps_tek < eps_pol) or (m > n);
{вывод полученных результатов}
IF (UpCase(Pr) = 'Y')
{обращение к ПП VIVOD2}
THEN VIVOD2(Xpol[n+1], F_up, eps_tek, Ypol[n+1], n, (m-1))
{обращение к ПП VIVOD1}
ELSE VIVOD1(Xpol[n+1], F_up, eps_tek, n, (m-1));
ReadKey;
END. {Конец основного блока}