Computational Science Lecture Notes: Derivatives, Integrals, Root Finding, and Python Code, Assignments of Physics

Download Computational Science Lecture Notes: Derivatives, Integrals, Root Finding, and Python Code and more Assignments Physics in PDF only on Docsity!

Lec1 - intro, basic tools

• Mechanics, syllabus
• General comments
• What is computational science?
• Derivatives, integrals and root finding
• Intro to Python

General comments

• Not programming course. Not traditional physics course.
• Topics drawn from many different areas of physics including recent developments eg. self-organized critical phenomena, black- hole orbits, fractals
• Brief discussion of programming issues - only as much as needed to do the science
• Python/VPython makes life easy. Don’t need to know much to do some quite com- plicated stuff (graphics built-in)

Examples

Type 1: throwing a ball into the air - time to reach ground? a)Newton’s laws, analytic solution (neglect air resistance, model as point particle) b)Real world: put in air resistance, spin of ball, can integrate the equations numerically with high accuracy.

Type 2: Show that water freezes when cooled. Features:

• Many d.o.f, complex interactions
• Qualitative change of state of system.

Errors: finite simulation time. Finite number of particles.

Derivatives

Want a numerical approx to (^) dxdf Why? Perhaps derivative too hard to compute analytically. Or perhaps f(x) produced by some other procedure at only a finite set of points. How? Taylor expand:

f (x + h) = f (x) + f 1 (x)h + f 2 (x)h^2 /2 +...

Many possibilities:

f 1 ∼ f (x + h) − f (x) h

• O(h)

or

f 1 ∼ f (x + h) − f (x − h) 2 h

• O(h^2 )

Improve?

f 1 ∼ f− 2 − f 2 − 8 f− 1 + 8f 1 12 h

• O(h^4 )

Python code

eg. df /dx for f = x^2 at x = 1

from math import * x=1. h=0. print (exp(x+h)-exp(x-h))/(2*h)

Improvements

Rather than type in the code from scratch each time I want to change h,x can create function

from math import * def deriv(x,h): print (exp(x+h)-exp(x-h))/(2*h)

Execute: eg deriv(1,0.1), deriv(1,0.01) etc

Better

Use quadratic approx to f in small interval

f (x) = f 0 + (f 1 − f− 1 ) 2 h

x+ f 1 − 2 f 0 + f− 1 2 h^2

x^2 +... ∫ (^) h −h

f (x)dx ∼ h 3

(f 1 + 4f 0 + f− 1 )

Leading to Simpson’s rule:

∫ (^) b

a

f dx = h 3

(f (a) + 4f (a + h) + 2f (a + 2h) +... +... 4 f (b − h) + f (b))

Python code

Eg. integrate ex^ over 0 → 1

def int(a,b,h): x=a integral=exp(a)h/2. while x<(b-h): x=x+h integral=integral+hexp(x)

integral=integral+exp(b)*h/2. print integral

• Note indentation - very important in Python.
• Note while loop.
• call: int(0,1,0.1)
• Place in separate file - called say int.py

Python code

Eg. solve ex^ − 0 .5 = 0 in range − 1 → 0

def root(a,b,tol): l=a u=b fl=exp(a)-0. fu=exp(b)-0. while (u-l) > tol: m=(u+l)/2. fm=exp(m)-0. if(fm*fu>0.0): u=m fu=fm else: l=m fl=fm

print m

Comments

• Note if/else structure
• if statement nested within while loop. In- dentation tells you which bits of code are associated with which control structures
• call with: root(0,1,0.01)

Summary

• General comments.
• Workhorse alg. for derivatives, integrals and roots
• Simple Python codes for doing this
• Simple assignment, arithmetic, I/O, while loops and if statements. Indentation
• Function definitions