How to use fmincon effctively to calculate time-changing control inputs?

I have a dynamic system which is described by differential equations:

dot_X = f(X, U, P) % dot_ means dx/dt

X is state vector, X = [x1, x2, x3]; U is control inputs vector, U = [u1, u2]; P is intermediate variables vector, P = [p1, p2], and both p1 and p2 are functions of X and U, p1 = g1(X, U), p2 = g2(X, U); f is complicated nonlinear function.

Initial condition of the system is: X(t=0) = 0, U is unknown and going to be optimized to satisfy my objective function.

My objective is to maximize x1 when t = 2s, that is obj = x1(t=2).

There are constrains on the intermediate variables p1 and p2, which I hope for every time instant, for example, t = 0.01, 0.02, 0.03, ... , 2.98, 2.99, 3, p1 and p2 both are always within [-3, 3]. I still didn't achieve writing this constraints yet. Other constrains are directly related to control inputs u1 and u2, which I set low and up boundary in my fmincon function already.

Without properly writing constrains now (I pasted another question "How to constrain a vector when using fmincon?" yesterday), I wrote my fmincon like below:

[ In_opt ] = fmincon( @( Inputs_vector )obj( other parameters_1, Inputs_vector), [Inputs_vector_0], [], [], [], [], [Inputs_vector_L],[Inputs_vector_U],... @( )ctr( ), optNLP);

For the optimization option "NLP", I set it default since I have little math knowledge on optimal control.

I set my simulation time from 0 to 2 seconds, and ode time interval is 0.02s. The results are not good because it vibrates and not physically nice neither. Then I tried to make the ode time interval smaller, like 0.005s, then the program becomes very very time consuming. To my understanding, when the time interval is set to 0.005s, there are 2/0.005 = 400 points in the simulation time stream which fmincon tries to optimize, which took a lot of effort. But I'm not sure whether my understanding of it is right.

Could you let me know whether my method using fmincon in this question is right? And if yes, how should I make the fmincon more efficient to calculate my optimal control inputs between [0s, 2s]? And if no, could you let me known other method which can deal with my question?

Thank you very much, I would really happy if I can solve this problem. Waiting for your reply.

Answer 1

[Alan Weiss]

I am sorry that you did not seem satisfied by the answers to your previous question. I believe that my answer was correct and reasonably complete. I hope that you will revisit it.

I noticed that you evaluated your ODE in a loop of time steps. I believe that this is a time-consuming mistake. Instead, evaluate ODE45 over the entire interval at once, yielding a sol solution structure. Then evaluate sol at the time points you want using deval.

Your nonlinear constraint should be a function that yields a vector such as the following.

 function [c,ceq] = mycon(x)
 ceq = [];
 % compute the solution to the ODE here
 % create a vector of times that you want to examine
 yvector = deval(sol,times); % solution vector
 c = yvector(:) - 3; % upper bounds
 c = [c;-3 -yvector(:)]; % lower bounds

I hope that this makes sense to you. If not, please ask a specific question about the specific step that does not work for you.

Alan Weiss

MATLAB mathematical toolbox documentation

[Tao]

Hi, Alan, I'm still trying to tweaking my code with your answer from the previous question, since I have some problem in model physical aspect where the K could be far larger than 3 or below -3 for x1 becomes high, I'm thinking how to constrain the K properly.

[Tao]

Hi, Alan. Could you elaborate more on "evaluate ODE45 over the entire interval at once"?

Since my objective is to maximize the state x1 at final time instant (t = 2s), I used ode45 to integrate my differential equations from t=0 to t=2, thus I could get x1 at t = 2s. Otherwise how should I get my objective function value and how to "evaluate ODE45 over the entire interval at once"?

I'm new to MATLAB, so thank you for your patience in advance.

[Tao]

I used a for-loop when integrating my differential equations:

Time = 0:0.02:2;

for i = 1:1:(length(Time)-1)

Tspan = [Time(i) Time(i+1)];

[T, X] = ode45(@(t,X) dynamics(...), Tspan, IC);

end

[Tao]

Also since P is depending on X and U, X and U are changing with time, so P is changing with time too. Essentially, the dynamics is not ode, but with time changing intermediate variable P, that's why I use for-loop together with ode45 to get my x1 at time t=2s.

[Alan Weiss]

Of course, I don't know as much about your problem as you do. But IF you have an initial guess at your control variables x0 THEN you should be able to integrate your differential equation from time 0 to time 2 in one go, without using a loop at all.

In more detail, maybe I don't understand what you are really trying to do, but if what you want to do is optimize something, then you first need to decide what your control variables are, meaning the variables that you can tweak in order to try to find an optimum. For fmincon to start you need to give an initial guess at the control variables, x0. For a given, fixed value of the control variables, you should be able to evaluate your nonlinear objective and constraint functions. It is entirely possible that x0 leads to an infeasible solution, meaning one where the constraints are violated. But that is OK, just let fmincon do its best to make the solution feasible later. And from what I understand, you can evaluate your objective and constraints by an integration from time 0 to time 2. I am suggesting that you do this integration in one go, without a loop, and use deval to evaluate your nonlinear constraints.

Alan Weiss

MATLAB mathematical toolbox documentation

[Tao]

Thanks, Alen, please give me a moment to think how to deal the code with your instruction on the ode45. I will let you know my progress in time, and really thanks.