csu-hmc / opty

A library for using direct collocation in the optimization of dynamic systems.
http://opty.readthedocs.io
Other
92 stars 20 forks source link

ackermann2010 example fails #20

Open moorepants opened 9 years ago

moorepants commented 9 years ago

The current state of affairs for the ackermann example is shown below.

Some things to try:

(inverted-pendulum-id)moorepants@moorepants-2170p:ackermann2010(master)$ ipython
Python 2.7.8 |Continuum Analytics, Inc.| (default, Aug 21 2014, 18:22:21) 
Type "copyright", "credits" or "license" for more information.

IPython 2.3.0 -- An enhanced Interactive Python.
Anaconda is brought to you by Continuum Analytics.
Please check out: http://continuum.io/thanks and https://binstar.org
?         -> Introduction and overview of IPython's features.
%quickref -> Quick reference.
help      -> Python's own help system.
object?   -> Details about 'object', use 'object??' for extra details.

In [1]: run ackermann2010.py 
Forming positions, velocities, accelerations and forces.
Initializing Kane's Method.
Forming Kane's Equations.
Printing ccode.
Computing the symbolic partial derivatives.
Printing ccode.

******************************************************************************
This program contains Ipopt, a library for large-scale nonlinear optimization.
 Ipopt is released as open source code under the Eclipse Public License (EPL).
         For more information visit http://projects.coin-or.org/Ipopt
******************************************************************************

This is Ipopt version 3.11, running with linear solver ma57.

Number of nonzeros in equality constraint Jacobian...:    44618
Number of nonzeros in inequality constraint Jacobian.:        0
Number of nonzeros in Lagrangian Hessian.............:        0

Reallocating memory for MA57: lfact (453385)
Reallocating memory for MA57: lfact (503593)
Total number of variables............................:     1440
                     variables with only lower bounds:        0
                variables with lower and upper bounds:     1440
                     variables with only upper bounds:        0
Total number of equality constraints.................:     1070
Total number of inequality constraints...............:        0
        inequality constraints with only lower bounds:        0
   inequality constraints with lower and upper bounds:        0
        inequality constraints with only upper bounds:        0

iter    objective    inf_pr   inf_du lg(mu)  ||d||  lg(rg) alpha_du alpha_pr  ls
   0  2.3742185e+07 2.72e+05 5.03e+00   0.0 0.00e+00    -  0.00e+00 0.00e+00   0
WARNING: Problem in step computation; switching to emergency mode.
   1r 2.3742185e+07 2.72e+05 9.99e+02   3.9 0.00e+00  20.0 0.00e+00 0.00e+00R  1
WARNING: Problem in step computation; switching to emergency mode.
Restoration phase is called at point that is almost feasible,
  with constraint violation 0.000000e+00. Abort.
Restoration phase in the restoration phase failed.

Number of Iterations....: 1

                                   (scaled)                 (unscaled)
Objective...............:   2.3742184587421797e+07    2.3742184587421797e+07
Dual infeasibility......:   5.0338983554237284e+00    5.0338983554237284e+00
Constraint violation....:   8.6174364489984910e+03    2.7227023016135045e+05
Complementarity.........:   5.9700000596999996e+02    5.9700000596999996e+02
Overall NLP error.......:   8.6174364489984910e+03    2.7227023016135045e+05

Number of objective function evaluations             = 2
Number of objective gradient evaluations             = 2
Number of equality constraint evaluations            = 2
Number of inequality constraint evaluations          = 0
Number of equality constraint Jacobian evaluations   = 2
Number of inequality constraint Jacobian evaluations = 0
Number of Lagrangian Hessian evaluations             = 0
Total CPU secs in IPOPT (w/o function evaluations)   =      2.998
Total CPU secs in NLP function evaluations           =      0.157

Timing Statistics:

OverallAlgorithm....................:      3.155 (sys:      7.384 wall:      3.090)
 PrintProblemStatistics.............:      0.000 (sys:      0.000 wall:      0.000)
 InitializeIterates.................:      0.332 (sys:      0.560 wall:      0.263)
 UpdateHessian......................:      0.000 (sys:      0.000 wall:      0.000)
 OutputIteration....................:      0.001 (sys:      0.000 wall:      0.000)
 UpdateBarrierParameter.............:      1.278 (sys:      2.345 wall:      0.952)
 ComputeSearchDirection.............:      0.000 (sys:      0.000 wall:      0.000)
 ComputeAcceptableTrialPoint........:      1.507 (sys:      4.392 wall:      1.835)
 AcceptTrialPoint...................:      0.000 (sys:      0.000 wall:      0.000)
 CheckConvergence...................:      0.036 (sys:      0.088 wall:      0.039)
PDSystemSolverTotal.................:      1.278 (sys:      2.345 wall:      0.952)
 PDSystemSolverSolveOnce............:      1.277 (sys:      2.345 wall:      0.952)
 ComputeResiduals...................:      0.000 (sys:      0.000 wall:      0.000)
 StdAugSystemSolverMultiSolve.......:      2.945 (sys:      6.892 wall:      2.859)
 LinearSystemScaling................:      0.000 (sys:      0.000 wall:      0.000)
 LinearSystemSymbolicFactorization..:      0.009 (sys:      0.000 wall:      0.002)
 LinearSystemFactorization..........:      2.875 (sys:      6.891 wall:      2.840)
 LinearSystemBackSolve..............:      0.000 (sys:      0.000 wall:      0.000)
 LinearSystemStructureConverter.....:      0.000 (sys:      0.000 wall:      0.000)
  LinearSystemStructureConverterInit:      0.000 (sys:      0.000 wall:      0.000)
QualityFunctionSearch...............:      0.000 (sys:      0.000 wall:      0.000)
TryCorrector........................:      0.000 (sys:      0.000 wall:      0.000)
Task1...............................:      0.000 (sys:      0.000 wall:      0.000)
Task2...............................:      0.000 (sys:      0.000 wall:      0.000)
Task3...............................:      0.000 (sys:      0.000 wall:      0.000)
Task4...............................:      0.000 (sys:      0.000 wall:      0.000)
Task5...............................:      0.000 (sys:      0.000 wall:      0.000)
Function Evaluations................:      0.157 (sys:      0.384 wall:      0.155)
 Objective function.................:      0.001 (sys:      0.000 wall:      0.000)
 Objective function gradient........:      0.004 (sys:      0.000 wall:      0.001)
 Equality constraints...............:      0.067 (sys:      0.192 wall:      0.075)
 Inequality constraints.............:      0.000 (sys:      0.000 wall:      0.000)
 Equality constraint Jacobian.......:      0.084 (sys:      0.192 wall:      0.078)
 Inequality constraint Jacobian.....:      0.000 (sys:      0.000 wall:      0.000)
 Lagrangian Hessian.................:      0.000 (sys:      0.000 wall:      0.000)

EXIT: Restoration Failed!

In [2]: pro
%profile  prob      property  

In [2]: prob.con(initial_guess)
Out[2]: 
array([ 914.98161943,  720.16048423, -974.78036008, ...,   -9.75610511,
         -3.85324192,   -8.17724079])

In [3]: prob.con_jac(initial_guess)
Out[3]: array([ 118.,    0.,    0., ...,    1.,    1.,    1.])

In [4]: prob.obj(initial_guess)
Out[4]: 150127264.38578781

In [5]: prob.obj_grad(initial_guess)
Out[5]: 
array([  0.        ,   0.        ,   0.        , ..., -10.56153419,
        -1.25221308,   4.8188792 ])

In [6]: %timeit prob.con_jac(initial_guess)
10 loops, best of 3: 21.8 ms per loop
moorepants commented 3 weeks ago

@Peter230655 the paper is this:

Ackermann, M., & van den Bogert, A. J. (2010). Optimality principles for model-based prediction of human gait. Journal of Biomechanics, 43(6), 1055–1060. https://doi.org/10.1016/j.jbiomech.2009.12.012

I'll email it to you.

I think the only thing that prevented me from solving this was a variable time step, which is now implemented.

I also built this package that implements a 2D walking model:

https://github.com/csu-hmc/gait2d/tree/master

So the model should already be sufficiently functional.