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

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

correction on rsedcloudco2/opacities/riceco2 computations, and more comments
File size: 12.7 KB
Line 
1      MODULE updatereffrad_mod
2     
3      IMPLICIT NONE
4     
5      CONTAINS
6     
7      SUBROUTINE updatereffrad(ngrid,nlayer,
8     &                rdust,rstormdust,rtopdust,rice,nuice,
9     &                reffrad,nueffrad, riceco2, nuiceco2,
10     &                pq,tauscaling,tau,pplay, pt)
11       USE updaterad, ONLY: updaterdust, updaterice_micro,
12     &                      updaterice_microco2, updaterice_typ
13       use tracer_mod, only: nqmx, igcm_dust_mass, igcm_dust_number,
14     &                       igcm_h2o_ice, igcm_ccn_mass, radius,
15     &                       igcm_co2_ice, nuiceco2_ref,
16     &                       igcm_ccnco2_number, igcm_ccnco2_mass,
17     &                       igcm_ccnco2_h2o_number,
18     &                       igcm_ccnco2_h2o_mass_ice,
19     &                       igcm_ccnco2_h2o_mass_ccn,
20     &                       igcm_ccn_number, nuice_ref, varian,
21     &                       ref_r0, igcm_dust_submicron,
22     &                       igcm_stormdust_mass,igcm_stormdust_number,
23     &                       igcm_topdust_mass,igcm_topdust_number,
24     &                       rho_ice
25       USE dimradmars_mod, only: nueffdust,naerkind,
26     &            name_iaer,
27     &            iaer_dust_conrath,iaer_dust_doubleq,
28     &            iaer_dust_submicron,iaer_h2o_ice,
29     &            iaer_stormdust_doubleq,iaer_topdust_doubleq
30       use dust_param_mod, only: doubleq, active
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
53      include "callkeys.h"
54
55c-----------------------------------------------------------------------
56c     Inputs/outputs:
57c     ------
58
59      INTEGER, INTENT(in) :: ngrid,nlayer
60c     Ice geometric mean radius (m)
61      REAL, INTENT(out) :: rice(ngrid,nlayer)
62c     Estimated effective variance of the size distribution (n.u.)
63      REAL, INTENT(out) :: nuice(ngrid,nlayer)
64c     Tracer mass mixing ratio (kg/kg)
65      REAL, INTENT(in) :: pq(ngrid,nlayer,nqmx)
66      REAL, INTENT(out) :: rdust(ngrid,nlayer) ! Dust geometric mean radius (m)
67      REAL, INTENT(out) :: rstormdust(ngrid,nlayer) ! Dust geometric mean radius (m)   
68      REAL, INTENT(out) :: rtopdust(ngrid,nlayer) ! Dust geometric mean radius (m)
69      REAL, INTENT(in) :: pplay(ngrid,nlayer) ! altitude at the middle of the layers
70      REAL, INTENT(in) :: tau(ngrid,naerkind)
71c     Aerosol effective radius used for radiative transfer (meter)
72      REAL, INTENT(out) :: reffrad(ngrid,nlayer,naerkind)
73c     Aerosol effective variance used for radiative transfer (n.u.)
74      REAL, INTENT(out) :: nueffrad(ngrid,nlayer,naerkind)
75      REAL, INTENT(in) :: tauscaling(ngrid)         ! Convertion factor for qccn and Nccn
76c     CO2 ice mean radius (m)
77      double precision, INTENT(out) :: riceco2(ngrid,nlayer) ! co2 ice radius
78      REAL, INTENT(out) :: nuiceco2(ngrid,nlayer)
79      REAL, INTENT(in) :: pt(ngrid,nlayer) ! temperature
80     
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
92      REAL, PARAMETER :: ccn0 = 1.3E8, threshold = 1e-30 ! limit value
93     
94c     For microphysics only:     
95      REAL Mo,No                       ! Mass and number of ccn
96      REAL rhocloud(ngrid,nlayer)  ! Cloud density (kg.m-3)
97c     For CO2 microphysics only:
98      REAL :: rhocloudco2(ngrid, nlayer) ! co2 cloud density
99
100      LOGICAL,SAVE :: firstcall=.true.
101      REAL Nccnco2, Qccnco2, Niceco2
102      REAL Nccnco2_h2o, Qccnco2_h2o
103      REAL CBRT
104      EXTERNAL CBRT
105
106
107!$OMP THREADPRIVATE(firstcall)
108
109c==================================================================
110c 1. Update radius from fields from dynamics or initial state
111c==================================================================
112
113c       1.1 Dust particles
114c       ------------------
115        IF (doubleq.AND.active) THEN
116          DO l=1,nlayer
117            DO ig=1, ngrid
118              call updaterdust(pq(ig,l,igcm_dust_mass),
119     &                         pq(ig,l,igcm_dust_number),rdust(ig,l))
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
129          ENDDO
130        ENDIF
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
142
143        ! updating radius of topdust particles
144        IF (topflows.AND.active) THEN
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
153       
154c       1.2 Water-ice particles
155c       -----------------------
156
157        IF (water.AND.activice) THEN
158         IF (microphys) THEN
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
164          IF (firstcall) THEN
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
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
178            ENDDO
179          ENDDO
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
188            ENDDO
189          ENDDO
190 
191        ENDIF ! of if microphys
192       ENDIF ! of if (water.AND.activice)
193
194
195c       1.3 CO2-ice particles
196
197        IF (co2clouds.AND.activeco2ice) THEN
198          DO l=1,nlayer
199            DO ig=1,ngrid
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)
205              Nccnco2_h2o = 0.
206              Qccnco2_h2o = 0.
207              if (co2useh2o) then
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
214                Nccnco2 = Nccnco2 - Nccnco2_h2o
215                Qccnco2 = Qccnco2 - Qccnco2_h2o
216                if (Nccnco2 <= 0) then
217                  Nccnco2 = threshold
218                  Qccnco2 = threshold
219                end if
220              end if
221              call updaterice_microco2(dble(Niceco2),
222     &                              dble(Qccnco2), dble(Nccnco2),
223     &                              dble(Qccnco2_h2o),
224     &                              dble(Nccnco2_h2o),
225     &                              pt(ig,l),
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
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
244              reffrad(ig,l,iaer) = rdust(ig,l) *
245     &          (1.e0 + nueffdust(ig,l))**2.5
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
254              reffrad(ig,l,iaer) = rdust(ig,l) * ref_r0
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
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:
279              reffrad(ig,l,iaer)=rice(ig,l)*(1.+nuice_ref)
280              nueffrad(ig,l,iaer)=nuice_ref
281            ENDDO
282          ENDDO
283c==================================================================
284        CASE("co2_ice") aerkind             ! CO2 ice crystals
285c==================================================================
286          DO l=1,nlayer
287            DO ig=1,ngrid
288              reffrad(ig,l,iaer)=real(riceco2(ig,l))*(1.+nuiceco2_ref)
289              nueffrad(ig,l,iaer)=nuiceco2_ref
290            ENDDO
291          ENDDO
292c==================================================================
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==================================================================
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==================================================================
313        END SELECT aerkind
314      ENDDO ! iaer (loop on aerosol kind)
315
316      END SUBROUTINE updatereffrad
317     
318      END MODULE updatereffrad_mod
Note: See TracBrowser for help on using the repository browser.