Unable to meet integration tolerances in fmincon to solve ODE optimization problem

[poster]

I have constrained ODE optimization problem to be solved using Matlab, I started by solving the ODE using ode15s which working will this initial values, also I had very good results without constraints using fminsearch, the problem started when I used fmincon it gave me these warnings:

Warning: Length of lower bounds is < length(x); filling in missing lower bounds with -Inf. 

> In checkbounds (line 33)

  In fmincon (line 318)

  In Optimization (line 64) 

Warning: Length of upper bounds is < length(x); filling in missing upper bounds with +Inf. 

> In checkbounds (line 47)

  In fmincon (line 318)

  In Optimization (line 64) 

Warning: Failure at t=2.340250e+01.  Unable to meet integration tolerances without reducing the step size below the smallest value

allowed (5.684342e-14) at time t. 

> In ode15s (line 668)

  In ObjFun (line 6)

  In barrier

  In fmincon (line 813)

  In Optimization (line 64) 

I tried to remove constraints, but nothing changed...

y1_max = 5; y2_max = 2.3; 

y01 = 1;    y02 = 1.6;



%Control function parameter

Umin = -0.3; Umax = -0.1;



%Creation the structure of model parameters

params.T0=T0; params.Tf=Tf;

params.y1d=y1d; params.y2d=y2d;

params.y01=y01; params.y02=y02;

params.y2_max=y2_max;



n=2;

U0=-0.2*ones(1,n);

params.n=n;

params.U=U0; params.Umax=Umax; params.Umin=Umin;

params.dt=(Tf-T0)/n;



%for initial value of optimization parameters

options = odeset('RelTol',1e-7,'AbsTol',1e-7);

U=U0; 



[t,y]=ode15s(@ODEsolver,[T0,Tf],[y01 y02],options,params);



%adding the ode solution as input prameters

params.t=t; params.y1= y(:,1); params.y2= y(:,2);



U0=-0.2*ones(1,n);

params.n=n; params.U=U0;

params.dt=(Tf-T0)/n;



A = [];

B = [];

Aq = []; Bq = [];

options1 = optimset('MaxIter',5000);

[x,fval,exitflag,output] = fmincon(@ObjFun,U0,A,B,Aq,Bq,Umin,Umax,IneqConst(U0,params),options1,params)

function y = ODEsolver(t,y,params)

    dt = params.dt;

    nu = floor(t/dt)+1;

    nu = min(nu,params.n-1);



    t1 = (nu-1)*dt;

    Ut = params.U(nu) + ((t-t1)/dt)*(params.U(nu+1) - params.U(nu));



    dy1a = ((0.63*y(1))/(1 + y(1)));

    dy1b = (Ut*y(1)*y(2))/(1+0.6*y(1));

    dy1 = dy1a + dy1b;



    dy2a = (0.15*y(2))/(0.05-0.2*y(2));

    dy2 = -0.2*y(2) + dy1b * dy2a;



    y = [dy1; dy2];

end

function ObjFun1 = ObjFun(u,params)

   % Calculating the value of the optimization criteria 

   params.U=u;



   options = odeset('RelTol',1e-8,'AbsTol',1e-8);

   [t,y]=ode15s(@ODEsolver,[0,100],[params.y01, 

   params.y02],options,params);



   ObjFun1 =trapz(t,(y(:,1)-params.y1d).^2);

end

function [c,Const] = IneqConst(u, params)

    params.U=u;

    options = odeset('RelTol',1e-8,'AbsTol',1e-8);

    [t,y]=ode15s(@ODEsolver,[0,100],[params.y01, params.y02],options,params);



    c =[];

    yCon = y(:,2)-params.y2_max;

    Const = trapz(t,(abs(yCon)+yCon).^2);

end

More answer...