离散信道容量迭代算法

哥花了一点小小时间写的，望各位不吝赐教！

//===============================================
//made by : 何盛昌
//number  : 496592219

//===============================================
#include <stdlib.h>
#include <iostream.h>
#include <math.h>

void main()
{
//===============================================
//确定输入及输出符号的个数
//因为一般先验概率是取等概率，所以为1/num_in
//===============================================
 int num_in,num_out;
 int i,j;
 cout<<"输入符号的个数为："<<endl;
 cin>>num_in;
 cout<<"输出符号的个数为："<<endl;
 cin>>num_out;
 cout<<"输入概率为："<<endl;
 float *P_in;
 P_in=new float[num_in];   //输入的符号的概率
 float *P_out;
 P_out=new float[num_out];  //输出的符号的概率
 for (i=0;i<num_in;i++)
 {
  P_in[i]=1.0/(float)num_in;
  cout<<P_in[i]<<" ";
 }
 cout<<endl;
//================================================
//再输入信道转移概率矩阵
//注意进行判断，某行如果总概率大于1，则出现错误应从新输入
//================================================
 cout<<"输入转移概率："<<endl;
 float **p_ji;                   //转移条件概率
 p_ji=new float *[num_in];
 for (i=0;i<num_in;i++)
 {
  p_ji[i]=new float[num_out];
 }
 
 for (i=0;i<num_in;i++)
 {
  for (j=0;j<num_out;j++)
  {
   cin>>p_ji[i][j];
  }
 }
 for (i=0;i<num_in;i++)
 {
  float validate=0.0;
  for (j=0;j<num_out;j++)
  {
   validate+=p_ji[i][j];
  }
  if (validate>1.000001||validate<0.999999)
  {
   cout<<"输入数据有误，请检查后再次输入。"<<endl;
   exit(-1);
  }
 }
 
 float **p_ij;     //反条件概率
 p_ij=new float *[num_in];
 for (i=0;i<num_in;i++)
 {
  p_ij[i]=new float[num_out];
 }
 float C_Pre,C;     //当前信道容量和前一次循环信道容量
 C=10.0;
 double Pe=0.000001;     //两次信道容量相差的阈值
 int r=0;      //迭代次数
 float *p_up;     //计算第r+1次循环输入分布p_in分子
 p_up=new float[num_in];
 float p_down;     //计算第r+1次循环输入分布p_in分母
 do
 {
  r++;
//================================================
//求第r次循环反条件概率p_ij
//================================================
  for (j=0;j<num_out;j++)
  {
   P_out[j]=0.0;
   for (i=0;i<num_in;i++)
   {
    P_out[j]+=p_ji[i][j]*P_in[i];     //
   }
   if (P_out[j]>=0.000001)
   {
    for (i=0;i<num_in;i++)
    {
     p_ij[i][j]=P_in[i]*p_ji[i][j]/P_out[j];
    }
   }
   else
   {
    for (i=0;i<num_in;i++)
    {
     p_ij[i][j]=0.0;
    }
   }
  }
//================================================
//求第r+1次循环输入分布p_in
//================================================
  p_down=0.0;
  for (i=0;i<num_in;i++)
  {
   p_up[i]=0.0;
   for (j=0;j<num_out;j++)
   {
    if (p_ij[i][j]>=0.000001)
    {
     p_up[i]+=p_ji[i][j]*log(p_ij[i][j])/log(2.0);
    }
   }
   p_up[i]=pow(2.0,p_up[i]);
   p_down+=p_up[i];
  }
  for (i=0;i<num_in;i++)
  {
   P_in[i]=p_up[i]/p_down;
  }
//================================================
//求C(r+1)
//================================================
  C_Pre=C;
  C=log(p_down)/log(2.0);
  cout<<"第"<<r<<"次的容量为："<<C<<endl;
 }
 while (fabs(C-C_Pre)/C>Pe);
 cout<<"迭代的次数为："<<r<<endl;
}