[SciPy-dev] scipy.optimize.nonlin.broyden2

Ondrej Certik ondrej@certik...
Wed Mar 19 17:18:43 CDT 2008


On Wed, Mar 19, 2008 at 9:30 PM, argriffi <argriffi@ncsu.edu> wrote:
> Hi Ondrej,
>
>  Running the attached script should demonstrate the problem I mentioned.
>  I removed dependencies on other modules I had written, but let me know
>  if you want it reduced even further.

Thanks, that is enough. I can reproduce the problem now:

------------------

test_solver_20 (__main__.TestCodonInversion) ... ok
test_solver_40 (__main__.TestCodonInversion) ... ERROR

======================================================================
ERROR: test_solver_40 (__main__.TestCodonInversion)
----------------------------------------------------------------------
Traceback (most recent call last):
  File "CodonInversion.py", line 181, in test_solver_40
    self._test_solver_helper(40)
  File "CodonInversion.py", line 173, in _test_solver_helper
    raise RuntimeError('failed on iteration %d of %d: %s' %
(len(objective_function.guesses), iterations, str(e)))
RuntimeError: failed on iteration 22 of 40: detected possibly invalid
input to the objective function: (nan, nan, nan)

----------------------------------------------------------------------
Ran 2 tests in 0.029s

FAILED (errors=1)

-------------------------


And the reason is simple. If you look at this code:


    for n in range(iter):
        deltaxm=-Gm*Fxm
        xm=xm+deltaxm
        Fxm1=myF(F,xm)
        deltaFxm=Fxm1-Fxm
        Fxm=Fxm1
        Gm=Gm+(deltaxm-Gm*deltaFxm)*deltaFxm.T/norm(deltaFxm)**2


And if you debug it, as I just did (for example using winpdb), you'll
notice, that deltaFxm gets zero, after a lot of iterations, because
you converged
to a solution already. Then norm(deltaFxm) is 0 and after dividing by
0, you get a NaN in Gm, thus you get NaN in deltaxm and xm and you get
NaN as an input into your function and here we are.

The solution is to check that deltaFxm is near zero and stop
iterating. Do you think you could please provide such a patch? :)

Also with a test case.

Thanks,
Ondrej

P.S. I CCed scipy-dev as I think you forgot to forward your email to
it. I hope it's ok.


More information about the Scipy-dev mailing list