The Gamma function is a complex special function such that $\Gamma(s) = (s-1)!$ for positive integers $s$.
{{{id=46| special.gamma(10) /// 362880.0 }}} {{{id=49| special.gamma(2+2j) /// (0.11229424234632605+0.32361288550192724j) }}} {{{id=47| complex_plot(special.gamma, (-2,2), (-2,2)) ///
}}}
All the scipy special functions can also be applied to numpy arrays (often very fast!).
{{{id=52| v = sp.r_[.5:4:10000j]; v /// array([ 0.5 , 0.50035004, 0.50070007, ..., 3.99929993, 3.99964996, 4. ]) }}} {{{id=51| special.gamma(v) /// array([ 1.77245385, 1.7712366 , 1.77002126, ..., 5.99472651, 5.99736257, 6. ]) }}}Matlab-style data plotting is supported through plt (remember, we did import matplotlib.pyplot as plt above).
{{{id=50| plt.clf() # clear the current figure plt.plot(v, special.gamma(v)) # plot(array of x vals, array of y values) plt.savefig('a.png') # makes plot get displayed /// }}}Here is an example from the Scipy tutorial that involves using a special Bessel function to model the vibrational modes of a thin drum head. Here we use Sage's integrated 3d plotting capabilities instead of matplotlib's pseudo-3d plotting.
{{{id=56| from scipy.special import jn, jn_zeros def drumhead_height(n, k, distance, angle, t): nth_zero = jn_zeros(n, k) # zeros of Bessel function return sp.cos(t)*sp.cos(n*angle)*jn(n, distance*nth_zero) # Now use a Sage function to create the 3d visualization pnt = (lambda r,theta: r*sp.cos(theta), lambda r,theta: r*sp.sin(theta), lambda r,theta: drumhead_height(1,1,r,theta,0.5)) G = parametric_plot3d(pnt, (0,1), (0,2*pi), plot_points=50) /// }}} {{{id=54| G.show() # default plotting requires java (e.g., not iphone/ipad/chrome friendly) /// }}} {{{id=58| G.show(viewer='tachyon') /// }}} {{{id=59| G.show(viewer='canvas3d') # works on everything (except Internet Explorer...) /// }}} {{{id=65| /// }}}Integration
The scipy.integrate module provides numerical integration and solving of ordinary differential equations.
{{{id=67| import scipy.integrate /// }}} {{{id=64| scipy.integrate? /// }}}Example: quad -- integrate a 1d functions over an interval
{{{id=63| from scipy.special import jv result = scipy.integrate.quad(lambda x: jv(2.5, x), 0, 4.5) result /// (1.1178179380783249, 7.8663171825372256e-09) }}} {{{id=69| plot(lambda x: jv(2.5, x), 0, 4.5, fill=True) ///
}}}
Optimization
The scipy.optimize module provides common optimization algorithms.
{{{id=95| import scipy.optimize /// }}} {{{id=94| scipy.optimize? /// }}}Example: Uses the Simplex Algorithm to minimize a function:
{{{id=96| def f(x): return x[0]*x[0] + x[1]*x[1] + 10 * cos(x[0])*cos(x[1]) A = scipy.optimize.fmin(f, (0,0)); A /// Optimization terminated successfully. Current function value: -1.808983 Iterations: 70 Function evaluations: 132 array([ 2.59576477e+00, -2.62130635e-05]) }}} {{{id=70| G = plot3d(lambda x,y: f((x,y)), (-1,4), (-1,1), opacity=.7) G += sphere((A[0],A[1],f(A)), color='red', size=.1) G.show() /// }}} {{{id=103| /// }}}Example: Least squares fit to noisy data.
{{{id=102| from scipy.optimize import leastsq from scipy import pi x = sp.r_[0:6e-2:30j] A = 10; k=1.0/3e-2; theta = pi/6 y_true = A * sp.sin(2*pi*k*x + theta) y_measured = y_true + 2*sp.random.randn(len(x)) # add random array def residuals(p, y, x): A, k, theta = p return y - A*sp.sin(2*pi*k*x + theta) def peval(x, p): return p[0]*sp.sin(2*pi*p[1]*x + p[2]) p0 = [8, 1/2.3e-2, pi/3] from scipy.optimize import leastsq plsq = leastsq(residuals, p0, args=(y_measured, x)) plt.clf() plt.plot(x, peval(x, plsq[0]), x, y_measured, 'o', x, y_true) plt.title('Least-squares fit to noisy data') plt.legend(['Fit', 'Noisy', 'True']) plt.savefig('a.png') /// }}} {{{id=101| /// }}} {{{id=98| /// }}} {{{id=73| /// }}}Interpolation
The scipy.interpolate package provides single and multivariate interpolation capabilities.
{{{id=109| import scipy.interpolate /// }}} {{{id=108| scipy.interpolate? /// }}}Example: 1-d interpolation
The function interp1d in scipy.interpolate creates an object that you can treat like a function that evaluate to the interpolation through a given list of x,y points in the plane (that partially define a *function*).
{{{id=107| random_vec = sp.random.randn(len(x)) @interact def demo(kind=['linear','nearest', 'zero', 'slinear', 'quadratic', 'cubic'], points=(50,(i for i in xrange(3,1000))), noisy=True): import scipy.interpolate # interpolate the function exp(-x/3): x = sp.r_[0:10:10j] y = sp.exp(-x/3.0) if noisy: y += 0.05*random_vec f = scipy.interpolate.interp1d(x, y, kind=kind) x_new = sp.r_[0:10:int(points)*1j] plt.clf() plt.plot(x, y, 'o', x_new, f(x_new), '-') plt.legend(['data', kind]) plt.axis([0,10,0,1]) plt.savefig('a.png') ///
|
|||||||||||||||||||
Fourier Transform
The scipy.fftpack package provides numerical Fast Fourier Transform algorithms.
{{{id=113| import scipy.fftpack /// }}} {{{id=112| scipy.fftpack? /// }}} {{{id=118| plt.savefig( /// }}} {{{id=116| @interact def demo(frequencies=input_box("1 2", type=str)): print frequencies, type(frequencies) frequencies = [eval(a) for a in frequencies.split()] x = sp.r_[0:4*sp.pi:200j] y = sum(sp.sin(f*x) for f in frequencies) z = sp.fftpack.fft(y) plt.clf() plt.plot(x, y) plt.savefig('a.png', dpi=40) plt.clf() plt.plot(sp.arange(len(z)), abs(z)) plt.savefig('b.png', dpi=40) ///
|
|||||||||||||
Signal Processing
The scipy.signal package provides a "signal processing toolbox", including filtering functions, and spline interpolation algorithms.
{{{id=114| import scipy.signal /// }}} {{{id=79| scipy.signal? /// }}} {{{id=82| /// }}}Linear Algebra
The scipy.linalg package provides algorithms for matrix and linear algebra, including a range of decomposition algorithms, special matrices, etc.
{{{id=121| import scipy.linalg /// }}} {{{id=84| scipy.linalg? /// }}} {{{id=87| /// }}}Statistics
The scipy.stats package has a large number of basic descriptive statistic routines.
{{{id=86| import scipy.stats /// }}} {{{id=119| scipy.stats? /// }}} {{{id=89| /// }}}Image Processing
The scipy.ndimage package provides general image processing and analysis algorithms.
{{{id=85| import scipy.ndimage /// }}} {{{id=120| scipy.ndimage? /// }}} {{{id=91| /// }}}IO Module
The scipy.io module lets you read and write various special file types, e.g., matlab matrices.
{{{id=74| import scipy.io /// }}} {{{id=93| scipy.io? /// }}}