source: LMDZ6/trunk/libf/phylmd/StratAer/aer_sedimnt.F90 @ 5077

Last change on this file since 5077 was 4950, checked in by lebasn, 7 months ago

StratAer?: New model version (microphysic, composition routine, code cleaning, new params...)

  • Property svn:keywords set to Id
File size: 5.5 KB
Line 
1!
2! $Id: aer_sedimnt.F90 4950 2024-05-22 13:16:36Z abarral $
3!
4SUBROUTINE AER_SEDIMNT(pdtphys, t_seri, pplay, paprs, tr_seri, dens_aer)
5
6!**** *AER_SEDIMNT* -  ROUTINE FOR PARAMETRIZATION OF AEROSOL SEDIMENTATION
7
8!      Christoph Kleinschmitt
9!      based on the sedimentation scheme of
10!      Olivier Boucher & Jean-Jacques Morcrette
11!      (following the ice sedimentation scheme of Adrian Tompkins)
12
13!**   INTERFACE.
14!     ----------
15!          *AER_SEDIMNT* IS CALLED FROM *traccoag_mod*.
16
17!-----------------------------------------------------------------------
18
19  USE phys_local_var_mod, ONLY: mdw, budg_sed_part, DENSO4, DENSO4B, f_r_wet, f_r_wetB, vsed_aer
20  USE strataer_local_var_mod, ONLY: flag_new_strat_compo
21  USE dimphy, ONLY : klon,klev
22  USE infotrac_phy
23  USE aerophys
24  USE YOMCST
25
26IMPLICIT NONE
27
28!-----------------------------------------------------------------------
29
30  ! transfer variables when calling this routine
31  REAL,INTENT(IN)                             :: pdtphys ! Pas d'integration pour la physique (seconde)
32  REAL,DIMENSION(klon,klev),INTENT(IN)        :: t_seri  ! Temperature
33  REAL,DIMENSION(klon,klev),INTENT(IN)        :: pplay   ! pression pour le mileu de chaque couche (en Pa)
34  REAL,DIMENSION(klon,klev+1),INTENT(IN)      :: paprs   ! pression pour chaque inter-couche (en Pa)
35  REAL,DIMENSION(klon,klev,nbtr),INTENT(INOUT):: tr_seri ! Concentration Traceur [U/KgA]
36  REAL,DIMENSION(klon,klev)                   :: dens_aer! density of aerosol particles [kg/m3 aerosol] with default H2SO4 mass
37
38  ! local variables in sedimentation routine
39  INTEGER :: JL,JK,nb   
40  REAL,DIMENSION(klon,klev)                   :: zvis     ! dynamic viscosity of air [kg/(m*s)]
41  REAL,DIMENSION(klon,klev)                   :: zlair    ! mean free path of air [m]
42  REAL                                        :: ZRHO     ! air density [kg/m^3]
43  REAL                                        :: ZGDP     ! =g/dp=1/(rho*dz)
44  REAL                                        :: ZDTGDP   ! =dt/(rho*dz)
45  REAL,DIMENSION(klon,nbtr_bin)               :: ZSEDFLX  ! sedimentation flux of tracer [U/(m^2*s)]
46  REAL,DIMENSION(nbtr_bin)                    :: ZAERONW  ! tracer concentration at current time step [U/KgA]
47  REAL,DIMENSION(klon,nbtr_bin)               :: ZAERONWM1! tracer concentration at preceding time step [U/KgA]
48  REAL,DIMENSION(klon,klev,nbtr_bin)          :: ZVAER    ! sedimentation velocity [m/s]
49  REAL,DIMENSION(nbtr_bin)                    :: ZSOLAERS ! sedimentation flux arriving from above [U/(m^2*s)]
50  REAL,DIMENSION(nbtr_bin)                    :: ZSOLAERB ! sedimentation flux leaving gridbox [U/(m^2*s)]
51  REAL,DIMENSION(klon,klev)                   :: m_sulf
52
53! dynamic viscosity of air (Pruppacher and Klett, 1978) [kg/(m*s)]
54WHERE (t_seri.GE.273.15)
55  zvis=(1.718 + 0.0049*(t_seri-273.15))*1.E-5
56  ELSEWHERE
57  zvis=(1.718 + 0.0049*(t_seri-273.15)-1.2E-05*(t_seri-273.15)**2)*1.E-5
58END WHERE
59
60! mean free path of air (Prupp. Klett) [m]
61zlair(:,:) = 0.066 *(1.01325E+5/pplay(:,:))*(t_seri(:,:)/293.15)*1.E-06
62
63!--initialisations of variables carried out from one layer to the next layer
64!--actually not needed if (JK>1) test is on
65DO JL=1,klon
66  DO nb=1,nbtr_bin
67    ZSEDFLX(JL,nb)=0.0
68    ZAERONWM1(JL,nb)=0.0
69  ENDDO
70ENDDO
71
72!--from top to bottom (!)
73DO JK=klev,1,-1
74  DO JL=1,klon
75    DO nb=1,nbtr_bin
76  !--initialisations
77      ZSOLAERS(nb)=0.0
78      ZSOLAERB(nb)=0.0
79      ZGDP=RG/(paprs(JL,JK)-paprs(JL,JK+1))
80      ZDTGDP=pdtphys*ZGDP
81
82  ! source from above
83      IF (JK<klev) THEN
84        ZSEDFLX(JL,nb)=ZSEDFLX(JL,nb)*ZAERONWM1(JL,nb) 
85        ZSOLAERS(nb)=ZSOLAERS(nb)+ZSEDFLX(JL,nb)*ZDTGDP
86      ENDIF
87
88  ! sink to next layer
89      ZRHO=pplay(JL,JK)/(RD*t_seri(JL,JK))
90
91      ! stokes-velocity with cunnigham slip- flow correction
92      IF(flag_new_strat_compo) THEN
93         ! stokes-velocity with cunnigham slip- flow correction
94         ZVAER(JL,JK,nb) = 2./9.*(DENSO4B(JL,JK,nb)*1000.-ZRHO)*RG/zvis(JL,JK)*(f_r_wetB(JL,JK,nb)*mdw(nb)/2.)**2.* &
95              (1.+ 2.*zlair(JL,JK)/(f_r_wetB(JL,JK,nb)*mdw(nb))*(1.257+0.4*EXP(-0.55*f_r_wetB(JL,JK,nb)*mdw(nb)/zlair(JL,JK))))
96      ELSE
97         ZVAER(JL,JK,nb) = 2./9.*(DENSO4(JL,JK)*1000.-ZRHO)*RG/zvis(JL,JK)*(f_r_wet(JL,JK)*mdw(nb)/2.)**2.* &
98              (1.+ 2.*zlair(JL,JK)/(f_r_wet(JL,JK)*mdw(nb))*(1.257+0.4*EXP(-0.55*f_r_wet(JL,JK)*mdw(nb)/zlair(JL,JK))))
99      ENDIF
100     
101      ZSEDFLX(JL,nb)=ZVAER(JL,JK,nb)*ZRHO
102      ZSOLAERB(nb)=ZSOLAERB(nb)+ZDTGDP*ZSEDFLX(JL,nb)
103
104  !---implicit solver
105      ZAERONW(nb)=(tr_seri(JL,JK,nb+nbtr_sulgas)+ZSOLAERS(nb))/(1.0+ZSOLAERB(nb))
106
107  !---new time-step AER variable needed for next layer
108      ZAERONWM1(JL,nb)=ZAERONW(nb)
109
110      tr_seri(JL,JK,nb+nbtr_sulgas)=ZAERONWM1(JL,nb)
111    ENDDO
112  ENDDO
113ENDDO
114
115!---sedimentation flux to the surface
116!---ZAERONWM1 now contains the surface concentration at the new timestep
117!---PFLUXAER in unit of xx m-2 s-1
118budg_sed_part(:)=0.0
119DO JL=1,klon
120  ZRHO=pplay(JL,1)/(RD*t_seri(JL,1))
121  DO nb=1,nbtr_bin
122    !compute budg_sed_part as sum over bins in kg(S)/m2/s
123    budg_sed_part(JL)=budg_sed_part(JL)+ZRHO*ZAERONWM1(JL,nb)*ZVAER(JL,1,nb)*(mSatom/mH2SO4mol) &
124                & *dens_aer_dry*4./3.*RPI*(mdw(nb)/2.)**3
125  ENDDO
126ENDDO
127
128vsed_aer(:,:)=0.0
129m_sulf(:,:)=0.0
130
131DO nb=1,nbtr_bin
132  !compute mass-weighted mean of sedimentation velocity
133  vsed_aer(:,:)=vsed_aer(:,:)+ZVAER(:,:,nb)*(mdw(nb)/2.)**3*MAX(1.e-30, tr_seri(:,:,nb+nbtr_sulgas))
134  m_sulf(:,:)=m_sulf(:,:)+(mdw(nb)/2.)**3*MAX(1.e-30, tr_seri(:,:,nb+nbtr_sulgas))
135ENDDO
136
137!divide by total aerosol mass in grid cell
138vsed_aer(:,:)=vsed_aer(:,:)/m_sulf(:,:)
139
140END SUBROUTINE AER_SEDIMNT
Note: See TracBrowser for help on using the repository browser.