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

Подпрограмма:  IIS4R (модуль IIS4R_p)

Назначение

Вычисление значений интерполяционного кубического сплайна в заданных точках.

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

Вычисляются значения кубического сплайна  S (x), соответствующего сетке  {xi} и интерполирующего заданные значения табличной функции  yi = y (xi),  i = 1, 2, ..., N, в точках  uk,  k = 1, 2, ..., M по формуле: 

   S (uk) = yj + c1 j ( uk - xj ) + c2 j ( uk - xj )2 + c3 j ( uk - xj )3  ,

где индекс  j такой, что  uk содержится в отрезке  [xj, xj + 1], а  c1 j, c2 j, c3 j,  j = 1, 2, ..., N - 1, - заданные коэффициенты кубического сплайна.

Дж.Алберг, Э.Нильсон, Дж.Уолш, Теория сплайнов и ее приложения, M., Мир, 1972.

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

procedure IIS4R(var X :Array of Real; var Y :Array of Real;
                NX :Integer; var C :Array of Real;
                var U :Array of Real; var S :Array of Real;
                M :Integer; var IERR :Integer);


   Параметры
X - вещественный вектоp длины  NX, содержащий заданные значения узлов сетки и упорядоченный так, что  X (1) < X (2) < ... < X (NX);
Y - вещественный вектоp длины  NX, содержащий заданные в узлах сетки значения интерполируемой функции;
NX - заданное число узлов сетки,  NX ≥ 2 (тип: целый);
C - вещественный двумерный массив размера  3 * (NX - 1), содержащий заданные коэффициенты интерполяционного кубического сплайна;
U - вещественный вектоp длины  M, содержащий заданные точки, в которых вычисляются значения интерполяционного кубического сплайна,  X (1) ≤ U (K) ≤ X (NX),   K = 1, 2, ..., M;
S - вещественный вектоp длины  M, содержащий вычисленные значения интерполяционного кубического сплайна;
M - заданное число точек, в которых вычисляются значения интерполяционного кубического сплайна (тип: целый);
IERR - целая переменная, служащая для сообщения об ошибках, обнаруженных в ходе работы подпрограммы; при этом:
IERR= 1 - когда одна или несколько заданных точек меньше  X (1);
IERR= 2 - когда одна или несколько заданных точек больше  X (NX);
IERR=65 - когда заданное число узлов сетки меньше 2;
IERR=66 - когда элементы вектоpа  X не упорядочены по возрастанию. 
   Версии:  нет 
   Вызываемые подпрограммы
UTI I10 - подпрограмма выдачи диагностических сообщений при работе подпрограммы I IS4R. 
   Замечания по использованию:  нет

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

Unit tiis4r_p;
interface
uses
SysUtils, Math, { Delphi }
Lstruct, Lfunc, UtRes_p, IIS3R_p, IIS4R_p;

function tiis4r: String;

implementation

function tiis4r: String;
var
NX,I,M,_i,IERR :Integer;
Y :Array [0..4] of Real;
BPAR :Array [0..3] of Real;
C :Array [0..11] of Real;
S :Array [0..7] of Real;
const
X :Array [0..4] of Real = ( 1.0,2.0,3.0,4.0,5.0 );
U :Array [0..7] of Real = ( 1.5,4.5,3.0,4.9,2.1,4.0,2.7,1.8 );
label
_5;
begin
Result := '';  { результат функции }
NX := 5;
for I:=1 to NX do
 begin
  Y[I-1] := ((X[I-1]-4.0)*X[I-1]-1.0)*X[I-1]+5.0;
_5:
 end;
BPAR[0] := 1.0;
BPAR[1] := 6.0/(X[1]-X[0])*((Y[1]-Y[0])/(X[1]-X[0])+6.0);
BPAR[2] := 1.0;
BPAR[3] := 6.0/(X[4]-X[3])*(34.0-(Y[4]-Y[3])/(X[4]-X[3]));
IIS3R(X,Y,NX,BPAR,C,IERR);
M := 8;
IIS4R(X,Y,NX,C,U,S,M,IERR);
Result := Result + Format(' %3d ',[IERR]) + #$0D#$0A;
Result := Result + Format('%s',[' S=' + #$0D#$0A]);
Result := Result + #$0D#$0A;
for _i:=0 to 7 do
 begin
  Result := Result + Format('%20.16f ',[S[_i]]);
  if ( ((_i+1) mod 3)=0 )
   then Result := Result + #$0D#$0A;
 end;
Result := Result + #$0D#$0A;
UtRes('tiis4r',Result);  { вывод результатов в файл tiis4r.res }
exit;
end;

end.


Результаты:  

      IERR =0  

      S(1)  =   -2.125         S(5)  =  -5.479
      S(2)  =   10.625         S(6)  =   1.000
      S(3)  =   -7.000         S(7)  =  -7.177
      S(4)  =   21.709         S(8)  =  -3.928