Changeset 1946 for trunk/LMDZ.TITAN/libf


Ignore:
Timestamp:
Jun 13, 2018, 4:18:12 PM (7 years ago)
Author:
jvatant
Message:

+ Fix a major typo in calchim, in the calculation of krate, leading to almost no photodisso ...
+ Also correct kedd profiles to match values from Yelle et al 08 ( slightly overestimated before )
--JVO

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/LMDZ.TITAN/libf/phytitan/calchim.F90

    r1928 r1946  
    1818  !     
    1919  !          + J. Vatant d'Ollone
    20   !               + 02/17 - adaptation for the issu du generic
     20  !               + 02/17 - adaptation for the new generic-forked physics
    2121  !               + 01/18 - 03/18 - Major transformations :
    2222  !                   - Upper chemistry fields are now stored in startfi
    2323  !                     and defined on a pressure grid from Vervack profile
    24   !                   - Altitudes above GCM top are computed using correct g(z)
    25   !                     to be coherent with obs and with dz=10km for UV processes
    26   !                     but this should be even better if all physiq "under" was doing this !
    2724  !                   - These modifs enables to run chemistry with others resolution than 32x48x55 !
    2825  !                   - Only the actinic fluxes are still read in a 49-lat input but interp. on lat grid
     
    307304     kedd(:) = 1.e3 ! Default value =/= zero
    308305
    309      ! NB : Eddy coeffs from 1D models (e.g. Lavvas et al 08) in altitude but they're rather linked to pressure
     306     ! NB : Eddy coeffs (e.g. Lavvas et al 08, Yelle et al 08) in altitude but they're rather linked to pressure
    310307     !      Below GCM top we have dynamic mixing and for levs < nld=klev-15 the chem. solver ignores diffusion
    311308
     
    316313        IF     ( logp.ge.1.0 .and. logp.le.4.0 ) THEN
    317314              kedd(l) = 10.**(6.0+1.3*(logp-1.0)/3.0)
    318      ! 2E7 at 600 km ~ 10-4 mbar
    319         ELSEIF ( logp.gt.4.0 .and. logp.le.6.0 ) THEN
    320               kedd(l) = 10.**(7.3+0.35*(logp-4.0))
    321      ! 1E8 above 900 km ~ 10-6 mbar
    322         ELSEIF ( logp.gt.6.0                   ) THEN
    323              kedd(l) = 1.e8
     315     ! 2E7 above 600 km ~ 10-4 mbar
     316        ELSEIF ( logp.gt.4.0                   ) THEN
     317             kedd(l) = 2.e7
    324318        ENDIF
    325319     ENDDO
     
    464458           DO i=1,nd_kim+1 ! nd_kim+1 is dissociation of N2 by GCR
    465459
    466               krpddec   =   (   krpd(idec+1,i,ialt  ,klat)   * (1.0-factalt)                   &
    467                               + krpd(idec+1,i,ialt+1,klat)   * factalt       ) * (1.0-factlat) &
    468                           * (   krpd(idec+1,i,ialt  ,klat+1) * (1.0-factalt)                   &
    469                               + krpd(idec+1,i,ialt+1,klat+1) * factalt       ) * factlat
     460                 krpddec   =   (   krpd(idec+1,i,ialt  ,klat)   * (1.0-factalt)                   &
     461                                 + krpd(idec+1,i,ialt+1,klat)   * factalt       ) * (1.0-factlat) &
     462                             + (   krpd(idec+1,i,ialt  ,klat+1) * (1.0-factalt)                   &
     463                                 + krpd(idec+1,i,ialt+1,klat+1) * factalt       ) * factlat
    470464
    471465              if      ( factdec.lt.0. ) then
    472466                 krpddecm1 =   (   krpd(idec  ,i,ialt  ,klat)   * (1.0-factalt)                   &
    473467                                 + krpd(idec  ,i,ialt+1,klat)   * factalt       ) * (1.0-factlat) &
    474                              * (   krpd(idec  ,i,ialt  ,klat+1) * (1.0-factalt)                   &
     468                             + (   krpd(idec  ,i,ialt  ,klat+1) * (1.0-factalt)                   &
    475469                                 + krpd(idec  ,i,ialt+1,klat+1) * factalt       ) * factlat
    476470                 krate(l,i) = krpddecm1 * abs(factdec) + krpddec   * ( 1.0 + factdec)
     
    478472                 krpddecp1 =   (   krpd(idec+2,i,ialt  ,klat)   * (1.0-factalt)                   &
    479473                                 + krpd(idec+2,i,ialt+1,klat)   * factalt       ) * (1.0-factlat) &
    480                              * (   krpd(idec+2,i,ialt  ,klat+1) * (1.0-factalt)                   &
     474                             + (   krpd(idec+2,i,ialt  ,klat+1) * (1.0-factalt)                   &
    481475                                 + krpd(idec+2,i,ialt+1,klat+1) * factalt       ) * factlat
    482476                 krate(l,i) = krpddecp1 * factdec      + krpddec   * ( 1.0 - factdec)
Note: See TracChangeset for help on using the changeset viewer.