Area difference in using trapz() through a loop and without a loop.

2 views (last 30 days)
clear
C1 = 0;
k1 = 1;
k2 = 0;
n1 = 0;
n2 = 0;
nm = 1;
D = 0;
O = 0;
gam = 1e-04.*k1;
alph = 0;
w = linspace(-5,5,1e+06);
delw = w(2)-w(1);
TA = (w+D)./(k1./2);
TB = (w-O)./(gam/2);
xa = 1-1i*TA;
xb = 1-1i*TB;
C = (2.*1i.*sqrt(C1))./(sqrt(gam).*((xa.*xb)+C1));
D = (2.*xa)./(sqrt(gam).*((xa.*xb)+C1));
MC = (abs(C)).^2;
MD = (abs(D)).^2;
Sbb = 2.*pi*(MC.*(n1+0.5)+MD.*(nm+0.5));
A = delw.*trapz(Sbb)./(2.*pi);
figure(1)
plot(w,Sbb);
set(gca,'FontSize',13)
Upon compiling the code above, one can observe that my area under the curve of Sbb versus w as given by A = delw.*trapz(Sbb)./(2.*pi) gives A = 9.4247. This setting is at C1 = 0. If I change my C1 to C1 = 1, my new A value gives A = 6.2834 (try it!). All is good.
Now, I know that A is dependent on C1, so I tried to loop C1 over a range of values like the following (it looks long but the code is almost exactly the same but with minor tweaks) to observe the trend:
clear
%C1 = 0:1:1;
C1 = linspace(0,500,1000);
k1 = 1;
k2 = 0;
n1 = 0;
n2 = 0;
nm = 1;
D = 0;
O = 0;
gam = 1e-04.*k1;
alph = 0;
w = linspace(-5,5,1e+06);
delw = w(2)-w(1);
for ii=1:length(C1)
TA = (w+D)./(k1./2);
TB = (w-O)./(gam/2);
xa = 1-1i*TA;
xb = 1-1i*TB;
C = (2.*1i.*sqrt(C1(ii)))./(sqrt(gam).*((xa.*xb)+C1(ii)));
D = (2.*xa)./(sqrt(gam).*((xa.*xb)+C1(ii)));
MC = (abs(C)).^2;
MD = (abs(D)).^2;
Sbb = 2.*pi.*(MC.*(n1+0.5)+MD.*(nm+0.5));
A = delw.*trapz(Sbb)./(2.*pi);
AA(ii)=A;
end
figure(2)
plot(C1,AA);
set(gca,'FontSize',13)
What I get is a bunch of zagged lines which does not seem to fit the bill. My first motivation was to switch the range of my C1 so that they are in increments of one unit and just up to 1 and recompile it:
C1 = 0:1:1;
%C1 = linspace(0,500,1000);
If I look over at my workspace, at C1 = 1, my AA(2) = 9.6485 which is clearly not 6.2834 as computed without the loop (but still at C1=1). I am absolutely convinced that 6.2834 is the correct value and 9.6485 isn't (Since A should decrease with increasing C1). What is the problem here? I would greatly appreciate any help with this.
Thanks in advance!
  2 Comments
John D'Errico
John D'Errico on 4 Nov 2017
Edited: John D'Errico on 4 Nov 2017
You accepted the answer the last time you asked virtually the identical question. What was wrong with that answer?
In fact, that answer was correct to your last question. That you do not understand the answer is no reason to ask the same question multiple times.
Alvin
Alvin on 4 Nov 2017
I did not say that anything is wrong with the answer to my last question. In fact the answer made sense. This is a different question looping over a particular parameter albeit in the same code structure as before.
I suggest you re-read my last post. If you've read it with proper scrutiny, you'd noticed the difference.

Sign in to comment.

Answers (0)

Tags

Community Treasure Hunt

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

Start Hunting!