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

回答(1 个)

David Hill
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 个评论
Deema42
Deema42 2021-3-23
Hello,
Thank you for your reply!
the difference between the relays is in the value of the channels: f,g, and h. These are generated randomely so I beleive tha rate for each relay will be different. This is how I differntaite between relays. Will you please suggest a solution for this?
Thank you in advance.
David Hill
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 CenterFile Exchange 中查找有关 Operating on Diagonal Matrices 的更多信息

产品

Community Treasure Hunt

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

Start Hunting!

Translated by