logo search
Лекция_3(Интерполяция)

3.2.12. Кубические сплайны

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

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

Например, непрерывная кусочно-линейная функция (ломаная) является сплайном первой степени с дефектом, равным единице, так непрерывна только сама функция (нулевая производная), а первая производная уже имеет разрыв.

На практике наиболее широкое распространение получили сплайны третьей степени, имеющие на непрерывную, по крайней мере, первую производную.

Величина называется наклоном сплайна в точке (узле) .

Кубический сплайн, заданный локально, – это интерполирующая функция в виде полинома третьей степени, вычисляемая по формуле

где и - первые производные :

; .

Производные локального сплайна могут задаваться тремя способами.

Способ 1. Производные и вычисляются с помощью формул численного дифференцирования по трём точкам:

Способ удобен тем, что для задания сплайна требуется вводить лишь ординаты (значения вычисляются программой). Для уменьшения времени многократных вычислений желательно предварительно вычислить массив и хранить его в памяти.

Способ 2. Значения (вычисленные отдельно или полученные из графиков как наклоны его в узлах) задаются непосредственно в виде массива .

Способы 1 и 2 называются локальными, поскольку с их помощью на каждом частичном отрезке сплайн строится отдельно (непосредственно по формуле). При этом тем не менее соблюдается непрерывность в узлах производной . Непрерывность же второй производной в узлах не гарантируется. Поэтому дефект такого сплайна равен двум. Это эрмитовы кубические сплайны, рассмотренные ранее.

Способ 3. (глобальный способ). Другой подход у выбору заключается в том, чтобы заставить функцию удовлетворять дополнительным условиям. Разумно требовать, чтобы эта функция была более гладкой, т. е. чтобы отдельные составляющие её полиномы при стыковке в узлах давали непрерывную вторую производную. Найдём вторые производные слева и справа от узла на отрезке :

Потребуем непрерывность второй производной в узлах :

при всех

и приходим к следующей системе линейных алгебраических уравнений относительно наклонов:

, . (*)

Так как неизвестных , то нужно задавать ещё два условия, которые называют краевыми (они обычно связаны с «крайними» значениями ).

Дадим три варианта краевых условий.

  1. Если известны , то задаём

.

  1. Производные аппроксимируем формулами численного дифференцирования третьего порядка с отбрасыванием остаточных членов:

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

Краевые условия можно комбинировать, т. е. в левом и правом узлах выбирать их независимо.

Система (*) при всех рассмотренных краевых условиях имеет единственное решение, для нахождения которого могут быть применены методы прогонки итераций.

Решая систему (*) при выбранных краевых условиях, находим наклоны во всех узлах. Затем задаём сплайн на каждом локальном отрезке. Построенный данным глобальным способом сплайн имеет дефект не больше единицы, так как этот сплайн обладает на отрезке непрерывной второй производной.

Пример.

Unit Splin;

{$N+ $E+}

interface

Uses Rtype;

Const Nmax = 50;

Type RealType = Extended;

Vect = Array [1..Nmax] of RealType;

Type SplineInterpolate = object

N1 : Integer;

X, Y, M : Vect;

Constructor Init;

procedure Input; Virtual;

procedure NDotSpline(NNN:Integer);

procedure CalculateSpline;

function Interpolate(XX:RealType) : RealType;

procedure Output; Virtual;

function GetX(I:Integer) : RealType;

function GetY(I:Integer) : RealType;

Destructor Done;

end; (* Object SplineInterpolate *)

implementation

constructor SplineInterpolate.Init;

Var i : Integer;

begin

for i:=1 to Nmax do X[i]:=0;

for i:=1 to Nmax do Y[i]:=0;

for i:=1 to Nmax do M[i]:=0;

end;

Procedure SplineInterpolate.Input;

Var i : Integer;

begin

WriteLn(' Кубическая сплайн интерполяция-экстраполяция...');

WriteLn(' ===============================================');

WriteLn;

WriteLn('Введите число точек интерполяции N : ');Read(N1);

for i:= 1 to N1 do

begin

WriteLn(' X[',I:2, '], Y[', I:2, '] = '); Read(X[i], Y[i]);

end;

end; (* SplineInterpolate.Input *)

Procedure SplineInterpolate.NDotSpline(NNN:Integer);

begin

N1:=NNN;

end; (* SplineInterpolate.NDotSpline *)

procedure SplineInterpolate.CalculateSpline;