| double mGRF(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) ) |
| { |
| fA=(*func)(A); |
| 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; |
| fB=(*func)(B); dB=((*func)(B+h)-(*func)(B))/h; |
| a =(fB-fA+A*dA-B*dB)/(dA-dB); // moving fixed point |
| fa=(*func)(a); da=((*func)(a+h)-(*func)(a))/h; |
| if ((a>A) && (a<B)) // bracketing moving point |
| if ((fa*fA<0) || (da*dA<0)) |
| B=a; |
| else |
| A=a; |
| } |
| else if ((fX*fB<0) || (dX*dB<0)) |
| { |
| A = X; |
| fA=(*func)(A); dA=((*func)(A+h)-(*func)(A))/h; |
| b =(fA-fB+B*dB-A*dA)/(dB-dA); // moving fixed point |
| fb=(*func)(b); db=((*func)(b+h)-(*func)(b))/h; |
| if ((b>A) && (b<B)) // bracketing moving point |
| if ((fb*fB<0) || (db*dB<0)) |
| A=b; |
| else |
| B=b; |
| } |
| |
| } |
