MatLab / Hermit

function [a,p] = Hermit(X,N)
%X - вектор узлов
%N - вектор кратностей узлов
%p - степень многочлена Эрмита

k = length(X);
n = sum(N);
p = n - 1;
F = fun(X);
m = max(N)-1;

%оптимальный шаг
x = X(1)-10^-3:10^-3:X(k)+10^-3;
f = fun(x);
d2f = diff(f,m)/(10^-3)^m;
h = m*sqrt(2^(-52)/max(abs(d2f)));

df = zeros(1,n-k);
A1 = ones(k,n);
A2 = zeros(n-k,n);
for i = 2:1:n
A1(:,i) = X'.^(i-1);
end
index = 1;
str = 1;
for i = 1:1:k
if (N(i) ~= 1)
for j = 2:1:N(i)
df(index) = DF_X(X(i),h,j-1);
index = index + 1;
st = 0;
for t = j:1:n
if (j == 2)
A2(str,t) = (st+1)*X(i)^st;
st = st + 1;
else
A2(str,t) = (st+1)*A2(str-1,t)/X(i);
st = st + 1;
end
end
str = str + 1;
end
end
end
A = [A1; A2];
H1 = [F'; df'];
a = A\H1;
end