Inconsistent looping issue with trapz()

 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)

Running the code above would give a straight line from 9.425 to 9.65. If that is the case, then based on my C1 array:

 C1 = 0:1:1

that at C1 = 0, A = 9.4247 and at C1 = 1, A = 9.6485. If I remove my for loop and set my C1 simply to be C1 = 0, I recover A = 9.4247. However, setting C1 = 1 returns A = 6.2834. This is problematic and I know that A = 6.2834 should be the correct answer. Furthermore I would like to loop this over a wide range of C1 values:

 C1 = linspace(0,500,1000)

But upon running it like this I get a lot of erratic spikes which are wrong and physically absurd. Does this have anything to do with how I'm using my trapz()? Would greatly appreciate help on the matter.

Thanks in advance!

Answer 1

[Star Strider]

The computations for ‘TA’, ‘Tb’, ‘xa’, and ‘xb’ don’t change w.r.t. the loop parameters, so it’s safe to evaluate them once before the loop.

Doing that seems to produce the result you want:

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;
xab = xa.*xb;
for ii=1:length(C1)
    C = (2.*1i.*sqrt(C1(ii)))./(sqrt(gam).*(xab+C1(ii)));
    D = 2*xa./(sqrt(gam).*(xab+C1(ii)));
    MC = (abs(C)).^2;
    MD = (abs(D)).^2;
    Sbb = 2.*pi.*(MC.*(n1+0.5)+MD.*(nm+0.5));
    A = trapz(w,Sbb)./(2.*pi);
    AA(ii)=A;
end

I evaluated this with:

C1 = linspace(0,500,1000);

and it seems to behave itself. It’s not obvious to me why putting those variable assignments before the loop also changes the integration results, since I don’t see anything in the loop that changes those variables as a function of the iteration variable.

That said, it also takes a while (on my laptop, I didn’t try it on my Ryzen 7 1800 desktop), so reducing the length of ‘w’ and possibly also ‘C1’ would be nice.

Check the code with my changes to see if it does what you want it to do.

[Alvin]

Thanks! That works! I for one do not know why it's behaving the way it is either given that w is independent of the loop so naturally all the other functions that are only dependent on w should not change as I loop over C1. But this is what I needed. Thanks again!

[Star Strider]

As always my pleasure!

If you find out what the problem was (I couldn’t, although I did my best to figure it out), please follow up.