Chemical & Biological Engineering CBE Computational Turbulent Reacting Flow CTRFlow

Publications Calendar Research Teaching Zoom LDS

Publications Calendar Research Teaching Zoom LDS

Python Arrays

Import Libraries

import numpy as np
import matplotlib.pyplot as plt
%matplotlib inline

Array basics were covered earlier

• Review those if needed
• See this page for more information on built-in Python Lists.

Numpy arrays

• Python lists are nice, but for numerical work, numpy arrays are much better.
• 90% of the work we will do in science, math, and engineering, is best done using numpy arrays.
• All of the above regarding Python lists is the same, except
  • Each element type needs to be the same.
  • The way we create arrays is different.

Some good links:

• Numpy for Matlab users
  • This too
• Tentative numpy tutorial
• Numpy quickstart tutorial
• Linear algebra
• Cookbook with good recipes

Length of a list

len function

a = np.array([1.1, 2.2, 3.3, 4.4])

Creating arrays

• array (we've already seen this)
• linspace, arange
• empty, zeros, ones, full

# linspace

# arange

# empty

# zeros

# ones

# full

Access elements

• use square brackets with index number
• array indexing starts at zero
• for 2D arrays (matrix) x[i,j] is row i, column j
• can access backwards0
  • x[-1] is the last element
  • x[-2] is the second to the last element, etc.
  • See the figure below. If we wrap the array into a circle. Then element 0 is at the top. We index by +1 in moving clockwise, and index by -1 in moving counter-clockwise.

x = np.arange(3,21,3)               # x is         [3 6 9 12 15 18]
print(x)

# Access elements with index values

# or use an "index array"

Array slices

• Can access multiple elements using syntax: a[s:e]
  • s is the start (optional)
  • e is the end (optional)
  • (stops short of the end)

x = np.arange(0,11,1.5)       

                                   # get all elements
                                   # get indices 2 to the end
                                   # get elements start to one before 2nd to last (stops short)
                                   # get indices 2 to 4 (stops short of 5)

Example

We have an array of times $t$, and we want the array of differences $dt$ between times.

• How big is $t$?
• How big will $dt$ be?
• Create $dt$ using a loop.
• Create $dt$ using slice-indexing.
  • Note, for arrays a and b, a-b is done element-wise.

t = np.array([ 
    0.01136646,  0.01655847,  0.02548247,  0.07698612,  0.17964401,
    0.19813924,  0.22419994,  0.28622096,  0.28884997,  0.3627257 ,
    0.42555072,  0.44965767,  0.47257237,  0.59335946,  0.73842352,
    0.85900369,  0.88851226,  0.89153229,  0.89355469,  0.97705246,
    ])

# loop

# slice

Array properties

• size and shape

x = np.array([[1,2,3,4],[7,8,9,10]])
print(x)

np.shape(x)            # rows x columns

np.size(x)             # total number of elements

len(x)                 # length of first dimension (rows)

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

Arithmetic

a = np.arange(9)
a

#------- THIS

b = 0.1 * a                 # multiplies each element of a by 0.1 and assigns whole resulting array to b

#------- NOT THIS

b = np.empty(len(a))
for i in range(len(b)):
    b[i] = 0.1*a[i]
    
#-------

print(b)

a+b                 # adds element-wise

a**2                # squares element-wise

BE CAREFUL about assignment

x = np.array([1,2,3,4]) 
y = x                    # note, y is now just a reference to x.
y[3] = 44                # if we change an element of y, it will change in x too!
print(x)                 # this is not usually what we want
print(y)

x = np.array([1,2,3,4])
y = x.copy()                 # now x and y are two different objects
y[3] = 44                    #   could also do y = np.array(x)
print(x)
print(y)

Basic Plots

x = np.linspace(0,20,1000)
y = np.exp(-x/5)*np.sin(x)

plt.plot(x,y);

Problem 1

• Make an array of $x$ from 0 to 5 with 1000 points.
• create the following arrays
  • $2x$
  • $x^2 -3x + 4$
  • $\sin(5x)$
• Plot each array versus $x$ on a single plot.

import numpy as np
import matplotlib.pyplot as plt
%matplotlib inline

Problem 2

• Ideal gas law $$\rho = \frac{MP}{RT}$$
• M=29 kg/kmol
• R=8314 J/kmol*K
• P=101325 Pa
• Find $\rho$ for $T=300,400,\ldots1000$ K
  • Store $\rho$ and $T$ as arrays
  • Plot $\rho$ versus $T$

Problem 3

• Redo Problem 2 by writing a function for $\rho(T)$ as you would have in a previous lecture.
  • Call the function with an array of $T$ and store the result as an array of $\rho$.
  • Plot $\rho$ versus $T$.

Problem 4

• You have an array of $x$ points $x=0.0, 0.1, \ldots, 50$.
• Create $x$
• Create $y=\sin(x)$
• At interior points, use a central difference approximation to find the derivative $$\left(\frac{dy}{dx}\right)_i \approx \frac{y_{i+1}-y_{i-1}}{2\Delta x},$$ where $\Delta x = x_1-x_0$.
• At edge point use a one-sided difference approximation, $$\left(\frac{dy}{dx}\right)_0 \approx \frac{y_{1}-y_{0}}{\Delta x},$$ $$\left(\frac{dy}{dx}\right)_{n-1} \approx \frac{y_{n-1}-y_{n-2}}{\Delta x}.$$
• Plot $dy/dx$ versus $x$ , and $\cos(x)$ versus $x$.

# y: * (*  *  *  *  *  *) *        # boundary, (interior), boundary
#    0  1  2  3  4  5  6  7
#                     n-2 n-1

Problem 5

• Solve the rate equation below to $t=5$ s using the Explicit Euler method $y_{i+1} = y_i + \Delta tf(y_i,t_i).$

$$\frac{dy}{dt} = -ky,$$ $$y_0=1,$$ $$k=2.$$

• Use $\Delta t = 0.5$.
• Store $t$ and $y$ in arrays.
• Plot $y$ versus $t$, and $y_0e^{-kt}$ versus $t$ on the same plot.

Problem 6

• The Haaland equation relates the friction factor $f$ in turbulent pipe flow to the Reynolds number $Re$. $$\frac{1}{\sqrt{f}} = -1.8\log_{10}\left[\left(\frac{\epsilon/D}{3.7}\right)^{1.11}+ \frac{6.9}{Re}\right].$$

• Write a function to compute $f(Re, \epsilon/D)$.

• Create an array of 100 points that is uniformly spaced on a log scale, for $1000\le Re\le 1\times 10^8$

• Plot $f$ versus $\log_{10}(Re)$.

• Use $\epsilon/D=0.001$.