1 / 232

Scientific Computing using Python

Scientific Computing using Python. Two day PET KY884 tutorial at HEAT Center, Aberdeen, MD Mon-Tue, July 20-21, 2009, 8:30am - 4:30pm Dr. Craig Rasmussen, and Dr. Sameer Shende info@paratools.com http://www.paratools.com/arl09. Tutorial Outline. Basic Python

csheryl
Download Presentation

Scientific Computing using Python

An Image/Link below is provided (as is) to download presentation Download Policy: Content on the Website is provided to you AS IS for your information and personal use and may not be sold / licensed / shared on other websites without getting consent from its author. Content is provided to you AS IS for your information and personal use only. Download presentation by click this link. While downloading, if for some reason you are not able to download a presentation, the publisher may have deleted the file from their server. During download, if you can't get a presentation, the file might be deleted by the publisher.

E N D

Presentation Transcript

  1. Scientific Computing using Python Two day PET KY884 tutorial at HEAT Center, Aberdeen, MD Mon-Tue, July 20-21, 2009, 8:30am - 4:30pm Dr. Craig Rasmussen, and Dr. Sameer Shende info@paratools.com http://www.paratools.com/arl09

  2. Tutorial Outline • Basic Python • IPython : Interactive Python • Advanced Python • NumPy : High performance arrays for Python • Matplotlib : Basic plotting tools for Python • MPI4py : Parallel programming with Python • F2py and SWIG : Language interoperability • Extra Credit • SciPy and SAGE : mathematical and scientific computing • Traits : Typing system for Python • Dune : A Python CCA-Compliant component framework • Portable performance evaluation using TAU • Labs

  3. Tutorial Goals • This tutorial is intended to introduce Python as a tool for high productivity scientific software development. • Today you should leave here with a better understanding of… • The basics of Python, particularly for scientific and numerical computing. • Toolkits and packages relevant to specific numerical tasks. • How Python is similar to tools like MATLAB or GNU Octave. • How Python might be used as a component architecture • …And most importantly, • Python makes scientific programming easy, quick, and fairly painless, leaving more time to think about science and not programming.

  4. SECTION 1 INTRODUCTION

  5. Interpreted and interactive Easy to learn and use Fun Free and portable Automatic garbage collection Object-oriented and Functional Truly Modular NumPy PySparse FFTW Plotting MPI4py Co-Array Python What Is Python? Python is an interpreted language that allows you to accomplish what you would with a compiled language, but without the complexity.

  6. $$ ipython Python 2.5.1 (… Feb 6 2009 …) Ipython 0.9.1 – An enhanced … '''a comment line …''' # another comment style # the IPython prompt In [1]: # the Python prompt, when native # python interpreter is run >>> # import a module >>> import math # what is math >>> type(math) <type 'module'> # what is in math >>> dir(math) ['__doc__',…, 'cos',…, pi, …] >>> cos(pi) NameError: name 'cos' is not defined # import into global namespace >>> from math import * >>> cos(pi) -1.0 Running Python

  7. # adding two values >>> 3 + 4 7 # setting a variable >>> a = 3 >>> a 3 # checking a variables type >>> type(a) <type 'int'> # an arbitrarily long integer >>> a = 1204386828483L >>> type(a) <type 'long'> # real numbers >>> b = 2.4/2 >>> print b 1.2 >>> type(b) <type 'float'> # complex numbers >>> c = 2 + 1.5j >>> c (2+1.5j) # multiplication >>> a = 3 >>> a*c (6+4.5j) Interactive Calculator

  8. # command line documentation $$ pydoc math Help on module math: >>> dir(math) ['__doc__', >>> math.__doc__ …mathematical functions defined… >>> help(math) Help on module math: >>> type(math) <type 'module'> # ipython documentation In[3]: math.<TAB> …math.pi math.sin math.sqrt In[4]: math? Type: module Base Class: <type 'module'> In[5]: import numpy In[6]: numpy?? Source:=== “““\ NumPy ========= Online Python Documentation

  9. Labs! Lab: Explore and Calculate

  10. # creating strings >>> s1 = "Hello " >>> s2 = 'world!' # string operations >>> s = s1 + s2 >>> print s Hello world! >>> 3*s1 'Hello Hello Hello ' >>> len(s) 12 # the string module >>> import string # split space delimited words >>> word_list = string.split(s) >>> print word_list ['Hello', 'world!'] >>> string.join(word_list) 'Hello world!' >>> string.replace(s,'world', 'class') 'Hello class!' Strings

  11. Labs! Lab: Strings

  12. # a tuple is a collection of obj >>> t = (44,) # length of one >>> t = (1,2,3) >>> print t (1,2,3) # accessing elements >>> t[0] 1 >>> t[1] = 22 TypeError: 'tuple' object does not support item assignment # a list is a mutable collection >>> l = [1,22,3,3,4,5] >>> l [1,22,3,3,4,5] >>> l[1] = 2 >>> l [1,2,3,3,4,5] >>> del l[2] >>> l [1,2,3,4,5] >>> len(l) 5 # in or not in >>> 4 in l True >>> 4 not in l False Tuples and Lists: sequence objects

  13. # negative indices count # backward from the end of # the list >>> l [1,2,3,4,5] >>> l[-1] 5 >>> l[-2] 4 >>> dir(list) [__add__, 'append', 'count', 'extend', 'index', 'insert', 'pop', 'remove', 'reverse', 'sort'] # what does count do? >>> list.count <method 'count' of 'list'…> >>> help(list.count) 'L.count(value) -> integer -- return number of occurrences of value' More on Lists

  14. >>> print l [1,2,3,4,5] # some ways to return entire # portion of the sequence >>> l[0:5] >>> l[0:] >>> l[:5] >>> l[:] [1,2,3,4,5] # middle three elements >>> l[1:4] >>> l[1:-1] >>> l[-4:-1] [2,3,4] # last two elements >>> l[3:] >>> l[-2:] [4,5] Slicing var[lower:upper] Slices extract a portion of a sequence (e.g., a list or a NumPy array). Mathematically the range is [lower, upper).

  15. # create data >>> pos = [1.0, 2.0, 3.0, 4.0, 5.0] >>> T = [9.9, 8.8. 7.7, 6.6, 5.5] # store data in a dictionary >>> data_dict = {'position': pos, 'temperature': T} # access elements >>> data_dict['position'] [1.0, 2.0, 3.0, 4.0, 5.0] Dictionaries: key/value pairs Dictionaries store key/value pairs. Indexing a dictionary by a key returns the value associate with it.

  16. Labs! Lab: Sequence Objects

  17. # if/elif/else example >>> print l [1,2,3,4,5] >>> if 3 in l: … print 'yes' … elif 3 not in l: … print 'no' … else: … print 'impossible!' … < hit return > yes # for loop examples >>> for i in range(1,3): print i … < hit return > 1 2 >>> for x in l: print x … < hit return > 1 … # while loop example >>> i = 1 >>> while i < 3: print i; i += 1 … < hit return > 1 2 If Statements and Loops

  18. # create a function in funcs.py def Celcius_to_F(T_C): T_F = (9./5.)*T_C + 32. return T_F ''' Note: indentation is used for scoping, no braces {} ''' # run from command line and # start up with created file $ python -i funcs.py >>> dir() ['Celcius_to_F', '__builtins__',… ' >>> Celsius_to_F = Celcius_to_F >>> Celsius_to_F <function Celsius_to_F at …> >>> Celsius_to_F(0) 32.0 >>> C = 100. >>> F = Celsius_to_F(C) >>> print F 212.0 Functions

  19. Labs! Lab: Functions

  20. # create a class in Complex.py class Complex: '''A simple Complex class''' def __init__(self, real, imag): '''Create and initialize''' self.real = real self.imag = imag def norm(self): '''Return the L2 Norm''' import math d = math.hypot(self.real,self.imag) return d #end class Complex # run from command line $ python -i Complex.py # help will display comments >>> help(Complex) Help on class Complex in module … # create a Complex object >>> c = Complex(3.0, -4.0) # print Complex attributes >>> c.real 3.0 >>> c.imag -4.0 # execute a Complex method >>> c.norm() 5.0 Classes

  21. Labs! Lab: Classes

  22. SECTION 2 Interactive Python IPython

  23. IPython Summary • An enhanced interactive Python shell • An architecture for interactive parallel computing • IPython contains • Object introspection • System shell access • Special interactive commands • Efficient environment for Python code development • Embeddable interpreter for your own programs • Inspired by Matlab • Interactive testing of threaded graphical toolkits

  24. $$ ipython -pylab IPython 0.9.1 -- An enhanced Interactive Python. ? -> Introduction and overview of IPython's features. %quickref -> Quick reference. help -> Python's own help system. object? -> Details about 'object'. ?object also works, ?? Prints # %fun_name are magic commands # get function info In [1]: %history? Print input history (_i<n> variables), with most recent last. In [2]: %history 1: #?%history 2: _ip.magic("history ") Running IPython

  25. # some shell commands are available In [27]: ls 01-Lab-Explore.ppt* 04-Lab-Functions.ppt* # TAB completion for more information about objects In [28]: %<TAB> %alias %autocall %autoindent %automagic %bg %bookmark %cd %clear %color_info %colors %cpaste %debug %dhist %dirs %doctest_mode # retrieve Out[] values In [29]: 4/2 Out[29]: 2 In [30]: b = Out[29] In [31]: print b 2 More IPython Commands

  26. # %run runs a Python script and loads its data into interactive # namespace; useful for programming In [32]: %run hello_script Hello # ! gives access to shell commands In [33]: !date Tue Jul 7 23:04:37 MDT 2009 # look at logfile (see %logstart and %logstop) In [34]: !cat ipython_log.py #log# Automatic Logger file. *** THIS MUST BE THE FIRST LINE *** #log# DO NOT CHANGE THIS LINE OR THE TWO BELOW #log# opts = Struct({'__allownew': True, 'logfile': 'ipython_log.py'}) #log# args = [] #log# It is safe to make manual edits below here. #log#----------------------------------------------------------------------- _ip.magic("run hello”) More IPython Commands

  27. Interactive Shell Recap • Object introspection (? and ??) • Searching in the local namespace (‘TAB’) • Numbered input/output prompts with command history • User-extensible ‘magic’ commands (‘%’) • Alias facility for defining your own system aliases • Complete system shell access • Background execution of Python commands in a separate thread • Expand python variables when calling the system shell • Filesystem navigation via a magic (‘%cd’) command • Bookmark with (‘%bookmark’) • A lightweight persistence framework via the (‘%store’) command • Automatic indentation (optional) • Macro system for quickly re-executing multiple lines of previous input • Session logging and restoring • Auto-parentheses (‘sin 3’) • Easy debugger access (%run –d) • Profiler support (%prun and %run –p)

  28. Labs! Lab: IPython Try out ipython commands as time allows

  29. SECTION 3 Advanced Python

  30. # The re module provides regular expression tools for advanced # string processing. >>> import re # Get a refresher on regular expressions >>> help(re) >>> help(re.findall) >>> help(re.sub) >>> re.findall(r'\bf[a-z]*', 'which foot or hand fell fastest') ['foot', 'fell', 'fastest’] >>> re.sub(r'(\b[a-z]+) \1', r'\1', 'cat in the the hat') 'cat in the hat' Regular Expressions

  31. Labs! Lab: Regular Expressions Try out the re module as time allows

  32. # a filter returns those items # for which the given function returns True >>> def f(x): return x < 3 >>> filter(f, [0,1,2,3,4,5,6,7]) [0, 1, 2] # map applies the given function to each item in a sequence >>> def square(x): return x*x >>> map(square, range(7)) [0, 1, 4, 9, 16, 25, 36] # lambda functions are small functions with no name (anonymous) >>> map(lambda x: x*x, range(7)) [0, 1, 4, 9, 16, 25, 36] Fun With Functions

  33. # reduce returns a single value by applying a binary function >>> reduce(lambda x,y: x+y, [0,1,2,3]) 6 # list comprehensions provide an easy way to create lists # [an expression followed by for then zero or more for or if] >>> vec = [2, 4, 6] >>> [3*x for x in vec] [6, 12, 18] >>> [3*x for x in vec if x > 3] [12, 18] >>> [x*y for x in vec for y in [3, 2, -1]] [6, 4, -2, 12, 8, -4, 18, 12, -6] More Fun With Functions

  34. Labs! Lab: Fun with Functions

  35. # dir(str) shows methods on str object # a string representation of a number >>> x = 3.25 >>> 'number is' + repr(x) 'number is3.25' # pad with zeros >>> '12'.zfill(5) '00012' # explicit formatting (Python 2.6) >>> 'The value of {0} is approximately {1:.3f}.'.format('PI', math.pi) The value of PI is approximately 3.142. Input/Output

  36. # file objects need to be opened # some modes - 'w' (write), 'r' (read), 'a' (append) # - 'r+' (read+write), 'rb', (read binary) >>> f = open('/tmp/workfile', 'w') >>> print f <open file '/tmp/workfile', mode 'w' at 80a0960> >>> help(f) >>> f.write('I want my binky!') >>> f.close() >>> f = open('/tmp/workfile', 'r+') >>> f.readline() 'I want my binky!' File I/O

  37. # file substitute.py import re fin = open('fadd.f90', 'r') p = re.compile('(subroutine)') try: while True: s = fin.readline() if s == "": break sout = p.sub('SUBROUTINE', s) print sout.replace('\n', "") # sys.stdout.write simpler except: print "Finished reading, file" # is this line reached? fin.close() Search and Replace

  38. class fibnum: def __init__(self): self.fn1 = 1 # f [n-1] self.fn2 = 1 # f [n-2] def next(self): # next() is the heart of any iterator oldfn2 = self.fn2 self.fn2 = self.fn1 self.fn1 = self.fn1 + oldfn2 return oldfn2 def __iter__(self): return self Iterators over Containers Interators require two methods: next() and __iter__() Fibonacci: f[n] = f[n-1] + f[n-2]; with f[0] = f[1] = 1

  39. # use Fibonacci iterator class >>> fromfibnum import * # construct a member of the class >>> f = fibnum() >>> l = [] >>> for i in f: l.append(i) if i > 20: break >>> l = [] [1, 1, 2, 3, 5, 8, 13, 21] # thanks to (and for more information on iterators): # http://heather.cs.ucdavis.edu/~matloff/Python/PyIterGen.pdf Iterators…

  40. Binary I/O Anticipating the next module NumPy (numerical arrays), you may want to look at the file PVReadBin.py to see how binary I/O is done in a practical application.

  41. Labs! Lab: Input/Output Try out file I/O as time allows

  42. SECTION 4 NUMERICAL PYTHON

  43. NumPy • Offers Matlab like capabilities within Python • Information • http://numpy.scipy.org/ • Download • http://sourceforge.net/projects/numpy/files/ • Numeric developers (initial coding Jim Hugunin) • Paul Dubouis • Travis Oliphant • Konrad Hinsen • Charles Waldman

  44. >>> from numpy import * >>> a = array([1.1, 2.2, 3.3]) >>> print a [ 1.1 2.2 3.3] # two-dimension array >>> b = array(([1,2,3],[4,5,6])) >>> print b [[1 2 3] [4 5 6]] >>> b.shape (2,3) >>> print ones((2,3), float) [[1. 1. 1.] [1. 1. 1.]] >>> print resize(b,(2,6)) [[1 2 3 4 5 6] [1 2 3 4 5 6]] >>> print reshape(b,(3,2)) [[1 2] [3 4] [5 6]] Creating Array: Basics

  45. # user reshape with range >>> a = reshape(range(12),(2,6)) >>> print a [[0 1 2 3 4 5] [6 7 8 9 10 11]] # set an entire row (or column) >>> a[0,:] = range(1,12,2) >>> print a [[1 3 5 7 9 11] [6 7 8 9 10 11]] >>> a = zeros([50,100]) # loop to set individual values >>> for i in range(50): … for j in range(100): … a[i,j] = i + j # call user function set(x,y) >>> shape = (50,100) >>> a = fromfunction(set, shape) # use scipy.io module to read # values from a file into an # array Creating Arrays: Strategies

  46. >>> a = arange(1,4); print a [1 2 3] # addition (element wise) >>> print 3 + a [4 5 6] # multiplication (element wise) >>> print 3*a [3 6 9] # it really is element wise >>> print a*a [1 4 9] # power: a**b -> power(a,b) >>> print a**a [1 4 27] # functions: sin(x), log(x), … >>> print sqrt(a*a) [1. 2. 3.] # comparison: ==, >, and, … >>> print a < a [False False False] # reductions >>> add.reduce(a) 6 Simple Array Operations

  47. >>> a = reshape(range(9),(3,3)) >>> print a [[0 1 2] [3 4 5] [6 7 8]] # second column >>> print a[:,1] [1 4 7] # last row >>> print a[-1,:] [6 7 8] # slices are references to # original memory, true for # all array/sequence assignment # work on the first row of a >>> b = a[0,:] >>> b[0] = 99 ; print b [99 1 2] # what is a[0,:] now? >>> print a[0,:] [99 1 2] Slicing Arrays

  48. >>> a = arange(10) >>> b = arange(10,20) # What will the following do? >>> a = a + b # Is the following different? >>> c = a + b >>> a = c # Does “a” reference old or new # memory? Answer, new memory! # Watch out for array # temporaries with large arrays! # Universal functions, ufuncs >>> type(add) <type 'numpy.ufunc'> # add is a binary operator >>> a = add(a,b) # in place operation >>> add(a,b,a) Array Temporaries and ufuncs

  49. >>> a = arange(1,11); print a [1 2 3 4 5 6 7 8 9 10] # create an index array >>> ind = [0, 5, 8] # take values from the array >>> print take(a,ind) >>> print a[ind] [1 6 9] # put values to the array >>> put(a,ind,[0,0,0]); print a >>> a[ind] = (0,0,0); print a [0 2 3 4 5 0 7 8 0 10] >>> a = reshape(range(9),(3,3)) >>> b = transpose(a); print b [[0 3 6] [1 4 7] [2 5 8]] >>> print diagonal(b) [0 4 8] >>> print trace(b) 12 >>> print where(b >= 3, 9, 0) [[0 9 9] [0 9 9] [0 9 9]] Array Functions

  50. Labs! Lab: NumPy Basics

More Related