新安江遗传算法c语言,基于遗传算法的新安江模型参数优化率定（四）
4.3.1新安江三水源模型
//新安江三水源模型.hios
#include?算法
#include?函數
#include?優化
#include?spa
const?intVariableNum = 11,?//待率定參數個數.net
M = 485,//其為降雨起止間時段數(應該為)blog
Mq = 14;//單位線中q[i]的個數，此課設中為ip
const?double? F = 686, dt = 6;//面積、系列數據的時間間隔get
const?double?FE =0.8;input
const?int? Days[4] = {30, 31, 30, 31};?//四至七月份每個月天數
const?doubleAveE[2][4] = {1.57, 2.29, 2.65, 3.41,0.97, 1.49, 1.71,???? 2.34};//四至七月份每個月的多年平均日蒸發量
double?K, WUM, WLM, C,?//蒸發(蒸散發能力折算系數、分層蓄水容量、深層蒸散發系數)
WM,?b,????//產流(蓄水容量、蓄水容量曲線指數)
SM,? EX,?KI,???//水源劃分(自由水蓄水容量、自由水蓄水容量曲線指數、自由水蓄水庫出流系數)
KKI, KKG; ????????????//匯流(壤中流、地下徑流消退系數)
doubleP[M], Ep[M], EU[M], EL[M], ED[M], E[M], PE[M], WU[M + 1], WL[M + 1], WD[M + 1],W[M], a[M], R[M];
doubleaf[M];//af[i]指產流面積比率(%),三水源劃分中須要的數據
doubleP1[M], P2[M], P3[M], Qo[ M ];?????//三個雨量站實測雨量，實測流域流量
doubleS[M], AU[M], RS[M], RI[M], RG[M];
doubleq[Mq], QS[Mq + M - 1],Qs[Mq + M - 1][M];
doubleQ[Mq + M - 1], QI[Mq + M - 1], QG[Mq + M - 1];
doubleSumQo, AveQo, S1, S2;
doubleU = F/3.6/dt;//折算系數
doubleDC;//肯定性系數
void?inputQ();//讀入原始數據
double?FuntionDC(doubleCanshu[]);//計算肯定性系數
void?outputQ(doubleCanshu[]);//輸出模擬流量
//**********讀入原始數據(函數)************
void?inputQ()
{
usingnamespace?std;
ifstream infile;//讀人三個雨量站實測流域流量
infile.open("infile_3P_Qo.txt");
for(int?i = 0; i < M; i++)
infile>>P1[i]>>P2[i]>>P3[i]>>Qo[i];
SumQo = 0, AveQo;
for?(int?i = 0; i < M; i++)
SumQo += Qo[i];
AveQo = SumQo/M;
S2 = 0;
for?(int?i = 0; i < M; i++)
S2 += pow(Qo[i] -AveQo, 2);
infile.close();
infile.open("infile_q.txt");
for?(int?i = 0; i < Mq; i++)
{//讀人時段單位線數據
infile>>q[i];
}
infile.close();
}
//**********計算肯定性系數(函數)************
doubleFuntionDC(double?Canshu[])
{
usingnamespace?std;
K = Canshu[10]; WM =Canshu[0]; WUM = Canshu[1]; WLM = Canshu[2]; C = Canshu[3];
b = Canshu[4];
SM = Canshu[5]; EX = Canshu[6];KI = Canshu[7];
KKI = Canshu[8]; KKG? = Canshu[9];
//******三層蒸發模式下的蓄滿產流模型(開始)
double?WDM =WM - WUM - WLM, KG = 0.7 - KI;
double?WMM =WM*(1 + b);
WU[0] = FE*WUM;
WL[0] = FE*WLM;
WD[0] = FE*WDM;
//********計算蒸發能力
intSumTime1 = 0, SumTime2 = 0;
for(int?j = 0; j < 4; j++)
{
SumTime1 =? SumTime2;
SumTime2+=4*Days[j];
if(SumTime2 > M) SumTime2 = M;
for(int?i = SumTime1; i < SumTime2; i++)
{
P[i] = (P1[i]+P2[i]+P3[i])/3;
if(P[i] < 3)??Ep[i] = AveE[0][j] * K;
else?? Ep[i] =AveE[1][j] * K;
}
}
for?(int?i = 0; i < M; i++)
{
W[i] = WU[i] + WL[i]+ WD[i];
if((1 -W[i]/WM)<0)
a[i] = WMM;
else?a[i] =WMM*(1 - pow(1 - W[i]/WM,1/(b + 1)));
if?( P[i] ==0)
{
if?(WU[i]
{
EU[i] = WU[i];
if?(WL[i]>= C*WLM)
{
EL[i] = (Ep[i] - EU[i])*WL[i]/WLM;
ED[i] = 0;
}
else
if?(WL[i]
{
EL[i] = WL[i];
ED[i] = (Ep[i] - EU[i])*C - EL[i];
}
else
{
EL[i] = (Ep[i] - EU[i])*C;
ED[i] =?0;
}
E[i] = EU[i] + EL[i] + ED[i];
WU[i + 1] = 0;
WL[i + 1] = WL[i] - EL[i];
WD[i + 1] = WD[i] - ED[i];
}
else
{
EL[i] = ED[i] = 0;
E[i] = EU[i] = Ep[i];
WU[i + 1] = WU[i] - Ep[i];
WL[i + 1] = WL[i];
WD[i + 1] = WD[i];
}
PE[i] = -E[i];???//PE為負值
R[i] = 0;
}
else
{//3
if?(P[i] +WU[i] < Ep[i])
{
EU[i] =P[i]?+ WU[i];
if?(WL[i] >=C*WLM)
{
EL[i] = (Ep[i] - EU[i])*WL[i]/WLM;
ED[i] = 0;
}
else
if?(WL[i]
{
EL[i] = WL[i];
ED[i] = (Ep[i] - EU[i])*C - EL[i];
}
else
{
EL[i] = (Ep[i] - EU[i])*C;
ED[i] =?0;
}
E[i] = EU[i] + EL[i] + ED[i];
PE[i] = P[i] - E[i];????//PE為負值
R[i] = 0;
WU[i + 1] = 0;
WL[i + 1] = WL[i] - EL[i];
WD[i + 1] = WD[i] - ED[i];
}
else
{
EL[i] = ED[i] = 0;
E[i] = EU[i] = Ep[i];
PE[i] = P[i] - E[i];
if?(P[i]
{
R[i] = 0;
WU[i + 1] = WU[i] +? P[i] - E[i];
WL[i + 1] = WL[i];
WD[i+ 1] = WD[i];
}//到此,因PE為負而無產流的狀況所有討論完畢
else
{//如下狀況出現產流，注意此時蒸發只發生在WU，即EL=ED=0
if?(a[i] +PE[i] <= WMM)
{
R[i] = PE[i] + W[i]- WM? + WM*(pow(1 - (a[i] + PE[i])/WMM,b + 1));
WU[i + 1] =?PE[i] + WU[i]- R[i];
WL[i + 1] = WL[i ];
WD[i + 1] = WD[i];
if?(WU[i + 1]> WUM)
{
WL[i+ 1] = WU[i + 1] - WUM + WL[i];
WU[i+ 1] =? WUM;
if?(WL[i + 1] > WLM)
{
WD[i+ 1] = WL[i + 1] - WLM + WD[i];
WL[i+ 1] = WLM;
if?(WD[i + 1] > WDM) WD[i + 1] = WDM;
}
}
}
else
{
R[i] = PE[i] + W[i]- WM;
WU[i + 1] = WUM;
WL[i + 1] = WLM;
WD[i + 1] = WDM;
}
}
}
}
if((a[i]+PE[i])>WMM)
af[i] = 1;
else
af[i] = 1 -pow(1 - (a[i]+PE[i])/WMM,b);
}
//**********三水源劃分(開始)
double?SMM =SM*(1 + EX);
double?S0 =FE*SM, af0 = 1 - pow(1 - a[0]/WMM,b);
//af0指W[0]對應的面積比率(%)
for?(int?i = 0; i < M; i++)
{
if(i == 0)???? S[i] = S0*af0/af[i];
else?S[i] =S[i - 1]*af[i - 1]/af[i]+( R[i - 1]- RS[i - 1] - RI[i - 1] - RG[i - 1])/af[i];
if(S[i]>SM)S[i] = SM;
AU[i] = SMM*(1 - pow(1 - S[i]/SM,1/(EX + 1)));
if(R[i] == 0)
RS[i]? =0;
else
{
if(AU[i] +PE[i] > SMM)
RS[i] =(S[i] + PE[i] - SM)*af[i];
else
RS[i] =(S[i] + PE[i] - SM + SM*(pow(1 - (AU[i] + PE[i])/SMM,EX + 1)))*af[i];
}
RI[i] = KI*(S[i] + (R[i] -RS[i])/af[i])*af[i];
RG[i]= KG/KI*RI[i];
}
//**********匯流(開始)
for(int? j = 0; j< M; j++)
for(int?i = 0; i < Mq + M - 1; i++)
if(i
else
{
if(i < j +Mq)? Qs[i][j] = RS[j]/10*q[i - j];
else?Qs[i][j]= 0;
}
for(int? i = 0; i< Mq + M - 1; i++)
{
QS[i] = 0;//必定要初始化
for(int?j = 0; j < M ; j++)
QS[i] += Qs[i][j];
}
QI[0] = RI[0]*(1 - KKI)*U;
QG[0] = RG[0]*(1 - KKG)*U;
for(int?i = 1; i < Mq + M - 1; i++)
if(i < M )
{
QI[i] = RI[i]*(1 - KKI)*U + QI[i - 1]*KKI;
QG[i] = RG[i]*(1 - KKG)*U + QG[i - 1]*KKG;
}
else
{
QI[i] = QI[i - 1]*KKI;
QG[i] = QG[i - 1]*KKG;
}
for(int? i = 0; i< Mq + M - 1; i++)
Q[i] = QS[i] + QI[i] + QG[i];
//*****肯定性系數計算(開始)
S1 = 0;
for?(int?i = 0; i < M; i++)
S1 += pow(Q[i] - Qo[i],2);
DC = 1 - S1/S2;
returnDC;
}
//**********輸出模擬流量過程(函數)************
voidoutputQ(double?Canshu[])
{
usingnamespace?std;
ofstream outfile;
outfile.open("outfile_Q.txt");
outfile<
<
for?(int?i = 0; i < M; i++)
outfile<
outfile.close();
}
4.3.2基因遺傳算法
//遺傳算法.h
#include?
#include?"新安江三水源模型.h"
const?intGenerationNum = 200,//最大演算世代次數
SumNum = 60,//當最優適應度重復次數超過此值時中止演算
IndividualNum =21,//該種群的個體數目
ChromosomeNum =11;//每一個個體的基因(待率定參數)數目
const?double?ChrTop[ChromosomeNum]//基因(待率定參數)的閾值上限
={200,20, 90, 0.20, 0.4, 25, 1.5, 0.7, 1.0, 1.0, 1.5},
ChrBottom[ChromosomeNum]//基因(待率定參數)的閾值下限
={120,10, 60, 0.09, 0.1,? 5,? 1.0, 0.0,?0.0, 0.0, 0.5},
Pc = 0.5,//個體的交叉率(crossoverrate)
PcChr = 0.7,//交叉對交叉的基因比率
PmInd = 0.7,//個體變異率(mutationrate)
PmChr = 0.5,//變異個體的基因變異率(mutationrate)
Bm = 4;//變異運算的形狀系數
intnc =? (int)((IndividualNum - 1)*Pc/2),//nc對個體的基因參與交叉
ncChr = (int) (ChromosomeNum*PcChr),//兩個體交叉的基因數
nmInd = (int) ((IndividualNum - 1)*PmInd),//nmInd個個體發生變異
nmChr = (int) (ChromosomeNum*PmChr),//個體的nmChr個基因發生變異
x, y,tx1,tx2, Best,Worst,//挑出最優及最差的個體
CountNum= 1;//計數最優適應度重復次數
doubleIndChr[IndividualNum][ChromosomeNum],//每代種群的全部個體的基因
Fitness[IndividualNum],//個體的適應度
BestFitness = 0,//備份最優個體的適應度
BestIndChr[ChromosomeNum],//備份最優個體的基因
SumFitness = 0,//累積適應度
SumPs[IndividualNum] ={0},?//累積選擇幾率
dSumPs = 0,//用來求累積選擇幾率的
r = 0,//偽隨機數,交叉變異時使用
rs[IndividualNum] ={ 0},//偽隨機數,選擇時使用
temp;//中間變量
void?YiChuanSuanFa()
{
usingnamespace?std;
ofstream outfile,outtext;
outfile.open("outfile_BestIndividual.txt");//寫出文件，用于繪制遺傳算法的進化過程
outfile<
<
<
<
<
//**********初始化
srand( (unsigned)time(NULL ) );
for(int?i=0; i
for?(int?j=0; j
IndChr[i][j] = (double)rand()/RAND_MAX*(ChrTop[j]-ChrBottom[j])+ChrBottom[j];
//**********世代更替演算(開始)
for(int?g = 1;?g< GenerationNum; g++)
{
//**********適應度(開始)
for(int?i = 0; i
Fitness[i]=? FuntionDC(IndChr[i]);
if(g == 1)?//計算初始化的最后一個個體的適應度
Fitness[IndividualNum- 1] =? FuntionDC(IndChr[IndividualNum -1]);
for(tx1 = 0; tx1 < IndividualNum; tx1++)
{//找出最優個體
Best = tx1;
for( tx2 = tx1+1; tx2< IndividualNum; tx2++)
if(Fitness[tx1]< Fitness[tx2])
{
Best= tx2;
break;
}
if(tx2 ==IndividualNum)?break;
}
for(tx1 = 0; tx1 < IndividualNum; tx1++)
{//找出最差個體
Worst = tx1;
for( tx2 = tx1+1; tx2< IndividualNum; tx2++)
if(Fitness[tx1]> Fitness[tx2])
{
Worst= tx2;
break;
}
if(tx2 ==IndividualNum)?break;
}
for?(int?k=0; k
{//將最優個體排至最后,
outfile<
temp =IndChr[IndividualNum - 1][k] ;
IndChr[IndividualNum- 1][k] = IndChr[Best][k];
IndChr[Best][k]= temp;
}//最優個體不參加選擇、交叉
outfile<
temp =? Fitness[IndividualNum - 1];
Fitness[IndividualNum -1] = Fitness[Best];
Fitness[Best] = temp;
SumFitness = 0;//初始化
for(int?i = 0; i < IndividualNum - 1; i++)
{
rs[i] = (double)rand()/RAND_MAX;
SumFitness +=Fitness[i];
}
outfile<
if(BestFitness == Fitness[Best])//進化停滯
{
if((++CountNum) >= SumNum)break;
}
else?//進化成功
{
BestFitness =Fitness[Best];
CountNum = 1;
}
//**********選擇(輪盤開始)
dSumPs = 0;
for(int?i = 0; i
{
SumPs[i] =Fitness[i]/SumFitness;
dSumPs +=SumPs[i];
SumPs[i] =dSumPs;
for(int?j = i+1; j< IndividualNum - 1; j++)
{//按升序排列隨機數
if(rs[j]
{
temp= rs[i];
rs[i]= rs[j];
rs[j]= temp;
}
}
}
for(int?i = 0; i < IndividualNum - 1; i++)
{//最優個體不參加選擇
for(int?j = 0 ; j< IndividualNum - 1; j++)
if?(SumPs[j] > rs[i])
{
for?(int?k=0;k
IndChr[i][k]= IndChr[j][k];
break;
}
}
//**********交叉(開始)************
for(int?i = 0; i < nc; i++)
{//隨機交叉個體對
x = y = 0;//最優個體不參加交叉
while(x == y)
{//+0.5四舍五入
x = (int)((double)rand()/RAND_MAX*(IndividualNum - 2) + 0.5);
y = (int)((double)rand()/RAND_MAX*(IndividualNum- 2) + 0.5);
}
for(int?j = 0; j
{//隨機交叉基因對
r = (double)rand()/RAND_MAX;
int?k = (int)((double)rand()/RAND_MAX*(ChromosomeNum - 1) + 0.5);
temp =IndChr[x][k] ;
IndChr[x][k]= r*temp+(1-r)*IndChr[y][k];
IndChr[y][k]= r*IndChr[y][k]+(1-r)*temp;
}
}
//**********變異(開始)************
double?t = g/GenerationNum;//變異運算的進化標記
for(int?i = 0; i < nmInd; i++)//隨機變異個體
{//最優個體只能進行優化變異
x = (int)((double)rand()/RAND_MAX*(IndividualNum- 1) + 0.5);
if(x== (IndividualNum - 1))
for?(int?k=0;k
BestIndChr[k]= IndChr[x][k];//備份最優基因
for(int?j = 0; j
{
int?k = (int)((double)rand()/RAND_MAX*(ChromosomeNum - 1) + 0.5);
r = (double)rand()/RAND_MAX;
if((rand()%2) ==0)
IndChr[x][k]+= (ChrTop[k] - IndChr[x][k])*pow(r*(1-t),Bm);
else
IndChr[x][k]-= (IndChr[x][k]? -ChrBottom[k])*pow(r*(1-t),Bm);
}
if(x== (IndividualNum - 1))
{//判斷最優基因變異后是否優化了
Fitness[x]=? FuntionDC(IndChr[x]);
if(Fitness[x] < BestFitness)//變異退化了
{
for?(int?k=0;k
IndChr[x][k] = BestIndChr[k];//換回最優基因
Fitness[x]= BestFitness;
}
else
{
BestFitness= Fitness[x];//變異優化了
CountNum= 0;
}
}
}
}
outfile.close();
outputQ(IndChr[Best]);//輸出模擬流量
}
//新安江模型參數率定.cpp
#include?"stdafx.h"
#include?"遺傳算法.h"
void?main()
{
inputQ();
YiChuanSuanFa();
}
文章轉自：http://blog.csdn.net/superwen_go/article/details/7669429
《新程序員》：云原生和全面數字化實踐50位技術專家共同創作，文字、視頻、音頻交互閱讀總結
以上是生活随笔為你收集整理的新安江遗传算法c语言,基于遗传算法的新安江模型参数优化率定（四）的全部內容，希望文章能夠幫你解決所遇到的問題。
- 上一篇： c语言表现一些简单的图片,C语言的一些简
- 下一篇： android 启动服务权限,andro