source: trunk/LMDZ.MARS/libf/phymars/gwprofil.F @ 308

Last change on this file since 308 was 38, checked in by emillour, 14 years ago

Ajout du modè Martien (mon LMDZ.MARS.BETA, du 28/01/2011) dans le rértoire mars, pour pouvoir suivre plus facilement les modifs.
EM

File size: 5.6 KB
Line 
1      SUBROUTINE GWPROFIL
2     *         ( klon, klev
3     *         , kgwd ,kdx  , ktest
4     *         , KKCRIT, KKCRITH, KCRIT ,  kkenvh, kknu,kknu2
5     *         , PAPHM1, PRHO   , PSTAB , PTFR , PVPH , PRI , PTAU
6     *         , ptauf ,pdmod   , pnu   , psig ,pgamma, pvar      )
7
8C**** *GWPROFIL*
9C
10C     PURPOSE.
11C     --------
12C
13C**   INTERFACE.
14C     ----------
15C          FROM *GWDRAG*
16C
17C        EXPLICIT ARGUMENTS :
18C        --------------------
19C     ==== INPUTS ===
20C     ==== OUTPUTS ===
21C
22C        IMPLICIT ARGUMENTS :   NONE
23C        --------------------
24C
25C     METHOD:
26C     -------
27C     THE STRESS PROFILE FOR GRAVITY WAVES IS COMPUTED AS FOLLOWS:
28C     IT IS CONSTANT (NO GWD) AT THE LEVELS BETWEEN THE GROUND
29C     AND THE TOP OF THE BLOCKED LAYER (KKENVH).
30C     IT DECREASES LINEARLY WITH HEIGHTS FROM THE TOP OF THE
31C     BLOCKED LAYER TO 3*VAROR (kKNU), TO SIMULATES LEE WAVES OR
32C     NONLINEAR GRAVITY WAVE BREAKING.
33C     ABOVE IT IS CONSTANT, EXCEPT WHEN THE WAVE ENCOUNTERS A CRITICAL
34C     LEVEL (KCRIT) OR WHEN IT BREAKS.
35C     
36C
37C
38C     EXTERNALS.
39C     ----------
40C
41C
42C     REFERENCE.
43C     ----------
44C
45C        SEE ECMWF RESEARCH DEPARTMENT DOCUMENTATION OF THE "I.F.S."
46C
47C     AUTHOR.
48C     -------
49C
50C     MODIFICATIONS.
51C     --------------
52C     PASSAGE OF THE NEW GWDRAG TO I.F.S. (F. LOTT, 22/11/93)
53C-----------------------------------------------------------------------
54      implicit none
55C
56
57C
58
59#include "dimensions.h"
60#include "dimphys.h"
61#include "dimradmars.h"
62      integer klon,klev,kidia,kfdia
63#include "yoegwd.h"
64
65C-----------------------------------------------------------------------
66C
67C*       0.1   ARGUMENTS
68C              ---------
69C
70      integer kgwd
71      INTEGER KKCRIT(NDLO2),KKCRITH(NDLO2),KCRIT(NDLO2)
72     *       ,kdx(NDLO2),ktest(NDLO2)
73     *       ,kkenvh(NDLO2),kknu(NDLO2),kknu2(NDLO2)
74C
75      REAL PAPHM1(NDLO2,klev+1), PSTAB(NDLO2,klev+1),
76     *     PRHO  (NDLO2,klev+1), PVPH (NDLO2,klev+1),
77     *     PRI   (NDLO2,klev+1), PTFR (NDLO2), PTAU(NDLO2,klev+1),
78     *     ptauf (NDLO2,klev+1)
79     
80      REAL pdmod (NDLO2) , pnu (NDLO2) , psig(NDLO2),
81     *     pgamma(NDLO2) , pvar(NDLO2)
82     
83C-----------------------------------------------------------------------
84C
85C*       0.2   LOCAL ARRAYS
86C              ------------
87C
88c   declarations pour 'implicit none"
89      real zsqr,zalfa,zriw,zalpha,zb,zdel,zdz2n,zdelp,zdelpt
90
91      integer ji,jk,jl,ilevh
92      REAL ZDZ2 (NDLO2,nlayermx) , ZNORM(NDLO2) , zoro(NDLO2)
93      REAL ZTAU (NDLO2,nlayermx+1)
94C
95C-----------------------------------------------------------------------
96C
97C*         1.    INITIALIZATION
98C                --------------
99
100
101      kidia=1
102      kfdia=klon
103
104 100  CONTINUE
105C
106C
107C*    COMPUTATIONAL CONSTANTS.
108C     ------------- ----------
109C
110      ilevh=KLEV/3
111C
112      DO 400 ji=1,kgwd
113      jl=kdx(ji)
114      Zoro(JL)=Psig(JL)*Pdmod(JL)/4./max(pvar(jl),1.0)
115      ZTAU(JL,KKNU(JL)+1)=PTAU(JL,KKNU(JL)+1)
116      ZTAU(JL,KLEV+1)=PTAU(JL,KLEV+1)
117  400 CONTINUE
118C
119      DO 430 JK=KLEV,2,-1
120C
121C
122C*         4.1    CONSTANT WAVE STRESS UNTIL TOP OF THE
123C                 BLOCKING LAYER.
124  410 CONTINUE
125C
126      DO 411 ji=1,kgwd
127      jl=kdx(ji)
128           IF(JK.GE.KKNU2(JL)) THEN
129           PTAU(JL,JK)=ZTAU(JL,KLEV+1)
130           ENDIF
131 411  CONTINUE             
132C
133C*         4.15   CONSTANT SHEAR STRESS UNTIL THE TOP OF THE
134C                 LOW LEVEL FLOW LAYER.
135 415  CONTINUE
136C       
137C
138C*         4.2    WAVE DISPLACEMENT AT NEXT LEVEL.
139C
140  420 CONTINUE
141C
142      DO 421 ji=1,kgwd
143      jl=kdx(ji)
144      IF(JK.LT.KKNU2(JL)) THEN
145      ZNORM(JL)=gkdrag*PRHO(JL,JK)*SQRT(PSTAB(JL,JK))*PVPH(JL,JK)
146     *                                                    *zoro(jl)
147      ZDZ2(JL,JK)=PTAU(JL,JK+1)/max(ZNORM(JL),gssec)
148      ENDIF
149  421 CONTINUE
150C
151C*         4.3    WAVE RICHARDSON NUMBER, NEW WAVE DISPLACEMENT
152C*                AND STRESS:  BREAKING EVALUATION AND CRITICAL
153C                 LEVEL
154C
155      DO 431 ji=1,kgwd
156      jl=kdx(ji)
157          IF(JK.LT.KKNU2(JL)) THEN
158          IF((PTAU(JL,JK+1).LT.GTSEC).OR.(JK.LE.KCRIT(JL))) THEN
159            PTAU(JL,JK)=0.0
160          ELSE
161               ZSQR=SQRT(PRI(JL,JK))
162               ZALFA=SQRT(PSTAB(JL,JK)*ZDZ2(JL,JK))/PVPH(JL,JK)
163               ZRIW=PRI(JL,JK)*(1.-ZALFA)/(1+ZALFA*ZSQR)**2
164               IF(ZRIW.LT.GRCRIT) THEN
165                 ZDEL=4./ZSQR/GRCRIT+1./GRCRIT**2+4./GRCRIT
166                 ZB=1./GRCRIT+2./ZSQR
167                 ZALPHA=0.5*(-ZB+SQRT(ZDEL))
168                 ZDZ2N=(PVPH(JL,JK)*ZALPHA)**2/PSTAB(JL,JK)
169                 PTAU(JL,JK)=ZNORM(JL)*ZDZ2N
170               ELSE
171                 PTAU(JL,JK)=ZNORM(JL)*ZDZ2(JL,JK)
172               ENDIF
173            PTAU(JL,JK)=MIN(PTAU(JL,JK),PTAU(JL,JK+1))
174          ENDIF
175          ENDIF
176  431 CONTINUE
177 
178  430 CONTINUE
179  440 CONTINUE
180 
181c     write(*,*) 'ptau'
182c     write(*,99) ((ji,ilevh,ptau(ji,ilevh),ji=1,NDLO2),
183c    .                  ilevh=1,nlayermx+1)
184 99   FORMAT(i3,i3,f15.5)
185
186
187C  REORGANISATION OF THE STRESS PROFILE
188C  IF BREAKING OCCURS AT LOW LEVEL:
189
190      DO 530 ji=1,kgwd
191      jl=kdx(ji)
192      ZTAU(JL,KKENVH(JL))=PTAU(JL,KKENVH(JL))
193      ZTAU(JL,KKCRITH(JL))=PTAU(JL,KKCRITH(JL))
194 530  CONTINUE     
195
196      DO 531 JK=1,KLEV
197     
198      DO 532 ji=1,kgwd
199      jl=kdx(ji)
200               
201         IF(JK.GT.KKCRITH(JL).AND.JK.LT.KKENVH(JL))THEN
202
203          ZDELP=PAPHM1(JL,JK)-PAPHM1(JL,KKENVH(JL))
204          ZDELPT=PAPHM1(JL,KKCRITH(JL))-PAPHM1(JL,KKENVH(JL))
205          PTAU(JL,JK)=ZTAU(JL,KKENVH(JL)) +
206     .                (ZTAU(JL,KKCRITH(JL))-ZTAU(JL,KKENVH(JL)) )*
207     .                ZDELP/ZDELPT
208     
209        ENDIF
210           
211 532  CONTINUE   
212 
213 531  CONTINUE       
214
215      RETURN
216      END
Note: See TracBrowser for help on using the repository browser.