Arrays have incompatible sizes for this operation. Error in (line 34), can you please tell me how to solve it?
1 次查看(过去 30 天)
显示 更早的评论
% Define problem parameters
L = 1; % size of domain
Ul = 1; % velocity of upper lid
Re = 1000; % Reynolds number
nu = Ul*L/Re; % viscosity
nx = 40; % number of cells in x direction
ny = 40; % number of cells in y direction
dx = L/nx; % cell size in x direction
dy = L/ny; % cell size in y direction
dt = 0.001; % time step
nt = 100; % number of time steps
% Create grid
x = linspace(0, L, nx+1);
y = linspace(0, L, ny+1);
[X, Y] = meshgrid(x, y);
% Set initial conditions
u = zeros(nx+1, ny+2);
v = zeros(nx+2, ny+1);
p = zeros(nx+2, ny+2);
% Define coefficients for finite volume discretization
ae = repmat(1/nu*dt/dx^2, nx+2, ny+2);
aw = repmat(1/nu*dt/dx^2, nx+2, ny+2);
an = repmat(1/nu*dt/dy^2, nx+2, ny+2);
as = repmat(1/nu*dt/dy^2, nx+2, ny+2);
ap = ae+aw+an+as+repmat(1/dt, nx+2, ny+2);
b = zeros(nx+2, ny+2);
% Main loop for time integration
for n = 1:nt
% Calculate intermediate velocities
ue = u(2:end-1, 2:end-1) + dt*(-(p(2:end-1, 3:end)-p(2:end-1, 2:end-1))/dy + nu*(u(3:end, 2:end-1)-2*u(2:end-1, 2:end-1)+u(1:end-2, 2:end-1))/dx^2);
uw = u(1:end-2, 2:end-1) + dt*(-(p(2:end-1, 2:end-1)-p(2:end-1, 1:end-2))/dy + nu*(u(2:end-1, 2:end-1)-2*u(1:end-2, 2:end-1)+u(1:end-3, 2:end-1))/dx^2);
vn = v(2:end-1, 2:end-1) + dt*(-(p(3:end, 2:end-1)-p(2:end-1, 2:end-1))/dx + nu*(v(2:end-1, 3:end)-2*v(2:end-1, 2:end-1)+v(2:end-1, 1:end-2))/dy^2);
vs = v(2:end-1, 1:end-2) + dt*(-(p(2:end-1, 2:end-1)-p(1:end-2, 2:end-1))/dx + nu*(v(2:end-1, 2:end-1)-2*v(2:end-1, 1:end-2)+v(2:end-1, 1:end-3))/dy^2);
end
% Impose boundary conditions
% Velocity boundary conditions
u(1,:) = Ul; % Top wall
v(1,:) = 0;
u(end,:) = 0; % Bottom wall
v(end,:) = 0;
u(:,1) = 0; % Left wall
v(:,1) = 0;
u(:,end) = 0; % Right wall
v(:,end) = 0;
% Pressure boundary conditions
% Top wall (Neumann boundary condition)
for i = 2:Nx-1
Ap =zeros (1, i);
Ap(i,end) = 0; % diagonal element
Aw =zeros (1, i);
Aw(i,end) = 0; % west coefficient
Ae = zeros (1, i);
Ae(i,end) = 0; % east coefficient
An = zeros (1, i);
An(i,end) = 1/dy; % north coefficient
As = zeros (1, i);
As(i,end) = -1/dy; % south coefficient
b(i,end) = 0; % right-hand side
end
% Bottom wall (Neumann boundary condition)
for i = 2:Nx-1
Ap(i,1) = 0; % diagonal element
Aw(i,1) = 0; % west coefficient
Ae(i,1) = 0; % east coefficient
An(i,1) = -1/dy; % north coefficient
As(i,1) = 1/dy; % south coefficient
b(i,1) = 0; % right-hand side
end
% Left wall (Neumann boundary condition)
for j = 2:Ny-1
Ap(1,j) = 1; % diagonal element
Aw(1,j) = 0; % west coefficient
Ae(1,j) = -1/dx; % east coefficient
An(1,j) = 0; % north coefficient
As(1,j) = 0; % south coefficient
b(1,j) = 0; % right-hand side
end
% Right wall (Dirichlet boundary condition)
for j = 2:Ny-1
Ap(end,j) = 1; % diagonal element
Aw(end,j) = -1/dx; % west coefficient
Ae(end,j) = 0; % east coefficient
An(end,j) = 0; % north coefficient
As(end,j) = 0; % south coefficient
b(end,j) = 0; % right-hand side
end
1 个评论
Dyuman Joshi
2023-3-21
编辑:Dyuman Joshi
2023-3-21
As p has different size than u and v, using 2:end-1 or 3:end results in a larger array for p as compared to u and v.
I suggest you to modify the indices accordingly.
L = 1; % size of domain
Ul = 1; % velocity of upper lid
Re = 1000; % Reynolds number
nu = Ul*L/Re; % viscosity
nx = 40; % number of cells in x direction
ny = 40; % number of cells in y direction
dx = L/nx; % cell size in x direction
dy = L/ny; % cell size in y direction
nt = 100; % number of time steps
% Set initial conditions
u = zeros(nx+1, ny+2);
v = zeros(nx+2, ny+1);
p = zeros(nx+2, ny+2);
%for n=1:nt
%Breaking down the formula of intermediate velocity
%ue
A1=-(p(2:end-1, 3:end)-p(2:end-1, 2:end-1))/dy
A2=nu*(u(3:end, 2:end-1)-2*u(2:end-1, 2:end-1)+u(1:end-2, 2:end-1))/dx^2
%end
回答(1 个)
Ashu
2023-3-21
Hi Nusrat
I understand that you are facing error with incompatible array sizes.
Few corrections that you need to consider are -
- Rename nx, ny to Nx, Ny respectively, because in the further code you are refering them with the later variable names.
Nx = 40; % number of cells in x direction
Ny = 40; % number of cells in y direction
- Check your array operations and their sizes in the 'Main operation loop'
Example analysis -
ue = u(2:end-1, 2:end-1) + dt*(-(p(2:end-1, 3:end) - p(2:end-1, 2:end-1))/dy + nu*(u(3:end, 2:end-1)-2*u(2:end-1, 2:end-1)+u(1:end-2, 2:end-1))/dx^2);
% size = 39 x 40 size = 40 x 40 size = 40 x 40
It's clear here that the size for first operand and the second operand don't match.
A simple correction can be like this
ue = u(2:end-1, 2:end-1) + dt*(-(p(2:end-2, 3:end) - p(2:end-2, 2:end-1))/dy + nu*(u(3:end, 2:end-1)-2*u(2:end-1, 2:end-1)+u(1:end-2, 2:end-1))/dx^2);
%------changes in the index----------%
These changes will depend upon your use case, so you can change accordingly.
0 个评论
另请参阅
类别
在 Help Center 和 File Exchange 中查找有关 Computational Fluid Dynamics (CFD) 的更多信息
产品
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!