Index: trunk/LMDZ.PLUTO/libf/phypluto/aeroptproperties.F90 =================================================================== --- trunk/LMDZ.PLUTO/libf/phypluto/aeroptproperties.F90 (revision 3334) +++ trunk/LMDZ.PLUTO/libf/phypluto/aeroptproperties.F90 (revision 3353) @@ -178,5 +178,5 @@ INTEGER :: ngrid,nlayer ! Aerosol effective radius used for radiative transfer (meter) - REAL,INTENT(IN) :: reffrad(ngrid,nlayer,naerkind) + REAL :: reffrad(ngrid,nlayer,naerkind) ! Aerosol effective variance used for radiative transfer (n.u.) REAL,INTENT(IN) :: nueffrad(ngrid,nlayer,naerkind) @@ -324,12 +324,19 @@ minrad=min(MINVAL(radiustab(1,1,1:nsize(1,1))),MINVAL(radiustab(1,2,1:nsize(1,2)))) maxrad=min(MAXVAL(radiustab(1,1,1:nsize(1,1))),MAXVAL(radiustab(1,2,1:nsize(1,2)))) - IF ((MINVAL(reffrad).LE.minrad).OR.(MAXVAL(reffrad).GE.maxrad)) then + IF ((MINVAL(reffrad).LT.minrad).OR.(MAXVAL(reffrad).GT.maxrad)) then WRITE(*,*) 'Warning: particle size in grid box #' - WRITE(*,*) ig,' is too large to be used by the ' - WRITE(*,*) 'radiative transfer; please extend the ' - WRITE(*,*) 'interpolation grid to larger grain sizes.' + WRITE(*,*) ig,' is larger than optical properties. ' + WRITE(*,*) 'reffrad=',MINVAL(reffrad),'-',MAXVAL(reffrad) WRITE(*,*) 'radiustab=',minrad,'-',maxrad - WRITE(*,*) 'reffrad=',MINVAL(reffrad),'-',MAXVAL(reffrad) - stop + + ! ensure reffrad is within bounds of radiustab + WHERE(reffrad.LT.minrad) + reffrad=minrad + ENDWHERE + WHERE(reffrad.GT.maxrad) + reffrad=maxrad + ENDWHERE + WRITE(*,*) 'Truncated reffrad within radiustab bounds:' + WRITE(*,*) 'New reffrad=',MINVAL(reffrad),'-',MAXVAL(reffrad) ENDIF Index: trunk/LMDZ.PLUTO/libf/phypluto/aerosol_mod.F90 =================================================================== --- trunk/LMDZ.PLUTO/libf/phypluto/aerosol_mod.F90 (revision 3334) +++ trunk/LMDZ.PLUTO/libf/phypluto/aerosol_mod.F90 (revision 3353) @@ -67,4 +67,5 @@ parameter(Nfine=701) character(len=100) :: file_path + character(len=100) :: file_name real,save :: levdat(Nfine),densdat(Nfine) @@ -74,5 +75,16 @@ IF (firstcall) then firstcall=.false. - file_path=trim(datadir)//'/haze_prop/'//hazemmr_file + + if (hazemmr_file/='None')then + file_name = hazemmr_file + print*, 'Read Haze MMR file: ',hazemmr_file + else if(hazedens_file/='None')then + file_name = hazedens_file + print*, 'Read Haze density: ',hazedens_file + else + STOP "No filename given for haze profile. Either set hazemmr_file or hazedens_file" + endif + + file_path=trim(datadir)//'/haze_prop/'//file_name open(224,file=file_path,form='formatted') do ifine=1,Nfine @@ -80,5 +92,5 @@ enddo close(224) - print*, 'TB22 read Haze MMR profile' + print*, 'Read Haze profile: ',file_path ENDIF @@ -91,5 +103,5 @@ ! --> kg m-3 --> kg/kg do iaer=1,naerkind - if(iaer.eq.iaero_haze.and.1.eq.2) then !TB22 activate/deactivate mmr or part density + if(iaer.eq.iaero_haze.and.hazedens_file/='None') then !AF24 activate/deactivate mmr or part density !print*, 'Haze profile is fixed' do ig=1,ngrid @@ -97,4 +109,5 @@ !from number density in cm-3 profmmr(ig,l)=profmmr(ig,l)*1.e6*4./3.*pi*reffrad(ig,l,iaer)**3*rho_q(i_haze)/(pplay(ig,l)/(r*pt(ig,l))) + ! print*, profmmr(ig,l) enddo enddo Index: trunk/LMDZ.PLUTO/libf/phypluto/calc_cpp_mugaz.F90 =================================================================== --- trunk/LMDZ.PLUTO/libf/phypluto/calc_cpp_mugaz.F90 (revision 3334) +++ trunk/LMDZ.PLUTO/libf/phypluto/calc_cpp_mugaz.F90 (revision 3353) @@ -11,5 +11,5 @@ ! ! Authors -! ------- +! ------- ! Robin Wordsworth (2009) ! A. Spiga: make the routine OK with latest changes in rcm1d @@ -25,5 +25,5 @@ character(len=20),parameter :: myname="calc_cpp_mugaz" - real cpp_c + real cpp_c real mugaz_c @@ -41,5 +41,5 @@ else ! all values at 300 K from Engineering Toolbox - if(igas.eq.igas_N2)then + if(igas.eq.igas_CO2)then mugaz_c = mugaz_c + 44.01*gfrac(igas) elseif(igas.eq.igas_N2)then @@ -59,5 +59,5 @@ elseif(igas.eq.igas_NH3)then mugaz_c = mugaz_c + 17.03*gfrac(igas) - elseif(igas.eq.igas_C2H6)then + elseif(igas.eq.igas_C2H6)then ! C2H6 http://encyclopedia.airliquide.com/Encyclopedia.asp?GasID=28 mugaz_c = mugaz_c + 30.07*gfrac(igas) @@ -92,7 +92,7 @@ else ! all values at 300 K from Engineering Toolbox - if(igas.eq.igas_N2)then + if(igas.eq.igas_CO2)then !cpp_c = cpp_c + 0.744*gfrac(igas) ! @ ~210 K (better for - !Mars conditions) + !Mars conditions) cpp_c = cpp_c + 0.846*gfrac(igas)*44.01/mugaz_c elseif(igas.eq.igas_N2)then Index: trunk/LMDZ.PLUTO/libf/phypluto/calc_rayleigh.F90 =================================================================== --- trunk/LMDZ.PLUTO/libf/phypluto/calc_rayleigh.F90 (revision 3334) +++ trunk/LMDZ.PLUTO/libf/phypluto/calc_rayleigh.F90 (revision 3353) @@ -2,14 +2,14 @@ !================================================================== -! +! ! Purpose ! ------- -! Average the Rayleigh scattering in each band, weighting the +! Average the Rayleigh scattering in each band, weighting the ! average by the blackbody function at temperature tstellar. -! Works for an arbitrary mix of gases (currently N2, N2 and +! Works for an arbitrary mix of gases (currently N2, N2 and ! H2 are possible). -! +! ! Authors -! ------- +! ------- ! Robin Wordsworth (2010) ! Jeremy Leconte (2012): Added option for variable gas. Improved water rayleigh (Bucholtz 1995). @@ -19,14 +19,14 @@ ! --------- ! setspv.F -! +! ! Calls ! ----- ! none -! +! !================================================================== use radinc_h, only: L_NSPECTV use radcommon_h, only: WAVEV, BWNV, DWNV, tstellar, tauray, taurayvar, scalep - use gases_h, only: ngasmx, vgas, gnom, gfrac, igas_N2, igas_H2, & + use gases_h, only: ngasmx, vgas, gnom, gfrac, igas_CO2, igas_H2, & igas_H2O, igas_He, igas_N2, igas_CH4, igas_CO use comcstfi_mod, only: g, mugaz, pi @@ -37,5 +37,5 @@ integer N,Nfine,ifine,igas parameter(Nfine=500.0) - real*8 :: P0 = 9.423D+6 ! Rayleigh scattering reference pressure in pascals. P0 = 93*101325 Pa + real*8 :: P0 = 9.423D+6 ! Rayleigh scattering reference pressure in pascals. P0 = 93*101325 Pa logical typeknown @@ -70,10 +70,10 @@ if((igas/=vgas).and.(gfrac(igas).lt.5.e-2))then print*,'Ignoring ',trim(gnom(igas)),' in Rayleigh scattering '// & - 'as its mixing ratio is less than 0.05.' + 'as its mixing ratio is less than 0.05.' ! ignore variable gas in Rayleigh calculation ! ignore gases of mixing ratio < 0.05 in Rayleigh calculation tauconsti(igas) = 0.0 else - if(igas.eq.igas_N2) then + if(igas.eq.igas_CO2) then ! Hansen 1974 : equation (2.32) tauconsti(igas) = (8.7/g)*1.527*scalep/P0 @@ -82,5 +82,5 @@ tauconsti(igas) = (9.81/g)*8.569E-3*scalep/(P0/93.0) elseif(igas.eq.igas_H2O)then - ! tauconsti(igas) = (10.0/g)*9.22E-3*scalep/101325.0 + ! tauconsti(igas) = (10.0/g)*9.22E-3*scalep/101325.0 tauconsti(igas) = 1.98017E-10*scalep/(g*mugaz*1E4) elseif(igas.eq.igas_H2)then @@ -103,5 +103,5 @@ ! 16.04*1.6726*1E-27 is CH4 molecular mass tauconsti(igas) = 24.*pi**3/(1E5/(1.380648813E-23*273.15))**2 * 4./9. * 1E24 * (1/( g *16.04*1.6726*1E-27)) * scalep - + elseif(igas.eq.igas_CO)then ! see above for explanation @@ -127,7 +127,7 @@ endif - + do N=1,L_NSPECTV - + tausum = 0.0 tauwei = 0.0 @@ -147,5 +147,5 @@ tauvari(igas) = 0.0 else - if(igas.eq.igas_N2)then + if(igas.eq.igas_CO2)then tauvari(igas) = (1.0+0.013/wl**2)/wl**4 elseif(igas.eq.igas_N2)then @@ -153,12 +153,12 @@ elseif(igas.eq.igas_H2O)then ! tauvari(igas) = 1.0/wl**4 ! to be improved... - tauvari(igas) = (5.7918E6/(2.38E2-1/wl**2)+1.679E5/(57.36E0-1/wl**2))**2/wl**4 + tauvari(igas) = (5.7918E6/(2.38E2-1/wl**2)+1.679E5/(57.36E0-1/wl**2))**2/wl**4 elseif(igas.eq.igas_H2)then - tauvari(igas) = 1.0/wl**4 + tauvari(igas) = 1.0/wl**4 ! below is exo_k formula : this tauvari is consistent with commented tauconsti for H2 above ! tauvari(igas) = (8.14E-49+1.28E-50/wl**2+1.61E-51/wl**4) / wl**4 elseif(igas.eq.igas_He)then - tauvari(igas) = 1.0/wl**4 + tauvari(igas) = 1.0/wl**4 elseif(igas.eq.igas_CH4)then tauvari(igas) = (4.6662E-4+4.02E-6/wl**2)**2 / wl**4 @@ -186,5 +186,5 @@ tauweivar=tauweivar+df tausumvar=tausumvar+tauvarvar*df - + enddo TAURAY(N)=tausum/tauwei Index: trunk/LMDZ.PLUTO/libf/phypluto/callcorrk_pluto_mod.F90 =================================================================== --- trunk/LMDZ.PLUTO/libf/phypluto/callcorrk_pluto_mod.F90 (revision 3334) +++ trunk/LMDZ.PLUTO/libf/phypluto/callcorrk_pluto_mod.F90 (revision 3353) @@ -345,6 +345,4 @@ CALL cooling_hcn_c2h2(ngrid,nlayer,pplay, & pt,dtlw_hcn_c2h2) - ! write(*,*) "Not supported yet" - STOP endif Index: trunk/LMDZ.PLUTO/libf/phypluto/callkeys_mod.F90 =================================================================== --- trunk/LMDZ.PLUTO/libf/phypluto/callkeys_mod.F90 (revision 3334) +++ trunk/LMDZ.PLUTO/libf/phypluto/callkeys_mod.F90 (revision 3353) @@ -197,5 +197,6 @@ logical,save :: oldplutocorrk !$OMP THREADPRIVATE(oldplutocorrk) - + logical,save :: oldplutosedim +!$OMP THREADPRIVATE(oldplutosedim) logical,save :: global1d real,save :: szangle Index: trunk/LMDZ.PLUTO/libf/phypluto/callsedim_pluto.F90 =================================================================== --- trunk/LMDZ.PLUTO/libf/phypluto/callsedim_pluto.F90 (revision 3353) +++ trunk/LMDZ.PLUTO/libf/phypluto/callsedim_pluto.F90 (revision 3353) @@ -0,0 +1,156 @@ + SUBROUTINE callsedim_pluto(ngrid,nlay, ptimestep, & + pplev,zlev,pt,pdt,rice_ch4,rice_co, & + pq, pdqfi, pdqsed,pdqs_sed,nq,pphi) + + use radinc_h, only : naerkind + use tracer_h, only: igcm_ch4_ice,igcm_co_ice,radius,rho_q + use comcstfi_mod, only: g + IMPLICIT NONE +!================================================================== +! +! Purpose +! ------- +! Calculates sedimentation of aerosols depending on their +! density and radius. +! +! Authors +! ------- +! F. Forget (1999) +! Tracer generalisation by E. Millour (2009) +! +!================================================================== + +!----------------------------------------------------------------------- +! declarations: +! ------------- + +! +! arguments: +! ---------- + + INTEGER ngrid ! number of horizontal grid points + INTEGER nlay ! number of atmospheric layers + REAL ptimestep ! physics time step (s) + REAL pplev(ngrid,nlay+1) ! pressure at inter-layers (Pa) + REAL pt(ngrid,nlay) ! temperature at mid-layer (K) + REAL pdt(ngrid,nlay) ! tendency on temperature + REAL zlev(ngrid,nlay+1) ! altitude at layer boundaries + REAL pphi(ngrid,nlay) ! geopotential + + +! Traceurs : + integer nq ! number of tracers + real pq(ngrid,nlay,nq) ! tracers (kg/kg) + real pdqfi(ngrid,nlay,nq) ! tendency before sedimentation (kg/kg.s-1) + real pdqsed(ngrid,nlay,nq) ! tendency due to sedimentation (kg/kg.s-1) + real pdqs_sed(ngrid,nq) ! flux at surface (kg.m-2.s-1) + +! local: +! ------ + + INTEGER l,ig, iq + + real zqi(ngrid,nlay) ! to locally store tracers + real zt(ngrid,nlay) ! to locally store temperature (K) + real masse (ngrid,nlay) ! Layer mass (kg.m-2) + real epaisseur (ngrid,nlay) ! Layer thickness (m) + real wq(ngrid,nlay+1) ! displaced tracer mass (kg.m-2) + real rfall_ch4(ngrid,nlay) + real rfall_co(ngrid,nlay) + real rice_ch4(ngrid,nlay) + real rice_co(ngrid,nlay) + + LOGICAL firstcall + SAVE firstcall + DATA firstcall/.true./ + +! ** un petit test de coherence +! -------------------------- + + IF (firstcall) THEN + IF(ngrid.NE.ngrid) THEN + PRINT*,'STOP dans callsedim' + PRINT*,'probleme de dimensions :' + PRINT*,'ngrid =',ngrid + PRINT*,'ngrid =',ngrid + STOP + ENDIF + + firstcall=.false. + ENDIF ! of IF (firstcall) + +!======================================================================= +! Preliminary calculation of the layer characteristics +! (mass (kg.m-2), thickness (m), etc. + + do l=1,nlay + do ig=1, ngrid + masse(ig,l)=(pplev(ig,l) - pplev(ig,l+1)) /g + epaisseur(ig,l)= zlev(ig,l+1) - zlev(ig,l) + zt(ig,l)=pt(ig,l)+pdt(ig,l)*ptimestep + end do + end do + + do iq=1,nq + if(radius(iq).gt.1.e-12) then ! no sedimentation for gases (defined by radius=0) +! Radius values are defined in initracer +! The value of q is updated after the other parameterisations + + do l=1,nlay + do ig=1,ngrid + ! store locally updated tracers + zqi(ig,l)=pq(ig,l,iq)+pdqfi(ig,l,iq)*ptimestep + +! cf sur Mars: +! On affecte un rayon moyen aux poussieres a chaque altitude du type : +! r(z)=r0*exp(-z/H) avec r0=0.8 micron et H=18 km. +! rfall(ig,l)=max( rice(ig,l)*1.5,rdust(ig,l) ) +! Pluton : choix de rfall a la place de radius + if (iq.eq.igcm_ch4_ice) then + ! TB: rfall_ch4(ig,l)=max( rice_ch4(ig,l)*1.5,2.e-7) + rfall_ch4(ig,l)=max( rice_ch4(ig,l),2.e-7) + rfall_ch4(ig,l)=min(rfall_ch4(ig,l),1.e-4) + endif + if (iq.eq.igcm_co_ice) then + rfall_co(ig,l)=max( rice_co(ig,l),2.e-7) + rfall_co(ig,l)=min(rfall_co(ig,l),1.e-4) + endif + enddo + enddo ! of do l=1,nlay + +!======================================================================= +! Calculate the transport due to sedimentation for each tracer + + if (iq.eq.igcm_ch4_ice) then + !if (iceparty.and.(iq.eq.igcm_ch4_ice)) then + call newsedim_pluto(ngrid,nlay,ngrid*nlay,ptimestep, & + pplev,masse,epaisseur,zt,rfall_ch4,rho_q(iq),zqi,wq,pphi) + else if (iq.eq.igcm_co_ice) then + call newsedim_pluto(ngrid,nlay,ngrid*nlay,ptimestep, & + pplev,masse,epaisseur,zt,rfall_co,rho_q(iq),zqi,wq,pphi) + else if ((radius(iq).gt.0.)) then ! haze tracers + call newsedim_pluto(ngrid,nlay,1,ptimestep, & + pplev,masse,epaisseur,zt,radius(iq),rho_q(iq),zqi,wq,pphi) + endif + +!======================================================================= +! Calculate the tendencies + + do ig=1,ngrid +! Ehouarn: with new way of tracking tracers by name, this is simply + pdqs_sed(ig,iq)=wq(ig,1)/ptimestep + end do + + DO l = 1, nlay + DO ig=1,ngrid + pdqsed(ig,l,iq)=(zqi(ig,l)- & + (pq(ig,l,iq) + pdqfi(ig,l,iq)*ptimestep))/ptimestep + ENDDO + ENDDO + + endif ! of if(radius(iq).gt.1.e-12) + enddo ! of do iq=1,nq + + RETURN + END + Index: trunk/LMDZ.PLUTO/libf/phypluto/datafile_mod.F90 =================================================================== --- trunk/LMDZ.PLUTO/libf/phypluto/datafile_mod.F90 (revision 3334) +++ trunk/LMDZ.PLUTO/libf/phypluto/datafile_mod.F90 (revision 3353) @@ -25,4 +25,6 @@ character(len=300),save :: hazerad_file character(len=300),save :: hazemmr_file + character(len=300),save :: hazedens_file + end module datafile_mod Index: trunk/LMDZ.PLUTO/libf/phypluto/inifis_mod.F90 =================================================================== --- trunk/LMDZ.PLUTO/libf/phypluto/inifis_mod.F90 (revision 3334) +++ trunk/LMDZ.PLUTO/libf/phypluto/inifis_mod.F90 (revision 3353) @@ -13,5 +13,5 @@ use radcommon_h, only: ini_radcommon_h use radii_mod, only: radfixed, Nmix_n2 - use datafile_mod, only: datadir, hazeprop_file, hazerad_file, hazemmr_file + use datafile_mod, only: datadir,hazeprop_file,hazerad_file,hazemmr_file,hazedens_file use comdiurn_h, only: sinlat, coslat, sinlon, coslon use comgeomfi_h, only: totarea, totarea_planet @@ -148,8 +148,11 @@ if (is_master) write(*,*) "Haze mmr datafile" - hazemmr_file="hazemmr.txt" ! default file + hazemmr_file="None" ! default file call getin_p("hazemmr_file",hazemmr_file) if (is_master) write(*,*) trim(rname)//" hazemmr_file = ",trim(hazemmr_file) - + if (is_master) write(*,*) "Haze density datafile" + hazedens_file="None" ! default file + call getin_p("hazedens_file",hazedens_file) + if (is_master) write(*,*) trim(rname)//" hazedens_file = ",trim(hazedens_file) if (is_master) write(*,*) trim(rname)//& @@ -312,16 +315,4 @@ call getin_p("UseTurbDiff",UseTurbDiff) if (is_master) write(*,*) trim(rname)//": UseTurbDiff = ",UseTurbDiff - - if (is_master) write(*,*) trim(rname)//": use vdifc from old Pluto ?" - oldplutovdifc=.true. ! default value - call getin_p("oldplutovdifc",oldplutovdifc) - if (is_master) write(*,*) trim(rname)//": oldplutovdifc = ",oldplutovdifc - - if (is_master) write(*,*) trim(rname)//& - ": call pluto.old correlated-k radiative transfer ?" - oldplutocorrk=.true. ! default value - call getin_p("oldplutocorrk",oldplutocorrk) - if (is_master) write(*,*) trim(rname)//": oldplutocorrk = ",oldplutocorrk - if (is_master) write(*,*) trim(rname)//": call convective adjustment ?" @@ -636,4 +627,22 @@ " thresh_non2 = ",thresh_non2 + ! use Pluto.old routines + if (is_master) write(*,*) trim(rname)//": use vdifc from old Pluto ?" + oldplutovdifc=.true. ! default value + call getin_p("oldplutovdifc",oldplutovdifc) + if (is_master) write(*,*) trim(rname)//": oldplutovdifc = ",oldplutovdifc + + if (is_master) write(*,*) trim(rname)//& + ": call pluto.old correlated-k radiative transfer ?" + oldplutocorrk=.true. ! default value + call getin_p("oldplutocorrk",oldplutocorrk) + if (is_master) write(*,*) trim(rname)//": oldplutocorrk = ",oldplutocorrk + + if (is_master) write(*,*) trim(rname)//& + ": call pluto.old sedimentation ?" + oldplutosedim=.true. ! default value + call getin_p("oldplutosedim",oldplutosedim) + if (is_master) write(*,*) trim(rname)//": oldplutosedim = ",oldplutosedim + !*************************************************************** !! Haze options @@ -707,4 +716,39 @@ if (is_master)write(*,*)trim(rname)//& "nb_monomer = ",nb_monomer + + if (is_master)write(*,*)trim(rname)//& + "fixed gaseous CH4 mixing ratio for rad transfer ?" + ch4fix=.false. ! default value + call getin_p("ch4fix",ch4fix) + if (is_master)write(*,*)trim(rname)//& + " ch4fix = ",ch4fix + if (is_master)write(*,*)trim(rname)//& + "fixed gaseous CH4 mixing ratio for rad transfer ?" + vmrch4fix=0.5 ! default value + call getin_p("vmrch4fix",vmrch4fix) + if (is_master)write(*,*)trim(rname)//& + " vmrch4fix = ",vmrch4fix + vmrch4_proffix=.false. ! default value + call getin_p("vmrch4_proffix",vmrch4_proffix) + if (is_master)write(*,*)trim(rname)//& + " vmrch4_proffix = ",vmrch4_proffix + + if (is_master)write(*,*)trim(rname)//& + "call specific cooling for radiative transfer ?" + cooling=.false. ! default value + call getin_p("cooling",cooling) + if (is_master)write(*,*)trim(rname)//& + " cooling = ",cooling + + if (is_master)write(*,*)trim(rname)//& + "NLTE correction for methane heating rates?" + nlte=.false. ! default value + call getin_p("nlte",nlte) + if (is_master)write(*,*)trim(rname)//& + " nlte = ",nlte + strobel=.false. ! default value + call getin_p("strobel",strobel) + if (is_master)write(*,*)trim(rname)//& + " strobel = ",strobel if (is_master)write(*,*)trim(rname)//& Index: trunk/LMDZ.PLUTO/libf/phypluto/physiq_mod.F90 =================================================================== --- trunk/LMDZ.PLUTO/libf/phypluto/physiq_mod.F90 (revision 3334) +++ trunk/LMDZ.PLUTO/libf/phypluto/physiq_mod.F90 (revision 3353) @@ -55,6 +55,7 @@ n2cond,nearn2cond,noseason_day,conservn2, & convergeps,kbo,triton,paleo,paleoyears, & - carbox, methane, oldplutovdifc, oldplutocorrk, & - aerohaze,haze_proffix,source_haze, & + carbox, methane,& + oldplutovdifc,oldplutocorrk,oldplutosedim, & + aerohaze,haze_proffix,source_haze,& season, sedimentation,generic_condensation, & specOLR, & @@ -979,4 +980,5 @@ cloudfrac,totcloudfrac,.false., & firstcall,lastcall) + albedo_equivalent(1:ngrid)=albedo(1:ngrid,1) else call callcorrk(ngrid,nlayer,pq,nq,qsurf, & @@ -1480,15 +1482,4 @@ end if - ! if (.not. water) then - ! Compute GCS (Generic Condensable Specie) cloud fraction. For now we can not have both water cloud fraction and GCS cloud fraction - ! Water is the priority - ! If you have set water and generic_condensation, then cloudfrac will be water cloudfrac - ! - ! If you have set generic_condensation (and not water) and you have set several GCS - ! then cloudfrac will be only the cloudfrac of the last generic tracer - ! (Because it is rewritten every tracer in the loop) - ! - ! Maybe one should create a cloudfrac_generic(ngrid,nlayer,nq) with 3 dimensions, the last one for tracers - ! Let's loop on tracers cloudfrac(:,:)=0.0 @@ -1507,6 +1498,4 @@ endif !generic_condensation - !Generic Rain !AF24: removed - ! Sedimentation. if (sedimentation) then @@ -1515,11 +1504,14 @@ zdqssed(1:ngrid,1:nq) = 0.0 - ! if(watertest)then !AF24: removed - - call callsedim(ngrid,nlayer,ptimestep, & + if (oldplutosedim)then + call callsedim_pluto(ngrid,nlayer,ptimestep, & + pplev,zzlev,pt,pdt,rice_ch4,rice_co, & + pq,pdq,zdqsed,zdqssed,nq,pphi) + + else + call callsedim(ngrid,nlayer,ptimestep, & pplev,zzlev,pt,pdt,pq,pdq, & zdqsed,zdqssed,nq) - - ! if(watertest)then !AF24: removed + endif ! Whether it falls as rain or snow depends only on the surface temperature @@ -1527,12 +1519,5 @@ dqsurf(1:ngrid,1:nq) = dqsurf(1:ngrid,1:nq) + zdqssed(1:ngrid,1:nq) -!! call writediagfi(ngrid,"callsedim_post_dqsurf"," "," ",2,dqsurf(1:ngrid,igcm_h2o_gas)) - - ! ! Test water conservation !AF24: removed - end if ! end of 'sedimentation' - -!! call writediagfi(ngrid,"mass_redist_pre_dqsurf"," "," ",2,dqsurf(1:ngrid,igcm_h2o_gas)) -!! call writediagfi(ngrid,"mass_redist_pre_qsurf"," "," ",2,qsurf(1:ngrid,igcm_h2o_gas)) ! --------------- @@ -2096,4 +2081,5 @@ !call wstats(ngrid,"ALB_1st","First Band Surface albedo"," ",2,albedo(:,1)) call wstats(ngrid,"p","Pressure","Pa",3,pplay) + call wstats(ngrid,"emis","Emissivity","",2,emis) call wstats(ngrid,"temp","Atmospheric temperature","K",3,zt) call wstats(ngrid,"u","Zonal (East-West) wind","m.s-1",3,zu) @@ -2157,4 +2143,5 @@ endif + call writediagfi(ngrid,"emis","Emissivity","",2,emis) call writediagfi(ngrid,"temp","temperature","K",3,zt) call writediagfi(ngrid,"teta","potential temperature","K",3,zh) Index: trunk/LMDZ.PLUTO/libf/phypluto/radii_mod.F90 =================================================================== --- trunk/LMDZ.PLUTO/libf/phypluto/radii_mod.F90 (revision 3334) +++ trunk/LMDZ.PLUTO/libf/phypluto/radii_mod.F90 (revision 3353) @@ -106,5 +106,5 @@ enddo close(223) - print*, 'TB22 READ HAZERAD' + print*, 'Read hazerad from ',file_path ENDIF