code: def _C(xi, N,L, alpha = lambda x: 2*x + 3./2):
      
      CC = nu.zeros((N,L), float) 
       
      for l in range(L):
          for n in range(N):
              a = alpha(l)
              if n==0:
                  CC[n][l] = 1.
                  continue 
              elif n==1: CC[n][l] = 2.*a*xi
              if n + 1 != N:
                  CC[n+1][l] = (n + 1.)**-1. * (2*(n + a)*xi*CC[n][l] - (n + 2*a - 1)*CC[n-1][l])
      return CC
docs: NAME:
         _C
      PURPOSE:
         Evaluate C_n,l (the Gegenbauer polynomial) for 0 <= l < L and 0<= n < N 
      INPUT:
         xi - radial transformed variable
         N - Size of the N dimension
         L - Size of the L dimension
         alpha = A lambda function of l. Default alpha = 2l + 3/2 
         
      OUTPUT:
         An LxN Gegenbauer Polynomial 
      HISTORY:
         2016-05-16 - Written - Aladdin