isa.cpp

一个模拟退火算法的改进程序

		}
		else
		{
			m2=m2+1;
			f=f+D_cb2-D_cb1;
		}
		D_cb1=D_cb2;
	}

	//求f均值
	f=f/m;

	//计算t0
	t0=f/float(log(float(m2/(m2*p-m1*(1-p)))));
	return t0;
}

//
//SA算法程序
//
void SA(float* codebook_SA)
{
	//
	//变量定义
	//
	int n;						//扰动计数器n
	int t;						//最优解连续扰动不变计数器
	int s;						//当前温度连续扰动不变计数器
	int rand1;
	int rand2;
	float CB_temp[CB_ORDER];
	float CB_percent_temp;
	float D_cb1;				//当前的失真量
	float D_cb2;				//扰动后的失真量
	float mindist_D_cb;			//最小的失真量
    float globmindist_D_cb;		//最小的失真量
	float chgprobability;		//概率变量
	float randomnum;			//产生的验证概率
	
	float Tm;					//当前温度
	
	float aver1_D_cb;			//当前温度的失真量的总和
	float aver2_D_cb;			//当前温度的失真量平方的总和
    float fc_D_cb;				//当前温度的失真量方差


	int i;						
	int j;
	int p;						
	int q;

	int k;						//flag


	//变量初始化
	//
	
	

	//设置初始温度值
	Tm=inti_t(codebook_SA,0.8,10000);
    
	//计算当前码序列的失真量
	D_cb1 = distortion(codebook, CB_percent);

    //初始化最小的失真量
	mindist_D_cb = 2000;
	globmindist_D_cb=2000;
    do
	{
		t = 0;
		while(t <= t_UNCHAGE)
		{
			//初始化温度的失真量的总和和失真量平方的总和
			aver1_D_cb=D_cb1;
			aver2_D_cb=D_cb1*D_cb1;
			//初始化s
			s = 0 ;

			//每次降温进行N次随机扰动
			for(i = 0; i < s_UNCHAGE; i++)
			{
				//随机产生两个交换的序号
				rand1 = random();
				rand2 = random();

				//交换序号所对应的矢量
				for(j = 0; j < CB_ORDER; j++)
				{
					CB_temp[j] = codebook[rand1][j];
					codebook[rand1][j] = codebook[rand2][j];
					codebook[rand2][j] = CB_temp[j];
				}
				//交换序号所对应矢量先验概率
				CB_percent_temp = CB_percent[rand1];
				CB_percent[rand1] = CB_percent[rand2];
				CB_percent[rand2] = CB_percent_temp;

				//计算新的失真量
				D_cb2 = distortion(codebook, CB_percent);

				//判断当前扰动是否被接受
				if(D_cb1 > D_cb2)
				{
					s=s+4;
					//更新当前最小的失真量
					D_cb1 = D_cb2;
					//计算当前温度的失真量的总和和失真量平方的总和
					aver1_D_cb=aver1_D_cb+D_cb1;
					aver2_D_cb=aver2_D_cb+D_cb1*D_cb1;

					if(mindist_D_cb > D_cb1)
					{
						//更新最小的失真量
						mindist_D_cb=D_cb1;

						//更新最小失真码本矢量
						for(p = 0; p < CB_SIZE; p++)
						{
							for(q = 0; q < CB_ORDER; q++)
							{
								mindist_codebook[p][q]=codebook[p][q];

							}

						}

						//更新最小失真码矢量的先验分布概率
						for(p = 0; p < CB_SIZE; p++)
						{
							mindist_CB_percent[p]=CB_percent[p];

						}

					}

					printf("temp=%f; dist=%f ; mindist=%f ; t=%d ;\n",Tm, D_cb1,mindist_D_cb,t);
				}
				else
				{

					//计算当前接受概率
					chgprobability = exp((D_cb1-D_cb2)/Tm);
					//产生随机数，0到1均匀分布
					randomnum = rand()*1.0/32767;
					if(randomnum > chgprobability)
					{
						s++;
						//还原原来的码书序号，交换序号所对应的矢量
						for(j = 0; j < CB_ORDER; j++)
						{
							CB_temp[j] = codebook[rand1][j];
							codebook[rand1][j] = codebook[rand2][j];
							codebook[rand2][j] = CB_temp[j];
						}
						//还原原来的码书先验概率，交换序号所对应矢量先验概率
						CB_percent_temp = CB_percent[rand1];
						CB_percent[rand1] = CB_percent[rand2];
						CB_percent[rand2] = CB_percent_temp;
								
					}
					else
					{
						if(D_cb1 != D_cb2)
						{
							s=s+2;
							//更新当前最小的失真量
							D_cb1 = D_cb2;

							//计算当前温度的失真量的总和和失真量平方的总和
							aver1_D_cb=aver1_D_cb+D_cb1;
							aver2_D_cb=aver2_D_cb+D_cb1*D_cb1;

						}
						else
						{
						
						}
						printf("temp=%f; dist=%f ; mindist=%f ; t=%d ;\n",Tm, D_cb1,mindist_D_cb,t);
					}			
				}
				if(s>N)
				{
					goto loop;
				}
			}
					
	loop:
			
			//计算当前温度的失真量方差
			fc_D_cb=float(sqrt(double(aver2_D_cb/N-(aver1_D_cb/N)*(aver1_D_cb/N))));
			//降温
			Tm=Tm/(1+Tm*float(log(1+1))/(3*fc_D_cb));

			if(mindist_D_cb == globmindist_D_cb)
			{
				t++;

			}
			else
			{
				t = 0;

			}

			globmindist_D_cb = mindist_D_cb;
	
		}
			//更新失真码本矢量
		for(p = 0; p < CB_SIZE; p++)
		{
			for(q = 0; q < CB_ORDER; q++)
			{
				codebook[p][q]=mindist_codebook[p][q];

			}

		}

		//更新失真码矢量的先验分布概率
		for(p = 0; p < CB_SIZE; p++)
		{
			CB_percent[p]=mindist_CB_percent[p];

		}


	}
	while(abs(D_cb2-mindist_D_cb)>df);

	
	//训练好的码书赋值
	for(i = 0; i < CB_SIZE; i++)
	{
		for(j = 0; j < CB_ORDER; j++)
		{
			*(codebook_SA+CB_ORDER*i+j) = mindist_codebook[i][j];
		
		}
	}







}