Monte carlo simulation code/ unreasonable results in averaging.
3 次查看(过去 30 天)
显示 更早的评论
Hello,
I am trying to write the following monte carlo simulation code. But, something is wrong and I really need help to fix it.
What I am trying to do is: for every value of k, I am calculating the rate for each one of the relays (MM) and I am selecting the one with the highest rate. For every value of k, I am caluculating evrything for N realizations and I am averaging the value at the end. However, I get the same value for all the values of k and even the value i get is unreasonable! I am sure that i am doing something wrong in the averaging line but I do not know how to fix it!
Any help is appreciated!
function [br]= Rate_optimum(k)
PL=3;
N=10^5;
MM=3;
si=1;
PT=10.^(5./10);
PR=10.^(0./10);
sign=10.^(0./10);sigd=10.^(0./10);sige=10.^(-7./10);
lamdap=0.1;
Rate_opt2 =0;
for kk=1:N
fprintf('Running %d per %d \n',kk,length(k));
w_opt=PR;
for bb = 1 : MM %relays
f= exprnd(1./lamdap);
h=exprnd(1./lamdap);
g= exprnd(1./lamdap);
alpha_opt=(k.*si.^(-1))./(PT.*f+(PR.*h.*(PT.*g+sign))+2.*sige);
SNR1_opt=PT.*f./(sige+sigd);
SNR2_opt=(w_opt.*h.*g.*PT)./(w_opt.*h.*sign+sige+sigd);
SNR_SD_opt=SNR1_opt+SNR2_opt;
Rate=0.5.*(1-alpha_opt).*log2(1+SNR_SD_opt);
if (Rate > Rate_opt2 )
Rate_opt2 = Rate ;
ID_max = bb % This is the selected relay from the set of relays we have
end
end
Rate_final=Rate_opt2;
Rate_opt2
end
br= Rate_final/N ;
br
end
%%%%
%The call function:
% %Test file: run this
% clear all
% close all
% k = 0:0.5:4;
% for kk=1:length(k)
% br(kk) = Rate_optimum(k(kk) );
% end
%
% semilogy(k,br,'o','LineWidth',1.5)
%
% grid on
0 个评论
回答(1 个)
David Hill
2021-3-23
Very confusing code.
function [br]= Rate_optimum(k)
PL=3;
N=10^5;
MM=3;
si=1;
PT=10.^(5./10);
PR=10.^(0./10);
sign=10.^(0./10);sigd=10.^(0./10);sige=10.^(-7./10);
lamdap=0.1;
Rate_opt2 =0;
for kk=1:N
fprintf('Running %d per %d \n',kk,length(k));
w_opt=PR;
for bb = 1 : MM %relays
f= exprnd(1./lamdap);
h=exprnd(1./lamdap);
g= exprnd(1./lamdap);
alpha_opt=(k.*si.^(-1))./(PT.*f+(PR.*h.*(PT.*g+sign))+2.*sige);%k is an array so alpha_out will also be an array
SNR1_opt=PT.*f./(sige+sigd);
SNR2_opt=(w_opt.*h.*g.*PT)./(w_opt.*h.*sign+sige+sigd);%everything else is just constant and can be taken outside the loops
SNR_SD_opt=SNR1_opt+SNR2_opt;
Rate=0.5.*(1-alpha_opt).*log2(1+SNR_SD_opt);%Rate is also be an array
if (Rate > Rate_opt2 ) %does not make sense comparing an array to a scalar
Rate_opt2 = Rate ;
ID_max = bb % This is the selected relay from the set of relays we have
end
end
Rate_final=Rate_opt2;
Rate_opt2
end
br= Rate_final/N ;
br
end
How are you differentiating between the different relays? I believe there should be three values for each of your constants (different values for each relay). You never index with bb; therefore, all the relays have the same constants and the rates will be the same. You should be able to do this without any loops.
2 个评论
David Hill
2021-3-23
I believe this is what you are looking for.
N=10^5;
MM=3;
si=1;
PT=10^(5/10);
PR=10^(0/10);
w_opt=PR;
sign=10^(0/10);sigd=10^(0/10);sige=10^(-7/10);
lamdap=0.1;
Rate_opt2 =0;
f=exprnd(1./lamdap,MM,N);
g=exprnd(1./lamdap,MM,N);
h=exprnd(1./lamdap,MM,N);
SNR1_opt=PT*f/(sige+sigd);%MMxNN matrix
SNR2_opt=(w_opt*h.*g*PT)./(w_opt*h*sign+sige+sigd);%MMxNN matrix
SNR_SD_opt=SNR1_opt+SNR2_opt;%MMxNN matrix
alpha_opt=si^(-1)./(PT*f+(PR*h.*(PT*g+sign))+2*sige);%changed alpha_opt so to be constant (MMxNN matrix)
Rate_final=zeros(size(h));
ID_max=zeros(size(h));
for K=1:length(k)%only need to loop the values of k
Rate=0.5*(1-k(K)*alpha_opt).*log2(1+SNR_SD_opt);%MMxNN matrix
[Rate_final(K,:)]=max(Rate);%1xNN array
end
[~,ID_max]=max(Rate);%1xNN array (ID_max does not change with k)
br=sum(Rate_final)/N;%1xNN array
另请参阅
类别
在 Help Center 和 File Exchange 中查找有关 Operating on Diagonal Matrices 的更多信息
产品
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!