subroutine srfmod(date,dtime,dx,dy,cellat,cellon,swe,snowrat, & snowfrc,snowalb,zenith,ctrns,ldrk,pwr,tempk, & temp0,press,press0,wind,depth,nspc,fsurf, & lrdlai,lai,solmas,vegmas,reemis) c use chmstry c c----CAMx v7Beta6 190902 c c----------------------------------------------------------------------- c Description: c This routine determines fraction of deposited material that is c sorbed to surface elements. Unsorbed fraction is available for c re-emission. Sorbed fraction undergoes chemistry on soil/vegetation c using information from the chemistry parameters input file. c c Copyright 1996 - 2018 c Ramboll c c Argument descriptions: c Inputs: c date I model date (YYJJJ) c dtime R time step (hours) c dx,dy R cell size (m) c cellat R latitude (deg) c cellon R longitude (deg) c swe R snow cover water equivalent (m) c snowrat R snow water accumulation rate (m/hr) c snowfrc R snow cover fraction c snowalb R snow albedo c zenith R zenith angle (deg) c ctrns R cloud transmissivity (unitless) c ldrk L darkness flag c pwr R rain water content (g/m3) c tempk R layer 1 temperature (K) c temp0 R surface temperature (K) c press R layer 1 pressure (mb) c press0 R surface pressure (mb) c wind R layer 1 wind (m/s) c depth R layer 1 thickness (m) c nspc I number of surface model species c fsurf R fraction LU c lrdlai L LAI inputs were read c lai R gridded LAI c solmas R Soil mass (mol/ha) c vegmas R Vegetation mass (mol/ha) c Outputs: c solmas R Soil mass (mol/ha) c vegmas R Vegetation mass (mol/ha) c reemis R Re-emission rate (mol/s) c c Called by: c CHEMDRIV c c Routines called: c CALDATE c MICROMET c VE_GAS c c----------------------------------------------------------------------- c LOG: c----------------------------------------------------------------------- c c 04/24/13 ---cemery--- Original development (from RTRAC) c 08/22/14 ---cemery--- Extended Zhang LU to surface model c 08/25/14 ---cemery--- Added snow effects to surface model c c----------------------------------------------------------------------- c Include files: c----------------------------------------------------------------------- c implicit none include 'camx.prm' include 'deposit.inc' c c----------------------------------------------------------------------- c Argument declarations: c----------------------------------------------------------------------- c integer date integer nspc real dtime real dx, dy real cellat, cellon real swe real snowrat real snowfrc real snowalb real zenith real ctrns real pwr real fsurf(NLU) real tempk, temp0 real press, press0 real wind real depth real lai real solmas(nspc) real vegmas(nspc) real reemis(nspc) logical ldrk logical lrdlai c c----------------------------------------------------------------------- c Local variables: c----------------------------------------------------------------------- c integer l integer isp integer m integer mbin, latbin integer isesn integer idate, iyear, month, iday integer nday(12) real fshads real fshadv real fzenith real deg2rad real zenang real totland real smkrat real smjrat real effrates real effratev real delcons real delconv real eps real rhoh2o real volrat real rain real snow real pbl,ustar,el,psih,wstar real ve real z0 real ssorb,slch real solemis,tmpsol real vegemis,tmpveg real smas,vmas real vemis(nspc) real ref_lai, lai_ref_intpl, rlai, lai_f logical lstable logical lsnow c data deg2rad /0.01745329/ data eps /1.0e-20/ data rhoh2o /1.e6/ ! water density (g/m3) data nday/31,28,31,30,31,30,31,31,30,31,30,31/ c c----------------------------------------------------------------------- c Entry point: c----------------------------------------------------------------------- c idate = date call caldate(idate) iyear = idate/10000 month = (idate - 10000*iyear)/100 if (mod(iyear,4).eq.0) nday(2)=29 iday = idate - 100*(idate/100) lsnow = .false. if (swe.ge.0.01) lsnow = .true. ! Snow effects occur when > 10 cm deep c c-----Prepare for LU-dependent surface roughness calculation c mbin = month if (cellat.lt.0.) then mbin = mod(month+6,12) if (mbin.eq.0) mbin = 12 endif latbin = 1 if (abs(cellat).gt.20.) then latbin = 2 elseif (abs(cellat).gt.35.) then latbin = 3 elseif (abs(cellat).gt.50.) then latbin = 4 elseif (abs(cellat).gt.75.) then latbin = 5 endif if ((cellat.gt.50. .and. cellat.lt.75.) .and. & (cellon.gt.-15. .and. cellon.lt.15.)) latbin = 3 isesn = iseason(latbin,mbin) c c-----Use input snow cover to set season, if specified c if (swe.ge.0.001) isesn = 4 ! Snow cover check > 1 cm c c-----Determine cell-average relative LAI c if (lrdlai .and. NLU.eq.NLUZ03) then ref_lai = 0. do m = 1,NLU lai_ref_intpl = lai_ref(m,mbin) + & float(min(nday(mbin),iday))/float(nday(mbin))* & (lai_ref(m,mbin+1) - lai_ref(m,mbin)) ref_lai = ref_lai + fsurf(m)*lai_ref_intpl enddo rlai = lai/(ref_lai + 1.e-10) endif c c-----Loop over land use; surface roughness for water is dependent on c wind speed c do l = 1,nspc vemis(l) = 0. enddo do 10 m = 1,NLU if (fsurf(m).lt.0.01) goto 10 c c-----Set surface roughness for Wesely (1989) scheme c if (NLU.eq.NLUW89) then z0 = z0lu(m,isesn) if (m.eq.7) z0 = amax1(z0,2.0e-6*wind**2.5) c c-----Set surface roughness and LAI for Zhang (2003) scheme c else lai_f = lai_ref(m,mbin) + & float(min(nday(mbin),iday))/float(nday(mbin))* & (lai_ref(m,mbin+1) - lai_ref(m,mbin)) if (lrdlai) then lai_f = lai_f*rlai lai_f = amin1(lai_ref(m,15),lai_f) lai_f = amax1(lai_ref(m,14),lai_f) endif if (m.eq.1 .or. m.eq.3) then z0 = 2.0e-6*wind**2.5 else if (z02(m).gt.z01(m)) then z0 = z01(m) + (lai_f - lai_ref(m,14))/ & (lai_ref(m,15) - lai_ref(m,14))* & (z02(m) - z01(m)) else z0 = z01(m) endif endif endif c c-----Get surface layer micrometeorological parameters for this cell and c landuse type c call micromet(tempk,temp0,press,press0,depth/2.,wind,z0,pbl, & ustar,el,psih,wstar,lstable) c c-----Loop over species, and calculate re-emission velocity for this cell, c landuse, and current species. Aggregated into cell-average re-emission c velocity. c do l = 1,nspc isp = idsmsp(l) if ((NLU.eq.NLUW89 .and. m.eq.7) .or. (NLU.eq.NLUZ03 .and. & (m.eq.1 .or. m.eq.3)))then ve = 0. else call ve_gas(z0,depth/2.,psih,ustar,diffrat(isp),ve) endif vemis(l) = vemis(l) + ve*fsurf(m) enddo 10 continue c c-----Calculate fraction of sorbed and un-sorbed material: c Un-sorbed fraction is available for re-emission c Sorbed fraction remains on surface and undergoes chemistry, leaching c and vegetative penetration c do l = 1,nspc ssorb = smssrb(l) if (lsnow) ssorb = smisrb(l) solemis = (vemis(l)/depth)*solmas(l)/(1. + ssorb) ! (mol/ha/s) vegemis = (vemis(l)/depth)*vegmas(l)/(1. + smvsrb(l)) ! (mol/ha/s) tmpsol = MAX(eps,solmas(l) - solemis*dtime*3600.) ! (mol/ha) tmpveg = MAX(eps,vegmas(l) - vegemis*dtime*3600.) ! (mol/ha) solemis = (solmas(l) - tmpsol)/(dtime*3600.) ! (mol/ha/s) vegemis = (vegmas(l) - tmpveg)/(dtime*3600.) ! (mol/ha/s) solmas(l) = tmpsol ! (mol/ha) vegmas(l) = tmpveg ! (mol/ha) reemis(l) = (solemis + vegemis)*dx*dy*1.e-4 ! (mol/s) enddo c c --- Calculate photolysis factors (shade, zenith, clouds, snow) c zenang = deg2rad*zenith fshads = 0.0 fshadv = 0.0 fzenith = 0.0 totland = 0.0 if (.not.ldrk) then fshadv = 0.5 do m = 1,NLU if (NLU.eq.NLUW89) then fshads = fshads + fsoil(m)*fsurf(m) else fshads = fshads + fsoilz(m)*fsurf(m) endif totland = totland + fsurf(m) enddo fshads = fshads/totland if (lsnow) then fshads = fshads + snowfrc*snowalb fshadv = fshadv + snowfrc*snowalb endif fzenith = COS(zenang) endif c c-----Chemical decay of precursor species to product species c do l = 1,nsmrxn isp = idsmpre(l) smkrat = smskrat(l) smjrat = smsjrat(l) if (lsnow) then smkrat = smikrat(l) smjrat = smijrat(l) endif effrates = smkrat + smjrat*fshads*fzenith*ctrns effratev = smvkrat(l) + smvjrat(l)*fshadv*fzenith*ctrns delcons = solmas(isp)*(1. - EXP(-effrates*dtime*60.)) delconv = vegmas(isp)*(1. - EXP(-effratev*dtime*60.)) smas = MAX(eps,solmas(isp) - delcons) vmas = MAX(eps,vegmas(isp) - delconv) delcons = solmas(isp) - smas delconv = vegmas(isp) - vmas solmas(isp) = smas vegmas(isp) = vmas if (idsmprd(l).ne.0) then isp = idsmprd(l) solmas(isp) = solmas(isp) + delcons vegmas(isp) = vegmas(isp) + delconv endif enddo c c --- Apply decay rates for leaching, snow melt, and plant penetration c Accelerate if raining (1 mm/hr = e-folding in 1 hr) c Accelerate if snowing (1 cm/hr = e-folding in 1 day) c rain = (pwr/rhoh2o/1.0e-7)**1.27 !rain rate (mm/hr) if (rain.lt.0.2) rain = 0. rain = rain/60. !mm/hr -> 1/min snow = 0. if (snowrat.gt.0.) & snow = 100.*10.*snowrat/(24.*60.) !m/hr SWE -> cm/hr depth -> 1/min do l = 1,nspc slch = smlch(l) if (lsnow) then slch = 0. if (temp0.gt.272.) slch = smmlt(l) endif delcons = EXP(-(slch + rain + snow)*dtime*60.) delconv = EXP(-(smpen(l) + rain + snow)*dtime*60.) solmas(l) = MAX(eps,solmas(l)*delcons) vegmas(l) = MAX(eps,vegmas(l)*delconv) enddo c c----------------------------------------------------------------------- c Return point: c----------------------------------------------------------------------- c return end