main.cpp

自己做的小波与神经网络组合方法实现负荷预测程序

        ///////////////////////矩阵求逆///////////
        ////////////////////////////////////////// 
  double tran[R_d3][N_d3];    //转置矩阵   
  double jz_ji[R_d3][R_d3];     //矩阵乘积矩阵tran[i][k]*yacobi_d3[k][j]
  double I[R_d3][R_d3];    //单位矩阵  
  double jz_he[R_d3][R_d3];  //矩阵求和矩阵jz_ji[i][j]+I[i][j]
  double a[R_d3][R_d3];       //求逆矩阵jz_he[R_d3][R_d3]
  double b[R_d3][R_d3];       //逆矩阵 
  double mulmatrix[R_d3][N_d3];   ////b[R][R]*tran[R][N]
  double v_x[N_d3][1];     //各训练样本误差
  double jz[R_d3][1];    //mulmatrix[i][k]*v_x[k][j]    
        
 qiu_detx_d3()
{
  int i,j,k;  
  int m;


             ///////////////////求矩阵转置tran[R][N]
  for(i=0;i<N_d3;i++)
  {
    for(j=0;j<R_d3;j++) 
    {
	 tran[j][i]=yacobi_d3[i][j];
	
	}
  }  

    
               ////////////实矩阵相乘tran[R][N]*yacobi[N][R]
  
   for(i=0;i<R_d3;i++)          //矩阵初始化
	{
	  for(j=0;j<R_d3;j++)
	   {
		jz_ji[i][j]=0;
	   }
    }
  
         ////////////矩阵相乘
  for(i=0;i<R_d3;i++)
   for(j=0;j<R_d3;j++)
	{
     for(k=0;k<N_d3;k++)
     jz_ji[i][j]+=tran[i][k]*yacobi_d3[k][j];
	}
 
	  

             ////////////求正定矩阵I[R][R]
  for(i=0;i<R_d3;i++)    
  {
   for(j=0;j<R_d3;j++)
     {
	   if(i==j)
	      I[i][j]=1;
       else
		   I[i][j]=0;  
      }
  }

  
  for(i=0;i<R_d3;i++)
   for(j=0;j<R_d3;j++)
    I[i][j]=niu_d3*I[i][j];

                        ///////实矩阵相加
  for(i=0;i<R_d3;i++)
   for(j=0;j<R_d3;j++)
    jz_he[i][j]=jz_ji[i][j]+I[i][j];

                         //////////求逆
  
  for (i=0;i<R_d3;i++)
    for (j=0;j<R_d3;j++)
      a[i][j]=jz_he[i][j];


  for(i=0;i<R_d3;i++)
   {
     for(j=0;j<R_d3;j++) b[i][j]=0.0;
     b[i][i]=1.0;
   }
 
    
                     //////////////矩阵求逆
  if(gjcpegd3( 1, a, b ) == 1)
   {
     printf(" The linear system has'nt solution!\n");
     printf(" Strike any key to exit!\n"); getchar(); exit(1);
   }

               
               
           /////////////矩阵乘积////////////// 
   for(i=0;i<R_d3;i++)     //矩阵初始化
	{
	  for(j=0;j<N_d3;j++)
	   {
		mulmatrix[i][j]=0;
	   }
    }
           
  for(i=0;i<R_d3;i++)
   for(j=0;j<N_d3;j++)
	{
     for(k=0;k<R_d3;k++)
     mulmatrix[i][j]+=b[i][k]*tran[k][j];
	}                         
  
                                           //求det_x            (3)
  
  for(m=0;m<N_d3;m++)
   {
	  for(i=0;i<1;i++)               //ON即输出为1 ,求v(x)
	  {   
		  v_x[m][i]=err_O_d3[m][i];           /////////////误差T[k]-O[k];
	  }
  }


     /////////////矩阵乘积////////////// 
   for(i=0;i<R_d3;i++)     //矩阵初始化
	{
	  for(j=0;j<1;j++)
	   {
		jz[i][j]=0;
	   }
    }
                       ///////矩阵相乘
  for(i=0;i<R_d3;i++)
   for(j=0;j<1;j++)
	{
     for(k=0;k<N_d3;k++)
      jz[i][j]+=mulmatrix[i][k]*v_x[k][j];
	}                
	    
  for(i=0;i<R_d3;i++)
	for(j=0;j<1;j++)
      det_x_d3[i]=-jz[i][j];    ////得到 det_x
 
  return 1; 

 }


   ////////////////////////////////////得到det_x,下一步修正权、阀值  
   ////////////////////////////////////
   ///////////////修正权阀值///////////
   ///////////////////////////////////
 JZ_quanfazhi_d3()
 {   
   int i,j,jj;
   jj=0;
    for(i=0;i<HN_d3;i++)                       
	 {
      for(j=0;j<IN_d3;j++)
		{ 
		  W_d3[i][j]=W_d3[i][j]+det_x_d3[jj];    ///W[HN][IN]; 输入层至隐层权值
	      jj=jj+1;
		}
      sita_d3[i]=sita_d3[i]+det_x_d3[jj];		//sita[HN]; 隐层的阈值
	  jj=jj+1;
     }
    for(j=0;j<ON_d3;j++)
	{
	 for(i=0;i<HN_d3;i++)
	  {
		V_d3[j][i]=V_d3[j][i]+det_x_d3[jj];     //V[ON][HN]; 隐层至输出层权值
        jj=jj+1;
	  }
      gama_d3[j]=gama_d3[j]+det_x_d3[jj];       //gama[ON]; 输出层的阈值               
	  jj=jj+1;
    } 
  return 1;
 }



double NW_d3[HN_d3][IN_d3]; //输入层至隐层权值
double NV_d3[ON_d3][HN_d3]; //隐层至输出层权值
double Nsita_d3[HN_d3]; //隐层的阈值
double Ngama_d3[ON_d3]; //输出层的阈值
void savequan_d3()  /////////////////保存权值和阈值
{
	
	for(int i=0;i<HN_d3;i++)
	{
		for(int j=0;j<IN_d3;j++)
		{
			NW_d3[i][j]=W_d3[i][j];
		}
	}
	for(int ii=0;ii<ON_d3;ii++)
	{
		for(int jj=0;jj<HN_d3;jj++)
		{
			NV_d3[ii][jj]=V_d3[ii][jj];
		}
	}
	
	for(int k=0;k<ON_d3;k++)
	{
		Ngama_d3[k]=gama_d3[k];  
	}
	
	for(int kk=0;kk<HN_d3;kk++)
	{
		Nsita_d3[kk]=sita_d3[kk];  
	}
	
}
 

void huanyuanquan_d3()  /////////////////还原权值和阈值
{
	
	for(int i=0;i<HN_d3;i++)
	{
		for(int j=0;j<IN_d3;j++)
		{
			W_d3[i][j]=NW_d3[i][j];
		}
	}
	for(int ii=0;ii<ON_d3;ii++)
	{
		for(int jj=0;jj<HN_d3;jj++)
		{
			V_d3[ii][jj]=NV_d3[ii][jj];
		}
	}
	
	for(int k=0;k<ON_d3;k++)
	{
		gama_d3[k]=Ngama_d3[k];  
	}
	
	for(int kk=0;kk<HN_d3;kk++)
	{
		sita_d3[kk]=Nsita_d3[kk];  
	}
	
}






////////////////////////////
////////d3神经网络预测//////
////////////////////////////
void netd3_1()
{
FILE *fp,*in;
double Pre_error;  //预定误差
double sum_err,sum_err1;   //定义误差平方和
int study,m,k,i,j;
niu_d3=0.01;
if ((fp=fopen("d3bp.txt","w"))==NULL)
{
printf("不能创建d3bp.txt文件!\n");
exit(1);
}
if ((in=fopen("d3in.txt","a+"))==NULL)
{
printf("不能打开d3in.txt文件!\n");
exit(1);
}

fscanf(in,"%lf",&Pre_error);
study=0;
for (m=0;m<N_d3;m++)
{ 
   for (i=0;i<IN_d3;i++)
	   fscanf(in,"%lf",&Study_Data_d3[m].input[i]);
}
 
for (m=0;m<N_d3;m++)
{
	for (k=0;k<ON_d3;k++)
	{
	  fscanf(in,"%lf",&Study_Data_d3[m].teach[k]);     
	} 
}

initial_d3(); //隐层、输出层权、阈值初始化 

for (m=0;m<N_d3;m++) 
 {
  input_P_d3(m); //输入第m个学习样本 
  input_T_d3(m);//输入第m个样本的教师信号 
  H_I_O_d3(); //第m个学习样本隐层各单元输入、输出值 
  O_I_O_d3(); //第m个学习样本输出层各单元输入、输出值   
  Err_yangben_d3(m);    ///第m个样本的误差平方和                
  Err_O_H_s_d3();    ///求敏感度s2
  Err_H_I_s_d3();     ///求敏感度s1
  qiuyacobi_d3(m);     ////求雅可比矩阵                    (2)
 }
  sum_err1=Err_Sum_d3();       //////所有样本误差平方和         (1)	   
 
  FILE *fp1;
  fp1=fopen("yacobi_d3.txt","w");
  for (i=0;i<N_d3;i++)
   { 
 	 for (j=0;j<R_d3;j++)
     fprintf(fp1,"%lf\t",yacobi_d3[i][j]);
     fprintf(fp1,"\n");
   }
  fclose(fp1);


  qiu_detx_d3();              //求det_x            (3) 

 JZ_quanfazhi_d3();  ///修正权、阀值          (4)  

 /////////////////////////////////////////////
 ////////////求相应的网络输出及相应误差/////// (4)
 /////////////////////////////////////////////

 for (m=0;m<N_d3;m++) 
 {  
  H_I_O_d3(); //第m个学习样本隐层各单元输入、输出值 
  O_I_O_d3(); //第m个学习样本输出层各单元输入、输出值 
  Err_yangben_d3(m);   ///第m个样本的误差平方和
 }
 sum_err=Err_Sum_d3();    //修正权、阀值后所有样本误差平方和

  
 while (sum_err>Pre_error)
 {   
     if(sum_err<sum_err1)
	   {         
	     niu_d3=niu_d3/beta_d3;              ///参数修正
		 
         for (m=0;m<N_d3;m++) 
		  {  			             
             H_I_O_d3(); //第m个学习样本隐层各单元输入、输出值 
             O_I_O_d3(); //第m个学习样本输出层各单元输入、输出值 
             Err_yangben_d3(m);   ///第m个样本的误差平方和
             Err_O_H_s_d3();    ///求敏感度s2
             Err_H_I_s_d3();     ///求敏感度s1            
             qiuyacobi_d3(m);     ////求雅可比矩阵                    (2)
		  }
          sum_err1=Err_Sum_d3();          //(1)


          qiu_detx_d3();              //求det_x            (3) 
 
          JZ_quanfazhi_d3();  ///修正权、阀值          (4)  

          for (m=0;m<N_d3;m++) 
		   {  
             H_I_O_d3(); //第m个学习样本隐层各单元输入、输出值 
             O_I_O_d3(); //第m个学习样本输出层各单元输入、输出值 
             Err_yangben_d3(m);   ///第m个样本的误差平方和
		   }
          sum_err=Err_Sum_d3();    //修正权、阀值后所有样本误差平方和
	 }
    if((sum_err>sum_err1)||(sum_err==sum_err1))
       { 		
         niu_d3=niu_d3*beta_d3;                //参数修正

         qiu_detx_d3();              //求det_x            (3) 

         JZ_quanfazhi_d3();  ///修正权、阀值              (4)

         for (m=0;m<N_d3;m++) 
		  {
            H_I_O_d3(); //第m个学习样本隐层各单元输入、输出值 
            O_I_O_d3(); //第m个学习样本输出层各单元输入、输出值 
            Err_yangben_d3(m);   ///第m个样本的误差平方和
		  }
         sum_err=Err_Sum_d3();
	   }
     ++study;
	 
     fprintf(fp,"网络已经学习了%d次，现在的全局误差为：%f\n",study,sum_err); 
 }

 savequan_d3();

}


////////////////////////////
////////d3神经网络预测//////
////////////////////////////
void netd3()
{
FILE *fp,*in;
double Pre_error;  //预定误差
double sum_err,sum_err1;   //定义误差平方和
int study,m,k,i,j;
niu_d3=0.01;
if ((fp=fopen("d3bp.txt","w"))==NULL)
{
printf("不能创建d3bp.txt文件!\n");
exit(1);
}
if ((in=fopen("d3in.txt","a+"))==NULL)
{
printf("不能打开d3in.txt文件!\n");
exit(1);
}

fscanf(in,"%lf",&Pre_error);
study=0;
for (m=0;m<N_d3;m++)
{ 
   for (i=0;i<IN_d3;i++)
	   fscanf(in,"%lf",&Study_Data_d3[m].input[i]);
}
 
for (m=0;m<N_d3;m++)
{
	for (k=0;k<ON_d3;k++)
	{
	  fscanf(in,"%lf",&Study_Data_d3[m].teach[k]);     
	} 
}

huanyuanquan_d3(); //隐层、输出层权、阈值初始化 

for (m=0;m<N_d3;m++) 
 {
  input_P_d3(m); //输入第m个学习样本 
  input_T_d3(m);//输入第m个样本的教师信号 
  H_I_O_d3(); //第m个学习样本隐层各单元输入、输出值 
  O_I_O_d3(); //第m个学习样本输出层各单元输入、输出值   
  Err_yangben_d3(m);    ///第m个样本的误差平方和                
  Err_O_H_s_d3();    ///求敏感度s2
  Err_H_I_s_d3();     ///求敏感度s1
  qiuyacobi_d3(m);     ////求雅可比矩阵                    (2)
 }
  sum_err1=Err_Sum_d3();       //////所有样本误差平方和         (1)	   
 
  FILE *fp1;
  fp1=fopen("yacobi_d3.txt","w");
  for (i=0;i<N_d3;i++)
   { 
 	 for (j=0;j<R_d3;j++)
     fprintf(fp1,"%lf\t",yacobi_d3[i][j]);
     fprintf(fp1,"\n");
   }
  fclose(fp1);


  qiu_detx_d3();              //求det_x            (3) 

 JZ_quanfazhi_d3();  ///修正权、阀值          (4)  

 /////////////////////////////////////////////
 ////////////求相应的网络输出及相应误差/////// (4)
 /////////////////////////////////////////////

 for (m=0;m<N_d3;m++) 
 {  
  H_I_O_d3(); //第m个学习样本隐层各单元输入、输出值 
  O_I_O_d3(); //第m个学习样本输出层各单元输入、输出值 
  Err_yangben_d3(m);   ///第m个样本的误差平方和
 }
 sum_err=Err_Sum_d3();    //修正权、阀值后所有样本误差平方和

  
 while (sum_err>Pre_error)
 {   
     if(sum_err<sum_err1)
	   {         
	     niu_d3=niu_d3/beta_d3;              ///参数修正
		 
         for (m=0;m<N_d3;m++) 
		  {  			             
             H_I_O_d3(); //第m个学习样本隐层各单元输入、输出值 
             O_I_O_d3(); //第m个学习样本输出层各单元输入、输出值 
             Err_yangben_d3(m);   ///第m个样本的误差平方和
             Err_O_H_s_d3();    ///求敏感度s2
             Err_H_I_s_d3();     ///求敏感度s1            
             qiuyacobi_d3(m);     ////求雅可比矩阵                    (2)
		  }
          sum_err1=Err_Sum_d3();          //(1)


          qiu_detx_d3();              //求det_x            (3) 
 
          JZ_quanfazhi_d3();  ///修正权、阀值          (4)  

          for (m=0;m<N_d3;m++) 
		   {  
             H_I_O_d3(); //第m个学习样本隐层各单元输入、输出值 
             O_I_O_d3(); //第m个学习样本输出层各单元输入、输出值 
             Err_yangben_d3(m);   ///第m个样本的误差平方和
		   }
          sum_err=Err_Sum_d3();    //修正权、阀值后所有样本误差平方和
	 }
    if((sum_err>sum_err1)||(sum_err==sum_err1))
       { 		
         niu_d3=niu_d3*beta_d3;                //参数修正

         qiu_detx_d3();              //求det_x            (3) 

         JZ_quanfazhi_d3();  ///修正权、阀值              (4)

         for (m=0;m<N_d3;m++) 
		  {
            H_I_O_d3(); //第m个学习样本隐层各单元输入、输出值 
            O_I_O_d3(); //第m个学习样本输出层各单元输入、输出值 
            Err_yangben_d3(m);   ///第m个样本的误差平方和
		  }
         sum_err=Err_Sum_d3();
	   }
     ++study;
	 
     fprintf(fp,"网络已经学习了%d次，现在的全局误差为：%f\n",study,sum_err); 
 }

}



///////////////////////////////////////////
///////////////////d3预测程序//////////////
/////////////////////////////////////////// 
 double net_d3[200];
 double tmp1_d3[100],tmp2_d3[100],tmp3_d3[100],tmp4_d3[100],tmp5_d3[100];  //存学习样本
 double tmp6_d3[100],tmp7_d3[100],tmp8_d3[100],tmp9_d3[100],tmp10_d3[100];  //
 double tmp11_d3[100],tmp12_d3[100],tmp13_d3[100],tmp14_d3[100],tmp15_d3[100];//
 double tmp16_d3[100],tmp17_d3[100],tmp18_d3[100],tmp19_d3[100],tmp20_d3[100];//
 double tmp21_d3[100],tmp22_d3[100],tmp23_d3[100],tmp24_d3[100],tmp25_d3[100];//
 double tmp26_d3[100],tmp27_d3[100],tmp28_d3[100],tmp29_d3[100],tmp30_d3[100];//存教师样本 
 double tmp31_d3[100],tmp32_d3[100],tmp33_d3[100],tmp34_d3[100],tmp35_d3[100];   //存预测样本数据
 
 double d3_temp[146];          
 double temp_d3,max_d3,min_d3;
 

 void yuce_d3_chushi()
 {
   
   int i;
   CString pday;
   int day;
   CString strday;
   int day1,day2,day3,day4,day5; 
   int yucetian;
  // cout<<"预测的日期(天)为:"<<endl;
   //cin>>yucetian;
   yucetian=30;  /////预测第30天的负荷序列  

  
  day1=yucetian-14;  ////////训练样本一预测天为29   
  day=day1;  //29
  pday="d3_day";
  strday.Format("%d",day);
  pday+=strday;