Question: odd fsolve instabilities

 

 

So many jags! Although it is not immediately clear from the above plot, the maximum solution and the instability both peak near 4.69

 

plot(P,4.6887..4.6892,smartview=false);

 

This is not a matter of precision -- the same results are achieved with plotting over 1000 points with Digits 100.

 

Why 4.69-ish?

 

plot(exp(-(1/100)*t^2)*cos(2*t), t = 4 .. 5, smartview = false);

 

Ah, because there is a minimum there, near 4.6889613332 . So we can imagine that secant-type methods might overshoot really badly when you are very close to the minimum. So the jumps very close to the mimumum are somewhat understandable.

 

What I do not understand, though, is the discontinuities further out from there, such as near 4.5.

 

plot(P, 4.515 .. 4.525, smartview = false);

The discontinuity is as large as any near the minimum, so something is leading it to overshoot by the same degree, and yet the slope is practically linear in the region:

 

plot(diff(exp(-(1/100)*t^2)*cos(2*t), t), t = 4.515 .. 4.525);

 

My speculation at the moment is that this has to do with the algorithm for selecting points "nearby" to evaluate function at. If trying to solve for P(x) involves evalating P(x-deltax) and P(x+deltax) then is would follow that if you happened to start at x-deltax then the forward projection from P(x-deltax) would be to P(x) leading to the same slope instabilities -- that there would be a kind of "ringing" near the minimum due to the algorithms for forward and backwards projection.

If my speculation is correct, then possibly a way to avoid this in the code would be to project ahead by a different distance than projecting back, that if it were evaluated at P(x-deltax_back) and P(x+deltax_forward) then when one were trying to solve at P(x-deltax_back) the projection ahead would be to P(x-deltax_back+deltax_forward) which would not project ahead to the minimum unless the two deltax were the same.