Monday 13 January 2020

Newton Method

The Newton method (or Newton-Raphson Method)is an algorithm that allows one to find the roots of an equation when solving the equation (f(x)=0) is not an option.

The general form of the algorithm is here:
$${b_{n} = b_{n-1}- {{f(x)}\over{f'(x)}}}$$

$$end:f(x) < err$$

It is simply an arbitrary starting value less the function of the value over the derivative of the value.The idea is that the result of the first iteration feeds into the subsequent iteration, and the process terminates when f(x) < err threshold.
We start with parameter b0 and an error threshold which remains constant.

We want to find the x intercept of the following function
$$f(x)=-x^3+2x^2-e^x$$



We will need its derivative
$${{df(x)\over{dx}} =-3x^2+4x-e^x}$$

The form of the algorithm will look as follows:
$${b_{1} = b_{0}- {{-x^3+2x^2-e^x}\over{-3x^2+4x-e^x}}}$$

We arbitrarily set b0 as 10 and the error threshold to .001

I will print the values out instead of computing the algorithm with MathJax.

b0: -6.47058884010279  fx: 354.649469421832  fxdx: -151.489463487155
b0: -4.12950543233229  fx: 104.509232053146  fxdx: -67.6925579111767
b0: -2.58562475655285  fx: 30.5816407170309  fxdx: -30.474214161426
b0: -1.58209959409794  fx: 8.76059234819741  fxdx: -14.0430588411127
b0: -0.95826169330173  fx: 2.33291052253774  fxdx: -6.97140224182995
b0: -0.623621624975194 fx: 0.484337035554275 fxdx: -4.19719802152472
b0: -0.508226298600212 fx: 0.046298105996643 fxdx: -3.4093487299717
b0: -0.494646547843232 fx: 0.000591761687825487 fxdx: -3.32239821278166
-0.4946465

The value of x where f(x) = 0 is ~ -.495. This can be seen in the above graph where the blue line intercepts the x axis.

R Script
library(ggplot2)

#Function
f1_x <- function(x) {
  return (-x^3+2*(x^2)-exp(x))
}

#Derivative
f1_dxdy <- function(x) {
  return (-3*(x^2)+(4*x)-exp(x))
}

#Load data
dfFX <- data.frame(x=seq(-2.5,2.5,.01),f1_x=f1_x(seq(-2.5,2.5,.01)),
                         f1_dxdy=f1_dxdy(seq(-2.5,2.5,.01)))

#plot the data
ggplot(data=dfFX) + 
  geom_line(aes(x=x, y=f1_x),size=1,color="blue")  + 
  ylim(c(-2.5,2.5)) +
  xlim(c(-2.5,2.5)) +
  geom_vline(xintercept=0) +
  geom_hline(yintercept=0) 

#The Newton Method - caps N
NewtonM <- function(b0,err,fx=f1_x,fxdx=f1_dxdy) {
  b0 <- b0-(fx(b0)/fxdx(b0))

  print(paste("b0:",b0,"fx:",fx(b0),"fxdx:",fxdx(b0)))
  
  if (abs(f1_x(b0)) <= err) {
    return(b0)
  }
  if (abs(f1_x(b0)) >= err) {
    return(newtonM(b0,err))
  }
}

No comments:

Post a Comment