Why I get index out of bounds error even though I debug the loop and its Ok

Hi all I try to optimize three variables using fmincon. Inside the objective function there is a loop which the error message point on it. I debug the loop and its work and its reach the end and go to fmincon to complete the check after that I get the error Kindly, can any on explain why I get this error. Regards

??? Attempted to access t02(206); index out of bounds because numel(t02)=205. % t02 is a time vector (206x1)

Error in ==> ObjFun_Zdd_F at 52
        ttt=t02(i+1)-t02(i);

Error in ==> @(Fkc)ObjFun_Zdd_F(Fkc)

Error in ==> E:\Program
Files\MATLAB\R2011a\toolbox\optim\optim\private\evalObjAndConstr.p>evalObjAndConstr
at 135

Error in ==> E:\Program
Files\MATLAB\R2011a\toolbox\optim\optim\sqpLineSearch.p>sqpLineSearch at 265

Error in ==> fmincon at 832
    [X,FVAL,EXITFLAG,OUTPUT,LAMBDA,GRAD,HESSIAN] =
    sqpLineSearch(funfcn,X,full(A),full(B),full(Aeq),full(Beq), ...

Error in ==> Fmain at 9
[Fkcopt,fval,exitflag,output] = fmincon(@(Fkc) ObjFun_Zdd_F(Fkc),Fkco,...

and this is fmincon

Fkco=[100,1000,500];
lb=[0,100,50];
ub=[3000,160000,50000];
options = optimset('Display','iter','TolFun',1e-8,'algorithm','sqp');
[Fkcopt,fval,exitflag,output] = fmincon(@(Fkc) ObjFun_Zdd_F(Fkc),Fkco,...
                                 [],[],[],[],lb,ub,[],options);

and the objective function is:

    function [val] = ObjFun_Zdd_F(Fkc)

    F=Fkc(1);
    kc(1)=Fkc(2);
    kc(2)=Fkc(3);
    load IRI_737b
    zdot = calc_zdot(road_x, road_z, 10); 
    zt=zdot.u;
    x00=[0,0,0,0];
    input=@(t)lookup_u(zdot,t);
    opt = odeset('RelTol', 1e-2, 'AbsTol', 1e-3);
    [t02,y] = ode23(@(t,x)f(t,x,F,kc, @(t)lookup_u(zdot,t)), [0 2], x00,opt);

    xx=size(t02,1);
    x_frame=y'; 
    FF=F * ones(1, xx);
    xu_frame=[x_frame;FF];

    val = 0;
    r1 = 1e+5
    r2 = 0.5;
    q = 1e-5;
    R = diag([r1, 0, 0, 0, q,r2]);
    xu_frame_zs = zeros(6,length(t02));
        for i=1:205
            ddx=ff(i,xu_frame(1:4,i),F,kc);
            Acc=ddx(4);
            xu_frame_zs(:,i) = [xu_frame(:,i);Acc];
            ttt=t02(i+1)-t02(i);      % here the error pointed
            val=val+ttt*xu_frame_zs(:,i)'*R*xu_frame_zs(:,i);
        end
   end

[Walter Roberson]

In other places in the code you are careful to check size(t02,1) or length(t02), but in that one place you assume that t02 is at least 206 long. Is there a reason you are not using xx-1 as the upper limit instead of hard-coding 205 ?

[Muna Shehan]

Thanks Walter for your replay.Actually the t02 length is 206 and because I get this error I try to enforce the exactly value of ( length(t02)-1).

[Walter Roberson]

Your code has

          for i=1:205

and you access up to i+1 which requires index 206, but your data is coming out length 205 in one of the ode calls.

Remember, when you only use a tspan range like [0 2] then the ode will decide by itself how many points to use, based upon what it decides it needs to meet the integration tolerances. When you use a slightly different input to it, it could decide that it needed more or fewer integration points.

[Muna Shehan]

Thanks Walter for your replay. The used points number that is decided by ode23 is 206 which is the length of t02, but I need to use these points inside the loop where the error occur. Inside the loop I need to determine ( ttt=t02(i+1)-t02), so I suppose the upper limit will be ( length(t02)-1).

[Walter Roberson]

Is there a reason you have coded it as

for i=1:205

instead of

for i=1:length(t02)-1

??

I suspect, by the way, that on the run that is failing, ode23 is deciding to return 205 timepoints instead of 206.

[Muna Shehan]

there is no reason I used 205 which I suppose is ( length(t02)-1). but what you mentioned is interesting because I suppose ode23 return 206 timepoints because the ode23 output is t02 with dimension( 206x1) and y with dimension ( 206x4). I can not understand this message error

??? Attempted to access t02(206); index out of bounds because numel(t02)=205.

because t02 dimension is ( 206x1)

[Walter Roberson]

At the command line, give the command

dbstop if error

and run the program. When it stops with that error message, show size(t02) and size(y) and look at xx, and between all of those figure out whether the ode routine returned 205 or 206 values.

[Muna Shehan]

I followed your instruction, and what you mention at the beginning is exactly what occur in my case. It is true size(t02) is 206 for the first fmincon iteration and since the ode23 is inside the objective function ObjFun_Zdd_F(Fkc) so size(t02) is changed based on what you explain before "the ode will decide by itself how many points to use, based upon what it decides it needs to meet the integration tolerances" which what happened in the next fmincon iteration. I have to put ode23 inside the objective function ObjFun_Zdd_F(Fkc) and even if I use length(t02)-1 as the loop upper limit I need to determine the output of lookup_u(zdot,t) which I did it before manually by using a breakpoint and debug the change in both t and u inside lookup_u(zdot,t) and re-write t and u in a vectors to filter it later based on the t02 values. After that I write filtered u values in a vector and use it later in the appropriate location inside ff(i,xu_frame(1:4,i),F,kc). Now I use length(t02)-1 as the upper limit and inside the loop I call lookup_u(zdot, t02(i)) and save its value in Zot02 then send this value to ff(i,xu_frame(1:4,i),F,kc,Zot02) as a passing parameter

Zot02= lookup_u(zdot, t02(i));
ddx=ff(i,xu_frame(1:4,i),F,kc,Zot02);

and lookup_u(tu, t) is

function u = lookup_u(tu, t)
i = find(tu.t <= t, 1, 'last');
[m,n] = size(tu.u);
if (m==1) || (n==1)
    if i<length(tu.t)
        u = (tu.u(i)*(tu.t(i+1) -t) + tu.u(i+1)*(t-tu.t(i)))/(tu.t(i+1)-tu.t(i));
    else
        u = tu.u(i); 
    end
else
    if i < length(tu.t)
        u = (tu.u(:,i)*(tu.t(i+1) -t) + ...
            tu.u(:,i+1)*(t-tu.t(i)))/(tu.t(i+1)-tu.t(i));
    else
        u = tu.u(:,i);
    end
end
end

This modification is allowed to calculate u corresponding to the t02 index, the code is run but it seem I still need to enhance fmincon results. Thanks Walter for all your replays because your comments are the reason which highlights my way to understand what happen and lead me to this conclusion.

                                                            Norm of First-order
   Iter F-count            f(x) Feasibility  Steplength        step  optimality
      0       4   4.771200e+002  0.000e+000                          4.350e+002
      1       8   3.767602e+002  0.000e+000  1.000e+000  2.081e+002  5.935e+003
      2      12   2.864263e+002  0.000e+000  1.000e+000  2.656e+002  3.344e+002
      3      16   2.097120e+002  0.000e+000  1.000e+000  4.654e+002  4.192e+003
      4      30   1.709489e+002  0.000e+000  2.825e-002  7.196e+002  1.239e+004
      5      58   1.696533e+002  0.000e+000  1.916e-004  4.805e+000  1.103e+002
      6      73   1.696308e+002  0.000e+000  1.977e-002  1.265e-001  1.106e+002
      7      77   1.696229e+002  0.000e+000  1.000e+000  3.172e-002  2.103e-001
      8      81   1.696210e+002  0.000e+000  1.000e+000  1.031e-002  1.301e-001
      9      86   1.696206e+002  0.000e+000  7.000e-001  1.277e-002  3.906e-002
     10      90   1.696206e+002  0.000e+000  1.000e+000  4.313e-003  3.354e-002
     11      94   1.696205e+002  0.000e+000  1.000e+000  1.991e-002  2.437e-002
     12      98   1.696203e+002  0.000e+000  1.000e+000  6.413e-002  3.515e-002
     13     104   1.696198e+002  0.000e+000  4.900e-001  1.593e-001  3.342e-001
     14     110   1.696197e+002  0.000e+000  4.900e-001  5.265e-002  6.015e-002
     15     114   1.696196e+002  0.000e+000  1.000e+000  2.605e-002  5.945e-002
     16     120   1.696193e+002  0.000e+000  4.900e-001  6.390e-002  3.182e-001
     17     129   1.696193e+002  0.000e+000  1.681e-001  4.015e-003  5.966e-002
     18     133   1.696192e+002  0.000e+000  1.000e+000  2.304e-002  9.673e-002
     19     137   1.696190e+002  0.000e+000  1.000e+000  1.701e-001  3.285e-001
     20     141   1.696186e+002  0.000e+000  1.000e+000  1.593e-001  3.340e-001
     21     145   1.696171e+002  0.000e+000  1.000e+000  8.108e-001  4.778e-001
     22     149   1.696056e+002  0.000e+000  1.000e+000  4.002e+000  4.248e-001
     23     153   1.687674e+002  0.000e+000  1.000e+000  4.070e+001  1.103e+001
     24     171   1.685538e+002  0.000e+000  6.782e-003  2.743e-001  2.253e+000
     25     176   1.685358e+002  0.000e+000  7.000e-001  5.257e-002  4.522e+000
     26     180   1.685354e+002  0.000e+000  1.000e+000  3.341e-002  1.945e+001
     27     184   1.685342e+002  0.000e+000  1.000e+000  2.767e-002  1.488e+000
     28     188   1.685340e+002  0.000e+000  1.000e+000  3.975e-003  1.894e-001
     29     192   1.685340e+002  0.000e+000  1.000e+000  2.070e-004  1.683e-001

Local minimum possible. Constraints satisfied.

fmincon stopped because the size of the current step is less than
the default value of the step size tolerance and constraints are 
satisfied to within the default value of the constraint tolerance.

<stopping criteria details>

>>

Answer 1

[Walter Roberson]

Poster had hard-coded the expected number of points returned by ode23s, but ode23s decides dynamically how many points to use and was not always using the same number of points. User has changed the code to not hard-code the number of points, solving the original problem.

[Muna Shehan]

what you say "'ode23s decides dynamically how many points to use and was not always using the same number of points" is true.The ode23s change the number of points at every optimization iteration which I did not know this information before this time. Thanks again