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

Last change on this file since 2616 was 2584, checked in by romain.vande, 3 years ago

Second stage of implementation of Open_MP in the physic.
Run with callrad=.true.
File size: 12.3 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
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
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 (slpwind.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              Nccnco2 = pq(ig,l,igcm_ccnco2_number)
201              Qccnco2 = pq(ig,l,igcm_ccnco2_mass)
202              Nccnco2_h2o = 0.
203              Qccnco2_h2o = 0.
204              if (co2useh2o) then
205                Nccnco2_h2o = pq(ig,l,igcm_ccnco2_h2o_number)
206                Qccnco2_h2o = pq(ig,l,igcm_ccnco2_h2o_mass_ice)
207     &                        + pq(ig,l,igcm_ccnco2_h2o_mass_ccn)
208                Nccnco2 = Nccnco2 - Nccnco2_h2o
209                Qccnco2 = Qccnco2 - Qccnco2_h2o
210              end if
211
212              call updaterice_microco2(dble(pq(ig,l,igcm_co2_ice)),
213     &                              dble(Qccnco2), dble(Nccnco2),
214     &                              dble(Qccnco2_h2o),
215     &                              dble(Nccnco2_h2o),
216     &                              pt(ig,l),
217     &                              tauscaling(ig),riceco2(ig,l),
218     &                              rhocloudco2(ig,l))
219              nuiceco2(ig,l) = nuiceco2_ref
220            END DO
221          ENDDO
222        ENDIF ! of if (co2clouds.AND.activeco2ice)
223
224c==================================================================
225c 2. Radius used in the radiative transfer code (reffrad)
226c==================================================================
227
228      DO iaer = 1, naerkind ! Loop on aerosol kind
229        aerkind: SELECT CASE (name_iaer(iaer))
230c==================================================================
231        CASE("dust_conrath") aerkind         ! Typical dust profile
232c==================================================================
233          DO l=1,nlayer
234            DO ig=1,ngrid
235              reffrad(ig,l,iaer) = rdust(ig,l) *
236     &          (1.e0 + nueffdust(ig,l))**2.5
237              nueffrad(ig,l,iaer) = nueffdust(ig,l)
238            ENDDO
239          ENDDO
240c==================================================================
241        CASE("dust_doubleq") aerkind! Two-moment scheme for dust
242c==================================================================
243          DO l=1,nlayer
244            DO ig=1,ngrid
245              reffrad(ig,l,iaer) = rdust(ig,l) * ref_r0
246              nueffrad(ig,l,iaer) = nueffdust(ig,l)
247            ENDDO
248          ENDDO
249c==================================================================
250        CASE("dust_submicron") aerkind   ! Small dust population
251c==================================================================
252          DO l=1,nlayer
253            DO ig=1,ngrid
254              reffrad(ig,l,iaer)=radius(igcm_dust_submicron)
255              nueffrad(ig,l,iaer)=0.03
256            ENDDO
257          ENDDO     
258c==================================================================
259        CASE("h2o_ice") aerkind             ! Water ice crystals
260c==================================================================
261          DO l=1,nlayer
262            DO ig=1,ngrid
263c             About reffice, do not confuse the mass mean radius
264c             (rayon moyen massique) and the number median radius
265c             (or geometric mean radius, rayon moyen géométrique).
266c             rice is a mass mean radius, whereas rdust
267c             is a geometric mean radius:
268c             number median rad = mass mean rad x exp(-1.5 sigma0^2)
269c             (Montmessin et al. 2004 paragraph 30). Therefore:
270              reffrad(ig,l,iaer)=rice(ig,l)*(1.+nuice_ref)
271              nueffrad(ig,l,iaer)=nuice_ref
272            ENDDO
273          ENDDO
274c==================================================================
275        CASE("co2_ice") aerkind             ! CO2 ice crystals
276c==================================================================
277          DO l=1,nlayer
278            DO ig=1,ngrid
279              reffrad(ig,l,iaer)=real(riceco2(ig,l))*(1.+nuiceco2_ref)
280              nueffrad(ig,l,iaer)=nuiceco2_ref
281            ENDDO
282          ENDDO
283c==================================================================
284        CASE("stormdust_doubleq") aerkind! Two-moment scheme for
285c       stormdust; same distribution than normal dust
286c==================================================================
287          DO l=1,nlayer
288            DO ig=1,ngrid
289              reffrad(ig,l,iaer) = rstormdust(ig,l) * ref_r0
290              nueffrad(ig,l,iaer) = nueffdust(ig,l)
291            ENDDO
292          ENDDO
293c==================================================================
294        CASE("topdust_doubleq") aerkind! MV18: Two-moment scheme for
295c       topdust; same distribution than normal dust
296c==================================================================
297          DO l=1,nlayer
298            DO ig=1,ngrid
299              reffrad(ig,l,iaer) = rtopdust(ig,l) * ref_r0
300              nueffrad(ig,l,iaer) = nueffdust(ig,l)
301            ENDDO
302          ENDDO
303c==================================================================
304        END SELECT aerkind
305      ENDDO ! iaer (loop on aerosol kind)
306
307      END SUBROUTINE updatereffrad
308     
309      END MODULE updatereffrad_mod
Note: See TracBrowser for help on using the repository browser.