[Scipy-svn] r7102 - trunk/scipy/special
scipy-svn at scipy.org
scipy-svn at scipy.org
Sat Jan 29 16:38:13 EST 2011
Author: warren.weckesser
Date: 2011-01-29 15:38:13 -0600 (Sat, 29 Jan 2011)
New Revision: 7102
Modified:
trunk/scipy/special/orthogonal.py
Log:
ENH: special: Don't use plain assert for argument validation.
Modified: trunk/scipy/special/orthogonal.py
===================================================================
--- trunk/scipy/special/orthogonal.py 2011-01-29 20:54:21 UTC (rev 7101)
+++ trunk/scipy/special/orthogonal.py 2011-01-29 21:38:13 UTC (rev 7102)
@@ -102,7 +102,7 @@
'eval_sh_legendre', 'eval_sh_chebyt', 'eval_sh_chebyu',
'eval_sh_jacobi', 'poch', 'binom']
-def poch(z,m):
+def poch(z, m):
"""Pochhammer symbol (z)_m = (z)(z+1)....(z+m-1) = gamma(z+m)/gamma(z)"""
return _gam(z+m) / _gam(z)
@@ -141,7 +141,7 @@
self.__dict__['_eval_func'] = lambda x: evf(x) * p
self.__dict__['normcoef'] *= p
-def gen_roots_and_weights(n,an_func,sqrt_bn_func,mu):
+def gen_roots_and_weights(n, an_func, sqrt_bn_func, mu):
"""[x,w] = gen_roots_and_weights(n,an_func,sqrt_bn_func,mu)
Returns the roots (x) of an nth order orthogonal polynomial,
@@ -168,7 +168,7 @@
return answer
# Jacobi Polynomials 1 P^(alpha,beta)_n(x)
-def j_roots(n,alpha,beta,mu=0):
+def j_roots(n, alpha, beta, mu=0):
"""[x,w] = j_roots(n,alpha,beta)
Returns the roots (x) of the nth order Jacobi polynomial, P^(alpha,beta)_n(x)
@@ -177,7 +177,8 @@
"""
if any(alpha <= -1) or any(beta <= -1):
raise ValueError("alpha and beta must be greater than -1.")
- assert(n>0), "n must be positive."
+ if n < 1:
+ raise ValueError("n must be positive.")
olderr = np.seterr(all='ignore')
try:
@@ -201,12 +202,14 @@
else:
return val
-def jacobi(n,alpha,beta,monic=0):
+def jacobi(n, alpha, beta, monic=0):
"""Returns the nth order Jacobi polynomial, P^(alpha,beta)_n(x)
orthogonal over [-1,1] with weighting function
(1-x)**alpha (1+x)**beta with alpha,beta > -1.
"""
- assert(n>=0), "n must be nonnegative"
+ if n < 0:
+ raise ValueError("n must be nonnegative.")
+
wfunc = lambda x: (1-x)**alpha * (1+x)**beta
if n==0:
return orthopoly1d([],[],1.0,1.0,wfunc,(-1,1),monic,
@@ -222,7 +225,7 @@
return p
# Jacobi Polynomials shifted G_n(p,q,x)
-def js_roots(n,p1,q1,mu=0):
+def js_roots(n, p1, q1, mu=0):
"""[x,w] = js_roots(n,p,q)
Returns the roots (x) of the nth order shifted Jacobi polynomial, G_n(p,q,x),
@@ -232,7 +235,7 @@
# from recurrence relation
if not ( any((p1 - q1) > -1) and any(q1 > 0) ):
raise ValueError("(p - q) > -1 and q > 0 please.")
- if (n <= 0):
+ if n <= 0:
raise ValueError("n must be positive.")
p,q = p1,q1
@@ -266,8 +269,9 @@
orthogonal over [0,1] with weighting function
(1-x)**(p-q) (x)**(q-1) with p>q-1 and q > 0.
"""
- if (n<0):
- raise ValueError("n must be nonnegative")
+ if n < 0:
+ raise ValueError("n must be nonnegative.")
+
wfunc = lambda x: (1.0-x)**(p-q) * (x)**(q-1.)
if n==0:
return orthopoly1d([],[],1.0,1.0,wfunc,(-1,1),monic,
@@ -283,7 +287,7 @@
return pp
# Generalized Laguerre L^(alpha)_n(x)
-def la_roots(n,alpha,mu=0):
+def la_roots(n, alpha, mu=0):
"""[x,w] = la_roots(n,alpha)
Returns the roots (x) of the nth order generalized (associated) Laguerre
@@ -292,7 +296,9 @@
"""
if not all(alpha > -1):
raise ValueError("alpha > -1")
- assert(n>0), "n must be positive."
+ if n < 1:
+ raise ValueError("n must be positive.")
+
(p,q) = (alpha,0.0)
sbn_La = lambda k: -sqrt(k*(k + p)) # from recurrence relation
an_La = lambda k: 2*k + p + 1
@@ -303,14 +309,16 @@
else:
return val
-def genlaguerre(n,alpha,monic=0):
+def genlaguerre(n, alpha, monic=0):
"""Returns the nth order generalized (associated) Laguerre polynomial,
L^(alpha)_n(x), orthogonal over [0,inf) with weighting function
exp(-x) x**alpha with alpha > -1
"""
if any(alpha <= -1):
raise ValueError("alpha must be > -1")
- assert(n>=0), "n must be nonnegative"
+ if n < 0:
+ raise ValueError("n must be nonnegative.")
+
if n==0: n1 = n+1
else: n1 = n
x,w,mu0 = la_roots(n1,alpha,mu=1)
@@ -323,7 +331,7 @@
return p
# Laguerre L_n(x)
-def l_roots(n,mu=0):
+def l_roots(n, mu=0):
"""[x,w] = l_roots(n)
Returns the roots (x) of the nth order Laguerre polynomial, L_n(x),
@@ -332,11 +340,13 @@
"""
return la_roots(n,0.0,mu=mu)
-def laguerre(n,monic=0):
+def laguerre(n, monic=0):
"""Return the nth order Laguerre polynoimal, L_n(x), orthogonal over
[0,inf) with weighting function exp(-x)
"""
- assert(n>=0), "n must be nonnegative"
+ if n < 0:
+ raise ValueError("n must be nonnegative.")
+
if n==0: n1 = n+1
else: n1 = n
x,w,mu0 = l_roots(n1,mu=1)
@@ -349,14 +359,16 @@
# Hermite 1 H_n(x)
-def h_roots(n,mu=0):
+def h_roots(n, mu=0):
"""[x,w] = h_roots(n)
Returns the roots (x) of the nth order Hermite polynomial,
H_n(x), and weights (w) to use in Gaussian Quadrature over
[-inf,inf] with weighting function exp(-x**2).
"""
- assert(n>0), "n must be positive."
+ if n < 1:
+ raise ValueError("n must be positive.")
+
sbn_H = lambda k: sqrt(k/2) # from recurrence relation
an_H = lambda k: 0*k
mu0 = sqrt(pi) # integral of weight over interval
@@ -366,11 +378,13 @@
else:
return val
-def hermite(n,monic=0):
+def hermite(n, monic=0):
"""Return the nth order Hermite polynomial, H_n(x), orthogonal over
(-inf,inf) with weighting function exp(-x**2)
"""
- assert(n>=0), "n must be nonnegative"
+ if n < 0:
+ raise ValueError("n must be nonnegative.")
+
if n==0: n1 = n+1
else: n1 = n
x,w,mu0 = h_roots(n1,mu=1)
@@ -383,14 +397,16 @@
return p
# Hermite 2 He_n(x)
-def he_roots(n,mu=0):
+def he_roots(n, mu=0):
"""[x,w] = he_roots(n)
Returns the roots (x) of the nth order Hermite polynomial,
He_n(x), and weights (w) to use in Gaussian Quadrature over
[-inf,inf] with weighting function exp(-(x/2)**2).
"""
- assert(n>0), "n must be positive."
+ if n < 1:
+ raise ValueError("n must be positive.")
+
sbn_He = lambda k: sqrt(k) # from recurrence relation
an_He = lambda k: 0*k
mu0 = sqrt(2*pi) # integral of weight over interval
@@ -400,11 +416,13 @@
else:
return val
-def hermitenorm(n,monic=0):
+def hermitenorm(n, monic=0):
"""Return the nth order normalized Hermite polynomial, He_n(x), orthogonal
over (-inf,inf) with weighting function exp(-(x/2)**2)
"""
- assert(n>=0), "n must be nonnegative"
+ if n < 0:
+ raise ValueError("n must be nonnegative.")
+
if n==0: n1 = n+1
else: n1 = n
x,w,mu0 = he_roots(n1,mu=1)
@@ -419,7 +437,7 @@
## The remainder of the polynomials can be derived from the ones above.
# Ultraspherical (Gegenbauer) C^(alpha)_n(x)
-def cg_roots(n,alpha,mu=0):
+def cg_roots(n, alpha, mu=0):
"""[x,w] = cg_roots(n,alpha)
Returns the roots (x) of the nth order Ultraspherical (Gegenbauer)
@@ -428,7 +446,7 @@
"""
return j_roots(n,alpha-0.5,alpha-0.5,mu=mu)
-def gegenbauer(n,alpha,monic=0):
+def gegenbauer(n, alpha, monic=0):
"""Return the nth order Gegenbauer (ultraspherical) polynomial,
C^(alpha)_n(x), orthogonal over [-1,1] with weighting function
(1-x**2)**(alpha-1/2) with alpha > -1/2
@@ -443,14 +461,16 @@
# Chebyshev of the first kind: T_n(x) = n! sqrt(pi) / _gam(n+1./2)* P^(-1/2,-1/2)_n(x)
# Computed anew.
-def t_roots(n,mu=0):
+def t_roots(n, mu=0):
"""[x,w] = t_roots(n)
Returns the roots (x) of the nth order Chebyshev (of the first kind)
polynomial, T_n(x), and weights (w) to use in Gaussian Quadrature
over [-1,1] with weighting function (1-x**2)**(-1/2).
"""
- assert(n>0), "n must be positive."
+ if n < 1:
+ raise ValueError("n must be positive.")
+
# from recurrence relation
sbn_J = lambda k: np.where(k==1,sqrt(2)/2.0,0.5)
an_J = lambda k: 0.0*k
@@ -462,11 +482,13 @@
else:
return val
-def chebyt(n,monic=0):
+def chebyt(n, monic=0):
"""Return nth order Chebyshev polynomial of first kind, Tn(x). Orthogonal
over [-1,1] with weight function (1-x**2)**(-1/2).
"""
- assert(n>=0), "n must be nonnegative"
+ if n < 0:
+ raise ValueError("n must be nonnegative.")
+
wfunc = lambda x: 1.0/sqrt(1-x*x)
if n==0:
return orthopoly1d([],[],pi,1.0,wfunc,(-1,1),monic,
@@ -481,7 +503,7 @@
# Chebyshev of the second kind
# U_n(x) = (n+1)! sqrt(pi) / (2*_gam(n+3./2)) * P^(1/2,1/2)_n(x)
-def u_roots(n,mu=0):
+def u_roots(n, mu=0):
"""[x,w] = u_roots(n)
Returns the roots (x) of the nth order Chebyshev (of the second kind)
@@ -490,7 +512,7 @@
"""
return j_roots(n,0.5,0.5,mu=mu)
-def chebyu(n,monic=0):
+def chebyu(n, monic=0):
"""Return nth order Chebyshev polynomial of second kind, Un(x). Orthogonal
over [-1,1] with weight function (1-x**2)**(1/2).
"""
@@ -502,7 +524,7 @@
return base
# Chebyshev of the first kind C_n(x)
-def c_roots(n,mu=0):
+def c_roots(n, mu=0):
"""[x,w] = c_roots(n)
Returns the roots (x) of the nth order Chebyshev (of the first kind)
@@ -516,11 +538,13 @@
[x,w] = j_roots(n,-0.5,-0.5,mu=0)
return [x*2,w]
-def chebyc(n,monic=0):
+def chebyc(n, monic=0):
"""Return nth order Chebyshev polynomial of first kind, Cn(x). Orthogonal
over [-2,2] with weight function (1-(x/2)**2)**(-1/2).
"""
- assert(n>=0), "n must be nonnegative"
+ if n < 0:
+ raise ValueError("n must be nonnegative.")
+
if n==0: n1 = n+1
else: n1 = n
x,w,mu0 = c_roots(n1,mu=1)
@@ -534,7 +558,7 @@
return p
# Chebyshev of the second kind S_n(x)
-def s_roots(n,mu=0):
+def s_roots(n, mu=0):
"""[x,w] = s_roots(n)
Returns the roots (x) of the nth order Chebyshev (of the second kind)
@@ -548,11 +572,13 @@
[x,w] = j_roots(n,0.5,0.5,mu=0)
return [x*2,w]
-def chebys(n,monic=0):
+def chebys(n, monic=0):
"""Return nth order Chebyshev polynomial of second kind, Sn(x). Orthogonal
over [-2,2] with weight function (1-(x/)**2)**(1/2).
"""
- assert(n>=0), "n must be nonnegative"
+ if n < 0:
+ raise ValueError("n must be nonnegative.")
+
if n==0: n1 = n+1
else: n1 = n
x,w,mu0 = s_roots(n1,mu=1)
@@ -567,7 +593,7 @@
return p
# Shifted Chebyshev of the first kind T^*_n(x)
-def ts_roots(n,mu=0):
+def ts_roots(n, mu=0):
"""[x,w] = ts_roots(n)
Returns the roots (x) of the nth order shifted Chebyshev (of the first kind)
@@ -576,7 +602,7 @@
"""
return js_roots(n,0.0,0.5,mu=mu)
-def sh_chebyt(n,monic=0):
+def sh_chebyt(n, monic=0):
"""Return nth order shifted Chebyshev polynomial of first kind, Tn(x).
Orthogonal over [0,1] with weight function (x-x**2)**(-1/2).
"""
@@ -592,7 +618,7 @@
# Shifted Chebyshev of the second kind U^*_n(x)
-def us_roots(n,mu=0):
+def us_roots(n, mu=0):
"""[x,w] = us_roots(n)
Returns the roots (x) of the nth order shifted Chebyshev (of the second kind)
@@ -601,7 +627,7 @@
"""
return js_roots(n,2.0,1.5,mu=mu)
-def sh_chebyu(n,monic=0):
+def sh_chebyu(n, monic=0):
"""Return nth order shifted Chebyshev polynomial of second kind, Un(x).
Orthogonal over [0,1] with weight function (x-x**2)**(1/2).
"""
@@ -612,7 +638,7 @@
return base
# Legendre
-def p_roots(n,mu=0):
+def p_roots(n, mu=0):
"""[x,w] = p_roots(n)
Returns the roots (x) of the nth order Legendre polynomial, P_n(x),
@@ -621,11 +647,13 @@
"""
return j_roots(n,0.0,0.0,mu=mu)
-def legendre(n,monic=0):
+def legendre(n, monic=0):
"""Returns the nth order Legendre polynomial, P_n(x), orthogonal over
[-1,1] with weight function 1.
"""
- assert(n>=0), "n must be nonnegative"
+ if n < 0:
+ raise ValueError("n must be nonnegative.")
+
if n==0: n1 = n+1
else: n1 = n
x,w,mu0 = p_roots(n1,mu=1)
@@ -637,7 +665,7 @@
return p
# Shifted Legendre P^*_n(x)
-def ps_roots(n,mu=0):
+def ps_roots(n, mu=0):
"""[x,w] = ps_roots(n)
Returns the roots (x) of the nth order shifted Legendre polynomial, P^*_n(x),
@@ -646,11 +674,13 @@
"""
return js_roots(n,1.0,1.0,mu=mu)
-def sh_legendre(n,monic=0):
+def sh_legendre(n, monic=0):
"""Returns the nth order shifted Legendre polynomial, P^*_n(x), orthogonal
over [0,1] with weighting function 1.
"""
- assert(n>=0), "n must be nonnegative"
+ if n < 0:
+ raise ValueError("n must be nonnegative.")
+
wfunc = lambda x: 0.0*x + 1.0
if n==0: return orthopoly1d([],[],1.0,1.0,wfunc,(0,1),monic,
lambda x: eval_sh_legendre(n,x))
More information about the Scipy-svn
mailing list