spsofc.cpp

应用随机微粒群算法求解一个六元的非线性方程组.

#include "stdio.h"
#include "stdlib.h"
#include "time.h"
#include "math.h"

//rand number defination
#define rdint(i)(rand()%(int)(i))       //rand integer
#define rdft()(float)((double)rdint(16384)/(16383.0)) //rand floating-point number
#define rnd(a,b)(rdint((int)(b)-(int(a)+1)+(int)(a))  //rand integer between a and b

//macro defination
#define POPSIZE 20
#define DIMENSION 6
#define POP 21
//This a program which uses stochatic particle swarm optimization to solve the nonlinear equations.
//global vector defination
float W=0.0;
float C1=1.8;
float C2=1.8;
float VMAX=2.0;
float XMIN=-2.0;
float XMAX=2.0;
float P[DIMENSION];
double PBEST;

struct indi                   //define individual structure
{
	float number[DIMENSION];
	float best[DIMENSION];
	double bestfitness;
	double fitness;
	float speed[DIMENSION];
}individual[POP];

void initiate(void);
void calculation(int number);
void globalbest(int number);
void localbest(int number);
void generate(int number);

void generate(int num)
{
   int j;
   for(j=0;j<DIMENSION;j++)
	  individual[num].number[j]=rdft()*(XMAX-XMIN)+XMIN;
   calculation(num);
   individual[num].bestfitness=individual[num].fitness;
   for(j=0;j<DIMENSION;j++)
	  individual[num].best[j]=individual[num].number[j];
}

//initialize 
void initiate()
{
	int i,j;
//	srand((int)time(0));

	for(i=0;i<POPSIZE;i++)
		for(j=0;j<DIMENSION;j++)
			individual[i].number[j]=rdft()*(XMAX-XMIN)+XMIN; //rand value of individual position
	for(i=0;i<POPSIZE;i++)
		for(j=0;j<DIMENSION;j++)
			individual[i].speed[j]=VMAX*rdft();   //rand value of individual speed 
	for(i=0;i<POPSIZE;i++)
		for(j=0;j<DIMENSION;j++)
			individual[i].best[j]=individual[i].number[j]; //the best individual position 
	for(i=0;i<POPSIZE;i++)
		calculation(i);
    for(i=0;i<POPSIZE;i++)
		individual[i].bestfitness=individual[i].fitness; //bestfitness of individual position
	globalbest(0);
}

//history optimal position modification of particle
void localbest(int number)
{
	int i;
	if(individual[number].bestfitness>individual[number].fitness)
		for(i=0;i<DIMENSION;i++)
			individual[number].best[i]=individual[number].number[i];
		individual[number].bestfitness=individual[number].fitness;
}

//history optimal position modification of swarm
void globalbest(int number)
{
	int i,j;
	float s=0;
	int flag=0;
	if(number==0)
	{
		s=individual[0].fitness;
		flag=0;
		for(i=1;i<POPSIZE;i++)
			if(individual[i].fitness<s)
			{
				s=individual[i].fitness;
				flag=i;
			}
        for(i=0;i<DIMENSION;i++)
            P[i]=individual[flag].number[i];
		    PBEST=individual[flag].fitness;
	}
    else
	{
		for(i=0;i<POPSIZE;i++)
			if(individual[i].bestfitness<PBEST)
			{
				for(j=0;j<DIMENSION;j++)
				    P[j]=individual[i].best[j];
				PBEST=individual[i].bestfitness;
			}
	}
}

//fitness function   (including the nonlinear equations)
void calculation(int num)
{
  double f1,f2,f3,f4,f5,f6;

  f1=individual[num].number[0]+individual[num].number[1]*individual[num].number[1]
		 *individual[num].number[3]*individual[num].number[5]/4.0+0.75;
  f2=individual[num].number[1]+0.405*exp(1.0+individual[num].number[0]*individual[num].number[1])
		-1.405;
  f3=individual[num].number[2]-individual[num].number[1]*individual[num].number[5]/4.0+1.5;

  f4=individual[num].number[3]-0.605*exp(1-individual[num].number[2]*individual[num].number[2])
		 -0.395;
  f5=individual[num].number[4]-individual[num].number[1]*individual[num].number[5]/2.0
		 +1.5;
  f6=individual[num].number[5]-individual[num].number[0]*individual[num].number[4];

  individual[num].fitness=fabs(f1)+fabs(f2)+fabs(f3)+fabs(f4)+fabs(f5)+fabs(f6);
}

void main()
{
	int i,j,k,t,m,flag=-1,fflag;
	srand(time(NULL));
    clock_t start, finish; 
	double duration;          //test time duration of spso
	float EE[31],EEcs[31];

	for(j=0;j<50;j++)
	{
		initiate();
		start = clock();
		for(i=0;i<25000;i++)
		{
		   flag=-1;
		   for(k=0;k<POPSIZE;k++)
		   {
			  if((individual[k].fitness==individual[k].bestfitness)&&(individual[k].fitness==PBEST))
			  {
				 flag=k;
				 break;
			  }
		   }
		   if(flag!=-1)
			  generate(flag);
		   for(k=0;k<POPSIZE;k++)
		   if(k!=flag)
		   {
			  for(t=0;t<DIMENSION;t++)
			  {
		   individual[k].speed[t]=W*individual[k].speed[t]+C1*rdft()*(individual[k].best[t]
			-individual[k].number[t])+C2*rdft()*(P[t]-individual[k].number[t]);
				 if(individual[k].speed[t]>VMAX)
					individual[k].speed[t]=VMAX;
			individual[k].number[t]=individual[k].number[t]+individual[k].speed[t];
				 if(individual[k].number[t]<XMIN)
					   individual[k].number[t]=XMIN;
				 if(individual[k].number[t]>XMAX)
					   individual[k].number[t]=XMAX;
			  }
			  calculation(k);
			  localbest(k);
		   }
		   globalbest(1);

		   for(m=0;m<DIMENSION;m++)
			  individual[POPSIZE].number[m]=P[m];

		}
		printf("%d %f\t",i,PBEST);
		for (j=0;j<DIMENSION;j++)
		   printf("%f\t",individual[POPSIZE].number[j]);

		finish = clock();
		duration = (double)(finish-start)/CLOCKS_PER_SEC;

		printf( "The time of one experiment is %f seconds\n", duration);
	}

}