Root Finding

Metode Fixed-Point

Metode fixed-point atau metode iterasi titik tetap merupakan metode penyelesaian persamaan non-linier dengan cara menyelesaikan setiap variabel \(x\) yang ada dalam suatu persamaan dengan sebagian yang lain sehingga diperoleh \(x=g(x)\) untuk masing-masing variabel \(x\).

Algoritma Metode Iterasi Titik Tetap

1. Definisikan \(f(x)\) dan \(g(x)\)
2. Tentukan nilai toleransi \(e\) dan iterasi maksimum \(N\)
3. Tentukan tebakan awal \(x_0\)
4. Untuk iterasi \(i=1\) hingga \(N\), Hitung \(f(xi)\)
5. Akar persamaan adalah \(x\) terakhir yang diperoleh

Jumlah iterasi akan bergantung dengan nilai tebakan awal yang kita berikan. Semakin dekat nilai tersebut dengan akar, semakin cepat nilai akar diperoleh.

Berikut ini user-defined function dari fungsi fixedpoint dengan beberapa argumen yaitu ftn yang berupa function di R, x0 merupakan nilai awal, tol merupakan nilai tolerasi dan max.iter adalah iterasi maksimum. tol dan max.iter masing-masing memiliki nilai default 1e-9 dan 100.

fixedpoint <- function(ftn, fx, x0, tol = 1e-9, max.iter = 100) {
  
  xold <- x0
  xnew <- ftn(xold)
  iter <- 1
  err <- abs(xnew-xold)
  
  # membuat tempat penyimpanan hasil
  all_res <- data.frame(
    iteration=iter, x_old=xold, x_new=xnew, `f(x)`=fx(xnew), error = err
  )
  
  while ((err > tol) && (iter < max.iter)) {
    # xold digunakan untuk menyimpan akar-akar persamaan dari iterasi sebelumnya
    xold <- xnew
    # xnew digunakan untuk menghitung akar-akar persamaan pada iterasi yang sedang berjalan
    xnew <- ftn(xold) #ingat ftn bentuknya fungsi
      
    # menghitung error
    err <- abs(xnew-xold) 
    
    # iter digunakan untuk menyimpan banyaknya iterasi
    iter <- iter + 1
    
    #menyimpan semua hasil pada iterasi tertentu
    all_res <- rbind(all_res,
      data.frame(
        iteration = iter, x_old = xold, x_new = xnew, `f(x)` = fx(xnew), error = err
      )
    )
  }
  
  # output bergantung pada kesuksesan algoritma
  if (err> tol) {
    cat("Algorithm failed to converge\n")
    return(list(all_result=all_res,x_roots=NULL))
  } else {
    cat("Algorithm converged\n")
    return(list(all_result=all_res,x_roots=xnew))
  }
}

Metode Newton-Raphson

Metode Newton-Raphson merupakan metode penyelesaian persamaan nonlinear dengan menggunakan pendekatan satu titik awal dan mendekatinya dengan memperhatikan slope atau gradien. Titik pendekatan ditanyakan pada persamaan:

\[ x_{n+1}=x_n-\frac{f(x_n)}{f'(x_n)} \]

Algoritma Metode Newton-Raphson

1. Definisikan \(f(x)\) dan \(f′(x)\)
2. Tentukan nilai toleransi \(e\) dan iterasi maksimum \(N\)
3. Tentukan tebakan awal \(x_0\)
4. Untuk iterasi \(i=1\) hingga \(N\) atau \(|f(x)| \ge e\), Hitung \(x\)
5. Akar persamaan adalah nilai \(x_i\) terakhir yang diperoleh

Berikut ini fungsi newtonraphson dengan beberapa argumen yaitu ftn yang berupa function di R, x0 merupakan nilai awal, tol merupakan nilai tolerasi dan max.iter adalah iterasi maksimum. tol dan max.iter masing-masing memiliki nilai default 1e-9 dan 100.

newtonraphson <- function(ftn, x0, tol = 1e-5, max.iter = 100) {
  x <- x0
  fx <- ftn(x)
  iter <- 0
  # menyimpan semua hasil
  all_res <- data.frame(
    iteration=iter, x_old=x0, x_new=x, `f(x)`=fx[[1]], error = abs(fx[[1]])
  )
  
  while ((abs(fx[1]) > tol) && (iter < max.iter)) {
    x0 <- x
    x <- x - fx[1]/fx[2]
    fx <- ftn(x)
    iter <- iter + 1
    
    #menyimpan semua hasil 
    all_res <- rbind(all_res,
      data.frame(
        iteration=iter, x_old=x0, x_new=x, `f(x)`=fx[[1]], error = abs(fx[[1]])
      )
    )
  }
  
  # output bergantung pada kesuksesan algoritma
  if (abs(fx[1]) > tol) {
    cat("Algorithm failed to converge\n")
    return(list(all_result=all_res,x_roots=NULL))
  } else {
    cat("Algorithm converged\n")
    return(list(all_result=all_res,x_roots=x))
  }
}

Metode Secant

Metode Secant merupakan perbaikan dari metode Newton-Raphson, dimana kemiringan dua titik dinyatakan secara diskrit dengan mengambil bentuk garis lurus yang melalui satu titik. Persamaan yang disajikan untuk metode Secant:

\[ x_{n+1}=x_n-f(x_n)\frac{f(x_n-x_{n+1})}{f(x_n)-f(x_{n+1})} \]

Algoritma Metode Secant

1. Definisikan \(f(x)\) dan \(f′(x)\)
2. Tentukan nilai toleransi \(e\) dan iterasi maksimum \(N\)
3. Tentukan tebakan awal \(x_0\) dan \(x_1\)
4. Hitung \(f(x_0)\) dan \(f(x_1)\)
5. Untuk iterasi \(i=1\) hingga \(N\) atau \(|f(x)| \ge e\), Hitung \(x\)
6. Akar persamaan adalah nilai \(x\) terakhir yang diperoleh

Fungsi secant memiliki beberapa argumen yaitu ftn yang berupa function di R, x0 dan x1 merupakan nilai awal, tol merupakan nilai tolerasi dan max.iter adalah iterasi maksimum. tol dan max.iter masing-masing memiliki nilai default 1e-9 dan 100.

secant <- function(ftn, x0, x1, tol = 1e-9, max.iter = 100) {
  xa <- x0
  fxa <- ftn(xa)
  xb <- x1
  fxb <- ftn(xb)
  iter <- 0
  
  # menyimpan semua hasil
  all_res <- data.frame(
    iteration=iter, xa=x0, xb=x1, `f(x_a)`=fxa, `f(x_b)`=fxb, error = abs(xb-xa)
  )
  
  while ((abs(xb-xa) > tol) && (iter < max.iter)) {
    xt <- fxb*((xb-xa)/(fxb-fxa))
    xc <- xb-xt
  
    iter <- iter + 1
  
    xa <- xb
    xb <- xc
    fxa <- ftn(xa)
    fxb <- ftn(xb)
    
    # menyimpan semua hasil
    all_res <- rbind(all_res,
      data.frame(
        iteration=iter, xa=xa, xb=xb, `f(x_a)`=fxa, `f(x_b)`=fxb, error = abs(xb-xa)
      )
    )
  }

  if (abs(xb-xa) > tol) {
    cat("Algorithm failed to converge\n")
    return(list(all_result=all_res,x_roots=NULL))
  } else {
    cat("Algorithm converged\n")
    return(list(all_result=all_res,x_roots=xb))
  }
}

Secara umum metode Secant menawarkan sejumlah keuntungan dibanding metode lainnya. Pertama, seperti metode Newton-Raphson dan tidak seperti metode tertutup lainnya, metode ini tidak memerlukan rentang pencarian akar penyelesaian. Kedua, tidak seperti metode Newton-Raphson, metode ini tidak memerlukan pencarian turunan pertama persamaan nonlinear secara analitik, dimana tidak dapat dilakukan otomasi pada setiap kasus.

Adapun kerugian dari metode ini adalah berpotensi menghasilkan hasil yang tidak konvergen sama seperti metode terbuka lainnya. Selain itu, kecepatan konvergensinya lebih lambat dibanding metode Newton-Raphson.

Metode Bisection

Prinsip metode bagi dua adalah mengurung akar fungsi pada interval \(x=[a,b]\) atau pada nilai \(x\) batas bawah \(a\) dan batas atas \(b\). Selanjutnya interval tersebut terus menerus dibagi 2 hingga sekecil mungkin, sehingga nilai hampiran yang dicari dapat ditentukan dengan tingkat toleransi tertentu.

Algoritma Metode Biseksi

1. Definisikan \(f(x)\)
2. Tentukan rentang untuk \(x\) yang berupa batas bawah \(a\) dan batas atas \(b\)
3. Tentukan nilai toleransi \(e\) dan iterasi maksimum \(N\)
4. Hitung \(f(a)\) dan \(f(b)\)
   
5. Hitung \(x=\frac{a+b}{2}\)
6. Hitung \(f(x)\)
7. Bila\(f(x).f(a)<0\), maka \(b=x\) dan \(f(b)=f(x)\). Bila tidak, \(a=x\) dan \(f(a)=f(x)\)
8. Bila \(|b−a|<e\) atau iterasi maksimum maka proses dihentikan dan didapatkan akar=\(x\), dan bila tidak ulangi langkah 6
9. Jika sudah diperoleh nilai dibawah nilai toleransi, nilai akar selanjutnya dihitung berdasarkan persamaan 5 dengan nilai \(a\) dan \(b\) merupakan nilai baru yang diperoleh dari proses iterasi.

Fungsi bisection memiliki beberapa argumen yaitu fx yang berupa function di R, xl dan xr merupakan nilai awal, tol merupakan nilai tolerasi dan max.iter adalah iterasi maksimum. tol dan max.iter masing-masing memiliki nilai default 1e-9 dan 100.

bisection <- function(f, xl, xr, tol = 1e-7,max.iter = 100) {
  # Jika tanda dari hasil perkalian  evaluasi fungsi pada titik xl dan xr bernilai positif (f(xl)*f(xr)>0) maka program dihentikan
  if (f(xl)*f(xr) > 0) {
    stop('signs of f(xl) and f(xr) must be differ')
  }
  
  # menyimpan semua hasil
  all_res <-data.frame(
    iteration=0, x_left=0, x_right=0, x_mid = 0, `f(xm)`=0, error = 1000
  )
  
  for (iter in 1:max.iter) {
    xm <- (xl + xr) / 2 # menghitung nilai tengah
    
    # menyimpan semua hasil
    all_res <-rbind(all_res,
      data.frame(
        iteration=iter, x_left=xl, x_right=xr, x_mid = xm, `f(xm)`=f(xm), error = abs(xr - xl)
      )
    )
    
    #Jika fungsinya sama dengan 0 di titik tengah atau titik tengah di bawah toleransi yang diinginkan, hentikan dan keluarkan nilai tengah sebagai akar
    if ((f(xm) == 0) || abs(xr - xl) < tol) {
      return(list(all_result=all_res,x_roots=xm))
    }
    
    # Jika diperlukan iterasi lain,
    # periksa tanda-tanda fungsi pada titik xm dan xl dan tetapkan kembali
    # xl atau xr sebagai titik tengah yang akan digunakan pada iterasi berikutnya.    
    ifelse(sign(f(xm)) == sign(f(xl)), xl <- xm,xr <- xm)
  }
}

Metode biseksi merupakan metode yang paling mudah dan paling sederhana dibanding metode lainnya. Adapun sifat metode ini antara lain:

• Konvergensi lambat
• Caranya mudah
• Tidak dapat digunakan untuk mencari akar imaginer
• Hanya dapat mencari satu akar pada satu siklus.

Soal 1

Misalkan diberikan suatu persamaan \(x^3-2x-5=0\).

1. Gambarkan fungsi tersebut!
2. Carilah akar-akar persamaan dengan menggunakan metode fixed-point, newton-raphson, secant dan bisection!
3. Buatlah animasi untuk metode newton-raphson dan bisection dengan bantuan package gganimate!

Jawaban 1

Mendefinisikan persamaan \(x^{3}-2x-5=0\) dalam bentuk fungsi matematika:

\[f(x)=x^{3}-2x-5\]

kemudian buat fungsi matematika tersebut ke dalam R.

fx <- function(x) {
  x^3 - 2 * x - 5
}

Khusus metode fixed-point, kita ubah persamaan \(x^{3}-2x-5=0\) menjadi \(x=g(x)\) dengan cara berikut:

\[x^{3}-2x-5=0\] \[2x=5-x^3\] \[x=\frac{5-x^3}{2}\] sehingga diperoleh

\[g(x)=\frac{5-x^3}{2}\]

gx <- function(x){
  5-x^{3}/(2)
}

mendefinisikan persamaan dalam bentuk fungsi untuk newton-raphson

fx_nr <- function(x) {
  rumus_fungsi <- function(x) x^3-2*x-5
  fungsi <- rumus_fungsi(x)
  # mencari turunan dengan menggunakan fungsi Deriv
  # dari package Deriv
  rumus_turunan <- Deriv::Deriv(fx)
  turunan <- rumus_turunan(x)
  return(c("fungsi"=fungsi,"turunan"=turunan))
}

1. Gambarkan fungsi tersebut!

curve(fx, xlim=c(-3,3), col='steelblue',lwd=2)
abline(h=0)
abline(v=0)

2. Carilah akar-akar persamaan dengan menggunakan metode fixed-point, newton-raphson, secant dan bisection!

#fixed-point
fixed_point_res <- fixedpoint(ftn = gx, fx = fx, x0 =0)

## Algorithm failed to converge

fixed_point_res

## $all_result
##     iteration          x_old          x_new           f.x.         error
## 1           1   0.000000e+00   5.000000e+00   1.100000e+02  5.000000e+00
## 2           2   5.000000e+00  -5.750000e+01  -1.899994e+05  6.250000e+01
## 3           3  -5.750000e+01   9.505969e+04   8.589921e+14  9.511719e+04
## 4           4   9.505969e+04  -4.294960e+14  -7.922777e+43  4.294960e+14
## 5           5  -4.294960e+14   3.961389e+43  6.216449e+130  3.961389e+43
## 6           6   3.961389e+43 -3.108224e+130           -Inf 3.108224e+130
## 7           7 -3.108224e+130            Inf            NaN           Inf
## 8           8            Inf           -Inf            NaN           Inf
## 9           9           -Inf            Inf            NaN           Inf
## 10         10            Inf           -Inf            NaN           Inf
## 11         11           -Inf            Inf            NaN           Inf
## 12         12            Inf           -Inf            NaN           Inf
## 13         13           -Inf            Inf            NaN           Inf
## 14         14            Inf           -Inf            NaN           Inf
## 15         15           -Inf            Inf            NaN           Inf
## 16         16            Inf           -Inf            NaN           Inf
## 17         17           -Inf            Inf            NaN           Inf
## 18         18            Inf           -Inf            NaN           Inf
## 19         19           -Inf            Inf            NaN           Inf
## 20         20            Inf           -Inf            NaN           Inf
## 21         21           -Inf            Inf            NaN           Inf
## 22         22            Inf           -Inf            NaN           Inf
## 23         23           -Inf            Inf            NaN           Inf
## 24         24            Inf           -Inf            NaN           Inf
## 25         25           -Inf            Inf            NaN           Inf
## 26         26            Inf           -Inf            NaN           Inf
## 27         27           -Inf            Inf            NaN           Inf
## 28         28            Inf           -Inf            NaN           Inf
## 29         29           -Inf            Inf            NaN           Inf
## 30         30            Inf           -Inf            NaN           Inf
## 31         31           -Inf            Inf            NaN           Inf
## 32         32            Inf           -Inf            NaN           Inf
## 33         33           -Inf            Inf            NaN           Inf
## 34         34            Inf           -Inf            NaN           Inf
## 35         35           -Inf            Inf            NaN           Inf
## 36         36            Inf           -Inf            NaN           Inf
## 37         37           -Inf            Inf            NaN           Inf
## 38         38            Inf           -Inf            NaN           Inf
## 39         39           -Inf            Inf            NaN           Inf
## 40         40            Inf           -Inf            NaN           Inf
## 41         41           -Inf            Inf            NaN           Inf
## 42         42            Inf           -Inf            NaN           Inf
## 43         43           -Inf            Inf            NaN           Inf
## 44         44            Inf           -Inf            NaN           Inf
## 45         45           -Inf            Inf            NaN           Inf
## 46         46            Inf           -Inf            NaN           Inf
## 47         47           -Inf            Inf            NaN           Inf
## 48         48            Inf           -Inf            NaN           Inf
## 49         49           -Inf            Inf            NaN           Inf
## 50         50            Inf           -Inf            NaN           Inf
## 51         51           -Inf            Inf            NaN           Inf
## 52         52            Inf           -Inf            NaN           Inf
## 53         53           -Inf            Inf            NaN           Inf
## 54         54            Inf           -Inf            NaN           Inf
## 55         55           -Inf            Inf            NaN           Inf
## 56         56            Inf           -Inf            NaN           Inf
## 57         57           -Inf            Inf            NaN           Inf
## 58         58            Inf           -Inf            NaN           Inf
## 59         59           -Inf            Inf            NaN           Inf
## 60         60            Inf           -Inf            NaN           Inf
## 61         61           -Inf            Inf            NaN           Inf
## 62         62            Inf           -Inf            NaN           Inf
## 63         63           -Inf            Inf            NaN           Inf
## 64         64            Inf           -Inf            NaN           Inf
## 65         65           -Inf            Inf            NaN           Inf
## 66         66            Inf           -Inf            NaN           Inf
## 67         67           -Inf            Inf            NaN           Inf
## 68         68            Inf           -Inf            NaN           Inf
## 69         69           -Inf            Inf            NaN           Inf
## 70         70            Inf           -Inf            NaN           Inf
## 71         71           -Inf            Inf            NaN           Inf
## 72         72            Inf           -Inf            NaN           Inf
## 73         73           -Inf            Inf            NaN           Inf
## 74         74            Inf           -Inf            NaN           Inf
## 75         75           -Inf            Inf            NaN           Inf
## 76         76            Inf           -Inf            NaN           Inf
## 77         77           -Inf            Inf            NaN           Inf
## 78         78            Inf           -Inf            NaN           Inf
## 79         79           -Inf            Inf            NaN           Inf
## 80         80            Inf           -Inf            NaN           Inf
## 81         81           -Inf            Inf            NaN           Inf
## 82         82            Inf           -Inf            NaN           Inf
## 83         83           -Inf            Inf            NaN           Inf
## 84         84            Inf           -Inf            NaN           Inf
## 85         85           -Inf            Inf            NaN           Inf
## 86         86            Inf           -Inf            NaN           Inf
## 87         87           -Inf            Inf            NaN           Inf
## 88         88            Inf           -Inf            NaN           Inf
## 89         89           -Inf            Inf            NaN           Inf
## 90         90            Inf           -Inf            NaN           Inf
## 91         91           -Inf            Inf            NaN           Inf
## 92         92            Inf           -Inf            NaN           Inf
## 93         93           -Inf            Inf            NaN           Inf
## 94         94            Inf           -Inf            NaN           Inf
## 95         95           -Inf            Inf            NaN           Inf
## 96         96            Inf           -Inf            NaN           Inf
## 97         97           -Inf            Inf            NaN           Inf
## 98         98            Inf           -Inf            NaN           Inf
## 99         99           -Inf            Inf            NaN           Inf
## 100       100            Inf           -Inf            NaN           Inf
## 
## $x_roots
## NULL

# newton raphson
nr_res <- newtonraphson(ftn=fx_nr, x0=0)

## Algorithm converged

nr_res

## $all_result
##          iteration      x_old      x_new          f.x.        error
## 1                0  0.0000000  0.0000000 -5.000000e+00 5.000000e+00
## fungsi           1  0.0000000 -2.5000000 -1.562500e+01 1.562500e+01
## fungsi1          2 -2.5000000 -1.5671642 -5.714632e+00 5.714632e+00
## fungsi2          3 -1.5671642 -0.5025924 -4.121770e+00 4.121770e+00
## fungsi3          4 -0.5025924 -3.8207065 -5.313249e+01 5.313249e+01
## fungsi4          5 -3.8207065 -2.5493934 -1.647076e+01 1.647076e+01
## fungsi5          6 -2.5493934 -1.6081115 -5.942390e+00 5.942390e+00
## fungsi6          7 -1.6081115 -0.5761004 -4.039002e+00 4.039002e+00
## fungsi7          8 -0.5761004 -4.5977096 -9.299526e+01 9.299526e+01
## fungsi8          9 -4.5977096 -3.0835431 -2.815198e+01 2.815198e+01
## fungsi9         10 -3.0835431 -2.0221943 -9.224909e+00 9.224909e+00
## fungsi10        11 -2.0221943 -1.1237641 -4.171613e+00 4.171613e+00
## fungsi11        12 -1.1237641  1.2086516 -5.651658e+00 5.651658e+00
## fungsi12        13  1.2086516  3.5807900  3.375152e+01 3.375152e+01
## fungsi13        14  3.5807900  2.6552332  8.409627e+00 8.409627e+00
## fungsi14        15  2.6552332  2.2161063  1.451367e+00 1.451367e+00
## fungsi15        16  2.2161063  2.1021250  8.489238e-02 8.489238e-02
## fungsi16        17  2.1021250  2.0945836  3.582354e-04 3.582354e-04
## fungsi17        18  2.0945836  2.0945515  6.472653e-09 6.472653e-09
## 
## $x_roots
##   fungsi 
## 2.094551

# secant
sec_res <- secant(fx, x0=1, x1=3)

## Algorithm converged

sec_res

## $all_result
##    iteration       xa       xb        f.x_a.        f.x_b.        error
## 1          0 1.000000 3.000000 -6.000000e+00  1.600000e+01 2.000000e+00
## 2          1 3.000000 1.545455  1.600000e+01 -4.399699e+00 1.454545e+00
## 3          2 1.545455 1.859163 -4.399699e+00 -2.292151e+00 3.137087e-01
## 4          3 1.859163 2.200350 -2.292151e+00  1.252384e+00 3.411868e-01
## 5          4 2.200350 2.079799  1.252384e+00 -1.632926e-01 1.205509e-01
## 6          5 2.079799 2.093704 -1.632926e-01 -9.451895e-03 1.390506e-02
## 7          6 2.093704 2.094559 -9.451895e-03  7.903548e-05 8.543201e-04
## 8          7 2.094559 2.094551  7.903548e-05 -3.771081e-08 7.084470e-06
## 9          8 2.094551 2.094551 -3.771081e-08 -1.501022e-13 3.378656e-09
## 10         9 2.094551 2.094551 -1.501022e-13 -8.881784e-16 1.332268e-14
## 
## $x_roots
## [1] 2.094551

# bisection
bis_res <- bisection(fx, xl=-3, xr=3)
bis_res

## $all_result
##    iteration    x_left  x_right    x_mid         f.xm.        error
## 1          0  0.000000 0.000000 0.000000  0.000000e+00 1.000000e+03
## 2          1 -3.000000 3.000000 0.000000 -5.000000e+00 6.000000e+00
## 3          2  0.000000 3.000000 1.500000 -4.625000e+00 3.000000e+00
## 4          3  1.500000 3.000000 2.250000  1.890625e+00 1.500000e+00
## 5          4  1.500000 2.250000 1.875000 -2.158203e+00 7.500000e-01
## 6          5  1.875000 2.250000 2.062500 -3.513184e-01 3.750000e-01
## 7          6  2.062500 2.250000 2.156250  7.127991e-01 1.875000e-01
## 8          7  2.062500 2.156250 2.109375  1.668358e-01 9.375000e-02
## 9          8  2.062500 2.109375 2.085938 -9.567881e-02 4.687500e-02
## 10         9  2.085938 2.109375 2.097656  3.471428e-02 2.343750e-02
## 11        10  2.085938 2.097656 2.091797 -3.069771e-02 1.171875e-02
## 12        11  2.091797 2.097656 2.094727  1.954348e-03 5.859375e-03
## 13        12  2.091797 2.094727 2.093262 -1.438516e-02 2.929688e-03
## 14        13  2.093262 2.094727 2.093994 -6.218774e-03 1.464844e-03
## 15        14  2.093994 2.094727 2.094360 -2.133056e-03 7.324219e-04
## 16        15  2.094360 2.094727 2.094543 -8.956468e-05 3.662109e-04
## 17        16  2.094543 2.094727 2.094635  9.323389e-04 1.831055e-04
## 18        17  2.094543 2.094635 2.094589  4.213739e-04 9.155273e-05
## 19        18  2.094543 2.094589 2.094566  1.659013e-04 4.577637e-05
## 20        19  2.094543 2.094566 2.094555  3.816751e-05 2.288818e-05
## 21        20  2.094543 2.094555 2.094549 -2.569879e-05 1.144409e-05
## 22        21  2.094549 2.094555 2.094552  6.234310e-06 5.722046e-06
## 23        22  2.094549 2.094552 2.094551 -9.732252e-06 2.861023e-06
## 24        23  2.094551 2.094552 2.094551 -1.748974e-06 1.430511e-06
## 25        24  2.094551 2.094552 2.094552  2.242667e-06 7.152557e-07
## 26        25  2.094551 2.094552 2.094552  2.468460e-07 3.576279e-07
## 27        26  2.094551 2.094552 2.094551 -7.510643e-07 1.788139e-07
## 28        27  2.094551 2.094552 2.094551 -2.521091e-07 8.940697e-08
## 
## $x_roots
## [1] 2.094551

3. Membuat animasi untuk metode newton-raphson dan bisection dengan bantuan package gganimate

Memanggil package yang dibutuhkan

Untuk membuat animasi dengan gganimate, terlebih dahulu buat plot dengan ggplot2

Animasi Newton Raphson

# menyimpan hasil newton rapshon kedalam objek baru
nr_data <- nr_res$all_result %>% 
# convert iteration ke dalam bentuk faktor untuk keperluan animasi
  mutate(iteration=as.factor(iteration))

#membuat plot dengan ggplot2 
p_nr <- ggplot(data=nr_data) +
  # membuat grafik fungsi
  stat_function(aes(x=seq(-3, 3, length.out=nrow(nr_data))), fun=fx , color="steelblue", size=1.25) +
  geom_hline(yintercept = 0) + # membuat sumbu y
  geom_vline(xintercept = 0) + # membuat sumbu x
  # membuat titik untuk setiap x_new
  geom_point(aes(x=x_new, y=fx(x_new))) +
  # membuat text untuk setiap x_new
  geom_text(aes(x=x_new, y=fx(x_new), label=str_c("Iteration ", iteration)), nudge_y = 4 ) +
  xlab("x") + ylab("f(x)") + theme_bw() +
  # membuat animasi dengan gganimate
  # animasi didasarkan pada kolom iteration
  transition_states(iteration, transition_length=3, state_length =3, wrap=F)

## Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
## ℹ Please use `linewidth` instead.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.

Jika ingin menampilkan langsung pada RStudio dapat menjalankan sintaks dibawah

animate(p_nr, nframes = 200, fps=5)

Jika ingin menyimpan hasilnya dalam bentuk .gif gunakan sintaks dibawah ini.

anim_save("img/nr_anim.gif", animate(p_nr,nframes = 200,fps=5))

Animasi Secant

#menyimpan hasil newton rapshon kedalam objek baru
sec_data <- sec_res$all_result %>%   
# convert iteration ke dalam bentuk faktor untuk keperluan animasi
  mutate(iteration=as.factor(iteration))

#membuat plot dengan ggplot2 
p_bis <- ggplot(data=sec_data) +
  # membuat grafik fungsi
  stat_function(aes(x=seq(-3, 3, length.out=nrow(sec_data))), fun=fx, color="steelblue", size=1.25 ) +
  geom_hline(yintercept = 0) + # membuat sumbu y
  geom_vline(xintercept = 0) + # membuat sumbu x
  # membuat titik untuk setiap x_left
  geom_point(aes(x=xa, y=fx(xa)), color="red" ) +
  # membuat titik untuk setiap x_right
  geom_point(aes(x=xb, y=fx(xb)), color="green") +
  # membuat text untuk setiap x_left
  geom_text(aes(x=xa, y=fx(xa), label=str_c("Iteration ", iteration)), nudge_y = 2, color="red") +
  # membuat text untuk setiap x_right
  geom_text(aes(x=xb, y=fx(xb), label=str_c("Iteration ", iteration)), nudge_y = 2, color="green") +
  xlab("x") + ylab("f(x)") + theme_bw() +
  # membuat animasi dengan gganimate
  # animasi didasarkan pada kolom iteration
  transition_states(iteration, transition_length=3,state_length =3, wrap=F)

anim_save("img/sec_anim.gif",animate(p_bis,nframes = 200,fps=5))

Animasi Bisection

#menyimpan hasil newton rapshon kedalam objek baru
bs_data <- bis_res$all_result %>%   
# convert iteration ke dalam bentuk faktor untuk keperluan animasi
  mutate(iteration=as.factor(iteration))

#membuat plot dengan ggplot2 
p_bis <- ggplot(data=bs_data) +
  # membuat grafik fungsi
  stat_function(aes(x=seq(-3, 3, length.out=nrow(bs_data))), fun=fx, color="steelblue", size=1.25 ) +
  geom_hline(yintercept = 0) + # membuat sumbu y
  geom_vline(xintercept = 0) + # membuat sumbu x
  # membuat titik untuk setiap x_left
  geom_point(aes(x=x_left, y=fx(x_left)), color="red" ) +
  # membuat titik untuk setiap x_right
  geom_point(aes(x=x_right, y=fx(x_right)), color="green") +
  # membuat titik untuk setiap x_mid
  geom_point(aes(x=x_mid, y=fx(x_mid)), color="black") +
  # membuat text untuk setiap x_left
  geom_text(aes(x=x_left, y=fx(x_left), label=str_c("Iteration ", iteration)), nudge_y = 2, color="red") +
  # membuat text untuk setiap x_right
  geom_text(aes(x=x_right, y=fx(x_right), label=str_c("Iteration ", iteration)), nudge_y = 2, color="green") +
  # membuat text untuk setiap x_mid
  geom_text(aes(x=x_mid, y=fx(x_mid), label=str_c("Iteration ", iteration)), nudge_y = 2, color="black") +
  xlab("x") + ylab("f(x)") + theme_bw() +
  # membuat animasi dengan gganimate
  # animasi didasarkan pada kolom iteration
  transition_states(iteration, transition_length=3,state_length =3, wrap=F)

anim_save("img/bis_anim.gif",animate(p_bis,nframes = 200,fps=5))

Soal 2

Misalkan diberikan suatu persamaan \(x^3-x-1=0\).

1. Gambarkan fungsi tersebut!
2. Carilah akar-akar persamaan dengan menggunakan metode fixed-point, newton-raphson, secant dan bisection!

Jawaban 2

Mendefinisikan persamaan \(x^3-x-1=0\) dalam bentuk fungsi matematika:

\[f(x)=x^3-x-1\]

Lalu masukan fungsi matematika tersebut menjadi fungsi dalam R.

fx <- function(x) {
  x^3 - x - 1
}

Khusus metode fixed-point, kita ubah persamaan \(x^3-x-1=0\) menjadi \(x=g(x)\) dengan cara berikut:

\[x^3-x-1=0\] \[x=x^3-1\]

sehingga diperoleh

\[g(x)=x^3-1\]

gx <- function(x){
  x^{3}-1
}

mendefinikan persamaan dalam bentuk fungsi untuk newton-raphson

fx_nr <- function(x) {
  rumus_fungsi <- function(x)  x^3 - x - 1
  fungsi <- rumus_fungsi(x)
  # mencari turunan dengan menggunakan fungsi Deriv
  # dari package Deriv
  rumus_turunan <- Deriv::Deriv(fx)
  turunan <- rumus_turunan(x)
  return(c("fungsi"=fungsi,"turunan"=turunan))
}

1. Gambarkan fungsi tersebut!

curve(fx,xlim=c(-3,3), col='steelblue',lwd=2)
abline(h=0)
abline(v=0)

2. Carilah akar-akar persamaan dengan menggunakan metode fixed-point, newton-raphson, secant dan bisection!

#fixed-point
fixed_point_res <- try(fixedpoint(ftn = gx, fx = fx, x0 =0))

## Error in while ((err > tol) && (iter < max.iter)) { : 
##   missing value where TRUE/FALSE needed

fixed_point_res

## [1] "Error in while ((err > tol) && (iter < max.iter)) { : \n  missing value where TRUE/FALSE needed\n"
## attr(,"class")
## [1] "try-error"
## attr(,"condition")
## <simpleError in while ((err > tol) && (iter < max.iter)) {    xold <- xnew    xnew <- ftn(xold)    err <- abs(xnew - xold)    iter <- iter + 1    all_res <- rbind(all_res, data.frame(iteration = iter, x_old = xold,         x_new = xnew, `f(x)` = fx(xnew), error = err))}: missing value where TRUE/FALSE needed>

# newton raphson
nr_res <-newtonraphson(ftn=fx_nr, x0=0)

## Algorithm converged

nr_res

## $all_result
##          iteration      x_old      x_new          f.x.        error
## 1                0  0.0000000  0.0000000 -1.000000e+00 1.000000e+00
## fungsi           1  0.0000000 -1.0000000 -1.000000e+00 1.000000e+00
## fungsi1          2 -1.0000000 -0.5000000 -6.250000e-01 6.250000e-01
## fungsi2          3 -0.5000000 -3.0000000 -2.500000e+01 2.500000e+01
## fungsi3          4 -3.0000000 -2.0384615 -7.432010e+00 7.432010e+00
## fungsi4          5 -2.0384615 -1.3902821 -2.296973e+00 2.296973e+00
## fungsi5          6 -1.3902821 -0.9116119 -8.459706e-01 8.459706e-01
## fungsi6          7 -0.9116119 -0.3450285 -6.960453e-01 6.960453e-01
## fungsi7          8 -0.3450285 -1.4277507 -2.482679e+00 2.482679e+00
## fungsi8          9 -1.4277507 -0.9424179 -8.945920e-01 8.945920e-01
## fungsi9         10 -0.9424179 -0.4049494 -6.614559e-01 6.614559e-01
## fungsi10        11 -0.4049494 -1.7069046 -4.266202e+00 4.266202e+00
## fungsi11        12 -1.7069046 -1.1557564 -1.388072e+00 1.388072e+00
## fungsi12        13 -1.1557564 -0.6941918 -6.403408e-01 6.403408e-01
## fungsi13        14 -0.6941918  0.7424943 -1.333159e+00 1.333159e+00
## fungsi14        15  0.7424943  2.7812959  1.773372e+01 1.773372e+01
## fungsi15        16  2.7812959  1.9827252  4.811763e+00 4.811763e+00
## fungsi16        17  1.9827252  1.5369274  1.093519e+00 1.093519e+00
## fungsi17        18  1.5369274  1.3572625  1.430341e-01 1.430341e-01
## fungsi18        19  1.3572625  1.3256631  4.034214e-03 4.034214e-03
## fungsi19        20  1.3256631  1.3247188  3.545493e-06 3.545493e-06
## 
## $x_roots
##   fungsi 
## 1.324719

# secant
sec_res <- secant(fx, x0=1, x1=3)

## Algorithm converged

sec_res

## $all_result
##    iteration       xa       xb        f.x_a.        f.x_b.        error
## 1          0 1.000000 3.000000 -1.000000e+00  2.300000e+01 2.000000e+00
## 2          1 3.000000 1.083333  2.300000e+01 -8.119213e-01 1.916667e+00
## 3          2 1.083333 1.148686 -8.119213e-01 -6.330171e-01 6.535308e-02
## 4          3 1.148686 1.379925 -6.330171e-01  2.477203e-01 2.312390e-01
## 5          4 1.379925 1.314886  2.477203e-01 -4.154635e-02 6.503936e-02
## 6          5 1.314886 1.324227 -4.154635e-02 -2.091088e-03 9.341372e-03
## 7          6 1.324227 1.324722 -2.091088e-03  1.930334e-05 4.950831e-04
## 8          7 1.324722 1.324718  1.930334e-05 -8.827270e-09 4.528428e-06
## 9          8 1.324718 1.324718 -8.827270e-09 -3.685940e-14 2.069869e-09
## 10         9 1.324718 1.324718 -3.685940e-14  2.220446e-16 8.659740e-15
## 
## $x_roots
## [1] 1.324718

# bisection
bis_res <- bisection(fx, xl=-3, xr=3)
bis_res

## $all_result
##    iteration    x_left  x_right    x_mid         f.xm.        error
## 1          0  0.000000 0.000000 0.000000  0.000000e+00 1.000000e+03
## 2          1 -3.000000 3.000000 0.000000 -1.000000e+00 6.000000e+00
## 3          2  0.000000 3.000000 1.500000  8.750000e-01 3.000000e+00
## 4          3  0.000000 1.500000 0.750000 -1.328125e+00 1.500000e+00
## 5          4  0.750000 1.500000 1.125000 -7.011719e-01 7.500000e-01
## 6          5  1.125000 1.500000 1.312500 -5.151367e-02 3.750000e-01
## 7          6  1.312500 1.500000 1.406250  3.746643e-01 1.875000e-01
## 8          7  1.312500 1.406250 1.359375  1.526146e-01 9.375000e-02
## 9          8  1.312500 1.359375 1.335938  4.834890e-02 4.687500e-02
## 10         9  1.312500 1.335938 1.324219 -2.127945e-03 2.343750e-02
## 11        10  1.324219 1.335938 1.330078  2.297349e-02 1.171875e-02
## 12        11  1.324219 1.330078 1.327148  1.038860e-02 5.859375e-03
## 13        12  1.324219 1.327148 1.325684  4.121792e-03 2.929688e-03
## 14        13  1.324219 1.325684 1.324951  9.947910e-04 1.464844e-03
## 15        14  1.324219 1.324951 1.324585 -5.671101e-04 7.324219e-04
## 16        15  1.324585 1.324951 1.324768  2.137072e-04 3.662109e-04
## 17        16  1.324585 1.324768 1.324677 -1.767348e-04 1.831055e-04
## 18        17  1.324677 1.324768 1.324722  1.847785e-05 9.155273e-05
## 19        18  1.324677 1.324722 1.324699 -7.913056e-05 4.577637e-05
## 20        19  1.324699 1.324722 1.324711 -3.032687e-05 2.288818e-05
## 21        20  1.324711 1.324722 1.324717 -5.924640e-06 1.144409e-05
## 22        21  1.324717 1.324722 1.324719  6.276573e-06 5.722046e-06
## 23        22  1.324717 1.324719 1.324718  1.759583e-07 2.861023e-06
## 24        23  1.324717 1.324718 1.324717 -2.874343e-06 1.430511e-06
## 25        24  1.324717 1.324718 1.324718 -1.349193e-06 7.152557e-07
## 26        25  1.324718 1.324718 1.324718 -5.866174e-07 3.576279e-07
## 27        26  1.324718 1.324718 1.324718 -2.053296e-07 1.788139e-07
## 28        27  1.324718 1.324718 1.324718 -1.468565e-08 8.940697e-08
## 
## $x_roots
## [1] 1.324718

set.seed(0)
x_values <- 1:30
r_values <- seq(0, 1, by=0.025)

x <- rep(x_values, length(r_values))
y <- unlist(map(r_values, ~runif(x_values, .)*x_values))
r <- rep(r_values, each=length(x_values))


data <- data.frame(x=x, y=y, r=r)

legend <- c("Linear Regression Line" = "blue", "SD Line" = "red")

anim_graph <- ggplot(data %>% 
                       group_by(r) %>% 
                       mutate(
                             sd_slope=sd(y)/sd(x), 
                             sd_intercept = mean(y)- (sd(y)/sd(x))*mean(x),
                             r_slope = r*sd(y)/sd(x),
                             r_intercept = mean(y)- (r*sd(y)/sd(x))*mean(x))
                      , aes(x, y)) + 
  geom_point(color = "purple4", alpha = 0.7, size = 3) + 
  #sd line
  geom_abline(aes(intercept=sd_intercept, slope=sd_slope), col="firebrick1") +
  #corr. coefficient line
  geom_abline(aes(intercept=r_intercept, slope=r_slope), col="royalblue2") +
  labs(title = "Showing Change in the Coefficient Correlation r",
    subtitle = 'r =  {round(frame_time, 2)}',
       x = 'x',
       y = 'y') +
  theme(panel.background = element_blank()) +
  transition_time(r) +
  ease_aes('linear')

anim_graph