help correct the error
5 views (last 30 days)
Show older comments
kindly help,
the error message reads 'index exceeds matrix dimensions
% Parameters
Lx = 1; % Length of domain in x-direction
Ly = 1; % Length of domain in y-direction
Lz = 1; % Length of domain in z-direction
Nx = 50; % Number of grid points in x-direction
Ny = 50; % Number of grid points in y-direction
Nz = 50; % Number of grid points in z-direction
T = 1; % Final time
alpha = 0.1; % Thermal diffusivity
dt = 0.001; % Time step size
% Grid spacings
dx = Lx/Nx;
dy = Ly/Ny;
dz = Lz/Nz;
% Create grid
x = linspace(0, Lx, Nx+1);
y = linspace(0, Ly, Ny+1);
z = linspace(0, Lz, Nz+1);
[X, Y, Z] = meshgrid(x, y, z);
% Initial condition
u0 = sin(pi*X).*sin(pi*Y).*sin(pi*Z);
% Initialize solution
u = u0;
% Time-stepping loop
t = 0;
while t < T
% Compute Laplacian using finite differences
d2u_dx2 = (-u(3:Nx+2, 2:Ny, 2:Nz) + 16*u(2:Nx+1, 2:Ny, 2:Nz) - 30*u(1:Nx, 2:Ny, 2:Nz) + 16*u(4:Nx-1, 2:Ny, 2:Nz) - u(5:Nx-2, 2:Ny, 2:Nz)) /12*dx^2; %HERE
d2u_dy2 = (-u(3:Nx, 2:Ny+2, 2:Nz) + 16*u(2:Nx, 2:Ny+1, 2:Nz) - 30*u(1:Nx, 2:Ny, 2:Nz) + 16*u(4:Nx, 2:Ny-1, 2:Nz) - u(5:Nx, 2:Ny-2, 2:Nz)) / 12*dy^2;
d2u_dz2 = (-u(3:Nx, 2:Ny, 2:Nz+2) + 16*u(2:Nx, 2:Ny, 2:Nz+1) - 30*u(1:Nx, 2:Ny, 2:Nz) + 16*u(4:Nx, 2:Ny, 2:Nz-1) - u(5:Nx, 2:Ny, 2:Nz-2)) / 12*dz^2;
% Update solution using the method of lines
u(2:Nx, 2:Ny, 2:Nz) = u(2:Nx, 2:Ny, 2:Nz) + alpha*dt*(d2u_dx2 + d2u_dy2 + d2u_dz2);
% Update time
t = t + dt;
end
% Plot final solution
figure;
slice(X, Y, Z, u, Lx/2, Ly/2, Lz/2);
xlabel('x');
ylabel('y');
zlabel('z');
colorbar;
0 Comments
Answers (1)
Walter Roberson
on 17 Oct 2023
Nx = 50; % Number of grid points in x-direction
Ny = 50; % Number of grid points in y-direction
Nz = 50; % Number of grid points in z-direction
50's.
x = linspace(0, Lx, Nx+1);
y = linspace(0, Ly, Ny+1);
z = linspace(0, Lz, Nz+1);
[X, Y, Z] = meshgrid(x, y, z);
50+1 by 50+1 by 50+1 --> 51 x 51 x 51
d2u_dx2 = (-u(3:Nx+2, 2:Ny, 2:Nz) + 16*u(2:Nx+1, 2:Ny, 2:Nz) - 30*u(1:Nx, 2:Ny, 2:Nz) + 16*u(4:Nx-1, 2:Ny, 2:Nz) - u(5:Nx-2, 2:Ny, 2:Nz)) /12*dx^2; %HERE
Nx+2 --> 50+2 --> 52. But the array is only 51.
0 Comments
See Also
Categories
Find more on Numeric Types in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!