Plot using while loop

Please how can I avoid these errors? The graph it doesn't appear (The script is bellow the graph )

close all

clear all

pi=8000;

ri=28;

rp=45;

rp2=35;

ro=80;

r=[ri:0.01:ro];

n=1;

while(r(n)<=rp2)

pa=0;

pb=pi;

Sr(n)= ((pa*ri^2 - pb*rp^2)/(rp^2 - ri^2)) + (((ri^2*rp^2)/r(n)^2)*((pa-pb)/(rp^2-ri^2)));

n=n+1;

end

%n=n-1;

while(r(n)>rp2 & r(n)<=rp)

pa=pi;

pb=0;

Sr(n)= ((pa*rp^2 - pb*ro^2)/(ro^2 - rp^2)) + (((rp^2*ro^2)/r(n)^2)*((pa-pb)/(ro^2-rp^2)));

n=n+1;

end

while (r(n)>rp & r(n)<=ro)

pa=pi;

pb=0;

Sr(n)= ((pa*ri^2 - pb*rp^2)/(rp^2 - ri^2)) + (((ri^2*rp^2)/r(n)^2)*((pa-pb)/(rp^2-ri^2)));

n=n+1;

end

plot(r',Sr','k')

xlabel('r')

ylabel('Sr')

[Mohammad Sami]

The way your code is written, it is possible for the value of n to exceed the length of the variable r.

At n = 5201, your last while loop does not terminate (80 > 45 & 80 <= 80), hence n gets incremented to 5202. However r only has 5201 elements. In the next iteration it gives you the error on index exceed number of array elements.

Also similarly for the plot, you need to make sure both r ad Sr will contain the same number of elements at the end of your code

Answer 1

[KALYAN ACHARJYA]

pi=8000;
ri=28;
rp=45;
rp2=35;
ro=80;
r=[ri:0.01:ro];
n=1;
while(r(n)<=rp2)
pa=0;
pb=pi;
Sr(n)= ((pa*ri^2 - pb*rp^2)/(rp^2 - ri^2)) + (((ri^2*rp^2)/r(n)^2)*((pa-pb)/(rp^2-ri^2)));
n=n+1;
end
%n=n-1;
n=1;
while r(n)>rp2 & r(n)<=rp
pa=pi;
pb=0;
Sr(n)= ((pa*rp^2 - pb*ro^2)/(ro^2 - rp^2)) + (((rp^2*ro^2)/r(n)^2)*((pa-pb)/(ro^2-rp^2)));
n=n+1;
end
while (r(n)>rp & r(n)<=ro)
  pa=pi;
  pb=0;
  Sr(n)= ((pa*ri^2 - pb*rp^2)/(rp^2 - ri^2)) + (((ri^2*rp^2)/r(n)^2)*((pa-pb)/(rp^2-ri^2)));
  n=n+1;
end
plot(r(1:length(Sr))',Sr','k')
xlabel('r')
ylabel('Sr')

Answer 2

[Allen]

Sami is correct that you are running n to a higher value than the total number of elements in r. I better way to approach this problem is to eliminate the loops. Using logical indexing is a cleaner and much more efficient way of doing this. The following code is a direct replacement to yours but with logical indexing.

close all
clear
Pi = 8000; % Change pi to Pi to avoid confusion with pi = 3.1415;
ri = 28;
rp = 45;
rp2 = 35;
ro = 80;
r = ri:0.01:ro;
Sr = NaN(size(r));

pa = 0;
pb = Pi;
idx = r<=rp2; % Logical index for all values of r that meet criteria.
Sr(idx) = (pa*ri^2-pb*rp^2)/(rp^2-ri^2)+ri^2*rp^2./r(idx).^2.*((pa-pb)/(rp^2-ri^2));

pa = Pi;
pb = 0;
idx = r>rp2 & r<=rp; % Logical index for all values of r that meet criteria.
Sr(idx) = (pa*rp^2-pb*ro^2)/(ro^2-rp^2)+ro^2*rp^2./r(idx).^2.*((pa-pb)/(ro^2-rp^2));

pa = Pi;
pb = 0;
idx = r>rp & r<=ro; % Logical index for all values of r that meet criteria.
Sr(idx) = (pa*ri^2-pb*rp^2)/(rp^2-ri^2)+ri^2*rp^2./r(idx).^2.*((pa-pb)/(rp^2-ri^2));

plot(r,Sr,'k')
xlabel('r')
ylabel('Sr')

Since the equation for Sr is similar, but with a change to the variables you could also make your code much more concise by incorporating an anonymous function for forming Sr.

close all
clear
Pi = 8000;
ri = 28;
rp = 45;
rp2 = 35;
ro = 80;
r = ri:0.01:ro;
Sr = NaN(size(r));

% FuncName = @(A,B,R1,R2,R) (A*R1^2-B*R2^2)/(R2^2-R1^2)+R1^2*R2^2./R.^2.*((A-B)/(R2^2-R1^2));
FuncName = @(A,B,R1,R2,R) (A*R1^2-B*R2^2+(A-B)*(R1*R2./R).^2)./(R2^2-R1^2); % More concise arrangement
Sr(r<=rp2) = FuncName(0,Pi,ri,rp,r(r<=rp2));
Sr(r>rp2 & r<=rp) = FuncName(Pi,0,rp,ro,r(r>rp2 & r<=rp));
Sr(r>rp & r<=ro) = FuncName(Pi,0,ri,rp,r(r>rp & r<=ro));

plot(r,Sr)
xlabel('r')
ylabel('Sr')