% Boyd & Vandenberghe "Convex Optimization"
% Joëlle Skaf - 04/24/08
%
% We consider a discrete memoryless communication channel, with input
% X(t) \in {1,...,n}, and output Y(t) \in {1,...,m}, for t = 1,2,...
% The relation between the input and output is given statistically:
%           p_ij = Prob(Y(t)=i|X(t)=j), i=1,...,m,  j=1,...,n
% The matrix P is called the channel transition matrix.
% The channel capacity C is given by
%           C = sup{ I(X;Y) | x >= 0, sum(x) = 1},
% I(X;Y) is the mutual information between X and Y, and it can be shown
% that:     I(X;Y) = c'*x - sum_{i=1}^m y_i*log_2(y_i)
% where     c_j = sum_{i=1}^m p_ij*log_2(p_ij), j=1,...,m

% Input data
rand('state', 0);
n = 15;
m = 10;
P = rand(m,n);
P = P./repmat(sum(P),m,1);
c = sum(P.*log2(P))';

% Channel capacity
cvx_begin
    variable x(n)
    y = P*x;
    maximize (c'*x + sum(entr(y))/log(2))
    x >= 0;
    sum(x) == 1;
cvx_end
C = cvx_optval;

% Results
display(['The channel capacity is: ' num2str(C) ' bits.'])

 
Successive approximation method to be employed.
   SDPT3 will be called several times to refine the solution.
   Original size: 45 variables, 21 equality constraints
   10 exponentials add 80 variables, 50 equality constraints
-----------------------------------------------------------------
 Cones  |             Errors              |
Mov/Act | Centering  Exp cone   Poly cone | Status
--------+---------------------------------+---------
 10/ 10 | 2.785e+00  4.909e-01  0.000e+00 | Solved
 10/ 10 | 8.556e-02  5.525e-04  0.000e+00 | Solved
 10/ 10 | 8.090e-03  4.960e-06  0.000e+00 | Solved
  8/ 10 | 6.944e-04  3.574e-08  0.000e+00 | Solved
  0/  5 | 6.072e-05  9.724e-11  0.000e+00 | Solved
-----------------------------------------------------------------
Status: Solved
Optimal value (cvx_optval): +0.296291
 
The channel capacity is: 0.29629 bits.