could you help me in this code
1 回表示 (過去 30 日間)
古いコメントを表示
% Define grid parameters
h = 1/4; % Grid spacing
x = -1:h:1; % Grid points in x direction
y = -1:h:1; % Grid points in y direction
n = length(x); % Number of grid points
% Initialize solution matrix
U = zeros(n,n);
% Set boundary conditions
U(1,:) = 0; % u(x,y)=0 along y=1
U(n,:) = 0; % u(x,y)=0 along y=-1
U(:,1) = 0.5*U(:,2); % u(x,y)=(1/2)u, x=-1,-1<y<1
U(:,n) = -0.5*U(:,n-1); % u(x,y)=-(1/2)u, x=1,-1<y<1
% Construct matrix A
B = [-6 2 0 0 0; 1 -6 1 0 0; 0 1 -6 1 0; 0 0 1 -6 1; 0 0 0 2 -6.25];
I = eye(5);
A = kron(eye(n-2), B) + kron(diag(ones(n-3,1), -1), I) + kron(diag(ones(n-3,1), 1), I);
% Construct right-hand side vector rhs from boundary conditions
rhs = zeros(n^2, 1);
rhs(1:n) = U(1,:)' / h^2;
rhs(end-n+1:end) = -U(n,:)' / h^2;
rhs(1:n:n*(n-2)+1) = rhs(1:n:n*(n-2)+1) + U(:,1) / h^2;
rhs(n:n:n^2-n+1) = rhs(n:n:n^2-n+1) - U(:,n) / h^2;
Arrays have incompatible sizes for this operation.
% Solve the system of equations using matrix inversion
U_vec = A \ rhs;
U = reshape(U_vec, [n, n]);
% Plot the solution as a 3-D mesh
figure;
mesh(x, y, U');
xlabel('x');
ylabel('y');
zlabel('u(x,y)');
title('Solution of the Poisson equation using central differences');
Error using \
0 件のコメント
回答 (1 件)
Walter Roberson
2023 年 5 月 4 日
You cannot add an 8 x 1 vector and a 9 x 1 vector.
% Define grid parameters
h = 1/4; % Grid spacing
x = -1:h:1; % Grid points in x direction
y = -1:h:1; % Grid points in y direction
n = length(x); % Number of grid points
% Initialize solution matrix
U = zeros(n,n);
% Set boundary conditions
U(1,:) = 0; % u(x,y)=0 along y=1
U(n,:) = 0; % u(x,y)=0 along y=-1
U(:,1) = 0.5*U(:,2); % u(x,y)=(1/2)u, x=-1,-1<y<1
U(:,n) = -0.5*U(:,n-1); % u(x,y)=-(1/2)u, x=1,-1<y<1
% Construct matrix A
B = [-6 2 0 0 0; 1 -6 1 0 0; 0 1 -6 1 0; 0 0 1 -6 1; 0 0 0 2 -6.25];
I = eye(5);
A = kron(eye(n-2), B) + kron(diag(ones(n-3,1), -1), I) + kron(diag(ones(n-3,1), 1), I);
% Construct right-hand side vector rhs from boundary conditions
rhs = zeros(n^2, 1);
rhs(1:n) = U(1,:)' / h^2;
rhs(end-n+1:end) = -U(n,:)' / h^2;
part1 = rhs(1:n:n*(n-2)+1);
part2 = U(:,1) / h^2;
size(part1)
size(part2)
rhs(1:n:n*(n-2)+1) = part1 + part2
rhs(n:n:n^2-n+1) = rhs(n:n:n^2-n+1) - U(:,n) / h^2;
2 件のコメント
Walter Roberson
2023 年 5 月 4 日
I do not know.
I would suggest to you that it might be easier to costruct rhs as an n x n 2D array and then reshape it to vector afterwards.
参考
カテゴリ
Help Center および File Exchange で Resizing and Reshaping Matrices についてさらに検索
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!