pso_binary.c

采用二进制对PSO进行编码

 cp1=cp0/z; // Decreasing. Equal to initial cp0 just after initialisation

 cp2=cp1; // Decreasing
 if(K>2)
  cp3=cp2/log(K-1);// Decreasing. Slightly smaller than cp2
 else cp3=cp2;

 

   n1=(int)(cp1*D);   // Decreasing.
   n2=(int)(cp2*D);// Decreasing. The more you trust p, the less you update it
   n3=(int)(cp3*D); //Decreasing.  The more you trust g, the less you update it

  Pp=P_m[s];Gp=P_m[g];

  for (d=0;d<D;d++)   Vp.x.b[d]=-1;  // Arbitrary value: means "non assigned"

  for (k=0;k<n1;k++) // Velocity. At most n1 non assigned components
  {
    d=alea_integer(0,D-1); Vp.x.b[d]=alea_integer(0,1);
  }
  for (k=0;k<n2;k++)  // Around p. At most n2 modified components
  {
    d=alea_integer(0,D-1); Pp.x.b[d]=1- Pp.x.b[d];
  }
    for (k=0;k<n3;k++) // Around g. At most n3 modified components
  {
    d=alea_integer(0,D-1); Gp.x.b[d]=1- Gp.x.b[d];
  } 
 
	for (d=0;d<D;d++) // Choose each component
	{
		bit1=Vp.x.b[d];  bit2=Pp.x.b[d];   bit3=Gp.x.b[d];

		//Majority choice
		if(bit1>=0) // If assigned velocity
		{ bit=bit1+bit2+bit3;  
		if(bit<2)  P[s].x.b[d]=0; else P[s].x.b[d]=1;
		}
		else // If non assigned velocity
		{
			 bit=bit2+bit3;
			if(bit==2) P[s].x.b[d]=1; else
			{if(bit==0) P[s].x.b[d]=0; else P[s].x.b[d]=alea_integer(0,1);}
		}
	}
break;


case 6:
// Method 6---------------------------------
  //cp2: Confidence coefficient for p
  //cp3: Confidence coefficient for g

cp0=1-0.5*log((double)D)/log((double)S);
 cp0=cp0+((double)K/(double)S)*log((double)D)/log((double)S);
if (K>1) cp0=cp0*log((double)K);

 z=(double)nb_eval*K/(double)S;
 if (z>1) z=log(z);

 cp1=cp0/z;
 cp2=cp1;


 if(K>2)
  cp3=cp1/log(K-1);// Decreasing. 
 else cp3=cp1;

   n2=(int)cp2*D;// Decreasing. The more you trust p, the less you update it
   n3=(int)cp3*D; //Decreasing.  The more you trust g, the less you update it

  Pp=P_m[s];Gp=P_m[g];

  for (k=0;k<n2;k++)  // Around p. At most n2 modified components
  {
    d=alea_integer(0,D-1); Pp.x.b[d]=1- Pp.x.b[d];
  }
    for (k=0;k<n3;k++) // Around g. At most n3 modified components
  {
    d=alea_integer(0,D-1); Gp.x.b[d]=1- Gp.x.b[d];
  } 
  
	for (d=0;d<D;d++) // Choose each component
	{
		 bit2=Pp.x.b[d];   bit3=Gp.x.b[d];


		//Majority choice
	//	if(bit2==bit3 && alea(0,1)<cp1) P[s].x.b[d]=bit2;
    if(bit2==bit3) P[s].x.b[d]=bit2;  
		else P[s].x.b[d]=alea_integer(0,1);
	}
break;

case 7:
// Method 7---------------------------------
// Like 5, but coefficients are modified only if there has been


// some improvement

  //cp2: Confidence coefficient for p
  //cp3: Confidence coefficient for g

if (error<error_prev) goto skip7;

	z7=z7+(double)K/(double)S;
	z=log(z7);
	cp1=cp0/z;

	cp2=cp1;

	if(K>2)

	cp3=cp1/log(K-1);// Decreasing. 
	else cp3=cp1;

   n2=(int)cp2*D;// Decreasing. The more you trust p, the less you update it
   n3=(int)cp3*D; //Decreasing.  The more you trust g, the less you update it
skip7:

  Pp=P_m[s];Gp=P_m[g];

  for (k=0;k<n2;k++)  // Around p. At most n2 modified components
  {
    d=alea_integer(0,D-1); Pp.x.b[d]=1- Pp.x.b[d];
  }
    for (k=0;k<n3;k++) // Around g. At most n3 modified components
  {
    d=alea_integer(0,D-1); Gp.x.b[d]=1- Gp.x.b[d];
  } 
  
	for (d=0;d<D;d++) // Choose each component
	{
		 bit2=Pp.x.b[d];   bit3=Gp.x.b[d];

		//Majority choice
		if(bit2==bit3) P[s].x.b[d]=bit2;	
		else P[s].x.b[d]=alea_integer(0,1);
	}
break;

case 71:
// Like 7, but with probabilistic majority choice
if (error<error_prev) goto skip71;

	z7=z7+(double)K/(double)S;
	z=log(z7);
	cp1=cp0/z;
	cp2=cp1;

	if(K>2)

	cp3=cp1/log(K-1);// Decreasing. 
	else cp3=cp1;

   n2=(int)cp2*D;// Decreasing. The more you trust p, the less you update it
   n3=(int)cp3*D; //Decreasing.  The more you trust g, the less you update it
skip71:

  Pp=P_m[s];Gp=P_m[g];

  for (k=0;k<n2;k++)  // Around p. At most n2 modified components
  {
    d=alea_integer(0,D-1); Pp.x.b[d]=1- Pp.x.b[d];
  }
    for (k=0;k<n3;k++) // Around g. At most n3 modified components
  {

    d=alea_integer(0,D-1); Gp.x.b[d]=1- Gp.x.b[d];
  }
  
	for (d=0;d<D;d++) // Choose each component
	{
		 bit2=Pp.x.b[d];   bit3=Gp.x.b[d];
		bit=bit2+bit3;
		//Probabilistic majority choice
		P[s].x.b[d]=alea_integer(0,1);
		pr=alea(0,1);
		// You may have 00 just by chance.
		if (bit==0 && pr>0.025*cp2*cp3) {P[s].x.b[d]=0;continue;}
		// You may have 11 just by chance
		if (bit==2 && pr>0.025*cp2*cp3) {P[s].x.b[d]=1;continue;}

//continue;
		// method 71a: to be more rigourous, you may also take into account
		// 01 and 01: 

		if (bit2==0 && bit3==1 && pr<0.5*cp2*(1-cp3)) {P[s].x.b[d]=1;continue;}
		if (bit2==0 && bit3==1 && pr<0.5*(1-cp2)*cp3) {P[s].x.b[d]=0;continue;}
		if (bit2==1 && bit3==0 && pr<0.5*(1-cp2)*cp3) {P[s].x.b[d]=1;continue;}
		if (bit2==1 && bit3==0 && pr<0.5*cp2*(1-cp3)) {P[s].x.b[d]=0;continue;}
		
	}
break;

case 72: // Mixing 7 and 11. Constant coefficients

    dist=(int)log(D);
	  
	Gp=P_m[g];
	r=alea(1,dist);

    for (k=0;k<r;k++) // Around g.
  {

    d=alea_integer(0,D-1); Gp.x.b[d]=1- Gp.x.b[d]; 
  } 
  
	if(s==g) { P[s]=Gp; break;}
 
	Pp=P_m[s];	
	r=alea(1,dist);
  for (k=0;k<r;k++)  // Around p.
  {
    d=alea_integer(0,D-1); Pp.x.b[d]=1- Pp.x.b[d];
  }

	for (d=0;d<D;d++) // Choose each component
	{
		 bit2=Pp.x.b[d];   bit3=Gp.x.b[d];

		//Majority choice
		if(bit2==bit3) P[s].x.b[d]=bit2;	
		else 
		{	
			P[s].x.b[d]=alea_integer(0,1);
	/*		
		// "a" variant. Trust g a bit more
			pr=alea(0,1+1/log(K-1)); 
			if (pr<=1) P[s].x.b[d]=bit2; else P[s].x.b[d]=bit3;
	*/	
		}
	}
break;

case 73:
if (error<error_prev) goto skip73;

	z7=z7+(double)K/(double)S;
	z=log(z7);
	cp1=cp0/z;
	cp2=cp1;

	if(K>2)

	cp3=cp1/log(K-1);// Decreasing. 
	else cp3=cp1;

skip73:

  Pp=P_m[s];Gp=P_m[g];
	dist=(int)cp2*D;
	r=alea(1,dist);
  for (k=0;k<r;k++)  // Around p.
  {
    d=alea_integer(0,D-1); Pp.x.b[d]=1- Pp.x.b[d];
  }
	dist=(int)cp3*D;
	r=alea(1,dist);
    for (k=0;k<r;k++) // Around g. 
  {
    d=alea_integer(0,D-1); Gp.x.b[d]=1- Gp.x.b[d];
  } 
  
	for (d=0;d<D;d++) // Choose each component
	{
		 bit2=Pp.x.b[d];   bit3=Gp.x.b[d];

		//Majority choice
			if(bit2==bit3) P[s].x.b[d]=bit2;	
		else   // P[s].x.b[d]=alea_integer(0,1);
		{	
			pr=alea(0,1+log(K-1));
			if (pr<=1) P[s].x.b[d]=bit2; else P[s].x.b[d]=bit3;
		}
		
	}
break;

case 8:
// Method 8---------------------------------

  Pp=P_m[s];Gp=P_m[g];


	for (d=0;d<D;d++) // Choose each component
	{
     bit2=Pp.x.b[d];   bit3=Gp.x.b[d];

    //Majority choice
		if(bit2==bit3) P[s].x.b[d]=bit2; else P[s].x.b[d]=alea_integer(0,1);
	}
break;

case 9:
// Method 9 ------------------------------------------------

//pr=1.0/pow(2,D);

c2=(double)nb_eval/(double)S;
c3=K*c2;
//cp2=pow(1-pr,c2);
//cp3=pow(1-pr,c3);
cp2=c2*log(1-pr);  // log(Failure proba)
cp3=c3*log(1-pr);

//printf("\n %f %f ",cp2,cp3);

if(s!=g)
{
	for (d=0;d<D;d++) // Choose each component
	{
		 bit2=P_m[s].x.b[d];   bit3=P_m[g].x.b[d];

     if (log(alea(0,1))>cp2) bit2=1-bit2;

     if(log(alea(0,1))>cp3) bit3=1-bit3;
		//Majority choice
    if(bit2==bit3) P[s].x.b[d]=bit2;
		else P[s].x.b[d]=alea_integer(0,1);
	}
}
else
{	
	for (d=0;d<D;d++) // Choose each component
	{
		bit3=P_m[g].x.b[d];
		if (alea(0,1)>cp3) bit3=1-bit3;
		P[s].x.b[d]=bit3;
	}
}
break;

case 10:
// Method 10 ------------------------------------------ 
cp0=0.5;
 cp1=cp0;
  cp2=(double)nb_eval*K/(double)S;
  cp2=1-cp1/log(cp2);
  cp3=cp2;

	for (d=0;d<D;d++)
	{
    if (deter==-1)
    {
      cp1=alea(cp0/2,cp0);
      cp2=(double)nb_eval*K/(double)S;
    cp2=1-cp1/log(cp2);
    cp3=cp2;
    }

    bit1=P[s].x.b[d];  bit2=P_m[s].x.b[d];   bit3=P_m[g].x.b[d];
    pr=alea(0,1); if (pr<cp1) choice[0]= bit1; else choice[0]=1-bit1;
    pr=alea(0,1); if (pr<cp2) choice[1]=bit2; else choice[1]=1-bit2;
    pr=alea(0,1); if (pr<cp3) choice[2]=bit3; else choice[2]=1-bit3;

    //Majority
    bit=choice[0]+choice[1]+choice[2];
    if(bit<2)  P[s].x.b[d]=0; else P[s].x.b[d]=1;
    V[s].b[d]=  P[s].x.b[d]-bit1; // -1 or 0 or 1. For information
	}

  break;

case 11:  // Pivot method, from an idea of P. Serra, 1997
          // (seriously) modified and adapted to binary problems by M. Clerc, 2004
case11:
	// Works pretty well on some problems .. and pretty bad on some others
  // Note that there is no velocity, and the current position is not used
  // It works better when information links are modified at random
  // if there has been no improvement after the previous iteration
  
 P[s]=P_m[g]; // Best previous position g of the best informant

 // (for D>=2) Simplified formula.  Note that "dist" is an integer
 //dist=log(D);  if (dist<3) dist=3;  // 11
 //dist=D*log(D)/S;// 11e
 //dist=1+(K-1)*Diam/(double)S; //if (dist<2)dist=2;   // 11f
 //printf("\n %i %i",nb_pivots, Diam);
 dist=(int)(1+Diam/(double)nb_pivots); if (dist<3) dist=3; //* 11g

//if(dist>3)printf("\n %i",dist); 
 if (dist>D) dist=D;  // For D=1.

  r=alea(1,dist); //11. Note that r is never equal to dist
 // r=1+alea_integer(1,dist-1);  //11b. Theoretically equivalent
 // r=1+alea_integer(1,dist); //11c. A bit larger
 // r=alea_integer(1,dist-1);  // 11d.   More uniform , on dist k values
                              // Sometimes far better but also sometimes far worse
  
 for (d=0;d<D;d++) distrib[d]=0; // Just for some stats
 
 for (k=0;k<r;k++) // Define a position around g
 //    Note that even if there at least two k values
 //  there is a small probability (1/D^k)
 // that just one bit is modified
  {
    d=alea_integer(0,D-1);
    P[s].x.b[d]=1- P[s].x.b[d];
  distrib[d]=1;
  }
   i=0;    for (d=0;d<D;d++) i=i+distrib[d];
   k_tot[i]=k_tot[i]+1; // Just for some statistics
 
 break;

case 12:  // Pivot method, variant "independent particles"
		// K is not taken into account
 Gp=P_m[s];

 dist=(int)log(D); 

  r=alea(1,dist);
 //r=1+alea_normal(0,dist);// Test

 for (k=0;k<r;k++)
  {
    d=alea_integer(0,D-1);
    Gp.x.b[d]=1- Gp.x.b[d]; // Around s
  }
    P[s]=Gp;
 break;


 case 13:  // Pivot method, variant "one random informer"
		// K is not taken into account
d=alea_integer(0,D-1);
if (total_error(P_m[d].f)<total_error(P_m[s].f)) g=d; else g=s;
 Gp=P_m[g];

  dist=(int)log(D); // c

  r=alea(1,dist);
 //r=1+alea_normal(0,dist);// Test

 for (k=0;k<r;k++)
  {
    d=alea_integer(0,D-1);
    Gp.x.b[d]=1- Gp.x.b[d]; // Around s
  }
    P[s]=Gp;
 break;


case 14: 

case14:
/*
 for (i=0;i<S;i++)  // For each particle ...
{  
    g_[i]=best_info(i,option_best); // .. find the best informer
} 
	
// Find the nearest one
 z=D+1;
 for (i=0;i<S;i++)
 {
	m=g_[i];
	if(m==g) continue;
		r=distance(g,m,-2);
		if(r<=0) continue;
			if(r<z) z=r;
 }
*/
/*

 z=D+1; // Find the nearest informer
 for (m=0;m<S;m++)
 {
	if(m==g) continue;
		r=distance(g,m,-2);
		if(r<=0) continue;
			if(r<z) z=r;
 }

*/

//----------- 
	// Good for M黨lenbein, Multimodal
	// Bad for Zebra3
 Gp=P_m[g];

  // A "center" between s and g
 for (d=0;d<D;d++) 
 {
	 if (P_m[s].x.b[d]!=P_m[g].x.b[d] && alea(0,1)<0.5) Gp.x.b[d]=P_m[s].x.b[d];
 }
z=distance(s,g,-1)+1;
r=alea(1,z);

 for (k=0;k<r;k++)
  {

    d=alea_integer(0,D-1);
    Gp.x.b[d]=1- Gp.x.b[d]; // Around g
  }
    P[s]=Gp;
 break;

case 15:
 if (error<error_prev) goto case14; else goto case11;


case 16: // Like 11, but with "taboo"
Gp=P_m[g];

  dist=(int)log(D);  // We suppose here D>=2
  r=alea(1,dist);

 for (k=0;k<r;k++)
  {
    d=alea_integer(0,D-1);
    if(B_mod[s].b[d]<1) // "Taboo". Modify only if it has not just been modified
      {
        Gp.x.b[d]=1- Gp.x.b[d]; // Around g
      }
  }
 
    for (d=0;d<D;d++)

    {
      if (Gp.x.b[d]!=P[s].x.b[d])  B_mod[s].b[d]=1; // Has been modified
          else B_mod[s].b[d]=0; // Has not been modified
       P[s].x.b[d] =Gp.x.b[d];
      }