Hi all,
I have encountered a peculiar problem when performing some simple arithmetic in a program I am developing. As the function of the program is pretty specific to my field I have included a smaller program which illustrates the same problem without the technicalities.
Given a function (in the mathematical sense of the word) e.g. f(x), and an interval on which the function is defined e.g. x1 to x2, the program divides this interval into n segments of equal width and searches for a zero crossing of the function, i.e. a point on the interval where there is a solution to the equation f(x) = 0. The program loops over each segment incrementing x from x1, evaluating f(x) as it goes. It stops when it detects a zero crossing and outputs the refined interval on which the solution lies. Here is the program...
//...problem code
#include <iostream>
#include <cmath>
#include <vector>
using namespace std;
int main(){
int n = 1000; // number of segments
double x, x1, x2, dx, fp, fc, xb1, xb2;
x1=-0.1; x2=0.1; dx=(x2-x1)/n; // inteval defined by x1 and x2, dx is width of segments
fp=erf(x=x1); // error function defined in cmath
for(int i=0; i<n; i++){ // loop over intervals
fc=erf(x+=dx); // determine incremented value of erf()
if(fc*fp<=0.0){ // sign change indicates root of function
xb1=x-dx; // record bound and break from loop
xb2=x;
break;
}
}
cout << xb1 << " " << xb2 << endl;
}In this case f(x) is the erf() function defined in cmath, this has the property erf(0)=0. The problem arises when n is defined such that the program encounters x=0 (as above), when this happens x-dx is not 0.0 as I would expect, but is a very small value of order 1e-16 and whose exact value is dependent on the value of n. Any help on this matter would be very much appreciated.
Best wishes
