Текст подпрограммы и версий
is03r_c.zip , is03d_c.zip 		Тексты тестовых примеров
tis03r_c.zip , tis03d_c.zip

Подпрограмма:  is03r_c

Назначение

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

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

Пусть в узлах  xi :  x1 < x2 < ...< xl , заданы значения табличной функции  gi . Строится сглаживающая сплайн - функция 

                 n
    f(x) =   ∑   ai Ni k (x) ,  определяемая условиями :
               i=1
      l
     ∑  wi ( gi - f( xi ) )2  -  min a  ;
    i=1

здесь  wi > 0 - заданные весовые коэффициенты, Ni k,  i = 1, ..., n - нормированные В - сплайны k - го порядка, соответствующие узлам  ti, ti + 1, ..., ti + k таким, что :  t1 ≤ t2 ≤ ...≤ tk < tk + 1 < ...< tn < tn + 1 ≤ ... ≤ tn + k. При этом требуется, чтобы  tk ≤x1 < xl ≤ tn + 1,  а  l ≥ n.

C. de Boor, Package for Calculating with B - splines, SIAM J. Numerical Analysis, 14(3), 1977, pp. 441-472.

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

    int is03r_c (integer *n, integer *k, integer *lx, real *t,
            real *x, real *g, real *w, real *a, real *r, real *r1, real *r2, 
            integer *ierr)

   Параметры
n - заданное число нормированных B - сплайнов (тип: целый);
k - заданный порядок B - сплайна (тип: целый);
lx - заданное число узлов аппроксимации, lx ≥ n (тип: целый);
t - вещественный вектоp длины n + k заданных значений узлов сплайна: t (1) ≤ t (2) ≤ ...≤ t (k) < t (k + 1) < ...< t (n + 1) ≤ t (n + 2) ≤ ...≤ t (n + k);
x - вещественный вектоp длины lx заданных значений узлов аппроксимации:  t (k) ≤ x (1) < x (2) < ... < x (lx) ≤ t (n + 1);
g - вещественный вектоp длины lx заданных значений сглаживаемой функции, g (i) = gi ,  i = 1, 2, ..., lx;
w - вещественный вектоp длины lx весовых коэффициентов w (i) = wi > 0 ,  i = 1, 2, ..., lx;
a - вещественный вектоp длины n, коэффициентов сглаживающего сплайна,  ai = a (i) ,  i = 1, 2, ..., n;
r - вещественный двумерный рабочий массив размера n на (2*k - 1);
r1 - вещественный рабочий вектоp длины k;
r2 - вещественный двумерный рабочий массив размера n на k;
ierr - целая переменная, служащая для сообщения об ошибках, обнаруженных в ходе работы подпрограммы; при этом:
ierr=65 - если матрица нормальной системы уpавнений для определения коэффициентов сглаживающего сплайна вырождена. 
   Версии
is03d_c - среднеквадратическое сглаживание дискретно заданной функции сплайном k - го порядка с повышенной точностью. Массивы t, x, g, w, a, r, r1, r2 имеют тип double. 
   Вызываемые подпрограммы
utis10_c - подпрограмма выдачи диагностических сообщений при работе подпрограммы is03r_c.
utis11_c - подпрограмма выдачи диагностических сообщений при работе подпрограммы is03d_c.
asb1r_c - решение системы линейных алгебраических уpавнений с ленточной матрицей, заданной в компактной форме, с выбором ведущего элемента по столбцу.
asb1d_c - решение системы линейных алгебраических уpавнений с ленточной матрицей, заданной в компактной форме с удвоенной точностью (выбор ведущего элемента по столбцу). 
   Замечания по использованию
  В подпрограммах is03r_c, is03d_c используются служебные подпрограммы i i21r1_c, i i21d1_c.

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

int main(void)
{
    /* Initialized data */
    static float t[9] = { 0.f,0.f,0.f,2.f,4.f,6.f,8.f,8.f,8.f };
    static float w[12] = { 1.f,1.f,1.f,1.f,1.f,1.f,1.f,1.f,1.f,1.f,1.f,1.f };

    /* Local variables */
    extern int is03r_c(int *, int *, int *, float *, float *, float *,
                       float *, float *, float *, float *, float *, int *);
    static int ierr;
    static float a[6], g[12];
    static int i__, k, n;
    static float r__[30] /* was [6][5] */, x[12], r1[3],
                  r2[18] /* was [6][3] */;
    static int lx;
    static float xx;
    int i__1, i;

    n = 6;
    lx = 12;
    k = 3;
    i__1 = lx;
    for (i__ = 1; i__ <= i__1; ++i__) {
        xx = (float) i__ * .5f + .5f;
        x[i__ - 1] = xx;
        g[i__ - 1] = xx * xx * xx;
/* l5: */
    }
    is03r_c(&n, &k, &lx, t, x, g, w, a, r__, r1, r2, &ierr);

    printf("\n %5i \n",ierr);
    for (i = 0; i <= 3; i += 3) {
        printf("\n %16.7e %16.7e %16.7e \n",a[i], a[i+1], a[i+2]);
    }
    return 0;
} /* main */


Результаты:

       ierr  =  0
       
       a  =  ( .290,  -2.236,  18.081,  110.,  321.72,  520.93 )