SUBROUTINE XAJMX( FILELEN ,& N ,& ! 单元出流数组大小 //输入变量 M ,& ! 无因次单位数组大小 //输入变量 PAR ,& ! //输入变量 AREA ,& ! 单元面积 //输入变量 UH ,& ! 无因次单位线 //输入变量 DT ,& ! 时段步长 //输入变量 P ,& ! 降雨系列 //输入变量 EP ,& ! 蒸发皿蒸发能力 //输入变量 W ,& ! 土壤含水层 1.上层 2.下层 3.深层 //输入变量 FR ,& ! 初始产流面积 //输入变量 S ,& ! 初始自由水深 //输入变量 QRSS0 ,& ! 初始壤中流流量 //输入变量 QRG0 ,& ! 初始地下水径流量 //输入变量 QR ) ! 单元出流 //输出变量 IMPLICIT NONE !///////////////////////////////////////变量声明////////////////////////////////////////////// INTEGER::N ! INTEGER::M ! REAL::PAR(13) ! 1.上层张力水容量wum 2.下层张力水容量wl 3.深层张力水容量wdm ! 4.蒸发能力折算系数KC.深层蒸发系数c 6.张力水蓄水容量系数b ! 7.不透水面积比率imp1 8.自由水蓄水容量sm 9.自由水蓄水容量指数ex !10.地下水出流系数kg 11.壤中流出流系数kss 12.地下水出流系数kkg !13.壤中流出流系数kkss REAL::AREA ! 单元面积 REAL::UH(M) ! 无因次单位线 REAL::DT ! 时段步长 REAL::P(N) ! 降雨系列 REAL::EP(N) ! 蒸发皿蒸发能力 REAL::QR(N) ! 单元出流 REAL::W(3) ! 土壤含水层 1.上层 2.下层 3.深层 REAL::FR ! 初始产流面积 REAL::S ! 初始自由水深 REAL::QRSS0 ! 初始壤中流流量 REAL::QRG0 ! 初始地下水径流量 INTEGER::D REAL::KSSD REAL::KGD REAL::E(3) REAL::WM(3) REAL::KC ! 蒸发能力折算系数 REAL::C ! 深层蒸发系数 REAL::B ! 张力水蓄水容量系数 REAL::IMP1 ! 不透水面积比率 REAL::SM ! 自由水蓄水容量 REAL::EX ! 自由水蓄水容量指数 REAL::KG ! 地下水出流系数 REAL::KSS ! 壤中流出流系数 REAL::KKG ! 地下水出流系数 REAL::KKSS ! 壤中流出流系数 !以下变量是原vb程序中未声明的变量 INTEGER::I ! 计数器 //临时变量 INTEGER::J ! 计数器 //临时变量 INTEGER::ICHECK ! 判断计算时段长度是否合适的识别码 INTEGER::NN REAL(KIND=8),PARAMETER::C5=5.000000000 REAL::U REAL::CI REAL::CG REAL::WM0 REAL::W0 REAL::PE REAL::R REAL::RIMP REAL::WMM REAL::A REAL::X REAL::RS REAL::RSS REAL::RGD REAL::RG REAL::SS REAL::Q REAL::KSSDD REAL::KGDD REAL::SMM REAL::SMMF REAL::SMF REAL::AU REAL::RSD REAL::RSSD REAL::QRS REAL::QRSS REAL::QRG REAL::QTR INTEGER:: FILELEN REAL::WTEMP(3) REAL::FRTEMP REAL::STEMP REAL::QRSS0TEMP REAL::QRG0TEMP !///////////////////////////////////////变量声明////////////////////////////////////////////// DO I = 1,3 WTEMP(I) = W(I) END DO FRTEMP = FR STEMP =S QRSS0TEMP=QRSS0 QRG0TEMP=QRG0 !///////////////////////////////////////计算区域////////////////////////////////////////////// D = 0 KSSD = 0.0 KGD = 0.0 do i = 1,3 E(I)=0.0 WM(I) =0.0 END DO KC = 0.0 C = 0.0 B = 0.0 IMP1 = 0.0 SM = 0.0 EX = 0.0 KG = 0.0 KSS = 0.0 KKG = 0.0 KKSS = 0.0 U = 0.0 CI = 0.0 CG = 0.0 WM0 = 0.0 W0 = 0.0 PE = 0.0 R = 0.0 RIMP = 0.0 WMM = 0.0 A = 0.0 X = 0.0 RS = 0.0 RSS = 0.0 RGD = 0.0 RG = 0.0 SS = 0.0 Q = 0.0 KSSDD = 0.0 KGDD = 0.0 SMM = 0.0 SMMF = 0.0 SMF = 0.0 AU = 0.0 RSD = 0.0 RSSD = 0.0 QRS = 0.0 QRSS = 0.0 QRG = 0.0 QTR = 0.0 ! 赋值 ICHECK = 1 DO I = 1, 3 WM(I) = PAR(I) END DO KC = PAR(4) C = PAR(5) B = PAR(6) IMP1 = PAR(7) SM = PAR(8) EX = PAR(9) KG = PAR(10) KSS = PAR(11) KKG = PAR(12) KKSS = PAR(13) DO I = 1, N QR(I) = 0.0 END DO U = AREA / 3.6 / DT IF( DT.LE.24.0 )THEN D = 24 / DT CI = KKSS ** (1.0 / D) CG = KKG ** (1.0 / D) KSSD = (1.0 - (1.0 - (KG + KSS)) ** (1.0 / D)) / (1.0 + & KG / KSS) KGD = KSSD * KG / KSS ELSE ICHECK = 0 END IF DO I = 1 , N IF (EP(I) .LT. 0.0)THEN EP(I) = 0.0 END IF IF (P(I) .LT. 0.0)THEN P(I) = 0.0 END IF EP(I) = EP(I) * KC WM0 = WM(1) + WM(2) + WM(3) W0 = W(1) + W(2) + W(3) PE = P(I) - EP(I) ! 赋初值 R = 0.0 RIMP = 0.0 IF(PE.GT.0.0)THEN WMM = (1.0 + B) * WM0 / (1.0 - IMP1) IF ((WM0 - W0).LE.0.0001)THEN A = WMM ELSE A = WMM * (1.0 - (1.0 - W0 / WM0) ** (1.0 / (1.0 + B))) END IF IF ((PE + A) .LT. WMM) THEN R = PE - WM0 + W0 + WM0 * ((1.0 - (PE + A) / WMM) & ** (1.0 + B)) ELSE R = PE - (WM0 - W0) END IF RIMP = PE * IMP1 END IF IF ((W(1) + P(I)) .GT. EP(I)) THEN E(1) = EP(I) E(2) = 0.0 E(3) = 0.0 ELSE E(1) = W(1) + P(I) E(2) = (EP(I) - E(1)) * W(2) / WM(2) IF (W(2) .LE.( C * WM(2))) THEN E(2) = C * (EP(I) - E(1)) E(3) = 0.0 IF (W(2) .GE. (C * (EP(I) - E(1)))) THEN E(2) = C * (EP(I) - E(1)) E(3) = 0.0 ELSE E(2) = W(2) E(3) = C * (EP(I) - E(1) - E(2)) END IF END IF END IF W(1) = W(1) + P(I) - R - E(1) W(2) = W(2) - E(2) W(3) = W(3) - E(3) IF (W(1) .GT. WM(1)) THEN W(2) = W(1) - WM(1) + W(2) W(1) = WM(1) IF (W(2) .GT. WM(2)) THEN W(3) = W(3) + W(2) - WM(2) W(2) = WM(2) END IF END IF X = FR IF (PE .LE.0.0) THEN RS = 0.0 RSS = S * KSSD * FR RGD = S * FR * KGD S = S - (RSS + RG) / FR ELSE FR = R / PE S = X * S / FR SS = S Q = R / FR NN = INT(Q / C5) + 1 ! 在vb程序中C5是双精度常数5 Q = Q / NN KSSDD = (1.0 - (1.0 - (KGD + KSSD)) ** (1.0 / NN)) /& (1.0 + KGD / KSSD) KGDD = KSSDD * KGD / KSSD RS = 0.0 RSS = 0.0 RG = 0.0 SMM = (1 + EX) * SM IF (EX .LT. 0.001) THEN SMMF = SMM ELSE SMMF = SMM * (1.0 - (1.0 - FR) ** (1.0 / EX)) END IF SMF = SMMF / (1.0 + EX) DO J = 1 , NN IF (S .GT. SMF) THEN S = SMF END IF AU = SMMF * (1.0 - (1.0 - S / SMF) ** (1.0 /& (1.0 + EX))) IF ((Q + AU) .LE. 0.0) THEN RSD = 0.0 RSSD = 0.0 RGD = 0.0 S = 0.0 ELSE IF ((Q + AU) .GE. SMMF) THEN RSD = (Q + S - SMF) * FR RSSD = SMF * KSSDD * FR RGD = SMF * FR * KGDD S = SMF - (KSSD + KGD) / FR ELSE IF ((Q + AU) .LT. SMMF) THEN RSD = (Q - SMF + S + SMF * (1.0 - (Q + AU)& / SMMF) ** (1.0 + EX)) * FR RSSD = (S + Q - RSD / FR) * KSSDD * FR RGD = (S + Q - RSD / FR) * KGDD * FR S = S + Q - (RSD + RSSD + RGD) / FR END IF RS = RS + RSD RSS = RSS + RSSD RG = RG + RGD END DO END IF RS = RS * (1.0 - IMP1) RSS = RSS * (1.0 - IMP1) RG = RG * (1.0 - IMP1) QRS = (RS + RIMP) * U QRSS = QRSS0 * CI + RSS * (1.0 - CI) * U QRG = QRG0 * CG + RG * (1.0 - CG) * U QTR = QRS + QRSS + QRG DO J = 1 , M IF ((I + J - 1) .LE. N) THEN QR(I + J - 1) = QR(I + J - 1) + QTR * UH(J) END IF END DO QRSS0 = QRSS QRG0 = QRG END DO DO I = 1,3 W(I) = WTEMP(I) END DO FR = FRTEMP S =STEMP QRSS0=QRSS0TEMP QRG0=QRG0TEMP END SUBROUTINE XAJMX