Cumtrapz with cell arrays in a double loop

2 次查看(过去 30 天)
clear
%g = 1; %Coupling strength
k = 20; %Decay linewidth for cavity
no = 0; %input photon number for cavity
nm = 40; %input phonon number for mechanical oscillator
Del = -10^2; %Delta, laser detuning
Ohm = 10^2; %Omega, natural frequency of mechanical oscillator
gam = 0.0032; %Decay linewidth for mechanical oscillator
w = (-1199:3:1200);
g = (0:1.25:999);
xat = 1 - (1i.*(w + Del).*(2./k));
xbt = 1 - (1i.*(w - Ohm).*(2./gam));
ca_c = cell(length(g),length(w));
da_c = cell(length(g),length(w));
Sbb_c = cell(length(g),length(w));
for u=1:length(g)
for v=1:length(w)
C = sqrt((4.*(g(u).^2))./(k.*gam)); %Cooperativity
Ca(u) = C;
c = (2.*1i.*Ca(u))./(sqrt(gam).*(xat(v).*xbt(v)+(Ca(u).^2))); %power spectrum coefficient
ca_c{u,v}= c; %c-cell
d = (2.*xat(v))./(sqrt(gam).*(xat(v).*xbt(v)+(Ca(u).^2))); %power spectrum coefficient
da_c{u,v} = d; %d-cell
Sbb = (abs(ca_c{u,v}).^2).*(no + 0.5)+(abs(da_c{u,v}).^2).*(nm + 0.5);
Sbb_c{u,v} = Sbb; %Sbb-cell
A = cumtrapz(w,Sbb_c{u,v})
end
end
The code above wouldn't compile because of the cumtrapz function. The error says ORDER contains an invalid permutation index. I've ask something somewhat similar on MathWorks a while ago but it was not answered. https://www.mathworks.com/matlabcentral/answers/350535-on-double-looping-and-confusion-with-cell-arrays-and-cumtrapz-function
Maybe my explanation was too detailed that people didn't want to get into the details. In this case I am simply trying to find the area under the curve using cumtrapz for Sbb against w. cumptrapz has to stay in the loop in order for me to obtain an array of areas for different g values. Much help would be appreciated. Thank you.

回答(1 个)

dpb
dpb 2017-7-30
编辑:dpb 2017-7-30
Well, I've no idea what you're really trying to do, either, but setting a breakpoint at the line containing cumtrapz we find the arguments to it are:
K>> whos w
Name Size Bytes Class Attributes
w 1x800 6400 double
K>> whos Sbb
Name Size Bytes Class Attributes
Sbb 1x1 8 double
K>> Sbb_c{u,v}
ans =
7.6805e-08
K>>
I'll admit the error message is more than a little confusing as to what it means, but the problem is you're trying to integrate 800 points in the X vector with a single element for the corresponding Y.
Not sure why all the huge cell arrays; looks like you could simply compute an array of Y of the length of X (w) to integrate.
But, that's why it crashes; there's a big mismatch in the two arrays that must have 1:1 consonance in size in order to work...
ADDENDUM
Try
...
Lg=length(g); Lw=length(w);
Sbb = zeros(Lg,Lw);
A=Sbb;
C=sqrt((4*(g.^2))/(k*gam));
for u=1:Lg
for v=1:Lw
Ca=(2*1i*C(u))/(sqrt(gam)*(xat(v)*xbt(v)+C(u)^2));
Da=(2*xat(v))/(sqrt(gam)*(xat(v)*xbt(v)+(C(u)^2)));
Sbb(u,v) = (abs(Ca)^2)*(no + 0.5)+(abs(Da)^2)*(nm + 0.5);
end
A(u,:) = cumtrapz(w,Sbb(u,:));
end
that eliminates the cell arrays and the unneeded storage for the intermediaries--I think it is what you were looking for. Note your C is independent of v so doesn't need to be in the innermost loop on v and I went ahead and moved it out and precomputed the vector, computing the two coefficients dependent upon it on the fly.
Then the integral can be computed on each vector at the end of the inner loop.
Unless there's other need for it later, there's no reason to save the full 2D array, either; simply the 1D vector.
The inner loops could also be vectorized and the Sbb array computed without any loops at all via expanding the two x vectors to a grid via meshgrid.
I believe this addresses the comments/questions raised in the various responses to me and Walter.
Walter, thnx for the digging on the error from cumtrapz -- I was somewhat confused; didn't realize it would take a 2nd-position single value as the DIM argument without the [] placeholders...not one of the more prevalent functions I've called over the years! :)
ADDENDUM 2:
It's not that hard to vectorize over v altho following is air code; check thoroughly...
sqrtgam=sqrt(gam); % take constant outside loops at initialization section
xabt=xat.*xbt; % ditto--this is a constant vector, too...
CSq=C(u).^2; % and so is this...
...
for u=1:Lg
Ca=2*1i*C(u)./(sqrtgam*(xabt+Csq(u)));
Da=2*xat./(sqrtgam*(xabt+Csq(u)));
Sbb(u,:) = (abs(Ca).^2)*(no + 0.5)+(abs(Da).^2)*(nm + 0.5);
A(u,:) = cumtrapz(w,Sbb(u,:));
end
  3 个评论
Walter Roberson
Walter Roberson 2017-7-30
You have
cumtrapz(w,Sbb_c{u,v})
but your Sbb_c{u,v} is a scalar.
When you call cumptraz() with two arguments, cumtrapz() looks at the size of the second argument. In the case where it is a scalar, it assumes that it is the dimension number. So you are effectively calling
cumtrapz(w, [], Sbb_c{u,v})
but 7.6805e-08 is not a valid dimension number.
Alvin
Alvin 2017-7-30
编辑:Alvin 2017-7-30
Thank you for the response! I think I have remedied that, my current issue is with the looping, and it is of the following:
Y = Sbb_c(100,:);
Z = cell2mat(Y);
A = cumtrapz(w,Z);
I added the three lines above and being outside the double loops and
plot(w,Z)
clearly works. What I'm doing here is picking the 100th row out of my 800x800 Sbb_c cell array and doing a cumtrapz() of that row with w. My wish is to loop it so that it includes all 800 rows (not just the 100th).

请先登录,再进行评论。

类别

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

Community Treasure Hunt

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

Start Hunting!

Translated by