clc;
clear;
close all;
M=1;%感知目标数
N=4;%天线数
K=3;%用户数
sigma_c=10^(-7);%通信接收机噪声
sigma_s=3.16*10^(-6);%回波信号接收机处噪声
P_max=1;
beta=0.025;%感知精度
R_min=0.3;%通信用户最小数据速率
alpha=0.01;%感知的往返路径损耗
L=30;%通信帧长
d=0.5;%相邻天线之间的归一化间隔
zeta=0;%信道误差上界
h_k = [1.01000000000000 - 0.01390100000000000i;-1.0151000000000000 + 0.0007500000000000i;1.61000000000000 - 0.01390100000000000i];
A = ones(3,3); %定义A（theta）
derivative_A = zeros(3,3);%A矩阵求导
%初始化泰勒参数
bar_x2_k = [0.1;0.1;0.1];
bar_x3_k = [0;0;0];
R0=[];
R1=[];
R2=[];
I=[];
Wk=[];
for i=1:1:10
cvx_begin
cvx_solver mosek
variable Wk(K,K,K) complex hermitian;
variable t;
variables phi(K) chi(K) r_k(K) lamda1(K) lamda2(K)
variables x1_k(K) x2_k(K) x3_k(K) z(K);
expressions B C D;
%定义sum_Wk
sum_Wk = 0;
for k = 1:K
sum_Wk = sum_Wk + Wk(:,:,k);
end
sum_Wi = 0;
for s = 1:1:K
if s ~= k
sum_Wi = sum_Wi + Wk(:,:,s);
end
end
%目标函数
maximize t
subject to
%约束（35）
0 <= real(trace(sum_Wk)) <= P_max;
%约束 （36）
for k=1:K
r_k(k) >= R_min;
end
%约束(37)
sum(r_k(k)) >= t;
%约束（38）
Wk(:,:,k) == hermitian_semidefinite(K);
%约束（26）
for k=1:1:K
D = [real(trace(derivative_A’derivative_A.sum_Wk))-(betasigma_s^2)/(2L*(abs(alpha))^2), real(trace(derivative_A’Asum_Wk));…
real(trace(A’derivative_Asum_Wk)), real(trace(A’Asum_Wk))] == hermitian_semidefinite(2);
end
for k=1:1:K
%约束（30）
%(19)的第一个式子
x1_k(k) + rel_entr(1,phi(k)) <= 0;
%a-b>=c
x1_k(k) - x2_k(k) >= x3_k(k);
%(19)的第二个式子bate
chi(k) <= exp( bar_x2_k(k) )( x2_k(k) - bar_x2_k(k) + 1 );
%(19)的第三个式子
z(k) - 1 <= exp( bar_x3_k(k) )( x3_k(k) - bar_x3_k(k) + 1 );
log(2)*r_k(k) + rel_entr(1,z(k)) <= 0;
end
%约束（32）（33）
for k=1:1:K
B = [real(Wk(:,:,k)+lamda1(k).eye(K)), real(Wk(:,:,k)h_k);…
real(h_k’ Wk(:,:,k)), real(h_k’ Wk(:,:,k)*h_k-phi(k)-lamda1(k)*zeta^2)]==hermitian_semidefinite(K+1);
end
for k=1:1:K
C = [real(-sum_Wi+lamda2(k).eye(K)), real(-sum_Wih_k);…
real(-h_k’*sum_Wi), real(-h_k’sum_Wih_k-sigma_c^2-lamda2(k)*zeta^2+chi(k))] == hermitian_semidefinite(K+1);
end
cvx_end
%bar_b_2=b2 bar_b_3=b3
for k=1:1:K
bar_x2_k(k) = x2_k(k);
bar_x3_k(k) = x3_k(k);
end
R0=[R0 sum(r_k)];
I=[I i];
end
figure;
plot(I,R0(:),‘K-d’)
xlabel(‘迭代次数’);
ylabel(‘系统总吞吐量(bps/Hz)’);
