Euler method

Method description

Euler's method numerically approximates solutions of first-order ordinary differential equations (ODEs) with a given initial value. It is a explicit method for solving initial value problems (IVPs), as described in Euler method

The ODE has to be provided in the folloving form:

{\frac {dy(t)}{dt}}=f(t,y(t))

with an initial value

y(t_{0})=y_{0}

To get a numeric solution, we replace the derivative on the LHS with a finite difference approximation:

{\frac {dy(t)}{dt}}\approx {\frac {y(t+h)-y(t)}{h}}

then solve for $y(t+h)$ :

y(t+h)\approx y(t)+h\,{\frac {dy(t)}{dt}}

which is the same as

y(t+h)\approx y(t)+h\,f(t,y(t))

The iterative solution rule is then:

y_{n+1}=y_{n}+h\,f(t_{n},y_{n})

$h$ is the step size, the most relevant parameter for accuracy of the solution. A smaller step size increases accuracy but also the computation cost, so it has always has to be hand-picked according to the problem at hand.

Example: Newton's Cooling Law

Newton's cooling law describes how an object of initial temperature $T(t_{0})=T_{0}$ cools down in an environment of temperature $T_{R}$ :

{\frac {dT(t)}{dt}}=-k\,\Delta T

or

{\frac {dT(t)}{dt}}=-k\,(T(t)-T_{R})

It says that the cooling rate ${\frac {dT(t)}{dt}}$ of the object is proportional to the current temperature difference $\Delta T=(T(t)-T_{R})$ to the surrounding environment.

The analytical solution, which we will compare to the numerical approximation, is

T(t)=T_{R}+(T_{0}-T_{R})\;e^{-kt}

Task

The task is to implement a routine of Euler's method and then to use it to solve the given example of Newton's cooling law with it for three different step sizes of 2 s, 5 s and 10 s and to compare with the analytical solution.

The initial temperature $T_{0}$ shall be 100 °C, the room temperature $T_{R}$ 20 °C, and the cooling constant $k$ 0.07. The time interval to calculate shall be from 0 s to 100 s.

An reference solution (Common Lisp) can be seen below. We see that bigger step sizes lead to reduced approximation accuracy.

Common Lisp

<lang lisp> (defconstant true 'cl:t) (shadow 't)

Approximates y(t) in y'(t)=f(t,y) with y(a)=y0 and t=a..b and the step size h.

(defun euler (f y0 a b h)

 (do ((t a  (incf t h))
      (y y0 (incf y (* h (funcall f t y)))))
     ((>= t b)
      'DONE)
     (format true "~6,3F  ~6,3F~%" t y)))

Example

Newton's cooling law, f(t,T) = -0.07*(T-20)

(defun newton-cooling (time T) (* -0.07 (- T 20)))

Generate data.

(euler #'newton-cooling 100 0 100 2) (euler #'newton-cooling 100 0 100 5) (euler #'newton-cooling 100 0 100 10)

Example output

0.000  100.000

10.000 44.000 20.000 27.200 30.000 22.160 40.000 20.648 50.000 20.194 60.000 20.058 70.000 20.017 80.000 20.005 90.000 20.002 DONE </lang> --Avi 12:33, 6 March 2011 (UTC)