Wingopt add actuator constraints

[avikde]

• add them
• do a small timestep simulation with the optimized knot points - would help the animation as well...
• as another optimization method, could use gradient descent in the objective+penalty (as before), followed by a QP solve for constraint satisfaction?

TODO:

• anim
• fix newton method - Hessian singular??
• inequality constraints
• analytical gradients? confirm that evaluating the autograd is the slowest
• periodicity constraint (not needed?)

For a lot of this, the Coros work in Geilinger seems quite suited.

• use of penalty based method
• handling inequality constraints using psi(x)
• solution method: newton's method http://www.stat.cmu.edu/~ryantibs/convexopt-S15/lectures/14-newton.pdf followed by backtracking line search https://en.wikipedia.org/wiki/Backtracking_line_search
• Should I write a whole solution method based on this? May make sense if it solves the problem...

Result

Alternate

• gauss newton for traj with dyn penalty 1e-3
• newton method for params with dyn penalty 1e3

[3.6155 0.0156 0.127  5.25   0.0027 0.0149] cost -132.9 -> -359.8
[0.     0.0293 0.1349 0.3796 0.0001 0.0575] cost 12390945.0 -> 9651572.7
[3.5883 0.0142 0.1241 5.2144 0.0007 0.0147] cost -336.1 -> -706.0
[0.     0.0294 0.1475 0.3781 0.0001 0.2003] cost 39636777.4 -> 13562140.0
[3.499  0.0143 0.1457 5.1834 0.0006 0.0149] cost -272.8 -> -360.1
[0.     0.0296 0.1307 0.3702 0.0001 0.0575] cost 10592822.2 -> 4832206.4
[3.5135 0.0151 0.1265 5.2246 0.0006 0.0151] cost -244.0 -> -286.1
[0.     0.0293 0.1313 0.3664 0.0001 0.0569] cost 7393967.0 -> 3723089.0
[array([0.005, 1.5  ]), array([0.0054, 1.3946]), array([0.0043, 1.0109]), array([0.0039, 0.8853]), array([0.0037, 0.8118])]

Observations:

• Looks like params converge
• the cost from the traj opts seem to not go monotonically down, so there is a problem caused by this alternation (especially with the different penalties).

Future work

• needing to switch penalties motivates using a better method - augmented lagrangian?
• timestep in the opt
• other TODOs from Noah
• Julia implementation?

[avikde]

With no input limit penalty

With penalty of 1e3

[avikde]

Using (traj, params) as the opt vec results in a garbage result

update the fact that this does not work naively is good. this means that being clever about design optimization might be needed. convex wrt avg performance?

[avikde]

yellow is only lift, green is lift and drag (+x direction)

[avikde]

Gauss newton implementation worked and reduced the run time from 15s to 7s. (first term now includes jacobian of r)

	Jr	J0	DJ0	D2J0	Step
Gauss Newton	3.2427	0.015	0.1125	4.0763	0.0279
Newton	0	0.0275	0.2949	13.5016	0.0314

Still, the most time-consuming calculation is the Hessian evaluation.

The remaining non-quadratic terms include the objective and the inequality constraints

[avikde]

With the new Gauss Newton solver the results are actually better!

However, there is still some kind of non-convexity with both G and cbar...

[3.3169 0.014  0.1152 4.9416 0.003  0.0146] cost -147.9 -> -221.2
[0.     0.0288 0.1281 0.352  0.     4.262 ] cost -219.2 -> -334.2
[3.3161 0.0139 0.1167 4.9468 0.0007 0.0147] cost -483.0 -> -1055.7
[0.     0.0306 0.1399 0.3667 0.     0.0279] cost -941.3 -> -26192793084794296.0
[array([0.005, 1.   ]), array([0.0108, 1.0001]), array([110784.6877,   1763.6069])]

[avikde]

confirmed: nonconvexity wrt params

Looking at the "slices" wrt the params:

Looks ~quadratic with T, ~linear with cbar?

[avikde]

Evaluating dydt with the initial params and the params after returns these at t=0

[   -765.5635      23.9391 -727110.7488    9699.3942] [-7.6556e+02  2.3939e+01  2.6673e+09 -5.6170e+02]

i.e. the dynamics constraint be being violated. Perhaps the penalty is too small?

update:

yes, that was the problem. However, now, it seems like the dynamics constraint is "too strong". It may not be possible to update the params without also updating the traj.

Also did not make progress when trying to update both at the same time.

[avikde]

Unfortunately the results are quite similar when optimizing just the traj starting from different params: there is a nonconvexity
Dyn penalty 1e-3

Higher dynamics penalty 1e3

1e6

conclusion

This is proving that #56 (comment) is actually wrong. The dynamics penalty is fine and was being satisfied. The problem is that this nonconvexity truly exists. Need to add constraints or missing some modeled effects.

[avikde]

This is the simplest representation of the problem. f() changes as a function of T, and hence should be causing dynamics violations, but J seems to be going down with higher T

[avikde]

Finally an optimal T!

Optimum agrees for different values of the penalty (good)

also optimum for ci

Wrt both:

Also holds on an optimized traj