Search code examples
cnumerical-methods

modifying regula falsi method to the secant method


I have implemented the regula falsi method. I am trying to modify it so it becomes the secant method. A pdf I read mentioned that it is essentially the same with just one change. Future guesses for my 'm' value should have a slightly different formula, instead of:

m = a - f(a) * ( (b-a)/( f(b)-f(a) ) );

it should be:

m = a - f(a) * ( (m-a)/( f(m)-f(a) ) );

but unfortunately it doesn't work (It never finds the root). What should I fix to get this into the secant method?

source code as follows:

#include <stdio.h>
#include <math.h>


void
secant(double a, double b, double e, double (*f)(double), int maxiter ) {
  double m, fm, fa, fb;
  int i;

  fa=(*f)(a);
  fb=(*f)(b);

  m = a - fa * ( (b-a)/( fb - fa ) );

  fm=(*f)(m);


   for(i=0; i<maxiter; i++) {
    if ( fabs(fm) <= e ) {                               
      printf("f(%f) = %f\n", m, fm);
      return;
    } else if ((fa*fm) < 0) {
      b=m;
      fb=fm;
    } else {
      a=m;
      fa=fm;
    }

     // the guess below works for regula falsi method:   
     // m = a - fa * ( (b-a)/(fb - fa)); 

     //this was supposed to be the change to turn this into the secant method 
     m = a - fa * ( (m-a)/(fm - fa) ); 

     fm=(*f)(m);
  }
}

int main(){
secant(1,4,0.0001,sin,500);
return 0;
}

Thanks in advance

EDIT: Ok after playing around with pen and paper I finally got it it wasnt a simple change as I initially thought:

void secant(double a, double b, double e, double (*f)(double), int maxiter ) {
  double m, fm, fa, fb;
  int i;
  fa=(*f)(a);
  fb=(*f)(b);

   for(i=0; i<maxiter; i++) {
     m = a - fa * ( (b-a)/(fb - fa) );
     fm=(*f)(m);
     if ( fabs(fm) <= e ) {                               
        printf("f(%f)=%f, iter: %d\n", m,fm,i);
         return;
     }
     a=b;
     b=m;
     fa=fb;
     fb=fm;
  }
}

Solution

  • It's easier for the secant method to not find the root. Are you sure it should find it?

    For testing, here is an example: http://www.mathcs.emory.edu/ccs/ccs315/ccs315/node18.html (example 4.7) You'd want to run that example ( f(x)=x^6-x-1 , x0=1 x1=2, root x=1.347)