-
Notifications
You must be signed in to change notification settings - Fork 190
/
diff_1st2nd_order.py
39 lines (35 loc) · 1.18 KB
/
diff_1st2nd_order.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
def first_order(f, x, h):
return (f(x+h) - f(x))/h
def second_order(f, x, h):
return (f(x+h) - f(x-h))/(2*h)
def derivative_on_mesh(formula, f, a, b, n):
"""
Differentiate f(x) at all internal points in a mesh
on [a,b] with n+1 equally spaced points.
The differentiation formula is given by formula(f, x, h).
"""
h = (b-a)/float(n)
x = linspace(a, b, n+1)
df = zeros(len(x))
for i in xrange(1, len(x)-1):
df[i] = formula(f, x[i], h)
# Return x and values at internal points only
return x[1:-1], df[1:-1]
def example(n):
a = 0; b = pi;
x, dF = derivative_on_mesh(first_order, sin, a, b, n)
x, dS = derivative_on_mesh(second_order, sin, a, b, n)
# Accurate plot of the exact derivative at internal points
h = (b-a)/float(n)
xfine = linspace(a+h, b-h, 1001)
exact = cos(xfine)
plot(x, dF, 'r-', x, dS, 'b-', xfine, exact, 'y-',
legend=('First-order derivative',
'Second-order derivative',
'Correct function'),
title='Approximate and correct discrete '\
'functions, n=%d' % n)
# Main program
from scitools.std import *
n = int(sys.argv[1])
example(n)