Why do I keep getting "Index exceeds matrix dimensions" while solving partial differential equation?

So I have a computational fluid dynamics problem to solve the generalized Burger equation. Using two different PDE schemes. But I keep getting this error "Index exceeds matrix dimensions." The thing that is puzzling me is I didn't specify any dimensions. So I don't know how to fix it.

for x = 1:X
    ub(1,x+1) = 0.5*(1+atan(250*(x-20)));
end
% Solving the Burger equation
for n = 1:T
    u(n+1,1) = 0;
    for j = 2:X-1

          Fnj1 = 0.5*u(n,j+1)*(1-u(n,j+1));
          Fnj = 0.5*u(n,j)*(1-u(n,j));
          ******Fbnj = 0.5*ub(n,j+1)*(1-ub(n+1,j));***************Here is where it says the error is
          Fbnj1 = 0.5*ub(n+1,j-1)*(1-ub(n+1,j-1));

          % Predictor
          ub(n+1,j) = u(n,j)-dt/dx*(Fnj1-Fnj)+r*(u(n,j+1)-...
              2*u(n,j)+u(n,j-1));

          % Corrector
          u(n+1,j) = 0.5*(u(n,j)+ub(n+1,j)-(dt/dx)*(Fbnj-...
              Fbnj1)+r*(ub(n+1,j+1)-2*ub(n+1,j)+ub(n+1,j-1)));

      end
      u(n+1,X) = 1;
  end

Answer 1

[KSSV]

Check line 9 u(n+1,1) = 0 and line 11 0.5*u(n,j+1). In line line 9 you said u is a row vector and in line you are treating u as a matrix. That's why there is a error. You have to initialize u properly.

[Ben McDaniel]

I'm sorry but I don't follow what you mean. I used this same notation on another project and it worked fine.

[KSSV]

You have u(n+1,1) and after you are accessing u(n,2), so obviously error pops out.

[Ben McDaniel]

But the error is with a ub variable not the u.

Answer 2

[Walter Roberson]

Follow the execution from the beginning. The first iteration of "for n", n will be assigned 1. You have

      u(n+1,1) = 0;

so that is u(2,1) = 0; and that will act to create the previously non-existent array u, implicitly defining u(1,1) = 0 (matrix extension fills empty slots with 0), and explicitly defining u(2,1) = 0

Then you start the "for j" loop, in which the j first gets assigned 2. You proceed from there to the next line which has

          Fnj1 = 0.5*u(n,j+1)*(1-u(n,j+1));

which uses u(1,3) because n = 1 and j = 2. But u(1,3) does not exist. Your code will fail before you reach the line you marked,

        Fbnj = 0.5*ub(n,j+1)*(1-ub(n+1,j));

From this we deduce that the code you posted is not your actual code, or possibly that you had a u variable left over from a previous execution, with some value in it that we cannot guess at.

[Walter Roberson]

Your loop

for x = 1:X
    ub(1,x+1) = 0.5*(1+atan(250*(x-20)));
end

starts with ub undefined, so the first iteration assigns to ub(1,2), implicitly assigning 0 to ub(1,1) because of matrix extension automatically filling with 0. The next iteration assigns to ub(1,3) and so on to ub(1,X+1) which is ub(1,42)

You start the for n loop with n = 1. You start the for j loop with j = 2. You reach

            Fb(n,j) = 0.5*ub(n+1,j)*(1-ub(n+1,j));

this calls upon ub(1+1,2) which is ub(2,2) . But you have not assigned anything to ub(2,:) yet. You do not assign anything to ub(2,:) until two lines further on,

            ub(n+1,j) = u(n,j)-dt/dx*(F(n,j+1)-F(n,j))+r*(u(n,j+1)-...
              2*u(n,j)+u(n,j-1));

which would assign to ub(2,2)

I would have thought it likely in a 2D PDE that you would be wanting to use ub(n,:) in your predictions of ub(n+1,:) but you do not do so.

Answer 3

[Ben McDaniel]

Okay. Here is the whole code because a similar version worked on a different problem.

clear variables;
clc;

T=18;       % Time      
X=41;       % # of mesh points
ne=3;
nu=0.6;     
dx=1;       % mesh spacing
c=1;
dt=0.5;   
v=(c*dt)/dx;
r=0.001*(dt/dx^2);

u = zeros(1,X);       % Creates 1x41 matrix for analytical solution
u_e = zeros(1,X);     % Creates 1x41 matrix for exact Sol

% Solving I.C. from x=1 to x=41

for x = 1:X
    u(1,x) = 0.5*(1+atan(250*(x-20))); 
end  

F = 0.5*u-0.5*u.^2;
Fb = zeros(1,X);

% Using predictor with IC

for x = 1:X
    ub(1,x+1) = 0.5*(1+atan(250*(x-20)));
end
% Solving the Burger equation 

for n = 1:T
    u(n+1,1) = 0;
    for j = 2:X-1        

          F(n,j+1) = 0.5*u(n,j+1)*(1-u(n,j+1));
          F(n,j) = 0.5*u(n,j)*(1-u(n,j));
          Fb(n,j) = 0.5*ub(n+1,j)*(1-ub(n+1,j));
          Fb(n,j-1) = 0.5*ub(n+1,j-1)*(1-ub(n+1,j-1));

          % Predictor
          ub(n+1,j) = u(n,j)-dt/dx*(F(n,j+1)-F(n,j))+r*(u(n,j+1)-...
              2*u(n,j)+u(n,j-1));

          % Corrector
          u(n+1,j) = 0.5*(u(n,j)+ub(n+1,j)-(dt/dx)*(Fbnj-...
              Fbnj1)+r*(ub(n+1,j+1)-2*ub(n+1,j)+ub(n+1,j-1)));

      end
      u(n+1,X) = 1;
  end

% Solving for the exact solution 

for n = 1:18
    for x = 1:X
        u_e(n,x) = 0.5*(1+atan(250*(x-20)));
    end
end

x = 1:X;

% Comparison Graph

plot(x,u(n+1,j),x,u_e(n,x)),title('MacCormack Method'),grid on;

[KSSV]

Attempted to access ub(2,2); index out of bounds because size(ub)=[1,42].

Error in Main (line 31) Fb(n,j) = 0.5*ub(n+1,j)*(1-ub(n+1,j));

[Ben McDaniel]

I managed to change some indices and get that equation to work, and I made a few changes but I still keep getting the "predictor" equation coming up with the matrix dimension error. How can I fix this? I really suck with matlab apparently. This problem has been driving me crazy.

[Walter Roberson]

Please post the current version of the code.

Answer 4

[Ben McDaniel]

Here is the current version.

clear variables;
clc;

T=18;       % Time      
X=41;       % # of mesh points    
dx=1;       % mesh spacing
c=1;
dt=1;   
v=(c*dt)/dx;
r=0.001*(dt/(dx^2));

u = zeros(T,X);       % Creates 18x41 matrix for analytical solution
u_e = zeros(T,X);     % Creates 18x41 matrix for exact Sol
ub = zeros(T,X);      % Creates 18x41 matrix for predictor

% Solving I.C. from x=1 to x=41

for x = 1:X
    u(1,x) = 0.5*(1+atan(250*(x-20))); 
end  

F = 0.5*u-0.5*u.^2;

% Using predictor with IC

for x = 1:X
    ub(1,x) = 0.5*(1+atan(250*(x-20)));
end
% Solving the Burger equation 

for n = 1:dt:T
    u(n,1) = 0;
    for j = 2:(X+1)

          F1 = 0.5*u(n,j+1)*(1-u(n,j+1));
          F = 0.5*u(n,j)*(1-u(n,j));
          % Predictor
          ub(n+1,j) = u(n,j)-(dt/dx)*(F1-F)+...
              r*(u(n,j+1)-2*u(n,j)+u(n,j-1));

          % Corrector
          Fb = 0.5*ub(n+1,j)*(1-ub(n+1,j));
          Fb1 = 0.5*ub(n+1,j-1)*(1-ub(n+1,j-1));
          u(n+1,j) = (0.5)*((u(n,j)+ub(n+1,j)-(dt/dx)*(Fb-Fb1)+...
              r*(ub(n+1,j+1)-(2*ub(n+1,j))+ub(n+1,j-1))));

      end
      u(n,X) = 1;
  end

% Solving for the exact solution 

for n = 1:T
    for x = 1:X
        u_e(n,x) = 0.5*(1+atan(250*(x-20)));
    end
end

x = 1:X;

% Comparison Graph

plot(x,u(5,j),x,u_e(5,x)),title('MacCormack Method'),grid on;

[Walter Roberson]

In the line

            F1 = 0.5*u(n,j+1)*(1-u(n,j+1));

your j can be as large as X+1 and then you add 1 to it, so you could be accessing a subscript as large as X+2

You initialized

u = zeros(T,X);       % Creates 18x41 matrix for analytical solution

which is only length X, and your loop

for x = 1:X
    u(1,x) = 0.5*(1+atan(250*(x-20))); 
end

never writes past X

Your corrector does have

          u(n+1,j) = (0.5)*((u(n,j)+ub(n+1,j)-(dt/dx)*(Fb-Fb1)+...
              r*(ub(n+1,j+1)-(2*ub(n+1,j))+ub(n+1,j-1))));

and since j can be as large as X+1 that could potentially extend u to X+1 columns, but it cannot extend it to X+2 columns . And, besides, you would already have failed out on the

            F1 = 0.5*u(n,j+1)*(1-u(n,j+1));

line when j is X, trying to access u(n,X+1) in the loop iteration before you could potentially extend u to width X+1.