Построение и работа с PPFORM

Этот пример показывает, как создать и работать с ppform сплайна в Curve Fitting Toolbox™.

Введение

(Одномерный) кусочный полином или стр для краткости характеризуется его последовательностью пропуска, breaks говорят, и его массив коэффициентов, coefs говорит локальной формы степени ее полиномиальных частей. Последовательность пропуска принята, чтобы строго увеличиться,

  breaks(1) < breaks(2) < ... < breaks(l+1),

с l количество полиномиальных частей, которые составляют стр. В фигуре ниже, breaks [0,1,4,6], следовательно l равняется 3.

В то время как эти полиномы могут иметь различные степени, они все зарегистрированы как полиномы того же порядка k, т.е. массив коэффициентов, coefs имеет размер [l,k] с coefs(j,:), содержащим коэффициенты k в локальной форме степени для j-th полиномиальная часть.

Вот пример, стр составили из трех квадратичных полиномов, т.е. l = k = 3. Пропуски отмечены красным.

sp = spmak([0 1 4 4 6],[2 -1]);
pp = fn2fm(sp,'pp') ;
breaks = fnbrk(pp,'b');
coefs = fnbrk(pp,'c');
coefs(3,[1 2]) = [0 1];
pp = ppmak(breaks,coefs,1);
fnplt(pp,[breaks(1)-1 breaks(2)],'g',1.8)
hold on
fnplt(pp, breaks([2 3]),'b',1.8)
fnplt(pp,[breaks(3),breaks(4)+1],'m',1.8)
lp1 = length(breaks);
xb = repmat(breaks,3,1);
yb = repmat([2;-2.2;NaN],1,lp1);
plot(xb(:),yb(:),'r')
axis([-1 7 -2.5 2.3])
hold off

Точное описание стр с точки зрения последовательности пропуска breaks и массив коэффициентов coefs

  pp(t) = polyval(coefs(j,:), t-breaks(j))
                              for breaks(j) <= t < breaks(j+1)

где, чтобы вспомнить,

  polyval(a,x) = a(1)*x^(k-1) + a(2)*x^(k-2) + ... + a(k)*x^0.

Для стр в фигуре выше, breaks(1) 0, и coefs(1,:) [-1/2 0 0], в то время как breaks(3) равняется 4, и coefs(3,:) [0 1 - 1]. Для t не в [breaks(1) .. breaks(l+1)), pp(t) задан путем расширения первой или последней полиномиальной части.

Стр обычно создается посредством процесса интерполяции или приближения. Но также возможно составить один в ppform с нуля, с помощью команды ppmak. Например, стр выше могут быть получены как

breaks = [0 1 4 6];
coefs = [1/2 0 0 -1/2 1 1/2 0 1 -1];
fn = ppmak(breaks,coefs)
fn = 

  struct with fields:

      form: 'pp'
    breaks: [0 1 4 6]
     coefs: [3x3 double]
    pieces: 3
     order: 3
       dim: 1

Это возвращается, в fn, полном описании этой функции стр в так называемой ppform.

Различные команды в Curve Fitting Toolbox могут работать с этой формой. Остающиеся разделы показывают некоторые примеры.

Работа на кусочных полиномах

Чтобы оценить стр, используйте команду fnval.

fnval(fn, -1:7)
ans =

  Columns 1 through 7

    0.5000         0    0.5000    1.0000    0.5000   -1.0000         0

  Columns 8 through 9

    1.0000    2.0000

Чтобы дифференцировать стр, используйте команду fnder.

dfn = fnder ( fn );
hold on
fnplt(dfn, 'jumps','y', 3)
hold off
h1 = findobj(gca,'Color','y');
legend(h1,{'First Derivative'},'location','SW')

Обратите внимание на то, что производная стр в качестве примера непрерывна в 1, но прерывиста в 4. Также обратите внимание, что по умолчанию fnplt строит ppform на своем основном интервале, т.е. на интервале [breaks(1) .. breaks(end)].

Можно также использовать fnder, чтобы взять вторую производную стр.

ddfn = fnder(fn, 2);
hold on
fnplt( ddfn ,'j', 'k', 1.6)
hold off
h2 = findobj(gca,'Color','k');
legend([h1 h2],{'First Derivative' 'Second Derivative'},'location','SW')

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

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

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

iddfn = fnint(ddfn);
hold on
fnplt(iddfn, 'b', .5)
hold off
h3 = findobj(gca,'Color','b', 'LineWidth',.5);
legend([h1 h2 h3],{'First Derivative' 'Second Derivative' ...
                   'Integral of Second Derivative'},'location','SW')

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

Можно получить части при помощи команды fnbrk. Например,

breaks = fnbrk(fn, 'breaks')
breaks =

     0     1     4     6

восстанавливает последовательность пропуска стр в fn, в то время как

piece2 = fnbrk(fn, 2);

восстанавливает вторую полиномиальную часть, как этот график подтверждает.

fnplt(pp,[breaks(1)-1 breaks(2)],'g',1.8)
hold on
fnplt(piece2, 'b', 2.5, breaks([2 3])+[-1 .5])
fnplt(pp,[breaks(3),breaks(4)+1],'m',1.8)
plot(xb(:),yb(:),'r')
title('The Polynomial that Supplies the Second Polynomial Piece')
hold off
axis([-1 7 -2.5 2.3])

Кусочные полиномы с векторным знаком

Стр могут также быть с векторным знаком, чтобы описать кривую, в с 2 пробелами или с 3 пробелами. В этом случае каждый локальный полиномиальный коэффициент является вектором, а не номером, но ничем иным об изменениях ppform. Существует одна дополнительная часть ppform, чтобы записать это, размерность ее цели.

Например, вот оцененные стр 2 векторов, описывающие модульный квадрат, когда его график показывает. Это - 2D кривая, следовательно ее размерность равняется 2.

square = ppmak(0:4, [1 0  0 1  -1 1  0 0 ; 0 0  1 0  0 1  -1 1]);
fnplt(square,'r',2)
axis([-.5 1.5 -.5 1.5])
axis equal
title('A Vector-Valued PP that Describes a Square')

Многомерные кусочные полиномы

Стр в Curve Fitting Toolbox могут также быть многомерными, а именно, продукт тензора одномерных функций стр. ppform таких многомерных стр незначительно более сложна с breaks теперь массив ячеек, содержащий последовательность пропуска для каждой переменной и coefs теперь многомерный массив. Намного более трудно составить неслучайное такая функция с нуля, таким образом, мы не попробуем это здесь, особенно поскольку тулбокс предназначается, чтобы помочь с конструкцией функций стр от некоторых условий о них. Например, сфера в этой фигуре создается при помощи csape.

x = 0:4;
y = -2:2;
s2 = 1/sqrt(2);
v = zeros(3,7,5);
v(3,:,:) = [0 1 s2 0 -s2 -1 0].'*[1 1 1 1 1];
v(2,:,:) = [1 0 s2 1 s2 0 -1].'*[0 1 0 -1 0];
v(1,:,:) = [1 0 s2 1 s2 0 -1].'*[1 0 -1 0 1];
sph = csape({x,y},v,{'clamped','periodic'});
fnplt(sph)
axis equal
axis off
title('A Sphere Described by a Bicubic 3-Vector-Valued Spline')

В то время как ppform сплайна эффективна для оценки, конструкция сплайна от некоторых данных обычно обрабатывается более эффективно первым определением его B-формы, т.е. его представления как линейная комбинация B-сплайнов. Для получения дополнительной информации смотрите Построение и работу с B-формой.

Для просмотра документации необходимо авторизоваться на сайте