Euler's method in python Euler's method in python numpy numpy

Euler's method in python


The formula you are trying to use is not Euler's method, but rather the exact value of e as n approaches infinity wiki,

$n = \lim_{n\to\infty} (1 + \frac{1}{n})^n$

Euler's method is used to solve first order differential equations.

Here are two guides that show how to implement Euler's method to solve a simple test function: beginner's guide and numerical ODE guide.

To answer the title of this post, rather than the question you are asking, I've used Euler's method to solve usual exponential decay:

$\frac{dN}{dt} = -\lambda N$

Which has the solution,

$N(t) = N_0 e^{-\lambda t}$

Code:

import numpy as npimport matplotlib.pyplot as pltfrom __future__ import division# Concentration over timeN = lambda t: N0 * np.exp(-k * t)# dN/dtdef dx_dt(x):    return -k * xk = .5h = 0.001N0 = 100.t = np.arange(0, 10, h)y = np.zeros(len(t))y[0] = N0for i in range(1, len(t)):    # Euler's method    y[i] = y[i-1] + dx_dt(y[i-1]) * hmax_error = abs(y-N(t)).max()print 'Max difference between the exact solution and Euler's approximation with step size h=0.001:'print '{0:.15}'.format(max_error)

Output:

Max difference between the exact solution and Euler's approximation with step size h=0.001:0.00919890254720457

Note: I'm not sure how to get LaTeX displaying properly.


Are you sure you are not trying to implement the Newton's method? Because Newton's method is used to approximate the roots.

In case you decide to go with Newton's method, here is a slightly changed version of your code that approximates the square-root of 2. You can change f(x) and fp(x) with the function and its derivative you use in your approximation to the thing you want.

import numpy as npdef f(x):    return x**2 - 2def fp(x):    return 2*xdef Newton(f, y0, N):    y = np.zeros(N+1)    y[0] = y0    for n in range(N):        y[n+1] = y[n] - f(y[n])/fp(y[n])    return yprint Newton(f, 1, 10)

gives

[ 1. 1.5 1.41666667 1.41421569 1.41421356 1.41421356 1.41421356 1.41421356 1.41421356 1.41421356 1.41421356]

which are the initial value and the first ten iterations to the square-root of two.

Besides this a big problem was the usage of ^ instead of ** for powers which is a legal but a totally different (bitwise) operation in python.