新安江模型程序C++代码.docx
《新安江模型程序C++代码.docx》由会员分享，可在线阅读，更多相关《新安江模型程序C++代码.docx（20页珍藏版）》请在冰豆网上搜索。
新安江模型程序C++代码
新安江模型程序C++代码
以下是类的声明：
classXinanjiangModel
{
private:
//FORCING
double*m_pP;//降水数据
double*m_pEm;//水面蒸发数据
//
longm_nSteps;//模型要运行的步长（一共m_nSteps步）
longsteps;
//OUTPUT
double*m_pR;//流域内每一步长的产流量（径流深度）
double*m_pRs;//每一步长的地表径流深（毫米）
double*m_pRi;//每一步长的壤中流深（毫米）
double*m_pRg;//每一步长的地下径流深（毫米）
double*m_pE;//每一步长的蒸发（毫米）
double*m_pQrs;//流域出口地表径流量
double*m_pQri;//流域出口壤中流径流流量
double*m_pQrg;//流域出口地下径流量
double*m_pQ;//流域出口的总流量
doublem_U;//for24h.U=A（km^2）/3.6/delta_t
//SOIL
double*m_pW;//流域内土壤湿度
double*m_pWu;//流域内上层土壤湿度
double*m_pWl;//流域内下层土壤适度
double*m_pWd;//流域内深层土壤湿度
doublem_Wum;//流域内上层土壤蓄水容量
doublem_Wlm;//流域内下层土壤蓄水容量
doublem_Wdm;//流域内深层土壤蓄水容量，WDM=WM-WUM-WLM
//EVAPORATION
double*m_pEu;//上层土壤蒸发量（毫米）
double*m_pEl;//下层土壤蒸发量（毫米）
double*m_pEd;//深层土壤蒸发量（毫米）
//runoff
double*RF;
//PARAMETER
doublem_Kc;//流域蒸散发能力与实测蒸散发值的比
doublem_IM;//不透水面积占全流域面积之比
doublem_B;//蓄水容量曲线的方次，小流域（几平方公里）B0.1左右
//中等面积（平方公里以内）.2~0.3，较大面积.3~0.4
doublem_WM;//流域平均蓄水容量（毫米）（WM=WUM+WLM+WDM）
doublem_C;//流域内深层土壤蒸发系数，江南湿润地区：
0.15-0.2，
//华北半湿润地区：
.09-0.12
doublem_SM;//自由水蓄水容量
doublem_EX;//自由水蓄水容量～面积分布曲线指数
doublem_KG;//地下水日出流系数
doublem_KI;//壤中流日出流系数
doublem_CG;//地下水消退系数
doublem_CI;//壤中流消退系数
double*m_UH;//单元流域上地面径流的单位线
doublem_WMM;//流域内最大蓄水容量
doublem_Area;//流域面积
intm_DeltaT;//每一步长的小时数
intm_PD;//给定数据，用以判断是否时行河道汇流计算
public:
XinanjiangModel（void）;
~XinanjiangModel（void）;
//初始化模型
voidInitModel（longnSteps,doubleArea,intDeltaT,intPD,char*ForcingFile）;
//设置模型参数
voidSetParameters（double*Params）;
//运行新安江模型
voidRunModel（void）;
//保存模拟结果到文件
voidSaveResults（char*FileName）;
//记录出流数据，用以作图分析
voidRunoff（char*runoff）;
private:
//进行汇流计算，将径流深度转换为流域出口的流量
voidRouting（void）;
};
以下是类的定义
#include"stdafx.h"
#include"xinanjiangmodel.h"
#include
#include
#include
usingnamespacestd;
#include"math.h"
#include"stdio.h"
#include"conio.h"
XinanjiangModel:
:
XinanjiangModel（void）
{
this->m_pP=NULL;
this->m_pEm=NULL;
this->m_pE=NULL;
this->m_pEd=NULL;
this->m_pEl=NULL;
this->m_pEu=NULL;
this->m_pW=NULL;
this->m_pWd=NULL;
this->m_pWl=NULL;
this->m_pWu=NULL;
this->m_pR=NULL;
this->m_pRg=NULL;
this->m_pRi=NULL;
this->m_pRs=NULL;
this->m_pQ=NULL;
this->m_pQrg=NULL;
this->m_pQri=NULL;
this->m_pQrs=NULL;
}
XinanjiangModel:
:
~XinanjiangModel（void）
{
delete[]this->m_pP;
delete[]this->m_pEm;
delete[]this->m_pE;
delete[]this->m_pEd;
delete[]this->m_pEl;
delete[]this->m_pEu;
delete[]this->m_pW;
delete[]this->m_pWd;
delete[]this->m_pWl;
delete[]this->m_pWu;
delete[]this->m_pR;
delete[]this->m_pRg;
delete[]this->m_pRi;
delete[]this->m_pRs;
delete[]this->m_pQ;
delete[]this->m_pQrg;
delete[]this->m_pQrs;
delete[]this->m_pQri;
}
//初始化模型
voidXinanjiangModel:
:
InitModel（longnSteps,doubleArea,intDeltaT,intPD,char*ForcingFile）
{
FILE*fp;
inti;
this->m_nSteps=nSteps;
this->steps=this->m_nSteps+18;
//驱动数据
this->m_pP=newdouble[this->steps];
this->m_pEm=newdouble[this->steps];
//模型输出，蒸散发项
this->m_pE=newdouble[this->steps];
this->m_pEd=newdouble[this->steps];
this->m_pEl=newdouble[this->steps];
this->m_pEu=newdouble[this->steps];
//模型输出，出流项，经过汇流的产流
this->m_pQrg=newdouble[this->steps];
this->m_pQrs=newdouble[this->steps];
this->m_pQri=newdouble[this->steps];
this->m_pQ=newdouble[this->steps];
//模型输出，产流项
this->m_pR=newdouble[this->steps];
this->m_pRg=newdouble[this->steps];
this->m_pRi=newdouble[this->steps];
this->m_pRs=newdouble[this->steps];
//模型状态量，土壤湿度
this->m_pW=newdouble[this->steps];
this->m_pWd=newdouble[this->steps];
this->m_pWl=newdouble[this->steps];
this->m_pWu=newdouble[this->steps];
//runoff值
this->RF=newdouble[this->steps];
for（i=0;isteps;i++）
{
//驱动数据
this->m_pP[i]=0.00;
this->m_pEm[i]=0.00;
//模型输出，蒸散发项
this->m_pE[i]=0.00;
this->m_pEd[i]=0.00;
this->m_pEl[i]=0.00;
this->m_pEu[i]=0.00;
//模型输出，出流项，经过汇流的产流
this->m_pQrg[i]=0.00;
this->m_pQrs[i]=0.00;
this->m_pQri[i]=0.00;
this->m_pQ[i]=0.00;
//模型输出，产流项
this->m_pR[i]=0.00;
this->m_pRg[i]=0.00;
this->m_pRi[i]=0.00;
this->m_pRs[i]=0.00;
//模型状态量，土壤湿度
this->m_pW[i]=0.00;
this->m_pWd[i]=0.00;
this->m_pWl[i]=0.00;
this->m_pWu[i]=0.00;
}
this->m_Area=Area;
this->m_DeltaT=DeltaT;
this->m_PD=PD;
this->m_U=this->m_Area/（3.6*this->m_DeltaT）;
//Forcing文件格式：
第一列：
降水（单位毫米）空格第二列水面蒸发（毫米）
if（（fp=fopen（ForcingFile,"r"））==NULL）
{printf（"Cannotopenforcingfile!
\n"）;return;}
for（i=0;im_nSteps;i++）
{fscanf（fp,"%lf%lf",&（this->m_pP[i]）,&（this->m_pEm[i]））;}
fclose（fp）;
}
//设置模型参数
voidXinanjiangModel:
:
SetParameters（double*Params）
{
this->m_Kc=Params[0];//
（1）流域蒸散发能力与实测水面蒸发之比
this->m_IM=Params[1];//
（2）流域不透水面积占全流域面积之比
this->m_B=Params[2];//（3）蓄水容量曲线的方次
this->m_Wum=Params[3];//（4）上层蓄水容量
this->m_Wlm=Params[4];//（5）下层蓄水容量
this->m_Wdm=Params[5];//（6）深层蓄水容量
this->m_C=Params[6];//（7）深层蒸散发系数
this->m_SM=Params[7];//（8）自由水蓄水容量
this->m_EX=Params[8];//（9）自由水蓄水容量～面积分布曲线指数
this->m_KG=Params[9];//（10）地下水日出流系数
this->m_KI=Params[10];//（11）壤中流日出流系数
this->m_CG=Params[11];//（12）地下水消退系数
this->m_CI=Params[12];//（13）壤中流消退系数
this->m_WM=this->m_Wum+this->m_Wlm+this->m_Wdm;
this->m_WMM=this->m_WM*（1.0+this->m_B）/（1.0-this->m_IM）;
}
//运行新安江模型
voidXinanjiangModel:
:
RunModel（void）
{
longi;
//模型的状态变量
doublePE;//>0时为净雨量;<0为蒸发不足量（mm）
doubleEp;//m_Kc*m_pEm[i]
doubleP;
doubleR;//产流深度，包括地表径流、壤中流和地下径流（mm）
doubleRB;//不透水面上产生的径流深度（mm）
doubleRG;//地下径流深度（mm）
doubleRI;//壤中流深度（mm）
doubleRS;//地表径流深（mm）
doubleA;//土壤湿度为W时土壤含水量折算成的径流深度（mm）
doubleE=0.0;//蒸散发（mm）
doubleEU=0.0;//上层土壤蒸散发量（mm）
doubleEL=0.0;//下层土壤蒸散发量（mm）
doubleED=0.0;//深层土壤蒸散发量（mm）
doubleS;
doubleFRo;
doubleFR;
doubleMS;
doubleAU;
doubleWU=5.0;//流域内上层土壤湿度
doubleWL=55.0;//流域内下层土壤适度
doubleWD=40.0;//流域内深层土壤湿度
doubleW=100.0;
doubleSo=5.0;
MS=m_SM*（1+m_EX）;
FRo=1-pow（（1-So/MS）,m_EX）;
for（i=0;im_nSteps;i++）
{
//——————蒸散发计算————————————//
RB=m_pP[i]*m_IM;//RB是降在不透水面的降雨量
P=m_pP[i]*（1-m_IM）;
Ep=m_Kc*m_pEm[i];
if（（WU+P）>=Ep）
{EU=Ep;EL=0;ED=0;}
elseif（（WU+P）{
EU=WU+P;
if（WL>=（m_C*m_Wlm））
{EL=（Ep-EU）*WL/m_Wlm;ED=0;}
elseif（（m_C*（Ep-EU））<=WL&&WL<（m_C*m_Wlm））
{EL=m_C*（Ep-EU）;ED=0;}
elseif（WL{EL=WL;ED=m_C*（Ep-EU）-EL;}
}
E=EU+EL+ED;
PE=P-E;
/*———————蒸散发计算结束———————————*/
//——————子流域产流量计算————————————//
if（PE<=0）
{R=0.00;W=W+PE;}
else
{
A=m_WMM*（1-pow（（1.0-W/m_WM）,1.0/（1+m_B）））;
//土壤湿度折算净雨量+降水后蒸发剩余雨量<流域内最大含水容量
if（（A+PE）m_WMM）
{
//流域内的产流深度计算
R=PE/*降水蒸发后的剩余量*/
+W/*流域内土壤湿度*/
+m_WM*pow（（1-（PE+A）/m_WMM）,（1+m_B））
-m_WM/*减去流域平均蓄水容量（m_WM:
参数）*/
+RB;/*不透水面上产生的径流*/
}
//土壤湿度折算净雨量+降水后蒸发剩余雨量<流域内最大含水容量
else
{
//流域内的产流深度计算
R=PE/*降水蒸发后的剩余量+W/*流域内土壤湿度*/
-m_WM/*减去流域平均蓄水容量*/
+RB;/*不透水面上产生的径流*/
}
}
//三层蓄水量的计算:
WU,WL,WD
if（WU+P-EU-R<=m_Wum）
{WU=WU+P-EU-R;WL=WL-EL;WD=WD–ED;}
else
{
WU=m_Wum;
if（WL-EL+（WU+P-EU-R-m_Wum）<=m_Wlm）
{
WL=WL–EL+（WU+P-EU-R-m_Wum）;
WD=WD-ED;
}
else
{
WL=m_Wlm;
if（WD-ED+WL-EL+（WU+P-EU-R-m_Wum）
-m_Wlm<=m_Wdm）
WD=WD-ED+WL-EL
+（WU+P-EU-R-m_Wum）-m_Wlm;
else
WD=m_Wdm;
}
}
W=WU+WL+WD;
////三水源划分汇流计算
if（PE>0）
{
FR=（R-RB）/PE;
AU=MS*（1-pow（（1-So*FRo/FR/m_SM）,1/（1+m_EX）））;
if（PE+AURS=FR*（PE+So*FRo/FR-m_SM+m_SM*pow（（1-（PE
+AU）/MS）,m_EX+1））;
elseif（PE+AU>=MS）
RS=FR*（PE+So*Fro/FR-m_SM）;
S=So*Fro/FR+（R–RS）/FR;
RI=m_KI*S*FR;
RG=m_KG*S*FR;
RS+=RB;
R=RS+RI+RG;
So=S*（1-m_KI-m_KG）;
FRo=FR;
}
else
{
S=So;
FR=1-pow（（1–S/MS）,m_EX）;
RI=0.00;
RG=0.00;
So=S*（1-m_KI-m_KG）;
RS=RB;
R=RS+RI+RG;
FRo=FR;
}
////三水源划分计算结束
/*以下部分是状态量：
总蒸发量、上、下和深层土壤的蒸发的保存*/
/*1*/this->m_pE[i]=E;//当前步长的蒸发（模型重要输出）
/*2*/this->m_pEu[i]=EU;//当前步长上层土壤蒸发
/*3*/this->m_pEl[i]=EL;//当前步长下层土壤蒸发
/*4*/this->m_pEd[i]=ED;//当前步长深层土壤蒸发
/*5*/this->m_pW[i]=W;//当前步长流域平均土壤含水量
/*6*/this->m_pWu[i]=WU;//当前步长流域上层土壤含水量
/*7*/this->m_pWl[i]=WL;//当前步长流域下层土壤含水量
/*8*/this->m_pWd[i]=WD;//当前步长流域深层土壤含水量
/*9*/this->m_pRg[i]=RG;//当前步长流域地下径流深度
/*10*/this->m_pRi[i]=RI;//当前步长流域壤中流深度
/*11*/this->m_pRs[i]=RS;//当前步长流域地表径流径流深度
/*12*/this->m_pR[i]=R;//当前步长的总产流径流深度
}
this->Routing（）;
}
//保存模拟结果到文件
voidXinanjiangModel:
:
SaveResults（char*FileName）
{
inti;
FILE*fp;
if（（fp=fopen（FileName,"w"））==NULL）
{
printf（"Cannotcreateoutputfile!
\n"）;
return;
}
fprintf（fp,"------------------------------------------------------------------\n"）;
fprintf（fp,"E（mm）EU（mm）EL（mm）ED（mm）W（mm）WU（mm）WL（mm）WD（m
升级会员 