Table of Contents

• 1  Numpy, Scipy and Matplotlib
• 2  NumPy
  • 2.1  Numpy arrays
    • 2.1.1  Memory layout - Python list vs numpy ndarray
    • 2.1.2  Creating arrays
  • 2.2  Numpy datatypes
  • 2.3  Working with numpy arrays
    • 2.3.1  numpy ufuncs
  • 2.4  Array indexing and slicing
    • 2.4.1  Array filtering
    • 2.4.2  Array IO
    • 2.4.3  Basic statistics
  • 2.5  Linear Algebra
    • • 2.5.0.1  Solving system of equations
    • 2.5.1  numpy - many more features

Numpy, Scipy and Matplotlib

There are numerous scientific programs, packages and libraries written in various programming languages:

Mathematica, Maple, Matlab, Root, R, Numerical Recipes, etc

For python the de-facto standard are the Numpy, Scipy and Matplotlib packages. They allow to combine the ease and flexibility of Python programming with efficient and powerful libraries, which provide state-of-the art functionality and fast execution also for large applications.

NumPy

Applications in Science and Physics often require very cpu-intensive operations, e.g:

• large equation systems
• matrices
• numerical integration
• fourier-transformation
• random number generation
• ...

Rich legacy of numerical libraries written e.g. in Fortran or C/C++ (e.g. Numerical Recipes). This is also used with numpy, which offers an interface to BLAS (Basic Linear Algebra Library) and ATLAS (Automatically Tuned Linear Algebra Software).

Moreover, numpy offers a very efficient Array-Datatype ndarray which provides a homogenous store for numerical data types and operations and functions for it. The backend is written in C and allows efficient vectorised or parallelized operations on large array structures.

import numpy as np

%%timeit
# std pythons lists
a = range(10000000)
b = range(10000000)
c = []
for i in range(len(a)):
  c.append(a[i] + b[i])

%%timeit
# same with numpy
import numpy as np
a = np.arange(10000000)
b = np.arange(10000000) 
c = a + b

Numpy arrays

The basis of numpy is the ndarray, this is a homogenous array especially for numerical data types.

• It supports many different kinds of integer, float and boolean data types
• but homogenous, all elements in array of same type
  • (not quite true, elements can also be Python objects and thereby any Python type, but that's exceptional case)
• not just a list, it contains additional information on shape
• many utility functions to reate, fill, store, extract, manipulate

# std python list very flexible
L = list(range(10)) # homogenous list
L

# but also arbitrary elements
L3 = [True, "2", 3.0, 4]
[type(item) for item in L3]

# numpy array from python list
np.array([1, 4, 2, 5, 3]) #  integer elements --> integer array

# numpy array from python list
np.array([1, 4.2, 2, 5, 3]) #  one float --> float array

# can also enforce data type
np.array([1, 4, 2, 5, 3],dtype='float32')

Memory layout - Python list vs numpy ndarray

• A python list is actually only a list of pointers/references to the actula objects
  • very flexible, but also inefficient when accessing large arrays
• in contrast, numpy ndarray contains directly all elements
  • efficient to process but needs to have specific data type

Creating arrays

# Create a length-10 integer array filled with zeros
np.zeros(10, dtype=int)

# Create a 3x5 floating-point array filled with ones
np.ones((3, 5), dtype=float)

# Create an array filled with a linear sequence
# Starting at 0, ending at 20, stepping by 2
# (this is similar to the built-in range() function)
np.arange(0, 20, 2)

# Create an array of five values evenly spaced between 0 and 1
np.linspace(0, 1, 11)

# Create an array of 10 values logarithmically spaced between 10**0 and 10**5
np.logspace(0,5,10)

# Create an array of 5 uniformly distributed
# random values between 0 and 1
np.random.random(50)

# Create a 3x3 array of uniformly distributed
# random values between 0 and 1
np.random.random((3, 7))

# Create an array of 10 normal distributed
# random values with mean 0 and std deviation 1
np.random.normal(0,1,10)

Numpy datatypes

In most cases Python/numpy automatically identifies suitable data-type, though sometimes it can be useful to explicitly specify the desired type.

Here the list of supported types:

Data type	Description
bool_	Boolean (True or False) stored as a byte
int_	Default integer type (same as C long; normally either int64 or int32)
intc	Identical to C int (normally int32 or int64)
intp	Integer used for indexing (same as C ssize_t; normally either int32 or int64)
int8	Byte (-128 to 127)
int16	Integer (-32768 to 32767)
int32	Integer (-2147483648 to 2147483647)
int64	Integer (-9223372036854775808 to 9223372036854775807)
uint8	Unsigned integer (0 to 255)
uint16	Unsigned integer (0 to 65535)
uint32	Unsigned integer (0 to 4294967295)
uint64	Unsigned integer (0 to 18446744073709551615)
float_	Shorthand for float64.
float16	Half precision float: sign bit, 5 bits exponent, 10 bits mantissa
float32	Single precision float: sign bit, 8 bits exponent, 23 bits mantissa
float64	Double precision float: sign bit, 11 bits exponent, 52 bits mantissa
complex_	Shorthand for complex128.
complex64	Complex number, represented by two 32-bit floats
complex128	Complex number, represented by two 64-bit floats

a=np.random.normal(0,1,10)
a.dtype   # data type available as field in ndarray

Working with numpy arrays

a=np.array([[1,2,3,4],[4,5,6,7],[9,10,11,12]]) 
print(a)

# many fields and functions in ndarray
[m for m in dir(a) if not m.startswith('__')]

a.shape # layout stored in shape

b=a.reshape(2,6) # can be changed

print(a.shape) # original unchanged
print(b.shape) 

b

a.T # transpose

-1*a # multiply with scalar

a+2.1 # add scalar

a+a # adding arrays element-by-element

a+b # error - requires same shape

a*a # multiplying array element-by-element

a.dot(a) # matrix multiplication -- error aligned shapes required

a.dot(a.reshape(4,3)) # matrix multiplication: 3x4 times 4x3 works ok

c=np.arange(4,8)
print(c, c.shape)

a.dot(c) # 3x4 matrix times 4-vec ok

numpy ufuncs

np.sqrt(a) # can apply std maths functions

a**0.5 # and operators

np.sin(a)

See here for detailed list of Numpy functions: https://scipy.github.io/old-wiki/pages/Numpy_Example_List_With_Doc.html

Array indexing and slicing

As for Python lists there are many sophisticated ways how to access or extract elements or sub-arrays

a=np.array([[1,2,3,4],[4,5,6,7],[9,10,11,12]])
a

a[1,2] # 2nd line, 3rd colun (index starts at 0)

a[1] # 2nd line

a[-1] #  last line

a[:2] # 1st and 2nd line

a[:,2] # 3rd column

a[:2,1:3] # ...

Array filtering

a=np.array(range(100))

np.where(a>50)

a>80 # returns boolean array if elements fullfill criteria

a[a>80] # sub-array fullfilling criteria

Array IO

#data = np.loadtxt('numbers.dat') # reads numbers from text file and stores in numpy array
# also works for url
data = np.loadtxt('http://www-static.etp.physik.uni-muenchen.de/kurs/Computing/python/source/numbers.dat')

print(data)

a=np.array(np.arange(1,10).reshape(3,3))
print('a=',a)

a.tofile('a.data')
b=np.fromfile('a.data', dtype=int) # watch out, shape info lost
print('b=',b)

Basic statistics

T=np.array([1.3,4.5,2.8,3.9])
print (T.mean(), T.std(), T.var()) # Mean, std deviation, Variance

T.min()

T=np.random.normal(10.,2.,1000) # 1000 normal distributed random numbers, mean 10, std 2
print (T.mean(), T.std(), T.var())

# correlation
T=np.array([1.3,4.5,2.8,3.9])
P=np.array([2.7,8.7,4.7,8.2])
print (np.corrcoef([T,P]))

Linear Algebra

a=np.arange(12).reshape(3,4)
print(a)

a.diagonal() # Diagonalwerte

# standard matrices
print(np.ones(3))
print(np.zeros(3))
print(np.eye(3))

Solving system of equations

import numpy as np
from numpy.linalg import inv
a=np.array([[3,1,5],[1,0,8],[2,1,4]])
inva=inv(a) # matrix inversion
print(a)
print(inva)
np.dot(a,inva) # should give unity matrix

from numpy.linalg import solve
print(a)
b=np.array([6,7,8])
x=solve(a,b) # solve equation a*x = b
print(x)
np.dot(a,x) # test: should give b

numpy - many more features

• advanced Linalg tools
• numerical methods, ...