Re: FindMinimum function does not seem to work for simple problems
- To: mathgroup at smc.vnet.net
- Subject: [mg128106] Re: FindMinimum function does not seem to work for simple problems
- From: vassilios.s.vassiliadis at gmail.com
- Date: Sun, 16 Sep 2012 03:23:54 -0400 (EDT)
- Delivered-to: l-mathgroup@mail-archive0.wolfram.com
- Delivered-to: l-mathgroup@wolfram.com
- Delivered-to: mathgroup-newout@smc.vnet.net
- Delivered-to: mathgroup-newsend@smc.vnet.net
- References: <k2s2qq$dq1$1@smc.vnet.net> <k2ub00$h9b$1@smc.vnet.net>
The problem arises from the discretization of a simple optimal control problem.
The description is as follows:
dv/dt = u
dx/dt = v
x(0)=0, v(0)=0,
x(tf)=300, v(tf)=0,
minimize tf.
-2<=u<=+1
The solution for this analytically is 30 time units as the minimum value of tf.
I used implicit Euler steps (backward finite differences) over a few steps to try to solve it with FindMinimum to test it (before I linked in larger optimal control problems). The performance is very very poor, unjustifiably so for a code that Wolfram advertise as "industrial strength". I tested this problem on CONOPT through GAMS and has no problems converging whatsoever for even very large numbers of discretization steps.
Here is the set of inputs to Mathematica. Please note that Mathematica 8.0 does not allow you to put a linear constraint with a single variable equated to a constant value (I note that others have reported this bug in forums such as this). I had to introduce a dummy variable "eps" equated to 1 indirectly to handle this problem.
=== BEGIN MODEL ==================
nsteps = 13;
vlist = Table[v[i], {i, 0, nsteps}];
xlist = Table[x[i], {i, 0, nsteps}];
hlist = Table[h[i], {i, 1, nsteps}];
ulist = Table[u[i], {i, 1, nsteps}];
varslist = Join[vlist, xlist, hlist, ulist, {eps}];
obj0 = Sum[h[i], {i, 1, nsteps}] + 10.^3*(eps - 1.0)^2;
odes = Flatten[
Table[{x[i] - x[i - 1] == v[i]*h[i],
v[i] - v[i - 1] == h[i]*u[i]}, {i, 1, nsteps}]];
ubounds = Flatten[Table[{u[i] <= 1.0, u[i] >= -2.0}, {i, 1, nsteps}]];
hbounds = Flatten[Table[{h[i] <= 100.0, h[i] >= 1.0}, {i, 1, nsteps}]];
(* these boundary conditions have a problem... equality of a single variable to a fixed value does not seem to work with FindMinimum as a valid constraint for some reason *)
boundary = {x[0]*eps == 0.0, v[0]*eps == 0.0, x[nsteps]*eps == 300.0, v[nsteps]*eps == 0.0} ;
constr0 =
Join[odes, ubounds, hbounds, boundary, {eps <= 1.000001, 0 <= eps}];
nvars0 = Dimensions[varslist][[1]];
vars0 = Transpose[{varslist, Table[1.0, {nvars0}]}]; vars0 = Join[Take[vars0, nvars0 - 1], {{eps, 1.0}}];
{objopt, solopt} = FindMinimum[
{obj0, constr0}, vars0, Method -> "InteriorPoint"
]
=== END MODEL ===================
The output that Mathematica gives after a very long run is the following:
FindMinimum::eit: The algorithm does not converge to the tolerance of 4.806217383937354`*^-6 in 500 iterations. The best estimated solution, with feasibility residual, KKT residual, or complementary residual of {0.000289013,0.0284091,2.28413*10^-11}, is returned. >>
I note that CONOPT gives a solution in a fraction of a second, e.g. milliseconds as it is really designed to solve problems with thousands of variables. If we believe what Wolfram are saying on their introduction to "Contrained Optimization", then the FindMinimum function implements a state-of-the-art Interior Point method and should be able to handle very easily large sparse problems,pretty much like the one I used to test it.
I note that the problem is essentially linear, were it not for bilinear terms in the equalities. Although this will give rise to many local solutions (the problem is nonconvex) it should not be hard to solve at all. This is proven at least by the use of one more commercial large scale solver like CONOPT...
Any ideas as to how it might be possible to still use FindMinimum are most welcome. Alternatively, I would be forced to consider how to link a serious solver like IPOPT from the COIN-OR project into Mathematica --- which I'd like to avoid and use totally native routines from Mathematica...
With sincere regards,
Vassilis