source: trunk/LMDZ.PLUTO/libf/phypluto/mp2m_calmufi.F90 @ 3966

Last change on this file since 3966 was 3962, checked in by debatzbr, 3 months ago

Pluto PCM: Calculates the latent heat released by condensation from the microphysical module.
BBT

File size: 12.2 KB
Line 
1MODULE mp2m_calmufi
2    use tracer_h
3    use comcstfi_mod, only : mugaz, cpp
4    use callkeys_mod, only : call_haze_prod_pCH4, haze_rho,&
5                             callmuclouds
6   
7    ! Microphysical model MP2M
8    use mp2m_intgcm
9    use mp2m_diagnostics
10   
11
12    implicit none
13   
14    !============================================================================
15    !
16    !     Purpose
17    !     -------
18    !     Interface subroutine to YAMMS model for LMD's PCM.
19    !
20    !     The subroutine computes the microphysics processes for a single vertical column.
21    !       - All input vectors are assumed to be defined from GROUND to TOP of the atmosphere.
22    !       - All output vectors are defined from GROUND to TOP of the atmosphere.
23    !       - Only tendencies are returned.
24    !
25    !     @important
26    !     The method assumes global initialization of YAMMS model (and extras) has been already
27    !     done elsewhere.
28    !
29    !     @warning
30    !     Microphysical tracers from physics must be in X/kg_of_air and convert into X/m2 for microphysics.
31    !
32    !     @warning
33    !     We suppose a given order of microphysical tracers in micro_indx(:):
34    !     1. mu_m0as, 2. mu_m3as, 3. mu_m0af, 4. mu_m3af.
35    !     If clouds are activated:
36    !     5. mu_m0ccn, 6. mu_m3ccn, 7(+). mu_m3ices.
37    !
38    !     Authors
39    !     -------
40    !     B. de Batz de Trenquelléon (11/2024)
41    !
42    !============================================================================
43
44    CONTAINS
45   
46    SUBROUTINE calmufi(dt, plev, zlev, play, zlay, g3d, temp, pq, zdqfi, zdqmufi_prod, zdqmufi, zdtcond)
47        !! Interface subroutine to YAMMS model for LMD PCM.
48        !!
49        !! The subroutine computes the microphysics processes for a single vertical column.
50        !!
51        REAL(kind=8), INTENT(IN) :: dt  !! Physics timestep (s).
52       
53        REAL(kind=8), DIMENSION(:,:), INTENT(IN) :: plev            ! Pressure levels (Pa).
54        REAL(kind=8), DIMENSION(:,:), INTENT(IN) :: zlev            ! Altitude levels (m).
55        REAL(kind=8), DIMENSION(:,:), INTENT(IN) :: play            ! Pressure layers (Pa).
56        REAL(kind=8), DIMENSION(:,:), INTENT(IN) :: zlay            ! Altitude at the center of each layer (m).
57        REAL(kind=8), DIMENSION(:,:), INTENT(IN) :: g3d             ! Latitude-Altitude depending gravitational acceleration (m.s-2).
58        REAL(kind=8), DIMENSION(:,:), INTENT(IN) :: temp            ! Temperature at the center of each layer (K).
59       
60        REAL(kind=8), DIMENSION(:,:,:), INTENT(IN)  :: pq           ! Tracers (X.kg-1).
61        REAL(kind=8), DIMENSION(:,:,:), INTENT(IN)  :: zdqfi        ! Tendency from former processes for tracers (X.kg-1.s-1).
62        REAL(kind=8), DIMENSION(:,:,:), INTENT(IN)  :: zdqmufi_prod ! Aerosols production tendency (kg/kg_of_air/s).
63        REAL(kind=8), DIMENSION(:,:,:), INTENT(OUT) :: zdqmufi      ! Microphysical tendency for tracers (X.m-2 --> X.kg-1.s-1).
64        REAL(kind=8), DIMENSION(:,:),   INTENT(OUT) :: zdtcond      ! Condensation heating rate (K.s-1).
65       
66        ! Local tracers:
67        !~~~~~~~~~~~~~~~
68        REAL(kind=8), DIMENSION(:,:,:), ALLOCATABLE :: zq    ! Local tracers updated from former processes (X.kg-1).
69       
70        ! Related to aerosols:
71        REAL(kind=8), DIMENSION(:), ALLOCATABLE :: m0as      ! 0th order moment of the spherical mode (m-2).
72        REAL(kind=8), DIMENSION(:), ALLOCATABLE :: m3as      ! 3rd order moment of the spherical mode (m3.m-2).
73        REAL(kind=8), DIMENSION(:), ALLOCATABLE :: m0af      ! 0th order moment of the fractal mode (m-2).
74        REAL(kind=8), DIMENSION(:), ALLOCATABLE :: m3af      ! 3rd order moment of the fractal mode (m3.m-2).
75
76        REAL(kind=8), DIMENSION(:), ALLOCATABLE :: m3as_prod ! Production of 3rd order moment of the spherical mode (m3.m-2).
77       
78        REAL(kind=8), DIMENSION(:), ALLOCATABLE   :: dm0as   ! Tendency of the 0th order moment of the spherical mode distribution (m-2).
79        REAL(kind=8), DIMENSION(:), ALLOCATABLE   :: dm3as   ! Tendency of the 3rd order moment of the spherical mode distribution (m3.m-2).
80        REAL(kind=8), DIMENSION(:), ALLOCATABLE   :: dm0af   ! Tendency of the 0th order moment of the fractal mode distribution (m-2).
81        REAL(kind=8), DIMENSION(:), ALLOCATABLE   :: dm3af   ! Tendency of the 3rd order moment of the fractal mode distribution (m3.m-2).
82
83        ! Related to clouds:
84        REAL(kind=8), DIMENSION(:), ALLOCATABLE   :: m0ccn   ! 0th order moment of the ccn distribution (m-2).
85        REAL(kind=8), DIMENSION(:), ALLOCATABLE   :: m3ccn   ! 3rd order moment of the ccn distribution (m3.m-2).
86        REAL(kind=8), DIMENSION(:,:), ALLOCATABLE :: m3ices  ! 3rd order moments of the ice components (m3.m-2).
87        REAL(kind=8), DIMENSION(:,:), ALLOCATABLE :: mugases ! Condensible species gas molar fraction (mol.mol-1).
88
89        REAL(kind=8), DIMENSION(:), ALLOCATABLE   :: dm0ccn   ! Tendency of the 0th order moment of the ccn distribution (m-2).
90        REAL(kind=8), DIMENSION(:), ALLOCATABLE   :: dm3ccn   ! Tendency of the 3rd order moment of the ccn distribution (m3.m-2).
91        REAL(kind=8), DIMENSION(:,:), ALLOCATABLE :: dm3ices  ! Tendencies of the 3rd order moments of each ice components (m3.m-2).
92        REAL(kind=8), DIMENSION(:,:), ALLOCATABLE :: dmugases ! Tendencies of each condensible gas species (mol.mol-1).
93        REAL(kind=8), DIMENSION(:), ALLOCATABLE   :: dtlc     ! Latent heat of condensation (J.kg-1).
94       
95        ! Local variables:
96        !~~~~~~~~~~~~~~~~~
97        REAL(kind=8), DIMENSION(:,:), ALLOCATABLE :: int2ext ! Conversion intensive to extensive (kg.m-2).
98        REAL(kind=8), DIMENSION(:), ALLOCATABLE   :: tmp
99        TYPE(error) :: err
100       
101        INTEGER           :: i,ilon,iq
102        INTEGER           :: nq,nlon,nlay
103        CHARACTER(len=10) :: tname
104       
105        ! Read size of arrays:
106        !~~~~~~~~~~~~~~~~~~~~~
107        nq    = size(pq,DIM=3)
108        nlon  = size(play,DIM=1)
109        nlay  = size(play,DIM=2)
110       
111        ! Allocate arrays:
112        !~~~~~~~~~~~~~~~~~
113        ALLOCATE(zq(nlon,nlay,nq))
114       
115        ALLOCATE(m0as(nlay))
116        ALLOCATE(m3as(nlay))
117        ALLOCATE(m0af(nlay))
118        ALLOCATE(m3af(nlay))
119        ALLOCATE(m3as_prod(nlay))
120        ALLOCATE(dm0as(nlay))
121        ALLOCATE(dm3as(nlay))
122        ALLOCATE(dm0af(nlay))
123        ALLOCATE(dm3af(nlay))
124
125        ALLOCATE(m0ccn(nlay))
126        ALLOCATE(m3ccn(nlay))
127        ALLOCATE(m3ices(nlay,nmicro_ices))
128        ALLOCATE(mugases(nlay,nmicro_ices))
129        ALLOCATE(dm0ccn(nlay))
130        ALLOCATE(dm3ccn(nlay))
131        ALLOCATE(dm3ices(nlay,nmicro_ices))
132        ALLOCATE(dmugases(nlay,nmicro_ices))
133        ALLOCATE(dtlc(nlay))
134
135        ALLOCATE(int2ext(nlon,nlay)) 
136       
137        !------------------
138        ! 1. Initialization
139        !------------------
140
141        ! Initialization of zdqmufi here since intent=out and no action performed on every tracers
142        zdqmufi(:,:,:) = 0.D0
143
144        ! Initialization of zdtcond here since intent=out
145        zdtcond(:,:) = 0.D0
146       
147        ! Initialize tracers updated with former processes from physics
148        zq(:,:,:) = pq(:,:,:) + zdqfi(:,:,:)*dt
149       
150        ! Loop on horizontal grid points
151        DO ilon = 1, nlon
152
153            ! Convert tracers to extensive [X.kg-1 --> X.m-2]
154            int2ext(ilon,:) = (plev(ilon,1:nlay)-plev(ilon,2:nlay+1)) / g3d(ilon,1:nlay)
155           
156            m0as(:) = zq(ilon,:,micro_indx(1)) * int2ext(ilon,:)
157            m3as(:) = zq(ilon,:,micro_indx(2)) * int2ext(ilon,:)
158           
159            m0af(:) = zq(ilon,:,micro_indx(3)) * int2ext(ilon,:)
160            m3af(:) = zq(ilon,:,micro_indx(4)) * int2ext(ilon,:)
161           
162            if (callmuclouds) then
163                m0ccn(:) = zq(ilon,:,micro_indx(5)) * int2ext(ilon,:)
164                m3ccn(:) = zq(ilon,:,micro_indx(6)) * int2ext(ilon,:)
165                do i = 1, nmicro_ices
166                    m3ices(:,i)  = zq(ilon,:,micro_ice_indx(i)) * int2ext(ilon,:)
167                    mugases(:,i) = zq(ilon,:,micro_gas_indx(i)) / (mmol(micro_gas_indx(i))/mugaz)
168                enddo
169            endif ! end of callmuclouds
170
171            ! Production of haze in the atmosphere by photolysis of CH4
172            if (call_haze_prod_pCH4) then
173                do iq = 1, nq
174                    tname = noms(iq)
175                    if (tname(1:4).eq."haze") then
176                        m3as_prod(:) = zdqmufi_prod(ilon,:,iq) * (int2ext(ilon,:) / haze_rho) * dt
177                    endif
178                enddo
179            else
180                m3as_prod(:) = 0.D0
181            endif
182           
183            ! Hackin the pressure level
184            tmp = plev(ilon,:)
185            if (tmp(nlay+1) == 0.0) then
186            tmp(nlay+1) = 2*tmp(nlay) - tmp(nlay-1)
187            endif
188       
189            ! Initialize YAMMS atmospheric column
190            err = mm_column_init(tmp,zlev(ilon,:),play(ilon,:),zlay(ilon,:),temp(ilon,:)) ; IF (err /= 0) call abort_program(err)
191            ! Initialize YAMMS aerosol moments column
192            err = mm_aerosols_init(m0as,m3as,m0af,m3af) ; IF (err /= 0) call abort_program(err)
193            ! Initialize YAMMS cloud moments column
194            if (callmuclouds) then
195                err = mm_clouds_init(m0ccn,m3ccn,m3ices,mugases) ; IF (err /= 0) call abort_program(err)
196            endif ! end of callmuclouds
197           
198            ! Initializes tendencies
199            dm0as(:)  = 0._mm_wp ; dm3as(:)  = 0._mm_wp ; dm0af(:)     = 0._mm_wp ; dm3af(:)      = 0._mm_wp
200            dm0ccn(:) = 0._mm_wp ; dm3ccn(:) = 0._mm_wp ; dm3ices(:,:) = 0._mm_wp ; dmugases(:,:) = 0._mm_wp
201            dtlc(:)   = 0._mm_wp
202           
203        !----------------------------
204        ! 2. Call microphysical model
205        !----------------------------
206           
207            ! Call microphysics
208            if (callmuclouds) then
209                if(.NOT.mm_muphys(m3as_prod,dm0as,dm3as,dm0af,dm3af,dm0ccn,dm3ccn,dm3ices,dmugases,dtlc)) then
210                    call abort_program(error("mm_muphys (clouds) aborted -> initialization not done !",-1))
211                endif
212            else
213                if (.NOT.mm_muphys(m3as_prod,dm0as,dm3as,dm0af,dm3af)) then
214                    call abort_program(error("mm_muphys (no cloud) aborted -> initialization not done !",-1))
215                endif
216            endif ! end of callmuclouds
217           
218            ! Save diagnostics
219            call mm_diagnostics(dt,mp2m_aer_s_prec(ilon),mp2m_aer_f_prec(ilon),  &
220                                mp2m_aer_s_w(ilon,:),mp2m_aer_f_w(ilon,:),       &
221                                mp2m_aer_s_flux(ilon,:),mp2m_aer_f_flux(ilon,:), &
222                                mp2m_rc_sph(ilon,:),mp2m_rc_fra(ilon,:),         &
223                                mp2m_ccn_prec(ilon),mp2m_ice_prec(ilon,:),       &
224                                mp2m_cld_w(ilon,:),mp2m_ccn_flux(ilon,:),        &
225                                mp2m_ice_fluxes(ilon,:,:),mp2m_rc_cld(ilon,:),   &
226                                mp2m_gas_sat(ilon,:,:),mp2m_nrate(ilon,:,:),mp2m_grate(ilon,:,:))
227       
228            ! Convert tracers back to intensives
229            zdqmufi(ilon,:,micro_indx(1)) = dm0as(:) / int2ext(ilon,:)
230            zdqmufi(ilon,:,micro_indx(2)) = dm3as(:) / int2ext(ilon,:)
231            zdqmufi(ilon,:,micro_indx(3)) = dm0af(:) / int2ext(ilon,:)
232            zdqmufi(ilon,:,micro_indx(4)) = dm3af(:) / int2ext(ilon,:)
233
234            if (callmuclouds) then
235                zdqmufi(ilon,:,micro_indx(5)) = dm0ccn(:) / int2ext(ilon,:)
236                zdqmufi(ilon,:,micro_indx(6)) = dm3ccn(:) / int2ext(ilon,:)
237                do i = 1, nmicro_ices
238                    zdqmufi(ilon,:,micro_ice_indx(i)) = dm3ices(:,i) / int2ext(ilon,:)
239                    zdqmufi(ilon,:,micro_gas_indx(i)) = dmugases(:,i) * (mmol(micro_gas_indx(i))/mugaz)
240                enddo
241
242                ! Compute condensation heating rate in K.s-1
243                zdtcond(ilon,:) = dtlc(:) / cpp / dt
244            endif ! End of callmuclouds
245           
246        END DO ! End loop on ilon
247
248        ! YAMMS gives a tendency which is integrated for all the timestep but in the GCM
249        ! we want to have routines spitting tendencies in s-1 -> let's divide !
250        zdqmufi(:,:,:) = zdqmufi(:,:,:) / dt
251       
252    END SUBROUTINE calmufi
253
254END MODULE mp2m_calmufi
Note: See TracBrowser for help on using the repository browser.