Bisection method. Vector<> global variable? function pointers?

evilcman (2)

I have this simple code, it should search for a root of function f with the bisection method:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
 
#include <iostream>
#include <vector>
#include <cstdlib>
#include <cmath>

using namespace std;

int n=30;

vector<double> lambda;

double f(double alpha)
{
  double sum;

  for(int i=0;i<n;i++){
    sum += lambda[i]/(lambda[i]+alpha);
  }

  return -0.1*alpha+sum;

}

double bisect(double * a, double * b , double (*f)(double)){
  double c;

  if(f(a[0]) * f(b[0])>=0.0){
    cerr << "bisect error: f(a)*f(b)>=0.0\n";
    return -123456789;
  }

  int itermax=20;
  double EPS=1e-5;
  cout << "iter\n";
  for(int i=0;i<itermax;i++)
  {
    c=((a[0])+(b[0]))/2.0;

    cout << " f(" << a[0] << ")= " << f(a[0])
         << " f(" << b[0] << ")= " << f(b[0])
         << " f(" << c << ")= " << f(c) << endl;

    if(abs(f(c))<EPS){ break;}
    else{
      if(f(a[0])*f(c)<0.0){
      b[0]=c;}
    else{
      a[0]=c;}
    }
 }

 return c;
}

int main(void){

for(int i=0;i<n;i++){
double l = 100((double)rand()/((double)RAND_MAX)-1.0);
lambda.push_back(exp(l));
cout << "lambda[" << i <<"]= " << exp(l) << endl;
}

double a=1e-5;
double b=1e5;
double x;
cout << "bisect\n";
x = bisect(&a,&b,&f);
cout << "\n\nsolution= " << x << endl;

return 0;

}


Here are two lines from the output:
f(24.4141)= 1.61311 f(48.8281)= -0.864139 f(36.6211)= 0.338324
f(36.6211)= -3.51471 f(48.8281)= -4.75959 f(42.7246)= -4.16748
It looks like the function f is not deterministic! Really strange.

If I change to code so that lambda is an additional argument of both the f and bisect functions than the program finds the root correctly. But, this is part of a bigger code, and it would really help me if i could use lambda as a global variable, and I am seriously out of ideas on why this doesn't work as I think it should.

I also tried checking if something changes lambda by writing to stdout the components at each call of f. They remain unchanged.

What the hell is the problem?
Last edited on
viliml (791)
try to replace line 14 with double sum=0.0;. I think that the hell is the problem.
evilcman (2)
You're right. I'm an idiot... :S
viliml (791)
:D
Topic archived. No new replies allowed.