模拟退火算法解决路径优化 的源代码.docx

模拟退火算法解决路径优化 的源代码.docx

《模拟退火算法解决路径优化 的源代码.docx》由会员分享，可在线阅读，更多相关《模拟退火算法解决路径优化 的源代码.docx（9页珍藏版）》请在冰点文库上搜索。

模拟退火算法解决路径优化的源代码

ÎÒÓÐsimulatedannealingwithmetropolies（MonteCarlo）×öµÄÒ»¸öÏîÄ¿µÄ´úÂ룬ÄãÒª¿´¿´Ã´£¿

voidanneal（intnparam,intnstep,intnstep_per_block,doublet0,

constdouble*param_in,

doublecost_in,double*params_out,double*cost_out）{

intnblock;

intstep;

intblock;

intnactive;

intrank;

intn_accepted=0;

inti,j,n;

doublecost_current,cost_trial;

int*param_index;

double*param_current;

double*param_trial;

double*Q;

double*S;

double*u;

double*dp;

double*A;

FILE*fp_log_file;

charfname[FILENAME_MAX];

doubletemp=t0;

doubletempmax=temp;

doubleebar,evar,emin,eta,specific_heat;

doubledelta;

doublechi=0.8;//Annealingschedule

doublechi_s=3.0;//Vanderbilt/Louie'growthfactor'

doublerm;

doubleroot3=sqrt（3.0）;

doublep=0.02/sqrt（3.0）;//maxsizeofannealingstep

param_current=newdouble[nparam];

param_trial=newdouble[nparam];

cost_current=cost_in;

MPI_Comm_rank（MPI_COMM_WORLD,&rank）;

sprintf（fname,"a_%4.4d.log",rank）;

fp_log_file=fopen（fname,"a"）;

if（fp_log_file==（FILE*）NULL）errorMessage（"fopen（log）failed\n"）;

//Workoutthenumberofactiveparameters,andsetupthe

//indextableoftheactiveparameters.

//Notethatthecompletearrayofparameters（param_trial）must

//beusedtoevaluatethecostfunction.

nactive=0;

for（n=0;n

param_current[n]=param_in[n];

param_trial[n]=param_in[n];

if（P.is_active[n]）nactive++;

}

param_index=newint[nactive];

i=0;

for（n=0;n

if（P.is_active[n]）param_index[i++]=n;

}

//InitialisethestepdistributionmatrixQ_ij

Q=newdouble[nactive*nactive];

S=newdouble[nactive*nactive];

u=newdouble[nactive];

dp=newdouble[nactive];

A=newdouble[nactive];

double*Qtmp;

Qtmp=newdouble[nactive*nactive];

for（i=0;i

for（j=0;j

delta=（i==j）;

Q[i*nactive+j]=p*delta*param_current[param_index[j]];

}

}

//carryoutannealingpoints

nblock=nstep/nstep_per_block;

rm=1.0/（double）nstep_per_block;

for（block=0;block

//Setthescheduleforthisblock,andinitialiseblockwisequantities.

//Wealsoensurethestepdistributionmatrixisdiagonal.

temp=chi*temp;

for（i=0;i

A[i]=0.0;

for（j=0;j

S[i*nactive+j]=0.0;

delta=（i==j）;

Q[i*nactive+j]*=delta;

}

}

ebar=0.0;

evar=0.0;

emin=cost_current;

for（i=0;i

printf（"Step:

%d%g\n",i,Q[i*nactive+i]）;

}

for（step=0;step

//Settherandomvectoru,andcomputethestepsizedp

for（i=0;i

u[i]=root3*（r_uniform（）*2.0-1.0）;

}

for（i=0;i

dp[i]=0.0;

for（j=0;j

dp[i]+=Q[i*nactive+j]*u[j];

}

}

for（i=0;i

n=param_index[i];

param_trial[n]=param_current[n]+dp[i];

if（param_trial[n]

if（param_trial[n]>P.max[n]）param_trial[n]=P.max[n];

}

//calculatenewcostfunctionscore

p_model->setParameters（param_trial）;

cost_trial=p_costWild->getCost（）;

cost_trial+=p_costLHY->getCost（）;

cost_trial+=p_costTOC1->getCost（）;

cost_trial+=p_costAPRR->getCost（）;

//Metropolis

delta=cost_trial-cost_current;

if（delta<0.0||r_uniform（）

for（n=0;n

param_current[n]=param_trial[n];

}

cost_current=cost_trial;

++n_accepted;

}

//'Energy'statistics

ebar+=cost_current;

evar+=cost_current*cost_current;

if（cost_current

//Pertimesteplog

fprintf（fp_log_file,"%6d%6d%10.4f%10.4f%10.4f%10.4f\n",

block,step,temp,

cost_current,cost_trial,

（float）n_accepted/（float）（block*nstep_per_block+step））;

//Accumulateaverage,measuredcovariance

for（i=0;i

A[i]+=param_current[param_index[i]];

for（j=0;j

S[i*nactive+j]

+=param_current[param_index[i]]*param_current[param_index[j]];

}

}

/*Nextstep*/

}

//Setthepreviousblockaverageandmeasuredcovariance

for（i=0;i

A[i]=rm*A[i];

}

for（i=0;i

for（j=0;j

S[i*nactive+j]=rm*S[i*nactive+j]-A[i]*A[j];

if（i==j）printf（"Average:

%d%g%g\n",i,A[i],S[i*nactive+j]）;

//Settheconvarienceforthenextiterations=6chi_sS/M

S[i*nactive+j]=6.0*chi_s*rm*S[i*nactive+j];

}

}

//Resetthestepdistributionmatrixforthenextblock

i=do_cholesky（nactive,S,Q）;

j=test_cholesky（nactive,S,Q）;

printf（"Cholesky%d%d\n",i,j）;

//Blockstatistics

ebar=rm*ebar;

evar=rm*evar;

specific_heat=（evar-ebar*ebar）/temp*temp;

eta=（ebar-emin）/ebar;

fprintf（fp_log_file,"%d%d%f%f%f%f%f%f\n",

block,nstep_per_block,temp,ebar,evar,emin,

specific_heat,eta）;

/*Nextblock*/

}

*cost_out=cost_current;

for（n=0;n

params_out[n]=param_current[n];

}

fclose（fp_log_file）;

deleteparam_index;

deleteparam_current;

deleteparam_trial;

deleteQ;

deleteu;

deletedp;

deleteS;

deleteA;

return;

}