Need help getting this code to plot

6 次查看(过去 30 天)
I want to plot the capacitance in this graph (cap) vs. the grid size to show how (in theory) it should get moe and more exact as grid size gets smaller. But I can't figure out how to do it. I tried putting it in the for loop, and outside the for loop, but have made no progress. Can someone help me figure out how to do this plot?
% --------------------------------------------------------------% Compute capacitance per unit length of
% a coaxial pair of rectangles
% -------------------------------------------------------------
function cap = capacitor(a, b, c, d, n, tol, rel)
% Arguments:
a = 0.015; %width of inner conductor
b = 0.015; %height of inner conductor
c = 0.02; %width of outer conductor
d = 0.02; %height of outer conductor
n = 10; %number of points in the x-direction (horizontal)
tol = 0.001; %relative tolerance for capacitance
rel = 1.7; %relaxation parameter
% (a good choice is 2-c/n, where c is about pi)
% Returns:
% cap = capacitance per unit length [pF/m]
% Make grids
h = 0.5*c/n; % Grid size
na = round(0.5*a/h); % Number of segments on ’a’
x = linspace(0,0.5*c,n+1); % Grid points along x-axis
m = round(0.5*d/h); % Number of segments on ’d’
mb = round(0.5*b/h); % Number of segments on ’b’
y = linspace(0,0.5*d,m+1); % Grid points along y-axis
% Initialize potential and mask array
f = zeros(n+1,m+1); % 2D-array with solution
mask = ones(n+1,m+1)*rel; % 2D-array with relaxation
% [mask(i,j) = 0 implies
% unchanged f(i,j)]
for i = 1:na+1;
for j = 1:mb+1;
mask(i,j) = 0;
f(i,j) = 1;
end
end
% Gauss Seidel iteration
oldcap = 0;
for iter = 1:1000000; % Maximum number of iterations
f = seidel(f,mask,n,m); % Perform Gauss-Seidel iteration
cap = gauss(n,m,h,f);% Compute the capacitance
if (abs(cap-oldcap)/cap<tol);
break % Stop if change in capacitance
% is sufficiently small
else
oldcap = cap; % Contiue until converged
end
end
str = sprintf('Number of iterations = %4i',iter); disp(str);
  1 个评论
CAM
CAM 2015-12-2
other relevant code:
function cap = gauss(n, m, h, f)
% Arguments:
% n = number of points in the x-direction (horizontal)
% m = number of points in the y-direction (vertical)
% h = cell size
% f = 2D-array with solution
% Returns:
% cap = capacitance per unit length [pF/m]
q=0;
for i = 1:n;
q=q+(f(i,m)+f(i+1,m))*0.5; % integrate along upper boundary
end
for j = 1:m;
q=q+(f(n,j)+f(n,j+1))*0.5; % integrate along right boundary
end
cap = q*4; % 4 quadrants
cap = cap*8.854187; % epsilon0*1e12 gives answer in pF/m
and
% --------------------------------------------------------------% Make one Seidel iteration
% -------------------------------------------------------------
function f = seidel(f, mask, n, m)
% Arguments:
% f = 2D-array with solution
% mask = 2D-array with relaxation
% n = number of points in the x-direction (horizontal)
% m = number of points in the y-direction (vertical)
% Returns:
% f = 2D-array with solution after Gauss-Seidel iteration
% Gauss seidel iteration
for i = 2:n
for j = 2:m
f(i,j) = f(i,j) + mask(i,j)* ...
(0.25*( f(i-1,j) + f(i+1,j) ...
+ f(i,j-1) + f(i,j+1)) - f(i,j));
end
end
% Symmetry on left boundary i-1 -> i+1
i=1;
for j = 2:m
f(i,j) = f(i,j) + mask(i,j)* ...
(0.25*( f(i+1,j) + f(i+1,j) ...
+ f(i,j-1) + f(i,j+1)) - f(i,j));
end
% Symmetry on lower boundary j-1 -> j+1
j=1;
for i = 2:n;
f(i,j) = f(i,j) + mask(i,j)* ...
(0.25*( f(i-1,j) + f(i+1,j) ...
+ f(i,j+1) + f(i,j+1)) - f(i,j));
end

请先登录,再进行评论。

采纳的回答

Ingrid
Ingrid 2015-12-2
it is not clear what you want to do but is it something like this?
function plotting(gridSizes)
cap = zeros(length(gridSizes),1);
for ii = 1:length(gridSizes)
cap(ii) = capacitor(gridSizes(ii));
end
plot(gridSizes,cap);
end
you could call this function with
gridSizes = 1:10;
plotting(gridSizes);

更多回答(0 个)

类别

Help CenterFile Exchange 中查找有关 Line Plots 的更多信息

Community Treasure Hunt

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

Start Hunting!

Translated by