1 \documentclass[12pt]{article}

     2 \usepackage{amsmath}

     3 \title{Python Workshop\\Problems and Exercises}

     4 \author{Asokan Pichai\\Prabhu Ramachandran}

     5 \begin{document}

     6 \maketitle

     7 

     8 \section{Matrices and Arrays \& 2D Plotting}

     9 \subsection{Matrices and Arrays}

    10 \begin{verbatim}

    11 # Simple array math example

    12 >>> import numpy as np

    13 >>> a = np.array([1,2,3,4])

    14 >>> b = np.arange(2,6)

    15 >>> b

    16 array([2,3,4,5])

    17 >>> a*2 + b + 1 # Basic math!

    18 array([5, 8, 11, 14])

    19 

    20 # Pi and e are defined.

    21 >>> x = np.linspace(0.0, 10.0, 1000)

    22 >>> x *= 2*np.pi/10 # inplace.

    23 # apply functions to array.

    24 >>> y = np.sin(x)

    25 >>> z = np.exp(y)

    26 

    27 >>> x = np.array([1., 2, 3, 4])

    28 >>> np.size(x)

    29 4

    30 >>> x.dtype # What is a.dtype?

    31 dtype('float64')

    32 >>> x.shape

    33 (4,)

    34 >>> print np.rank(x), x.itemsize

    35 1 8

    36 >>> x[0] = 10

    37 >>> print x[0], x[-1]

    38 10.0 4.0

    39 

    40 >>> a = np.array([[ 0, 1, 2, 3],

    41 ...            [10,11,12,13]])

    42 >>> a.shape # (rows, columns)

    43 (2, 4)

    44 # Accessing and setting values

    45 >>> a[1,3] 

    46 13

    47 >>> a[1,3] = -1

    48 >>> a[1] # The second row

    49 array([10,11,12,-1])

    50 

    51 >>> b=np.array([[0,2,4,2],[1,2,3,4]])

    52 >>> np.add(a,b,a)

    53 >>> np.sum(x,axis=1)

    54 

    55 >>> np.greater(a,4)

    56 >>> np.sqrt(a)

    57 

    58 >>> np.array([2,3,4])  

    59 array([2, 3, 4])

    60 

    61 >>> np.linspace(0, 2, 4)   

    62 array([0.,0.6666667,1.3333333,2.])

    63 

    64 >>>np.ones([2,2])

    65 array([[ 1.,  1.],

    66      [ 1.,  1.]])

    67 

    68 >>>a = np.array([[1,2,3],[4,5,6]])

    69 >>>np.ones_like(a)

    70 array([[1, 1, 1],

    71        [1, 1, 1]])

    72 

    73 >>> a = np.array([[1,2,3], [4,5,6], 

    74                [7,8,9]])

    75 >>> a[0,1:3]

    76 array([2, 3])

    77 >>> a[1:,1:]

    78 array([[5, 6],

    79        [8, 9]])

    80 >>> a[:,2]

    81 array([3, 6, 9])

    82 >>> a[...,2]

    83 array([3, 6, 9])

    84 

    85 >>> a[0::2,0::2]

    86 array([[1, 3],

    87        [7, 9]])

    88 # Slices are references to the 

    89 # same memory!

    90 

    91 >>> np.random.rand(3,2)

    92 array([[ 0.96276665,  0.77174861],

    93        [ 0.35138557,  0.61462271],

    94        [ 0.16789255,  0.43848811]])

    95 >>> np.random.randint(1,100)

    96 42

    97 \end{verbatim}

    98 

    99 \subsection{Problem Set}

   100 \begin{verbatim}

   101     >>> from scipy import misc

   102     >>> A=misc.imread(name)

   103     >>> misc.imshow(A)

   104 \end{verbatim}

   105 \begin{enumerate}

   106   \item Convert an RGB image to Grayscale. $ Y = 0.5R + 0.25G + 0.25B $.

   107   \item Scale the image to 50\%.

   108   \item Introduce some random noise.

   109   \item Smooth the image using a mean filter.

   110   \\\small{Each element in the array is replaced by mean of all the neighbouring elements}

   111   \\\small{How fast does your code run?}

   112 \end{enumerate}

   113 

   114 \subsection{2D Plotting}

   115 \begin{verbatim}

   116 $ ipython -pylab

   117 >>> x = linspace(0, 2*pi, 1000)

   118 >>> plot(x, sin(x)) 

   119 >>> plot(x, sin(x), 'ro')

   120 >>> xlabel(r'$\chi$', color='g')

   121 # LaTeX markup!

   122 >>> ylabel(r'sin($\chi$)', color='r')

   123 >>> title('Simple figure', fontsize=20)

   124 >>> savefig('/tmp/test.eps')

   125 

   126 # Set properties of objects:

   127 >>> l, = plot(x, sin(x))

   128 # Why "l,"?

   129 >>> setp(l, linewidth=2.0, color='r')

   130 >>> l.set_linewidth(2.0)

   131 >>> draw() # Redraw.

   132 >>> setp(l) # Print properties.

   133 >>> clf() # Clear figure.

   134 >>> close() # Close figure.

   135 

   136 >>> w = arange(-2,2,.1)

   137 >>> plot(w,exp(-(w*w))*cos)

   138 >>> ylabel('$f(\omega)$')

   139 >>> xlabel('$\omega$')

   140 >>> title(r"$f(\omega)=e^{-\omega^2}

   141             cos({\omega^2})$")

   142 >>> annotate('maxima',xy=(0, 1), 

   143              xytext=(1, 0.8), 

   144              arrowprops=dict(

   145              facecolor='black', 

   146              shrink=0.05))

   147 

   148 >>> x = linspace(0, 2*np.pi, 1000)

   149 >>> plot(x, cos(5*x), 'r--', 

   150          label='cosine')

   151 >>> plot(x, sin(5*x), 'g--', 

   152          label='sine')

   153 >>> legend() 

   154 # Or use:

   155 >>> legend(['cosine', 'sine'])

   156 

   157 >>> figure(1)

   158 >>> plot(x, sin(x))

   159 >>> figure(2)

   160 >>> plot(x, tanh(x))

   161 >>> figure(1)

   162 >>> title('Easy as 1,2,3')

   163 

   164 \end{verbatim}

   165 

   166 \subsection{Problem Set}

   167   \begin{enumerate}

   168     \item Write a function that plots any regular n-gon given n.

   169     \item Consider the logistic map, $f(x) = kx(1-x)$, plot it for

   170           $k=2.5, 3.5$ and $4$ in the same plot.

   171 

   172   \item Consider the iteration $x_{n+1} = f(x_n)$ where $f(x) =

   173         kx(1-x)$.  Plot the successive iterates of this process.

   174   \item Plot this using a cobweb plot as follows:

   175   \begin{enumerate}

   176     \item Start at $(x_0, 0)$

   177     \item Draw line to $(x_i, f(x_i))$; 

   178     \item Set $x_{i+1} = f(x_i)$

   179     \item Draw line to $(x_i, x_i)$

   180     \item Repeat from 2 for as long as you want 

   181   \end{enumerate}

   182 \end{enumerate}

   183 

   184 \section{Advanced Numpy}

   185 

   186 \subsection{Broadcasting}

   187 \begin{verbatim}

   188 >>> a = np.arange(4)

   189 >>> b = np.arange(5)

   190 >>> a+b

   191 >>> a+3

   192 >>> c=np.array([3])

   193 >>> a+c

   194 >>> b+c

   195 

   196 >>> a = np.arange(4)

   197 >>> a+3

   198 array([3, 4, 5, 6])

   199 

   200 >>> x = np.ones((3, 5))

   201 >>> y = np.ones(8)

   202 >>> (x[..., None] + y).shape

   203 (3, 5, 8)

   204 

   205 \end{verbatim}

   206 

   207 \subsection{Copies \& Views}

   208 \begin{verbatim}

   209 >>> a = np.array([[1,2,3],[4,5,6]])

   210 >>> b = a

   211 >>> b is a

   212 >>> b[0,0]=0; print a

   213 >>> c = a.view()

   214 >>> c is a

   215 >>> c.base is a

   216 >>> c.flags.owndata

   217 >>> d = a.copy()

   218 >>> d.base is a

   219 >>> d.flags.owndata

   220 

   221 >>> a = np.arange(1,9)

   222 >>> a.shape=3,3

   223 >>> b = a[0,1:3]

   224 >>> c = a[0::2,0::2]

   225 >>> a.flags.owndata

   226 >>> b.flags.owndata

   227 >>> b.base

   228 >>> c.base is a

   229 

   230 >>> b = a[np.array([0,1,2])]

   231 array([[1, 2, 3],

   232        [4, 5, 6],

   233        [7, 8, 9]])

   234 >>> b.flags.owndata

   235 >>> abool=np.greater(a,2)

   236 >>> c = a[abool]

   237 >>> c.flags.owndata

   238 

   239 \end{verbatim}

   240 

   241 \section{Scipy}

   242 \subsection{Linear Algebra}

   243 \begin{verbatim}

   244 >>> import scipy as sp

   245 >>> from scipy import linalg

   246 >>> A=sp.mat(np.arange(1,10))

   247 >>> A.shape=3,3

   248 >>> linalg.inv(A)

   249 >>> linalg.det(A)

   250 >>> linalg.norm(A)

   251 >>> linalg.expm(A) #logm

   252 >>> linalg.sinm(A) #cosm, tanm, ...

   253 

   254 >>> A = sp.mat(np.arange(1,10))

   255 >>> A.shape=3,3

   256 >>> linalg.lu(A)

   257 >>> linalg.eig(A)

   258 >>> linalg.eigvals(A)

   259 \end{verbatim}

   260 

   261 \subsection{Solving Linear Equations}

   262 

   263 \begin{align*}

   264   3x + 2y - z  & = 1 \\

   265   2x - 2y + 4z  & = -2 \\

   266   -x + \frac{1}{2}y -z & = 0

   267 \end{align*}

   268 To Solve this, 

   269 \begin{verbatim}

   270 >>> A = sp.mat([[3,2,-1],[2,-2,4]

   271               ,[-1,1/2,-1]])

   272 >>> B = sp.mat([[1],[-2],[0]])

   273 >>> linalg.solve(A,B)

   274 \end{verbatim}

   275 

   276 \subsection{Integrate}

   277 \subsubsection{Quadrature}

   278 Calculate the area under $(sin(x) + x^2)$ in the range $(0,1)$

   279 \begin{verbatim}

   280 >>> def f(x):

   281         return np.sin(x)+x**2

   282 >>> integrate.quad(f, 0, 1)

   283 \end{verbatim}

   284 

   285 \subsubsection{ODE Integration}

   286 Numerically solve ODEs\\

   287 \begin{align*}

   288 \frac{dx}{dt} &=-e^{-t}x^2\\ 

   289          x(0) &=2    

   290 \end{align*}

   291 \begin{verbatim}

   292 >>> def dx_dt(x,t):

   293 ...     return -np.exp(-t)*x**2

   294 >>> x=integrate.odeint(dx_dt, 2, t)

   295 >>> plt.plot(x,t)

   296 \end{verbatim}

   297 

   298 \subsection{Interpolation}

   299 \subsubsection{1D Interpolation}

   300 \begin{verbatim}

   301 >>> from scipy import interpolate

   302 >>> interpolate.interp1d?

   303 >>> x = np.arange(0,2*np.pi,np.pi/4)

   304 >>> y = np.sin(x)

   305 >>> fl = interpolate.interp1d(

   306             x,y,kind='linear')

   307 >>> fc = interpolate.interp1d(

   308              x,y,kind='cubic')

   309 >>> fl(np.pi/3)

   310 >>> fc(np.pi/3)

   311 \end{verbatim}

   312 

   313 \subsubsection{Splines}

   314 Plot the Cubic Spline of $sin(x)$

   315 \begin{verbatim}

   316 >>> x = np.arange(0,2*np.pi,np.pi/4)

   317 >>> y = np.sin(x)

   318 >>> tck = interpolate.splrep(x,y)

   319 >>> X = np.arange(0,2*np.pi,np.pi/50)

   320 >>> Y = interpolate.splev(X,tck,der=0)

   321 >>> plt.plot(x,y,'o',x,y,X,Y)

   322 >>> plt.show()

   323 \end{verbatim}

   324 

   325 \subsection{Signal \& Image Processing}

   326 Applying a simple median filter

   327 \begin{verbatim}

   328 >>> from scipy import signal, ndimage

   329 >>> from scipy import lena

   330 >>> A=lena().astype('float32')

   331 >>> B=signal.medfilt2d(A)

   332 >>> imshow(B)

   333 \end{verbatim}

   334 Zooming an array - uses spline interpolation

   335 \begin{verbatim}

   336 >>> b=ndimage.zoom(A,0.5)

   337 >>> imshow(b)

   338 \end{verbatim}

   339 

   340 \section{3D Data Visualization}

   341 \subsection{Using mlab}

   342 \begin{verbatim}

   343 >>> from enthought.mayavi import mlab

   344 

   345 >>> mlab.test_<TAB>

   346 >>> mlab.test_contour3d()

   347 >>> mlab.test_contour3d??

   348 \end{verbatim}

   349 

   350 \subsubsection{Plotting Functions}

   351 \begin{verbatim}

   352 >>> from numpy import *

   353 >>> t = linspace(0, 2*pi, 50)

   354 >>> u = cos(t)*pi

   355 >>> x, y, z = sin(u), cos(u), sin(t)

   356 

   357 >>> x = mgrid[-3:3:100j,-3:3:100j]

   358 >>> z = sin(x*x + y*y)

   359 >>> mlab.surf(x, y, z)

   360 \end{verbatim}

   361 

   362 \subsubsection{Large 2D Data}

   363 \begin{verbatim}

   364 >>> mlab.mesh(x, y, z)

   365 

   366 >>> phi, theta = numpy.mgrid[0:pi:20j, 

   367 ...                         0:2*pi:20j]

   368 >>> x = sin(phi)*cos(theta)

   369 >>> y = sin(phi)*sin(theta)

   370 >>> z = cos(phi)

   371 >>> mlab.mesh(x, y, z, 

   372 ...           representation=

   373 ...           'wireframe')

   374 \end{verbatim}

   375 

   376 \subsubsection{Large 3D Data}

   377 \begin{verbatim}

   378 >>> x, y, z = ogrid[-5:5:64j, 

   379 ...                -5:5:64j, 

   380 ...                -5:5:64j]

   381 >>> mlab.contour3d(x*x*0.5 + y*y + 

   382                    z*z*2)

   383 

   384 >>> mlab.test_quiver3d()

   385 \end{verbatim}

   386 

   387 \subsection{Motivational Problem}

   388 Atmospheric data of temperature over the surface of the earth. Let temperature ($T$) vary linearly with height ($z$)\\

   389 $T = 288.15 - 6.5z$

   390 

   391 \begin{verbatim}

   392 lat = linspace(-89, 89, 37)

   393 lon = linspace(0, 360, 37)

   394 z = linspace(0, 100, 11)

   395 \end{verbatim}

   396 

   397 \begin{verbatim}

   398 x, y, z = mgrid[0:360:37j,-89:89:37j,

   399                 0:100:11j]

   400 t = 288.15 - 6.5*z

   401 mlab.contour3d(x, y, z, t)

   402 mlab.outline()

   403 mlab.colorbar()

   404 \end{verbatim}

   405 

   406 \subsection{Lorenz equation}

   407 \begin{eqnarray*}

   408 \frac{d x}{dt} &=& s (y-x)\\

   409 \frac{d y}{d t} &=& rx -y -xz\\

   410 \frac{d z}{d t} &=& xy - bz\\

   411 \end{eqnarray*}

   412 

   413 Let $s=10,$

   414 $r=28,$ 

   415 $b=8./3.$

   416 

   417 \begin{verbatim}

   418 x, y, z = mgrid[-50:50:20j,-50:50:20j,

   419                 -10:60:20j]

   420 

   421 \end{verbatim}

   422 

   423 \end{document}