discrete time optimization with non-linear constraints

Hello,

I have an optimization problem to solve with non-linear constraints. It is a control theory based discrete time model (which i feel fules out using fmincon) over a time horizon say N seconds.

I found a few old posts similar to this, but none I felt had clear answers.

I am going nuts trying to find a way to implement it in matlab. Which Matlab tool would be best suited in this case?

Thanks,

italic EDIT: the system is continuous, but we analyse it in a discrete time domain. Thus, variables have discrete values. There are n_v entities and each of the entity has each parameter described below:

variables: p,v,u; size(p)=size(v)=size(u)=(1,N) vectors

Obj. fn: minimize sum(u(1,:))

constraints:-

Upper and lower bounds like:

p_min<= p(1,:)<= p_max

v_min<= v(1,:)<= v_max

u_min<= u(1,:)<= u_max

linear equalities:

p(n+1)=p(n)+ v(n)t + 0.5u(n)t^2;

v(n+1)=v(n)+u(n)t;

u(n+1)= (s_star / (p1-p2))^2 ; p1 and p2 correspond to entity 1 and 2's p parameter

non linear equalities:

s_star= v1(n) * (v1(n)-v2(n)) ; where v1 and v2 are parameters of entities 1 and 2 respectively at nth instant. * this is the non-linear equality/constraint.

s_star2 = v1(n) * (v1(n)-v1(n-1));

Note: we use either s_star or s_star2 in our code

initial conditions:

p(1)=p_init;

v(1)=vel_init;

end conditions:

p(N)<= p_max;

v(N)=0;

I can give more information if required.

thanks,

[Matt J]

If there's nothing in the Control Systems Toolbox, then write down the mathematical details of the problem and we can say more. The fact that it is discrete time doesn't automatically rule out fmincon. The real question is whether the unknown variables are discrete-valued or continuous-valued.

[Torsten]

You will have to describe the problem in more detail in order to get help.

Best wishes

Torsten.

Answer 1

[Tohru Kikawada]

I'm not sure what you want to get, but I tried to create an example. Hope this helps.

function [x,f,eflag,outpt] = myModelParamsSolver

xLast = []; % Last place computeall was called
myf = []; % Use for objective at xLast
myc = []; % Use for nonlinear inequality constraint
myceq = []; % Use for nonlinear equality constraint
N = 10;
x0 = rand(1,2+2*N);
fun = @objfun; % the objective function, nested below
cfun = @constr; % the constraint function, nested below
options = optimoptions('fmincon');
lb = zeros(1,2+2*N);
ub = ones(1,2+2*N);
% Call fmincon
[x,f,eflag,outpt] = fmincon(fun,x0,[],[],[],[],lb,ub,cfun,options);

      function y = objfun(x)
          if ~isequal(x,xLast) % Check if computation is necessary
              [myf,myc,myceq] = myDiscreteTimeModel(x);
              xLast = x;
          end
          % Now compute objective function
          y = myf;
      end

      function [c,ceq] = constr(x)
          if ~isequal(x,xLast) % Check if computation is necessary
              [myf,myc,myceq] = myDiscreteTimeModel(x);
              xLast = x;
          end
          % Now compute constraint functions
          c = myc; % In this case, the computation is trivial
          ceq = myceq;
      end

      function [f,c,ceq] = myDiscreteTimeModel(x)
          p = zeros(1,N);
          v = p;
          u = p;
          Ts = 0.01;
          p_max = 100;
          p_init = 10;
          vel_init = 1;

          % initial conditions:
          p(1)=p_init;
          v(1)=vel_init;
          p1 = x(1);
          p2 = x(2);
          v1 = x(3:3+N-1);
          v2 = x(3+N:end);

          for n = 1:N
              t = N*Ts;

              % non linear equalities:
              s_star= v1(n) * (v1(n)-v2(n));
              % where v1 and v2 are parameters of entities 1 and 2 respectively at nth instant.
              % * this is the non-linear equality/constraint.
              %s_star2 = v1(n) * (v1(n)-v1(n-1));
              % Note: we use either s_star or s_star2 in our code

              % linear equalities:
              p(n+1)= p(n)+ v(n)*t + 0.5*u(n)*t^2;
              v(n+1)= v(n)+ u(n)*t;
              u(n+1)= (s_star / (p1-p2))^2 ; % p1 and p2 correspond to entity 1 and 2's p parameter
          end

          % end conditions:
          ceq = v(N);
          c = p(N)-p_max;
          f = sum(u(:));
      end
  end