source: trunk/LMDZ.MARS/libf/phymars/updatereffrad_mod.F

Last change on this file was 2660, checked in by cmathe, 3 years ago

correction on rsedcloudco2/opacities/riceco2 computations, and more comments

File size: 12.7 KB
RevLine 
[1969]1      MODULE updatereffrad_mod
2     
3      IMPLICIT NONE
4     
5      CONTAINS
6     
[38]7      SUBROUTINE updatereffrad(ngrid,nlayer,
[2199]8     &                rdust,rstormdust,rtopdust,rice,nuice,
[2447]9     &                reffrad,nueffrad, riceco2, nuiceco2,
[2494]10     &                pq,tauscaling,tau,pplay, pt)
[1969]11       USE updaterad, ONLY: updaterdust, updaterice_micro,
[2447]12     &                      updaterice_microco2, updaterice_typ
[1036]13       use tracer_mod, only: nqmx, igcm_dust_mass, igcm_dust_number,
14     &                       igcm_h2o_ice, igcm_ccn_mass, radius,
[2447]15     &                       igcm_co2_ice, nuiceco2_ref,
16     &                       igcm_ccnco2_number, igcm_ccnco2_mass,
[2562]17     &                       igcm_ccnco2_h2o_number,
18     &                       igcm_ccnco2_h2o_mass_ice,
19     &                       igcm_ccnco2_h2o_mass_ccn,
[1036]20     &                       igcm_ccn_number, nuice_ref, varian,
[1974]21     &                       ref_r0, igcm_dust_submicron,
[2199]22     &                       igcm_stormdust_mass,igcm_stormdust_number,
[2494]23     &                       igcm_topdust_mass,igcm_topdust_number,
24     &                       rho_ice
[1246]25       USE dimradmars_mod, only: nueffdust,naerkind,
26     &            name_iaer,
27     &            iaer_dust_conrath,iaer_dust_doubleq,
[1974]28     &            iaer_dust_submicron,iaer_h2o_ice,
[2199]29     &            iaer_stormdust_doubleq,iaer_topdust_doubleq
[2409]30       use dust_param_mod, only: doubleq, active
[38]31       IMPLICIT NONE
32c=======================================================================
33c   subject:
34c   --------
35c   Subroutine designed to update the aerosol size distribution used by
36c     the radiative transfer scheme. This size distribution is assumed
37c     to be a log-normal distribution, with effective radius "reffrad" and
38c     variance "nueffrad".
39c   At firstcall, "rice" and "nuice" are not known, because
40c     the H2O ice microphysical scheme is called after the radiative
41c     transfer in physiq.F. That's why we assess the size of the
42c     water-ice particles at firstcall (see part 1.2 below).
43c
44c   author:   
45c   ------
46c   J.-B. Madeleine (2009-2010)
47c
48c=======================================================================
49c
50c    Declarations :
51c    -------------
52c
[1969]53      include "callkeys.h"
[38]54
55c-----------------------------------------------------------------------
[1974]56c     Inputs/outputs:
[38]57c     ------
58
[1974]59      INTEGER, INTENT(in) :: ngrid,nlayer
[38]60c     Ice geometric mean radius (m)
[1974]61      REAL, INTENT(out) :: rice(ngrid,nlayer)
[38]62c     Estimated effective variance of the size distribution (n.u.)
[1974]63      REAL, INTENT(out) :: nuice(ngrid,nlayer)
[38]64c     Tracer mass mixing ratio (kg/kg)
[1974]65      REAL, INTENT(in) :: pq(ngrid,nlayer,nqmx)
66      REAL, INTENT(out) :: rdust(ngrid,nlayer) ! Dust geometric mean radius (m)
[2199]67      REAL, INTENT(out) :: rstormdust(ngrid,nlayer) ! Dust geometric mean radius (m)   
68      REAL, INTENT(out) :: rtopdust(ngrid,nlayer) ! Dust geometric mean radius (m)
[1974]69      REAL, INTENT(in) :: pplay(ngrid,nlayer) ! altitude at the middle of the layers
70      REAL, INTENT(in) :: tau(ngrid,naerkind)
[38]71c     Aerosol effective radius used for radiative transfer (meter)
[1974]72      REAL, INTENT(out) :: reffrad(ngrid,nlayer,naerkind)
[38]73c     Aerosol effective variance used for radiative transfer (n.u.)
[1974]74      REAL, INTENT(out) :: nueffrad(ngrid,nlayer,naerkind)
75      REAL, INTENT(in) :: tauscaling(ngrid)         ! Convertion factor for qccn and Nccn
[2447]76c     CO2 ice mean radius (m)
[2459]77      double precision, INTENT(out) :: riceco2(ngrid,nlayer) ! co2 ice radius
[2447]78      REAL, INTENT(out) :: nuiceco2(ngrid,nlayer)
[2494]79      REAL, INTENT(in) :: pt(ngrid,nlayer) ! temperature
[2447]80     
[38]81c     Local variables:
82c     ---------------
83
84      INTEGER :: ig,l          ! 3D grid indices
85      INTEGER :: iaer          ! Aerosol index
86
87c     Number of cloud condensation nuclei near the surface
88c     (only used at firstcall). This value is taken from
89c     Montmessin et al. 2004 JGR 109 E10004 p5 (2E6 part m-3), and
90c     converted to part kg-1 using a typical atmospheric density.
91
[2660]92      REAL, PARAMETER :: ccn0 = 1.3E8, threshold = 1e-30 ! limit value
[629]93     
94c     For microphysics only:     
95      REAL Mo,No                       ! Mass and number of ccn
[1047]96      REAL rhocloud(ngrid,nlayer)  ! Cloud density (kg.m-3)
[2447]97c     For CO2 microphysics only:
98      REAL :: rhocloudco2(ngrid, nlayer) ! co2 cloud density
[38]99
[1224]100      LOGICAL,SAVE :: firstcall=.true.
[2660]101      REAL Nccnco2, Qccnco2, Niceco2
[2562]102      REAL Nccnco2_h2o, Qccnco2_h2o
[38]103      REAL CBRT
104      EXTERNAL CBRT
105
[2584]106
107!$OMP THREADPRIVATE(firstcall)
108
[38]109c==================================================================
[629]110c 1. Update radius from fields from dynamics or initial state
111c==================================================================
[38]112
[358]113c       1.1 Dust particles
114c       ------------------
115        IF (doubleq.AND.active) THEN
116          DO l=1,nlayer
117            DO ig=1, ngrid
[744]118              call updaterdust(pq(ig,l,igcm_dust_mass),
119     &                         pq(ig,l,igcm_dust_number),rdust(ig,l))
[358]120              nueffdust(ig,l) = exp(varian**2.)-1.
121             ENDDO
122           ENDDO
123        ELSE
124          DO l=1,nlayer
125            DO ig=1, ngrid
126              rdust(ig,l) = 0.8E-6
127              nueffdust(ig,l) = 0.3
128            ENDDO
[38]129          ENDDO
[358]130        ENDIF
[1974]131
132        ! updating radius of stormdust particles
133        IF (rdstorm.AND.active) THEN
134          DO l=1,nlayer
135            DO ig=1, ngrid
136              call updaterdust(pq(ig,l,igcm_stormdust_mass),
137     &                 pq(ig,l,igcm_stormdust_number),rstormdust(ig,l))
138              nueffdust(ig,l) = exp(varian**2.)-1.
139             ENDDO
140           ENDDO
141        ENDIF
[2199]142
143        ! updating radius of topdust particles
[2628]144        IF (topflows.AND.active) THEN
[2199]145          DO l=1,nlayer
146            DO ig=1, ngrid
147              call updaterdust(pq(ig,l,igcm_topdust_mass),
148     &                 pq(ig,l,igcm_topdust_number),rtopdust(ig,l))
149              nueffdust(ig,l) = exp(varian**2.)-1.
150             ENDDO
151           ENDDO
152        ENDIF
[629]153       
[358]154c       1.2 Water-ice particles
155c       -----------------------
[744]156
157        IF (water.AND.activice) THEN
158         IF (microphys) THEN
[1208]159       
160c    At firstcall, the true number and true mass of cloud condensation nuclei are not known.
161c    Indeed it is scaled on the prescribed dust opacity via a 'tauscaling' coefficient
162c    computed after radiative transfer. If tauscaling is not in startfi, we make an assumption for its value.
163
[744]164          IF (firstcall) THEN
[1974]165            !IF (minval(tauscaling).lt.0) tauscaling(:) = 1.e-3 ! default value when non-read in startfi is -1
166            !IF (freedust)                tauscaling(:) = 1.    ! if freedust, enforce no rescaling at all
[1208]167            firstcall = .false.
168          ENDIF
169 
170          DO l=1,nlayer
171            DO ig=1,ngrid
172              call updaterice_micro(pq(ig,l,igcm_h2o_ice),
173     &                              pq(ig,l,igcm_ccn_mass),
174     &                              pq(ig,l,igcm_ccn_number),
175     &                              tauscaling(ig),rice(ig,l),
176     &                              rhocloud(ig,l))
177              nuice(ig,l) = nuice_ref
[358]178            ENDDO
[1208]179          ENDDO
[744]180         
181        ELSE ! if not microphys
182         
183          DO l=1,nlayer
184            DO ig=1,ngrid   
185              call updaterice_typ(pq(ig,l,igcm_h2o_ice),
186     &                          tau(ig,1),pplay(ig,l),rice(ig,l))
187              nuice(ig,l) = nuice_ref
[629]188            ENDDO
[744]189          ENDDO
190 
191        ENDIF ! of if microphys
192       ENDIF ! of if (water.AND.activice)
[38]193
[2447]194
[2494]195c       1.3 CO2-ice particles
[2447]196
197        IF (co2clouds.AND.activeco2ice) THEN
198          DO l=1,nlayer
199            DO ig=1,ngrid
[2660]200              Niceco2 = max(pq(ig,l,igcm_co2_ice), threshold)
201              Nccnco2 = max(pq(ig,l,igcm_ccnco2_number),
202     &                              threshold)
203              Qccnco2 = max(pq(ig,l,igcm_ccnco2_mass),
204     &                              threshold)
[2562]205              Nccnco2_h2o = 0.
206              Qccnco2_h2o = 0.
207              if (co2useh2o) then
[2660]208                Nccnco2_h2o = max(pq(ig,l,igcm_ccnco2_h2o_number),
209     &                              threshold)
210                Qccnco2_h2o = max(pq(ig,l,igcm_ccnco2_h2o_mass_ice)
211     &                        + pq(ig,l,igcm_ccnco2_h2o_mass_ccn),
212     &                              threshold)
213
[2562]214                Nccnco2 = Nccnco2 - Nccnco2_h2o
215                Qccnco2 = Qccnco2 - Qccnco2_h2o
[2660]216                if (Nccnco2 <= 0) then
217                  Nccnco2 = threshold
218                  Qccnco2 = threshold
219                end if
[2562]220              end if
[2660]221              call updaterice_microco2(dble(Niceco2),
[2562]222     &                              dble(Qccnco2), dble(Nccnco2),
223     &                              dble(Qccnco2_h2o),
224     &                              dble(Nccnco2_h2o),
[2494]225     &                              pt(ig,l),
[2447]226     &                              tauscaling(ig),riceco2(ig,l),
227     &                              rhocloudco2(ig,l))
228              nuiceco2(ig,l) = nuiceco2_ref
229            END DO
230          ENDDO
231        ENDIF ! of if (co2clouds.AND.activeco2ice)
232
[38]233c==================================================================
234c 2. Radius used in the radiative transfer code (reffrad)
235c==================================================================
236
237      DO iaer = 1, naerkind ! Loop on aerosol kind
238        aerkind: SELECT CASE (name_iaer(iaer))
239c==================================================================
240        CASE("dust_conrath") aerkind         ! Typical dust profile
241c==================================================================
242          DO l=1,nlayer
243            DO ig=1,ngrid
[358]244              reffrad(ig,l,iaer) = rdust(ig,l) *
245     &          (1.e0 + nueffdust(ig,l))**2.5
[38]246              nueffrad(ig,l,iaer) = nueffdust(ig,l)
247            ENDDO
248          ENDDO
249c==================================================================
250        CASE("dust_doubleq") aerkind! Two-moment scheme for dust
251c==================================================================
252          DO l=1,nlayer
253            DO ig=1,ngrid
[358]254              reffrad(ig,l,iaer) = rdust(ig,l) * ref_r0
[38]255              nueffrad(ig,l,iaer) = nueffdust(ig,l)
256            ENDDO
257          ENDDO
258c==================================================================
259        CASE("dust_submicron") aerkind   ! Small dust population
260c==================================================================
261          DO l=1,nlayer
262            DO ig=1,ngrid
263              reffrad(ig,l,iaer)=radius(igcm_dust_submicron)
264              nueffrad(ig,l,iaer)=0.03
265            ENDDO
266          ENDDO     
267c==================================================================
268        CASE("h2o_ice") aerkind             ! Water ice crystals
269c==================================================================
270          DO l=1,nlayer
271            DO ig=1,ngrid
[358]272c             About reffice, do not confuse the mass mean radius
273c             (rayon moyen massique) and the number median radius
274c             (or geometric mean radius, rayon moyen géométrique).
275c             rice is a mass mean radius, whereas rdust
276c             is a geometric mean radius:
277c             number median rad = mass mean rad x exp(-1.5 sigma0^2)
278c             (Montmessin et al. 2004 paragraph 30). Therefore:
[38]279              reffrad(ig,l,iaer)=rice(ig,l)*(1.+nuice_ref)
280              nueffrad(ig,l,iaer)=nuice_ref
281            ENDDO
282          ENDDO
283c==================================================================
[2447]284        CASE("co2_ice") aerkind             ! CO2 ice crystals
285c==================================================================
286          DO l=1,nlayer
287            DO ig=1,ngrid
[2494]288              reffrad(ig,l,iaer)=real(riceco2(ig,l))*(1.+nuiceco2_ref)
[2447]289              nueffrad(ig,l,iaer)=nuiceco2_ref
290            ENDDO
291          ENDDO
292c==================================================================
[1974]293        CASE("stormdust_doubleq") aerkind! Two-moment scheme for
294c       stormdust; same distribution than normal dust
295c==================================================================
296          DO l=1,nlayer
297            DO ig=1,ngrid
298              reffrad(ig,l,iaer) = rstormdust(ig,l) * ref_r0
299              nueffrad(ig,l,iaer) = nueffdust(ig,l)
300            ENDDO
301          ENDDO
302c==================================================================
[2199]303        CASE("topdust_doubleq") aerkind! MV18: Two-moment scheme for
304c       topdust; same distribution than normal dust
305c==================================================================
306          DO l=1,nlayer
307            DO ig=1,ngrid
308              reffrad(ig,l,iaer) = rtopdust(ig,l) * ref_r0
309              nueffrad(ig,l,iaer) = nueffdust(ig,l)
310            ENDDO
311          ENDDO
312c==================================================================
[38]313        END SELECT aerkind
314      ENDDO ! iaer (loop on aerosol kind)
315
[1969]316      END SUBROUTINE updatereffrad
317     
318      END MODULE updatereffrad_mod
Note: See TracBrowser for help on using the repository browser.