Ignore:
Timestamp:
Mar 15, 2018, 8:07:03 PM (6 years ago)
Author:
oboucher
Message:

Implementing the MACspV2 aerosol plume climatology
which can be called by setting flag_aerosol=7
and aerosols1980.nc pointing to aerosols.nat.nc.
Also requires the MACv2.0-SP_v1.nc file.

File:
1 edited

Legend:

Unmodified
Added
Removed

  • LMDZ6/trunk/libf/phylmd/newmicro.F90

    r3265 r3274  
    11! $Id$
    22
    3 SUBROUTINE newmicro(ok_cdnc, bl95_b0, bl95_b1, paprs, pplay, t, pqlwp, pclc, &
     3SUBROUTINE newmicro(flag_aerosol, ok_cdnc, bl95_b0, bl95_b1, paprs, pplay, t, pqlwp, pclc, &
    44    pcltau, pclemi, pch, pcl, pcm, pct, pctlwp, xflwp, xfiwp, xflwc, xfiwc, &
    55    mass_solu_aero, mass_solu_aero_pi, pcldtaupi, re, fl, reliq, reice, &
    … …  
    88  USE dimphy
    99  USE phys_local_var_mod, ONLY: scdnc, cldncl, reffclwtop, lcc, reffclws, &
    10       reffclwc, cldnvi, lcc3d, lcc3dcon, lcc3dstra, icc3dcon, icc3dstra, zfice
     10      reffclwc, cldnvi, lcc3d, lcc3dcon, lcc3dstra, icc3dcon, icc3dstra,  &
     11      zfice, dNovrN
    1112  USE phys_state_var_mod, ONLY: rnebcon, clwcon
    1213  USE icefrac_lsc_mod ! computes ice fraction (JBM 3/14)
    … …  
    141142  ! within the grid cell)
    142143
     144  INTEGER flag_aerosol
    143145  LOGICAL ok_cdnc
    144146  REAL bl95_b0, bl95_b1 ! Parameter in B&L 95-Formula
    … …  
    216218        xflwc(i, k) = (1.-zfice(i,k))*pqlwp(i, k)
    217219        xfiwc(i, k) = zfice(i, k)*pqlwp(i, k)
    218       END DO
    219     END DO
     220      ENDDO
     221    ENDDO
    220222  ELSE ! of IF (iflag_t_glace.EQ.0)
    221223    DO k = 1, klev
    … …  
    234236        xflwc(i, k) = (1.-zfice(i,k))*pqlwp(i, k)
    235237        xfiwc(i, k) = zfice(i, k)*pqlwp(i, k)
    236       END DO
    237     END DO
     238      ENDDO
     239    ENDDO
    238240  ENDIF
    239241
    … …  
    244246    DO k = 1, klev
    245247      DO i = 1, klon
    246 
    247248        ! Formula "D" of Boucher and Lohmann, Tellus, 1995
    248249        ! Cloud droplet number concentration (CDNC) is restricted
    249250        ! to be within [20, 1000 cm^3]
    250 
    251         ! --present-day case
    252         cdnc(i, k) = 10.**(bl95_b0+bl95_b1*log(max(mass_solu_aero(i,k), &
    253           1.E-4))/log(10.))*1.E6 !-m-3
    254         cdnc(i, k) = min(1000.E6, max(cdnc_min_m3,cdnc(i,k)))
    255251
    256252        ! --pre-industrial case
    … …  
    258254          1.E-4))/log(10.))*1.E6 !-m-3
    259255        cdnc_pi(i, k) = min(1000.E6, max(cdnc_min_m3,cdnc_pi(i,k)))
     256
     257      ENDDO
     258    ENDDO
     259
     260    !--flag_aerosol=7 => MACv2SP climatology
     261    !--in this case there is an enhancement factor
     262    IF (flag_aerosol .EQ. 7) THEN
     263
     264      !--present-day
     265      DO k = 1, klev
     266        DO i = 1, klon
     267          cdnc(i, k) = cdnc_pi(i,k)*dNovrN(i)
     268        ENDDO
     269      ENDDO
     270
     271    !--standard case
     272    ELSE
     273
     274      DO k = 1, klev
     275        DO i = 1, klon
     276
     277          ! Formula "D" of Boucher and Lohmann, Tellus, 1995
     278          ! Cloud droplet number concentration (CDNC) is restricted
     279          ! to be within [20, 1000 cm^3]
     280
     281          ! --present-day case
     282          cdnc(i, k) = 10.**(bl95_b0+bl95_b1*log(max(mass_solu_aero(i,k), &
     283            1.E-4))/log(10.))*1.E6 !-m-3
     284          cdnc(i, k) = min(1000.E6, max(cdnc_min_m3,cdnc(i,k)))
     285
     286        ENDDO
     287      ENDDO
     288
     289    ENDIF !--flag_aerosol
     290
     291    !--computing cloud droplet size
     292    DO k = 1, klev
     293      DO i = 1, klon
    260294
    261295        ! --present-day case
    … …  
    292326            zfiwp_var*(3.448E-03+2.431/rei)
    293327
    294         END IF
    295 
    296       END DO
    297     END DO
     328        ENDIF
     329
     330      ENDDO
     331    ENDDO
    298332
    299333  ELSE !--not ok_cdnc
    … …  
    305339        rad_chaud(i, k) = rad_chau2
    306340        rad_chaud_pi(i, k) = rad_chau2
    307       END DO
    308     END DO
     341      ENDDO
     342    ENDDO
    309343    DO k = min(3, klev) + 1, klev
    310344      DO i = 1, klon
    311345        rad_chaud(i, k) = rad_chau1
    312346        rad_chaud_pi(i, k) = rad_chau1
    313       END DO
    314     END DO
    315 
    316   END IF !--ok_cdnc
     347      ENDDO
     348    ENDDO
     349
     350  ENDIF !--ok_cdnc
    317351
    318352  ! --computation of cloud optical depth and emissivity
    … …  
    389423        pclemi(i, k) = 1.0 - exp(-coef_chau*zflwp_var-df*k_ice*zfiwp_var)
    390424
    391       END IF
     425      ENDIF
    392426
    393427      reice(i, k) = rei
    … …  
    396430      xfiwp(i) = xfiwp(i) + xfiwc(i, k)*rhodz(i, k)
    397431
    398     END DO
    399   END DO
     432    ENDDO
     433  ENDDO
    400434
    401435  ! --if cloud droplet radius is fixed, then pcldtaupi=pcltau
    … …  
    406440        pcldtaupi(i, k) = pcltau(i, k)
    407441        reice_pi(i, k) = reice(i, k)
    408       END DO
    409     END DO
    410   END IF
     442      ENDDO
     443    ENDDO
     444  ENDIF
    411445
    412446  DO k = 1, klev
    … …  
    415449      reliq_pi(i, k) = rad_chaud_pi(i, k)
    416450      reice_pi(i, k) = reice(i, k)
    417     END DO
    418   END DO
     451    ENDDO
     452  ENDDO
    419453
    420454  ! COMPUTE CLOUD LIQUID PATH AND TOTAL CLOUDINESS
    … …  
    432466    pcl(i) = 1.0
    433467    pctlwp(i) = 0.0
    434   END DO
     468  ENDDO
    435469
    436470  ! --calculation of liquid water path
    … …  
    439473    DO i = 1, klon
    440474      pctlwp(i) = pctlwp(i) + pqlwp(i, k)*rhodz(i, k)
    441     END DO
    442   END DO
     475    ENDDO
     476  ENDDO
    443477
    444478  ! --calculation of cloud properties with cloud overlap
    … …  
    462496            (i),kind=8),1.-zepsec))
    463497          zcloudl(i) = pclc(i, k)
    464         END IF
     498        ENDIF
    465499        zcloud(i) = pclc(i, k)
    466       END DO
    467     END DO
     500      ENDDO
     501    ENDDO
    468502  ELSE IF (novlp==2) THEN
    469503    DO k = klev, 1, -1
    … …  
    477511        ELSE IF (paprs(i,k)>=prlmc) THEN
    478512          pcl(i) = min(pclc(i,k), pcl(i))
    479         END IF
    480       END DO
    481     END DO
     513        ENDIF
     514      ENDDO
     515    ENDDO
    482516  ELSE IF (novlp==3) THEN
    483517    DO k = klev, 1, -1
    … …  
    491525        ELSE IF (paprs(i,k)>=prlmc) THEN
    492526          pcl(i) = pcl(i)*(1.0-pclc(i,k))
    493         END IF
    494       END DO
    495     END DO
    496   END IF
     527        ENDIF
     528      ENDDO
     529    ENDDO
     530  ENDIF
    497531
    498532  DO i = 1, klon
    … …  
    500534    pcm(i) = 1. - pcm(i)
    501535    pcl(i) = 1. - pcl(i)
    502   END DO
     536  ENDDO
    503537
    504538  ! ========================================================
    … …  
    521555        ELSE
    522556          lcc3d(i, k) = pclc(i, k)*phase3d(i, k)
    523         END IF
     557        ENDIF
    524558        scdnc(i, k) = lcc3d(i, k)*cdnc(i, k) ! m-3
    525       END DO
    526     END DO
     559      ENDDO
     560    ENDDO
    527561
    528562    DO i = 1, klon
    … …  
    532566      IF (novlp.EQ.3 .OR. novlp.EQ.1) tcc(i) = 1.
    533567      IF (novlp.EQ.2) tcc(i) = 0.
    534     END DO
     568    ENDDO
    535569
    536570    DO i = 1, klon
    … …  
    546580              WRITE (*, *) 'Hypothese de recouvrement: MAXIMUM'
    547581              first = .FALSE.
    548             END IF
     582            ENDIF
    549583            flag_max = -1.
    550584            ftmp(i) = max(tcc(i), pclc(i,k))
    551           END IF
     585          ENDIF
    552586
    553587          IF (novlp.EQ.3) THEN
    … …  
    555589              WRITE (*, *) 'Hypothese de recouvrement: RANDOM'
    556590              first = .FALSE.
    557             END IF
     591            ENDIF
    558592            flag_max = 1.
    559593            ftmp(i) = tcc(i)*(1-pclc(i,k))
    560           END IF
     594          ENDIF
    561595
    562596          IF (novlp.EQ.1) THEN
    … …  
    566600                &                                          RANDOM'
    567601              first = .FALSE.
    568             END IF
     602            ENDIF
    569603            flag_max = 1.
    570604            ftmp(i) = tcc(i)*(1.-max(pclc(i,k),pclc(i,k+1)))/(1.-min(pclc(i, &
    571605              k+1),1.-thres_neb))
    572           END IF
     606          ENDIF
    573607          ! Effective radius of cloud droplet at top of cloud (m)
    574608          reffclwtop(i) = reffclwtop(i) + rad_chaud(i, k)*1.0E-06*phase3d(i, &
    … …  
    582616          tcc(i) = ftmp(i)
    583617
    584         END IF ! is there a visible, not-too-small cloud?
    585       END DO ! loop over k
     618        ENDIF ! is there a visible, not-too-small cloud?
     619      ENDDO ! loop over k
    586620
    587621      IF (novlp.EQ.3 .OR. novlp.EQ.1) tcc(i) = 1. - tcc(i)
    588622
    589     END DO ! loop over i
     623    ENDDO ! loop over i
    590624
    591625    ! ! Convective and Stratiform Cloud Droplet Effective Radius (REFFCLWC
    … …  
    615649        reffclws(i, k) = radius
    616650        reffclws(i, k) = reffclws(i, k)*lcc3dstra(i, k)
    617       END DO !klev
    618     END DO !klon
     651      ENDDO !klev
     652    ENDDO !klon
    619653
    620654    ! Column Integrated Cloud Droplet Number (CLDNVI) : variable 2D
    … …  
    628662        lcc_integrat(i) = lcc_integrat(i) + lcc3d(i, k)*dh(i, k)
    629663        height(i) = height(i) + dh(i, k)
    630       END DO ! klev
     664      ENDDO ! klev
    631665      lcc_integrat(i) = lcc_integrat(i)/height(i)
    632666      IF (lcc_integrat(i)<=1.0E-03) THEN
    … …  
    634668      ELSE
    635669        cldnvi(i) = cldnvi(i)*lcc(i)/lcc_integrat(i)
    636       END IF
    637     END DO ! klon
     670      ENDIF
     671    ENDDO ! klon
    638672
    639673    DO i = 1, klon
    … …  
    649683        IF (icc3dstra(i,k)<=0.0) icc3dstra(i, k) = 0.0
    650684!FC (CAUSES)
    651       END DO
     685      ENDDO
    652686      IF (reffclwtop(i)<=0.0) reffclwtop(i) = 0.0
    653687      IF (cldncl(i)<=0.0) cldncl(i) = 0.0
    654688      IF (cldnvi(i)<=0.0) cldnvi(i) = 0.0
    655689      IF (lcc(i)<=0.0) lcc(i) = 0.0
    656     END DO
    657 
    658   END IF !ok_cdnc
     690    ENDDO
     691
     692  ENDIF !ok_cdnc
    659693
    660694  first=.false. !to be sure
Note: See TracChangeset for help on using the changeset viewer.