source: LMDZ6/trunk/libf/phylmd/lmdz_wake_popdyn_3.f90 @ 5804

Last change on this file since 5804 was 5804, checked in by fhourdin, 3 months ago

Séparation de lmdz_wake en petits fichiers (JYG&FH)
File size: 13.0 KB
RevLine 
[5804]1MODULE lmdz_wake_popdyn_3
2PUBLIC wake_popdyn_3
3CONTAINS
4
5    SUBROUTINE wake_popdyn_3 ( klon, klev, phys_sub, wk_adv, dtimesub, wgen, &
6                             wdensmin, &
7                             sigmaw, asigmaw, wdens, awdens, &                       !! state variables
8                             gfl, agfl, cstar, cin, wape, &
9                             rad_wk, arad_wk, irad_wk, &
10                             d_sigmaw, d_asigmaw, d_wdens, d_awdens, &               !! tendencies
11                             d_sig_gen, d_sig_death, d_sig_col, d_sig_spread, d_sig_bnd, &
12                             d_asig_death, d_asig_aicol, d_asig_iicol, d_asig_spread, d_asig_bnd, &
13                             d_dens_gen, d_dens_death, d_dens_col, d_dens_bnd, &
14                             d_adens_death, d_adens_icol, d_adens_acol, d_adens_bnd )
15                             
16                                             
17
18         USE lmdz_wake_ini, ONLY: CPPKEY_IOPHYS_WK
19  USE lmdz_wake_ini , ONLY : wake_ini
20  USE lmdz_wake_ini , ONLY : prt_level,RG
21  USE lmdz_wake_ini , ONLY : stark, wdens_ref
22  USE lmdz_wake_ini , ONLY : tau_cv, rzero, aa0
23!!  USE lmdz_wake_ini , ONLY : iflag_wk_pop_dyn, wdensmin
24  USE lmdz_wake_ini , ONLY : iflag_wk_pop_dyn
25  USE lmdz_wake_ini , ONLY : sigmad, cstart, sigmaw_max
26  USE lmdz_wake_ini , ONLY : smallestreal
27 
28IMPLICIT NONE
29
30  INTEGER, INTENT(IN)                                   :: klon,klev
31  LOGICAL,                          INTENT(IN)          :: phys_sub
32  LOGICAL, DIMENSION (klon),        INTENT(IN)          :: wk_adv
33  REAL,                             INTENT(IN)          :: dtimesub
34  REAL,                             INTENT(IN)          :: wdensmin
35  REAL, DIMENSION (klon),           INTENT(IN)          :: wgen      !! B = birth rate of wakes
36  REAL, DIMENSION (klon),           INTENT(INOUT)       :: sigmaw    !! sigma = fractional area of wakes
37  REAL, DIMENSION (klon),           INTENT(INOUT)       :: asigmaw   !! sigma = fractional area of active wakes
38  REAL, DIMENSION (klon),           INTENT(INOUT)       :: wdens     !! D = number of wakes per unit area
39  REAL, DIMENSION (klon),           INTENT(INOUT)       :: awdens    !! A = number of active wakes per unit area
40  REAL, DIMENSION (klon),           INTENT(IN)          :: cstar     !! C* = spreading velocity of wakes
41  REAL, DIMENSION (klon),           INTENT(IN)          :: cin, wape  ! RM : A Faire disparaitre
42
43! 
44  REAL, DIMENSION (klon),           INTENT(OUT)         :: rad_wk    !! r = mean wake radius
45  REAL, DIMENSION (klon),           INTENT(OUT)         :: arad_wk    !! r_A = wake radius of active wakes
46  REAL, DIMENSION (klon),           INTENT(OUT)         :: irad_wk    !! r_I = wake radius of inactive wakes
47  REAL, DIMENSION (klon),           INTENT(OUT)         :: gfl       !! Lg = gust front length per unit area
48  REAL, DIMENSION (klon),           INTENT(OUT)         :: agfl      !! LgA = gust front length of active wakes
49                                                                     !!  per unit area
50  REAL, DIMENSION (klon),           INTENT(OUT)         :: d_sigmaw, d_asigmaw, d_wdens, d_awdens
51  ! Some components of the tendencies of state variables 
52  REAL, DIMENSION (klon),           INTENT(OUT)         :: d_sig_gen, d_sig_death, d_sig_col, d_sig_spread, d_sig_bnd
53  REAL, DIMENSION (klon),           INTENT(OUT)         :: d_asig_death, d_asig_aicol, d_asig_iicol, d_asig_spread, d_asig_bnd
54  REAL, DIMENSION (klon),           INTENT(OUT)         :: d_dens_gen, d_dens_death, d_dens_col, d_dens_bnd
55  REAL, DIMENSION (klon),           INTENT(OUT)         :: d_adens_death, d_adens_acol, d_adens_icol, d_adens_bnd
56
57
58!! internal variables
59 
60  INTEGER                                               :: i, k
61  REAL, DIMENSION (klon)                                :: iwdens, isigmaw !! inactive wake density and fractional area
62!!  REAL, DIMENSION (klon)                                :: d_arad, d_irad
63  REAL, DIMENSION (klon)                                :: igfl            !! LgI = gust front length of inactive wakes
64                                                                           !!  per unit area
65  REAL, DIMENSION (klon)                                :: s_wk            !! mean area of individual wakes
66  REAL, DIMENSION (klon)                                :: as_wk           !! mean area of individual active wakes
67  REAL, DIMENSION (klon)                                :: is_wk           !! mean area of individual inactive wakes
68  REAL, DIMENSION (klon)                                :: tau_wk_inv      !! tau = life time of wakes
69  REAL, DIMENSION (klon)                                :: tau_prime       !! tau_prime = life time of actives wakes
70  REAL                                                  :: d_wdens_targ, d_sigmaw_targ
71
72
73!! Equations
74!! ---------
75!! Gust fronts:
76!! Lg_A = 2 pi r_A A
77!! Lg_I = 2 pi r_I I
78!! Lg   = 2 pi r   D
79!!
80!! Areas:
81!! s = pi r^2
82!! s_A = pi r_A^2
83!! s_I = pi r_I^2
84!!
85!! Life expectancy:
86!! tau_I = 3 C* ((C*/C*t)^3/2 - 1) / r_I
87!!
88!! Time deratives:
89!! dD/dt = B - (D-A)/tau_I - 2 Lg C* D
90!! dA/dt = B - A/tau_A     + 2 Lg_I C* (D-A) - 2 Lg_A C* A
91!! dsigma/dt = B a0 - sigma_I/tau_I + Lg C* - 2 Lg_I C* (D-A) (2 s_I - a0)
92!! dsigma_A/dt = B a0 - sigma_A/tau_A + Lg_A C* + (Lg_A I + Lg_I A) C* s_I + 2 Lg_I C* I a0
93!!
94
95! Initialization
96 tau_wk_inv(:) = 0.
97! Initialization of output variables
98 rad_wk(:) = 0.
99 arad_wk(:) = 0.
100 irad_wk(:) = 0.
101 gfl(:) = 0.
102 agfl(:) = 0.
103!
104 d_wdens(:) = 0.
105 d_awdens(:) = 0.
106 d_sigmaw(:) = 0.
107 d_asigmaw (:) = 0.
108!
109 d_sig_gen(:) = 0.
110 d_sig_death(:) = 0.
111 d_sig_col(:) = 0.
112 d_sig_spread(:) = 0.
113 d_sig_bnd(:) = 0.
114 d_asig_death(:) = 0.
115 d_asig_aicol(:) = 0.
116 d_asig_iicol(:) = 0.
117 d_asig_spread(:) = 0.
118 d_asig_bnd(:) = 0.
119 d_dens_gen(:) = 0.
120 d_dens_death(:) = 0.
121 d_dens_col(:) = 0.
122 d_dens_bnd(:) = 0.
123 d_adens_death(:) = 0.
124 d_adens_icol(:) = 0.
125 d_adens_acol(:) = 0.
126 d_adens_bnd(:) = 0.
127
128
129      DO i = 1, klon
130        IF (wk_adv(i)) THEN
131         iwdens(i) = wdens(i) - awdens(i)
132         isigmaw(i) = sigmaw(i) - asigmaw(i)
133!
134         arad_wk(i) = max( sqrt(asigmaw(i)/(3.14*awdens(i))) , rzero)
135         irad_wk(i) = max( sqrt((sigmaw(i)-asigmaw(i))/  &
136                           (3.14*max(smallestreal,(wdens(i)-awdens(i))))), rzero)
137         rad_wk(i) = (awdens(i)*arad_wk(i)+(wdens(i)-awdens(i))*irad_wk(i))/wdens(i)
138!
139         s_wk(i) = 3.14*rad_wk(i)**2
140         as_wk(i) = 3.14*arad_wk(i)**2
141         is_wk(i) = 3.14*irad_wk(i)**2
142!
143         gfl(i)  = 2.*sqrt(3.14*wdens(i)*sigmaw(i))
144         agfl(i) = 2.*sqrt(3.14*awdens(i)*asigmaw(i))
145         igfl(i) = gfl(i) - agfl(i)
146        ENDIF
147      ENDDO
148
149
150      DO i = 1, klon
151        IF (wk_adv(i)) THEN
152!  print *,'ZZZZpopdyn3 wgen(1) ',wgen(1)
153!  print *,'ZZZZpopdyn3 cstar(1) ',cstar(1)
154!  print *,'ZZZZpopdyn3 isigmaw(1) ',isigmaw(1)
155!  print *,'ZZZZpopdyn3 gfl(1) ',gfl(1)
156!!          tau_wk_inv(i) = max( (3.*cstar(i))/(irad_wk(i)*((cstar(i)/cstart)**1.5 - 1)), 0.)
157          tau_wk_inv(i) = min(max( (3.*cstar(i))/(irad_wk(i)*((cstar(i)/cstart)**1.5 - 1)), 0.), 1./dtimesub)
158          tau_prime(i) = tau_cv
159
160          d_sig_gen(i) = wgen(i)*aa0
161          d_sig_death(i) = - isigmaw(i)*tau_wk_inv(i)
162          d_sig_col(i) = - 2.*igfl(i)*cstar(i)*iwdens(i)*(2.*is_wk(i)-aa0)
163          d_sig_spread(i) = gfl(i)*cstar(i)
164!
165          d_sig_gen(i) =  d_sig_gen(i)*dtimesub
166          d_sig_death(i) = d_sig_death(i)*dtimesub
167          d_sig_col(i) =  d_sig_col(i)*dtimesub
168          d_sig_spread(i) =  d_sig_spread(i)*dtimesub
169          d_sigmaw(i) =  d_sig_gen(i) + d_sig_death(i) + d_sig_col(i) + d_sig_spread(i)
170         
171          d_sigmaw_targ = max(d_sigmaw(i), sigmad-sigmaw(i))
172!!          d_sig_bnd(i) = d_sig_bnd(i) + d_sigmaw_targ - d_sigmaw(i)
173!!          d_sig_bnd_provis(i) = d_sigmaw_targ - d_sigmaw(i)
174          d_sig_bnd(i) = d_sigmaw_targ - d_sigmaw(i)
175          d_sigmaw(i) = d_sigmaw_targ
176!!          d_sigmaw(i) = max(d_sigmaw(i), sigmad-sigmaw(i))
177          d_asig_death(i) = - asigmaw(i)/tau_prime(i)
178!!  Bug : factor 2 omitted by mistake (bug found by Lamine Thiam)
179!!          d_asig_aicol(i) = (agfl(i)*iwdens(i) + igfl(i)*awdens(i))*cstar(i)*is_wk(i)
180          d_asig_aicol(i) = 2.*(agfl(i)*iwdens(i) + igfl(i)*awdens(i))*cstar(i)*is_wk(i)
181          d_asig_iicol(i) = 2.*igfl(i)*cstar(i)*iwdens(i)*aa0
182          d_asig_spread(i) = agfl(i)*cstar(i)
183!
184          d_asig_death(i) = d_asig_death(i)*dtimesub
185          d_asig_aicol(i) =  d_asig_aicol(i)*dtimesub
186          d_asig_iicol(i) =  d_asig_iicol(i)*dtimesub
187          d_asig_spread(i) =  d_asig_spread(i)*dtimesub
188          d_asigmaw(i) =  d_sig_gen(i) + d_asig_death(i) + d_asig_aicol(i) + d_asig_iicol(i) + d_asig_spread(i)
189!
190          d_sigmaw_targ = min(max(d_asigmaw(i),-asigmaw(i)), sigmaw(i)-asigmaw(i))
191!!          d_dens_bnd(i) = d_dens_bnd(i) + d_sigmaw_targ - d_sigmaw(i)
192          d_asig_bnd(i) = d_sigmaw_targ - d_asigmaw(i)
193          d_asigmaw(i) = d_sigmaw_targ
194          d_dens_gen(i) = wgen(i)
195          d_dens_death(i) = - iwdens(i)*tau_wk_inv(i)
196          d_dens_col(i) =  - 2.*gfl(i)*cstar(i)*wdens(i)
197!
198          d_dens_gen(i) =  d_dens_gen(i)*dtimesub
199          d_dens_death(i) = d_dens_death(i)*dtimesub
200          d_dens_col(i) =  d_dens_col(i)*dtimesub
201          d_wdens(i) = d_dens_gen(i) + d_dens_death(i) + d_dens_col(i)
202!!
203          d_wdens_targ = max(d_wdens(i), wdensmin-wdens(i))
204!!          d_dens_bnd(i) = d_dens_bnd(i) + d_wdens_targ - d_wdens(i)
205          d_dens_bnd(i) = d_wdens_targ - d_wdens(i)
206          d_wdens(i) = d_wdens_targ
207          d_adens_death(i) = -awdens(i)/tau_prime(i)
208          d_adens_icol(i) =   2.*igfl(i)*cstar(i)*iwdens(i)
209          d_adens_acol(i)  = - 2.*agfl(i)*cstar(i)*awdens(i)
210!
211          d_adens_death(i) =  d_adens_death(i)*dtimesub
212          d_adens_icol(i) =   d_adens_icol(i)*dtimesub
213          d_adens_acol(i)  =   d_adens_acol(i)*dtimesub
214          d_awdens(i) =   d_dens_gen(i) + d_adens_death(i) + d_adens_icol(i) + d_adens_acol(i)     
215          d_wdens_targ = min(max(d_awdens(i),-awdens(i)), wdens(i)-awdens(i))
216!!          d_dens_bnd(i) = d_dens_bnd(i) + d_wdens_targ - d_wdens(i)
217          d_adens_bnd(i) = d_wdens_targ - d_awdens(i)
218          d_awdens(i) = d_wdens_targ
219
220!!          d_irad(i) = (d_sigmaw(i)-d_asigmaw(i)-isigmaw(i)*(d_wdens(i)-awdens(i))/iwdens(i)) / &
221!!                      max(smallestreal,(2.*3.14*iwdens(i)*irad_wk(i)))
222!!          d_arad(i) = (d_asigmaw(i)-asigmaw(i)*d_awdens(i)/awdens(i)) / &
223!!                      max(smallestreal,(2.*3.14*awdens(i)*arad_wk(i)))
224!!          d_irad(i) = d_irad(i)*dtimesub
225!!          d_arad(i) = d_arad(i)*dtimesub
226!!          call iophys_ecrit('d_irad',1,'d_irad','m',d_irad)
227!!          call iophys_ecrit('d_airad',1,'d_arad','m',d_arad)
228!!
229        ENDIF
230      ENDDO
231IF (CPPKEY_IOPHYS_WK) THEN
232    IF (phys_sub) THEN
233       call iophys_ecrit('d_sigmaw0',1,'d_sigmaw0','',d_sigmaw)
234!
235       call iophys_ecrit('cstar',1,'cstar','',cstar)
236       call iophys_ecrit('wgen_pd3',1,'wgen_popdyn3','',wgen)
237       call iophys_ecrit('tauwk_inv',1,'tau_wk_inv','',tau_wk_inv)
238       call iophys_ecrit('d_sigmaw',1,'d_sigmaw','',d_sigmaw)
239       call iophys_ecrit('d_sig_gen',1,'d_sig_gen','',d_sig_gen)
240       call iophys_ecrit('d_sig_death',1,'d_sig_death','',d_sig_death)
241       call iophys_ecrit('d_sig_col',1,'d_sig_col','',d_sig_col)
242       call iophys_ecrit('d_sig_spread',1,'d_sig_spread','',d_sig_spread)
243       call iophys_ecrit('d_sig_bnd',1,'d_sig_bnd','',d_sig_bnd)
244!
245       call iophys_ecrit('d_asigmaw0',1,'d_asigmaw0','',d_asigmaw)
246!
247       call iophys_ecrit('d_asigmaw',1,'d_asigmaw','',d_asigmaw)
248       call iophys_ecrit('d_asig_death',1,'d_asig_death','',d_asig_death)
249       call iophys_ecrit('d_asig_aicol',1,'d_asig_aicol','',d_asig_aicol)
250       call iophys_ecrit('d_asig_iicol',1,'d_asig_iicol','',d_asig_iicol)
251       call iophys_ecrit('d_asig_spread',1,'d_asig_spread','',d_asig_spread)
252       call iophys_ecrit('d_asig_bnd',1,'d_asig_bnd','',d_asig_bnd)
253!
254       call iophys_ecrit('d_wdens',1,'d_wdens','',d_wdens)
255       call iophys_ecrit('d_dens_gen',1,'d_dens_gen','',d_dens_gen)
256       call iophys_ecrit('d_dens_death',1,'d_dens_death','',d_dens_death)
257       call iophys_ecrit('d_dens_col',1,'d_dens_col','',d_dens_col)
258! 
259       call iophys_ecrit('d_awdens',1,'d_awdens','',d_awdens)
260       call iophys_ecrit('d_adens_death',1,'d_adens_death','',d_adens_death)
261       call iophys_ecrit('d_adens_icol',1,'d_adens_icol','',d_adens_icol)
262       call iophys_ecrit('d_adens_acol',1,'d_adens_acol','',d_adens_acol)
263    ENDIF
264END IF
265
266
267      IF (prt_level >= 10) THEN
268        print *,'wake, cstar(1), cstar(1)/cstart, rad_wk(1), tau_wk_inv(1), gfl(1) ', &
269                       cstar(1), cstar(1)/cstart, rad_wk(1), tau_wk_inv(1), gfl(1)
270        print *,'wake, wdens(1), awdens(1), d_awdens(1) ', &
271                       wdens(1), awdens(1), d_awdens(1)
272        print *,'wake, d_sig_gen(1), d_sig_death(1), d_sig_col(1), d_sigmaw(1) ', &
273                       d_sig_gen(1), d_sig_death(1), d_sig_col(1), d_sigmaw(1)
274      ENDIF
275sigmaw=sigmaw+d_sigmaw
276asigmaw=asigmaw+d_asigmaw
277wdens=wdens+d_wdens
278awdens=awdens+d_awdens
279
280    RETURN
281    END SUBROUTINE wake_popdyn_3 
282END MODULE lmdz_wake_popdyn_3
Note: See TracBrowser for help on using the repository browser.