信道容量迭代算法

采用Matlab编程实现，程序代码如下：

%n代表行，m代表列
%p表示输入概率分布
%q表示输出概率分布
%C表示信道容量
%e为很小的正数
%k表示迭代次数
function [] = Capacity(P)
[n,m] = size(P);
p=ones(1,n)/n; p=p';
C=0; k=1;e=10^-5;
Cl=0;Cu=0;
while 1
    %(1)求q
    for j=1:m
        q(j,1) = 0;
        for i =1:n
            q(j,1) = p(i,1)*P(i,j)+q(j,1);
        end
    end
    %(2)求a
    for i=1:n
        lna(i,1)=0;
        for j=1:m
            middle = P(i,j)/q(j,1) 
            if middle ~= 0
            lna(i,1) = P(i,j)*log(middle)+lna(i,1);
            end
        end
        a(i,1) = exp(lna(i,1));
    end
    %(3)计算u
    u = 0;
    for i=1:n
        u = a(i,1)*p(i,1)+u;
    end
    %(4)计算Cl
    Cl = log2(u);
    Cu = log2(max(a));
    %(5)比较Cu和Cl
    if (Cu-Cl)<e
        break;
    end
    for i=1:n
        p(i,1) = p(i,1)*a(i,1)/u;
        k = k+1;
    end
end
    %(6)输出信道容量C 迭代次数k 输入概率分布p
    C=Cl
    p
    k