Proper set up of fmincon including ODEs

I am actually trying to understand how the fmincon function is used in a little more complex way(including solving ODE). So I thought of a little example for which I imagine a vehicle on a 2D object, like a table. So the vehicle can move in x and in y direction and it is supposed to "fly" through a given point(xf,yf) on the table which I defined as my equality constrained. The inequality constraints are just a region on the table where the vehicle is not allowed to fly through, so that it doesn't just fly in a straight line through xf,yf. The vehicle should not leave the table so I defined the upper and lower bounds fir the tables geometry. Also right now they include Acceleration in x and y direction and the limits of the vehicle mass(so that the mass is never negative).

My vehicle has a thruster, propellant and therefore a propellant flow.

My objective of this little example is to write a code that returns me the optimal trajectory(so the alpha is the control variable), so that the vehicle reaches xf,yf with the max. propellant(mp) possible. Obviously this is a Maximization problem, that is why I used mp=-....

The problem right now is that I don't understand how my Optimizer can return an optimal value for mp at the final time(And not just the highest possible mass, which is obviously at the beginning) .

Are the other variables I have(alpha,x,y,ax,ay) all to be defined as one control variable in form of a vector and if so how and when do I need to pass them to my Optimizer, because again I just want to output an optimal mpfinal and then plot the flight trajectory and the change of alpha in another step out of the results.

I attached my recent code

Hope someone can help

[Torsten]

Maybe someone can help if you include a mathematical problem formulation:

min :

s.t.

...

[Brendan Görres]

my equations:

alpha=atan(y/x);
x''=thrust*cos(alpha)/m;
y''=thrust*sin(alpha)/m;
m'=-propflow%propflow=const      

In my code I transformed the 2nd order equations to first order equations, so that de45 can solve it.

Then my objective function is:

mfinal=-m0+(propflow*sqrt(2*X(1)*X(5)/(thrust*cos(alpha))));

where

X=[x;ax;y;ay;m];

My linear constraints are:

if y>=5e3
    x>7e3;
else
    x>0;
end

That is the restricted area my vehicle should not pass.(I just realized that I formulated them as nonlinear constraints in the attached code, but I am going to change that)

What I am trying to do is to get the highest possible mass at the final time.

I tried to change the code so that the initial guess for the code is just

m0=68e3;

And I changed lb uand ub to:

lb=0;
ub=70e3;

so these are the bounds for my mass.

But now the problem is, that my A,b,Aeq,beq define my constraints for the position of the vehicle(x and y) and the limits are only for the mass that should be optimized. And I don't know, how to pass these information to fmincon, because my constraints are for other variables than my ub and lb. So that would lead to a different array size.

And furthermore I don't know how my fmincon can optimize only the mass at a specific time and not just find the maximum, which is obvously at the beginning.

[Torsten]

I don't understand what your control variables are (i.e. the variables you can vary over time as to satisfy your constraints and to influence your objective). Are thrust and/or propflow functions of time ?

[Brendan Görres]

Thrust and propflow are constants, but I can change the direction of the thrust by changing the angle. therefore my control variable should be alpha.

[Torsten]

But you can't control alpha since it is set to atan(y/x), and x and y are determined by your ODEs.

[Walter Roberson]

Your angle would change with time in a nonsmooth way in order to accomodate the geometry changes along the path. The ode* routines require smooth changes. Therefore every time you change angle nonsmoothly you need to terminate the ode and restart it from that position.

I have a suspicion that the solution will correspond to a spline fit .

[Brendan Görres]

@Torsten

But if I leave alpha as a variable without the function

alpha=atan(y/x);

the ode45 cannot solve the equations. Or if I input alpha to fmincon as X0, how could I pass my constraints(that are only on x and y and not on alpha) to fmincon, because then the matrix for the linear inequality constraints has 2 rows and 2 columns while the X0 would just be alpha(0) and that is not solvable.

And that is exactly my problem. I do not know what Inputs I should give fmincon in this example, because I couldn't find anything similar while browsing..

[Torsten]

Optimal control problems contain functions as unknown objects, not single values.

So alpha = (alpha(1),alpha(2),...,alpha(n)) would be your X vector for fmincon to optimized. Here, alpha(i) corresponds to alpha at a certain time t_i between tstart and tend.

To use this vector alpha within your ODEs, you can use interpolation to the time t the integrator has reached.

Using this vector alpha, you get certain derived quantities, e.g. (x(1),...x(n)) and (y(1),...,y(n)). Constraints on these values can be set in function "nonlcon".

I strongly recommend that you consult a book on optimal control.