function [ddT,dT,T,x] = Pseudo(N,a,b)
//N is the degree of the highest order polynomial used
//[a b] is the computational domain with a and b both being finite
//ddT_mn=(T''_n(x_m)),  dT_mn=(T'_n(x_m)),  T_mn=(T_n(x_m))
//x=-cos(pi*(0:N)'/N)


// Output variables initialisation (not found in input variables)
ddT=[];
dT=[];
T=[];
x=[];


// Display mode
mode(0);

// Display warning for floating point exception
ieee(1);

n = mtlb_imp(0,N);x = -cos((%pi*n')/N);xsub = x(mtlb_imp(2,N));T = cos(acos(x)*n);Usub = diag(1 ./sqrt(1-xsub .^2))*sin(acos(xsub)*n);
dT = (2/mtlb_s(b,a))*[-((-1) .^n) .*(n .^2);Usub*diag(n);n .^2];
ddT = (4/(mtlb_s(b,a)^2))*[((((-1) .^n) .*(n .^2)) .*(n .^2-1))/3;(diag(1 ./(1-xsub .^2))*mtlb_s(diag(xsub)*Usub,T(mtlb_imp(2,N),:)*diag(n)))*diag(n);((n .^2) .*(n .^2-1))/3];
x = mtlb_a((mtlb_s(b,a)/2)*x,mtlb_a(a,b)/2);
endfunction
