Changeset 667 for palm/trunk/SOURCE/init_grid.f90
 Timestamp:
 Dec 23, 2010 12:06:00 PM (12 years ago)
 Location:
 palm/trunk/SOURCE
 Files:

 2 edited
Legend:
 Unmodified
 Added
 Removed

palm/trunk/SOURCE

Property
svn:mergeinfo
set to
(toggle deleted branches)
/palm/branches/suehring 423666 /palm/branches/letzel/masked_output/SOURCE 296409

Property
svn:mergeinfo
set to
(toggle deleted branches)

palm/trunk/SOURCE/init_grid.f90
r559 r667 4 4 ! Current revisions: 5 5 !  6 ! 6 ! Definition of new array bounds nxlg, nxrg, nysg, nyng on each PE. 7 ! Furthermore the allocation of arrays and steering of loops is done with these 8 ! parameters. Call of exchange_horiz are modified. 9 ! In case of dirichlet bounday condition at the bottom zu(0)=0.0 10 ! dzu_mg has to be set explicitly for a equally spaced grid near bottom. 11 ! ddzu_pres added to use a equally spaced grid near bottom. 7 12 ! 8 13 ! Former revisions: … … 76 81 77 82 REAL, DIMENSION(:,:,:), ALLOCATABLE :: distance 78 83 84 ! 85 ! Computation of the array bounds. 86 nxlg = nxl  nbgp 87 nxrg = nxr + nbgp 88 nysg = nys  nbgp 89 nyng = nyn + nbgp 79 90 ! 80 91 ! Allocate grid arrays 81 92 ALLOCATE( ddzu(1:nzt+1), ddzw(1:nzt+1), dd2zu(1:nzt), dzu(1:nzt+1), & 82 dzw(1:nzt+1), l_grid(1:nzt), zu( 0:nzt+1), zw(0:nzt+1) )93 dzw(1:nzt+1), l_grid(1:nzt), zu(nzb:nzt+1), zw(nzb:nzt+1) ) 83 94 84 95 ! … … 97 108 ! 98 109 ! Grid for atmosphere with surface at z=0 (k=0, wgrid). 99 ! Since the wlevel lies on the surface, the first ulevel (staggered!)100 ! lies below the surface (used for "mirror" boundary condition).101 110 ! The first ulevel above the surface corresponds to the top of the 102 111 ! Prandtllayer. 103 zu(0) =  dz * 0.5 112 113 IF ( ibc_uv_b == 0 .OR. ibc_uv_b == 2 ) THEN 114 zu(0) = 0.0 115 ! zu(0) =  dz * 0.5 116 ELSE 117 zu(0) =  dz * 0.5 118 ENDIF 104 119 zu(1) = dz * 0.5 105 120 … … 173 188 dd2zu(k) = 1.0 / ( dzu(k) + dzu(k+1) ) 174 189 ENDDO 190 191 ! 192 ! In case of FFT method or SOR swap inverse grid lenght ddzu to ddzu_fft and 193 ! modify the lowest entry. This is necessary for atmosphere runs, because 194 ! zu(0) and so ddzu(1) changed. Accompanied with this modified arrays 195 ! pressure solver uses wrong values and this causes kinks in the profiles 196 ! of turbulent quantities. 197 IF ( psolver /= 'multigrid' ) THEN 198 ALLOCATE( ddzu_pres(1:nzt+1) ) 199 ddzu_pres = ddzu 200 IF( .NOT. ocean ) ddzu_pres(1) = ddzu_pres(2) 201 ENDIF 175 202 176 203 ! … … 187 214 188 215 dzu_mg(:,maximum_grid_level) = dzu 216 ! 217 ! To ensure a equally spaced grid. For ocean runs this is not necessary, 218 ! because zu(0) does not changed so far. Also this would cause errors 219 ! if a vertical stretching for ocean runs is used. 220 IF ( .NOT. ocean ) dzu_mg(1,maximum_grid_level) = dzu(2) 189 221 dzw_mg(:,maximum_grid_level) = dzw 190 222 nzt_l = nzt … … 239 271 ! the flag arrays needed for the multigrid method 240 272 gls = 2**( maximum_grid_level ) 273 241 274 ALLOCATE( corner_nl(nys:nyn,nxl:nxr), corner_nr(nys:nyn,nxl:nxr), & 242 275 corner_sl(nys:nyn,nxl:nxr), corner_sr(nys:nyn,nxl:nxr), & 243 nzb_local(gls:ny+gls,gls:nx+gls), nzb_tmp(1:ny+1,1:nx+1), & 276 nzb_local(gls:ny+gls,gls:nx+gls), & 277 nzb_tmp(nbgp:ny+nbgp,nbgp:nx+nbgp), & 244 278 wall_l(nys:nyn,nxl:nxr), wall_n(nys:nyn,nxl:nxr), & 245 279 wall_r(nys:nyn,nxl:nxr), wall_s(nys:nyn,nxl:nxr) ) 246 ALLOCATE( fwxm(nys1:nyn+1,nxl1:nxr+1), fwxp(nys1:nyn+1,nxl1:nxr+1), & 247 fwym(nys1:nyn+1,nxl1:nxr+1), fwyp(nys1:nyn+1,nxl1:nxr+1), & 248 fxm(nys1:nyn+1,nxl1:nxr+1), fxp(nys1:nyn+1,nxl1:nxr+1), & 249 fym(nys1:nyn+1,nxl1:nxr+1), fyp(nys1:nyn+1,nxl1:nxr+1), & 250 nzb_s_inner(nys1:nyn+1,nxl1:nxr+1), & 251 nzb_s_outer(nys1:nyn+1,nxl1:nxr+1), & 252 nzb_u_inner(nys1:nyn+1,nxl1:nxr+1), & 253 nzb_u_outer(nys1:nyn+1,nxl1:nxr+1), & 254 nzb_v_inner(nys1:nyn+1,nxl1:nxr+1), & 255 nzb_v_outer(nys1:nyn+1,nxl1:nxr+1), & 256 nzb_w_inner(nys1:nyn+1,nxl1:nxr+1), & 257 nzb_w_outer(nys1:nyn+1,nxl1:nxr+1), & 258 nzb_diff_s_inner(nys1:nyn+1,nxl1:nxr+1), & 259 nzb_diff_s_outer(nys1:nyn+1,nxl1:nxr+1), & 260 nzb_diff_u(nys1:nyn+1,nxl1:nxr+1), & 261 nzb_diff_v(nys1:nyn+1,nxl1:nxr+1), & 262 nzb_2d(nys1:nyn+1,nxl1:nxr+1), & 263 wall_e_x(nys1:nyn+1,nxl1:nxr+1), & 264 wall_e_y(nys1:nyn+1,nxl1:nxr+1), & 265 wall_u(nys1:nyn+1,nxl1:nxr+1), & 266 wall_v(nys1:nyn+1,nxl1:nxr+1), & 267 wall_w_x(nys1:nyn+1,nxl1:nxr+1), & 268 wall_w_y(nys1:nyn+1,nxl1:nxr+1) ) 269 270 ALLOCATE( l_wall(nzb:nzt+1,nys1:nyn+1,nxl1:nxr+1) ) 280 ALLOCATE( fwxm(nysg:nyng,nxlg:nxrg), fwxp(nysg:nyng,nxlg:nxrg), & 281 fwym(nysg:nyng,nxlg:nxrg), fwyp(nysg:nyng,nxlg:nxrg), & 282 fxm(nysg:nyng,nxlg:nxrg), fxp(nysg:nyng,nxlg:nxrg), & 283 fym(nysg:nyng,nxlg:nxrg), fyp(nysg:nyng,nxlg:nxrg), & 284 nzb_s_inner(nysg:nyng,nxlg:nxrg), & 285 nzb_s_outer(nysg:nyng,nxlg:nxrg), & 286 nzb_u_inner(nysg:nyng,nxlg:nxrg), & 287 nzb_u_outer(nysg:nyng,nxlg:nxrg), & 288 nzb_v_inner(nysg:nyng,nxlg:nxrg), & 289 nzb_v_outer(nysg:nyng,nxlg:nxrg), & 290 nzb_w_inner(nysg:nyng,nxlg:nxrg), & 291 nzb_w_outer(nysg:nyng,nxlg:nxrg), & 292 nzb_diff_s_inner(nysg:nyng,nxlg:nxrg), & 293 nzb_diff_s_outer(nysg:nyng,nxlg:nxrg), & 294 nzb_diff_u(nysg:nyng,nxlg:nxrg), & 295 nzb_diff_v(nysg:nyng,nxlg:nxrg), & 296 nzb_2d(nysg:nyng,nxlg:nxrg), & 297 wall_e_x(nysg:nyng,nxlg:nxrg), & 298 wall_e_y(nysg:nyng,nxlg:nxrg), & 299 wall_u(nysg:nyng,nxlg:nxrg), & 300 wall_v(nysg:nyng,nxlg:nxrg), & 301 wall_w_x(nysg:nyng,nxlg:nxrg), & 302 wall_w_y(nysg:nyng,nxlg:nxrg) ) 303 304 305 306 ALLOCATE( l_wall(nzb:nzt+1,nysg:nyng,nxlg:nxrg) ) 271 307 272 308 nzb_s_inner = nzb; nzb_s_outer = nzb … … 327 363 vertical_influence(0) = vertical_influence(1) 328 364 329 DO i = nxl 1, nxr+1330 DO j = nys 1, nyn+1365 DO i = nxlg, nxrg 366 DO j = nysg, nyng 331 367 DO k = nzb_s_inner(j,i) + 1, & 332 368 nzb_s_inner(j,i) + vertical_influence(nzb_s_inner(j,i)) … … 477 513 nzb_local(:,gls:1) = nzb_local(:,nxgls+1:nx) 478 514 nzb_local(:,nx+1:nx+gls) = nzb_local(:,0:gls1) 515 516 479 517 480 518 GOTO 12 … … 588 626 ! nzb_s_outer: 589 627 ! extend nzb_local east/westwards first, then north/southwards 590 nzb_tmp = nzb_local( 1:ny+1,1:nx+1)628 nzb_tmp = nzb_local(nbgp:ny+nbgp,nbgp:nx+nbgp) 591 629 DO j = 1, ny + 1 592 630 DO i = 0, nx … … 620 658 ! nzb_u_inner: 621 659 ! extend nzb_local rightwards only 622 nzb_tmp = nzb_local( 1:ny+1,1:nx+1)660 nzb_tmp = nzb_local(nbgp:ny+nbgp,nbgp:nx+nbgp) 623 661 DO j = 1, ny + 1 624 662 DO i = 0, nx + 1 … … 626 664 ENDDO 627 665 ENDDO 628 nzb_u_inner = nzb_tmp(nys 1:nyn+1,nxl1:nxr+1)666 nzb_u_inner = nzb_tmp(nysg:nyng,nxlg:nxrg) 629 667 630 668 ! … … 652 690 ! nzb_v_inner: 653 691 ! extend nzb_local northwards only 654 nzb_tmp = nzb_local( 1:ny+1,1:nx+1)692 nzb_tmp = nzb_local(nbgp:ny+nbgp,nbgp:nx+nbgp) 655 693 DO i = 1, nx + 1 656 694 DO j = 0, ny + 1 … … 658 696 ENDDO 659 697 ENDDO 660 nzb_v_inner = nzb_tmp(nys 1:nyn+1,nxl1:nxr+1)698 nzb_v_inner = nzb_tmp(nysnbgp:nyn+nbgp,nxlnbgp:nxr+nbgp) 661 699 662 700 ! … … 1096 1134 ! 1097 1135 ! Need to set lateral boundary conditions for l_wall 1098 CALL exchange_horiz( l_wall ) 1136 1137 CALL exchange_horiz( l_wall, nbgp ) 1099 1138 1100 1139 DEALLOCATE( corner_nl, corner_nr, corner_sl, corner_sr, nzb_local, &
Note: See TracChangeset
for help on using the changeset viewer.