subroutine prodhaze(ngrid,nlayer,nq,ptimestep, & pplev,pq,pdq,pdist_sol,mu0,declin,zdqfasthaze,zdqsfasthaze,gradflux, & flym_bot,flym_sol_bot,flym_ipm_bot,flym_sol,flym_ipm) use radinc_h, only : naerkind use comgeomfi_h implicit none !================================================================== ! Purpose ! ------- ! Fast production of haze in the atmosphere, direct accumulation to surface ! ! Authors ! ------- ! Tanguy Bertrand and Francois Forget (2014) ! !================================================================== #include "dimensions.h" #include "dimphys.h" #include "comcstfi.h" #include "surfdat.h" #include "comvert.h" #include "callkeys.h" #include "tracer.h" !----------------------------------------------------------------------- ! Arguments INTEGER ngrid, nlayer,nq LOGICAL firstcall SAVE firstcall DATA firstcall/.true./ REAL ptimestep REAL pplev(ngrid,nlayer+1) REAL,INTENT(IN) :: pq(ngrid,nlayer,nq) REAL,INTENT(IN) :: pdq(ngrid,nlayer,nq) REAL,INTENT(IN) :: pdist_sol ! distance SUN-pluto in AU REAL,INTENT(IN) :: mu0(ngrid) ! cosinus of solar incident flux REAL,INTENT(IN) :: declin ! declination REAL,INTENT(OUT) :: zdqfasthaze(ngrid,nq) ! Final tendancy ch4 REAL,INTENT(OUT) :: zdqsfasthaze(ngrid) ! Final tendancy haze surf REAL, INTENT(OUT) :: gradflux(ngrid) ! gradient flux Lyman alpha ph.m-2.s-1 REAL, INTENT(OUT) :: flym_bot(ngrid) ! flux Lyman alpha bottom ph.m-2.s-1 REAL, INTENT(OUT) :: flym_sol_bot(ngrid) ! Solar flux Lyman alpha ph.m-2.s-1 reaching the surface REAL, INTENT(OUT) :: flym_ipm_bot(ngrid) ! IPM (Interplanetary) flux Lyman alpha ph.m-2.s-1 reaching the surface REAL, INTENT(OUT) :: flym_sol(ngrid) ! Incident Solar flux Lyman alpha ph.m-2.s-1 REAL, INTENT(OUT) :: flym_ipm(ngrid) ! Incident IPM (Interplanetary) flux Lyman alpha ph.m-2.s-1 !----------------------------------------------------------------------- ! Local variables INTEGER l,ig,iq REAL zq_ch4(ngrid) REAL zq_prec(ngrid) REAL tauch4 REAL sigch4 ! aborption cross section ch4 cm-2 mol-1 REAL flym_sol_earth ! initial flux Earth ph.m-2.s-1 REAL flym_sol_pluto ! initial Lyman alpha SOLAR flux Pluto ph.m-2.s-1 REAL mu_ipm(ngrid) ! local Mean incident flux for IPM Lyman alpha photons REAL kch4 ! fraction of Carbon from ch4 directly dissociated by prec_haze REAL ncratio ! ration Nitrogen/Carbon observed in tholins REAL avogadro ! avogadro constant CHARACTER(len=10) :: tname REAL pagg,pform REAL tconv REAL longit REAL latit REAL valmin REAL valmax REAL valmin_dl REAL puis REAL dist REAL zdqprec(ngrid) REAL zdqch4(ngrid) REAL zdqconv_prec(ngrid) !----------------------------------------------------------------------- !******** INPUT avogadro = 6.022141e23 sigch4 = 1.85e-17 ! aborption cross section ch4 cm-2 mol-1 !! Initial Solar flux Lyman alpha on Earth (1AU) flym_sol_earth=5.e15 ! ph.m-2.s-1 --> 5e+11 ph.cm-2.s-1 !! Parameter of conversion precurseur to haze : here accelaeration/deceleration of haze producition kch4=1. ncratio=0.5 ! boost for haze considering nitrogen contribution ! ratio n/c : =0.25 if there is 1N for 3C pagg=0.05 ! Minimum pressure required for aggregation pform=0.006 ! Minimum pressure required for haze formation tconv=1.e7 ! (s) time needed to convert gas precursors into aerosol !******** Incident Solar flux !! Initial Solar flux Lyman alpha on Pluto flym_sol_pluto=0.25*flym_sol_earth/(pdist_sol*pdist_sol) ! ph.m-2.s-1 ! Decreasing by 12% the initial IPM flux to account for Interplanetary H absorption: ! Fig 3 Gladstone et al 2014 : Lyalpha at Pluto flym_sol_pluto=flym_sol_pluto*(1.-0.12*pdist_sol/32.91) !******** Incident IPM (Interplanetary) flux Lyman alpha ph.m-2.s-1 ! fit Fig. 4 in Randall Gladstone et al. (2015) ! -- New version : Integration over semi sphere of Gladstone data. ! Fit of results if (diurnal) then ! 1) get longitude/latitude (rad) of anti-subsolar point (max de mu0 - 180) longit=long((MAXLOC(mu0,DIM=1,MASK=mu0.GT.0.9))) ! rad IF (longit.GE.0) THEN longit=longit-pi ELSE longit=longit+pi ENDIF latit=-declin ! rad ! 2) Define input parameter for the fit valmin=48.74e10 ! minimum value of ipm flux in ph/m2/s valmax=115.014e10 ! maximal value of ipm flux in ph/m2/s valmin_dl=74.5e10 ! daylight minimum value puis=3.5 ! 3) Loop for each location and fit DO ig=1,ngrid !!! Solar flux flym_sol(ig)=flym_sol_pluto*mu0(ig) !!! IPM flux ! calculation of cosinus of incident angle for IPM flux mu_ipm(ig) = max(mu0(ig), 0.5) IF (mu0(ig).LT.1.e-4) THEN dist=acos(sin(lati(ig))*sin(latit)+cos(lati(ig))* & cos(latit)*cos(longit-long(ig)))*180/pi IF (dist.gt.90) dist=90 flym_ipm(ig)=(valmin_dl-valmin)/(90.**puis)*(dist)**puis & +valmin ELSE flym_ipm(ig)= mu0(ig)*(valmax-valmin_dl)+valmin_dl !! IPM flux must be weighted by solar flux flym_ipm(ig)=flym_ipm(ig)*flym_sol_pluto/1.0156e12 ! Value in 2015 for flym_sol_pluto ENDIF ENDDO ! ngrid ELSE ! diurnal=false DO ig=1,ngrid flym_sol(ig)=flym_sol_pluto*mu0(ig) flym_ipm(ig)=flym_sol(ig)*1.15 mu_ipm(ig) = max(mu0(ig), 0.5) ENDDO ENDIF ! diurnal IF (firstcall) then write(*,*) 'Fast haze parameters:' write(*,*) 'kch4, ncratio = ',kch4,ncratio firstcall=.false. ENDIF ! note: mu0 = cos(solar zenith angle) ! max(mu0) = declin !************************** !******** LOOP ************ !************************** zdqfasthaze(:,:)=0. zdqsfasthaze(:)=0. zdqch4(:)=0. zdqprec(:)=0. DO ig=1,ngrid !! Update of tracer ch4 and prec_haze zq_ch4(ig)=pq(ig,1,igcm_ch4_gas)+pdq(ig,1,igcm_ch4_gas)*ptimestep zq_prec(ig)=pq(ig,1,igcm_prec_haze)+pdq(ig,1,igcm_prec_haze)*ptimestep !! Security as we loose CH4 molecules if (zq_ch4(ig).lt.0.) then zq_ch4(ig)=0. endif !! Calculation of CH4 optical depth in Lyman alpha for whole atm tauch4=sigch4*1.e-4*avogadro*(zq_ch4(ig)/(mmol(igcm_ch4_gas)*1.e-3))*(pplev(ig,1))/g !! Calculation of Flux reaching the surface if (mu0(ig).gt.1.e-6) then flym_sol_bot(ig)=flym_sol(ig)*exp(-tauch4/(mu0(ig))) endif flym_ipm_bot(ig)=flym_ipm(ig)*exp(-tauch4/mu_ipm(ig)) !! Total flux reaching the surface and flux absorbed by CH4 gradflux(ig)=flym_sol(ig)-flym_sol_bot(ig) + flym_ipm(ig)-flym_ipm_bot(ig) flym_bot(ig)=flym_ipm_bot(ig)+flym_sol_bot(ig) !! Tendancy on ch4 considering 1 photon destroys 1 ch4 by photolysis zdqch4(ig)=-mmol(igcm_ch4_gas)*1.e-3/avogadro*gradflux(ig)*g/(pplev(ig,1)) !! Tendency of precursors created zdqprec(ig)=-zdqch4(ig) !! update precurseur zq_prec zq_prec(ig)=zq_prec(ig)+zdqprec(ig)*ptimestep !! If pressure > pagg, aggregation is possible !! If pressure > pform, haze formation is possible if (pplev(ig,1).gt.pform) then ! New tendancy for prec_haze, turning into haze aerosols zdqconv_prec(ig) = -zq_prec(ig)*(1.-exp(-ptimestep/tconv))/ptimestep else zdqconv_prec(ig) = 0. endif !! Final tendencies DO iq=1,nq tname=noms(iq) if (tname(1:4).eq."haze") then ! Tendency haze in atmosphere kg/kg zdqfasthaze(ig,iq) = -zdqconv_prec(ig)*kch4*(1.+ncratio*14./12.) ! Tendancy for haze directly deposited at the surface kg/m2 zdqsfasthaze(ig) = zdqfasthaze(ig,iq)*pplev(ig,1)/g else if (noms(iq).eq."prec_haze") then zdqfasthaze(ig,igcm_prec_haze)= zdqprec(ig) + zdqconv_prec(ig) else if (noms(iq).eq."ch4_gas".and.(.not.hazeconservch4)) then ! Tendancy for haze directly deposited at the surface kg/m2 zdqfasthaze(ig,igcm_ch4_gas)= zdqch4(ig) endif ENDDO ENDDO ! ngrid end subroutine prodhaze