Текст подпрограммы и версий
ia80r_p.zip 		Тексты тестовых примеров
tia80r_p.zip

Подпрограмма:  IA80R (модуль IA80R_p)

Назначение

Рациональная Чебышевская аппроксимация вещественной функции на конечном интервале.

Математическое описание

IA80R вычисляет лучшую в смысле Чебышева рациональную аппроксимацию с некоторым весом вещественной непрерывной на отрезке  [a, b] функции  F (X) : 

 R*[PHI(X)]  =  ( P1 + P2* PHI(X) + ... + PN+1* PHI(X)N ) /
                       / ( Q1 + Q2* PHI(X) + ... + QM+1* PHI(X)M ) ,

причем  R* [ PHI (X)] обладает тем свойством, что на ней достигается минимум 

     max  | R[PHI(X)] - F(X) |  /  | G(X) | .
    [A, B]

Минимизация проводится на классе рациональных функций  R [ РНI (Х)], имеющих степенью числителя  N, а степенью знаменателя  M.

W.J.Cody, W.Fraser, J.F.Hart, Rational Chebyshev Approximation Using Linear Equations, Numerishe Mathematik, 12 (4), 1968.

Использование 

procedure IA80R(F :Func_F1; PHI :Func_F1; G :Func_F1; A :Real; B :Real;
                EPS :Real; N :Integer; M :Integer;
                var P :Array of Real; var Q :Array of Real;
                var RAB :Array of Real; var IERR :Integer);


   Параметры
F - имя вещественной подпрограммы - функции вычисления функции  F (X);
PHI - имя вещественной подпрограммы - функции вычисления функции  PHI (X);
G - имя вещественной подпрограммы - функции вычисления функции  G (X);
A, B - заданные нижняя и верхняя границы отрезка, на котоpом производится аппроксимация (тип: вещественный);
EPS - заданная точность вычисления аппроксимации (тип: вещественный);
N, M - заданные степень числителя и знаменателя (тип: целый);
P - вещественный вектоp длины  N + 1, в который помещаются коэффициенты числителя;
Q - вещественный вектоp длины  M + 1, в который помещаются коэффициенты знаменателя;
RAB - вещественный вектоp длины  (N + M + 8) * (N + M + 2) , используемый в подпрограмме как рабочий; в  RAB (1) помещается вычисленное значение наилучшего приближения;
IERR - целая переменная, служащая для сообщения об ошибках, обнаруженных в ходе работы подпрограммы; при этом:
IERR=65 - когда заданная точность не может быть достигнута;
IERR=66 - когда алгоритм не сходится; возможно, что знаменатель имеет корень на заданном отрезке;
IERR=67 - когда сходимость не может быть достигнута в пределах 20 итераций;
IERR=68 - когда алгоритм не сходится из - за сингулярности матрицы системы уравнений для определения  P и  Q;
IERR=69 - когда получена хорошая аппроксимация, но не обязательно наилучшая в смысле Чебышева. 
   Версии:  нет 
   Вызываемые подпрограммы
UTIA10 - подпрограмма выдачи диагностических сообщений при работе подпрограммы IA80R. 
   Замечания по использованию
 

Если  IERR = 65, то числитель и знаменатель имеют общие множители. B этом случае значения  N и  M следует уменьшить на 1.

Bо всех других случаях, когда  IERR ≠ 0, следует попытаться или уменьшить степени числителя и знаменателя, или изменить пределы отрезка, или уменьшить точность вычислений.

Пример использования

Пусть требуется аппроксимировать функцию  √ x  на  [0.5, 2.]  следующим образом: 

    x / 4.  +  ( P1 + P2* x ) / ( Q1 + Q2* x ) .


Unit TIA80R_p;
interface
uses
SysUtils, Math, { Delphi }
Lstruct, Lfunc, UtRes_p, FFIA80R_p, FPIA80R_p, FGIA80R_p, IA80R_p;

function TIA80R: String;

implementation

function TIA80R: String;
var
N,M,_i,IERR :Integer;
A,B,EPS :Real;
P :Array [0..1] of Real;
Q :Array [0..1] of Real;
RАВ :Array [0..39] of Real;
begin
Result := '';  { результат функции }
for _i:=0 to 39 do  //обнуление рабочего массива
 RAB[_i] := 0.0;
A := 0.5;
B := 2.0;
EPS := 1.E-2;
N := 1;
M := 1;
IA80R(FFIA80R,FPIA80R,FGIA80R,A,B,EPS,N,M,P,Q,RAB,IERR);
Result := Result + Format(' %20.16f %20.16f %20.16f  %5d %5d ',
 [A,B,EPS,N,M]) + #$0D#$0A;
Result := Result + #$0D#$0A;
for _i:=0 to 1 do
 begin
  Result := Result + Format('%20.16f ',[P[_i]]);
  if ( ((_i+1) mod 2)=0 )
   then Result := Result + #$0D#$0A;
 end;
Result := Result + #$0D#$0A;
Result := Result + #$0D#$0A;
for _i:=0 to 1 do
 begin
  Result := Result + Format('%20.16f ',[Q[_i]]);
  if ( ((_i+1) mod 2)=0 )
   then Result := Result + #$0D#$0A;
 end;
Result := Result + #$0D#$0A;
Result := Result + Format(' %20.16f  %20.16f  %5d ',
 [RAB[0],RAB[1],IERR]) + #$0D#$0A;
UtRes('TIA80R',Result);  { вывод результатов в файл TIA80R.res }
exit;
end;

end.

Результаты:

       RAB (1)  =  -0.000319 ;    IERR  =  0 ;

       P (1)  =  0.232604 ,    P (2)  =  1.32410 ;

       Q (1)  =  1. ,    Q (2)  =  1.074789