Ignore:
Timestamp:
Jul 18, 2012, 10:54:57 AM (13 years ago)
Author:
jleconte
Message:

18/07/2012 == JL

  • New water cycle scheme:
    • largescale now in F90. Robustness increased by i) including evap inside largescale ii) computing the

condensed water amount iteratively

  • same improvements in moistadj.
  • Water thermodynamical data and saturation curves centralized in module watercommn_h
    • The saturation curves used are now Tetens formula as they are analyticaly inversible (Ts(P)-> Ps(T)).

New saturation curve yields very good agreement with the former one.

  • Saturation curves are now generalized for arbitrary water amount (not just q<<1)
  • The old watersat should be removed soon.
  • The effect of water vapor on total (surface) pressure can be taken into account by setting

mass_redistrib=.true. in callphys.def (routine mass_redistribution inspired from co2_condense in martian
model but with a different scheme as many routines evaporate/condense water vapor).

  • New cloud and precipitation scheme (JL + BC):
    • The default recovery assumption for computing the total cloud fraction has been changed (total random gave too

large cloud fractions). See totalcloudfrac.F90 for details and to change this.

  • Totalcloudfraction now set the total cloud fraction to the fraction of the

optically thickest cloud and totalcloudfrac is thus called in aeropacity.

  • Only the total cloud fraction is used to compute optical depth in aeropacity (no more effective

optical depth with exponential formula).

  • 4 precipitation schemes are now available (see rain.F90 for details). The choice can be made using precip_scheme

in callphys.def. Usage of the more physically based model of Boucher et al 95 (precip_scheme=4) is recommended.
default behavior is set to the former "simple scheme" (precip_scheme=1).

  • See rain.f90 to determine the parameter to be defined in callphys.def as a function of the precipitation scheme used.
  • Physiq.F90 now written in a matricial (more F90) way.
  • Radii (H2O and CO2 cloud particles, aerosols, duts, ...) calculations now centralized in module radii_mod.F90

and work with the new aerosol scheme implemented by Laura K. Some inconsistency may remain in callsedim.


Implementation compiled with ifort and pgf90.
gcm.e runs in Earth and Early Mars case with CO2 and H2O cycle + dust.

File:
1 edited

Legend:

Unmodified
Added
Removed

  • trunk/LMDZ.GENERIC/libf/phystd/callcorrk.F90

    r726 r728  
    1414      use ioipsl_getincom
    1515      use gases_h
     16      use radii_mod, only : su_aer_radii,co2_reffrad,h2o_reffrad,dust_reffrad,h2so4_reffrad
    1617      use aerosol_mod
    1718
    … …  
    7475
    7576      REAL reffrad(ngrid,nlayer,naerkind)
    76       REAL nueffrad(ngrid,nlayer,naerkind)
     77      REAL, SAVE :: nueffrad(ngridmx,nlayermx,naerkind)
    7778
    7879!     OUTPUT
    … …  
    171172      real pqtest(ngridmx,nlayer,nq)
    172173
    173 !     are particle radii fixed?
    174       logical, save ::  radfixed
    175174      real maxrad, minrad
    176 
    177 !     Local water cloud optical properties
    178       real, parameter ::  rad_froid=35.0e-6
    179       real, parameter ::  rad_chaud=10.0e-6
    180       real, parameter ::  coef_chaud=0.13
    181       real, parameter ::  coef_froid=0.09
    182       real zfice
    183      
     175           
    184176      real CBRT
    185177      external CBRT
    … …  
    201193      qsiaer(:,:,:)=0.0
    202194      giaer(:,:,:) =0.0
    203       radfixed=.false.
    204195
    205196      if(firstcall) then
    206          if(kastprof)then
    207             radfixed=.true.
    208          endif 
    209197
    210198         call system('rm -f surf_vals_long.out')
    … …  
    212200!--------------------------------------------------
    213201!     Effective radius and variance of the aerosols
    214 
    215 
    216 !     these values will change once the microphysics gets to work
    217 !     UNLESS tracer=.false., in which case we should be working with
    218 !     a fixed aerosol layer, and be able to define reffrad in a
    219 !     .def file. To be improved!
    220 !     Generalization of aerosols added by LK
    221             do iaer = 1,naerkind
    222               if (iaer.eq.iaero_co2) then ! active CO2 aerosols 
    223                 do l=1,nlayer
    224                    do ig=1,ngrid
    225                      reffrad(ig,l,iaer)  = 1.e-4
    226                      nueffrad(ig,l,iaer) = 0.1
    227                    enddo
    228                 enddo
    229               endif
    230               if (iaer.eq.iaero_h2o) then ! active H2O aerosols
    231                 do l=1,nlayer
    232                    do ig=1,ngrid
    233                      reffrad(ig,l,iaer)  = 1.e-5
    234                      nueffrad(ig,l,iaer) = 0.1
    235                    enddo
    236                 enddo
    237               endif
    238               if(iaer.eq.iaero_dust)then ! active dust
    239                 do l=1,nlayer
    240                    do ig=1,ngrid
    241                      reffrad(ig,l,iaer)  = 1.e-5
    242                      nueffrad(ig,l,iaer) = 0.1
    243                    enddo
    244                 enddo
    245               endif
    246               if(iaer.eq.iaero_h2so4)then ! Active h2so4 aerosols
    247                 do l=1,nlayer
    248                    do ig=1,ngrid
    249                      reffrad(ig,l,iaer)  = 1.e-6
    250                      nueffrad(ig,l,iaer) = 0.1
    251                    enddo
    252                 enddo
    253               endif
    254             enddo
    255 
     202         if(naerkind.gt.4)then
     203            print*,'Code not general enough to deal with naerkind > 4 yet.'
     204            call abort
     205         endif
     206         call su_aer_radii(reffrad,nueffrad)
     207         
     208         
     209
     210!--------------------------------------------------
     211!     set up correlated k
    256212         print*, "callcorrk: Correlated-k data base folder:",trim(datadir)
    257213         call getin("corrkdir",corrkdir)
    … …  
    305261!     Initialization on every call   
    306262
    307       do l=1,nlayer
    308          do ig=1,ngrid
    309             do iaer=1,naerkind
    310                nueffrad(ig,l,iaer) = 0.1 ! stays at 0.1
    311             enddo
    312          enddo
    313       enddo
    314 
    315 
     263!--------------------------------------------------
     264!     Effective radius and variance of the aerosols
    316265      do iaer=1,naerkind
    317266
    318        if(iaer.eq.iaero_co2)then
    319          if(radfixed)then
    320            do l=1,nlayer
    321             do ig=1,ngrid
    322                reffrad(ig,l,iaer)    = 5.e-5 ! CO2 ice
    323             enddo
    324            enddo
    325            print*,'CO2 ice particle size =',reffrad(1,1,iaer)/1.e-6,'um'
    326          else
    327            maxrad=0.0
    328            !averad=0.0
    329            minrad=1.0
    330            do l=1,nlayer
    331 
    332              !masse = (pplev(ig,l) - pplev(ig,l+1))/g
    333 
    334             do ig=1,ngrid
    335                if(tracer.and.igcm_co2_ice.gt.0)then
    336 
    337                   if(igcm_co2_ice.lt.1)then
    338                      print*,'Tracers but no CO2 ice still'
    339                      print*,'seems to be a problem...'
    340                      print*,'Aborting in callcorrk.'
    341                      stop
    342                   endif
    343                   reffrad(ig,l,iaer) = CBRT( 3*pq(ig,l,igcm_co2_ice)/ &
    344                        (4*Nmix_co2*pi*rho_co2) )
    345                endif
    346                reffrad(ig,l,iaer) = max(reffrad(ig,l,iaer),1.e-6)
    347                reffrad(ig,l,iaer) = min(reffrad(ig,l,iaer),500.e-6)
    348                !averad = averad + reffrad(ig,l,iaer)*area(ig)
    349                maxrad = max(reffrad(ig,l,iaer),maxrad)
    350                minrad = min(reffrad(ig,l,iaer),minrad)
    351             enddo
    352            enddo
    353            if(igcm_co2_ice.gt.0)then
    354              print*,'Max. CO2 ice particle size = ',maxrad/1.e-6,' um'
    355              print*,'Min. CO2 ice particle size = ',minrad/1.e-6,' um'
    356            endif
    357          endif
    358        endif
    359 
    360 
    361        if(iaer.eq.iaero_h2o)then
    362           if(radfixed) then
    363             do l=1,nlayer
    364                do ig=1,ngrid
    365                   !reffrad(ig,l,iaer) = 2.e-5 ! H2O ice
    366                   reffrad(ig,l,iaer) = rad_chaud ! H2O ice
    367 
    368                   zfice = 1.0 - (pt(ig,l)-T_h2O_ice_clouds) / (T_h2O_ice_liq-T_h2O_ice_clouds)
    369                   zfice = MIN(MAX(zfice,0.0),1.0)
    370                   reffrad(ig,l,2)= rad_chaud * (1.-zfice) + rad_froid * zfice
    371                   nueffrad(ig,l,2) = coef_chaud * (1.-zfice) + coef_froid * zfice
    372                enddo
    373             enddo
    374             print*,'H2O ice particle size=',reffrad(1,1,iaer)/1.e-6,'um'
    375           else
    376              if(water)then
    377                maxrad=0.0
    378                minrad=1.0
    379                do l=1,nlayer
    380                 do ig=1,ngrid
    381                   reffrad(ig,l,iaer) = CBRT( 3*pq(ig,l,igcm_h2o_ice)/ &
    382                        (4*Nmix_h2o*pi*rho_ice) )
    383                   reffrad(ig,l,iaer) = max(reffrad(ig,l,iaer),1.e-6)
    384                   reffrad(ig,l,iaer) = min(reffrad(ig,l,iaer),100.e-6)
    385 
    386                   maxrad = max(reffrad(ig,l,iaer),maxrad)
    387                   minrad = min(reffrad(ig,l,iaer),minrad)
    388                 enddo
    389                enddo
    390                print*,'Max. H2O ice particle size = ',maxrad/1.e-6,' um'
    391                print*,'Min. H2O ice particle size = ',minrad/1.e-6,' um'
    392              endif
    393          endif
    394        endif
    395 
     267         if ((iaer.eq.iaero_co2).and.tracer.and.(igcm_co2_ice.gt.0)) then ! treat condensed co2 particles.
     268            call co2_reffrad(pq,reffrad)
     269            print*,'Max. CO2 ice particle size = ',maxval(reffrad(1:ngridmx,1:nlayermx,iaer))/1.e-6,' um'
     270            print*,'Min. CO2 ice particle size = ',minval(reffrad(1:ngridmx,1:nlayermx,iaer))/1.e-6,' um'
     271         end if
     272         if ((iaer.eq.iaero_h2o).and.water) then ! treat condensed water particles. to be generalized for other aerosols
     273            call h2o_reffrad(pq,pt,reffrad,nueffrad)
     274            print*,'Max. H2O cloud particle size = ',maxval(reffrad(1:ngridmx,1:nlayermx,iaer))/1.e-6,' um'
     275            print*,'Min. H2O cloud particle size = ',minval(reffrad(1:ngridmx,1:nlayermx,iaer))/1.e-6,' um'
     276         endif
    396277         if(iaer.eq.iaero_dust)then
    397             do l=1,nlayer
    398                do ig=1,ngrid
    399                   reffrad(ig,l,iaer) = 2.e-6 ! dust
    400                enddo
    401             enddo
     278            call dust_reffrad(reffrad)
    402279            print*,'Dust particle size = ',reffrad(1,1,iaer)/1.e-6,' um'
    403280         endif
    404 
    405281         if(iaer.eq.iaero_h2so4)then
    406             do l=1,nlayer
    407                do ig=1,ngrid
    408                   reffrad(ig,l,iaer) = 1.e-6 ! h2so4
    409                enddo
    410             enddo
     282            call h2so4_reffrad(reffrad)
    411283            print*,'H2SO4 particle size =',reffrad(1,1,iaer)/1.e-6,' um'
    412284         endif
    413       enddo ! do iaer=1,naerkind
    414          ! added by LK
     285      end do !iaer=1,naerkind
    415286
    416287
    … …  
    641512      if(.not.kastprof)then
    642513      ! IMPORTANT: Now convert from kg/kg to mol/mol
    643       do k=1,L_LEVELS
    644          qvar(k) = qvar(k)/epsi
    645       end do
     514         do k=1,L_LEVELS
     515            qvar(k) = qvar(k)/(epsi+qvar(k)*(1.-epsi))
     516         end do
    646517      end if
    647518
    … …  
    764635            print*,'Maximum temperature is outside the radiative'
    765636            print*,'transfer kmatrix bounds, exiting.'
    766             !print*,'WARNING, OVERRIDING FOR TEST'
    767             call abort
     637            print*,'level,grid,T',k,ig,tlevrad(k)
     638            print*,'WARNING, OVERRIDING FOR TEST'
     639            !call abort
     640            tlevrad(k)=tgasmax
    768641         endif
    769642      enddo
Note: See TracChangeset for help on using the changeset viewer.