| #include <math.h> |
| #define ACY 1.0E–6 accuracy for estimates, function, and finite difference values |
| #define h 1.0E–7 infinitesimal for finite diferences |
| double GRF(float (*func)(double), double A, double B, int n) |
| { |
| double X, fX, dX, A, fA, dA, B, fB, dB, a, fa, da, b, fb, db; |
| while ( (fabs(B-A)>ACY) && (–n>0) ) |
| { |
| double fA=(*func)(A); |
| double dA=((*func)(A+h)-(*func)(A))/h; |
| if ((fabs(fA)< ACY) || (fabs(dA)< ACY)) // zero or extremum |
| return A; |
| fB=(*func)(B); |
| dB=((*func)(B+h)-(*func)(B))/h; |
| if ((fabs(fB)< ACY) || (fabs(dB)< ACY)) // zero or extremum |
| return B; |
| X=A+fabs(fA)/(fabs(fA)+fabs(fB))*(B-A); |
| fX=(*func)(X); |
| dX=((*func)(X+h)-(*func)(X))/h; |
| if ((fabs(fX)< ACY) || (fabs(dX)< ACY)) // zero or extremum |
| return X; |
| if ((fA*fX<0) || (dA*dX<0)) |
| B = X; |
| else |
| if ((fX*fB<0) || (dX*dB<0)) |
| A = X; |
| } |
| } |
