Euler’s Method

Introduction to Data Science Ordinary Differential Equations Euler’s Method

In this section we’ll take a brief look at a fairly simple method for approximating solutions to differential equations. We derive the formulas used by Euler’s Method and give a brief discussion of the errors in the approximations of the solutions.

Let’s recall the very basic aspect of Differential Equation.

We know, $\frac{\text{d}y(x)}{\text{d}x}=f(x,y)$

It can be expressed as, $\frac{\text{d}y(x)}{\text{d}x}\sim \frac {y(x+h)-y(x)}{h}$

Combining both the equations, we get, $\frac{\text{d}y(x)}{\text{d}x}\sim \frac {y(x+h)-y(x)}{h}=f(x,y)$

Summarizing the equation as, $y(x+h)= y(x)+h.f(x,y)$

Thus, we arrive at, $y_{n+1}= y_{n}+h.f(x_{n},y_{n})$ .

where $h$ is the step size (the smaller, the better), $f(x,y)$ is right side of the differential equation.

Example 1: $\frac{\text{d}y}{\text{d}x}=2y$ with the initial condition as $y(0)=1$ .

Here, we have two approaches to find the solution of the Differential Equation, which are as follows:

1) Analytical Solution: This is the exact solution of an ordinary differential equation.

Integrating both sides, $\int_{}{}\frac{\text{d}y}{ y}=\int_{}{}{d}x$

We get, $\ln y = x + \ln C$

Therefore, $y(x) = e^{x}$

2) Numerical Solution: This is just an approximation of the solution.

We require this formula, $y(x+h)= y(x)+h.f(x,y)$

along with initial conditions, $x_{0}=0$ and $y_{0}=1$

Let’s continue with Example 1, to find the Numerical Solution.

Example 1: $\frac{\text{d}y}{\text{d}x}=2y$ with the initial condition as $y(0)=1$ and step-size $h=2$ .

We know the formula, $y_{n+1}= y_{n}+h.f(x_{n},y_{n})$ .

In this equation we put the values of $x_{n}$ at interval of $2$ because of Step-Size ( $h=2$ ) and obtain the value of $y_{n}.$ By initial conditions, we know $x_{0}=0$ and $y_{0}=1$ .

Here, we are assuming the slope to be same between two consecutive points of $x$ and $y$ , whereas in the exact case the slope is changing constantly.

Example 2: Find approximate value of ‘ $y$ ‘ at $x=1$ in five steps by taking h=0.2 for the DE $\frac{\text{d}y}{\text{d}x}= x+y$ where $y(0)=1$ .

SOLUTION: $\frac{\text{d}y}{\text{d}x}=x+y$ and $y(0)=1$

Hence, $f(x,y)=x+y$ and $x_{0} = 0$ , $y_{0} = 1$

We know the formula, $y_{n+1}= y_{n}+h.f(x_{n},y_{n})$

Therefore, equation can be modified as $y_{n+1}= y_{n}+h.(x_{n}+y_{n})$

Calculations required is put into a tabular form as shown below in Table 2.

The approximate value of $y$ at $x=1$ , i.e $y(1)= 2.976664$ .

Code:

				
					// Euler's Method
// Solving dy/dx = x + y by Euler's method

%use s2

// Defining function f(x,y) for this case
val f: BivariateRealFunction = object : AbstractBivariateRealFunction() {
    override fun evaluate(x: Double, y: Double): Double {
        return x+y
    }
}

val x_at_0 = 0.0
val y_at_0 = 1.0
val h = 0.2

// Evaluation of h*f(x(n),y(n))
val e_term: BivariateRealFunction = object : AbstractBivariateRealFunction() {
    override fun evaluate(h: Double, f_at_n: Double): Double {
        return h*f_at_n
    }
}

// y(n+1) = y(n) + h*f(x(n),y(n))      ..........Euler's Method 
val Y_at_n_plus_1: BivariateRealFunction = object : AbstractBivariateRealFunction() {
    override fun evaluate(y_at_n: Double, e_term: Double): Double {
        return y_at_n + e_term
    }
}

val h_term: BivariateRealFunction = object : AbstractBivariateRealFunction() {
    override fun evaluate(h: Double, t: Double): Double {
        return h*t
    }
}

val xn: BivariateRealFunction = object : AbstractBivariateRealFunction() {
    override fun evaluate(x_at_0: Double, h_term: Double): Double {
        return x_at_0 + h_term
    }
}

println("Solution of Differential Equation by Euler's Method")
val steps = mutableListOf(0.0, 1.0, 2.0, 3.0, 4.0)    //no. of steps = 5
var y_n = y_at_0
for (index in steps){
    val h_n = h_term.evaluate(h,index) 
    val x_n = xn.evaluate(x_at_0,h_n)
    val f_n = f.evaluate(x_n,y_n)
    val e_n = e_term.evaluate(h,f_n)
    val y = Y_at_n_plus_1.evaluate(y_n,e_n)
    y_n = y
    println("Value of y at n=%f: %f".format(index,y))
}
println("The approximate value of y at x=1 is %f".format(y_n))

Output: