БЧА НИВЦ МГУ. IP04R_P. Полиномиальная интерполяция

Документ взят из кэша поисковой машины. Адрес оригинального документа : http://num-anal.srcc.msu.su/lib_na/cat/i_htm_p/ip04r_p.htm
Дата изменения: Fri Nov 27 12:06:06 2015
Дата индексирования: Sun Apr 10 01:15:59 2016
Кодировка: Windows-1251

Текст подпрограммы и версий
ip04r_p.zip , ip04e_p.zip

Тексты тестовых примеров
tip04r_p.zip , tip04e_p.zip

Подпрограмма: IP04R (модуль IP04R_p)

Назначение

Бикубическая интерполяция функции двух переменных, заданной вместе с частными производными первого порядка и смешанной производной второго порядка в вершинах элементарного сеточного прямоугольника.

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

Пусть X1L, X1U, X2L и X2U следующим образом определяют вершины элементарного сеточного прямоугольника:

        (X1L,X2U)   |                           |    (X1U,X2U)
      ______________________________________
                             |  4                   3  |
                             |       (X1,X2)       |

                             |  1                   2  |
      ______________________________________
        (X1L,X2L)   |                            |    (X1U,X2L)

Эти вершины пронумерованы на рисунке цифрами, начиная с нижней левой вершины против часовой стрелки. Пусть в каждой из вершин известны значения интерполируемой функции

        Y  =  f(x₁, x₂) ,
 обеих ее частных производных первого порядка:
        Y1  =  ∂f / ∂x₁ ,   Y2  =  ∂f / ∂x₂
 и смешанной производной   Y12  =  ∂² f / ∂x₁ ∂x₂ .

Подпрограмма IP04R при помощи бикубической интерполяции вычисляет в точке (X1,X2), заданной внутри исходного элементарного сеточного прямоугольника, значение функции

        AY = f (X1, X2) ,
 обеих ее частных производных первого порядка:
      AY1  =  ∂ f (X1, X2) / ∂x₁ ,   AY2  =  ∂ f (X1, X2) / ∂x₂
 и смешанной производной   AY12  =  ∂² f (X1, X2) / ∂x₁ ∂x₂.

Используются следующие формулы:

                          ₄     ₄                                               
           AY  =    ∑    ∑   c_{i j}  t^i-1 u^j-1 ;
                         ⁱ⁼¹  ^j=1
                          ₄     ₄
          AY1  =   ∑    ∑   (i - 1) c_{i j}  t^i-2 u^j-1 ;
                         ⁱ⁼¹  ^j=1
                          ₄     ₄
          AY2  =   ∑    ∑   (j - 1) c_{i j}  t^i-1 u^j-2 ;
                         ⁱ⁼¹  ^j=1
                          ₄     ₄
        AY12  =   ∑    ∑   (i - 1) (j - 1) c_{i j}  t^i-2 u^j-2 ;
                         ⁱ⁼¹  ^j=1

Здесь t = (X1 - X1L) / (X1U - X1L), U = (X2 - X2L) / (X2U - X2L). Коэффициенты c_i j ( i, j = 1, 2, 3, 4 ) вычисляются в подпрограмме IP04R.

Н.С.Бахвалов. Численные методы. Изд - во "Наука", 1973.

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

procedure IP04R(var Y :Array of Real; var Y1 :Array of Real;
                var Y2 :Array of Real; var Y12 :Array of Real;
                X1L :Real; X1U :Real; X2L :Real; X2U :Real; X1 :Real;
                X2 :Real; var AY :Real; var AY1 :Real; var AY2 :Real;
                var AY12 :Real; var C :Array of Real);

   Параметры

Y, Y1 - Y2, Y12	вещественные векторы длины 4, содержащие соответственно значения интерполируемой функции, обеих ее частных производных первого порядка и смешанной производной в вершинах заданного элементарного сеточного прямоугольника, пронумерованных начиная с нижней левой вершины против часовой стрелки;
X1L - X1U	левая и правая координаты вершин элементарного сеточного прямоугольника по первой переменной (тип: вещественный);
X2L - X2U	нижняя и верхняя координаты вершин элементарного сеточного прямоугольника по второй переменной (тип: вещественный);
X1, X2 -	координаты точки, заданной внутри исходного элементарного сеточного прямоугольника, в которой выполняется интерполяция (тип: вещественный);
AY, AY1 - AY2, AY12	вещественные переменные, значения которых полагаются равными соответственно значениям в точке (X1, X2) интерполируемой функции, ее частных производных первого порядка и смешанной производной;
C -	вещественный двумерный массив размеров 4 на 4, в котором помещаются коэффициенты, используемые при бикубической интерполяции.

   Версии

IP04E -

бикубическая интерполяция функции двух переменных в режиме расширенной (Extended) точности; при этом все формальные параметры должны иметь тип Extended.

   Вызываемые подпрограммы:  нет

   Замечания по использованию:  нет

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

Unit tip04r_p;
interface
uses
SysUtils, Math, { Delphi }
Lstruct, Lfunc, UtRes_p, IP04R_p;

function tip04r: String;

implementation

function tip04r: String;
var
_i :Integer;
X1L,X1U,X2L,X2U,X1,X2,TY,TY1,TY2,TY12,AY1,AY2,AY12,AY :Real;
Y :Array [0..3] of Real;
Y1 :Array [0..3] of Real;
Y2 :Array [0..3] of Real;
Y12 :Array [0..3] of Real;
C :Array [0..15] of Real;
begin
Result := '';  { результат функции }
X1L := 0.0;
X1U := 0.1;
X2L := 0.0;
X2U := 0.15;
Y[0] := Sin(X1L)*Sin(X2L);
Y[1] := Sin(X1U)*Sin(X2L);
Y[2] := Sin(X1U)*Sin(X2U);
Y[3] := Sin(X1L)*Sin(X2U);
Y1[0] := Cos(X1L)*Sin(X2L);
Y1[1] := Cos(X1U)*Sin(X2L);
Y1[2] := Cos(X1U)*Sin(X2U);
Y1[3] := Cos(X1L)*Sin(X2U);
Y2[0] := Sin(X1L)*Cos(X2L);
Y2[1] := Sin(X1U)*Cos(X2L);
Y2[2] := Sin(X1U)*Cos(X2U);
Y2[3] := Sin(X1L)*Cos(X2U);
Y12[0] := Cos(X1L)*Cos(X2L);
Y12[1] := Cos(X1U)*Cos(X2L);
Y12[2] := Cos(X1U)*Cos(X2U);
Y12[3] := Cos(X1L)*Cos(X2U);
X1 := 0.06;
X2 := 0.07;
IP04R(Y,Y1,Y2,Y12,X1L,X1U,X2L,X2U,X1,X2,
     AY,AY1,AY2,AY12,C);
Result := Result + #$0D#$0A;
for _i:=0 to 15 do
 begin
  Result := Result + Format('%20.16f ',[C[_i]]);
  if ( ((_i+1) mod 3)=0 )
   then Result := Result + #$0D#$0A;
 end;
Result := Result + #$0D#$0A;
TY := Sin(X1)*Sin(X2);
TY1 := Cos(X1)*Sin(X2);
TY2 := Sin(X1)*Cos(X2);
TY12 := Cos(X1)*Cos(X2);
Result := Result + Format(' %20.16f  %20.16f  %20.16f  %20.16f  %20.16f  %20.16f  %20.16f %20.16f ',
 [AY,AY1,AY2,AY12,TY,TY1,TY2,TY12]) + #$0D#$0A;
UtRes('tip04r',Result);  { вывод результатов в файл tip04r.res }
exit;
end;

end.

Результаты:
        
           AY = 0.4194o5E-02 
         AY1 = 0.698169E-01 
         AY2 = 0.598171E-01 
       AY12 = 0.995755