Changeset 272
- Timestamp:
- Aug 26, 2011, 10:40:39 AM (13 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/LMDZ.MARS/libf/phymars/surflayer_interpol.F
r268 r272 42 42 REAL, INTENT(IN) :: wmax(ngrid) 43 43 REAL, INTENT(IN) :: pts(ngrid),ph(ngrid,nlay) 44 REAL, INTENT(IN OUT) :: z_out ! output height (in m above surface)44 REAL, INTENT(IN) :: z_out ! output height (in m above surface) 45 45 REAL, INTENT(OUT) :: Teta_out(ngrid),u_out(ngrid)! interpolated fields at z_out : potential temperature and norm(uv) 46 46 REAL, INTENT(OUT) :: ustar(ngrid), tstar(ngrid) ! u* and teta* … … 57 57 ! Local(2): 58 58 ! --------- 59 REAL zout 59 60 60 61 REAL rib(ngrid) ! Bulk Richardson number … … 84 85 ! couche de surface: 85 86 ! ------------------ 86 87 zout=z_out 87 88 tstar(:)=0. 88 89 ustar(:)=0. … … 184 185 endif 185 186 186 reynolds(ig)=karman*sqrt(fm(ig))*sqrt(pu(ig,1)**2 + pv(ig,1)**2) 187 reynolds(ig)=karman*sqrt(fm(ig)) 188 & *sqrt(zu2(ig)) 189 ! & *sqrt(pu(ig,1)**2 + pv(ig,1)**2) 187 190 & *pz0(ig)/(log(z1z0)*nu) 188 191 pz0tcomp(ig)=pz0(ig)*exp(-karman*7.3* … … 244 247 245 248 DO ig=1,ngrid 246 IF(z _out .ge. pz(ig,1)) THEN247 z _out=1.249 IF(zout .ge. pz(ig,1)) THEN 250 zout=1. 248 251 print*, 'Warning, z_out should be less than the first 249 & 252 & vertical grid-point' 250 253 print*, 'z1 =',pz(ig,1) 251 254 print*, 'z_out=',z_out 252 255 print*, 'z_out has been changed to 1m 253 & 256 & and computation is resumed' 254 257 ENDIF 255 258 256 IF(z _out .lt. pz0(ig)) THEN257 z _out=1.258 print*, 'Warning, z_out should be more than the roughness259 & 260 print*, 'z0 =',pz0 (ig)259 IF(zout .lt. pz0tcomp(ig)) THEN 260 zout=pz0tcomp(ig) 261 print*, 'Warning, z_out should be more than the thermal 262 & roughness length' 263 print*, 'z0 =',pz0tcomp(ig) 261 264 print*, 'z_out=',z_out 262 print*, 'z_out has been changed to z0 263 & 265 print*, 'z_out has been changed to z0t 266 & and computation is resumed' 264 267 ENDIF 265 268 ENDDO 266 269 267 print*, 'interpolation of u and teta at z_out=',z _out270 print*, 'interpolation of u and teta at z_out=',zout 268 271 269 272 DO ig=1,ngrid 270 273 IF (L_mo(ig) .gt. 0.) THEN 271 u_out(ig)=(ustar(ig)/karman)*log(z _out/pz0(ig)) +272 & 5.*(ustar(ig)/(karman*L_mo(ig)))*(z _out-pz0(ig))274 u_out(ig)=(ustar(ig)/karman)*log(zout/pz0(ig)) + 275 & 5.*(ustar(ig)/(karman*L_mo(ig)))*(zout-pz0(ig)) 273 276 Teta_out(ig)=pts(ig)+(tstar(ig)/(prandtl(ig)*karman)) 274 & *log(z _out/pz0tcomp(ig)) +277 & *log(zout/pz0tcomp(ig)) + 275 278 & 5.*(tstar(ig)/(prandtl(ig)*karman*L_mo(ig))) 276 & *(z _out-pz0tcomp(ig))279 & *(zout-pz0tcomp(ig)) 277 280 ELSEIF (L_mo(ig) .lt. 0.) THEN 278 281 … … 280 283 281 284 u_out(ig)=(ustar(ig)/karman)*( 282 & 2.*atan((1.-16.*z _out/L_mo(ig))**0.25)285 & 2.*atan((1.-16.*zout/L_mo(ig))**0.25) 283 286 & -2.*atan((1.-16.*pz0(ig)/L_mo(ig))**0.25) 284 & -2.*log(1.+(1.-16.*z _out/L_mo(ig))**0.25)287 & -2.*log(1.+(1.-16.*zout/L_mo(ig))**0.25) 285 288 & +2.*log(1.+(1.-16.*pz0(ig)/L_mo(ig))**0.25) 286 & - log(1.+sqrt(1.-16.*z _out/L_mo(ig)))289 & - log(1.+sqrt(1.-16.*zout/L_mo(ig))) 287 290 & + log(1.+sqrt(1.-16.*pz0(ig)/L_mo(ig))) 288 & + log(z _out/pz0(ig))291 & + log(zout/pz0(ig)) 289 292 & ) 290 293 291 294 Teta_out(ig)=pts(ig)+(tstar(ig)/(prandtl(ig)*karman))*( 292 295 & 2.*log(1.+sqrt(1.-16.*pz0tcomp(ig)/L_mo(ig))) 293 & -2.*log(1.+sqrt(1.-16.*z _out/L_mo(ig)))294 & + log(z _out/pz0tcomp(ig))296 & -2.*log(1.+sqrt(1.-16.*zout/L_mo(ig))) 297 & + log(zout/pz0tcomp(ig)) 295 298 & ) 296 299 … … 305 308 306 309 u_out(ig)=(ustar(ig)/karman)*( 307 & (4./L_mo(ig))*(z _out-pz0(ig))308 & + (20./(L_mo(ig))**2)*(z _out**2-pz0(ig)**2)309 & + (160./(L_mo(ig))**3)*(z _out**3-pz0(ig)**3)310 & + log(z _out/pz0(ig))310 & (4./L_mo(ig))*(zout-pz0(ig)) 311 & + (20./(L_mo(ig))**2)*(zout**2-pz0(ig)**2) 312 & + (160./(L_mo(ig))**3)*(zout**3-pz0(ig)**3) 313 & + log(zout/pz0(ig)) 311 314 & ) 312 315 313 316 Teta_out(ig)=pts(ig)+(tstar(ig)/(prandtl(ig)*karman))*( 314 & (8./L_mo(ig))*(z _out-pz0tcomp(ig))315 & + (48./(L_mo(ig))**2)*(z _out**2-pz0tcomp(ig)**2)316 & + (1280./(3.*(L_mo(ig))**3))*(z _out**3-pz0tcomp(ig)**3)317 & + log(z _out/pz0tcomp(ig))317 & (8./L_mo(ig))*(zout-pz0tcomp(ig)) 318 & + (48./(L_mo(ig))**2)*(zout**2-pz0tcomp(ig)**2) 319 & + (1280./(3.*(L_mo(ig))**3))*(zout**3-pz0tcomp(ig)**3) 320 & + log(zout/pz0tcomp(ig)) 318 321 & ) 319 322 … … 327 330 ! Usefull diagnostics for the interpolation routine : 328 331 329 call WRITEDIAGFI(ngrid,'L',330 & 'Monin Obukhov length','m',331 & 2,L_mo)332 call WRITEDIAGFI(ngrid,'z0T',333 & 'thermal roughness length','m',334 & 2,pz0tcomp)335 call WRITEDIAGFI(ngrid,'z0',336 & 'roughness length','m',337 & 2,pz0)332 ! call WRITEDIAGFI(ngrid,'L', 333 ! & 'Monin Obukhov length','m', 334 ! & 2,L_mo) 335 ! call WRITEDIAGFI(ngrid,'z0T', 336 ! & 'thermal roughness length','m', 337 ! & 2,pz0tcomp) 338 ! call WRITEDIAGFI(ngrid,'z0', 339 ! & 'roughness length','m', 340 ! & 2,pz0) 338 341 339 342
Note: See TracChangeset
for help on using the changeset viewer.