gusucode.com > m4k20功率分配图 mimo功率优化,OFDM典型仿真程序，贪婪算法仿真程序 > code/mimo功率优化/opt_mutual_mimo.m

    clear all;
close all;
clc;

disp('Initializing...');

%--------------------PARAMETERS--------------------
% change these parameters to run different simulations
N = 4;
M = 4;
NSAMPLE = 100;

SNR_dB = [-10:5:15]; % SNR in dB
K = 5;  % channel K factor

random_data = 0; % generate random channel mean and correlation matrices
spec_case = 0; % 1 is no mean, 2 is no tx corr

% setting up various file names
if spec_case == 1
  suffix = 'corr';
elseif spec_case == 2
  suffix = 'mean';
else
  suffix = 'both';
end


%--------------------PROGRAM--------------------
I = eye(N);
SNR = 10.^(SNR_dB/10);

  Rt_half = randn(N,N) + j*randn(N,N)%;
  Rt = Rt_half*Rt_half'%;
  Rt = Rt/(K+1)%;
  
  H0 = randn(M,N) + j*randn(M,N)%;
  H0 = H0/norm(H0, 'fro')*sqrt(K*M*N/(K+1))%;  
end

if (spec_case == 2)
  % check for mean only case
  Rt = I*real(trace(Rt))/N; % this is to retain the same K factor
elseif (spec_case == 1)
  % check for correlation only case
  Rt = Rt/trace(Rt)*N; % full correlation
  H0 = zeros(M,N);
end

Rt_sqrt = sqrtm(Rt);


S0 = (H0'*H0 + M*Rt);
[U lambda] = eig(S0);
lambda = real(diag(lambda));


lambda_equal = ones(N,1)/N;  

disp('Calculating Jensen''s upper bound on average mutual information...');
for k = 1:length(SNR_dB)
  disp(['  SNR = ', num2str(SNR_dB(k)), 'dB']);  
  lambda_jen_wf(:,k) = waterfill(SNR(k), lambda);
  mutual_jen_wf(k) = ...
      - mutual_mimo_func_diag(M,N,H0*U,U'*Rt_sqrt*U,...
                              SNR(k),lambda_jen_wf(:,k),NSAMPLE,I);
  mutual_equal(k) = ...
      - mutual_mimo_func_diag(M,N,H0,Rt_sqrt,...
                              SNR(k),lambda_equal,NSAMPLE,I); 
end

figure;
plot(SNR_dB, mutual_jen_wf/log(2), 'b+-');
hold on;
plot(SNR_dB, mutual_equal/log(2), 'r:s');
xlabel('SNR(dB)');
ylabel('频谱效率(bps/Hz)');
legend('功率优化算法', '等功率分配算法',2);


figure;
plot(SNR_dB, lambda_equal*ones(1, length(SNR_dB)), 'r:s');
hold on;
plot(SNR_dB, lambda_jen_wf, 'b+-');
xlabel('SNR(dB)');
ylabel('功率分配');
legend('等功率分配算法','功率优化算法',2);