source: LMDZ5/trunk/libf/phylmd/alpale.F90 @ 4947

Last change on this file since 4947 was 2554, checked in by jyg, 8 years ago

1/ Correcting a bug in alpale.F90: now ale_wk and
alp_wk are arguments of subroutine alpale;
their values are present in the output files.

2/ Minor changes to alpale_th.F90: (i) the
(commented) part relative to TKE transport is put
back in "physiq"; (ii) the part dealing with
iflag_coupl=2 terminates with a STOP.

File size: 8.1 KB
Line 
1SUBROUTINE alpale ( debut, itap, dtime, paprs, omega, t_seri,   &
2                    alp_offset, it_wape_prescr,  wape_prescr, fip_prescr, &
3                    ale_bl_prescr, alp_bl_prescr, &
4                    wake_pe, wake_fip,  &
5                    Ale_bl, Ale_bl_trig, Alp_bl, &
6                    Ale, Alp, Ale_wake, Alp_wake )
7
8! **************************************************************
9! *
10! ALPALE                                                       *
11! *
12! *
13! written by   : Jean-Yves Grandpeix, 12/05/2016              *
14! modified by :                                               *
15! **************************************************************
16
17  USE dimphy
18  USE ioipsl_getin_p_mod, ONLY : getin_p
19  USE print_control_mod, ONLY: mydebug=>debug , lunout, prt_level
20  USE phys_local_var_mod, ONLY: zw2       ! Variables internes non sauvegardees de la physique
21!
22  IMPLICIT NONE
23
24!================================================================
25! Auteur(s)   : Jean-Yves Grandpeix, 12/05/2016
26! Objet : Sums up all contributions to Ale and Alp
27!================================================================
28
29! Input arguments
30!----------------
31  LOGICAL, INTENT(IN)                                        :: debut
32  INTEGER, INTENT(IN)                                        :: itap
33  REAL, INTENT(IN)                                           :: dtime
34  INTEGER, INTENT(IN)                                        :: it_wape_prescr
35  REAL, INTENT(IN)                                           :: wape_prescr, fip_prescr
36  REAL, INTENT(IN)                                           :: Ale_bl_prescr, Alp_bl_prescr
37  REAL, INTENT(IN)                                           :: alp_offset
38  REAL, DIMENSION(klon,klev+1), INTENT(IN)                   :: paprs
39  REAL, DIMENSION(klon,klev), INTENT(IN)                     :: t_seri
40  REAL, DIMENSION(klon,klev), INTENT(IN)                     :: omega
41  REAL, DIMENSION(klon), INTENT(IN)                          :: wake_pe, wake_fip
42  REAL, DIMENSION(klon), INTENT(IN)                          :: Ale_bl, Ale_bl_trig, Alp_bl
43
44
45! Output arguments
46!----------------
47  REAL, DIMENSION(klon), INTENT(OUT)                         :: Ale, Alp
48  REAL, DIMENSION(klon), INTENT(OUT)                         :: Ale_wake, Alp_wake
49
50  include "thermcell.h"
51  include "YOMCST.h"
52  include "YOETHF.h"
53
54! Local variables
55!----------------
56  INTEGER                                                    :: i, k
57  REAL, DIMENSION(klon)                                      :: www
58  REAL, SAVE                                                 :: ale_max=1000.
59  REAL, SAVE                                                 :: alp_max=2.
60  CHARACTER*20 modname
61  CHARACTER*80 abort_message
62
63
64    !$OMP THREADPRIVATE(ale_max,alp_max)
65
66       ! Calcul de l'energie disponible ALE (J/kg) et de la puissance
67       ! disponible ALP (W/m2) pour le soulevement des particules dans
68       ! le modele convectif
69       !
70       do i = 1,klon
71          ALE(i) = 0.
72          ALP(i) = 0.
73       enddo
74       !
75       !calcul de ale_wake et alp_wake
76       if (iflag_wake>=1) then
77          if (itap .le. it_wape_prescr) then
78             do i = 1,klon
79                ale_wake(i) = wape_prescr
80                alp_wake(i) = fip_prescr
81             enddo
82          else
83             do i = 1,klon
84                !jyg  ALE=WAPE au lieu de ALE = 1/2 Cstar**2
85                !cc           ale_wake(i) = 0.5*wake_cstar(i)**2
86                ale_wake(i) = wake_pe(i)
87                alp_wake(i) = wake_fip(i)
88             enddo
89          endif
90       else
91          do i = 1,klon
92             ale_wake(i) = 0.
93             alp_wake(i) = 0.
94          enddo
95       endif
96       !combinaison avec ale et alp de couche limite: constantes si pas
97       !de couplage, valeurs calculees dans le thermique sinon
98       if (iflag_coupl.eq.0) then
99          if (debut.and.prt_level.gt.9) &
100               WRITE(lunout,*)'ALE et ALP imposes'
101          do i = 1,klon
102             !on ne couple que ale
103             !           ALE(i) = max(ale_wake(i),Ale_bl(i))
104             ALE(i) = max(ale_wake(i),ale_bl_prescr)
105             !on ne couple que alp
106             !           ALP(i) = alp_wake(i) + Alp_bl(i)
107             ALP(i) = alp_wake(i) + alp_bl_prescr
108          enddo
109       else
110          IF(prt_level>9)WRITE(lunout,*)'ALE et ALP couples au thermique'
111          !         do i = 1,klon
112          !             ALE(i) = max(ale_wake(i),Ale_bl(i))
113          ! avant        ALP(i) = alp_wake(i) + Alp_bl(i)
114          !             ALP(i) = alp_wake(i) + Alp_bl(i) + alp_offset ! modif sb
115          !         write(20,*)'ALE',ALE(i),Ale_bl(i),ale_wake(i)
116          !         write(21,*)'ALP',ALP(i),Alp_bl(i),alp_wake(i)
117          !         enddo
118
119          ! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
120          ! Modif FH 2010/04/27. Sans doute temporaire.
121          ! Deux options pour le alp_offset : constant si >?? 0 ou
122          ! proportionnel ??a w si <0
123          ! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
124          ! Estimation d'une vitesse verticale effective pour ALP
125          if (1==0) THEN
126             www(1:klon)=0.
127             do k=2,klev-1
128                do i=1,klon
129                   www(i)=max(www(i),-omega(i,k)*RD*t_seri(i,k) &
130                        /(RG*paprs(i,k)) *zw2(i,k)*zw2(i,k))
131                   ! if (paprs(i,k)>pbase(i)) then
132                   ! calcul approche de la vitesse verticale en m/s
133                   !  www(i)=max(www(i),-omega(i,k)*RD*temp(i,k)/(RG*paprs(i,k))
134                   !             endif
135                   !   Le 0.1 est en gros H / ps = 1e4 / 1e5
136                enddo
137             enddo
138             do i=1,klon
139                if (www(i)>0. .and. ale_bl(i)>0. ) www(i)=www(i)/ale_bl(i)
140             enddo
141          ENDIF
142
143
144          do i = 1,klon
145             ALE(i) = max(ale_wake(i),Ale_bl(i))
146             !cc nrlmd le 10/04/2012----------Stochastic triggering------------
147             if (iflag_trig_bl.ge.1) then
148                ALE(i) = max(ale_wake(i),Ale_bl_trig(i))
149             endif
150             !cc fin nrlmd le 10/04/2012
151             if (alp_offset>=0.) then
152                ALP(i) = alp_wake(i) + Alp_bl(i) + alp_offset ! modif sb
153             else
154                abort_message ='Ne pas passer la car www non calcule'
155                CALL abort_physic (modname,abort_message,1)
156
157                ! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
158                !                                _                  _
159                ! Ajout d'une composante 3 * A * w w'2 a w'3 avec
160                ! w=www : w max sous pbase ou A est la fraction
161                ! couverte par les ascendances w' on utilise le fait
162                ! que A * w'3 = ALP et donc A * w'2 ~ ALP / sqrt(ALE)
163                ! (on ajoute 0.1 pour les singularites)
164                ALP(i)=alp_wake(i)*(1.+3.*www(i)/( sqrt(ale_wake(i))+0.1) ) &
165                     +alp_bl(i)  *(1.+3.*www(i)/( sqrt(ale_bl(i))  +0.1) )
166                !    ALP(i)=alp_wake(i)+Alp_bl(i)+alp_offset*min(omega(i,6),0.)
167                !             if (alp(i)<0.) then
168                !                print*,'ALP ',alp(i),alp_wake(i) &
169                !                     ,Alp_bl(i),alp_offset*min(omega(i,6),0.)
170                !             endif
171             endif
172          enddo
173          ! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
174
175       endif
176       do i=1,klon
177          if (alp(i)>alp_max) then
178             IF(prt_level>9)WRITE(lunout,*)                             &
179                  'WARNING SUPER ALP (seuil=',alp_max, &
180                  '): i, alp, alp_wake,ale',i,alp(i),alp_wake(i),ale(i)
181             alp(i)=alp_max
182          endif
183          if (ale(i)>ale_max) then
184             IF(prt_level>9)WRITE(lunout,*)                             &
185                  'WARNING SUPER ALE (seuil=',ale_max, &
186                  '): i, alp, alp_wake,ale',i,ale(i),ale_wake(i),alp(i)
187             ale(i)=ale_max
188          endif
189       enddo
190
191       !fin calcul ale et alp
192       !=======================================================================
193
194
195  RETURN
196  END
197
Note: See TracBrowser for help on using the repository browser.