richard
User
Posts: 23
Joined: 8 months ago

Hi,
I am triyng to calculate the route wich takes less time to go from point A to point B in the presence of a constant flow (I. E. a simple version of Zermelo's navigation problem) using the GAMS software. However, if I put both points on a straight line and make the constant flow equal to zero, I don't get the results that I should.

Code: Select all

``````    \$set n 10
set j /0*%n%/;
sets
jlast(j)
jnotlast(j);
jlast(j)\$(ord(j)=card(j))=yes;
jnotlast(j)=not jlast(j);
scalar

v aeroplane velocity /10/

u wind velocity in the x axis /0/

n number of intervals /%n%/

density density /1.225/
variable
gamma(j),
x(j),
y(j),
time,
objective;

positive variable
t(j)

step;

equation
diffx(j),
diffy(j),
obj;

* diffx[j]\$(jnotlast(j)).. x[j+1]-x[j] =e=0.5*(t(j+1)-t(j))*(v*cos(gamma(j+1)) + v*cos(gamma(j)) );
* diffy[j]\$(jnotlast(j)).. y[j+1]-y[j] =e=0.5*(t(j+1)-t(j))*(v*sin(gamma(j+1)) + v*sin(gamma(j)) );

diffx[j]\$(jnotlast(j)).. x[j+1]-x[j] =e=0.5*step*(v*cos(gamma(j+1))-u + v*cos(gamma(j))-u );
diffy[j]\$(jnotlast(j)).. y[j+1]-y[j] =e=0.5*step*(v*sin(gamma(j+1)) + v*sin(gamma(j)) );

obj.. time =e= n*step;

x.fx('0') = 1.0e-12;
x.fx('%n%') = 1.0e-12;
y.fx('0') = 1.0e-12;
y.fx('%n%') = 10;
t.fx('0')= 1.0e-12;

y.up(j) = 10;
y.lo(j) = 1.0e-12;
t.lo(j)=1.0e-12;
gamma.up(j)=pi;
gamma.lo(j)=0;

model brahstron1 /all/;

option
nlp=ipopt;
solve brahstron1 using nlp minimize time;
``````
In this code I have defined two diferential equations, for x'=dx/dt and y'=dy/dt using the trapezoid rule. Here, gamma is our control variable and that's why it's free. Step is the difference between the times j and j+1, that is why time is n (number of intervals)*step(size of interval) and is what I seek to reduce.

Below this one can find the restrictions, these include the initial and final points for x and y, also the initial t=0.

I would spect as a result that the angle gamma became pi/2 and y going between 0 and 10, however, I find weird results, as Y being 5 on all the steps.

I have looked at this code for a while and I cannot find where my mistake is.

Can somebody please tell me what am I doing wrong?
Thank you very much for reading.