analysis: analysis::hermite_numerov::hermite_gauss

def analysis.hermite_numerov.hermite_gauss_integral	(	m,
		n,
		x0,
		x1,
		stepsize = `0.01`
	)

Definition at line 85 of file hermite_numerov.py.

00085                                                        :
00086     "compute integral(h[m](x}*h[n](x)*exp(-x^2),{x,x0,x1})"
00087     if m>n:
00088         m,n = n,m #always make m<=n
00089     xmax=max(abs(x0), abs(x1))
00090     if not integral_cache.has_key((m,n)) or integral_cache[(m,n)][0] < xmax:
00091         #just use Simpson's rule for this
00092         #first, load up cache
00093         hermite_gauss(m,xmax)
00094         hermite_gauss(n,xmax)
00095         xm, xa1, ya1, y2a1 = hermite_cache[m]
00096         xm, xa2, ya2, y2a2 = hermite_cache[n]
00097 
00098         stepcount=int(2.0*xmax/stepsize)
00099         if stepcount%2==0:
00100             stepcount=stepcount+1 #always odd for Simpson
00101         stepsize=2*xmax/(stepcount-1) #exact step size
00102         xlist=Numeric.array(range(stepcount),Numeric.Float)*stepsize - xmax
00103         y=spline.splint(xa1, ya1, y2a1, xlist)*spline.splint(xa2, ya2, y2a2, xlist)
00104         hmn2=spline.spline(xlist, y)
00105         yint=Numeric.cumsum(y[0:-2:2]+y[2::2]+4.0*y[1::2])*(stepsize/3.0)
00106         yint=Numeric.concatenate( ( (0,), yint )) #first element of integral is 0
00107         yi2=spline.spline(xlist[::2], yint)
00108         integral_cache[(m,n)]=(xmax, xlist[::2], yint, yi2, stepsize, xlist, y, hmn2)
00109 
00110     xmax, xa, ya, y2a, stepjunk, xjunk, yjunk, hmn2junk=integral_cache[(m,n)]
00111     iv1, iv2=spline.splint(xa, ya, y2a, (x0, x1) )
00112     return iv2-iv1
00113 
def hermite_gauss_matrix_element(m, n, x0, x1, f, stepsize=0.01):