📄 main.c
字号:
for (s = 0; s < R.SW.S; s++)
{ if(alea(0,1)<param.w0Rate) R.SW.w[s]=param.w[0]; else {R.SW.w[s]= (alea(0,1)+alea(0,1))*(param.w[1]-param.w[0])+2*param.w[0]-param.w[1]; } R.SW.c[s]=param.c1; } break; case 3: // Random w (uniform) for all particles zz=param.w[0]*0.5; zzz=zz+0.5; for (s = 0; s < R.SW.S; s++)
{ R.SW.w[s]=alea(zz,zzz); R.SW.c[s]=param.c1; //R.SW.c[s]=pow(R.SW.w[s]+1,2)*0.5; } break; case 4: // Forced "scout" particle(s) // m=sqrt(pb.SS.D); m=1; for (s = 0; s < m; s++)
{ R.SW.w[s]=alea(param.w[0],param.w[1]); R.SW.c[s]=param.c1; // R.SW.w[s]= (alea(0,1)+alea(0,1))*(param.w[1]-param.w[0])+2*param.w[0]-param.w[1]; // R.SW.c[s]=pow(R.SW.w[s]+1,2)*0.5; } for (s = m; s < R.SW.S; s++)
{ R.SW.w[s]=param.w[0]; R.SW.c[s]=param.c1; } break; }
// Take quantisation into account
R.SW.X[s] = quantis (R.SW.X[s], pb.SS);
// First evaluations
for (s = 0; s < R.SW.S; s++)
{
R.SW.X[s].f = perf (R.SW.X[s], pb);
R.SW.P[s] = R.SW.X[s]; // Best position = current one
R.SW.P[s].improved = 0; // No improvement
}
// If the number max of evaluations is smaller than
// the swarm size, just keep evalMax particles, and finish
if (R.SW.S>pb.evalMax) R.SW.S=pb.evalMax;
R.nEval = R.SW.S;
// Find the best
// Remember that f is fabs(fitness-objective)
// We want to minimise it
R.SW.best = 0;
errorPrev = fabs(pb.epsilon-R.SW.P[R.SW.best].f);
for (s = 1; s < R.SW.S; s++)
{
zz=fabs(pb.epsilon-R.SW.P[s].f);
if (zz < errorPrev)
{
R.SW.best = s;
errorPrev=zz;
}
}
/*
// Display the best printf( " Best value after init. %f ", errorPrev );
printf( "\n Position :\n" );
for ( d = 0; d < SS.D; d++ ) printf( " %f", R.SW.P[R.SW.best].x[d] );
*/
initLinks = 1; // So that information links will beinitialized
// Note: It is also a flag saying "No improvement"
noStop=0;
iter=0; iterBegin=0;
rho=param.exploit2;
explRate=param.exploit1;
// ---------------------------------------------- ITERATIONS
while (noStop == 0)
{
iter=iter+1;
if (initLinks==1) // Random topology
{
// Who informs who, at random
for (s = 0; s < R.SW.S; s++)
{
for (m = 0; m < R.SW.S; m++)
{
if (alea (0, 1)<p) LINKS[m][s] = 1; // Probabilistic method
else LINKS[m][s] = 0;
}
}
// Each particle informs itself
for (m = 0; m < R.SW.S; m++)
{ LINKS[m][m] = 1;
}
}
// Prepare the move
//printf("\nIteration %i",iter);
for (s = 0; s < R.SW.S; s++)
{
index[s]=s;
R.SW.Pp[s] = R.SW.P[s]; // Save the previous best, for information
R.SW.lBest[s]=0; // Initialise the flags (for information)
R.SW.X[s].forced=0; // Non "forced" move (for the moment)
R.SW.P[s].forced=0;
}
explForced=0; // In order the evaluate the exploitation rate
switch (param.randOrder)
{
default:
break;
case 1: //Define a random permutation
length=R.SW.S;
for (s=0;s<length;s++) indexTemp[s]=index[s];
for (s=0;s<R.SW.S;s++)
{
rank=alea_integer(0,length-1);
index[s]=indexTemp[rank];
if (rank<length-1) // Compact
{
for (t=rank;t<length;t++)
indexTemp[t]=indexTemp[t+1];
}
length=length-1;
}
break;
}
for (s0 = 0; s0 < R.SW.S; s0++) // For each particle ...
{
s=index[s0];
// ... find the first informant
s1 = 0;
while (LINKS[s1][s] == 0) s1++;
if (s1 >= R.SW.S) s1 = s; // If no informant at all, just himself
// Find the best informant
g[s] = s1;
for (m = s1; m < R.SW.S; m++)
{
if (LINKS[m][s] == 1 && R.SW.P[m].f < R.SW.P[g[s]].f)
g[s] = m;
}
R.SW.lBest[g[s]]=1; // Particle g[s] is a best informant
}
// Find the best of the best of all informants
sBest=0;
for (s = 1; s < R.SW.S; s++)
{
if(R.SW.P[g[s]].f < R.SW.P[g[sBest]].f) sBest=s;
}
// Define the exploitation area
// (when localSearch=0, this is just for information)
ex=exploitArea(R.SW, pb.SS,rho);
// For information, compute the size of the best local area
explBestSize=1;
for(d=0;d<pb.SS.D;d++)
explBestSize=
explBestSize*(ex.exI[g[sBest]].max[d]-ex.exI[g[sBest]].min[d]);
fprintf(f_expl, "explBestSize %f ",explBestSize);
// When the position is "forced", theoretically speaking the velocity
// of the particle should be also modified
// However, in practice, it apears it does not make any difference
switch(param.localSearch)
{
case 0:
break;
case 1: // At random in some local areas
for (s = 0; s < R.SW.S; s++) // For each particle ...
{
if(alea(0,1)>param.exploit1) break;
explForced=explForced+1;
for (d = 0; d < pb.SS.D; d++)
{
local.x[d]=ex.exI[g[s]].min[d]+
alea(0,1)*(ex.exI[g[s]].max[d]-ex.exI[g[s]].min[d]);
}
local = quantis (local, pb.SS);
local.f =perf(local,pb);
R.nEval = R.nEval + 1;
if(local.f<R.SW.X[s].f)
{
R.SW.X[s]=local;
R.SW.X[s].forced=1;
}
}
case 2: // At random in the best local area
for (s = 0; s < R.SW.S; s++) // For each particle ...
{
if(alea(0,1)>param.exploit1) break;// ... may or may not be forced
explForced=explForced+1;
for (d = 0; d < pb.SS.D; d++)
{
local.x[d]=ex.exI[g[sBest]].min[d]+
alea(0,1)*(ex.exI[g[sBest]].max[d]-ex.exI[g[sBest]].min[d]);
}
local = quantis (local, pb.SS);
local.f =perf(local,pb);
R.nEval = R.nEval + 1;
if(local.f<R.SW.X[s].f)
{
R.SW.X[s]=local;
R.SW.X[s].forced=1;
}
}
break;
case 3: // In the middle of some local areas
for (s = 0; s < R.SW.S; s++) // For each particle ...
{
if(alea(0,1)>param.exploit1) break;// ... may or may not be forced
if(R.SW.P[g[s]].forced==1) break; // No need to choose the
//middle several times!
explForced=explForced+1;
for (d = 0; d < pb.SS.D; d++)
{
local.x[d]=
0.5*(ex.exI[g[s]].max[d]+ex.exI[g[s]].min[d]);
}
local = quantis (local, pb.SS);
local.f =perf(local,pb);
R.nEval = R.nEval + 1;
if(local.f<R.SW.X[s].f)
{
R.SW.X[s]=local;
R.SW.X[s].forced=1;
R.SW.P[g[s]].forced=1;
}
}
break;
case 4: // In the middle of the best local exploitation area
s=g[sBest];
if(R.SW.P[g[s]].forced==1) break; // Already done
explForced=explForced+1;
for (d = 0; d < pb.SS.D; d++)
{
local.x[d]=
0.5*(ex.exI[g[sBest]].max[d]+ex.exI[g[sBest]].min[d]);
}
local = quantis (R.SW.X[s], pb.SS);
local.f =perf(R.SW.X[s],pb);
R.nEval = R.nEval + 1;
if(local.f<R.SW.X[s].f)
{
R.SW.X[s]=local;
R.SW.X[s].forced=1;
}
break;
case 5: // Gaussian in the best LEA
s=g[sBest];
if(R.SW.P[g[s]].forced==1) break; // Already done
explForced=explForced+1;
for (d = 0; d < pb.SS.D; d++)
{
zz=0.5*(ex.exI[g[sBest]].max[d]-ex.exI[g[sBest]].min[d]);
zzz=alea_normal(0,1.48*zzz); // For a fifty-fifty probability distribution
if(fabs(zzz)>zz)
{ if (zzz>0) zzz=zz; else zzz=-zz;};
local.x[d]=
zzz+0.5*(ex.exI[g[sBest]].max[d]+ex.exI[g[sBest]].min[d]);
}
local = quantis (R.SW.X[s], pb.SS);
local.f =perf(R.SW.X[s],pb);
R.nEval = R.nEval + 1;
if(local.f<R.SW.X[s].f)
{
R.SW.X[s]=local;
R.SW.X[s].forced=1;
}
break;
}
//** Move, for non "forced" particles =====================
for (s0 = 0; s0 < R.SW.S; s0++) // For each particle ...
{
s=index[s0];
if(R.SW.X[s].forced==1) goto update_best;
// Exploration tendency for (d = 0; d < pb.SS.D; d++)
{
R.SW.V[s].v[d]=R.SW.w[s]*R.SW.V[s].v[d];
}
// Prepare Exploitation tendency
for (d = 0; d < pb.SS.D; d++)
{
PX.v[d]= R.SW.P[s].x[d] - R.SW.X[s].x[d];
}
PX.size=pb.SS.D;
if(g[s]!=s)
{
for (d = 0; d < pb.SS.D; d++)
{
GX.v[d]= R.SW.P[g[s]].x[d] - R.SW.X[s].x[d];
}
GX.size=pb.SS.D;
}
// Exploitation tendencies
//c1=param.c1; c2=param.c2; c1=R.SW.c[s]; c2=c1;
switch(R.SW.R[s])
{
case 0:
for (d = 0; d < pb.SS.D; d++)
{
R.SW.V[s].v[d]=R.SW.V[s].v[d] +
+ c1*alea(0, 1)*PX.v[d];
}
if (g[s]!=s) // If the best neighbour is not the particle itself
{
for (d = 0; d < pb.SS.D; d++)
{
R.SW.V[s].v[d]=R.SW.V[s].v[d]
+ c2*alea(0,1) * GX.v[d];
}
}
break;
case 1:
for (d = 0; d < pb.SS.D; d++)
{
R.SW.V[s].v[d]=R.SW.V[s].v[d] +
+ c1*alea_normal(alea(0,1), param.sigma)*PX.v[d];
}
if (g[s]!=s)
{
for (d = 0; d < pb.SS.D; d++)
{
R.SW.V[s].v[d]=R.SW.V[s].v[d]
+ c2*alea_normal(alea(0,1), param.sigma) * GX.v[d];
}
}
break; case 2: for (d = 0; d < pb.SS.D; d++)
{
R.SW.V[s].v[d]=R.SW.V[s].v[d] + c1*alea(0, 1)*PX.v[d];
}
if (g[s]!=s) // If the best neighbour is not the particle itself
{
for (d = 0; d < pb.SS.D; d++)
{
zz=fabs(GX.v[d])/(pb.SS.max[d]-pb.SS.min[d]); zzz=(1-zz); R.SW.V[s].v[d]=R.SW.V[s].v[d]
+ c2*alea(0,1) *pow(zzz,distPow) *GX.v[d];
}
}
}
// Update the position
for (d = 0; d < pb.SS.D; d++)
{
R.SW.X[s].x[d] = R.SW.X[s].x[d] + R.SW.V[s].v[d];
}
if (R.nEval >= pb.evalMax) goto end;
// --------------------------
noEval = 1;
// Quantisation
R.SW.X[s] = quantis (R.SW.X[s], pb.SS);
switch (param.clamping)
{
case 0: // No clamping AND no evaluation
outside = 0;
for (d = 0; d < pb.SS.D; d++)
{
if (R.SW.X[s].x[d] < pb.SS.min[d] || R.SW.X[s].x[d] > pb.SS.max[d])
outside++;
}
if (outside == 0) // If inside, the position is evaluated
{
R.SW.X[s].f =
perf (R.SW.X[s], pb);
R.nEval = R.nEval + 1;
}
break;
case 1: // Set to the bounds, and v to zero outside=0;
for (d = 0; d < pb.SS.D; d++)
{
if (R.SW.X[s].x[d] < pb.SS.minS[d])
{
R.SW.X[s].x[d] = pb.SS.minS[d];
R.SW.V[s].v[d] = 0; if(pb.SS.minS[d]<pb.SS.min[d]-zero) outside++;
}
if (R.SW.X[s].x[d] > pb.SS.maxS[d])
{
R.SW.X[s].x[d] = pb.SS.maxS[d];
R.SW.V[s].v[d] = 0; if(pb.SS.maxS[d]>pb.SS.max[d]+zero) outside++;
}
}
if(outside==0) {
R.SW.X[s].f =perf(R.SW.X[s],pb);
R.nEval = R.nEval + 1; }
break;
}
update_best:
// ... update the best previous position
if (R.SW.X[s].f < R.SW.P[s].f) // Improvement
{
R.SW.P[s] = R.SW.X[s];
// ... update the best of the bests
if (R.SW.P[s].f < R.SW.P[R.SW.best].f)
{
R.SW.best = s;
}
}
} // End of "for (s0=0 ... "
explRate0=explRate;
explRate=explNb(R.SW,ex); // Exploitation rate
// WARNING: count only "successful" moves
/*
// The following count takes into account all particles that have been
// forced into a local area, even if the final position
// has not been kept (unsuccessful)
explRate=(double)explForced/R.SW.S;
fprintf(f_expl,"%f\n",explRate);
*/
// Modify the future size of the local area
// in order to "balance" exploitation/exploration
switch(param.explOption)
{
default: // Keep constant *
break;
case 1:
rho=pow(pow(rho,pb.SS.D)+explRate-explRate0,1./pb.SS.D);
break;
case 2:
if(explRate-explRate0>0)
{
rho=rho*pow( 1-(explRate-explRate0)/(1-explRate0),1./pb.SS.D);
}
else
{
rho=pow(pow(rho,pb.SS.D) +
(explRate-explRate0)*(pow(rho,pb.SS.D)-1)/explRate0,1./pb.SS.D);
}
break;
case 3:
if(alea(0,1)<0.5) rho=param.exploit2;
else rho=1./param.S;
break;
case 4: // Random Uniform
rho=alea(0,param.exploit2);
break;
// case 5: // Random Gauss
// rho=fabs(alea_normal(1./param.S,param.exploit2));
rho=fabs(alea_normal(0,param.exploit2));
// rho=fabs(alea_normal(1./param.S,1));
if (rho>1) rho=1;
break;
}
//printf("\nExploitation rate %f",explRate);
//printf(" => rho %f",rho);
// Check if finished
error = R.SW.P[R.SW.best].f;
if (error < errorPrev) // Improvement
{
initLinks = 0;
} else // No improvement
{
initLinks = 1; // Information links will be reinitialized
} // Prepare the curve "success rate" vs FEmax if(R.nEval>FEmax[FEclass]) { if(error<=pb.epsilon) successFlag[run][FEclass]=1; else successFlag[run][FEclass]=0;//printf("\n %i flag %i",FEmax[FEclass],successFlag[run][FEclass]); FEclass=FEclass+1; }
errorPrev = error;
end:
switch (param.stop)
{
case 0:
if (error > pb.epsilon && R.nEval < pb.evalMax) noStop = 0; // Won't stop
else noStop = 1; // Will stop
break;
case 1:
if (R.nEval < pb.evalMax) noStop = 0; // Won't stop
else noStop = 1; // Will stop
break;
}
} // End of "while nostop ...
// printf( "\n and the winner is ... %i", R.SW.best );
// fprintf( f_stag, "\nEND" );
R.error = error;
return R;
}
// ===========================================================
double alea (double a, double b)
{ // random number (uniform distribution) in [a b] // randOption is a global parameter (see also param.rand)
double r; if(randOption==0) r=a+(double)rand_kiss()*(b-a)/RAND_MAX_KISS; else r =a+ (double) rand ()*(b-a)/RAND_MAX;
return r;
}
// ===========================================================int alea_integer (int a, int b)
{ // Integer random number in [a b]
int ir; double r;
r = alea (0, 1); ir = (int) (a + r * (b + 1 - a));
if (ir > b) ir = b;
return ir;
}
// ===========================================================
double alea_normal (double mean, double std_dev)
{
/* Use the polar form of the Box-Muller transformation to obtain a pseudo
random number from a Gaussian distribution */
double x1, x2, w, y1;
// double y2;
do
{
x1 = 2.0 * alea (0, 1) - 1.0;
x2 = 2.0 * alea (0, 1) - 1.0;
w = x1 * x1 + x2 * x2;
} while (w >= 1.0);
w = sqrt (-2.0 * log (w) / w); y1 = x1 * w; // y2 = x2 * w;
if(alea(0,1)<0.5) y1=-y1;
y1 = y1 * std_dev + mean;
return y1;
}
/*
// =============================================================
struct velocity aleaVector(int D,double coeff)
{
struct velocity V;
int d;
int i;
int K=2; // 1 => uniform distribution in a hypercube
// 2 => "triangle" distribution
double rnd;
V.size=D;
for (d=0;d<D;d++)
{
rnd=0;
for (i=1;i<=K;i++) rnd=rnd+alea(0,1);
V.v[d]=rnd*coeff/K;
}
return V;
}
*/
//================================================
static int compareLBest (void const *a, void const *b)
{
// For qsort (increasing order of any list of real numbers)
struct lBest const *pa = a;
struct lBest const *pb = b;
// NOTE: needs the global variable dLBest
if(pa->x[dLBest] > pb->x[dLBest]) return 1;
if(pa->x[dLBest] < pb->x[dLBest]) return -1;
return 0;
}
//============================================================
struct exploit exploitArea(struct swarm SW, struct SS SS, double rho)
{
// Define the exploitation area
int d;
struct exploit ex; // local search space (D, max, min, quantum)
// around each local best
int i;
struct lBest lB[S_max+2];
int nbLBest; // Number of different local bests
//double rho; // Define what "around" means. Should be <= 0.5
int s;
//printf("\nExploitation area with rho=%f",rho);
/*
For each local best, define the local exploitation area (LEA)
(hyperparallelepid "around" it.
*/
// Build the list of local bests (index, coordinates) + min + max
i=1;
for(s=0;s<SW.S;s++)
{
if(SW.lBest[s]==0) continue;
lB[i].rank=s;
for (d=0;d<SS.D;d++) lB[i].x[d]=SW.Pp[s].x[d];
i=i+1;
}
lB[0].rank=-1; // "Left" border point
lB[i].rank=-2; // "Right" border point
for (d=0;d<SS.D;d++)
{
lB[0].x[d]=SS.min[d];
lB[i].x[d]=SS.max[d];
}
nbLBest=i-2; // Number of local bests
// For each dimension
// Sort the list
// Define the intervals
for(d=0;d<SS.D;d++)
{
// Sort the list by ascending order
dLBest=d; // Global variable. Needed for compareLBest
qsort(lB,nbLBest+2,sizeof(lB[0]),compareLBest);
// For each local best
// add the interval in dimension d that defines the LEA
for (i=0;i<SW.S;i++)
{
ex.rank[i]=0; // Means N/A
ex.exI[i].D=SS.D; // For future use (local search)
ex.exI[i].q.size=SS.D; // "
}
for (i=1;i<nbLBest+2;i++)
{
ex.rank[lB[i].rank]=1;
ex.exI[lB[i].rank].min[d]=lB[i].x[d]-rho*(lB[i].x[d]-lB[i-1].x[d]);
ex.exI[lB[i].rank].max[d]=lB[i].x[d]+rho*(lB[i+1].x[d]-lB[i].x[d]);
ex.exI[lB[i].rank].q.q[d]=SS.q.q[d];
}
}
return ex;
}
//==================================================
double explNb(struct swarm SW, struct exploit ex )
{
// For each position
// check if it is inside a LEA
int d;
int Dim;
int exNb=0;
double rate;
int i;
int inExArea;
int s;
Dim=ex.exI[0].D;
for (s=0;s<SW.S;s++) // For each position
{
for(i=0;i<SW.S;i++) // For each "exploitation area" (hyperparallelepid)
{
if(ex.rank[i]<=0) continue;
inExArea=1;
for(d=0;d<Dim;d++) // For each "exploitation interval"
{
if(SW.X[s].x[d]<=ex.exI[i].min[d] || SW.X[s].x[d]>= ex.exI[i].max[d])
{inExArea=0; break;}
}
if(inExArea==1)
{
exNb=exNb+1;
break; // Counted just once
}
}
}
rate=(double)exNb/SW.S;
fprintf(f_expl,"%f\n",rate); // Rate of particles that "exploit"
return rate;
}
//==================================================
struct position hammersley(struct SS SS, int N, int deb, int s)
{
struct position x;
double a;
int d,D;
int dd;
int kp;
int pp;
D=SS.D;
// A sequence of prime numbers. You may modify it
int P[100]={2,3,5,7,11, 13, 17, 19, 23, 29,
31, 37, 41, 43, 47, 53, 59, 61, 67, 71,
73,79,83,89,97,101,103,107,109,113,
127,131,137,139,149,151,157,163,167,173,
179, 181,191,193,197,199,211,223,227,229,
233,239,241,251,257,263,269,271,277,281,
283,293,307,311,313,317,331,337,347,349,
353,359,367,373,379,383,389,397,401,409,
419,421,431,433,439,443,449,457,461,463,
467,479,487,491,499,503,509,521,523,541
};
double z;
if(D>100) { printf("\nDimension D = %i is too high\n",D); ERROR("You must increase the list of prime numbers"); }
x.x[deb]=SS.min[deb]+(SS.max[0]-SS.min[0])*(double)s/(N-1);
for(d=1;d<D;d++) // for each dimension
{ // compute the position
dd=deb+d; if(dd>=D) dd=0;
pp=P[dd];
kp=s+1;
z=0;
while (kp>0)
{
a=kp % P[dd]; // kp modulo P[dd]
z=z+(double)a/pp;
kp=kp/P[dd]; // Integer division
pp=pp*P[dd];
}
x.x[dd]=SS.min[dd]+z*(SS.max[dd]-SS.min[dd]);
}
/*
fprintf(f_run,"\n");
for(d=0;d<D;d++)
fprintf(f_run,"%f ",(double)x.x[d]);
*/
return x;
}
//================================================== KISS
/*
A good pseudo-random numbers generator
The idea is to use simple, fast, individually promising
generators to get a composite that will be fast, easy to code
have a very long period and pass all the tests put to it.
The three components of KISS are
x(n)=a*x(n-1)+1 mod 2^32
y(n)=y(n-1)(I+L^13)(I+R^17)(I+L^5),
z(n)=2*z(n-1)+z(n-2) +carry mod 2^32
The y's are a shift register sequence on 32bit binary vectors
period 2^32-1;
The z's are a simple multiply-with-carry sequence with period
2^63+2^32-1. The period of KISS is thus
2^32*(2^32-1)*(2^63+2^32-1) > 2^127
*/
static ulong kiss_x = 1;
static ulong kiss_y = 2;
static ulong kiss_z = 4;
static ulong kiss_w = 8;
static ulong kiss_carry = 0;
static ulong kiss_k;
static ulong kiss_m;
void seed_rand_kiss(ulong seed)
{
kiss_x = seed | 1;
kiss_y = seed | 2;
kiss_z = seed | 4;
kiss_w = seed | 8;
kiss_carry = 0;
}
ulong rand_kiss()
{
kiss_x = kiss_x * 69069 + 1;
kiss_y ^= kiss_y << 13;
kiss_y ^= kiss_y >> 17;
kiss_y ^= kiss_y << 5;
kiss_k = (kiss_z >> 2) + (kiss_w >> 3) + (kiss_carry >> 2);
kiss_m = kiss_w + kiss_w + kiss_z + kiss_carry;
kiss_z = kiss_w;
kiss_w = kiss_m;
kiss_carry = kiss_k >> 30;
return kiss_x + kiss_y + kiss_w;
}
/*
// ===========================================================
double normL (struct velocity v,double L)
{ // L-norm of a vector
int d;
double n;
n = 0;
for (d = 0; d < v.size; d++) n = n + pow(fabs(v.v[d]),L);
n = pow (n, 1/L);
return n;
} */
// ===========================================================
double distanceL (struct position x1, struct position x2,double L)
{ // Distance between two positions
// L = 2 => Euclidean
int d;
double n;
n = 0;
for (d = 0; d < x1.size; d++) n = n + pow (fabs(x1.x[d] - x2.x[d]), L);
n = pow (n, 1/L); return n;
}
// ===========================================================
int sign (double x)
{
//if (x == 0) return 0; if (x <= 0) return -1;
return 1;
}
// ===========================================================
struct position quantis (struct position x, struct SS SS)
{
/* Quantisatition of a position
Only values like x+k*q (k integer) are admissible */
int d; double qd; struct position quantx;
quantx = x;
for (d = 0; d < x.size; d++) { qd = SS.q.q[d];
if (qd > zero) // Note that qd can't be < 0
{
quantx.x[d] = qd * floor (0.5 + x.x[d] / qd);
}
}
return quantx;
}
// ===========================================================
double perf (struct position x, struct problem pb)
{ // Evaluate the fitness value for the particle of rank s
double c; int d; double dist,dMax;
double DD; int dim;
int k; int n,n1; double f, f1,f2; double p;
double s11, s12, s21, s22;
double sum1,sum2;
double t0, tt, t1; double x1,x2,x3,x4; double xd;
struct position xs; double y;
double z,z2;
#include "cec2005data.c"
// Data for center-bias test function (code 6)double radius=10; // int shift=1; // Set to zero in order to have the center on {0, ...,0}, to 1 elsedouble center[30]={80, -20, 1.9231600e+001, -2.5231000e+000, 7.0433800e+001,
4.7177400e+001, -7.8358000e+000, -8.6669300e+001, 5.7853200e+001, -9.9533000e+000,
2.0777800e+001, 5.2548600e+001, 7.5926300e+001, 4.2877300e+001, -5.8272000e+001,
-1.6972800e+001, 7.8384500e+001, 7.5042700e+001, -1.6151300e+001, 7.0856900e+001,
-7.9579500e+001, -2.6483700e+001, 5.6369900e+001, -8.8224900e+001, -6.4999600e+001,
-5.3502200e+001, -5.4230000e+001, 1.8682600e+001, -4.1006100e+001, -5.4213400e+001};struct position summit; // Variables specific to Coil compressing spring static double Fmax=1000.0; static double Fp=300;⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -