I am trying to use the Crank Nicolson method for solving second order equation. But there is something wrong with the code. I am not very good at Matlab.

          A = [ 0 1; -1 -0.1];
          h = 0.1;
          t = 0:h:10^2;
          y0 = [-pi;1];
          n=t/h;
          I = eye(2,2);
          y = zeros(2,length(t));
          y(:,1) = y0;

         b = (I - h/2.*A)^(-1).*(I + h/2.*A);

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

          y(i) = b^(i).*y0;

        end
        plot(t,y)

Answer 1

[Walter Roberson]

Your b is 2 x 2, so b^(i) is 2 x 2. You then use .* to multiply by a 2 x 1 vector. That would fail on versions before R2016b but as of R2016b the 2 x 1 vector would be replicated so that the computation would be b^(i) .* [y0, y0] which would produce a 2 x 2 result. You then try to store the 2 x 42 result in a single location.

I suspect you want b^(i) * y0 instead of b^(i) .* y0 . With b^(i) * y0 that would be a matrix multiplication between a 2 x 2 and a 2 x 1, which is well defined and gives a 2 x 1 result. But then you still try to store the vector into a single location.