Index exceeds matrix dimensions

Hi everyone,

I am trying to model a problem for FTCS and I am using this code:

% --- Define constants and initial condition

L =1; % length of domain in x direction

tmax= 0.8; % end time

nx =10; % number of nodes in x direction

nt =10; % number of time steps

dx = L/(nx-1);

dt = tmax/(nt-1);

r = 0.5; r2 = 1 - 2*r;

% --- Loop over time steps

t = 0;

u =100; % initial condition

for m=1:nt

uold = u; % prepare for next step

t = t + dt;

for i=2:nx-1

u(i) = r*uold(i-1) + r2*uold(i) + r*uold(i+1);

end

end

But I couldn't figure out whats the mistake here, I am new to this so anyhelp is much appreciated.

Answer 1

[Reshma Nerella]

Hi,

In the code, you are assigning a value to variable 'u'

u =100 % 1x1 Double

Implies you can only index till 1(since we have only one element).

for i=2:nx-1
u(i) = r*uold(i-1) + r2*uold(i) + r*uold(i+1); %You are indexing 'u' from 2, beyond dimensions.
end

If you want to find 'u' matrix, consider preallocating it for speed and assign required values to it.

u = zeros(1,nx-1); %creates a matrix of size 1 x (nx-1) with zeros

For the initial condition, if you want first element of 'u' to be 100, You can assign it this way.

u(1) = 100;

You may get the issue with 'uold'

uold = u; % assigning a single value

for i=2:nx-1
u(i) = r*uold(i-1) + r2*uold(i) + r*uold(i+1); %Indexing it beyond dimensions
end

Hope this helps!

Answer 2

[Harsha]

Thank you very much for your reply,

I actually modified the code to the follwoing:

c>> n = 10;

Lx = 1;

dx = Lx/(n-1);

x = 0:dx:Lx;

% 2. Parameters for the t vector

m = 80;

tf = 0.8;

dt = tf/(m-1);

t = 0:dt:tf;

% 3. Other parameters needed for the solution

% The value of alpha

Fo = .5; % a mutiplicative constant that should be < = 1/2

% to insure stability

% Initial and boundary conditions

f = @(x)100;

g1 = @(t)0;

g2 = @(t)0;

u = zeros(n,m);

u(2:n,1) = f(x(2:n)); % Put in the initial condition starting from

% 2 to n-1 since f(0) = 0 and f(N) = 1

u(1,:) = g1(t); % The boundary conditions, g1 and g2 at

u(n,:) = g2(t); % x = 0 and x = 1

% Implementation of the explicit method

for k = 2:m-1 % Time Loop

for i= 2:n-1 % Space Loop

u(i,k+1) = Fo*(u(i-1,k)+u(i+1,k))+(1-2*Fo)*u(i,k);

end

end

figure

hold all

for i=1:4:numel(t)

plot(x,u(:,i),'linewidth',1.5,'DisplayName',sprintf('t = %1.2f',t(i))); fprintf(x,u(:,i),'linewidth',1.5,'DisplayName',sprintf('t = %1.2f',t(i)));

end

a = ylabel('Heat');

set(a,'Fontsize',14);

a = xlabel('x');

set(a,'Fontsize',14);

a=title(['Using The Explicit Method - Fo =' num2str(Fo)]);

set(a,'Fontsize',16);

legend('-DynamicLegend','location','bestoutside');

grid;

figure

[X, T] = meshgrid(x,t);

s2 = mesh(X',T',u);

title(['3-D plot of the 1D Heat Equation using the Explicit Method - Fo =' num2str(Fo)])

set(s2,'FaceColor',[0 0 1],'edgecolor',[0 0 0],'LineStyle','--');

a = title('Exact solution of the 1D Diffusivity Equation');

set(a,'fontsize',14);

a = xlabel('x');

set(a,'fontsize',20);

a = ylabel('y');

which I got it from here: https://www.mathworks.com/matlabcentral/answers/404194-question-on-heat-equation-1d-forward-in-time-centered-in-space

Now my only problem is to get the values as a table like the attached pic, any idea for doing that?

Edit I also keep getting zeros for all t .

Very grateful for your help!

Answer 3

[VBBV]

n = 10;
Lx = 1;
dx = Lx/(n-1);
x = 0:dx:Lx;
% 2. Parameters for the t vector
m = 80; 
tf = 0.8;
dt = tf/(m-1);
t = 0:dt:tf;
% 3. Other parameters needed for the solution
% The value of alpha
Fo = .5;   % a mutiplicative constant that should be < = 1/2
            % to insure stability
% Initial and boundary conditions
f  = 100;   
g1 = 0;           
g2 = 0;
u  =  zeros(n,m);
% u(2:n-1,1) =  f;    % Put in the initial condition starting from
                          % 2 to n-1  since f(0) = 0 and f(N) = 1 
u(1,:)   =  g1;       % The boundary conditions, g1 and g2 at
u(n,:)   =  g2;       % x = 0 and x = 1
u(2:n-1,:) =  f;  
% Implementation of the explicit method
for k = 2:m-1   % Time Loop
      for i= 2:n-1  % Space Loop
%         u(i,k+1) = Fo*(u(i-1,k)+u(i+1,k))+(1-2*Fo)*u(i,k);
        u(i,k+1) = (Fo*(u(i-1,k)+u(i+1,k))+2*u(i,k))/4;
      end  
end
figure
hold all 
for i=1:2:numel(t)
    plot(x,u(:,i),'linewidth',1.5,'DisplayName',sprintf('t = %1.2f',t(i))); 
    hold on
    %fprintf(x,u(:,i),'linewidth',1.5,'DisplayName',sprintf('t = %1.2f',t(i)));
end
a = ylabel('Heat');
set(a,'Fontsize',10);
a = xlabel('x');
set(a,'Fontsize',10);
a=title(['Using The Explicit Method - Fo =' num2str(Fo)]);
set(a,'Fontsize',10);
legend('-DynamicLegend','location','bestoutside');
grid;

figure
[X, T] = meshgrid(x,t); 
s2 = contourf(x,t,u.',0:5:100);
colorbar
shading interp
axis([0 1 0 0.12])
title(['3-D plot of the 1D Heat Equation using the Explicit Method - Fo =' num2str(Fo)])
xlabel('Bar length')
ylabel('Time step')

You need to change the boundary conditions to get the right solution.

Also your finite difference scheme application seems to be incorrect