I'm solving using modefied euler method, and then finding Jacobian matrix (Fi_x). after running the code it gives me jacobian matrix 6x4. but it should be 6x6 because I have here n=6 and jacobian matrix should be in 6x6 dimension.

   kmax=(tc-t0)/h;
         Xd1=[x0];
         Tm1=[t0];
 for k=1:kmax
         t1=k*h;Tm1=[Tm1 t1];
         xd=modeu(Fdfpek,t1,x0,h);
         xp1=x0+eps*v1;
         Fi2=modeu(Fdfpek,t1,xp1,h);
    Fi1=xd;
    Fix1=(1./eps)*(Fi2-Fi1);
    %-----------------------
    xp2=x0+eps*v2;
    Fi2=modeu(Fdfpek,t1,xp2,h);
    Fi1=xd;
    % second column of Fi_x
    Fix2=(1./eps)*(Fi2-Fi1);
    Fi_x1=[Fix1 Fix2];
    Fi_xM1(:,:,k)=Fi_x1;
    %-----------------------
         Xd1=[Xd1;  xd'];
         x0=xd';

 end
         x0c =Xd1(kmax,:);
         Y=abs(Yaf); th=angle(Yaf);
         nmax=(tf-tc)/h;
         Xc1=[x0c];
        Tm2=[tc];
 for n=1:nmax
         t2=n*h+tc;Tm2=[Tm2 t2];
         xa=modeu(Fafpek,t2,x0c,h);
           xp11=x0c+eps*v3;
         Fi2=modeu(Fafpek,t2,xp11,h);
    Fi1=xa;
    Fix1=(1./eps)*(Fi2-Fi1);
    %-----------------------
    xp22=x0c+eps*v4;
    Fi2=modeu(Fafpek,t2,xp22,h);
    Fi1=xa;
    % second column of Fi_x
    Fix2=(1./eps)*(Fi2-Fi1);
    Fi_x2=[Fix1 Fix2];
    Fi_xM2(:,:,n)=Fi_x2;
    %-----------------------
      x1(n+1)=xa(1);
      x2(n+1)=xa(2);
       Xc1=[Xc1; xa'];
        x0c=xa';

 end
 Fi_x=[Fi_xM1(:,:,k), Fi_xM2(:,:,n)];

.......................................................................
        t =[Tm1'; Tm2']; x = [Xd1; Xc1];

[ammar samii]

The solution here divided into tow steps: from t0-tc and from tc-tf. then the total solution added into one (which is x as shown in the final row)