% 
% 	All subroutines in here are for [0,1]
% 

N = 10
h = pi/sqrt(N)


function result = ln(x)
  result = log(x);
end

function result = sinc(x)
  result = x;% init vector
  small = find(abs(x)<0.004);
  large = find(abs(x)>=0.0039);
  result(small) = 1-1/6.*pi.^2*x(small).^2+1/120.*pi.^4.*x(small).^4;
  result(large) = sin(pi.*x(large))./(pi.*x(large));
end


function result = sinc_p(x)
  result = x;% init vector
  small = find(abs(x)<0.004);
  large = find(abs(x)>=0.0039);
  result(small)=-1./3.*pi.^2.*x(small)+1./30.*pi.^4.*x(small).^3-1./840.*pi.^6.*x(small).^5;
  result(large)=cos(pi.*x(large))./x(large)-sin(pi.*x(large))./pi./x(large).^2;
end

function result = sinc_pp(x)
  result = x;% init vector
  small = find(abs(x)<0.004);
  large = find(abs(x)>=0.0039);
  result(small)=-1./3.*pi.^2+1./10.*pi.^4.*x(small).^2-1./168.*pi.^6.*x(small).^4;
  result(large)=-sin(pi.*x(large)).*pi./x(large)-2.*cos(pi.*x(large))./x(large).^2+2.*sin(pi.*x(large))./pi./x(large).^3;
end


%
% 	xk = psi(k*h), called as xk(k, h) 
% 
function result = xk(k, h)
  result = exp(h.*k)./(1+exp(h.*k));
end

% 
% 	The S(k,1) o phi(x) * g(x) function, with standard phi and g
% 	
function result = skh(k,h,x)
  result=sinc((ln(x./(1-x))-k.*h)./h).*x.*(1-x);
end

%
% 	The first derivative of the S(k,1) o phi(x) * g(x) function,
% 	with standard phi and g
%
function result = skh_p(k,h,x)
  result = \
  sinc_p(-(-ln(-x./(-1+x))+k.*h)./h)./h-sinc(-(-ln(-x./(-1+x))+k.*h)./h).*(-1+2.*x);
end

% 
% 	The second derivative of the S(k,1) o phi(x) * g(x) function,
% 	with standard phi and g
% 
function result = skh_pp(k,h,x)
  result = \
  (-1+2.*x).*sinc_p(-(-ln(-x./(-1+x))+k.*h)./h)./(-1+x)./x./h-sinc_pp(-(-ln(-x./(-1+x))+k.*h)./h)./h.^2./(-1+x)./x-2.*sinc(-(-ln(-x./(-1+x))+k.*h)./h);
end


%============================================================================= 
% 	The spline (0,0,1,0,0,1):
%============================================================================= 
% 
function result = rig_val(vec)
  result = 3*vec.^2-2*vec.^3;
end
%	and its first derivative
function result = rig_val_p(x)
  result = 6*x-6*x.^2;
end
%	and its second derivative
function result = rig_val_pp(x)
  result = 6-12*x;
end



%============================================================================= 
% 	The spline (0,1,0,0,0,1):
%============================================================================= 
% 	
function result = lef_der(x)
  result = x-2.*x.^2+x.^3;
end

%	and its first derivative
function result = lef_der_p(x)
  result = 1-4.*x+3.*x.^2;
end

%	and its second derivative
function result = lef_der_pp(x)
  result = -4+6.*x;
end




%============================================================================= 
% 	The spline (0,0,0,1,0,1):
%============================================================================= 
% 
function result = rig_der(vec)
  result = -vec.^2+vec.^3;
end
