source: LMDZ6/trunk/libf/phylmd/lmdz_wake_popdyn_2.f90 @ 5878

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

Séparation de lmdz_wake en petits fichiers (JYG&FH)
File size: 7.2 KB
Line 
1MODULE lmdz_wake_popdyn_2
2PUBLIC wake_popdyn_2
3CONTAINS
4
5    SUBROUTINE wake_popdyn_2 ( klon, klev, wk_adv, dtimesub, wgen, &
6                             wdensmin, &
7                             sigmaw, wdens, awdens, &   !! states variables
8                             gfl, cstar, cin, wape, rad_wk, &
9                             d_sigmaw, d_wdens, d_awdens, &  !! tendences
10                             cont_fact, &
11                             d_sig_gen, d_sig_death, d_sig_col, d_sig_spread, d_sig_bnd, &
12                             d_dens_gen, d_dens_death, d_dens_col, d_dens_bnd, &
13                             d_adens_death, d_adens_icol, d_adens_acol, d_adens_bnd )
14                             
15                                             
16
17  USE lmdz_wake_ini , ONLY : wake_ini
18  USE lmdz_wake_ini , ONLY : prt_level,RG
19  USE lmdz_wake_ini , ONLY : stark, wdens_ref
20  USE lmdz_wake_ini , ONLY : tau_cv, rzero, aa0
21!!  USE lmdz_wake_ini , ONLY : iflag_wk_pop_dyn, wdensmin
22  USE lmdz_wake_ini , ONLY : iflag_wk_pop_dyn
23  USE lmdz_wake_ini , ONLY : sigmad, cstart, sigmaw_max
24 
25IMPLICIT NONE
26
27  INTEGER, INTENT(IN)                                   :: klon,klev
28  LOGICAL, DIMENSION (klon),        INTENT(IN)          :: wk_adv
29  REAL,                             INTENT(IN)          :: dtimesub
30  REAL,                             INTENT(IN)          :: wdensmin
31  REAL, DIMENSION (klon),           INTENT(IN)          :: wgen      !! B = birth rate of wakes
32  REAL, DIMENSION (klon),           INTENT(INOUT)       :: sigmaw    !! sigma = fractional area of wakes
33  REAL, DIMENSION (klon),           INTENT(INOUT)       :: wdens     !! D = number of wakes per unit area
34  REAL, DIMENSION (klon),           INTENT(INOUT)       :: awdens    !! A = number of active wakes per unit area
35  REAL, DIMENSION (klon),           INTENT(IN)          :: cstar     !! C* = spreading velocity of wakes
36  REAL, DIMENSION (klon),           INTENT(IN)          :: cin, wape  ! RM : A Faire disparaitre
37
38! 
39  REAL, DIMENSION (klon),           INTENT(OUT)         :: rad_wk    !! r = wake radius
40  REAL, DIMENSION (klon),           INTENT(OUT)         :: gfl       !! Lg = gust front lenght per unit area
41  REAL, DIMENSION (klon),           INTENT(OUT)         :: d_sigmaw, d_wdens, d_awdens
42  REAL, DIMENSION (klon),           INTENT(OUT)         :: cont_fact  !! RM facteur de contact = 2 pi * rad * C*
43  ! Some components of the tendencies of state variables 
44  REAL, DIMENSION (klon),           INTENT(OUT)         :: d_sig_gen, d_sig_death, d_sig_col, d_sig_spread, d_sig_bnd
45  REAL, DIMENSION (klon),           INTENT(OUT)         :: d_dens_gen, d_dens_death, d_dens_col, d_dens_bnd
46  REAL, DIMENSION (klon),           INTENT(OUT)         :: d_adens_death, d_adens_icol, d_adens_acol, d_adens_bnd
47
48
49!! internal variables
50 
51  INTEGER                                               :: i, k
52  REAL, DIMENSION (klon)                                :: tau_wk_inv      !! tau = life time of wakes
53  REAL                                                  :: tau_wk_inv_min
54  REAL, DIMENSION (klon)                                :: tau_prime       !! tau_prime = life time of actives wakes
55  REAL                                                  :: d_wdens_targ, d_sigmaw_targ
56 
57
58!! Equations
59!! dD/dt = B - (D-A)/tau - f D^2
60!! dA/dt = B - A/tau_prime + f (D-A)^2 - f A^2
61!! dsigma/dt = B a0 - sigma/D (D-A)/tau + Lg C* - f (D-A)^2 (sigma/D-a0)
62!!
63!! f = 2 (B (a0-sigma/D) + Lg C*) / (2 (D-A)^2 (2 sigma/D-a0) + D (1-2 sigma))
64
65
66      DO i = 1, klon
67        IF (wk_adv(i)) THEN
68          rad_wk(i) = max( sqrt(sigmaw(i)/(3.14*wdens(i))) , rzero)
69          gfl(i)  = 2.*sqrt(3.14*wdens(i)*sigmaw(i))
70        END IF
71      END DO
72
73
74      DO i = 1, klon
75        IF (wk_adv(i)) THEN
76!!          tau_wk(i) = max(rad_wk(i)/(3.*cstar(i))*((cstar(i)/cstart)**1.5 - 1), 100.)
77          tau_wk_inv(i) = max( (3.*cstar(i))/(rad_wk(i)*((cstar(i)/cstart)**1.5 - 1)), 0.)
78          tau_wk_inv_min = min(tau_wk_inv(i), 1./dtimesub)
79          tau_prime(i) = tau_cv
80!!          cont_fact(i) = 2.*(wgen(i)*(aa0-sigmaw(i)/wdens(i)) + gfl(i)*cstar(i)) / &
81!!                             (2.*(wdens(i)-awdens(i))**2*(2.*sigmaw(i)/wdens(i) - aa0) + wdens(i)*(1.-2.*sigmaw(i)))
82!!          cont_fact(i) = 2.*3.14*rad_wk(i)*cstar(i)     ! bug
83!!          cont_fact(i) = 4.*3.14*rad_wk(i)*cstar(i)
84          cont_fact(i) = 2.*gfl(i)*cstar(i)/wdens(i)
85
86          d_sig_gen(i) = wgen(i)*aa0
87          d_sig_death(i) = - sigmaw(i)*(1.-awdens(i)/wdens(i))*tau_wk_inv_min
88          d_sig_col(i) = - cont_fact(i)*(wdens(i)-awdens(i))**2*(2.*sigmaw(i)/wdens(i)-aa0)
89          d_sig_spread(i) = gfl(i)*cstar(i)
90!
91          d_sig_gen(i) =  d_sig_gen(i)*dtimesub
92          d_sig_death(i) = d_sig_death(i)*dtimesub
93          d_sig_col(i) =  d_sig_col(i)*dtimesub
94          d_sig_spread(i) =  d_sig_spread(i)*dtimesub
95          d_sigmaw(i) =  d_sig_gen(i) + d_sig_death(i) + d_sig_col(i) + d_sig_spread(i)
96
97         
98          d_sigmaw_targ = max(d_sigmaw(i), sigmad-sigmaw(i))
99!!          d_sig_bnd(i) = d_sig_bnd(i) + d_sigmaw_targ - d_sigmaw(i)
100!!          d_sig_bnd_provis(i) = d_sigmaw_targ - d_sigmaw(i)
101          d_sig_bnd(i) = d_sigmaw_targ - d_sigmaw(i)
102          d_sigmaw(i) = d_sigmaw_targ
103!!          d_sigmaw(i) = max(d_sigmaw(i), sigmad-sigmaw(i))
104         
105         
106          d_dens_gen(i) = wgen(i)
107          d_dens_death(i) = - (wdens(i)-awdens(i))*tau_wk_inv_min
108          d_dens_col(i) =  - cont_fact(i)*wdens(i)**2
109!
110          d_dens_gen(i) =  d_dens_gen(i)*dtimesub
111          d_dens_death(i) = d_dens_death(i)*dtimesub
112          d_dens_col(i) =  d_dens_col(i)*dtimesub
113          d_wdens(i) = d_dens_gen(i) + d_dens_death(i) + d_dens_col(i)
114
115
116          d_adens_death(i) = -awdens(i)/tau_prime(i)
117          d_adens_icol(i) = cont_fact(i)*(wdens(i)-awdens(i))**2
118          d_adens_acol(i) = - cont_fact(i)*awdens(i)**2
119!
120          d_adens_death(i) =  d_adens_death(i)*dtimesub
121          d_adens_icol(i) =   d_adens_icol(i)*dtimesub
122          d_adens_acol(i) =   d_adens_acol(i)*dtimesub
123          d_awdens(i) =   d_dens_gen(i) + d_adens_death(i) + d_adens_icol(i) + d_adens_acol(i)     
124           
125!!
126          d_wdens_targ = max(d_wdens(i), wdensmin-wdens(i))
127!!          d_dens_bnd(i) = d_dens_bnd(i) + d_wdens_targ - d_wdens(i)
128          d_dens_bnd(i) = d_wdens_targ - d_wdens(i)
129          d_wdens(i) = d_wdens_targ
130         
131          d_wdens_targ = min(max(d_awdens(i),-awdens(i)), wdens(i)-awdens(i))
132!!          d_dens_bnd(i) = d_dens_bnd(i) + d_wdens_targ - d_wdens(i)
133          d_adens_bnd(i) = d_wdens_targ - d_awdens(i)
134          d_awdens(i) = d_wdens_targ
135
136
137
138        ENDIF
139      ENDDO
140
141      IF (prt_level >= 10) THEN
142        print *,'wake, cstar(1), cstar(1)/cstart, rad_wk(1), tau_wk_inv(1), cont_fact(1) ', &
143                       cstar(1), cstar(1)/cstart, rad_wk(1), tau_wk_inv(1), cont_fact(1)
144        print *,'wake, wdens(1), awdens(1), d_awdens(1) ', &
145                       wdens(1), awdens(1), d_awdens(1)
146        print *,'wake, d_sig_gen(1), d_sig_death(1), d_sig_col(1), d_sigmaw(1) ', &
147                       d_sig_gen(1), d_sig_death(1), d_sig_col(1), d_sigmaw(1)
148      ENDIF
149sigmaw=sigmaw+d_sigmaw
150wdens=wdens+d_wdens
151awdens=awdens+d_awdens
152
153    RETURN
154    END SUBROUTINE wake_popdyn_2 
155END MODULE lmdz_wake_popdyn_2
Note: See TracBrowser for help on using the repository browser.