UPP  V11.0.0
 All Data Structures Files Functions Pages
CALTAU.f
Go to the documentation of this file.
1 
27 
28  SUBROUTINE caltau(TAUX,TAUY)
29 
30 !
31 !
32  use vrbls3d, only: zint, pmid, q, t, uh, vh, el_pbl, zmid
33  use vrbls2d, only: z0, uz0, vz0
34  use masks, only: lmh
35  use params_mod, only: d00, d50, h1, d608, rd, d25
36  use ctlblk_mod, only: jsta_2l, jend_2u, lm, jsta, jend, spval, jsta_m,&
37  jm, im, jend_m, ista, iend, ista_m, iend_m, ista_2l, iend_2u
38  use gridspec_mod, only: gridtype
39 
40  implicit none
41 !
42 ! DECLARE VARIABLES.
43  INTEGER, dimension(4) :: kk(4)
44  INTEGER, dimension(jm) :: ive, ivw
45  REAL, dimension(ista:iend,jsta:jend), intent(inout) :: taux, tauy
46  REAL, ALLOCATABLE :: el(:,:,:)
47  REAL, dimension(ista:iend,jsta:jend) :: egridu,egridv,egrid4,egrid5, el0
48  REAL uz0v,vz0v
49  CHARACTER*1 agrid
50  integer i,j,lmhk,ie,iw,ii,jj
51  real dz,rdz,rsfc,tv,rho,ulmh,vlmh,deludz,delvdz,elsqr,zint1, &
52  zint2,z0v,psfc,tvv,qvv,elv,elv1,elv2
53 !
54 !********************************************************************
55 ! START CALTAU HERE.
56 !
57  ALLOCATE (el(ista_2l:iend_2u,jsta_2l:jend_2u,lm))
58 !
59 ! COMPUTE MASTER LENGTH SCALE.
60 !
61 ! CALL CLMAX(EL0,EGRIDU,EGRIDV,EGRID4,EGRID5)
62 ! CALL MIXLEN(EL0,EL)
63 !
64 ! INITIALIZE OUTPUT AND WORK ARRAY TO ZERO.
65 !
66  DO j=jsta,jend
67  DO i=ista,iend
68  egridu(i,j) = d00
69  egridv(i,j) = d00
70  taux(i,j) = spval
71  tauy(i,j) = spval
72  ENDDO
73  ENDDO
74 !
75 ! COMPUTE SURFACE LAYER U AND V WIND STRESSES.
76 !
77 ! ASSUME THAT U AND V HAVE UPDATED HALOS
78 !
79  IF(gridtype == 'A')THEN
80  CALL clmax(el0,egridu,egridv,egrid4,egrid5)
81  CALL mixlen(el0,el)
82 
83  DO j=jsta,jend
84  DO i=ista,iend
85 !
86  lmhk = nint(lmh(i,j))
87  IF(el(i,j,lmhk-1)<spval.and.z0(i,j)<spval.and. &
88  uz0(i,j)<spval.and.vz0(i,j)<spval)THEN
89 !
90 ! COMPUTE THICKNESS OF LAYER AT MASS POINT.
91 !
92  dz = d50*(zint(i,j,lmhk)-zint(i,j,lmhk+1))
93  dz = dz-z0(i,j)
94  rdz = 1./dz
95 !
96 ! COMPUTE REPRESENTATIVE AIR DENSITY.
97 !
98  psfc = pmid(i,j,lmhk)
99  tv = (h1+d608*q(i,j,lmhk))*t(i,j,lmhk)
100  rho = psfc/(rd*tv)
101 !
102 ! COMPUTE A MEAN MASS POINT WIND IN THE
103 ! FIRST ATMOSPHERIC ETA LAYER.
104 !
105  ulmh = uh(i,j,lmhk)
106  vlmh = vh(i,j,lmhk)
107 !
108 ! COMPUTE WIND SHEAR COMPONENTS ACROSS LAYER.
109 !
110  deludz = (ulmh-uz0(i,j))*rdz
111  delvdz = (vlmh-vz0(i,j))*rdz
112 !
113 ! COMPUTE U (EGRIDU) AND V (EGRIDV) WIND STRESSES.
114 !
115  elsqr = el(i,j,lmhk-1)*el(i,j,lmhk-1)
116  taux(i,j) = rho*elsqr*deludz*deludz
117  tauy(i,j) = rho*elsqr*delvdz*delvdz
118  ELSE
119  taux(i,j) = spval
120  tauy(i,j) = spval
121  ENDIF
122 
123 !
124  END DO
125  END DO
126  ELSE IF(gridtype == 'E')THEN
127  call exch(zint(1,jsta_2l,lm))
128  call exch(zint(1,jsta_2l,lm+1))
129  call exch(z0(1,jsta_2l))
130  call exch(pmid(1,jsta_2l,lm))
131  call exch(t(1,jsta_2l,lm))
132  call exch(q(1,jsta_2l,lm))
133  call exch(el_pbl(1,jsta_2l,lm))
134  call exch(el_pbl(1,jsta_2l,lm-1))
135 
136  DO j=jsta_m,jend_m
137  ive(j)=mod(j,2)
138  ivw(j)=ive(j)-1
139  ENDDO
140 
141  DO j=jsta_m,jend_m
142  DO i=ista_m,iend_m
143 !
144  lmhk = nint(lmh(i,j))
145  ie=i+ive(j)
146  iw=i+ivw(j)
147  zint1=(zint(iw,j,lmhk)+zint(ie,j,lmhk) &
148  +zint(i,j+1,lmhk)+zint(i,j-1,lmhk))*d25
149  zint2=(zint(iw,j,lmhk+1)+zint(ie,j,lmhk+1) &
150  +zint(i,j+1,lmhk+1)+zint(i,j-1,lmhk+1))*d25
151  dz = d50*(zint1-zint2)
152  z0v=(z0(iw,j)+z0(ie,j)+z0(i,j+1)+z0(i,j-1))*d25
153  dz = dz-z0v
154  rdz = 1./dz
155 !
156 ! COMPUTE REPRESENTATIVE AIR DENSITY.
157 !
158  psfc = (pmid(iw,j,lmhk)+pmid(ie,j,lmhk) &
159  +pmid(i,j+1,lmhk)+pmid(i,j-1,lmhk))*d25
160  tvv = (t(iw,j,lmhk)+t(ie,j,lmhk) &
161  +t(i,j+1,lmhk)+t(i,j-1,lmhk))*d25
162  qvv = (q(iw,j,lmhk)+q(ie,j,lmhk) &
163  +q(i,j+1,lmhk)+q(i,j-1,lmhk))*d25
164  tv = (h1+d608*qvv)*tvv
165  rho = psfc/(rd*tv)
166 
167 ! COMPUTE WIND SHEAR COMPONENTS ACROSS LAYER.
168 !
169  deludz = (uh(i,j,lmhk)-uz0(i,j))*rdz
170  delvdz = (vh(i,j,lmhk)-vz0(i,j))*rdz
171 
172 ! COMPUTE U (EGRIDU) AND V (EGRIDV) WIND STRESSES.
173 !
174  elv1=(el_pbl(iw,j,lmhk)+el_pbl(ie,j,lmhk) &
175  +el_pbl(i,j+1,lmhk)+el_pbl(i,j-1,lmhk))*d25
176  elv2=(el_pbl(iw,j,lmhk-1)+el_pbl(ie,j,lmhk-1) &
177  +el_pbl(i,j+1,lmhk-1)+el_pbl(i,j-1,lmhk-1))*d25
178  elv=(elv1+elv2)/2.0 ! EL is defined at the bottom of layer
179  elsqr =elv*elv
180  taux(i,j)=rho*elsqr*deludz*deludz
181  tauy(i,j)=rho*elsqr*delvdz*delvdz
182 ! ii=im/2
183 ! jj=(jsta+jend)/2
184 ! if(i==ii.and.j==jj)print*,'sample tau'
185 ! & ,RHO,ELSQR,DELUDZ,DELVDZ
186  END DO
187  END DO
188  ELSE IF(gridtype == 'B')THEN
189 ! PUT TAUX AND TAUY ON MASS POINTS
190  call exch(vh(1,jsta_2l,lm))
191  DO j=jsta_m,jend_m
192  DO i=ista_m,iend_m
193 !
194  lmhk = nint(lmh(i,j))
195 !
196 ! COMPUTE THICKNESS OF LAYER AT MASS POINT.
197 !
198 ! DZ = D50*(ZINT(I,J,LMHK)-ZINT(I,J,LMHK+1))
199 ! DZ = ZMID(I,J,LMHK)-Z0(I,J)
200  dz=zmid(i,j,lmhk)-(z0(i,j)+zint(i,j,lmhk+1))
201  if(dz==0.0)dz=0.2
202  rdz = 1./dz
203 !
204 ! COMPUTE REPRESENTATIVE AIR DENSITY.
205 !
206  psfc = pmid(i,j,lmhk)
207  tv = (h1+d608*q(i,j,lmhk))*t(i,j,lmhk)
208  rho = psfc/(rd*tv)
209 !
210 ! PUT U AND V ONTO MASS POINTS
211 !
212  ulmh = 0.5*(uh(i-1,j,lmhk)+uh(i,j,lmhk))
213  vlmh = 0.5*(vh(i,j-1,lmhk)+vh(i,j,lmhk))
214 !
215 ! COMPUTE WIND SHEAR COMPONENTS ACROSS LAYER.
216 !
217  deludz = (ulmh-uz0(i,j))*rdz
218  delvdz = (vlmh-vz0(i,j))*rdz
219 !
220 ! COMPUTE U (EGRIDU) AND V (EGRIDV) WIND STRESSES.
221 !
222  elv=0.5*(el_pbl(i,j,lmhk)+el_pbl(i,j,lmhk-1))
223  elsqr = elv*elv
224  taux(i,j) = rho*elsqr*deludz*deludz
225 ! if(TAUX(I,J)>1.0e2)print*,'Debug TAUX= ',i,j, &
226 ! ELV,ULMH,UZ0(I,J),ZMID(I,J,LMHK),Z0(I,J),RDZ,TAUX(I,J),zint(i,j,lm+1)
227  tauy(i,j) = rho*elsqr*delvdz*delvdz
228 
229  END DO
230  END DO
231  END IF
232 !
233  DEALLOCATE(el)
234 ! END OF ROUTINE.
235 !
236  RETURN
237  END
masks
Definition: MASKS_mod.f:1
params_mod
Definition: params.F:1