Tighten tolerance for first iteration in solve1()
Currently solve1() uses the same tolerance level (tolf) to decide to start Newton iterations or to terminate them. I suggest the the tolerance above which to start Newton iterations should be tighter.
It is mainly useful for the simulation of purely backward models. When the model is stationary, we expect the trajectory to converge asymptotically towards the steady state. This occurs as long as at least one iteration Newton takes place in each period. If the tolerance criteria is too lax, the simulation remains constant from too early a point.
I suggest to use the square of tolf at the start of the computation
For reference, csolve() always performs one Newton iteration. Matlab's fsolve() and lmmcp use a single tolerance level.