[6610] | 1 | !$Id: dustsource.F90 104 2008-12-23 10:28:51Z acosce $ |
---|
| 2 | !! ========================================================================= |
---|
| 3 | !! INCA - INteraction with Chemistry and Aerosols |
---|
| 4 | !! |
---|
| 5 | !! Copyright Laboratoire des Sciences du Climat et de l'Environnement (LSCE) |
---|
| 6 | !! Unite mixte CEA-CNRS-UVSQ |
---|
| 7 | !! |
---|
| 8 | !! Contributors to this INCA subroutine: |
---|
| 9 | !! |
---|
| 10 | !! Michael Schulz, LSCE, Michael.Schulz@cea.fr |
---|
| 11 | !! |
---|
| 12 | !! Anne Cozic, LSCE, anne.cozic@cea.fr |
---|
| 13 | !! Yann Meurdesoif, LSCE, yann.meurdesoif@cea.fr |
---|
| 14 | !! |
---|
| 15 | !! Olivier Boucher, LMD |
---|
| 16 | !! introduction of a Weibull subgrid scale wind distribution |
---|
| 17 | !! |
---|
| 18 | !! This software is a computer program whose purpose is to simulate the |
---|
| 19 | !! atmospheric gas phase and aerosol composition. The model is designed to be |
---|
| 20 | !! used within a transport model or a general circulation model. This version |
---|
| 21 | !! of INCA was designed to be coupled to the LMDz GCM. LMDz-INCA accounts |
---|
| 22 | !! for emissions, transport (resolved and sub-grid scale), photochemical |
---|
| 23 | !! transformations, and scavenging (dry deposition and washout) of chemical |
---|
| 24 | !! species and aerosols interactively in the GCM. Several versions of the INCA |
---|
| 25 | !! model are currently used depending on the envisaged applications with the |
---|
| 26 | !! chemistry-climate model. |
---|
| 27 | !! |
---|
| 28 | !! This software is governed by the CeCILL license under French law and |
---|
| 29 | !! abiding by the rules of distribution of free software. You can use, |
---|
| 30 | !! modify and/ or redistribute the software under the terms of the CeCILL |
---|
| 31 | !! license as circulated by CEA, CNRS and INRIA at the following URL |
---|
| 32 | !! "http://www.cecill.info". |
---|
| 33 | !! |
---|
| 34 | !! As a counterpart to the access to the source code and rights to copy, |
---|
| 35 | !! modify and redistribute granted by the license, users are provided only |
---|
| 36 | !! with a limited warranty and the software's author, the holder of the |
---|
| 37 | !! economic rights, and the successive licensors have only limited |
---|
| 38 | !! liability. |
---|
| 39 | !! |
---|
| 40 | !! In this respect, the user's attention is drawn to the risks associated |
---|
| 41 | !! with loading, using, modifying and/or developing or reproducing the |
---|
| 42 | !! software by the user in light of its specific status of free software, |
---|
| 43 | !! that may mean that it is complicated to manipulate, and that also |
---|
| 44 | !! therefore means that it is reserved for developers and experienced |
---|
| 45 | !! professionals having in-depth computer knowledge. Users are therefore |
---|
| 46 | !! encouraged to load and test the software's suitability as regards their |
---|
| 47 | !! requirements in conditions enabling the security of their systems and/or |
---|
| 48 | !! data to be ensured and, more generally, to use and operate it in the |
---|
| 49 | !! same conditions as regards security. |
---|
| 50 | !! |
---|
| 51 | !! The fact that you are presently reading this means that you have had |
---|
| 52 | !! knowledge of the CeCILL license and that you accept its terms. |
---|
| 53 | !! ========================================================================= |
---|
| 54 | |
---|
| 55 | #include <inca_define.h> |
---|
| 56 | |
---|
| 57 | |
---|
| 58 | #ifdef AER |
---|
| 59 | SUBROUTINE DUSTSOURCE(delt) |
---|
| 60 | ! ----------------------------------------------------------------------- |
---|
| 61 | ! |
---|
| 62 | ! Author |
---|
| 63 | ! ------ |
---|
| 64 | ! Michael Schulz |
---|
| 65 | ! / Laboratoire des Sciences du Climat et de l'Environnement / Saclay |
---|
| 66 | ! 9. June 2002 |
---|
| 67 | ! |
---|
| 68 | ! purpose |
---|
| 69 | ! ------- |
---|
| 70 | ! compute interactively emission flux of mineral dust |
---|
| 71 | ! |
---|
| 72 | ! interface |
---|
| 73 | ! --------- |
---|
| 74 | ! input |
---|
| 75 | ! delt time step duration [sec] |
---|
| 76 | ! |
---|
| 77 | ! from aerosol modules |
---|
| 78 | ! zspeed wind speed at 10 m [m/s] |
---|
| 79 | ! landmask sea land mask [] |
---|
| 80 | ! ztempc surface air temperature [°C] |
---|
| 81 | ! zprecip precipitation amount [kg m-2 s-1] |
---|
| 82 | ! zprecipinsoil proxy for soil humidity [mm] |
---|
| 83 | ! rop particle density |
---|
| 84 | ! rhv source strength factor [kg s2 m-5] |
---|
| 85 | ! wth threshold velocity [m/s] |
---|
| 86 | ! cly clay content [%] |
---|
| 87 | ! |
---|
| 88 | ! srcsigma sigma of aerosol mode |
---|
| 89 | ! srcsigmaln |
---|
| 90 | ! |
---|
| 91 | ! variables saved |
---|
| 92 | ! entered input fields read flag |
---|
| 93 | ! |
---|
| 94 | ! local variables |
---|
| 95 | ! fdist number/mass factor [#/kg] |
---|
| 96 | ! nbdays soil drying rate [days/mm] |
---|
| 97 | ! mask dust emission area mask [] |
---|
| 98 | ! |
---|
| 99 | ! output |
---|
| 100 | ! eflux emission of dust per tracer [kg m-2 s-1] |
---|
| 101 | ! |
---|
| 102 | ! method |
---|
| 103 | ! ------ |
---|
| 104 | ! |
---|
| 105 | ! The dust aerosol emitted is described as a lognormal size |
---|
| 106 | ! distribution with a mass median radius of 2.5 um and a standard |
---|
| 107 | ! deviation of s = 2. Such a source size distribution was found to |
---|
| 108 | ! provide a consistent mineral aerosol transport model, which |
---|
| 109 | ! provides enough mineral dust for long range transport from the |
---|
| 110 | ! Sahara to the Atlantic and the Mediterranean. Dust plumes |
---|
| 111 | ! simulated with such a model provide aerosol optical thickness |
---|
| 112 | ! values very close to those observed with satellites [Schulz et |
---|
| 113 | ! al., 1998] [Guelle et al., 2000]. Using a fixed lognormal |
---|
| 114 | ! distribution at the source implies that number concentrations can |
---|
| 115 | ! be derived from mass emission flux and imposed source size |
---|
| 116 | ! distribution [Schulz et al., 1998]. |
---|
| 117 | ! |
---|
| 118 | ! While surface wind speed is produced interactively or read from |
---|
| 119 | ! ECMWF winds (see routine aerosol_meteo_calc), the source function |
---|
| 120 | ! incorporates also spatial information on three parameters a) the |
---|
| 121 | ! source strength factor (rhv); b) the threshold velocity (wth) and |
---|
| 122 | ! c) the clay fraction of the uppermost soil layer (cly). |
---|
| 123 | ! |
---|
| 124 | ! eflux = max(rhv*(landwind-wth)*landwind**2,0.)*landmask |
---|
| 125 | ! |
---|
| 126 | ! This work is based on an offline dust source developped by |
---|
| 127 | ! [Claquin, 1999]. He used a global TM3 dust simulation with a |
---|
| 128 | ! dustsource in all arid regions to obtain an approximate dust |
---|
| 129 | ! field. In a second step he compared regional averages of TOMS |
---|
| 130 | ! observations of aerosol index to this simulation, to obtain a |
---|
| 131 | ! source strength factor for a given region [see also Balkanski et |
---|
| 132 | ! al. [submitted , 2002]. The threshold velocity was obtained by |
---|
| 133 | ! comparing the FAO soil database to the maps of Marticorena [1995] |
---|
| 134 | ! to derive by this procedure a soil specific threshold velocity for |
---|
| 135 | ! desert soils. Combined with the high resolution FAO dataset this |
---|
| 136 | ! results in a global map of threshold velocities for dust rise in |
---|
| 137 | ! arid and semi-arid regions. The offline source formulation was |
---|
| 138 | ! converted to an interactive dust source for this work [Schulz und |
---|
| 139 | ! Timmreck, in prep] including prognostic variables produced by the |
---|
| 140 | ! GCM, such as wind speed, precipitation and surface temperature. |
---|
| 141 | ! |
---|
| 142 | ! Dust emission is thought to be inhibited, when the surface soil is |
---|
| 143 | ! wettened or frozen. Drying of desert soils can be parameterised as |
---|
| 144 | ! a function of recent accumulated precipitation, surface |
---|
| 145 | ! temperature. The clay content of the soil is a key parameter for |
---|
| 146 | ! this process, since clay retains water for a longer period. While |
---|
| 147 | ! Claquin [1999] used an offline calculation to define the periods |
---|
| 148 | ! when the soil is available for emission of dust, we prefer an |
---|
| 149 | ! interactive approach here. Soil wettening and drying rate are |
---|
| 150 | ! computed from GCM calculated preceipitation and surface |
---|
| 151 | ! temperatures. The clay content is assembled from the FAO soil |
---|
| 152 | ! data. During freezing conditions, no drying of the soil and no |
---|
| 153 | ! emissions are assumed. The parameterisation uses the conditions |
---|
| 154 | ! given by Claquin [1999] and computes an evaporation rate acting on |
---|
| 155 | ! the accumulated precipitation amount [Schulz und Timmreck, in |
---|
| 156 | ! prep]. This might overestimate the time which is needed to dry off |
---|
| 157 | ! the top soil, since no run-off or losses to the ground water |
---|
| 158 | ! reservoirs are assumed. |
---|
| 159 | ! (parameterised M. Schulz LSCE 13.2.2001) |
---|
| 160 | ! ----------------------------------------------------------------------- |
---|
| 161 | |
---|
| 162 | |
---|
| 163 | USE SPECIES_NAMES |
---|
| 164 | USE AEROSOL_METEO, only : landwind,ztempc,landmask,zprecip,zprecipinsoil |
---|
| 165 | USE AEROSOL_MOD, ONLY : srcsigma,srcsigmaln,rop,rhv,wth,cly,cimode |
---|
| 166 | USE SFLX, ONLY : eflux |
---|
| 167 | USE AEROSOL_DIAG, ONLY : DUSTE |
---|
| 168 | USE INCA_DIM |
---|
| 169 | USE MOD_INCA_PARA |
---|
| 170 | USE CHEM_CONTROLS, ONLY : lafin |
---|
| 171 | USE CONST_MOD, ONLY : pi |
---|
| 172 | USE PRINT_INCA |
---|
| 173 | USE PARAM_CHEM |
---|
| 174 | |
---|
| 175 | IMPLICIT NONE |
---|
| 176 | |
---|
| 177 | |
---|
| 178 | REAL, INTENT(in) :: delt |
---|
| 179 | |
---|
| 180 | |
---|
| 181 | ! source mass median diameter [m] |
---|
| 182 | REAL :: srcmmd_CIDUSTM = 2.5e-6 |
---|
| 183 | |
---|
| 184 | ! temp variables |
---|
| 185 | REAL :: fdist ! Number/Mass factor to compute number flux |
---|
| 186 | REAL :: nbday(PLON) |
---|
| 187 | REAL :: clyfac(PLON) |
---|
| 188 | REAL :: avgdryrate(PLON),drying(PLON) |
---|
| 189 | INTEGER :: mask(PLON) |
---|
| 190 | REAL :: efluxECP_MPI(PLON_MPI) |
---|
| 191 | REAL :: efluxECP_MPI2(PLON_MPI) |
---|
| 192 | REAL :: efluxECP(PLON) |
---|
| 193 | REAL,ALLOCATABLE :: wth_glo(:) |
---|
| 194 | REAL,ALLOCATABLE :: rhv_glo(:) |
---|
| 195 | REAL,ALLOCATABLE :: cly_glo(:) |
---|
| 196 | REAL,ALLOCATABLE, SAVE :: zprecipinsoil_glo(:) |
---|
| 197 | !$OMP THREADPRIVATE(zprecipinsoil_glo) |
---|
| 198 | INTEGER :: i,j |
---|
| 199 | LOGICAL, SAVE :: entered = .FALSE. |
---|
| 200 | LOGICAL, SAVE :: ECM_HR = .TRUE. ! Interpolation from 1°x1° |
---|
| 201 | !$OMP THREADPRIVATE (entered, ECM_HR) |
---|
| 202 | |
---|
| 203 | ! Weibull parameters and variables |
---|
| 204 | !--number of bins in weibull distribution (nwb=1 reduces to having no distribution) |
---|
| 205 | INTEGER, PARAMETER :: nwb=12 |
---|
| 206 | INTEGER kwb |
---|
| 207 | REAL auxreal, wind10ms, zlandwind, weilambda, pdfu, probu, pdfcum |
---|
| 208 | REAL, PARAMETER :: rpi = 4.*ATAN(1.) |
---|
| 209 | |
---|
| 210 | |
---|
| 211 | mask(:)=0 |
---|
| 212 | |
---|
| 213 | |
---|
| 214 | IF( .not. LMDZ_10m_winds ) THEN |
---|
| 215 | !$OMP MASTER |
---|
| 216 | CALL DUSTECMWF_XIOS(efluxECP_MPI) |
---|
| 217 | !$OMP END MASTER |
---|
| 218 | |
---|
| 219 | CALL scatter_omp(efluxECP_MPI,efluxECP) |
---|
| 220 | |
---|
| 221 | ENDIF |
---|
| 222 | |
---|
| 223 | ! read of geographical source info only once |
---|
| 224 | IF( .NOT. entered ) THEN |
---|
| 225 | |
---|
| 226 | CALL read_threshold_map() |
---|
| 227 | ALLOCATE(zprecipinsoil_glo(nbp_glo)) |
---|
| 228 | |
---|
| 229 | !$OMP MASTER |
---|
| 230 | IF (is_mpi_root) THEN |
---|
| 231 | OPEN(53,file='precipinsoil.dat',status='old',form='formatted',err=999) |
---|
| 232 | READ(53,'(G18.10)') (zprecipinsoil_glo(i),i=1,nbp_glo) |
---|
| 233 | GOTO 1000 |
---|
| 234 | 999 zprecipinsoil_glo=0. |
---|
| 235 | 1000 CLOSE(53) |
---|
| 236 | ENDIF |
---|
| 237 | !$OMP END MASTER |
---|
| 238 | |
---|
| 239 | CALL scatter(zprecipinsoil_glo,zprecipinsoil) |
---|
| 240 | DEALLOCATE(zprecipinsoil_glo) |
---|
| 241 | |
---|
| 242 | entered = .TRUE. |
---|
| 243 | |
---|
| 244 | END IF |
---|
| 245 | |
---|
| 246 | WHERE (abs(cly).LT.9990. .AND. landmask.LT.0.5 .AND. abs(wth).LT.9990.) |
---|
| 247 | mask = 1 |
---|
| 248 | ENDWHERE |
---|
| 249 | |
---|
| 250 | ! precipitation is accumulated everywhere, where dust emission is in |
---|
| 251 | ! principal possible |
---|
| 252 | zprecipinsoil = (zprecipinsoil + zprecip*delt) * mask |
---|
| 253 | |
---|
| 254 | |
---|
| 255 | ! surface soil wettening computations |
---|
| 256 | ! old function, did accumulate zprecipinsoil at the edges of the arid zones |
---|
| 257 | ! thus avoiding dust emissions there |
---|
| 258 | ! where (ztempc(:,1).gt.0. .and. mask.eq.1) |
---|
| 259 | ! clyfac=MIN(16.,cly*0.4+8.) |
---|
| 260 | ! nbday=MIN(16.,CLYFAC/CLYFAC**((ZTEMPC(:,1)-10.)/15.)) |
---|
| 261 | ! zprecipinsoil=max(0.,zprecipinsoil - delt/86400./nbday) |
---|
| 262 | ! endwhere |
---|
| 263 | ! |
---|
| 264 | ! new soil wettening formulation f(clay content and soil temperature) |
---|
| 265 | ! Dependent on clay content and interpolated to obtain |
---|
| 266 | ! clyfac, which is the maximim number of dry days after a rain event and corresponds to a |
---|
| 267 | ! maximum amount of water hold in top soil [mm] at the lowest temperatures [5°C]. |
---|
| 268 | ! The avg drying rate scales the drying to keep the area an arid area with an aridity limit |
---|
| 269 | ! assumed to be 300 mm precipitation per year. |
---|
| 270 | ! More rapid drying only depends on temperature and |
---|
| 271 | ! fits the 'dry days table' of Claquin to a function with two fitting parameters. |
---|
| 272 | |
---|
| 273 | avgdryrate=300./365.*delt/86400. ! [mm/timestep] average drying rate at lowest temperature |
---|
| 274 | WHERE (mask.EQ.1) |
---|
| 275 | clyfac=MIN(16.,cly*0.4+8.) ! [days] number of days before dry |
---|
| 276 | ! = max amount of water hold in top soil [mm] |
---|
| 277 | drying=avgdryrate*EXP(0.03905491*EXP(0.17446*ztempc(:,1))) ! [mm] |
---|
| 278 | zprecipinsoil=MIN(MAX(0.,zprecipinsoil - drying),clyfac) ! [mm] remaining water in top soil |
---|
| 279 | ENDWHERE |
---|
| 280 | |
---|
| 281 | IF (lafin) THEN |
---|
| 282 | ALLOCATE(zprecipinsoil_glo(nbp_glo)) |
---|
| 283 | zprecipinsoil_glo(:) = 0. |
---|
| 284 | CALL gather(zprecipinsoil,zprecipinsoil_glo) |
---|
| 285 | !$OMP MASTER |
---|
| 286 | IF (is_mpi_root) THEN |
---|
| 287 | OPEN(54,file='reprecipinsoil.dat',status='new',form='formatted') |
---|
| 288 | WRITE(54,'(G18.10)') (zprecipinsoil_glo(i),i=1,nbp_glo) |
---|
| 289 | CLOSE(54) |
---|
| 290 | ENDIF |
---|
| 291 | !$OMP END MASTER |
---|
| 292 | DEALLOCATE(zprecipinsoil_glo) |
---|
| 293 | ENDIF |
---|
| 294 | |
---|
| 295 | ! conversion of mass median diameter**3 to diameter of average mass**3, |
---|
| 296 | ! => exp(1.5*srcsigmaln(asmode)**2)**3 = exp(4.5*srcsigmaln(asmode)**2) |
---|
| 297 | ! Source size distribution is assumed to be constant |
---|
| 298 | fdist=1./pi*6./rop(id_CIDUSTM) & |
---|
| 299 | /srcmmd_CIDUSTM**3 *EXP(4.5*srcsigmaln(cimode)**2) |
---|
| 300 | |
---|
| 301 | |
---|
| 302 | |
---|
| 303 | !****************source function |
---|
| 304 | !--------unit kg m**-2 s-1 |
---|
| 305 | ! computation of flux |
---|
| 306 | IF( LMDZ_10m_winds ) THEN |
---|
| 307 | WHERE (ztempc(:,1).GT.0. .AND. mask.EQ.1 .AND. zprecipinsoil.LT.1e-10) |
---|
| 308 | eflux(:,id_CIDUSTM) = MAX(rhv*(landwind-wth)*landwind**2,0.) |
---|
| 309 | ENDWHERE |
---|
| 310 | ENDIF |
---|
| 311 | |
---|
| 312 | |
---|
| 313 | IF( (.not. LMDZ_10m_winds) .and. ECM_HR ) THEN |
---|
| 314 | WHERE (ztempc(:,1).GT.0. .AND. mask.EQ.1 .AND. zprecipinsoil.LT.1e-10) |
---|
| 315 | DUSTE = MAX(rhv*(landwind-wth)*landwind**2,0.) |
---|
| 316 | eflux(:,id_CIDUSTM) = efluxECP |
---|
| 317 | ENDWHERE |
---|
| 318 | |
---|
| 319 | endif |
---|
| 320 | |
---|
| 321 | |
---|
| 322 | ! and for the number flux, allowing other aerosol species sources to add |
---|
| 323 | ! to the number flux |
---|
| 324 | eflux(:,id_CIN) = eflux(:,id_CIN)+ eflux(:,id_CIDUSTM)*fdist |
---|
| 325 | |
---|
| 326 | |
---|
| 327 | RETURN |
---|
| 328 | END SUBROUTINE dustsource |
---|
| 329 | |
---|
| 330 | #endif |
---|