unrecognized function or variable in parfor

11 次查看（过去 30 天）
Maurilio Matracia
Maurilio Matracia 2021-9-17
评论： Image Analyst 2021-9-20
Hello everyone!
I am running the following code and getting this error:
Error using Coverage_chetlur>(parfor supply)
Unrecognized function or variable 'ZZ'.
Error in Coverage_chetlur (line 173)
parfor B = 1:2
Could you please help me fixing it? I am not sure why it appears since I am defining ZZ after having initialized the parfor loop.
PS: Please note that the code works perfectly if I use the for loop instead of the parfor one.
tic
Pi = 3.14159 ;
nIntervals = 3 ;
nPoints = 3 ;
%% PARAMETERS
lambdaT = 0.5e-5 ;
alpha = [3.5 , 3.5] ;
m = [2 , 1] ;
iterations = 1e0 ;
p= [1.5, 10] ;
NA = 2 % Better not to exceed NA=20
Rd = 5e2 ;
R_T = Rd + 19e3 ;
infty = 1e9 ;
zero = 1e-6 ;
tau = 0.31623 ;
sigma_n2 = 1e-12 ;
R_U = linspace(0,1, 10) * Rd ; R_U = 0.7 * Rd ;
R_Plus = Rd+R_U-zero ;
R_Minus = Rd-R_U+zero ;
H = 1e2 ;
urbL = 3 ;
eta = [0.9772 0.007943 1 ; 0.7943 0.01 1 ; 0.6918 0.005 1 ; 0.5888 0.0004 1] ;
eta = [0.9772 1 ; 0.7943 1 ; 0.6918 1 ; 0.5888 1] ; eta(:,2) = eta(:,1) ; %%% NOT CONSIDERING NLoS
Sa = [4.88 9.6117 12.0810 27.304]' ; Sb = [0.429 0.1581 0.1139 0.0797]' ; Sa = Sa(urbL) ; Sb = Sb(urbL) ;
xi = eta(urbL,:) .* p ;
D{1,1} = @(z) sqrt(z.^2+H^2) ;
D{1,2} = @(z) (xi(2)/xi(1)).^(1/alpha(2)) .* ( sqrt(z.^2+H^2).^(alpha(1)/alpha(2)) ) ;
D{2,1} = @(z) (xi(1)/xi(2)).^(1/alpha(1)) .* z.^(alpha(2)/alpha(1)) .* (z>D{1,2}(H)) ...
+ H .* (z<= D{1,2}(H)) ;
D{2,2} = @(z) z ;
for B = [1,2]
Z{B,1} = @(z) sqrt(D{B,1}(z).^2-H^2) ;
Z{B,2} = @(z) D{B,2}(z) ;
end
for u = 1:length(R_U)
Ru = R_U(u) ; Rplus = Rd+Ru-zero ; Rminus = Rd-Ru+zero ;
for it = 1 : iterations
thetaA = 2*pi*rand(NA,1) ;
rA = sqrt( rand(NA,1)*(Rd^2 - 0^2) + 0^2 ) ;
coordA = rA .* [cos(thetaA), sin(thetaA)] ;
Zs{1} = sqrt( (Ru-coordA(:, 1)).^2 + coordA(:, 2).^2 ) ;
ED{1,it} = sqrt(Zs{1}.^2 + H^2) ;
EDmin{1}(it) = min(ED{1,it}) ;
coord{1} = coordA ;
nT = poissrnd (lambdaT * pi * (R_T^2-Rd^2)) ;
thetaT = rand(nT,1) * 2*pi ;
rT = sqrt(rand(nT,1)*(R_T^2-Rd^2) + Rd^2) ;
coord{2} = [ rT.*cos(thetaT) , rT.*sin(thetaT) ] ;
zT = sqrt( (Ru-coord{2}(:, 1)).^2 + coord{2}(:, 2).^2 ) ;
ED{2,it} = sqrt( (Ru-coord{2}(:, 1)).^2 + coord{2}(:, 2).^2 ) ; Zs{2} = ED{2,it} ;
EDmin{2}(it) = min(ED{2,it}) ;
end
for B = [1,2]
mu{B} = @(z) m(B).* tau./xi(B).* D{B,B}(z).^alpha(B) ;
end
%% COVERAGE SIMULATION
for it = 1:iterations
for B = [1,2]
pTAG(B) = xi(B) * EDmin{B}(it)^-alpha(B) ; % Conditional received power
Hc{B} = gamrnd( m(B),1/m(B), 1,length(ED{B,it}) ) ;
end
tag(it) = 1*(pTAG(1)>pTAG(2)) + 2*(pTAG(2)>pTAG(1)) ;
receivedPower = Hc{tag(it)}( find(ED{tag(it), it}==EDmin{tag(it)}(it)) ) * xi(tag(it))*EDmin{tag(it)}(it)^-alpha(tag(it)) ;
interference(it) = sum(Hc{1}'.*xi(1).*ED{1,it}.^-alpha(1)) + sum(Hc{2}'.*xi(2).*ED{2,it}.^-alpha(2)) - receivedPower ;
SINR(it) = receivedPower/(interference(it)+sigma_n2) ;
end
Pc_s(u) = sum(SINR >= tau) / iterations ;
%% ANALYSIS
Betai = @(z) real( (z>0).*acos( (z.^2 + Ru.^2 - Rd.^2) ./ ( (z>0).*2.*z.*Ru + (z==0)) ) ) ;
dBetai = @(z) (z>0 ).*(Ru.^2 - z.^2 - Rd.^2) ./ ( (z>0) .* (z .* sqrt(4*z.^2 .* Ru.^2 - (z.^2 + Ru.^2 - Rd.^2).^2)) + (z==0) ) ;
zX = @(beta) Ru.*cos(beta) + sqrt(Rd.^2 - Ru.^2.*sin(beta).^2) ;
rW = @(w,beta) sqrt(Ru^2 + w.^2 - 2*Ru.*w.*cos(beta)) ;
A0 = Pi* Rd^2 ;
IND{1} = @(z) z <= Rminus ;
IND{2} = @(z) (Rminus < z) & (z < Rplus) ;
IND{3} = @(z) z >= Rplus ;
Sigma_i{1} = @(z) pi*z.^2 ;
Sigma_i{2} = @(z) Betai(z).*z.^2 + integral(@(beta) zX(beta).^2, Betai(z),pi,'ArrayValued',1, 'RelTol',1e-4, 'AbsTol',1e-6) ;
dSigma_i{1} = @(z) 2*pi*z ;
dSigma_i{2} = @(z) 2*z.*Betai(z) + z.^2 .* dBetai(z) - zX(Betai(z)).^2 .* dBetai(z) ;
Sigma = @(z) Sigma_i{1}(z).*IND{1}(z) + Sigma_i{2}(z).*IND{2}(z) + (A0+zero) .* IND{3}(z) ;
dSigma = @(z) dSigma_i{1}(z).*IND{1}(z)+dSigma_i{2}(z).*IND{2}(z) ;
F_OmegaA= @(w) Sigma(w) / (pi*Rd^2) ; Fbar_OmegaA = @(w) 1-F_OmegaA(w) ; f_OmegaA = @(w) dSigma(w) / A0 ;
Fbar_Z{1} = @(z_) arrayfun(@(z) Fbar_OmegaA(z).^NA, z_) ;
Fbar_Z{2} = @(z_) arrayfun(@(z) exp( - lambdaT .* integral(@(beta) integral(@(w) (rW(w,beta)>Rd) .* w, 0, z, 'arrayvalued',1,'RelTol',1e-5, 'AbsTol',1e-8), 0,2*pi, 'arrayvalued',1,'RelTol',1e-5, 'AbsTol',1e-8) ) , z_) ;
F_Z{1} = @(z_) arrayfun(@(z) 1-Fbar_Z{1}(z), z_) ;
f_Z{1} = @(z_) arrayfun(@(z) NA * Fbar_OmegaA(z).^(NA-1) .* f_OmegaA(z) , z_) ;
F_Z{2} =@(z_) arrayfun(@(z) 1-Fbar_Z{2}(z), z_) ;
f_Z{2} = @(z_) arrayfun(@(z) lambdaT .* Fbar_Z{2}(z) .* z .* integral(@(beta) (rW(z,beta)>Rd), 0,2*pi,'arrayvalued',1,'RelTol',1e-5, 'AbsTol',1e-8), z_) ;
%% LAPLACE TRANSFORMS
f_hatOmegaA = @(w,z) dSigma(w) ./ (A0-Sigma(z)) ;
I_{1} = @(s,z) ( m(1) ./ (m(1) + xi(1) * s .* D{1,1}(z).^-alpha(1)) ).^m(1) ;
I_{2} = @(s,w) 1 - ( m(2) / (m(2)+xi(2) * s .* D{2,2}(w).^-alpha(2)) ).^m(2) ;
for B = [1,2]
hatUps{B,1} = @(s_,z_) arrayfun(@(s,z) integral(@(w) I_{1}(s,w) .* f_hatOmegaA(w,Z{B,1}(z)), Z{B,1}(z), Rplus,'arrayvalued',1, 'RelTol',1e-4, 'AbsTol',1e-8), s_,z_) ;
LapBC{B,1} = @(s_,z_) arrayfun(@(s,z) hatUps{B,1}(s,z).^(NA-(B~=2)), s_,z_) ;
LapBC{B,2} = @(s_,z_) arrayfun(@(s,z) exp( -lambdaT * integral(@(beta) integral(@(w) (rW(w,beta) > Rd) .* I_{2}(s,w) .* w, Z{B,2}(z),infty, ...
'arrayvalued',1, 'RelTol',1e-4, 'AbsTol',1e-8), 0,2*pi,'arrayvalued',1, 'RelTol',1e-4, 'AbsTol',1e-8) ) , s_,z_) ;
LapIII{B,2} = @(s,z) exp(-2*pi*lambdaT * integral(@(x) I_{2}(s,x) .* x, Z{B,2}(z),infty,'arrayvalued',1) ) ;
LapII{B,2} = @(s,z) LapIII{B,2}(s,z) .* exp( lambdaT * integral(@(beta) integral(@(x) I_{2}(s,x) .* x, Z{B,2}(z),zX(beta),'ArrayValued',1), -Betai(Z{B,2}(z)),Betai(Z{B,2}(z)),'ArrayValued',1) ) ;
% LapI{B,2} = @(s,z) LapIII{B,2}(s,z) .* exp( lambdaT * integral(@(beta) integral(@(x) I_{2}(s,x) .* x, Z{B,2}(z),zX(beta),'ArrayValued',1), -pi,pi,'ArrayValued',1) ) ;
LapBC{B,2} = @(s,z) LapIII{B,2}(s,z).*(Z{B,2}(z)>=Rplus) + LapII{B,2}(s,z).*( (Rminus<=Z{B,2}(z)) & (Z{B,2}(z)<Rplus) ) ; %+ LapI{B,2}(s,z).*(Z{B,2}(z)<=Rd-Ru(u)) ;
LapB{B} = @(s_,z_) arrayfun(@(s,z) LapBC{B,1}(s,z) .* LapBC{B,2}(s,z) , s_,z_) ;
end
%% ASSOCIATION PROBABILITIES
R_B = [zero,Rplus ; zero,Rplus; Rminus,infty] ;
a{1} = @(z) Fbar_Z{2}(Z{1,2}(z)) ;
a{2} = @(z) Fbar_Z{1}(Z{2,1}(z)) ;
for B = [1,2]
A(B) = integral(@(z) f_Z{B}(z) .* a{B}(z), R_B(B,1), R_B(B,2), 'arrayvalued',1,'RelTol',1e-4, 'AbsTol',1e-8) ;
%% COVERAGE PROBABILITY
Pcond{B} = @(z_) arrayfun(@(z) 0, z_) ;
epsilon(B) = (factorial(m(B)))^(-1/m(B)) ;
for k = 1:m(B)
Pcond{B} = @(z_) arrayfun(@(z) Pcond{B}(z) + nchoosek(m(B),k) * (-1)^(k+1) * exp(-k * epsilon(B) * mu{B}(z) * sigma_n2) .* LapB{B}(k * epsilon(B) *mu{B}(z),z), z_) ;
end
end
%% ANALYTICAL INTEGRATION
Pc_a(u) = 0 ;
%% MONTECARLO INTEGRATION
for B = 1:2
zAxis{B} = logspace( log10(R_B(B,1)), log10(R_B(B,2)), nIntervals) ;
prod{B} = zeros(1, nPoints) ;
F_zAxis{B} = F_Z{B}(zAxis{B}) ;
parfor nP = 1 : nPoints
nP
U = rand ;
diff = abs(U - F_zAxis{B}) ;
mindiff = min(diff) ;
index = find(diff==min(diff)) ;
ZZ{B}(nP) = zAxis{B}(index(1)) ; % We use index(1) because there might be many minimum values
prod{B}(nP) = a{B}( ZZ{B}(nP) ) .* Pcond{B}( ZZ{B}(nP) ) ;
end
pc_a{B} = sum(prod{B}) / nPoints ;
end
Pc_a(u) = pc_a{1} + pc_a{2}
end
toc

请先登录，再回答此问题。

回答（2 个）

Raymond Norris
Raymond Norris 2021-9-17
I don't have much time to look at this closer, so I'll give a couple of quick thoughts
  1. For starters, I would suggest going through and addressing the Code Analyzer prompts (the orange markings). There are quite a few cases of preallocation that ought to be done.
  2. Without running your code, I can't tell for sure, but I don't think you need to use cell arrays as much as you are.
  3. To address the parfor, I would suggest refactoring your code as such
%% MONTECARLO INTEGRATION
prod = zeros(2, nPoints);
for B = 1:2
zAxis{B} = logspace( log10(R_B(B,1)), log10(R_B(B,2)), nIntervals) ;
F_zAxis{B} = F_Z{B}(zAxis{B}) ;
parfor nP = 1 : nPoints
prod(B,nP) = unit_of_work(F_zAxis,zAxis,a,Pcond,B,nP) ;
end
pc_a{B} = sum(prod{B}) / nPoints ;
end
Pc_a(u) = pc_a{1} + pc_a{2}
Notice a couple of things
  • I'm switching prod from a cell array to a numeric array
  • I'm preallocating prod
  • I'm putting everything in the parfor into a subfunction (see below)
  • I didn't change how B is indexed (e.g. prod{B}), you'll need to do that if you stick with numeric arrays
  • I would suggest a different variable name then prod, as this referrs to the built-in function for multiplication of elements of a variable
unit_of_work is defined as follows
function result = unit_of_work(F_zAxis,zAxis,a,Pcond,B,nP)
nP
U = rand ;
diff = abs(U - F_zAxis{B}) ;
index = find(diff==min(diff)) ;
ZZ{B}(nP) = zAxis{B}(index(1)) ; % We use index(1) because there might be many minimum values
result = a{B}( ZZ{B}(nP) ) .* Pcond{B}( ZZ{B}(nP) ) ;
end
You might have to clean something up a bit (i.e. cell array to numeric), but this ought to get you started.
  2 个评论
Maurilio Matracia
Maurilio Matracia 2021-9-19
Thank you very much. Actually I just had to initialize ZZ and zAxis to make it work, do you know why?
Image Analyst
Image Analyst 2021-9-20
How would it know what ZZ is otherwise? How would it suppose to know that ZZ is a cell array with at least B cells? And if you're going to assign the nP index to the contents of the B'th cell, then the B'th cell must have something in it already.
See the FAQ to get a good intuitive feel for what a cell array is.
https://matlab.fandom.com/wiki/FAQ#What_is_a_cell_array.3F

请先登录，再进行评论。

Image Analyst
Image Analyst 2021-9-17
Evidently, it never gets to the line
ZZ{B}(nP) = zAxis{B}(index(1)) ; % We use index(1) because there might be many minimum values
Set a break point there and see if it does. Otherwise set other breakpoints and follow the execution path of your program by stepping through it:
Debugging in MATLAB | Doug's MATLAB Video Tutorials

类别

在 Help Center 和 File Exchange 中查找有关 Matrix Indexing 的更多信息

标签

Community Treasure Hunt

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

Start Hunting!

Translated by