Pella86

bisect Error Analysis

Feb 4th, 2019
153
0
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
C++ 2.43 KB | None | 0 0
  1. #include <iostream>
  2. #include <cmath>
  3. #include <iomanip>
  4. #include <stdexcept>
  5.  
  6. using namespace std;
  7.  
  8. // This function if is 0 is the solutioon of sqrt(2) - n = 0
  9. double f(double n){
  10.     return n*n - 2;
  11.     //return n*n*n*n - n*n;
  12. }
  13.  
  14. // function that returns the "sign" of the value -1, 0 for 0 or, 1 for positive
  15. template <typename T>
  16. int sgn(T val) {
  17.     return (T(0) < val) - (val < T(0));
  18. }
  19.  
  20. // contains the result of the function bysect, the iteration count and the
  21. // resulting number to calculate the error
  22. struct BisectResult{
  23.     int itr_count;
  24.     double result;
  25. };
  26.  
  27. // type of function pointer f(x)
  28. typedef double (*fptr_t) (double);
  29.  
  30. // calculates how many iterations to reach the given tolerance
  31. BisectResult bisect(double lower, double upper, double tolerance, fptr_t func){
  32.     int itr_max = 10000;
  33.     int itr_count = 0;
  34.     BisectResult bs_res;
  35.  
  36.     double half = (upper + lower) / 2.0;
  37.     double res;
  38.  
  39.     do{
  40.         res = func(half);
  41.         bs_res.result = half;
  42.         bs_res.itr_count = itr_count;
  43.  
  44.         if( sgn<double>(res) == sgn<double>(f(lower)) ){
  45.             lower = half;
  46.         }
  47.         else{
  48.             upper = half;
  49.         }
  50.  
  51.         half = (lower + upper) / 2;
  52.  
  53.         ++itr_count;
  54.  
  55.         if(itr_count >= itr_max){
  56.             throw runtime_error("Function bisect: Max iteration reached without convergence");
  57.         }
  58.  
  59.     } while(abs(res) > tolerance);
  60.  
  61.     return bs_res;
  62. }
  63.  
  64. // calculate the error given the expected result
  65. // 0 as expected result is a pain tho
  66. double calc_error(double expected_result, double experimental_result){
  67.     double err = (expected_result - experimental_result) / expected_result;
  68.     err = err * 100;
  69.     return err;
  70.  
  71. }
  72.  
  73.  
  74. int main()
  75. {
  76.     //------------------ user inputs --------------------------------------
  77.  
  78.     // define the search interval
  79.     double lower = 0, upper = 2.0;
  80.  
  81.     //------------------ algorithm starts --------------------------------------
  82.     BisectResult br;
  83.     for(double tolerance = 0.1; tolerance >= 0.00000001; tolerance/=10){
  84.             cout << "------------------------------" << endl;
  85.             cout << "tolerance: " << tolerance << endl;
  86.             br = bisect(lower, upper, tolerance, f);
  87.             cout << "iterations: " << br.itr_count << " result: " << br.result << endl;
  88.             cout << "error: " << calc_error(sqrt(2), br.result) << "%" << endl;
  89.     }
  90.  
  91.     return 0;
  92. }
Advertisement
Add Comment
Please, Sign In to add comment