| 1 | ! Copyright (2013-2015,2017,2022-2023) Université de Reims Champagne-Ardenne |
|---|
| 2 | ! Contributors : J. Burgalat (GSMA, URCA), B. de Batz de Trenquelléon (GSMA, URCA) |
|---|
| 3 | ! email of the author : jeremie.burgalat@univ-reims.fr |
|---|
| 4 | ! |
|---|
| 5 | ! This software is a computer program whose purpose is to compute |
|---|
| 6 | ! microphysics processes using a two-moments scheme. |
|---|
| 7 | ! |
|---|
| 8 | ! This library is governed by the CeCILL-B license under French law and |
|---|
| 9 | ! abiding by the rules of distribution of free software. You can use, |
|---|
| 10 | ! modify and/ or redistribute the software under the terms of the CeCILL-B |
|---|
| 11 | ! license as circulated by CEA, CNRS and INRIA at the following URL |
|---|
| 12 | ! "http://www.cecill.info". |
|---|
| 13 | ! |
|---|
| 14 | ! As a counterpart to the access to the source code and rights to copy, |
|---|
| 15 | ! modify and redistribute granted by the license, users are provided only |
|---|
| 16 | ! with a limited warranty and the software's author, the holder of the |
|---|
| 17 | ! economic rights, and the successive licensors have only limited |
|---|
| 18 | ! liability. |
|---|
| 19 | ! |
|---|
| 20 | ! In this respect, the user's attention is drawn to the risks associated |
|---|
| 21 | ! with loading, using, modifying and/or developing or reproducing the |
|---|
| 22 | ! software by the user in light of its specific status of free software, |
|---|
| 23 | ! that may mean that it is complicated to manipulate, and that also |
|---|
| 24 | ! therefore means that it is reserved for developers and experienced |
|---|
| 25 | ! professionals having in-depth computer knowledge. Users are therefore |
|---|
| 26 | ! encouraged to load and test the software's suitability as regards their |
|---|
| 27 | ! requirements in conditions enabling the security of their systems and/or |
|---|
| 28 | ! data to be ensured and, more generally, to use and operate it in the |
|---|
| 29 | ! same conditions as regards security. |
|---|
| 30 | ! |
|---|
| 31 | ! The fact that you are presently reading this means that you have had |
|---|
| 32 | ! knowledge of the CeCILL-B license and that you accept its terms. |
|---|
| 33 | |
|---|
| 34 | !! file: mm_microphysic.f90 |
|---|
| 35 | !! brief: Microphysic processes interface module. |
|---|
| 36 | !! author: J. Burgalat |
|---|
| 37 | !! date: 2013-2015,2017,2022-2023 |
|---|
| 38 | !! modifications: B. de Batz de Trenquelléon |
|---|
| 39 | |
|---|
| 40 | MODULE MM_MICROPHYSIC |
|---|
| 41 | !! Microphysic processes interface module. |
|---|
| 42 | USE MM_MPREC |
|---|
| 43 | USE MM_GLOBALS |
|---|
| 44 | USE MM_HAZE |
|---|
| 45 | USE MM_CLOUDS |
|---|
| 46 | USE MM_METHODS |
|---|
| 47 | IMPLICIT NONE |
|---|
| 48 | |
|---|
| 49 | PRIVATE |
|---|
| 50 | |
|---|
| 51 | PUBLIC :: mm_muphys, mm_diagnostics, mm_get_radii |
|---|
| 52 | |
|---|
| 53 | !! Interface to main microphysics subroutine. |
|---|
| 54 | !! |
|---|
| 55 | !! The interface computes calls either all the microphysics processes ([[mm_microphysic(module):muphys_all(function)]] |
|---|
| 56 | !! or only aerosols microphysics ([[mm_microphysic(module):muphys_nocld(function)]]) in a single call. |
|---|
| 57 | INTERFACE mm_muphys |
|---|
| 58 | MODULE PROCEDURE muphys_all, muphys_nocld |
|---|
| 59 | END INTERFACE mm_muphys |
|---|
| 60 | |
|---|
| 61 | CONTAINS |
|---|
| 62 | |
|---|
| 63 | |
|---|
| 64 | |
|---|
| 65 | FUNCTION muphys_all(dm0a_s,dm3a_s,dm0a_f,dm3a_f,dm0n,dm3n,dm3i,dgazs) RESULT(ret) |
|---|
| 66 | !! Compute the evolution of moments tracers through haze and clouds microphysics processes. |
|---|
| 67 | !! |
|---|
| 68 | !! This method computes the evolution of all the microphysics tracers, given under the form of moments |
|---|
| 69 | !! (and molar fraction for cloud condensible species) during a time step. |
|---|
| 70 | !! |
|---|
| 71 | !! The method requires that global variables of the model (i.e. variables declared in [[mm_globals(module)]] |
|---|
| 72 | !! module) are initialized/updated correctly (see [[mm_globals(module):mm_global_init(interface)]], |
|---|
| 73 | !! [[mm_globals(module):mm_column_init(function)]],[[mm_globals(module):mm_aerosols_init(function)]] and |
|---|
| 74 | !! [[mm_globals(module):mm_clouds_init(function)]]). |
|---|
| 75 | !! |
|---|
| 76 | !! The tendencies returned by the method are defined on the vertical __layers__ of the model from the __GROUND__ to |
|---|
| 77 | !! the __TOP__ of the atmosphere. They should be added to the input variables used in the initialization methods |
|---|
| 78 | !! before the latter are called to initialize a new step. |
|---|
| 79 | !! @note |
|---|
| 80 | !! __dm3i__ and __dgazs__ are 2D-arrays with vertical __layers__ in the 1st dimension and the number of |
|---|
| 81 | !! ice components in the 2nd. They share the same _species_ indexing. |
|---|
| 82 | !! |
|---|
| 83 | !! It should be a 2D-array with the vertical layers in first dimension and the number of ice components in the second. |
|---|
| 84 | REAL(kind=mm_wp), INTENT(out), DIMENSION(:) :: dm0a_s |
|---|
| 85 | !! Tendency of the 0th order moment of the spherical mode distribution (\(m^{-2}\)). |
|---|
| 86 | REAL(kind=mm_wp), INTENT(out), DIMENSION(:) :: dm3a_s |
|---|
| 87 | !! Tendency of the 3rd order moment of the spherical mode distribution (\(m^{3}.m^{-2}\)). |
|---|
| 88 | REAL(kind=mm_wp), INTENT(out), DIMENSION(:) :: dm0a_f |
|---|
| 89 | !! Tendency of the 0th order moment of the fractal mode distribution (\(m^{-2}\)). |
|---|
| 90 | REAL(kind=mm_wp), INTENT(out), DIMENSION(:) :: dm3a_f |
|---|
| 91 | !! Tendency of the 3rd order moment of the fractal mode distribution (\(m^{3}.m^{-2}\)). |
|---|
| 92 | REAL(kind=mm_wp), INTENT(out), DIMENSION(:) :: dm0n |
|---|
| 93 | !! Tendency of the 0th order moment of the _CCN_ distribution (\(m^{-2}\)). |
|---|
| 94 | REAL(kind=mm_wp), INTENT(out), DIMENSION(:) :: dm3n |
|---|
| 95 | !! Tendency of the 3rd order moment of the _CCN_ distribution (\(m^{3}.m^{-2}\)). |
|---|
| 96 | REAL(kind=mm_wp), INTENT(out), DIMENSION(:,:) :: dm3i |
|---|
| 97 | !! Tendencies of the 3rd order moments of each ice components (\(m^{3}.m^{-2}\)). |
|---|
| 98 | REAL(kind=mm_wp), INTENT(out), DIMENSION(:,:) :: dgazs |
|---|
| 99 | !! Tendencies of each condensible gaz species (\(mol.mol^{-1}\)). |
|---|
| 100 | LOGICAL :: ret |
|---|
| 101 | !! .true. on success (i.e. model has been initialized at least once previously), .false. otherwise. |
|---|
| 102 | REAL(kind=mm_wp), DIMENSION(SIZE(dm0a_s)) :: zdm0a_f,zdm3a_f |
|---|
| 103 | INTEGER :: i |
|---|
| 104 | ! Checks initialization |
|---|
| 105 | ret = (mm_ini_col.AND.mm_ini_aer.AND.(.NOT.mm_w_clouds.OR.mm_ini_cld)) |
|---|
| 106 | IF (.NOT.ret) RETURN |
|---|
| 107 | ! Calls haze microphysics (-> m-3) |
|---|
| 108 | call mm_haze_microphysics(dm0a_s,dm3a_s,dm0a_f,dm3a_f) |
|---|
| 109 | IF (mm_w_clouds) THEN |
|---|
| 110 | ! Calls cloud microphysics (-> m-3) |
|---|
| 111 | call mm_cloud_microphysics(zdm0a_f,zdm3a_f,dm0n,dm3n,dm3i,dgazs) |
|---|
| 112 | ! add temporary aerosols tendencies (-> m-3) |
|---|
| 113 | dm0a_f = dm0a_f + zdm0a_f ; dm3a_f = dm3a_f + zdm3a_f |
|---|
| 114 | ! reverse directly clouds tendencies (-> m-2) |
|---|
| 115 | dm0n = dm0n(mm_nla:1:-1) * mm_dzlev(mm_nla:1:-1) |
|---|
| 116 | dm3n = dm3n(mm_nla:1:-1) * mm_dzlev(mm_nla:1:-1) |
|---|
| 117 | DO i=1,mm_nesp |
|---|
| 118 | dm3i(:,i) = dm3i(mm_nla:1:-1,i) * mm_dzlev(mm_nla:1:-1) |
|---|
| 119 | dgazs(:,i) = dgazs(mm_nla:1:-1,i) |
|---|
| 120 | ENDDO |
|---|
| 121 | ELSE |
|---|
| 122 | dm0n = 0._mm_wp ; dm3n = 0._mm_wp ; dm3i = 0._mm_wp ; dgazs = 0._mm_wp |
|---|
| 123 | ENDIF |
|---|
| 124 | ! multiply by altitude thickness and reverse vectors so they go from ground to top :) |
|---|
| 125 | dm0a_s = dm0a_s(mm_nla:1:-1) * mm_dzlev(mm_nla:1:-1) |
|---|
| 126 | dm3a_s = dm3a_s(mm_nla:1:-1) * mm_dzlev(mm_nla:1:-1) |
|---|
| 127 | dm0a_f = dm0a_f(mm_nla:1:-1) * mm_dzlev(mm_nla:1:-1) |
|---|
| 128 | dm3a_f = dm3a_f(mm_nla:1:-1) * mm_dzlev(mm_nla:1:-1) |
|---|
| 129 | RETURN |
|---|
| 130 | END FUNCTION muphys_all |
|---|
| 131 | |
|---|
| 132 | FUNCTION muphys_nocld(dm0a_s,dm3a_s,dm0a_f,dm3a_f) RESULT(ret) |
|---|
| 133 | !! Compute the evolution of moments tracers through haze microphysics only. |
|---|
| 134 | !! |
|---|
| 135 | !! This method is a __light__ version of [[mm_microphysic(module):muphys_all(function)]] where |
|---|
| 136 | !! only haze microphysics is computed and its tendencies returned. |
|---|
| 137 | !! |
|---|
| 138 | !! The method has the same requirements and remarks than [[mm_microphysic(module):muphys_all(function)]]. |
|---|
| 139 | REAL(kind=mm_wp), INTENT(out), DIMENSION(:) :: dm0a_s |
|---|
| 140 | !! Tendency of the 0th order moment of the spherical mode distribution (\(m^{-2}\)). |
|---|
| 141 | REAL(kind=mm_wp), INTENT(out), DIMENSION(:) :: dm3a_s |
|---|
| 142 | !! Tendency of the 3rd order moment of the spherical mode distribution (\(m^{3}.m^{-2}\)). |
|---|
| 143 | REAL(kind=mm_wp), INTENT(out), DIMENSION(:) :: dm0a_f |
|---|
| 144 | !! Tendency of the 0th order moment of the fractal mode distribution (\(m^{-2}\)). |
|---|
| 145 | REAL(kind=mm_wp), INTENT(out), DIMENSION(:) :: dm3a_f |
|---|
| 146 | !! Tendency of the 3rd order moment of the fractal mode distribution (\(m^{3}.m^{-2}\)). |
|---|
| 147 | LOGICAL :: ret |
|---|
| 148 | !! .true. on succes (i.e. model has been initialized at least once previously), .false. otherwise. |
|---|
| 149 | ret = (mm_ini_col.AND.mm_ini_aer) |
|---|
| 150 | IF (.NOT.ret) RETURN |
|---|
| 151 | IF (mm_w_clouds.AND.mm_debug) THEN |
|---|
| 152 | WRITE(*,'(a)') "WARNING: clouds microphysics enabled but will not be & |
|---|
| 153 | &computed... (wrong interface)" |
|---|
| 154 | ENDIF |
|---|
| 155 | ! Calls haze microphysics |
|---|
| 156 | call mm_haze_microphysics(dm0a_s,dm3a_s,dm0a_f,dm3a_f) |
|---|
| 157 | ! reverse vectors so they go from ground to top :) |
|---|
| 158 | dm0a_s = dm0a_s(mm_nla:1:-1) * mm_dzlev(mm_nla:1:-1) |
|---|
| 159 | dm3a_s = dm3a_s(mm_nla:1:-1) * mm_dzlev(mm_nla:1:-1) |
|---|
| 160 | dm0a_f = dm0a_f(mm_nla:1:-1) * mm_dzlev(mm_nla:1:-1) |
|---|
| 161 | dm3a_f = dm3a_f(mm_nla:1:-1) * mm_dzlev(mm_nla:1:-1) |
|---|
| 162 | RETURN |
|---|
| 163 | END FUNCTION muphys_nocld |
|---|
| 164 | |
|---|
| 165 | SUBROUTINE mm_diagnostics(dt,aer_prec,aer_s_w,aer_f_w,aer_s_flux,aer_f_flux,ccn_prec,ccn_w,ccn_flux,ice_prec,ice_fluxes,gazs_sat) |
|---|
| 166 | !! Get various diagnostic fields of the microphysics. |
|---|
| 167 | !! |
|---|
| 168 | !! The current diagnostics saved during microphysic computation are : |
|---|
| 169 | !! |
|---|
| 170 | !! - Mass fluxes (aerosols -both mode-, CCN and ices) |
|---|
| 171 | !! - Settling velocity (aerosols -total-, CCN and ices) |
|---|
| 172 | !! - Precipitations (aerosols -total-, CCN and ices) |
|---|
| 173 | !! - condensible gazes saturation ratio |
|---|
| 174 | !! |
|---|
| 175 | !! @note |
|---|
| 176 | !! Fluxes values are always negative as they account for sedimentation fluxes. They are set as |
|---|
| 177 | !! vector (for aerosols and CCN) or 2D-array (with the vertical structure in the first dimension |
|---|
| 178 | !! and number of species in the second, for ice) and are ordered from __GROUND__ to __TOP__. |
|---|
| 179 | !! |
|---|
| 180 | !! @note |
|---|
| 181 | !! Precipitations are always positive and defined in meters. For ice, it is set as a vector with |
|---|
| 182 | !! the precipitations of each cloud ice components. |
|---|
| 183 | !! |
|---|
| 184 | !! @note |
|---|
| 185 | !! __ccnprec__, __iceprec__, __icefluxes__ and __gazsat__ are always set to 0 if clouds |
|---|
| 186 | !! microphysics is disabled (see [[mm_globals(module):mm_w_clouds(variable)]] documentation). |
|---|
| 187 | REAL(kind=8), INTENT(IN) :: dt !! Physics timestep (s). |
|---|
| 188 | REAL(kind=mm_wp), INTENT(out), OPTIONAL :: aer_prec !! Aerosols precipitations (both modes) (kg.m-2.s-1). |
|---|
| 189 | REAL(kind=mm_wp), INTENT(out), OPTIONAL :: ccn_prec !! CCN precipitations (kg.m-2.s-1). |
|---|
| 190 | REAL(kind=mm_wp), INTENT(out), OPTIONAL, DIMENSION(:) :: aer_s_w !! Spherical aerosol settling velocity (\(m.s^{-1}\)). |
|---|
| 191 | REAL(kind=mm_wp), INTENT(out), OPTIONAL, DIMENSION(:) :: aer_f_w !! Fractal aerosol settling velocity (\(m.s^{-1}\)). |
|---|
| 192 | REAL(kind=mm_wp), INTENT(out), OPTIONAL, DIMENSION(:) :: ccn_w !! CCN settling velocity (\(m.s^{-1}\)). |
|---|
| 193 | REAL(kind=mm_wp), INTENT(out), OPTIONAL, DIMENSION(:) :: aer_s_flux !! Spherical aerosol mass flux (\(kg.m^{-2}.s^{-1}\)). |
|---|
| 194 | REAL(kind=mm_wp), INTENT(out), OPTIONAL, DIMENSION(:) :: aer_f_flux !! Fractal aerosol mass flux (\(kg.m^{-2}.s^{-1}\)). |
|---|
| 195 | REAL(kind=mm_wp), INTENT(out), OPTIONAL, DIMENSION(:) :: ccn_flux !! CCN mass flux (\(kg.m^{-2}.s^{-1}\)). |
|---|
| 196 | REAL(kind=mm_wp), INTENT(out), OPTIONAL, DIMENSION(:,:) :: ice_fluxes !! Ice sedimentation fluxes (\(kg.m^{-2}.s^{-1}\)). |
|---|
| 197 | REAL(kind=mm_wp), INTENT(out), OPTIONAL, DIMENSION(:,:) :: gazs_sat !! Condensible gaz saturation ratios (--). |
|---|
| 198 | REAL(kind=mm_wp), INTENT(out), OPTIONAL, DIMENSION(:) :: ice_prec !! Ice precipitations (kg.m-2.s-1). |
|---|
| 199 | |
|---|
| 200 | IF (PRESENT(aer_prec)) aer_prec = ABS(mm_aer_prec) / dt |
|---|
| 201 | IF (PRESENT(aer_s_w)) aer_s_w = -mm_m3as_vsed(mm_nla:1:-1) |
|---|
| 202 | IF (PRESENT(aer_f_w)) aer_f_w = -mm_m3af_vsed(mm_nla:1:-1) |
|---|
| 203 | IF (PRESENT(aer_s_flux)) aer_s_flux = -mm_aer_s_flux(mm_nla:1:-1) |
|---|
| 204 | IF (PRESENT(aer_f_flux)) aer_f_flux = -mm_aer_f_flux(mm_nla:1:-1) |
|---|
| 205 | |
|---|
| 206 | IF (mm_w_clouds) THEN |
|---|
| 207 | IF (PRESENT(ccn_prec)) ccn_prec = ABS(mm_ccn_prec) / dt |
|---|
| 208 | IF (PRESENT(ice_prec)) ice_prec = ABS(mm_ice_prec) / dt |
|---|
| 209 | IF (PRESENT(ccn_w)) ccn_w = mm_ccn_vsed(mm_nla:1:-1) |
|---|
| 210 | IF (PRESENT(ccn_flux)) ccn_flux = mm_ccn_flux(mm_nla:1:-1) |
|---|
| 211 | IF (PRESENT(ice_fluxes)) ice_fluxes = mm_ice_fluxes(mm_nla:1:-1,:) |
|---|
| 212 | IF (PRESENT(gazs_sat)) gazs_sat = mm_gazs_sat(mm_nla:1:-1,:) |
|---|
| 213 | ELSE |
|---|
| 214 | IF (PRESENT(ccn_prec)) ccn_prec = 0._mm_wp |
|---|
| 215 | IF (PRESENT(ice_prec)) ice_prec = 0._mm_wp |
|---|
| 216 | IF (PRESENT(ccn_w)) ccn_w = 0._mm_wp |
|---|
| 217 | IF (PRESENT(ccn_flux)) ccn_flux = 0._mm_wp |
|---|
| 218 | IF (PRESENT(ice_fluxes)) ice_fluxes = 0._mm_wp |
|---|
| 219 | IF (PRESENT(gazs_sat)) gazs_sat = 0._mm_wp |
|---|
| 220 | ENDIF |
|---|
| 221 | END SUBROUTINE mm_diagnostics |
|---|
| 222 | |
|---|
| 223 | SUBROUTINE mm_get_radii(rcsph,rcfra,rccld) |
|---|
| 224 | !! Get characteristic radii of microphysical tracers on the vertical grid. |
|---|
| 225 | REAL(kind=mm_wp), INTENT(out), DIMENSION(:), OPTIONAL :: rcsph !! Spherical mode characteristic radius |
|---|
| 226 | REAL(kind=mm_wp), INTENT(out), DIMENSION(:), OPTIONAL :: rcfra !! Fractal mode characteristic radius |
|---|
| 227 | REAL(kind=mm_wp), INTENT(out), DIMENSION(:), OPTIONAL :: rccld !! Cloud drops mean radius |
|---|
| 228 | IF (mm_ini_aer) THEN |
|---|
| 229 | IF (PRESENT(rcsph)) rcsph = mm_rcs(mm_nla:1:-1) |
|---|
| 230 | IF (PRESENT(rcfra)) rcfra = mm_rcf(mm_nla:1:-1) |
|---|
| 231 | ELSE |
|---|
| 232 | IF (PRESENT(rcsph)) rcsph = 0._mm_wp |
|---|
| 233 | IF (PRESENT(rcfra)) rcfra = 0._mm_wp |
|---|
| 234 | ENDIF |
|---|
| 235 | IF (PRESENT(rccld)) THEN |
|---|
| 236 | IF (mm_w_clouds.AND.mm_ini_cld) THEN |
|---|
| 237 | rccld = mm_drad(mm_nla:1:-1) |
|---|
| 238 | ELSE |
|---|
| 239 | rccld = 0._mm_wp |
|---|
| 240 | ENDIF |
|---|
| 241 | ENDIF |
|---|
| 242 | END SUBROUTINE mm_get_radii |
|---|
| 243 | |
|---|
| 244 | END MODULE MM_MICROPHYSIC |
|---|
| 245 | |
|---|
