📄 standard_pso_2007.c
字号:
printf ("\n Eval. (mean)= %f", evalMean);
printf ("\n Error (mean) = %f", errorMean);
// Variance variance = 0;
for (run = 0; run < runMax; run++) variance = variance + pow (errorMeanBest[run] - errorMean, 2);
variance = sqrt (variance / runMax);
printf ("\n Std. dev. %f", variance);
printf("\n Log_progress (mean) = %f", logProgressMean); // Success rate and minimum value
printf("\n Failure(s) %i",nFailure); successRate = 100 * (1 - nFailure / (double) runMax);
printf ("\n Success rate = %.2f%%", successRate);
if (run > 1) printf ("\n Best min value = %f", errorMin); // Save /* fprintf(f_synth,"\n"); for (d=0;d<SS.D;d++) fprintf(f_synth,"%f ",
pb.offset[d]); fprintf(f_synth," %f %f %f %.0f%% %f",errorMean,variance,errorMin,
successRate,evalMean); fprintf( f_synth, "\n%f %f %f %f %.0f%% %f ", shift, errorMean, variance, errorMin, successRate, evalMean );
*/ fprintf (f_synth, "\n");
fprintf (f_synth, "%f %f %.0f%% %f ", errorMean, variance, successRate, evalMean); return 0; // End of main program
}// ===============================================================// PSOstruct result PSO (struct param param, struct problem pb) {
struct velocity aleaV;
int d;
double error;
double errorPrev;
struct velocity expt1,expt2;
int g;
struct velocity GX;
int index[S_max], indexTemp[S_max];
int initLinks; // Flag to (re)init or not the information links
int iter; // Iteration number (time step)
int iterBegin;
int length;
int LINKS[S_max][S_max]; // Information links int m;
int noEval;
double normPX, normGX;
int noStop;
int outside;
double p;
struct velocity PX;
struct result R;
int rank;
struct matrix RotatePX;
struct matrix RotateGX;
int s0, s,s1;
int t;
double zz;
aleaV.size=pb.SS.D;
RotatePX.size=pb.SS.D;
RotateGX.size=pb.SS.D;
// ----------------------------------------------------- // INITIALISATION
p=param.p; // Probability threshold for random topology R.SW.S = param.S; // Size of the current swarm // Position and velocity for (s = 0; s < R.SW.S; s++) {
R.SW.X[s].size = pb.SS.D; R.SW.V[s].size = pb.SS.D;
for (d = 0; d < pb.SS.D; d++) { R.SW.X[s].x[d] = alea (pb.SS.min[d], pb.SS.max[d]);
}
for (d = 0; d < pb.SS.D; d++)
{
R.SW.V[s].v[d] =
(alea( pb.SS.min[d], pb.SS.max[d] ) - R.SW.X[s].x[d])/2;
}
}
// 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.function,pb.SS);
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 R.SW.best = 0;
switch (param.stop)
{
default:
errorPrev = fabs(pb.epsilon-R.SW.P[R.SW.best].f);
break;
case 2:
errorPrev=distanceL(R.SW.P[R.SW.best],pb.solution,2);
break;
}
for (s = 1; s < R.SW.S; s++) {
switch (param.stop)
{
default:
zz=fabs(pb.epsilon-R.SW.P[s].f);
if (zz < errorPrev)
R.SW.best = s;
errorPrev=zz;
break;
case 2:
zz=distanceL(R.SW.P[R.SW.best],pb.solution,2);
if (zz<errorPrev)
R.SW.best = s;
errorPrev=zz;
break;
}
}
/* // 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; // ---------------------------------------------- ITERATIONS
iter=0; iterBegin=0; 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;
}
} // The swarm MOVES //printf("\nIteration %i",iter);
for (s = 0; s < R.SW.S; s++) index[s]=s;
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;
// Find the best informant
g = s1;
for (m = s1; m < R.SW.S; m++)
{
if (LINKS[m][s] == 1 && R.SW.P[m].f < R.SW.P[g].f)
g = m;
}
//.. compute the new velocity, and move
// Exploration tendency
for (d = 0; d < pb.SS.D; d++)
{
R.SW.V[s].v[d]=param.w *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)
{
for (d = 0; d < pb.SS.D; d++)
{
GX.v[d]= R.SW.P[g].x[d] - R.SW.X[s].x[d];
}
GX.size=pb.SS.D;
}
// Option "non sentivity to rotation"
if (param.rotation>0)
{
normPX=normL(PX,2);
if (g!=s) normGX=normL(GX,2);
if(normPX>0)
{
RotatePX=matrixRotation(PX);
}
if(g!= s && normGX>0)
{
RotateGX=matrixRotation(GX);
}
}
// Exploitation tendencies
switch (param.rotation)
{
default:
for (d = 0; d < pb.SS.D; d++)
{
R.SW.V[s].v[d]=R.SW.V[s].v[d] +
+ alea(0, param.c)*PX.v[d];
}
if (g!=s)
{
for (d = 0; d < pb.SS.D; d++)
{
R.SW.V[s].v[d]=R.SW.V[s].v[d]
+ alea(0,param.c) * GX.v[d];
}
}
break;
case 1:
// First exploitation tendency
if(normPX>0)
{
zz=param.c*normPX/sqrtD;
aleaV=aleaVector(pb.SS.D, zz);
expt1=matrixVectProduct(RotatePX,aleaV);
for (d = 0; d < pb.SS.D; d++)
{
R.SW.V[s].v[d]=R.SW.V[s].v[d]+expt1.v[d];
}
}
// Second exploitation tendency
if(g!=s && normGX>0)
{
zz=param.c*normGX/sqrtD;
aleaV=aleaVector(pb.SS.D, zz);
expt2=matrixVectProduct(RotateGX,aleaV);
for (d = 0; d < pb.SS.D; d++)
{
R.SW.V[s].v[d]=R.SW.V[s].v[d]+expt2.v[d];
}
}
break;
}
// 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.function, pb.SS);
R.nEval = R.nEval + 1;
}
break;
case 1: // Set to the bounds, and v to zero
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.min[d];
R.SW.V[s].v[d] = 0;
}
if (R.SW.X[s].x[d] > pb.SS.max[d])
{
R.SW.X[s].x[d] = pb.SS.max[d];
R.SW.V[s].v[d] = 0;
}
}
R.SW.X[s].f =perf(R.SW.X[s],pb.function, pb.SS);
R.nEval = R.nEval + 1;
break;
} // ... 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 (s=0 ... "
// Check if finished
switch (param.stop)
{
default:
error = R.SW.P[R.SW.best].f;
break;
case 2:
error=distanceL(R.SW.P[R.SW.best],pb.solution,2);
break;
} error= fabs(error - pb.epsilon);
if (error < errorPrev) // Improvement {
initLinks = 0;
} else // No improvement {
initLinks = 1; // Information links will be reinitialized
}
errorPrev = error; end:
switch (param.stop) {
case 0:
case 2:
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 ); ⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -