Changeset 517 for LMDZ.3.3


Ignore:
Timestamp:
Apr 16, 2004, 5:43:38 PM (20 years ago)
Author:
lmdzadmin
Message:

Inclusion des modifications de O. Boucher et de J. Quaas pour le calcul des
premiers effets directs et indirects dus aux aerosols
LF

Location:
LMDZ.3.3/branches/rel-LF/libf/phylmd
Files:
8 edited

Legend:

Unmodified
Added
Removed

  • LMDZ.3.3/branches/rel-LF/libf/phylmd/conf_phys.F90

    r498 r517  
    55  subroutine conf_phys(ocean, ok_veget, ok_journe, ok_mensuel, ok_instan, &
    66 &                     fact_cldcon, facttemps,ok_newmicro,iflag_cldcon, &
    7  &                     ratqsbas,ratqshaut,if_ebil)
     7 &                     ratqsbas,ratqshaut,if_ebil, &
     8 &                     ok_ade, ok_aie, &
     9 &                     bl95_b0, bl95_b1)
    810
    911   use IOIPSL
    … …  
    2931! ok_mensuel: sorties mensuelles
    3032! ok_instan:  sorties instantanees
    31 
     33! ok_ade, ok_aie: apply or not aerosol direct and indirect effects
     34! bl95_b*: parameters in the formula to link CDNC to aerosol mass conc
     35!
    3236
    3337
    … …  
    3539  character (len = 6)  :: ocean
    3640  logical              :: ok_veget, ok_newmicro
    37   logical              :: ok_journe, ok_mensuel, ok_instan
     41  logical              :: ok_journe, ok_mensuel, ok_instan       
     42  LOGICAL              :: ok_ade, ok_aie
     43  REAL                 :: bl95_b0, bl95_b1
    3844  real                 :: fact_cldcon, facttemps,ratqsbas,ratqshaut
    3945  integer              :: iflag_cldcon, if_ebil
    … …  
    8995  ok_instan = .false.
    9096  call getin('OK_instan', ok_instan)
     97!
     98!Config Key  = ok_ade
     99!Config Desc = Aerosol direct effect or not?
     100!Config Def  = .false.
     101!Config Help = Used in radlwsw.F
     102!
     103  ok_ade = .false.
     104  call getin('ok_ade', ok_ade)
     105
     106!
     107!Config Key  = ok_aie
     108!Config Desc = Aerosol indirect effect or not?
     109!Config Def  = .false.
     110!Config Help = Used in nuage.F and radlwsw.F
     111!
     112  ok_aie = .false.
     113  call getin('ok_aie', ok_aie)
     114
     115!
     116!Config Key  = bl95_b0
     117!Config Desc = Parameter in CDNC-maer link (Boucher&Lohmann 1995)
     118!Config Def  = .false.
     119!Config Help = Used in nuage.F
     120!
     121  bl95_b0 = 2.
     122  call getin('bl95_b0', bl95_b0)
     123
     124!Config Key  = bl95_b1
     125!Config Desc = Parameter in CDNC-maer link (Boucher&Lohmann 1995)
     126!Config Def  = .false.
     127!Config Help = Used in nuage.F
     128!
     129  bl95_b1 = 0.2
     130  call getin('bl95_b1', bl95_b1)
     131
     132!
    91133!
    92134!Config Key  = if_ebil
    … …  
    554596  write(numout,*)' ksta_ter = ',ksta_ter
    555597  write(numout,*)' ok_kzmin = ',ok_kzmin
     598  write(numout,*)' ok_ade = ',ok_ade
     599  write(numout,*)' ok_aie = ',ok_aie
     600  write(numout,*)' bl95_b0 = ',bl95_b0
     601  write(numout,*)' bl95_b1 = ',bl95_b1
    556602  write(numout,*)' lev_histhf = ',lev_histhf
    557603  write(numout,*)' lev_histday = ',lev_histday

  • LMDZ.3.3/branches/rel-LF/libf/phylmd/fisrtilp_tr.F

    r230 r517  
    55     s                   pfrac_impa, pfrac_nucl, pfrac_1nucl,
    66     s                   frac_impa, frac_nucl,
    7      s                   prfl, psfl)
     7     s                   prfl, psfl,
     8     s                   RHcl) ! relative humidity in clear sky (needed for aer optical properties; aeropt.F)
    89
    910c
    … …  
    3738      REAL prfl(klon,klev+1) ! flux d'eau precipitante aux interfaces (kg/m2/s)
    3839      REAL psfl(klon,klev+1) ! flux d'eau precipitante aux interfaces (kg/m2/s)
     40     
     41Cjq   For aerosol opt properties needed (see aeropt.F)
     42      REAL RHcl(klon,klev)
     43     
    3944cAA
    4045c Coeffients de fraction lessivee : pour OFF-LINE
    … …  
    290295            rneb(i,k) = MAX(0.0,MIN(1.0,rneb(i,k)))
    291296            zcond(i) = MAX(0.0,zqn(i)-zqs(i))*rneb(i,k)/(1.+zdqs(i))
     297           
     298c--Olivier
     299            RHcl(i,k)=(zqs(i)+zq(i)-zdelq)/2./zqs(i)
     300            IF (rneb(i,k) .LE. 0.0) RHcl(i,k)=zq(i)/zqs(i)
     301            IF (rneb(i,k) .GE. 1.0) RHcl(i,k)=1.0
     302c--fin
     303           
    292304         ENDDO
    293305      ELSE

  • LMDZ.3.3/branches/rel-LF/libf/phylmd/ini_histmth.h

    r511 r517  
    277277     .                "ave(X)", zsto,zout)
    278278c
     279c Effets des aerosols
     280c
     281c     IF (ok_ade.OR.ok_aie) THEN
     282         CALL histdef(nid_mth, "topsad", "ADE at TOA", "W/m2",
     283     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
     284     .                "ave(X)", zsto,zout)
     285c
     286         CALL histdef(nid_mth, "solsad", "ADE at sfc", "W/m2",
     287     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
     288     .                "ave(X)", zsto,zout)
     289c
     290         CALL histdef(nid_mth, "topsai", "AIE at TOA", "W/m2",
     291     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
     292     .                "ave(X)", zsto,zout)
     293c
     294         CALL histdef(nid_mth, "solsai", "AIE at sfc", "W/m2",
     295     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
     296     .                "ave(X)", zsto,zout)
     297c     endif
     298c
     299
     300c
    279301c          CALL histdef(nid_mth, "frtu", "Zonal wind stress", "Pa",
    280302c    .                iim,jjmp1,nhori, 1,1,1, -99, 32,
    … …  
    693715     .                "ave(X)", zsto,zout)
    694716         ENDIF
    695 C
     717c
     718c Effets des aerosols
     719c
     720c     IF (ok_ade.OR.ok_aie) THEN
     721         CALL histdef(nid_mth, "re", "CDR", "um",
     722     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
     723     .                "ave(X)", zsto,zout)
     724c
     725         CALL histdef(nid_mth, "redenom", "CDR denominator", "-",
     726     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
     727     .                "ave(X)", zsto,zout)
     728c
     729         CALL histdef(nid_mth, "tau", "cloud opt thickness", "-",
     730     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
     731     .                "ave(X)", zsto,zout)
     732c
     733         CALL histdef(nid_mth, "taupi", "cloud opt thickn. (pi)", "-",
     734     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
     735     .                "ave(X)", zsto,zout)
     736c     endif
     737c
     738         CALL histdef(nid_mth, "ozone", "Ozone concentration", "-",
     739     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
     740     .                "ave(X)", zsto,zout)
     741c
    696742         if (nqmax.GE.3) THEN
    697743         DO iq=1,nqmax-2

  • LMDZ.3.3/branches/rel-LF/libf/phylmd/newmicro.F

    r486 r517  
    11      SUBROUTINE newmicro (paprs, pplay,ok_newmicro,
    22     .                  t, pqlwp, pclc, pcltau, pclemi,
    3 cIM    .                  pch, pcl, pcm, pct, pctlwp)
    43     .                  pch, pcl, pcm, pct, pctlwp,
    5      .                  xflwp, xfiwp, xflwc, xfiwc)
    6 
     4     s                  xflwp, xfiwp, xflwc, xfiwc,
     5     e                  ok_aie,
     6     e                  sulfate, sulfate_pi,
     7     e                  bl95_b0, bl95_b1,
     8     s                  cldtaupi, re, fl)
    79      IMPLICIT none
    810c======================================================================
    … …  
    1517c pclc----input-R-couverture nuageuse pour le rayonnement (0 a 1)
    1618c
     19c ok_aie--input-L-apply aerosol indirect effect or not
     20c sulfate-input-R-sulfate aerosol mass concentration [um/m^3]
     21c sulfate_pi-input-R-dito, pre-industrial value
     22c bl95_b0-input-R-a parameter, may be varied for tests (s-sea, l-land)
     23c bl95_b1-input-R-a parameter, may be varied for tests (    -"-      )
     24c     
     25c cldtaupi-output-R-pre-industrial value of cloud optical thickness,
     26c                   needed for the diagnostics of the aerosol indirect
     27c                   radiative forcing (see radlwsw)
     28c re------output-R-Cloud droplet effective radius multiplied by fl [um]
     29c fl------output-R-Denominator to re, introduced to avoid problems in
     30c                  the averaging of the output. fl is the fraction of liquid
     31c                  water clouds within a grid cell           
    1732c pcltau--output-R-epaisseur optique des nuages
    1833c pclemi--output-R-emissivite des nuages (0 a 1)
    … …  
    6681      parameter (DF=1.66) ! diffusivity factor
    6782c sb --
    68 
     83cjq for the aerosol indirect effect
     84cjq introduced by Johannes Quaas (quaas@lmd.jussieu.fr), 27/11/2003
     85cjq     
     86      LOGICAL ok_aie            ! Apply AIE or not?
     87      LOGICAL ok_a1lwpdep       ! a1 LWP dependent?
     88     
     89      REAL sulfate(klon, klev)  ! sulfate aerosol mass concentration [ug m-3]
     90      REAL cdnc(klon, klev)     ! cloud droplet number concentration [m-3]
     91      REAL re(klon, klev)       ! cloud droplet effective radius [um]
     92      REAL sulfate_pi(klon, klev)  ! sulfate aerosol mass concentration [ug m-3] (pre-industrial value)
     93      REAL cdnc_pi(klon, klev)     ! cloud droplet number concentration [m-3] (pi value)
     94      REAL re_pi(klon, klev)       ! cloud droplet effective radius [um] (pi value)
     95     
     96      REAL fl(klon, klev)       ! xliq * rneb (denominator to re; fraction of liquid water clouds within the grid cell)
     97     
     98      REAL bl95_b0, bl95_b1     ! Parameter in B&L 95-Formula
     99     
     100      REAL cldtaupi(klon, klev) ! pre-industrial cloud opt thickness for diag
     101cjq-end   
    69102c
    70103c Calculer l'epaisseur optique et l'emmissivite des nuages
    … …  
    119152
    120153c for liquid water clouds:
     154         IF (ok_aie) THEN
     155            ! Formula "D" of Boucher and Lohmann, Tellus, 1995
     156            !             
     157            cdnc(i,k) = 10.**(bl95_b0+bl95_b1*
     158     .           log(MAX(sulfate(i,k),1.e-4))/log(10.))*1.e6 !-m-3
     159            ! Cloud droplet number concentration (CDNC) is restricted
     160            ! to be within [20, 1000 cm^3]
     161            !
     162            cdnc(i,k)=MIN(1000.e6,MAX(20.e6,cdnc(i,k)))
     163            !
     164            !
     165            cdnc_pi(i,k) = 10.**(bl95_b0+bl95_b1*
     166     .           log(MAX(sulfate_pi(i,k),1.e-4))/log(10.))*1.e6 !-m-3
     167            cdnc_pi(i,k)=MIN(1000.e6,MAX(20.e6,cdnc_pi(i,k)))
     168            !           
     169            !
     170            ! air density: pplay(i,k) / (RD * zT(i,k))
     171            ! factor 1.1: derive effective radius from volume-mean radius
     172            ! factor 1000 is the water density
     173            ! _chaud means that this is the CDR for liquid water clouds
     174            !
     175            rad_chaud =
     176     .           1.1 * ( (pqlwp(i,k) * pplay(i,k) / (RD * T(i,k)) ) 
     177     .               / (4./3. * RPI * 1000. * cdnc(i,k)) )**(1./3.)
     178            !
     179            ! Convert to um. CDR shall be at least 3 um.
     180            !
     181c           rad_chaud = MAX(rad_chaud*1.e6, 3.)
     182            rad_chaud = MAX(rad_chaud*1.e6, 5.)
     183           
     184            ! Pre-industrial cloud opt thickness
     185            !
     186            ! "radius" is calculated as rad_chaud above (plus the
     187            ! ice cloud contribution) but using cdnc_pi instead of
     188            ! cdnc.
     189            radius =
     190     .           1.1 * ( (pqlwp(i,k) * pplay(i,k) / (RD * T(i,k)) ) 
     191     .               / (4./3. * RPI * 1000. * cdnc_pi(i,k)) )**(1./3.)
     192            radius = MAX(radius*1.e6, 3.)
     193           
     194            tc = t(i,k)-273.15
     195            rei = 0.71*tc + 61.29
     196            if (tc.le.-81.4) rei = 3.5
     197            if (zflwp(i).eq.0.) radius = 1.
     198            if (zfiwp(i).eq.0. .or. rei.le.0.) rei = 1.
     199            cldtaupi(i,k) = 3.0/2.0 * zflwp(i) / radius
     200     .             + zfiwp(i) * (3.448e-03  + 2.431/rei)
     201         ENDIF                  ! ok_aie
     202         ! For output diagnostics
     203         !
     204         ! Cloud droplet effective radius [um]
     205         !
     206         ! we multiply here with f * xl (fraction of liquid water
     207         ! clouds in the grid cell) to avoid problems in the
     208         ! averaging of the output.
     209         ! In the output of IOIPSL, derive the real cloud droplet
     210         ! effective radius as re/fl
     211         !
     212         fl(i,k) = pclc(i,k)*(1.-zfice)           
     213         re(i,k) = rad_chaud*fl(i,k)
     214           
     215c-jq end         
     216         
    121217         rel = rad_chaud
    122 
    123218c for ice clouds: as a function of the ambiant temperature
    124219c [formula used by Iacobellis and Somerville (2000), with an
    … …  
    156251         IF (lo) pcltau(i,k) = 0.0
    157252         IF (lo) pclemi(i,k) = 0.0
     253         
     254         IF (lo) cldtaupi(i,k) = 0.0
     255         IF (.NOT.ok_aie) cldtaupi(i,k)=pcltau(i,k)           
    158256      ENDDO
    159257      ENDDO

  • LMDZ.3.3/branches/rel-LF/libf/phylmd/nuage.F

    r390 r517  
    11      SUBROUTINE nuage (paprs, pplay,
    22     .                  t, pqlwp, pclc, pcltau, pclemi,
    3      .                  pch, pcl, pcm, pct, pctlwp)
     3     .                  pch, pcl, pcm, pct, pctlwp,
     4     e                  ok_aie,
     5     e                  sulfate, sulfate_pi,
     6     e                  bl95_b0, bl95_b1,
     7     s                  cldtaupi, re, fl)
    48      IMPLICIT none
    59c======================================================================
    … …  
    1115c pqlwp---input-R-eau liquide nuageuse dans l'atmosphere (kg/kg)
    1216c pclc----input-R-couverture nuageuse pour le rayonnement (0 a 1)
     17c ok_aie--input-L-apply aerosol indirect effect or not
     18c sulfate-input-R-sulfate aerosol mass concentration [um/m^3]
     19c sulfate_pi-input-R-dito, pre-industrial value
     20c bl95_b0-input-R-a parameter, may be varied for tests (s-sea, l-land)
     21c bl95_b1-input-R-a parameter, may be varied for tests (    -"-      )
     22c     
     23c cldtaupi-output-R-pre-industrial value of cloud optical thickness,
     24c                   needed for the diagnostics of the aerosol indirect
     25c                   radiative forcing (see radlwsw)
     26c re------output-R-Cloud droplet effective radius multiplied by fl [um]
     27c fl------output-R-Denominator to re, introduced to avoid problems in
     28c                  the averaging of the output. fl is the fraction of liquid
     29c                  water clouds within a grid cell     
    1330c
    1431c pcltau--output-R-epaisseur optique des nuages
    … …  
    2037#include "dimensions.h"
    2138#include "dimphy.h"
    22 #include "nuage.h"
    2339      REAL paprs(klon,klev+1), pplay(klon,klev)
    2440      REAL t(klon,klev)
    … …  
    3854      REAL zflwp, zradef, zfice, zmsac
    3955c
    40       REAL radius, rad_chaud
    41 ccc      PARAMETER (rad_chau1=13.0, rad_chau2=9.0, rad_froid=35.0)
     56      REAL radius, rad_froid, rad_chaud, rad_chau1, rad_chau2
     57      PARAMETER (rad_chau1=13.0, rad_chau2=9.0, rad_froid=35.0)
    4258ccc      PARAMETER (rad_chaud=15.0, rad_froid=35.0)
    4359c sintex initial      PARAMETER (rad_chaud=10.0, rad_froid=30.0)
    … …  
    4864      INTEGER nexpo ! exponentiel pour glace/eau
    4965      PARAMETER (nexpo=6)
     66     
     67cjq for the aerosol indirect effect
     68cjq introduced by Johannes Quaas (quaas@lmd.jussieu.fr), 27/11/2003
     69cjq     
     70      LOGICAL ok_aie            ! Apply AIE or not?
     71     
     72      REAL sulfate(klon, klev)  ! sulfate aerosol mass concentration [ug m-3]
     73      REAL cdnc(klon, klev)     ! cloud droplet number concentration [m-3]
     74      REAL re(klon, klev)       ! cloud droplet effective radius [um]
     75      REAL sulfate_pi(klon, klev)  ! sulfate aerosol mass concentration [ug m-3] (pre-industrial value)
     76      REAL cdnc_pi(klon, klev)     ! cloud droplet number concentration [m-3] (pi value)
     77      REAL re_pi(klon, klev)       ! cloud droplet effective radius [um] (pi value)
     78     
     79      REAL fl(klon, klev)  ! xliq * rneb (denominator to re; fraction of liquid water clouds within the grid cell)
     80     
     81      REAL bl95_b0, bl95_b1     ! Parameter in B&L 95-Formula
     82     
     83      REAL cldtaupi(klon, klev) ! pre-industrial cloud opt thickness for diag
     84cjq-end     
     85     
    5086ccc      PARAMETER (nexpo=1)
    5187c
    … …  
    5692         rad_chaud = rad_chau1
    5793         IF (k.LE.3) rad_chaud = rad_chau2
     94           
    5895         pclc(i,k) = MAX(pclc(i,k), seuil_neb)
    5996         zflwp = 1000.*pqlwp(i,k)/RG/pclc(i,k)
    … …  
    6299         zfice = MIN(MAX(zfice,0.0),1.0)
    63100         zfice = zfice**nexpo
     101         
     102         IF (ok_aie) THEN
     103            ! Formula "D" of Boucher and Lohmann, Tellus, 1995
     104            !             
     105            cdnc(i,k) = 10.**(bl95_b0+bl95_b1*
     106     .           log(MAX(sulfate(i,k),1.e-4))/log(10.))*1.e6 !-m-3
     107            ! Cloud droplet number concentration (CDNC) is restricted
     108            ! to be within [20, 1000 cm^3]
     109            !
     110            cdnc(i,k)=MIN(1000.e6,MAX(20.e6,cdnc(i,k)))
     111            cdnc_pi(i,k) = 10.**(bl95_b0+bl95_b1*
     112     .           log(MAX(sulfate_pi(i,k),1.e-4))/log(10.))*1.e6 !-m-3
     113            cdnc_pi(i,k)=MIN(1000.e6,MAX(20.e6,cdnc_pi(i,k)))
     114            !           
     115            !
     116            ! air density: pplay(i,k) / (RD * zT(i,k))
     117            ! factor 1.1: derive effective radius from volume-mean radius
     118            ! factor 1000 is the water density
     119            ! _chaud means that this is the CDR for liquid water clouds
     120            !
     121            rad_chaud =
     122     .           1.1 * ( (pqlwp(i,k) * pplay(i,k) / (RD * T(i,k)) ) 
     123     .               / (4./3. * RPI * 1000. * cdnc(i,k)) )**(1./3.)
     124            !
     125            ! Convert to um. CDR shall be at least 3 um.
     126            !
     127            rad_chaud = MAX(rad_chaud*1.e6, 3.)
     128           
     129            ! For output diagnostics
     130            !
     131            ! Cloud droplet effective radius [um]
     132            !
     133            ! we multiply here with f * xl (fraction of liquid water
     134            ! clouds in the grid cell) to avoid problems in the
     135            ! averaging of the output.
     136            ! In the output of IOIPSL, derive the real cloud droplet
     137            ! effective radius as re/fl
     138            !
     139            fl(i,k) = pclc(i,k)*(1.-zfice)           
     140            re(i,k) = rad_chaud*fl(i,k)
     141           
     142            ! Pre-industrial cloud opt thickness
     143            !
     144            ! "radius" is calculated as rad_chaud above (plus the
     145            ! ice cloud contribution) but using cdnc_pi instead of
     146            ! cdnc.
     147            radius = MAX(1.1e6 * ( (pqlwp(i,k)*pplay(i,k)/(RD*T(i,k))) 
     148     .                / (4./3.*RPI*1000.*cdnc_pi(i,k)) )**(1./3.),
     149     .               3.) * (1.-zfice) + rad_froid * zfice           
     150            cldtaupi(i,k) = 3.0/2.0 * zflwp / radius
     151     .           
     152         ENDIF                  ! ok_aie
     153         
    64154         radius = rad_chaud * (1.-zfice) + rad_froid * zfice
    65155         coef = coef_chau * (1.-zfice) + coef_froi * zfice
    … …  
    70160         IF (lo) pcltau(i,k) = 0.0
    71161         IF (lo) pclemi(i,k) = 0.0
     162         
     163         IF (.NOT.ok_aie) cldtaupi(i,k)=pcltau(i,k)           
    72164      ENDDO
    73165      ENDDO

  • LMDZ.3.3/branches/rel-LF/libf/phylmd/physiq.F

    r515 r517  
    1 c
     1C
    22c $Header$
    33c
    … …  
    143143      LOGICAL ok_region ! sortir le fichier regional
    144144      PARAMETER (ok_region=.FALSE.)
     145c
     146c
     147      LOGICAL ok_polder ! sortir échantillonné de manière POLDER
     148      save ok_polder
    145149c======================================================================
    146150c
    … …  
    877881      REAL zx_lon(iim,jjmp1), zx_lat(iim,jjmp1)
    878882c
    879       INTEGER nid_day, nid_mth, nid_ins, nid_nmc
    880       SAVE nid_day, nid_mth, nid_ins, nid_nmc
     883      INTEGER nid_day, nid_mth, nid_ins, nid_nmc, nid_pol
     884      SAVE nid_day, nid_mth, nid_ins, nid_nmc, nid_pol
    881885c
    882886      INTEGER nhori, nvert
    … …  
    928932      REAL zu10m(klon), zv10m(klon)           !vents a 10m moyennes s/1 maille
    929933      CHARACTER*40 t2mincels, t2maxcels       !t2m min., t2m max
     934cjq   Aerosol effects (Johannes Quaas, 27/11/2003)
     935      REAL sulfate(klon, klev) ! SO4 aerosol concentration [ug/m3]
     936      REAL sulfate_pi(klon, klev) ! SO4 aerosol concentration [ug/m3] (pre-industrial value)
     937      SAVE sulfate_pi
     938
     939      REAL cldtaupi(klon,klev)  ! Cloud optical thickness for pre-industrial (pi) aerosols
     940
     941      REAL re(klon, klev)       ! Cloud droplet effective radius
     942      REAL fl(klon, klev)  ! denominator of re
     943
     944      REAL re_top(klon), fl_top(klon) ! CDR at top of liquid water clouds
     945
     946      ! Aerosol optical properties
     947      REAL tau_ae(klon,klev,2), piz_ae(klon,klev,2)
     948      REAL cg_ae(klon,klev,2)
     949
     950      REAL topswad(klon), solswad(klon) ! Aerosol direct effect.
     951      ! ok_ade=T -ADE=topswad-topsw
     952
     953      REAL topswai(klon), solswai(klon) ! Aerosol indirect effect.
     954      ! ok_aie=T ->
     955      !        ok_ade=T -AIE=topswai-topswad
     956      !        ok_ade=F -AIE=topswai-topsw
     957
     958      ! For POLDER swath
     959      INTEGER pyr, pmo, pday    ! Year, month and day
     960      INTEGER poldermask(klon)  ! POLDER swath mask (0 or 1)
     961
     962      REAL aerindex(klon)       ! POLDER aerosol index
     963     
     964      ! Parameters
     965      LOGICAL ok_ade, ok_aie    ! Apply aerosol (in)direct effects or not
     966      REAL bl95_b0, bl95_b1   ! Parameter in Boucher and Lohmann (1995)
     967cjq-end
    930968c
    931969c Declaration des constantes et des fonctions thermodynamiques
    … …  
    9911029         call conf_phys(ocean, ok_veget, ok_journe, ok_mensuel,
    9921030     .                  ok_instan, fact_cldcon, facttemps,ok_newmicro,
    993      .                  iflag_cldcon,ratqsbas,ratqshaut, if_ebil)
     1031     .                  iflag_cldcon,ratqsbas,ratqshaut, if_ebil,
     1032     .                  ok_ade, ok_aie,
     1033     .                  bl95_b0, bl95_b1)
    9941034cIM  .                  , RI0)
    9951035
    … …  
    11451185
    11461186#undef histmthNMC
    1147 #define histmthNMC
     1187cccccccc#define histmthNMC
    11481188#ifdef histmthNMC
    11491189#include "ini_histmthNMC.h"
    … …  
    11511191
    11521192#include "ini_histins.h"
     1193#include "ini_histpol.h"
    11531194
    11541195#ifdef histREGDYN
    … …  
    19121953     .            t_seri, convliq, convfra, dtau_c, dem_c,
    19131954     .            cldh_c, cldl_c, cldm_c, cldt_c, cldq_c,
    1914      .            flwp_c, fiwp_c, flwc_c, fiwc_c)
     1955     .            flwp_c, fiwp_c, flwc_c, fiwc_c,
     1956     e            ok_aie,
     1957     e            sulfate, sulfate_pi,
     1958     e            bl95_b0, bl95_b1,
     1959     s            cldtaupi, re, fl)
    19151960c
    19161961cIM calcul tau. emi nuages startiformes
    … …  
    19181963     .            t_seri, cldliq, cldfra, dtau_s, dem_s,
    19191964     .            cldh_s, cldl_s, cldm_s, cldt_s, cldq_s,
    1920      .            flwp_s, fiwp_s, flwc_s, fiwc_s)
     1965     .            flwp_s, fiwp_s, flwc_s, fiwc_s,
     1966     e            ok_aie,
     1967     e            sulfate, sulfate_pi,
     1968     e            bl95_b0, bl95_b1,
     1969     s            cldtaupi, re, fl)
    19211970c
    19221971      cldtot(:,:)=min(max(cldfra(:,:),rnebcon(:,:)),1.)
    … …  
    21852234      ENDDO
    21862235      ENDDO
    2187 c
     2236cjq - introduce the aerosol direct and first indirect radiative forcings
     2237cjq - Johannes Quaas, 27/11/2003 (quaas@lmd.jussieu.fr)
     2238      IF (ok_ade.OR.ok_aie) THEN
     2239         ! Get sulfate aerosol distribution
     2240         CALL readsulfate(rjourvrai, debut, sulfate)
     2241         CALL readsulfate_preind(rjourvrai, debut, sulfate_pi)
     2242
     2243         ! Calculate aerosol optical properties (Olivier Boucher)
     2244         CALL aeropt(pplay, paprs, t_seri, sulfate, rhcl,
     2245     .        tau_ae, piz_ae, cg_ae, aerindex)
     2246      ENDIF
     2247
     2248c     
    21882249c Calculer les parametres optiques des nuages et quelques
    21892250c parametres pour diagnostiques:
    … …  
    21932254     .            t_seri, cldliq, cldfra, cldtau, cldemi,
    21942255     .            cldh, cldl, cldm, cldt, cldq,
    2195      .            flwp, fiwp, flwc, fiwc)
     2256     .            flwp, fiwp, flwc, fiwc,
     2257     e            ok_aie,
     2258     e            sulfate, sulfate_pi,
     2259     e            bl95_b0, bl95_b1,
     2260     s            cldtaupi, re, fl)
    21962261      else
    21972262      CALL nuage (paprs, pplay,
    21982263     .            t_seri, cldliq, cldfra, cldtau, cldemi,
    2199      .            cldh, cldl, cldm, cldt, cldq)
     2264     .            cldh, cldl, cldm, cldt, cldq,
     2265     e            ok_aie,
     2266     e            sulfate, sulfate_pi,
     2267     e            bl95_b0, bl95_b1,
     2268     s            cldtaupi, re, fl)
     2269     
    22002270      endif
    22012271c
    … …  
    22332303     s             topsw0,toplw0,solsw0,sollw0,
    22342304     s             lwdn0, lwdn, lwup0, lwup,
    2235      s             swdn0, swdn, swup0, swup     )
     2305     s             swdn0, swdn, swup0, swup,
     2306     e             ok_ade, ok_aie, ! new for aerosol radiative effects
     2307     e             tau_ae, piz_ae, cg_ae, ! ="=
     2308     s             topswad, solswad, ! ="=
     2309     e             cldtaupi, ! ="=
     2310     s             topswai, solswai) ! ="=
    22362311      itaprad = 0
    22372312      ENDIF

  • LMDZ.3.3/branches/rel-LF/libf/phylmd/radlwsw.F

    r503 r517  
    1 cIM   SUBROUTINE radlwsw(dist, rmu0, fract, co2_ppm, solaire,
    21      SUBROUTINE radlwsw(dist, rmu0, fract,
    32     .                  paprs, pplay,tsol,albedo, alblw, t,q,wo,
    4      .                  cldfra, cldemi, cldtau,
     3     .                  cldfra, cldemi, cldtaupd,
    54     .                  heat,heat0,cool,cool0,radsol,albpla,
    65     .                  topsw,toplw,solsw,sollw,
    76     .                  sollwdown,
    8 cIM  .                  sollwdown, sollwdownclr,
    9 cIM  .                  toplwdown, toplwdownclr,
    107     .                  topsw0,toplw0,solsw0,sollw0,
    11 cIM BEG
    128     .                  lwdn0, lwdn, lwup0, lwup,
    13 cIM END
    14      .                  swdn0, swdn, swup0, swup    )
     9     .                  swdn0, swdn, swup0, swup,
     10     .                  ok_ade, ok_aie,
     11     .                  tau_ae, piz_ae, cg_ae,
     12     .                  topswad, solswad,
     13     .                  cldtaupi, topswai, solswai)
     14c     
    1515      IMPLICIT none
    1616c======================================================================
    … …  
    3131c wo-------input-R- contenu en ozone (en cm.atm)
    3232c cldfra---input-R- fraction nuageuse (entre 0 et 1)
    33 c cldtau---input-R- epaisseur optique des nuages dans le visible
     33c cldtaupd---input-R- epaisseur optique des nuages dans le visible (present-day value)
    3434c cldemi---input-R- emissivite des nuages dans l'IR (entre 0 et 1)
     35c ok_ade---input-L- apply the Aerosol Direct Effect or not?
     36c ok_aie---input-L- apply the Aerosol Indirect Effect or not?
     37c tau_ae, piz_ae, cg_ae-input-R- aerosol optical properties (calculated in aeropt.F)
     38c cldtaupi-input-R- epaisseur optique des nuages dans le visible
     39c                   calculated for pre-industrial (pi) aerosol concentrations, i.e. with smaller
     40c                   droplet concentration, thus larger droplets, thus generally cdltaupi cldtaupd
     41c                   it is needed for the diagnostics of the aerosol indirect radiative forcing     
    3542c
    3643c heat-----output-R- echauffement atmospherique (visible) (K/jour)
    … …  
    4249c solsw----output-R- flux solaire net a la surface
    4350c sollw----output-R- ray. IR montant a la surface
     51c solswad---output-R- ray. solaire net absorbe a la surface (aerosol dir)
     52c topswad---output-R- ray. solaire absorbe au sommet de l'atm. (aerosol dir)
     53c solswai---output-R- ray. solaire net absorbe a la surface (aerosol ind)
     54c topswai---output-R- ray. solaire absorbe au sommet de l'atm. (aerosol ind)
     55c
     56c ATTENTION: swai and swad have to be interpreted in the following manner:
     57c ---------
     58c ok_ade=F & ok_aie=F -both are zero
     59c ok_ade=T & ok_aie=F -aerosol direct forcing is F_{AD} = topsw-topswad
     60c                        indirect is zero
     61c ok_ade=F & ok_aie=T -aerosol indirect forcing is F_{AI} = topsw-topswai
     62c                        direct is zero
     63c ok_ade=T & ok_aie=T -aerosol indirect forcing is F_{AI} = topsw-topswai
     64c                        aerosol direct forcing is F_{AD} = topswai-topswad
     65c
     66     
    4467c======================================================================
    4568#include "dimensions.h"
    … …  
    5679      real albedo(klon), alblw(klon), tsol(klon)
    5780      real t(klon,klev), q(klon,klev), wo(klon,klev)
    58       real cldfra(klon,klev), cldemi(klon,klev), cldtau(klon,klev)
     81      real cldfra(klon,klev), cldemi(klon,klev), cldtaupd(klon,klev)
    5982      real heat(klon,klev), cool(klon,klev)
    6083      real heat0(klon,klev), cool0(klon,klev)
    … …  
    123146      REAL lwup(klon,kflev+1),lwup0(klon,kflev+1)
    124147cIM END
    125 c---------------------------------------------------------------
     148c-OB
     149cjq the following quantities are needed for the aerosol radiative forcings
     150
     151      real topswad(klon), solswad(klon) ! output: aerosol direct forcing at TOA and surface
     152      real topswai(klon), solswai(klon) ! output: aerosol indirect forcing atTOA and surface
     153      real tau_ae(klon,klev,2), piz_ae(klon,klev,2), cg_ae(klon,klev,2) ! aerosol optical properties (see aeropt.F)
     154      real cldtaupi(klon,klev)  ! cloud optical thickness for pre-industrial aerosol concentrations
     155                                ! (i.e., with a smaller droplet concentrationand thus larger droplet radii)
     156      logical ok_ade, ok_aie    ! switches whether to use aerosol direct (indirect) effects or not
     157      real*8 tauae(kdlon,kflev,2) ! aer opt properties
     158      real*8 pizae(kdlon,kflev,2)
     159      real*8 cgae(kdlon,kflev,2)
     160      REAL*8 PTAUA(kdlon,2,kflev) ! present-day value of cloud opt thickness (PTAU is pre-industrial value), local use
     161      REAL*8 POMEGAA(kdlon,2,kflev) ! dito for single scatt albedo
     162      REAL*8 ztopswad(kdlon), zsolswad(kdlon) ! Aerosol direct forcing at TOAand surface
     163      REAL*8 ztopswai(kdlon), zsolswai(kdlon) ! dito, indirect
     164cjq-end
     165     
     166c
     167c-------------------------------------------
    126168      nb_gr = klon / kdlon
    127169      IF (nb_gr*kdlon .NE. klon) THEN
    … …  
    202244         PCLDLU(i,k) = cldfra(iof+i,k)*cldemi(iof+i,k)
    203245         PCLDSW(i,k) = cldfra(iof+i,k)
    204          PTAU(i,1,k) = MAX(cldtau(iof+i,k), 1.0e-05)! 1e-12 serait instable
    205          PTAU(i,2,k) = MAX(cldtau(iof+i,k), 1.0e-05)! pour 32-bit machines
     246         PTAU(i,1,k) = MAX(cldtaupi(iof+i,k), 1.0e-05)! 1e-12 serait instable
     247         PTAU(i,2,k) = MAX(cldtaupi(iof+i,k), 1.0e-05)! pour 32-bit machines
    206248         POMEGA(i,1,k) = 0.9999 - 5.0e-04 * EXP(-0.5 * PTAU(i,1,k))
    207249         POMEGA(i,2,k) = 0.9988 - 2.5e-03 * EXP(-0.05 * PTAU(i,2,k))
    208250         PCG(i,1,k) = 0.865
    209251         PCG(i,2,k) = 0.910
     252c-OB
     253cjq Introduced for aerosol indirect forcings.
     254cjq The following values use the cloud optical thickness calculated from
     255cjq present-day aerosol concentrations whereas the quantities without the
     256cjq "A" at the end are for pre-industial (natural-only) aerosol concentrations
     257cjq
     258         PTAUA(i,1,k) = MAX(cldtaupd(iof+i,k), 1.0e-05)! 1e-12 serait instable
     259         PTAUA(i,2,k) = MAX(cldtaupd(iof+i,k), 1.0e-05)! pour 32-bit machines
     260         POMEGAA(i,1,k) = 0.9999 - 5.0e-04 * EXP(-0.5 * PTAUA(i,1,k))
     261         POMEGAA(i,2,k) = 0.9988 - 2.5e-03 * EXP(-0.05 * PTAUA(i,2,k))
     262cjq-end
    210263      ENDDO
    211264      ENDDO
    … …  
    222275         PAER(i,k,kk) = 1.0E-15
    223276      ENDDO
     277      ENDDO
     278      ENDDO
     279c-OB
     280      DO k = 1, kflev
     281      DO i = 1, kdlon
     282        tauae(i,k,1)=tau_ae(iof+i,k,1)
     283        pizae(i,k,1)=piz_ae(iof+i,k,1)
     284        cgae(i,k,1) =cg_ae(iof+i,k,1)
     285        tauae(i,k,2)=tau_ae(iof+i,k,2)
     286        pizae(i,k,2)=piz_ae(iof+i,k,2)
     287        cgae(i,k,2) =cg_ae(iof+i,k,2)
    224288      ENDDO
    225289      ENDDO
    … …  
    247311     S        zheat, zheat0,
    248312     S        zalbpla,ztopsw,zsolsw,ztopsw0,zsolsw0,
    249      S        ZFSUP,ZFSDN,ZFSUP0,ZFSDN0)
     313     S        ZFSUP,ZFSDN,ZFSUP0,ZFSDN0,
     314     S        tauae, pizae, cgae, ! aerosol optical properties
     315     s        PTAUA, POMEGAA,
     316     s        ztopswad,zsolswad,ztopswai,zsolswai, ! diagnosed aerosol forcing
     317     J        ok_ade, ok_aie) ! apply aerosol effects or not?
     318
    250319c======================================================================
    251320      DO i = 1, kdlon
    … …  
    292361c        swup  ( iof+i,2)   = ZFSUP  ( i,kflev + 1 )
    293362      ENDDO
     363cjq-transform the aerosol forcings, if they have
     364cjq to be calculated
     365      IF (ok_ade) THEN
     366      DO i = 1, kdlon
     367         topswad(iof+i) = ztopswad(i)
     368         solswad(iof+i) = zsolswad(i)
     369      ENDDO
     370      ELSE
     371      DO i = 1, kdlon
     372         topswad(iof+i) = 0.0
     373         solswad(iof+i) = 0.0
     374      ENDDO
     375      ENDIF
     376      IF (ok_aie) THEN
     377      DO i = 1, kdlon
     378         topswai(iof+i) = ztopswai(i)
     379         solswai(iof+i) = zsolswai(i)
     380      ENDDO
     381      ELSE
     382      DO i = 1, kdlon
     383         topswai(iof+i) = 0.0
     384         solswai(iof+i) = 0.0
     385      ENDDO
     386      ENDIF
     387cjq-end
    294388      DO k = 1, kflev
    295389c      DO i = 1, kdlon
    … …  
    321415     S              PHEAT, PHEAT0,
    322416     S              PALBPLA,PTOPSW,PSOLSW,PTOPSW0,PSOLSW0,
    323      S              ZFSUP,ZFSDN,ZFSUP0,ZFSDN0)
     417     S              ZFSUP,ZFSDN,ZFSUP0,ZFSDN0,
     418     S              tauae, pizae, cgae,
     419     s              PTAUA, POMEGAA,
     420     S              PTOPSWAD,PSOLSWAD,PTOPSWAI,PSOLSWAI,
     421     J              ok_ade, ok_aie )
     422     
    324423      IMPLICIT none
    325424
    … …  
    358457C        ORIGINAL : 89-07-14
    359458C        95-01-01   J.-J. MORCRETTE  Direct/Diffuse Albedo
     459c        03-11-27   J. QUAAS Introduce aerosol forcings (based on BOUCHER)
    360460C     ------------------------------------------------------------------
    361461C
    … …  
    426526      DATA itapsw /0/
    427527      DATA appel1er /.TRUE./
     528cjq-Introduced for aerosol forcings
     529      real*8 flag_aer
     530      logical ok_ade, ok_aie    ! use aerosol forcings or not?
     531      real*8 tauae(kdlon,kflev,2)  ! aerosol optical properties
     532      real*8 pizae(kdlon,kflev,2)  ! (see aeropt.F)
     533      real*8 cgae(kdlon,kflev,2)   ! -"-
     534      REAL*8 PTAUA(KDLON,2,KFLEV)    ! CLOUD OPTICAL THICKNESS (pre-industrial value)
     535      REAL*8 POMEGAA(KDLON,2,KFLEV)  ! SINGLE SCATTERING ALBEDO
     536      REAL*8 PTOPSWAD(KDLON)     ! SHORTWAVE FLUX AT T.O.A.(+AEROSOL DIR)
     537      REAL*8 PSOLSWAD(KDLON)     ! SHORTWAVE FLUX AT SURFACE(+AEROSOL DIR)
     538      REAL*8 PTOPSWAI(KDLON)     ! SHORTWAVE FLUX AT T.O.A.(+AEROSOL IND)
     539      REAL*8 PSOLSWAI(KDLON)     ! SHORTWAVE FLUX AT SURFACE(+AEROSOL IND)
     540cjq - Fluxes including aerosol effects
     541      REAL*8 ZFSUPAD(KDLON,KFLEV+1)
     542      REAL*8 ZFSDNAD(KDLON,KFLEV+1)
     543      REAL*8 ZFSUPAI(KDLON,KFLEV+1)
     544      REAL*8 ZFSDNAI(KDLON,KFLEV+1)
     545      SAVE ZFSUPAD, ZFSDNAD, ZFSUPAI, ZFSDNAI ! aerosol fluxes
     546cjq-end
     547     
    428548c
    429549      IF (appel1er) THEN
    … …  
    451571      INU = 1
    452572      CALL SW1S(INU,
    453      S     PAER, PALBD, PALBP, PCG, ZCLD, ZCLEAR, ZCLDSW0,
     573     S     PAER, flag_aer, tauae, pizae, cgae,
     574     S     PALBD, PALBP, PCG, ZCLD, ZCLEAR, ZCLDSW0,
    454575     S     ZDSIG, POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD,
    455576     S     ZFD, ZFU)
    456577      INU = 2
    457578      CALL SW2S(INU,
    458      S     PAER, ZAKI, PALBD, PALBP, PCG, ZCLD, ZCLEAR, ZCLDSW0,
     579     S     PAER, flag_aer, tauae, pizae, cgae,
     580     S     ZAKI, PALBD, PALBP, PCG, ZCLD, ZCLEAR, ZCLDSW0,
    459581     S     ZDSIG, POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD,
    460582     S     PWV, PQS,
    … …  
    466588      ENDDO
    467589      ENDDO
    468 c cloudy-sky:
    469 cIM ctes ds clesphys.h   CALL SWU(PSCT,RCO2,PCLDSW,PPMB,PPSOL,
     590     
     591      flag_aer=0.0
    470592      CALL SWU(PSCT,PCLDSW,PPMB,PPSOL,
    471593     S         PRMU0,PFRAC,PTAVE,PWV,
    … …  
    473595      INU = 1
    474596      CALL SW1S(INU,
    475      S     PAER, PALBD, PALBP, PCG, ZCLD, ZCLEAR, PCLDSW,
     597     S     PAER, flag_aer, tauae, pizae, cgae,
     598     S     PALBD, PALBP, PCG, ZCLD, ZCLEAR, PCLDSW,
    476599     S     ZDSIG, POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD,
    477600     S     ZFD, ZFU)
    478601      INU = 2
    479602      CALL SW2S(INU,
    480      S     PAER, ZAKI, PALBD, PALBP, PCG, ZCLD, ZCLEAR, PCLDSW,
     603     S     PAER, flag_aer, tauae, pizae, cgae,
     604     S     ZAKI, PALBD, PALBP, PCG, ZCLD, ZCLEAR, PCLDSW,
    481605     S     ZDSIG, POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD,
    482606     S     PWV, PQS,
    483607     S    ZFDOWN, ZFUP)
     608
     609c cloudy-sky:
     610     
    484611      DO JK = 1 , KFLEV+1
    485612      DO JL = 1, KDLON
    … …  
    488615      ENDDO
    489616      ENDDO
     617     
     618c     
     619      IF (ok_ade) THEN
    490620c
     621c cloudy-sky + aerosol dir OB
     622      flag_aer=1.0
     623      CALL SWU(PSCT,PCLDSW,PPMB,PPSOL,
     624     S         PRMU0,PFRAC,PTAVE,PWV,
     625     S         ZAKI,ZCLD,ZCLEAR,ZDSIG,ZFACT,ZRMU,ZSEC,ZUD)
     626      INU = 1
     627      CALL SW1S(INU,
     628     S     PAER, flag_aer, tauae, pizae, cgae,
     629     S     PALBD, PALBP, PCG, ZCLD, ZCLEAR, PCLDSW,
     630     S     ZDSIG, POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD,
     631     S     ZFD, ZFU)
     632      INU = 2
     633      CALL SW2S(INU,
     634     S     PAER, flag_aer, tauae, pizae, cgae,
     635     S     ZAKI, PALBD, PALBP, PCG, ZCLD, ZCLEAR, PCLDSW,
     636     S     ZDSIG, POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD,
     637     S     PWV, PQS,
     638     S    ZFDOWN, ZFUP)
     639      DO JK = 1 , KFLEV+1
     640      DO JL = 1, KDLON
     641         ZFSUPAD(JL,JK) = ZFSUP(JL,JK)
     642         ZFSDNAD(JL,JK) = ZFSDN(JL,JK)
     643         ZFSUP(JL,JK) = (ZFUP(JL,JK)   + ZFU(JL,JK)) * ZFACT(JL)
     644         ZFSDN(JL,JK) = (ZFDOWN(JL,JK) + ZFD(JL,JK)) * ZFACT(JL)
     645      ENDDO
     646      ENDDO
     647     
     648      ENDIF ! ok_ade
     649     
     650      IF (ok_aie) THEN
     651         
     652cjq   cloudy-sky + aerosol direct + aerosol indirect
     653      flag_aer=1.0
     654      CALL SWU(PSCT,PCLDSW,PPMB,PPSOL,
     655     S         PRMU0,PFRAC,PTAVE,PWV,
     656     S         ZAKI,ZCLD,ZCLEAR,ZDSIG,ZFACT,ZRMU,ZSEC,ZUD)
     657      INU = 1
     658      CALL SW1S(INU,
     659     S     PAER, flag_aer, tauae, pizae, cgae,
     660     S     PALBD, PALBP, PCG, ZCLD, ZCLEAR, PCLDSW,
     661     S     ZDSIG, POMEGAA, ZOZ, ZRMU, ZSEC, PTAUA, ZUD,
     662     S     ZFD, ZFU)
     663      INU = 2
     664      CALL SW2S(INU,
     665     S     PAER, flag_aer, tauae, pizae, cgae,
     666     S     ZAKI, PALBD, PALBP, PCG, ZCLD, ZCLEAR, PCLDSW,
     667     S     ZDSIG, POMEGAA, ZOZ, ZRMU, ZSEC, PTAUA, ZUD,
     668     S     PWV, PQS,
     669     S    ZFDOWN, ZFUP)
     670      DO JK = 1 , KFLEV+1
     671      DO JL = 1, KDLON
     672         ZFSUPAI(JL,JK) = ZFSUP(JL,JK)
     673         ZFSDNAI(JL,JK) = ZFSDN(JL,JK)         
     674         ZFSUP(JL,JK) = (ZFUP(JL,JK)   + ZFU(JL,JK)) * ZFACT(JL)
     675         ZFSDN(JL,JK) = (ZFDOWN(JL,JK) + ZFD(JL,JK)) * ZFACT(JL)
     676      ENDDO
     677      ENDDO
     678      ENDIF ! ok_aie     
     679cjq -end
     680     
    491681      itapsw = 0
    492682      ENDIF
    … …  
    512702         PSOLSW0(i) = ZFSDN0(i,1) - ZFSUP0(i,1)
    513703         PTOPSW0(i) = ZFSDN0(i,KFLEV+1) - ZFSUP0(i,KFLEV+1)
     704c-OB
     705         PSOLSWAD(i) = ZFSDNAD(i,1) - ZFSUPAD(i,1)
     706         PTOPSWAD(i) = ZFSDNAD(i,KFLEV+1) - ZFSUPAD(i,KFLEV+1)
     707c
     708         PSOLSWAI(i) = ZFSDNAI(i,1) - ZFSUPAI(i,1)
     709         PTOPSWAI(i) = ZFSDNAI(i,KFLEV+1) - ZFSUPAI(i,KFLEV+1)
     710c-fin
    514711      ENDDO
    515712C
    … …  
    707904      END
    708905      SUBROUTINE SW1S ( KNU
    709      S  ,  PAER  , PALBD , PALBP, PCG  , PCLD , PCLEAR, PCLDSW
     906     S  ,  PAER  , flag_aer, tauae, pizae, cgae
     907     S  ,  PALBD , PALBP, PCG  , PCLD , PCLEAR, PCLDSW
    710908     S  ,  PDSIG , POMEGA, POZ  , PRMU , PSEC , PTAU  , PUD 
    711909     S  ,  PFD   , PFU)
    … …  
    748946C
    749947      INTEGER KNU
     948c-OB
     949      real*8 flag_aer
     950      real*8 tauae(kdlon,kflev,2)
     951      real*8 pizae(kdlon,kflev,2)
     952      real*8 cgae(kdlon,kflev,2)
    750953      REAL*8 PAER(KDLON,KFLEV,5)
    751954      REAL*8 PALBD(KDLON,2)
    … …  
    8391042C
    8401043      CALL SWCLR ( KNU
    841      S  , PAER   , PALBP , PDSIG , ZRAYL, PSEC
     1044     S  , PAER   , flag_aer, tauae, pizae, cgae
     1045     S  , PALBP  , PDSIG , ZRAYL, PSEC
    8421046     S  , ZCGAZ  , ZPIZAZ, ZRAY1 , ZRAY2, ZREFZ, ZRJ0
    8431047     S  , ZRK0   , ZRMU0 , ZTAUAZ, ZTRA1, ZTRA2)
    … …  
    9391143      END
    9401144      SUBROUTINE SW2S ( KNU
    941      S  ,  PAER  ,PAKI, PALBD, PALBP, PCG   , PCLD, PCLEAR, PCLDSW
     1145     S  ,  PAER  , flag_aer, tauae, pizae, cgae
     1146     S  ,  PAKI, PALBD, PALBP, PCG   , PCLD, PCLEAR, PCLDSW
    9421147     S  ,  PDSIG ,POMEGA,POZ , PRMU , PSEC  , PTAU
    9431148     S  ,  PUD   ,PWV , PQS
    … …  
    9861191C
    9871192      INTEGER KNU
     1193c-OB
     1194      real*8 flag_aer
     1195      real*8 tauae(kdlon,kflev,2)
     1196      real*8 pizae(kdlon,kflev,2)
     1197      real*8 cgae(kdlon,kflev,2)
    9881198      REAL*8 PAER(KDLON,KFLEV,5)
    9891199      REAL*8 PAKI(KDLON,2)
    … …  
    11071317C
    11081318      CALL SWCLR ( KNU
    1109      S  , PAER   , PALBP , PDSIG , ZRAYL, PSEC
     1319     S  , PAER   , flag_aer, tauae, pizae, cgae
     1320     S  , PALBP  , PDSIG , ZRAYL, PSEC
    11101321     S  , ZCGAZ  , ZPIZAZ, ZRAY1 , ZRAY2, ZREFZ, ZRJ0
    11111322     S  , ZRK0   , ZRMU0 , ZTAUAZ, ZTRA1, ZTRA2)
    … …  
    14791690      END
    14801691      SUBROUTINE SWCLR  ( KNU
    1481      S  , PAER  , PALBP , PDSIG , PRAYL , PSEC
     1692     S  , PAER  , flag_aer, tauae, pizae, cgae
     1693     S  , PALBP , PDSIG , PRAYL , PSEC
    14821694     S  , PCGAZ , PPIZAZ, PRAY1 , PRAY2 , PREFZ , PRJ 
    14831695     S  , PRK   , PRMU0 , PTAUAZ, PTRA1 , PTRA2                   )
    … …  
    15121724C
    15131725      INTEGER KNU
     1726c-OB
     1727      real*8 flag_aer
     1728      real*8 tauae(kdlon,kflev,2)
     1729      real*8 pizae(kdlon,kflev,2)
     1730      real*8 cgae(kdlon,kflev,2)
    15141731      REAL*8 PAER(KDLON,KFLEV,5)
    15151732      REAL*8 PALBP(KDLON,2)
    … …  
    15761793C
    15771794      DO 108 JK = 1 , KFLEV
    1578       DO 104 JL = 1, KDLON
    1579       PCGAZ(JL,JK) = 0.
    1580       PPIZAZ(JL,JK) =  0.
    1581       PTAUAZ(JL,JK) = 0.
    1582  104  CONTINUE
    1583       DO 106 JAE=1,5
     1795c-OB
     1796c      DO 104 JL = 1, KDLON
     1797c      PCGAZ(JL,JK) = 0.
     1798c      PPIZAZ(JL,JK) =  0.
     1799c      PTAUAZ(JL,JK) = 0.
     1800c 104  CONTINUE
     1801c-OB
     1802c      DO 106 JAE=1,5
     1803c      DO 105 JL = 1, KDLON
     1804c      PTAUAZ(JL,JK)=PTAUAZ(JL,JK)
     1805c     S        +PAER(JL,JK,JAE)*TAUA(KNU,JAE)
     1806c      PPIZAZ(JL,JK)=PPIZAZ(JL,JK)+PAER(JL,JK,JAE)
     1807c     S        * TAUA(KNU,JAE)*RPIZA(KNU,JAE)
     1808c      PCGAZ(JL,JK) =  PCGAZ(JL,JK) +PAER(JL,JK,JAE)
     1809c     S        * TAUA(KNU,JAE)*RPIZA(KNU,JAE)*RCGA(KNU,JAE)
     1810c 105  CONTINUE
     1811c 106  CONTINUE
     1812c-OB
    15841813      DO 105 JL = 1, KDLON
    1585       PTAUAZ(JL,JK)=PTAUAZ(JL,JK)
    1586      S        +PAER(JL,JK,JAE)*TAUA(KNU,JAE)
    1587       PPIZAZ(JL,JK)=PPIZAZ(JL,JK)+PAER(JL,JK,JAE)
    1588      S        * TAUA(KNU,JAE)*RPIZA(KNU,JAE)
    1589       PCGAZ(JL,JK) =  PCGAZ(JL,JK) +PAER(JL,JK,JAE)
    1590      S        * TAUA(KNU,JAE)*RPIZA(KNU,JAE)*RCGA(KNU,JAE)
     1814      PTAUAZ(JL,JK)=flag_aer * tauae(JL,JK,KNU)
     1815      PPIZAZ(JL,JK)=flag_aer * pizae(JL,JK,KNU)
     1816      PCGAZ (JL,JK)=flag_aer * cgae(JL,JK,KNU)
    15911817 105  CONTINUE
    1592  106  CONTINUE
    1593 C
     1818C
     1819      IF (flag_aer.GT.0) THEN
     1820c-OB
    15941821      DO 107 JL = 1, KDLON
    1595       IF (KAER.NE.0) THEN
    1596          PCGAZ(JL,JK)=PCGAZ(JL,JK)/PPIZAZ(JL,JK)
    1597          PPIZAZ(JL,JK)=PPIZAZ(JL,JK)/PTAUAZ(JL,JK)
     1822c         PCGAZ(JL,JK)=PCGAZ(JL,JK)/PPIZAZ(JL,JK)
     1823c         PPIZAZ(JL,JK)=PPIZAZ(JL,JK)/PTAUAZ(JL,JK)
    15981824         ZTRAY = PRAYL(JL) * PDSIG(JL,JK)
    15991825         ZRATIO = ZTRAY / (ZTRAY + PTAUAZ(JL,JK))
    … …  
    16041830         PPIZAZ(JL,JK) =ZRATIO+(1.-ZRATIO)*PPIZAZ(JL,JK)*(1.-ZFF)
    16051831     S                       / (1. - PPIZAZ(JL,JK) * ZFF)
     1832 107  CONTINUE
    16061833      ELSE
     1834      DO JL = 1, KDLON
    16071835         ZTRAY = PRAYL(JL) * PDSIG(JL,JK)
    16081836         PTAUAZ(JL,JK) = ZTRAY
    16091837         PCGAZ(JL,JK) = 0.
    16101838         PPIZAZ(JL,JK) = 1.-REPSCT
    1611       END IF
    1612  107  CONTINUE
     1839      END DO
     1840      END IF   ! check flag_aer
     1841c     107  CONTINUE
    16131842c      PRINT 9107,JK,((PAER(JL,JK,JAE),JAE=1,5)
    16141843c     $ ,PTAUAZ(JL,JK),PPIZAZ(JL,JK),PCGAZ(JL,JK),JL=1,KDLON)

  • LMDZ.3.3/branches/rel-LF/libf/phylmd/write_histmth.h

    r506 r517  
    262262     $               iim*jjmp1,ndex2d)
    263263c
     264c
     265c effets des aerosols
     266c
     267c     IF (ok_ade.OR.ok_aie) THEN
     268      zx_tmp_fi2d(1:klon) = topswai(1:klon) - topswad(1:klon)
     269c      CALL gr_fi_ecrit(1, klon,iim,jjmp1, topswad,zx_tmp_2d)
     270      CALL gr_fi_ecrit(1, klon,iim,jjmp1,zx_tmp_fi2d ,zx_tmp_2d)
     271      CALL histwrite(nid_mth,"topsad",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
     272c
     273      zx_tmp_fi2d(1:klon) = solswai(1:klon) - solswad(1:klon)
     274c      CALL gr_fi_ecrit(1, klon,iim,jjmp1, solswad,zx_tmp_2d)
     275      CALL gr_fi_ecrit(1, klon,iim,jjmp1,zx_tmp_fi2d ,zx_tmp_2d)
     276      CALL histwrite(nid_mth,"solsad",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
     277c
     278      zx_tmp_fi2d(1:klon) = topsw(1:klon) - topswai(1:klon)
     279c      CALL gr_fi_ecrit(1, klon,iim,jjmp1, topswai,zx_tmp_2d)
     280      CALL gr_fi_ecrit(1, klon,iim,jjmp1,zx_tmp_fi2d ,zx_tmp_2d)
     281      CALL histwrite(nid_mth,"topsai",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
     282c
     283      zx_tmp_fi2d(1:klon) = solsw(1:klon) - solswai(1:klon)
     284c      CALL gr_fi_ecrit(1, klon,iim,jjmp1, solswai,zx_tmp_2d)
     285      CALL gr_fi_ecrit(1, klon,iim,jjmp1,zx_tmp_fi2d ,zx_tmp_2d)
     286      CALL histwrite(nid_mth,"solsai",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
     287c     endif
     288c
    264289      CALL gr_fi_ecrit(1, klon,iim,jjmp1, bils,zx_tmp_2d)
    265290      CALL histwrite(nid_mth,"bils",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
    … …  
    712737      ENDIF
    713738C
     739c
     740c effets des aerosols
     741c
     742c     IF (ok_ade.OR.ok_aie) THEN
     743      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, re, zx_tmp_3d)
     744      CALL histwrite(nid_mth,"re",itau_w,zx_tmp_3d,
     745     .                                   iim*jjmp1*klev,ndex3d)
     746c
     747      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, fl, zx_tmp_3d)
     748      CALL histwrite(nid_mth,"redenom",itau_w,zx_tmp_3d,
     749     .                                   iim*jjmp1*klev,ndex3d)
     750c
     751      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, cldtau, zx_tmp_3d)
     752      CALL histwrite(nid_mth,"tau",itau_w,zx_tmp_3d,
     753     .                                   iim*jjmp1*klev,ndex3d)
     754c
     755      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, cldtaupi, zx_tmp_3d)
     756      CALL histwrite(nid_mth,"taupi",itau_w,zx_tmp_3d,
     757     .                                   iim*jjmp1*klev,ndex3d)
     758c     endif
     759c
    714760      IF (nqmax.GE.3) THEN
    715761      DO iq=1,nqmax-2
Note: See TracChangeset for help on using the changeset viewer.