2D cylindrical boundriondition problem

2 次查看(过去 30 天)
function cyl3
close all;
L = 1; % total height of the coaxial cable
hr = 0.01;
hz = 0.01; % these should be the r coordinates of grid
vr = 0.8:hr:0.2;
vz = 0:hr:10; % these should be the z coordinates of grid
nr = length(vr);
nz = length(vz);
n = nr*nz;
A = zeros(n);
b = zeros(n,1);
for jz = 1:nz
for ir = 1:nr
i = (jz-1) * nr + ir; % the global index of the ir,jz vertex
r = vr(ir); % the geometrical r coordinate
z = vz(jz); % the geometrical z coordinate
if (r == 0.8) || (r == 0.2) || (z == 0) ||(z == 10)
% implement boundary condition for r_a ?
A(i,i) = 1;
b(i) = 0;
else
% implement the main equation here
A(i,i) = -2/hr^2 - 2/hz^2;
A(i,i+1) =1/hr^2 + (1/(2*r)*hr) ;
A(i,i-1) = 1/hr^2 - (1/(2*r)*hr);
A(i,i+nr) = 1/hz^2;
A(i,i-nr) = 1/hz^2;
end
end
end
% solve the problem
u = A\b;
s = reshape(u,nr,nz);
% ss = zeros(nz,nr);
% for jz = 1:nz
% for ir = 1:nr
% i = (jz-1) * nr + ir;
% ss(jz,i) = u(i);
% end
%end
surf(s)
end

回答(1 个)

Alan Stevens
Alan Stevens 2021-1-8
It looks like
vr = 0.8:hr:0.2;
is the cause of the problem. This leads to vr being empty!

类别

Help CenterFile Exchange 中查找有关 Geodesy and Mapping 的更多信息

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by