intlag_cgl

 INTLAG_CGL  Computes matrix "a" to evaluate 1D Lagrange interpolant at CGL 

    (Chebyshev Gauss Lobatto)  nodes in [-1,1] at another mesh x_new in [-1,1]
             (formula (2.4.30)  pag. 88 CHQZ2)

  [a]=intlag_cgl(x_cgl, x_new) 

 Input: x_cgl = array of np_cgl CGL nodes in [-1,1] (ordered from left to
                   right)
        x_new = array of np_new another set in [-1,1] (ordered from left to
                   right)

 Output: a = matrix of size (np_new,np_cgl) 


       If u_cgl is the array with evaluations  of a function u at nodes
          x_cgl, it holds:  u_new = a * u_cgl, i.e. u_new=u(x_new)

 Reference: CHQZ2 = C. Canuto, M.Y. Hussaini, A. Quarteroni, T.A. Zang,
                    "Spectral Methods. Fundamentals in Single Domains"
                    Springer Verlag, Berlin Heidelberg New York, 2006.

0001 function [a]=intlag_cgl(x_cgl, x_new);
0002 % INTLAG_CGL  Computes matrix "a" to evaluate 1D Lagrange interpolant at CGL
0003 %
0004 %    (Chebyshev Gauss Lobatto)  nodes in [-1,1] at another mesh x_new in [-1,1]
0005 %             (formula (2.4.30)  pag. 88 CHQZ2)
0006 %
0007 %  [a]=intlag_cgl(x_cgl, x_new)
0008 %
0009 % Input: x_cgl = array of np_cgl CGL nodes in [-1,1] (ordered from left to
0010 %                   right)
0011 %        x_new = array of np_new another set in [-1,1] (ordered from left to
0012 %                   right)
0013 %
0014 % Output: a = matrix of size (np_new,np_cgl)
0015 %
0016 %
0017 %       If u_cgl is the array with evaluations  of a function u at nodes
0018 %          x_cgl, it holds:  u_new = a * u_cgl, i.e. u_new=u(x_new)
0019 %
0020 % Reference: CHQZ2 = C. Canuto, M.Y. Hussaini, A. Quarteroni, T.A. Zang,
0021 %                    "Spectral Methods. Fundamentals in Single Domains"
0022 %                    Springer Verlag, Berlin Heidelberg New York, 2006.
0023 
0024 %   Written by Paola Gervasio
0025 %   $Date: 2007/04/01$
0026 
0027 x_cgl=-x_cgl;
0028 np_cgl=length(x_cgl);
0029 np_new=length(x_new);
0030 n=np_cgl-1;
0031 a=zeros(np_new,np_cgl);
0032 
0033 [pn,pn1]=jacobi_eval(x_new,n,-0.5,-0.5);
0034 pn1(:,1)=pn1(:,1)*2^(2*n)*(prod(1:n))^2/(prod(1:2*n));
0035 c=[2;ones(n-1,1);2];
0036 
0037 for i=1:np_new
0038 for l=1:np_cgl
0039 if abs(x_cgl(l)-x_new(i))>1.e-14
0040 a(i,np_cgl+1-l)=(-1)^(l)*(1-x_new(i)^2)*pn1(i,1)/(c(l)*n*n*(x_new(i)-x_cgl(l)));
0041 else
0042 a(i,np_cgl+1-l)=1;   
0043 end
0044 end
0045 end
0046 
0047 return