pnleg1

 PNLEG1    Evaluates the first derivative of Legendre polynomial
           of degree n 

    [p1,p]=pnleg1(x,n)  evaluates (L_n)'(x), L_n(x)
    at the node(s) x, using the three term relation  (2.3.3), pag. 75 CHQZ2.

 Input: x = scalar or column array
        n = degree of Legendre polynomial

 Output: p1 = scalar or column array (same dimension as x)
             with the evaluation of (L_n)'(x)
         p  = scalar or column array (same dimension as x)
             with the evaluation of L_n(x)

 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 [p1,p] = pnleg1 (x, n) 
0002 % PNLEG1    Evaluates the first derivative of Legendre polynomial
0003 %           of degree n
0004 %
0005 %    [p1,p]=pnleg1(x,n)  evaluates (L_n)'(x), L_n(x)
0006 %    at the node(s) x, using the three term relation  (2.3.3), pag. 75 CHQZ2.
0007 %
0008 % Input: x = scalar or column array
0009 %        n = degree of Legendre polynomial
0010 %
0011 % Output: p1 = scalar or column array (same dimension as x)
0012 %             with the evaluation of (L_n)'(x)
0013 %         p  = scalar or column array (same dimension as x)
0014 %             with the evaluation of L_n(x)
0015 %
0016 % Reference: CHQZ2 = C. Canuto, M.Y. Hussaini, A. Quarteroni, T.A. Zang,
0017 %                    "Spectral Methods. Fundamentals in Single Domains"
0018 %                    Springer Verlag, Berlin Heidelberg New York, 2006.
0019 
0020 %   Written by Paola Gervasio
0021 %   $Date: 2007/04/01$
0022 
0023 
0024 nn=size(x);
0025 
0026 if nn==1
0027     % x is a scalar
0028     p=0; p1=0;
0029 if n==0
0030 p=1;p1=0;
0031 elseif n==1
0032 p=x;p1=1;
0033 else
0034 p1=0;p2=0;p11=0;p21=0;p3=0;p31=0;
0035 p1=1;  p2=x; 
0036 p11=0; p21=1;
0037 for k =1:n-1
0038    duekp1=2*k+1;kp1=1/(k+1);
0039    p3=(duekp1*x*p2-k*p1)*kp1; 
0040    p31=(duekp1*(x*p21+p2)-k*p11)*kp1; 
0041    p11=p21; p21=p31; 
0042    p1=p2; p2=p3; 
0043 end 
0044 p1=p31; p=p3;
0045 end
0046 
0047 else
0048 x=x(:);
0049     % x is an array
0050 p=zeros(nn); p1=p;
0051 if n==0
0052 p=ones(nn);p1=zeros(nn);
0053 elseif n==1
0054 p=x;p1=ones(nn);
0055 else
0056 p1=p;p2=p;p11=p;p21=p;p3=p;p31=p;
0057 p1=1;  p2=x; 
0058 p11=0; p21=1;
0059 for k =1:n-1
0060    duekp1=2*k+1;kp1=1/(k+1);
0061    p3=(duekp1*x.*p2-k*p1)*kp1; 
0062    p31=(duekp1*(x.*p21+p2)-k*p11)*kp1; 
0063    p11=p21; p21=p31; 
0064    p1=p2; p2=p3; 
0065 end 
0066 p1=p31; p=p3;
0067 end
0068 end
0069 return