nn-2.cpp

三个输出的神经元网络拟合如下函数（前面上传的文件搞错了）。 y1=x1ln(x2) + x2ln(x3)；y2=x3ln(x4) + x4ln(x5)；y3=x5ln(x6) + x6ln(x1)

// Neural Networks
// Written by Microsoft Visual C++
// Copyright by UTLab at Tsinghua University

#include "math.h"
#include "stdlib.h"
#include "stdio.h"
#include "conio.h"
#include "UTLab.h"

#define N 6     // number of input neurons (decision variables)
#define H 15     // number of hidden neurons
#define O 3     // number of output neurons
#define D 3000  // number of training data
#define SIMU_PRINT_NUMBER 100

static double BOUND[N+O+1][2], WH[H+1][N+1], WO[O+1][H+1];
static void   Simu(double x[D+1][N+1], double y[D+1][O+1]);
static void   Input_Output(double x[N+1], double y[O+1]);
static void   Train(double original_x[D+1][N+1], double original_y[D+1][O+1],
			  double standard_x[D+1][N+1], double standard_y[D+1][O+1]);
static void   Standarlize(double original_x[D+1][N+1], double original_y[D+1][O+1],
			  double standard_x[D+1][N+1], double standard_y[D+1][O+1]);
static void   BPA(double standard_x[D+1][N+1], double standard_y[D+1][O+1]);
static void   Init_NN(void);
static void   NN(double x[N+1], double y[O+1]);

static void Input_Output(double x[N+1], double y[O+1])
{
   int i;
   
   for(i=1; i<=N; i++) x[i]=myu(1,5);

   y[1]=x[1]*log(x[2])+x[2]*log(x[3]);
   y[2]=x[3]*log(x[4])+x[4]*log(x[5]);
   y[3]=x[5]*log(x[6])+x[6]*log(x[1]);
}


int main( void )
{
   double x[D+1][N+1], y[D+1][O+1], standard_x[D+1][N+1], standard_y[D+1][O+1];

   srand(100);
   Simu(x,y);
   Standarlize(x,y,standard_x,standard_y);
   Train(x,y,standard_x,standard_y);

  
   return 1;
}

static void BPA(double x[D+1][N+1],double y[D+1][O+1])
{	
	
	double E0 = 0.001;
    double EE0 = 0.001;
	double TOLERANCE=0.0001;// *((N+1+O)*H+O);
	int    i,j,k,step=0;
	double ALPHA=0.05;
	double BETA=4./3;
	double ETA=0.5;
	double MU=1.;
	double DETA_WH[H+1][N+1],DETA_WO[O+1][H+1];
	double E_WO[O+1],E_WH[H+1],h[H+1],NN_y[O+1];
	double e[O+1],lambda,E,EE,TD,temp;
	
	Init_NN();
    for( i=1; i<=H; i++ )
		for( j=0; j<=N; j++ )
			DETA_WH[i][j] = 0.;
	for( i=1; i<=O; i++ )
		for( j=0; j<=H; j++ )
			DETA_WO[i][j] = 0.;

	lambda=1.;
	do
	{	step++; 
	    E=0.;
		TD=0.;
		EE=0.;
		for( k=1; k<=D; k++)
		{	/*calculate the output of the hidden neurons*/
			h[0]=1.;
			for( i=1; i<=H; i++ )
			{	h[i] = WH[i][0];
				for( j=1; j<=N; j++ )
					h[i] += WH[i][j] * x[k][j];
				h[i] = tanh(h[i]);
			}

			/*calculate the output of the outlayer neurons*/
			for( i=1; i<=O; i++ )
			{	NN_y[i] = WO[i][0];
				for( j=1; j<=H; j++)
					NN_y[i] += WO[i][j] * h[j];
			}

			/*adjust the weights of WO:*/
			for( i=1; i<=O; i++ )
			{	e[i] = y[k][i] - NN_y[i];
				E_WO[i] = lambda*e[i] + (1-lambda)*tanh(BETA*e[i]);
				for( j=0; j<=H; j++ )
				{	DETA_WO[i][j] = ETA*DETA_WO[i][j] + ALPHA*E_WO[i]*h[j];
					TD += fabs(DETA_WO[i][j]);
					WO[i][j] += DETA_WO[i][j];
				}
			}

			/*adjust the weights of WH:*/
			for( j=1; j<=H; j++ )
			{	for( i=1,temp=0.; i<=O; i++ )
					temp += E_WO[i] * WO[i][j];
				E_WH[j] = ( 1-h[j]*h[j])*temp;
				DETA_WH[j][0] = ETA*DETA_WH[j][0] +ALPHA*E_WH[j];
				TD += fabs(DETA_WH[j][0]);
				WH[j][0] += DETA_WH[j][0];
				for( i=1; i<=N; i++ )
				{   DETA_WH[j][i]= ETA*DETA_WH[j][i] +ALPHA*E_WH[j]*x[k][i];
					TD += fabs(DETA_WH[j][i]);
					WH[j][i] += DETA_WH[j][i];
				}
			}

			/*calculate the total errors*/
			for( i=1; i<=H; i++ )
			{	h[i] = WH[i][0];
				for( j=1; j<=N; j++ )
					h[i] += WH[i][j]*x[k][j];
				h[i] = tanh(h[i]);
			}
			for( i=1; i<=O; i++ )
			{	NN_y[i] = WO[i][0];
				for( j=1; j<=H; j++)
					NN_y[i] += WO[i][j] * h[j];
			}
			for( i=1; i<=O; i++ )
			{	e[i] = y[k][i] - NN_y[i];
				EE += fabs(e[i])/D;
				E += 0.5*e[i]*e[i];
			}
		}
		lambda = exp( -MU/(E*E));
		TD=TD/D;
		if(step%100==0) printf("\nBPA_%d  (Std)Sum-squared-E=%f  Average-E=%f  Gradient=%f;",step,E,EE,TD);
	} while ( step<200 && TD>TOLERANCE && E>E0 && EE>EE0);
	printf("\nBPA_%d  (Std)Sum-squared-E=%f  Average-E=%f  Gradient=%f;",step,E,EE,TD);
}

static void Simu(double x[D+1][N+1], double y[D+1][O+1])
{
   int	i, j;
   double xx[N+1], yy[O+1];
   float read;
   FILE  *fp;

   if((fp=fopen("TRAIN2.dat","r"))==NULL) {
     mark:
     fp=fopen("TRAIN2.dat","w");
     for(i=1; i<=D; i++) {
       Input_Output(xx, yy);
       if(i%SIMU_PRINT_NUMBER==0) 
		  printf("\n%d input-output Data have been generated",i);
       for(j=1; j<=N; j++) {
          x[i][j]=xx[j];
          fprintf(fp,"%3.6f ",xx[j]);
       }
       for(j=1; j<=O; j++) {
          y[i][j]=yy[j];
          fprintf(fp,"%3.6f ",yy[j]);
       }
	   fprintf(fp,"\n");
     }
     fclose(fp);
   }
   else {
     printf("The training data <<TRAIN2.dat>> already exists. Renew it? (y/n)");
     if(getche()=='y') goto mark;
     for(i=1; i<=D; i++) {
        for(j=1; j<=N; j++) {
          fscanf(fp,"%f",&read);
          x[i][j]=read;
        }
        for(j=1; j<=O; j++) {
          fscanf(fp,"%f",&read);
          y[i][j]=read;
        } 
     }
   }
   fclose(fp);
}


static void Train(double x[D+1][N+1], double y[D+1][O+1],double standard_x[D+1][N+1], double standard_y[D+1][O+1])
{
   int i, j, k;
   double E1, E2, a, xx[N+1], yy[O+1];
   float read;
   FILE *fp;

   if((fp=fopen("WEIGHT2.dat","r"))==NULL) {
        mark:
        BPA(standard_x,standard_y);
        E1=0;
        E2=0;
        for(i=1; i<=D; i++) {
           for(j=1; j<=N; j++) xx[j]=x[i][j];
           NN(xx,yy);
           a=0;
           for(j=1; j<=O; j++) a += (yy[j]-y[i][j])*(yy[j]-y[i][j]);
           E1=E1+a;
           E2=E2+sqrt(a);
        }
        E1=E1/2;
        E2=E2/D;
        printf("\n\nSum-Squared Error = %3.4f",E1);
        printf("\nAverage Error = %3.4f\n\n",E2);
        for(i=1; i<=18; i++){
           k=(int)myu(1,D);
           for(j=1; j<=N; j++) xx[j]=x[k][j];
           NN(xx,yy);
           printf("Sample %2d: y=(",i);
		   for(j=1; j<=O; j++) {
			   if(j<O) printf("%3.4f,",y[k][j]);
			   else printf("%3.4f",y[k][j]);
		   }
		   printf(") Dif.=(");
           for(j=1; j<=O; j++) {
			   if(j<O) printf("%3.4f,",yy[j]-y[k][j]);
			   else printf("%3.4f",yy[j]-y[k][j]);
		   }
           printf(")\n");
        }
        printf("\nAre you satisfactory?(y/n)");
        if(getche()=='n') goto mark;
		printf("\n");

        fp=fopen("WEIGHT2.dat","w");
        for(i=1; i<=H; i++) {
           for(j=0; j<=N; j++) fprintf(fp,"%3.4f ",WH[i][j]);
           fprintf(fp,"\n");
        }
        for(i=1; i<=O; i++) {
          for(j=0; j<=H; j++) fprintf(fp,"%3.4f ",WO[i][j]);
          fprintf(fp,"\n");
        }
        fclose(fp);

   }
   else {
     printf("\nThe weight data <<WEIGHT2.dat>> already exists. Renew it? (y/n)");
     if(getche()=='y') goto mark;
     for(i=1; i<=H; i++)
        for(j=0; j<=N; j++) {
 	     fscanf(fp, "%f", &read);
        WH[i][j]=read;
     }
	  for(i=1; i<=O; i++)
		  for(j=0; j<=H; j++) {
			  fscanf(fp,"%f", &read);
           WO[i][j]=read;
        }
     fclose(fp);
   }
}



static void Standarlize(double x[D+1][N+1],double y[D+1][O+1],double xx[D+1][N+1],double yy[D+1][O+1])
{  
	int i,k;
	for(i=1; i<=N; i++) {
		BOUND[i][0]=x[1][i];
		BOUND[i][1]=x[1][i];
	}
	for( i=1; i<=O; i++ ) {
		BOUND[N+i][0]=y[1][i];
		BOUND[N+i][1]=y[1][i];
	}
	for( k=2; k<=D; k++ ) {
		for( i=1; i<=N; i++ )
		{   if(BOUND[i][0]>x[k][i])	
				BOUND[i][0]=x[k][i];
			if(BOUND[i][1]<x[k][i])	
				BOUND[i][1]=x[k][i];
		}
		for( i=1; i<=O; i++ )
		{   if(BOUND[N+i][0]>y[k][i])	
				BOUND[N+i][0]=y[k][i];
			if(BOUND[N+i][1]<y[k][i])	
				BOUND[N+i][1]=y[k][i];
		}
	}
	for( i=1; i<=N+O; i++ )
		BOUND[i][1]=BOUND[i][1] - BOUND[i][0];
	for( k=1; k<=D; k++ )
	{	for( i=1; i<=N; i++ )
			xx[k][i]=(x[k][i]-BOUND[i][0])/BOUND[i][1];
		for( i=1; i<=O; i++ )
			yy[k][i]=(y[k][i]-BOUND[N+i][0])/BOUND[N+i][1];
	}
}

static void NN( double x[], double y[] )
{	
	double h[H+1];
	int i,j;
	for( i=1; i<=H; i++ )
	{	h[i] = WH[i][0];
		for( j=1; j<=N; j++ )
			h[i] += WH[i][j] * (x[j]-BOUND[j][0])/BOUND[j][1];
		h[i] = tanh( h[i] );
	}
	for( i=1; i<=O; i++ )
	{	y[i] = WO[i][0];
		for( j=1; j<=H; j++)
			y[i] += WO[i][j] * h[j];
		y[i] = y[i]*BOUND[N+i][1] + BOUND[N+i][0];
	}
}

static void Init_NN(void)
{
   int i, j;
   for(i=1; i<=H; i++)
      for(j=0; j<=N; j++)
         WH[i][j]= myu(-0.3,0.3);
   for(i=1; i<=O; i++)
      for(j=0; j<=H; j++)
         WO[i][j]= myu(-0.3,0.3);
}