1 | MODULE gwprofil_mod |
---|
2 | |
---|
3 | IMPLICIT NONE |
---|
4 | |
---|
5 | CONTAINS |
---|
6 | |
---|
7 | SUBROUTINE GWPROFIL( ngrid, nlayer, kgwd ,kdx, ktest, & |
---|
8 | IKCRITH, ICRIT, IKENVH, IKNU,IKNU2,pplev, & |
---|
9 | ZRHO, BV, ZVPH, PRI,zdmod, psig , pvar, & |
---|
10 | !not used variables |
---|
11 | !IKCRIT,ZTAUf,ZTFR,znu,pgam, |
---|
12 | ! in-output |
---|
13 | ZTAU ) |
---|
14 | |
---|
15 | !--------------------------------------------------------------------------- |
---|
16 | ! THe stress profile for gravity waves is computed as follows: |
---|
17 | ! It is constant (no GWD) at the levels between the ground and |
---|
18 | ! the top of the blocked layer (IKENVH). It decreases linearly |
---|
19 | ! with heights from the top of the blocked layer to 3*varor(IKNU) |
---|
20 | ! to simulates lee waves or nonlinear gravity wave breaking. |
---|
21 | ! Above it is constant, except when the wave encounters a critical |
---|
22 | ! level(ICRIT) or when it breaks.! |
---|
23 | ! PASSAGE OF THE NEW GWDRAG TO I.F.S. (F. LOTT, 22/11/93) |
---|
24 | ! REFERENCE. |
---|
25 | ! SEE ECMWF RESEARCH DEPARTMENT DOCUMENTATION OF THE "I.F.S." |
---|
26 | ! MODIFICATIONS. |
---|
27 | ! Rewrited by J.Liu 03/03/2022 |
---|
28 | !--------------------------------------------------------------------------- |
---|
29 | |
---|
30 | use dimradmars_mod, only: ndomainsz |
---|
31 | use yoegwd_h, only: gkdrag, grcrit, gssec, gtsec |
---|
32 | |
---|
33 | implicit none |
---|
34 | |
---|
35 | ! 0. DECLARATIONS: |
---|
36 | ! 0.1 ARGUMENTS |
---|
37 | integer, intent(in):: ngrid ! Number of atmospheric columns |
---|
38 | integer, intent(in):: nlayer ! Number of atmospheric layers |
---|
39 | INTEGER, intent(in) :: kgwd ! Number of points at which the scheme is called |
---|
40 | ! INTEGER, INTENT(IN) :: IKCRIT(ndomainsz ) ! Top of low level flow height (not used) |
---|
41 | INTEGER, INTENT(IN) :: IKCRITH(ndomainsz ) ! dynamical mixing height for the breaking of gravity waves |
---|
42 | INTEGER, INTENT(IN) :: ICRIT(ndomainsz ) ! Critical layer where orographic GW breaks |
---|
43 | INTEGER, INTENT(IN) :: kdx(ndomainsz ) ! Points at which to call the scheme |
---|
44 | INTEGER, INTENT(IN) :: ktest(ndomainsz ) ! Map of calling points |
---|
45 | INTEGER, INTENT(IN) :: IKENVH(ndomainsz ) ! Top of the blocked layer |
---|
46 | INTEGER, INTENT(IN) :: IKNU(ndomainsz ) ! 4*pvar layer |
---|
47 | INTEGER, INTENT(IN) :: IKNU2(ndomainsz ) ! 3*pvar layer |
---|
48 | |
---|
49 | REAL, INTENT(IN):: pplev(ndomainsz ,nlayer+1) ! Pressure at 1/2 levels(Pa) |
---|
50 | REAL, INTENT(IN):: BV(ndomainsz ,nlayer+1) ! Brunt–Väisälä frequency at 1/2 level (output from OROSETUP input for GWSTRESS and GWPROFIL) |
---|
51 | REAL, INTENT(IN):: ZRHO (ndomainsz ,nlayer+1) ! Atmospheric density at 1/2 level |
---|
52 | REAL, INTENT(IN):: ZVPH (ndomainsz ,nlayer+1) ! SUB-GRID SCALE STANDARD DEVIATION at 1/2 level |
---|
53 | REAL, INTENT(IN):: PRI (ndomainsz ,nlayer+1) ! Mean flow richardson number at 1/2 level |
---|
54 | ! REAL, INTENT(IN):: ZTFR (ndomainsz ) ! not used |
---|
55 | ! REAL, INTENT(IN):: ZTAUf (ndomainsz ,nlayer+1) ! not used |
---|
56 | ! |
---|
57 | REAL, INTENT(IN):: zdmod (ndomainsz ) |
---|
58 | ! REAL, INTENT(IN):: znu (ndomainsz ) ! not used |
---|
59 | REAL, INTENT(IN):: psig(ndomainsz ) ! Sub-grid scale slope |
---|
60 | ! REAL, INTENT(IN):: pgam(ndomainsz ) ! Sub-grid scale anisotropy(not used) |
---|
61 | REAL, INTENT(IN):: pvar(ndomainsz ) ! Sub-grid scale standard deviation |
---|
62 | REAL, INTENT(inout):: ZTAU(ndomainsz ,nlayer+1) ! Gravity waves induced stress |
---|
63 | |
---|
64 | ! 0.2 LOCAL ARRAYS |
---|
65 | real zsqr,zalfa,zriw,zalpha,zb,zdel,zdz2n,zdelp,zdelpt |
---|
66 | integer ji,jk,jl,ilevh |
---|
67 | REAL ZDZ2 (ndomainsz ,nlayer) , ZNORM(ndomainsz ) , zoro(ndomainsz ) |
---|
68 | REAL Z_TAU (ndomainsz ,nlayer+1) |
---|
69 | |
---|
70 | !----------------------------------------------------------------------- |
---|
71 | !* 1. Initialization |
---|
72 | !----------------------------------------------------------------------- |
---|
73 | ! kidia=1 |
---|
74 | ! kfdia=ngrid |
---|
75 | 100 CONTINUE |
---|
76 | ! 1.1 Computional constants . |
---|
77 | ilevh=nlayer/3 |
---|
78 | |
---|
79 | DO ji=1,kgwd |
---|
80 | jl=kdx(ji) |
---|
81 | Zoro(JL)=psig(JL)*zdmod(JL)/4./max(pvar(jl),1.0) |
---|
82 | Z_TAU(JL,IKNU(JL)+1)=ZTAU(JL,IKNU(JL)+1) |
---|
83 | Z_TAU(JL,nlayer+1)=ZTAU(JL,nlayer+1) |
---|
84 | ENDDO |
---|
85 | |
---|
86 | ! all start here with this long loop |
---|
87 | DO JK=nlayer,2,-1 |
---|
88 | ! 4.1 Constant wave stress until top of the blocking layer |
---|
89 | 410 CONTINUE |
---|
90 | DO ji=1,kgwd |
---|
91 | jl=kdx(ji) |
---|
92 | IF(JK.GE.IKNU2(JL)) THEN |
---|
93 | ZTAU(JL,JK)=Z_TAU(JL,nlayer+1) |
---|
94 | ENDIF |
---|
95 | ENDDO |
---|
96 | |
---|
97 | !* 4.15 Constant shear stress until the top of the low level flow layer |
---|
98 | 415 CONTINUE |
---|
99 | ! 4.2 Wave displacement at next level. |
---|
100 | 420 CONTINUE |
---|
101 | |
---|
102 | DO ji=1,kgwd |
---|
103 | jl=kdx(ji) |
---|
104 | IF(JK.LT.IKNU2(JL)) THEN |
---|
105 | ZNORM(JL)=gkdrag*ZRHO(JL,JK)*SQRT(BV(JL,JK))*ZVPH(JL,JK)*zoro(jl) |
---|
106 | ZDZ2(JL,JK)=ZTAU(JL,JK+1)/max(ZNORM(JL),gssec) |
---|
107 | ENDIF |
---|
108 | ENDDO |
---|
109 | |
---|
110 | ! 4.3 Wave Richardson number, new wave displacement and stress: |
---|
111 | ! Breaking evaluation and critical level C |
---|
112 | DO ji=1,kgwd |
---|
113 | jl=kdx(ji) |
---|
114 | IF(JK.LT.IKNU2(JL)) THEN |
---|
115 | IF((ZTAU(JL,JK+1).LT.GTSEC).OR.(JK.LE.ICRIT(JL))) THEN |
---|
116 | ZTAU(JL,JK)=0.0 |
---|
117 | ELSE |
---|
118 | ZSQR=SQRT(PRI(JL,JK)) |
---|
119 | ZALFA=SQRT(BV(JL,JK)*ZDZ2(JL,JK))/ZVPH(JL,JK) |
---|
120 | ZRIW=PRI(JL,JK)*(1.-ZALFA)/(1+ZALFA*ZSQR)**2 |
---|
121 | IF(ZRIW.LT.GRCRIT) THEN |
---|
122 | ZDEL=4./ZSQR/GRCRIT+1./GRCRIT**2+4./GRCRIT |
---|
123 | ZB=1./GRCRIT+2./ZSQR |
---|
124 | ZALPHA=0.5*(-ZB+SQRT(ZDEL)) |
---|
125 | ZDZ2N=(ZVPH(JL,JK)*ZALPHA)**2/BV(JL,JK) |
---|
126 | ZTAU(JL,JK)=ZNORM(JL)*ZDZ2N |
---|
127 | ELSE |
---|
128 | ZTAU(JL,JK)=ZNORM(JL)*ZDZ2(JL,JK) |
---|
129 | ENDIF |
---|
130 | ZTAU(JL,JK)=MIN(ZTAU(JL,JK),ZTAU(JL,JK+1)) |
---|
131 | ENDIF |
---|
132 | ENDIF |
---|
133 | ENDDO |
---|
134 | end DO !JK=nlayer,2,-1 |
---|
135 | !------------------------------------------------------------------------------------------- |
---|
136 | |
---|
137 | 440 CONTINUE |
---|
138 | |
---|
139 | ! write(*,*) 'ZTAU' |
---|
140 | ! write(*,99) ((ji,ilevh,ZTAU(ji,ilevh),ji=1,ndomainsz ), |
---|
141 | ! . ilevh=1,nlayer+1) |
---|
142 | 99 FORMAT(i3,i3,f15.5) ! not used |
---|
143 | |
---|
144 | ! Proganisation of the stress profile if breaking occurs at low level |
---|
145 | DO ji=1,kgwd |
---|
146 | jl=kdx(ji) |
---|
147 | Z_TAU(JL,IKENVH(JL))=ZTAU(JL,IKENVH(JL)) |
---|
148 | Z_TAU(JL,IKCRITH(JL))=ZTAU(JL,IKCRITH(JL)) |
---|
149 | ENDDO |
---|
150 | |
---|
151 | DO JK=1,nlayer |
---|
152 | DO ji=1,kgwd |
---|
153 | jl=kdx(ji) |
---|
154 | IF(JK.GT.IKCRITH(JL).AND.JK.LT.IKENVH(JL))THEN |
---|
155 | ZDELP=pplev(JL,JK)-pplev(JL,IKENVH(JL)) |
---|
156 | ZDELPT=pplev(JL,IKCRITH(JL))-pplev(JL,IKENVH(JL)) |
---|
157 | ZTAU(JL,JK)=Z_TAU(JL,IKENVH(JL))+(Z_TAU(JL,IKCRITH(JL))-Z_TAU(JL,IKENVH(JL)))*ZDELP/ZDELPT |
---|
158 | ENDIF |
---|
159 | ENDDO |
---|
160 | end DO |
---|
161 | |
---|
162 | END SUBROUTINE GWPROFIL |
---|
163 | |
---|
164 | END MODULE gwprofil_mod |
---|