新手小白，求助各位大神一份关于模拟退火算法的仿真计算代码，具体内容帖子里写了，不胜感激

内容是这样的：
课题是关于宽带光纤拉曼放大器（RFA）的优化设计，现在想选用的算法是基于遗传算法的模拟退火算法，但是不会进行编程仿真，找遍论文也没有找到相关的源代码，所以想求各位大神帮帮忙，帮忙写下仿真的代码，matlab或者C++6.0均可，不胜感激，目标函数、数值计算参数和相关的算法流程下面有写：
目标函数：
数值计算参数：
算法流程：
下面有一份做好的源代码，但是错误较多，不能运行，希望可以给各位大神一个参考，希望诸位大神能解我燃眉之急，不胜感激！
#include<iostream.h>                double interpl(int t,double u);
#include<math.h>                void main()
#include<stdio.h>                {
#include<stdlib.h>                    int ml=100;
const L=50000;                        double r=0.95;
const m=4;                        double stepfactor=0.002;
const n=40;                        double double T=100;
const h=1000;                        double q[8],x0[2*m],x2[2*m];
const 1=L/h+1;                        double array[8];
double function(double[]);                int i,j,k;
double raman_gain(double fP,double fS);            double y1,y2;
double p[8][2*m]={{1423,1446,1458,1500,2887,3643,3910,13524},{1424,1440,1458,1496,3000,3643,3910,13520},{1424,1440,1460,1480,2900,3645,3920,12000},{1424,1444,1464,1486,3000,3650,4000,12000},{1425,1445,1460,1495,2887,3643,3906,13254},{1425,1448,1456,1494,3100,3710,3950,13500},{1426,1448,1468,1488,2880,3500,3890,12900},{1426,1448,1470,1490,2980,3700,4000,13670}};
double c;                    {  for(j=0;j<2*m;j++)
int a=0;                        array[j]=p[i][j];
double a1,a2;                        q[i]=function(array);}
do{for(i=0;i<8;i++)
    double ma=(q[0]+q[1]+q[2]+q[3]+q[4]+q[5]+q[6]+q[7])/8;
    for(i=0;i<7;i++)
    { if((q[i]-ma)*(q[i+1]-ma)<0{  for(j=0;j<2*m;j++){
        y1=p[i][j]; y2=p[i+1][j]; p[i][j]=0.9*y1+0.1*y2; p[i+1][j]=0.9*y2+0.1*y1;
}}
    else if((q[i]-ma>0)
    {    for(j=0;j<2*m;j++)   p[i][j]=p[i][j]*(1+0.004*((rand()/(double)RAND_MAX-0.5));}}
if((q[7]-ma)>0)
{    for(j=0;j<2*m;j++)  p[7][j]=p[7][j]*(1+stepfactor*((rand()/(double)RAND_MAX)-0.5));}
for(i=0;i<8;i++)                {  if(q[i]<b)
{for(fj=0;j<2*m;j++)            {for(j=0;j<8;j++)
array[j]=p[i][j];            {x0[j]=p[i][j];}
q[i]=function(array);}            b=q[i];
double b=q[0];                d=i;}}
int d=0;            for(j=0;j<8;j++)
for(i=0;i<8;i++)            {p[d][j]=0;}
c=b;
for(i=0;i<m1;i++)
{ for(j=0;j<8;j++)
x2[j]=x0[j]*(1+stepfactor*((rand()/(double)RAND_MAX)-0.5));
a1=function(x0);    a2=function(x2);
if(a2<a1)
{ for(k=0;k<8;k++)  x0[k]=x2[k];
c=a2; i=m1;}
else
{ if(exp(-(a2-a1)/T)>(rand()/(double)RAND_MAX))
{ for(k=0;k<8;k++)
    x0[k]=x2[k];}}
cout<<c<<'\n';
a++;
cout<<a<<'\n';}
T=T*r;
for(i=0;i<8;i++)
    cout<<x0[i]<<'\n';
for(i=0;i<8;i++)
{    for(j=0;j<2*m;j++)
{    if(p[i][j]!=0)
        j=2*m;
    else{    p[i][j]=x0[j];}}}
}while(c>0.025);
return;}
double function(double p[])
{ int i,j;  double F[n],S[m],Slamda[n],fS[n],fP[m];
  for(i=0;i<n;i++)
{ Slamda[i]=1530+i*2;  fS[i]=3e17/Slamda[i];}
  for(i=0;i<m;i++)
{ fP[i]=3e8/p[i]/1e-9;  S[i]=p[i+m];}
double gPS[m][n];
for(i=0;i<m;i++)
{  for(j=0;j<m;j++)  gPS[i][j]=raman_gain(fP[i],fS[j])/50e-12/2;}
for(i=0;i<n;i++)
{  double k=0;
   for(j=0;j<m;j++)   k=k+4.343*gPS[j][i]*S[j];
  F[i]=(fabs(k-20))/20;}
double a=0;
for(i=0;i<n;i++)
{ if(F[i]>a)
    a=F[i];}
return(a);}