|
| 1 | + subroutine orgncswat2(iwave) |
| 2 | + |
| 3 | +!! ~ ~ ~ PURPOSE ~ ~ ~ |
| 4 | +!! this subroutine calculates the amount of organic nitrogen removed in |
| 5 | +!! surface runoff - when using CSWAT==2 it |
| 6 | + |
| 7 | + |
| 8 | +!! ~ ~ ~ INCOMING VARIABLES ~ ~ ~ |
| 9 | +!! name |units |definition |
| 10 | +!! ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ |
| 11 | +!! da_ha |ha |area of watershed in hectares |
| 12 | +!! enratio |none |enrichment ratio calculated for day in HRU |
| 13 | +!! erorgn(:) |none |organic N enrichment ratio, if left blank |
| 14 | +!! |the model will calculate for every event |
| 15 | +!! ihru |none |HRU number |
| 16 | +!! iwave |none |flag to differentiate calculation of HRU and |
| 17 | +!! |subbasin sediment calculation |
| 18 | +!! |iwave = 0 for HRU |
| 19 | +!! |iwave = subbasin # for subbasin |
| 20 | +!! sedyld(:) |metric tons |daily soil loss caused by water erosion in |
| 21 | +!! |HRU |
| 22 | +!! sol_bd(:,:) |Mg/m**3 |bulk density of the soil |
| 23 | +!! sol_z(:,:) |mm |depth to bottom of soil layer |
| 24 | +!! sub_bd(:) |Mg/m^3 |bulk density in subbasin first soil layer |
| 25 | +!! sub_fr(:) |none |fraction of watershed area in subbasin |
| 26 | +!! sub_orgn(:) |kg N/ha |amount of nitrogen stored in all organic |
| 27 | +!! sedc_d(:) |kg C/ha |amount of C lost with sediment |
| 28 | +!! |
| 29 | +!! |
| 30 | +!! |
| 31 | +!! |pools |
| 32 | +!! ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ |
| 33 | + |
| 34 | +!! ~ ~ ~ OUTGOING VARIABLES ~ ~ ~ |
| 35 | +!! name |units |definition |
| 36 | +!! ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ |
| 37 | +!! sedorgn(:) |kg N/ha |amount of organic nitrogen in surface runoff |
| 38 | +!! |in HRU for the day |
| 39 | +!! ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ |
| 40 | + |
| 41 | +!! ~ ~ ~ LOCAL DEFINITIONS ~ ~ ~ |
| 42 | +!! name |units |definition |
| 43 | +!! ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ |
| 44 | +!! conc | |concentration of organic N in soil |
| 45 | +!! er |none |enrichment ratio |
| 46 | +!! j |none |HRU number |
| 47 | +!! wt1 |none |conversion factor (mg/kg => kg/ha) |
| 48 | +!! xx |kg N/ha |amount of organic N in first soil layer |
| 49 | +!! ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ |
| 50 | + |
| 51 | +!! ~ ~ ~ ~ ~ ~ END SPECIFICATIONS ~ ~ ~ ~ ~ ~ |
| 52 | + |
| 53 | + use parm |
| 54 | + |
| 55 | + integer, intent (in) :: iwave |
| 56 | + integer :: j |
| 57 | + real :: xx, wt1, er, conc |
| 58 | + real :: sol_mass, QBC, VBC, YBC, YOC, YW, TOT, YEW, X1, PRMT_21, PRMT_44 |
| 59 | + real :: DK, V, X3, CO, CS, perc_clyr, latc_clyr |
| 60 | + integer :: k |
| 61 | + latc_clyr = 0. |
| 62 | + |
| 63 | + j = 0 |
| 64 | + j = ihru |
| 65 | + |
| 66 | + !!for debug purpose by zhang |
| 67 | + !if (iyr == 1991 .and. i==235) then |
| 68 | + ! write(*,*) 'stop' |
| 69 | + !end if |
| 70 | + |
| 71 | + xx = 0. |
| 72 | + wt1 = 0. !! conversion factor |
| 73 | + er = 0. !! enrichment ratio |
| 74 | + if (iwave <= 0) then |
| 75 | + !! HRU calculations |
| 76 | + !xx = sol_n(1,j) + sol_fon(1,j) + sol_mn(1,j) |
| 77 | + xx = sol_LSN(1,j)+sol_LMN(1,j)+sol_HPN(1,j)+sol_HSN(1,j) !+sol_BMN(1,j) |
| 78 | + !wt = sol_bd(1,j) * sol_z(1,j) * 10. (tons/ha) |
| 79 | + !wt1 = wt/1000 |
| 80 | + wt1 = sol_bd(1,j) * sol_z(1,j) / 100. |
| 81 | + |
| 82 | + if (erorgn(j) > .001) then |
| 83 | + er = erorgn(j) |
| 84 | + else |
| 85 | + er = enratio |
| 86 | + end if |
| 87 | + |
| 88 | + else |
| 89 | + !! subbasin calculations |
| 90 | + xx = sub_orgn(iwave) |
| 91 | + wt1 = sub_bd(iwave) * sol_z(1,j) / 100. |
| 92 | + |
| 93 | + er = enratio |
| 94 | + end if |
| 95 | + |
| 96 | + conc = 0. |
| 97 | + conc = xx * er / wt1 |
| 98 | + |
| 99 | + if (iwave <= 0) then |
| 100 | + !! HRU calculations |
| 101 | + sedorgn(j) = .001 * conc * sedyld(j) / hru_ha(j) |
| 102 | + else |
| 103 | + !! subbasin calculations |
| 104 | + sedorgn(j) = .001 * conc * sedyld(j) / (da_ha * sub_fr(iwave)) |
| 105 | + end if |
| 106 | + |
| 107 | + !! update soil nitrogen pools only for HRU calculations |
| 108 | + if (iwave <= 0 .and. xx > 1.e-6) then |
| 109 | + xx1 = (1. - sedorgn(j) / xx) |
| 110 | + |
| 111 | + !!add by zhang to update soil nitrogen pools |
| 112 | + |
| 113 | + sol_LSN(1,j) = sol_LSN(1,j) * xx1 |
| 114 | + sol_LMN(1,j) = sol_LMN(1,j) * xx1 |
| 115 | + sol_HPN(1,j) = sol_HPN(1,j) * xx1 |
| 116 | + sol_HSN(1,j) = sol_HSN(1,j) * xx1 |
| 117 | + !sol_BMN(1,j) = sol_BMN(1,j) * xx1 |
| 118 | + end if |
| 119 | + |
| 120 | + !return |
| 121 | + |
| 122 | + !Calculate runoff and leached C&N from micro-biomass |
| 123 | + latc_clyr = 0. |
| 124 | + sol_mass = 0. |
| 125 | + !kg/ha |
| 126 | + sol_mass = (sol_z(1,j) / 1000.) * 10000. * sol_bd(1,j)* 1000. * (1- sol_rock(1,j) / 100.) |
| 127 | + |
| 128 | + |
| 129 | + QBC=0. !c loss with runoff or lateral flow |
| 130 | + VBC=0. !c los with vertical flow |
| 131 | + YBC=0. !BMC LOSS WITH SEDIMENT |
| 132 | + YOC=0. !Organic C loss with sediment |
| 133 | + YW=0. !YW = WIND EROSION (T/HA) |
| 134 | + TOT=sol_HPC(1,j)+sol_HSC(1,j)+sol_LMC(1,j)+sol_LSC(1,j) !Total organic carbon in layer 1 |
| 135 | + !YEW = MIN(er*(sedyld(j)/hru_ha(j)+YW/hru_ha(j))/(sol_mass/1000.),.9) |
| 136 | + ! Not sure whether should consider enrichment ratio or not! |
| 137 | + YEW = MIN((sedyld(j)/hru_ha(j)+YW/hru_ha(j))/(sol_mass/1000.),.9) !fraction of soil erosion of total soil mass |
| 138 | + X1=1.-YEW |
| 139 | + !YEW=MIN(ER*(YSD(NDRV)+YW)/WT(LD1),.9) |
| 140 | + !ER enrichment ratio |
| 141 | + !YSD water erosion |
| 142 | + !YW wind erosion |
| 143 | + YOC=YEW*TOT |
| 144 | + sol_HSC(1,j)=sol_HSC(1,j)*X1 |
| 145 | + sol_HPC(1,j)=sol_HPC(1,j)*X1 |
| 146 | + sol_LS(1,j)=sol_LS(1,j)*X1 |
| 147 | + sol_LM(1,j)=sol_LM(1,j)*X1 |
| 148 | + sol_LSL(1,j)=sol_LSL(1,j)*X1 |
| 149 | + sol_LSC(1,j)=sol_LSC(1,j)*X1 |
| 150 | + sol_LMC(1,j)=sol_LMC(1,j)*X1 |
| 151 | + sol_LSLC(1,j)=sol_LSLC(1,j)*X1 |
| 152 | + sol_LSLNC(1,j)=sol_LSC(1,j)-sol_LSLC(1,j) |
| 153 | + if (surfq(j) > 0) then |
| 154 | + !write(*,*) 'stop' |
| 155 | + end if |
| 156 | + IF(sol_BMC(1,j)>.01) THEN |
| 157 | + PRMT_21 = 0. !KOC FOR CARBON LOSS IN WATER AND SEDIMENT(500._1500.) KD = KOC * C |
| 158 | + PRMT_21 = 1000. |
| 159 | + sol_WOC(1,j) = sol_LSC(1,j)+sol_LMC(1,j)+sol_HPC(1,j)+sol_HSC(1,j)+sol_BMC(1,j) |
| 160 | + DK=.0001*PRMT_21*sol_WOC(1,j) |
| 161 | + !X1=PO(LD1)-S15(LD1) |
| 162 | + X1 = sol_por(1,j)*sol_z(1,j)-sol_wpmm(1,j) !mm |
| 163 | + IF (X1 <= 0.) THEN |
| 164 | + X1 = 0.01 |
| 165 | + END IF |
| 166 | + XX=X1+DK |
| 167 | + !V=QD+Y4 |
| 168 | + V = surfq(j) + sol_prk(1,j) + flat(1,j) |
| 169 | + !QD surface runoff |
| 170 | + X3=0. |
| 171 | + IF(V>1.E-10)THEN |
| 172 | + X3=sol_BMC(1,j)*(1.-EXP(-V/XX)) !loss of biomass C |
| 173 | + PRMT_44 = 0. !RATIO OF SOLUBLE C CONCENTRATION IN RUNOFF TO PERCOLATE(0.1_1.) |
| 174 | + PRMT_44 = .5 |
| 175 | + CO=X3/(sol_prk(1,j) + PRMT_44*(surfq(j)+flat(1,j))) !CS is the horizontal concentration |
| 176 | + CS=PRMT_44*CO !CO is the vertical concentration |
| 177 | + VBC=CO*(sol_prk(1,j)) !!! sol_prk(:,:) |mm H2O |percolation from soil layer on current day |
| 178 | + sol_BMC(1,j)=sol_BMC(1,j)-X3 |
| 179 | + QBC=CS*(surfq(j)+flat(1,j)) |
| 180 | + ! COMPUTE WBMC LOSS WITH SEDIMENT |
| 181 | + IF(YEW>0.)THEN |
| 182 | + CS=DK*sol_BMC(1,j)/XX |
| 183 | + YBC=YEW*CS |
| 184 | + END IF |
| 185 | + END IF |
| 186 | + END IF |
| 187 | + |
| 188 | + sol_BMC(1,j)=sol_BMC(1,j)-YBC |
| 189 | + surfqc_d(j) = QBC*(surfq(j)/(surfq(j)+flat(1,j)+1.e-6)) |
| 190 | + |
| 191 | + sol_latc(1,j) = QBC*(flat(1,j)/(surfq(j)+flat(1,j)+1.e-6)) |
| 192 | + sol_percc(1,j) = VBC |
| 193 | + sedc_d(j) = YOC + YBC |
| 194 | + |
| 195 | + latc_clyr = latc_clyr + sol_latc(1,j) |
| 196 | + DO k=2,sol_nly(j) |
| 197 | + if (sol_prk(k,j) > 0 .and. k == sol_nly(j)) then |
| 198 | + !write (*,*) 'stop' |
| 199 | + end if |
| 200 | + sol_thick = 0. |
| 201 | + sol_thick = sol_z(k,j)-sol_z(k-1,j) |
| 202 | + sol_WOC(k,j) = sol_LSC(k,j)+sol_LMC(k,j)+sol_HPC(k,j)+sol_HSC(k,j) |
| 203 | + Y1=sol_BMC(k,j)+VBC |
| 204 | + VBC=0. |
| 205 | + IF(Y1>=.01)THEN |
| 206 | + V=sol_prk(k,j) + flat(k,j) |
| 207 | + IF(V>0.)VBC=Y1*(1.-EXP(-V/(sol_por(k,j)*sol_thick-sol_wpmm(k,j)+.0001*PRMT_21*sol_WOC(k,j)))) |
| 208 | + END IF |
| 209 | + sol_latc(k,j) = VBC*(flat(k,j)/(sol_prk(k,j) + flat(k,j)+1.e-6)) |
| 210 | + sol_percc(k,j) = VBC-sol_latc(k,j) |
| 211 | + sol_BMC(k,j)=Y1-VBC |
| 212 | + |
| 213 | + !! calculate nitrate in percolate |
| 214 | + !perc_clyr = 0. |
| 215 | + perc_clyr = perc_clyr + sol_percc(k,j) |
| 216 | + |
| 217 | + latc_clyr = latc_clyr + sol_latc(k,j) |
| 218 | + END DO |
| 219 | + |
| 220 | + latc_d(j) = latc_clyr |
| 221 | + percc_d(j) = perc_clyr |
| 222 | + |
| 223 | + |
| 224 | + return |
| 225 | + end |
0 commit comments