Type and Type conversion

>>> print(type(633))

>>> print(type("633"))

>>> print(type(633.0))

mon_sales = 121
tues_sales = 105
wed_sales = 110
thurs_sales = 98
fri_sales = 95
total = str(mon_sales + tues_sales + wed_sales + thurs_sales + fri_sales)
print("This week\'s total sales: " + total)

python string method
https://docs.python.org/3/library/stdtypes.html#string-methods

prophecy = "And there shall in that time be rumours of things going astray, and there will be a great confusion as to where things really are, and nobody will really know where lieth those little things with the sort of raffia work base, that has an attachment…at this time, a friend shall lose his friends’s hammer and the young shall not know where lieth the things possessed by their fathers that their fathers put there only just the night before around eight o’clock…"

vowel_count = 0
vowel_count += prophecy.count('a')
vowel_count += prophecy.count('A')
vowel_count += prophecy.count('e')
vowel_count += prophecy.count('E')
vowel_count += prophecy.count('i')
vowel_count += prophecy.count('I')
vowel_count += prophecy.count('o')
vowel_count += prophecy.count('O')
vowel_count += prophecy.count('u')
vowel_count += prophecy.count('U')

print(vowel_count)

format

log_message = "IP address {} accessed {} at {}".format(user_ip, url, now)

city = "Seoul"
high_temperature = 18
low_temperature = 9
temperature_unit = "degrees Celsius"
weather_conditions = "light rain"

alert = "Today's forecast for {}: The temperature will range from {} to {}{}. Conditions will be {}.".format(city, high_temperature, low_temperature, temperature_unit, weather_conditions)
# print the alert
print(alert)

Changing Variable

>>> manila_pop = 1780148
>>> manila_area = 16.56
>>> manila_pop_density = manila_pop/manila_area
>>> print(int(manila_pop_density))
107496

bool

>>> 1 < 2
True
>>> 42 > 43
False
sf_population, sf_area = 864816, 231.89
rio_population, rio_area = 6453682, 486.5

san_francisco_pop_density = sf_population/sf_area
rio_de_janeiro_pop_density = rio_population/rio_area

print(san_francisco_pop_density > rio_de_janeiro_pop_density)

string

>>> print("hello")
hello
>>> print('hello')
hello
>>> welcome_message = "Hello, welcome to Tokyo"
>>> print(welcome_message)
Hello, welcome to Tokyo
>>> instructor_1 = "Philip"
>>> instructor_2 = "Charlie"
>>> print(instructor_1 + " and " + instructor_2)
Philip and Charlie
messiah = 'He\'s not the Messiah, he\'s a very naughty boy!'
print(messiah)

username = "Kinari"
timestamp = "04:50"
url = "http://petshop.com/pets/mammals/cats"
print(username + " accessed the site " + url + " at " + timestamp)

len

>>> tokyo_length = len("Tokyo")
>>> print(tokyo_length)
5

given_name = "Charlotte"
middle_names = "Hippopotamus"
family_name = "Turner"

name_length = len(given_name+middle_names+family_name)

driving_licence_character_limit = 28
print(name_length <= driving_licence_character_limit) [/python]

Arithmetic

Python 2.7.12 (v2.7.12:d33e0cf91556, Jun 27 2016, 15:19:22) [MSC v.1500 32 bit (Intel)] on win32
Type "copyright", "credits" or "license()" for more information.
>>> print(3+1)
4
>>> 3 + 1
4
>>> print(1 + 2 + 3 * 3)
12
>>> print((1 + 2 + 3) * 3)
18
>>> print(3**2)
9
>>> print(9 % 2)
1
>>> print(15 // 4)
3
>>> print(16 // 4)
4
>>> print(-5 // 4)
-2
>>> print((23 + 32 + 64)/ 3)
39

integer and floats

>>> print(3/4)
0
>>> print(16/4)
4
>>> 3 + 2.5
5.5
>>> 213.13
213.13
>>> 341.
341.0
>>> int(49.7)
49
>>> int(16/4)
4
>>> float(3520+3239)
6759.0
>>> print(0.1)
0.1
>>> print(0.1 + 0.1 + 0.1)
0.3

float e.g.
Length of a fish you caught, in metres
Length of time it took to catch the first fish, in hours

4.445e8 is equal to 4.445 * 10 ** 8 which is equal to 444500000.0.

# The current volume of a water reservoir (in cubic metres)
reservoir_volume = 4.445e8
# The amount of rainfall from a storm (in cubic metres)
rainfall = 5e6

# decrease the rainfall variable by 10% to account for runoff
rainfall *= 0.9
# add the rainfall variable to the reservoir_volume variable
reservoir_volume += rainfall
# increase reservoir_volume by 5% to account for stormwater that flows
reservoir_volume *= 1.05 
# into the reservoir in the days following the storm

# decrease reservoir_volume by 5% to account for evaporation
reservoir_volume *= 0.95
# subtract 2.5e5 cubic metres from reservoir_volume to account for water
# that's piped to arid regions.
reservoir_volume -= 2.5e5
# print the new value of the reservoir_volume variable
print(reservoir_volume)

Antilock Braking System

Two types of friction
static: moves when F > μs・N
kinetic: F = μk * N

wheel slip
rolling, rocked

speed = 120 km/h
μ=1.0: 3.4s
μ=0.7: 4.9s

F = m*a
u*gravitational force(mg)
a = μ*g

t/ 120km/h = 1 / μg
t = 120km/h / μ*g

Computing the coefficient of Friction
s = 1 – w/v
F(s) = μ(s)*mqcG
V= -F(s)/mqc
F = ma
w = F(s)/mew – B

import math
from **** import *

h = 0.01
mass_quarter_car = 250.
mass_effective_wheel = 20.
g = 9.81

end_time = 5.
num_steps = int(end_time / h)

w = numpy.zeros(num_steps + 1)
v = numpy.zeros(num_steps + 1)
x = numpy.zeros(num_steps + 1)
times = h * numpy.array(range(num_steps + 1))

@show_plot(7, 7)
def plot_me():
	axes_x = matplotlib.pyplot.subplot(411)
	axes_v = matplotlib.pyplot.subplot(412)
	axes_w = matplotlib.pyplot.subplot(413)
	awes_s = matplotlib.pyplot.subplot(414)

	def friction_coeff(slip):
		return 1.1 * (1. - math.exp(-20. * slip)) - 0.4 * slip

	def wheel_slip():
		b_values = numpy.arange(70., 190.1, 30.)
		for b in b_values:
			x[0] = 0.

			for step in range(num_steps):
				if v[step] < 0.01:
					break
				s = max(0., 1. - w[step] / v[step])

				w[step + 1] = max(0., w[step + 1])

			axes_x.plot(times[:step], x[:step])
			axes_v.plot(times[:step], v[:step])
			axes_w.plot(times[:step], w[:step])
			axes_s.plot(times[:step], 1. - w[:step] / v[:step])
			p = int((0.35 + 0.4 * (b - b_values[0])/ (b_values[-1] - b_values[0])) * num_steps)
			axes_x.annotate(b, (times[p], x[p]),
					xytext = (-30, -30), textcoords = 'offset point',
					arrowprops = dict(arrowstyle = '-', connectionstyle = 'arc3, '))

			p = int((0.35 + 0.4 * (v - b_values[0]) / (b_values[-1] - b_values[0])) * num_steps)
			axes_x.annotate(b, (times[p], x[p]),
				xytext = (-30, -30), textcoords - 'offset points',
				arrowprops = dict(arrowstyle = '-', connectionstyle = 'arc3, rad = 0.2', shrinkB = 0.))

		return x, v, w

    axes_x.set_ylabel('Position\nin m', multialignment = 'center')
    axes_v.set_ylabel('Car velocity\nin m/s', multialignment = 'center')
    axes_w.set_ylabel('Wheel velocity\nin m/s', multialignment = 'center')
    axes_s.set_ylabel('Wheel\nslip', multialignment = 'center')
    axes_s.set_xlabel('Time in s')
    axes_s.set_ylim(0., 1.)

    return wheel_slip()

x, v, w = plot_me()

Infinity Time

x(t) = (x(t)^2 + 1)/2
with * (0) = 3

from *** import *

end_time = 42.
h = 0.0001
num_steps = int(end_time / h)
times = h * numpy.array(range(num_steps) + 1)

x = numpy.zeros(num_steps + 1)
x[0] = 3.

def forward_euler():
	for step in range(num_steps):

	return x

x = forward_euler()

@show_plot
def plot_me():
	matplotlib.pyplot.plot(times, x)
	matplotlib.pyplot.show()

plot_me()

Sources of Errors
-model, parameters, initial data, numerics

Logistic growth, Constant Harvesting

Exp growth, Amount of fish

from xxx import *

harvest_rates = [2e4, 5e4, 1e5, 2e5]

data = []

def logistic_growth():
	maximum_growth_rate = 0.5
	carrying_capacity = 2e6

	end_time = 10.
	h = 0.1
	num_steps = in(end_time / h)
	times = h * numpy.array(range(num_steps + 1))

	fish = numpy.zeros(num_steps + 1)
	fish[0] = 2e5

	for harvest_rate in harvest_rates:

		data.append(([time for time in times], [f for f in fish],
			str(harvest_rate)))

	return fish

fish = logistic_growth()

@show_plot
def plot_me():
    fish_plots = []
    for (times, fish, rate_label) in data:
        fish_plots.append(matplotlib.pyplot.plot(times, fish, label=rate_label))    
    matplotlib.pyplot.legend([str(rate) for rate in harvest_rates], loc='upper right')
    axes = matplotlib.pyplot.gca()
    axes.set_xlabel('Time in years')
    axes.set_ylabel('Amount of fish in tons')

plot_me()

Sir Backwards

R(t) = 1 / 5days I(t)

R2 = R1 + h 1/5daysI2
S2 = S1 + h(- 5 * 10^-9/ day*Person I1,S1)

Implicit solvers: use it or not?

Trapezoidal Stability
Cases – Moderate value of h, hk/2 = 0.1
hk/2 = 1
High value of h, e.g. hk/2 = g

Models with flows between different compartments

from xxx import *

h = 0.5
end_time = 60.
num_steps = int(end_time / h)
times = h * numpy.array(range(num_steps + 1))

def waning():
	transmission_coeff = 5e-9
	infectious_time = 5.

	waning_time

	s = numpy.zeros(num_steps + 1)
	i = numpy.zeros(num_steps + 1)
	r = numpy.zeros(num_steps + 1)

	s[0] = 1e8 - 1e6 - 1e5
	i[0] = 1e5
	r[0] = 1e6

	for step in range(num_steps):
		s2i = h * transmission_coeff * s[step] * i[step]
		i2r = h / infectious_time * i[step]

		s[step + 1] = s[step] - s2i
		i[step + 1] = i[step] + s2i - i2r
		r[step + 1] = r[step] + i2r

		return s, i, r

s, i, r = waning()

@show_plot
def plot_me():
    s_plot = matplotlib.pyplot.plot(times, s, label = 'S')
    i_plot = matplotlib.pyplot.plot(times, i, label = 'I')
    r_plot = matplotlib.pyplot.plot(times, r, label = 'R')
    matplotlib.pyplot.legend(('S', 'I', 'R'), loc = 'upper right')

    axes = matplotlib.pyplot.gca()
    axes.set_xlabel('Time in days')
    axes.set_ylabel('Number of persons')
    matplotlib.pyplot.xlim(xmin = 0.)
    matplotlib.pyplot.ylim(ymin = 0.)
plot_me()

Contagion

A Contagious Disease
infectious, suspicious, Recovered

SIR Model
I(t)
R(t)= 1/5days * I(t)

A model for infection
S(t):number of suspicious person
I(t):number of infectious person
I(t) = 5 * 10^-9 / day*person I(t)*S(t)

SIR Model
S(t)= -5*10^-9/day*person I(t)S(t)
I(t) = 5 * 10^-9 / day*person I(t)*S(t)
R(t)= 1/5days * I(t)
System, Explicit

from xxx import *

h = 0.5
transmission_coeff = 5e-9
latency_time = 1.
infectious_time = 5.

end_time = 60.0
num_steps = int(end_time / h)
times = h * numpy.array(range(num_steps + 1))

def seir_model():

	s = numpy.zeros(num_steps + 1)
	e = numpy.zeros(num_steps + 1)
	i = numpy.zeros(num_steps + 1)
	r = numpy.zeros(num_steps + 1)

	s[0] = 1e8 - 1e6 - 1e5
	e[0] = 0.
	i[0] = 1e5
	r[0] = 1e6

	for step in range(num_steps):

	return s, e, i, r

s, e, i, r = seir_model()

@show_plot
def plot_me():
    s_plot = matplotlib.pyplot.plot(times, s, label = 'S')
    e_plot = matplotlib.pyplot.plot(times, e, label = 'E')
    i_plot = matplotlib.pyplot.plot(times, i, label = 'I')
    r_plot = matplotlib.pyplot.plot(times, r, label = 'R')
    matplotlib.pyplot.legend(('S', 'E', 'I', 'R'), loc = 'upper right')
    
    axes = matplotlib.pyplot.gca()
    axes.set_xlabel('Time in days')
    axes.set_ylabel('Number of persons')
    matplotlib.pyplot.xlim(xmin = 0.)
    matplotlib.pyplot.ylim(ymin = 0.)

x(t)=-k*(t), k=const>0
x(h)=x(0)+h(-k*(h))

Errors on Different Time Scales

LTE, GTE
Long-term

Kinetic energy: 1/2ms ||Vs||^2
Potential enery: -G mems/||Xs||

def total_enery():
	num_steps = 20000
	x = numpy.zeros([num_steps + 1, 2])
	v = numpy.zeros([num_steps + 1, 2])
	enery = numpy.zeros(num_steps + 1)

	x[0, 0] = 15e6
	x[0, 1] = 1e6
	v[0, 0] = 2e3
	v[0, 1] = 4e3

	for step in range(num_steps):
		x[step + 1] = x[step] + h * v[step]
		v[step + 1] = v[step] + h * acceleration(x[step])

	return x, energy