[Ipopt] Infeasible Starting Point

[Darcy Allison]

  Here is my AMPL input file. It is a simple orbit optimization. The 
file was written to find the optimal transfer between two circular 
orbits which should be the Hohmann transfer. It is very sensative to the 
initial guess for the eccentricity parameter, ecc. It seems like if it 
is below a certain value (around 0.7) it will not converge. However, if 
above the value it converges very quickly and accurately. The code below 
does not start in a feasible region (ecc is too small). If you change 
the initial ecc to something in the feasible region (0.76 >= ecc <= 1) 
it will converge.

I would like to use Ipopt for much harder problems but I need a code 
that can also handle infeasible starting points because I will not know 
beforehand where the feasible space is for these harder problems.

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

#instruct ampl that the IPOPT executable is the solver - make sure IPOPT 
is in the path
option solver ipopt;

#declare design variables and their bounds
var ecc >= 0, <= 1;
var p >= 0;
param r1 = 6531000.0;
param r2 = 48571000.0;
param mu = 398600500000000.0;

#specify the objective function delv=delv1+delv2
#delv1 = v1^2 + vc1^2 - 2*vc1*sqrt(mu*p)/r1
#delv2 = v2^2 + vc2^2 - 2*vc2*sqrt(mu*p)/r2
minimize obj:
             (mu*((2.0/r1) + ((ecc^2-1.0)/p)) + mu/r1 - 
2.0*((mu/r1)^0.5)*(((mu*p)^0.5)/r1))^0.5 + (mu*((2.0/r2) + 
((ecc^2-1.0)/p)) + mu/r2 - 2.0*((mu/r2)^0.5)*(((mu*p)^0.5)/r2))^0.5;

#specify the constraints
s.t.
         rp:
             r1*(1+ecc) >= p;

         ra:
             r2*(1-ecc) <= p;

#specify the starting point
let ecc  := 0.55;
let p     := 10000000.0;

#instruct ampl to solve the problem
solve;

#print out the solution
display ecc;        #eccentricity
display p;            #semilatus rectum
display r1;            #radius of first circular orbit
display r2;            #radius of second circular orbit
display mu;            #Earth's gravitational constant
display p/(1-ecc^2);#semi-major axis
display (mu*((2.0/r1) + ((ecc^2-1.0)/p)))^0.5;                        
                             #v1
display (mu*((2.0/r2) + ((ecc^2-1.0)/p)))^0.5;                        
                             #v2
display (mu*((2.0/r1) + ((ecc^2-1.0)/p)) + mu/r1 - 
2.0*((mu/r1)^0.5)*(((mu*p)^0.5)/r1))^0.5;    #delta-v1
display (mu*((2.0/r2) + ((ecc^2-1.0)/p)) + mu/r2 - 
2.0*((mu/r2)^0.5)*(((mu*p)^0.5)/r2))^0.5;    #delta-v2
#now display total delta-v which is the sum of delta-v1 and delta-v2
display (mu*((2.0/r1) + ((ecc^2-1.0)/p)) + mu/r1 - 
2.0*((mu/r1)^0.5)*(((mu*p)^0.5)/r1))^0.5 + (mu*((2.0/r2) + 
((ecc^2-1.0)/p)) + mu/r2 - 2.0*((mu/r2)^0.5)*(((mu*p)^0.5)/r2))^0.5;
display 3.14159265359*(((r1+r2)^3)/(8*mu))^0.5;                    #time 
to complete the transfer

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

Regards.
Darcy


On 3/18/2011 11:01 PM, Ashutosh Mahajan wrote:
> Dear Darcy,
>
> It may be more helpful to describe the problem you are solving and how you
> input it to Ipopt. Do you call it through AMPL, GAMS or the C++ interface?
> Are there any points where the functions become infinite or nan (like 1/x)?
>
> --
> regards
> Ashutosh Mahajan
> http://coral.ie.lehigh.edu/~asm4
>
> On Thu, Mar 17, 2011 at 10:24:06AM -0400, Darcy Allison wrote:
>> Hi. I recently used Ipopt for a nonlinear optimization problem and
>> discovered that if I did not start in the feasible region that it
>> crashed and would not converge to the optimum. However, if I started
>> in the feasible region it would find the optimum just fine. I am
>> wondering if this deficiency is currently being worked on or not?
>> Thanks.
>>
>> _______________________________________________
>> Ipopt mailing list
>> Ipopt at list.coin-or.org
>> http://list.coin-or.org/mailman/listinfo/ipopt