Changeset 114 for palm/trunk
- Timestamp:
- Oct 10, 2007 12:03:15 AM (17 years ago)
- Location:
- palm/trunk
- Files:
-
- 12 edited
-
DOC/app/chapter_4.1.html (modified) (3 diffs)
-
SCRIPTS/mrun (modified) (1 diff)
-
SOURCE/CURRENT_MODIFICATIONS (modified) (2 diffs)
-
SOURCE/advec_particles.f90 (modified) (3 diffs)
-
SOURCE/check_parameters.f90 (modified) (2 diffs)
-
SOURCE/init_grid.f90 (modified) (9 diffs)
-
SOURCE/init_pegrid.f90 (modified) (2 diffs)
-
SOURCE/modules.f90 (modified) (2 diffs)
-
SOURCE/palm.f90 (modified) (1 diff)
-
SOURCE/poismg.f90 (modified) (20 diffs)
-
SOURCE/user_interface.f90 (modified) (4 diffs)
-
UTIL/combine_plot_fields.f90 (modified) (3 diffs, 1 prop)
Legend:
- Unmodified
- Added
- Removed
-
palm/trunk/DOC/app/chapter_4.1.html
r111 r114 16 16 17 17 <title>PALM chapter 4.1</title></head> 18 19 18 <body> 20 19 … … 9260 9259 9261 9260 In case of ocean runs (see <a href="chapter_4.1.html#ocean">ocean</a>), 9262 this parameter gives the velocity valueat the sea surface, which is9261 this parameter gives the geostrophic velocity value (i.e. the pressure gradient) at the sea surface, which is 9263 9262 at k=nzt. The profile is then constructed from the surface down to the 9264 9263 bottom of the model.<br> … … 10168 10167 10169 10168 In case of ocean runs (see <a href="chapter_4.1.html#ocean">ocean</a>), 10170 this parameter gives the velocity valueat the sea surface, which is10169 this parameter gives the geostrophic velocity value (i.e. the pressure gradient) at the sea surface, which is 10171 10170 at k=nzt. The profile is then constructed from the surface down to the 10172 10171 bottom of the model.</td> -
palm/trunk/SCRIPTS/mrun
r109 r114 200 200 read_from_config="" 201 201 restart_run=false 202 return_addres=$(nslookup `hostname` 2>&1 | grep "Address:" | tail -1 | awk '{print $2}')202 # return_addres=$(nslookup `hostname` 2>&1 | grep "Address:" | tail -1 | awk '{print $2}') 203 203 if [[ $return_addres = 130.75.105.158 ]] 204 204 then -
palm/trunk/SOURCE/CURRENT_MODIFICATIONS
r110 r114 1 1 New: 2 2 --- 3 4 Pressure boundary conditions for vertical walls added to the multigrid solver. 5 They are applied using new wall flag arrays (wall_flags_..) which are defined 6 for each grid level. New argument gls added to routine user_init_grid 7 (user_interface). 8 9 init_grid, init_pegrid, modules, user_interface 3 10 4 11 … … 9 16 Errors: 10 17 ------ 18 19 Bugfix: pleft/pright changed to pnorth/psouth in sendrecv of particle tail 20 numbers along y (advec_particles) 21 22 Bugfix: model_string needed a default value (combine_plot_fields) 23 24 advec_particles, combine_plot_fields -
palm/trunk/SOURCE/advec_particles.f90
r110 r114 4 4 ! Actual revisions: 5 5 ! ----------------- 6 ! Bugfix: pleft/pright changed to pnorth/psouth in sendrecv of particle tail 7 ! numbers along y 6 8 ! TEST: PRINT statements on unit 9 (commented out) 7 9 ! … … 2773 2775 IF ( use_particle_tails ) THEN 2774 2776 2775 CALL MPI_SENDRECV( trspt_count, 1, MPI_INTEGER, p left,0, &2776 trnpt_count_recv, 1, MPI_INTEGER, p right, 0, &2777 CALL MPI_SENDRECV( trspt_count, 1, MPI_INTEGER, psouth, 0, & 2778 trnpt_count_recv, 1, MPI_INTEGER, pnorth, 0, & 2777 2779 comm2d, status, ierr ) 2778 2780 … … 2850 2852 IF ( use_particle_tails ) THEN 2851 2853 2852 CALL MPI_SENDRECV( trnpt_count, 1, MPI_INTEGER, p right, 0, &2853 trspt_count_recv, 1, MPI_INTEGER, p left,0, &2854 CALL MPI_SENDRECV( trnpt_count, 1, MPI_INTEGER, pnorth, 0, & 2855 trspt_count_recv, 1, MPI_INTEGER, psouth, 0, & 2854 2856 comm2d, status, ierr ) 2855 2857 -
palm/trunk/SOURCE/check_parameters.f90
r110 r114 4 4 ! Actual revisions: 5 5 ! ----------------- 6 ! 6 ! Multigrid solver allows topography 7 7 ! 8 8 ! Former revisions: … … 286 286 WRITE( action, '(A,A)' ) 'timestep_scheme = ', timestep_scheme 287 287 ENDIF 288 IF ( psolver == ' multigrid' .OR. psolver == 'sor' ) THEN288 IF ( psolver == 'sor' ) THEN 289 289 WRITE( action, '(A,A)' ) 'psolver = ', psolver 290 290 ENDIF -
palm/trunk/SOURCE/init_grid.f90
r98 r114 4 4 ! Actual revisions: 5 5 ! ----------------- 6 ! 6 ! Calculation of wall flag arrays 7 7 ! 8 8 ! Former revisions: … … 44 44 IMPLICIT NONE 45 45 46 INTEGER :: bh, blx, bly, bxl, bxr, byn, bys, i, i_center, j, j_center, k, &47 l, nzb_si, nzt_l, vi46 INTEGER :: bh, blx, bly, bxl, bxr, byn, bys, gls, i, inc, i_center, j, & 47 j_center, k, l, nxl_l, nxr_l, nyn_l, nys_l, nzb_si, nzt_l, vi 48 48 49 49 INTEGER, DIMENSION(:), ALLOCATABLE :: vertical_influence … … 217 217 ! 218 218 !-- Allocate outer and inner index arrays for topography and set 219 !-- defaults 220 ALLOCATE( corner_nl(nys:nyn,nxl:nxr), corner_nr(nys:nyn,nxl:nxr), & 221 corner_sl(nys:nyn,nxl:nxr), corner_sr(nys:nyn,nxl:nxr), & 222 nzb_local(-1:ny+1,-1:nx+1), nzb_tmp(-1:ny+1,-1:nx+1), & 223 wall_l(nys:nyn,nxl:nxr), wall_n(nys:nyn,nxl:nxr), & 219 !-- defaults. 220 !-- nzb_local has to contain additional layers of ghost points for calculating 221 !-- the flag arrays needed for the multigrid method 222 gls = 2**( maximum_grid_level ) 223 ALLOCATE( corner_nl(nys:nyn,nxl:nxr), corner_nr(nys:nyn,nxl:nxr), & 224 corner_sl(nys:nyn,nxl:nxr), corner_sr(nys:nyn,nxl:nxr), & 225 nzb_local(-gls:ny+gls,-gls:nx+gls), nzb_tmp(-1:ny+1,-1:nx+1), & 226 wall_l(nys:nyn,nxl:nxr), wall_n(nys:nyn,nxl:nxr), & 224 227 wall_r(nys:nyn,nxl:nxr), wall_s(nys:nyn,nxl:nxr) ) 225 228 ALLOCATE( fwxm(nys-1:nyn+1,nxl-1:nxr+1), fwxp(nys-1:nyn+1,nxl-1:nxr+1), & … … 388 391 ENDDO 389 392 ENDDO 390 nzb_local(-1,0:nx) = nzb_local(ny,0:nx) 391 nzb_local(ny+1,0:nx) = nzb_local(0,0:nx) 392 nzb_local(:,-1) = nzb_local(:,nx) 393 nzb_local(:,nx+1) = nzb_local(:,0) 393 ! 394 !-- Add cyclic boundaries (additional layers are for calculating flag 395 !-- arrays needed for the multigrid sover) 396 nzb_local(-gls:-1,0:nx) = nzb_local(ny-gls+1:ny,0:nx) 397 nzb_local(ny+1:ny+gls,0:nx) = nzb_local(0:gls-1,0:nx) 398 nzb_local(:,-gls:-1) = nzb_local(:,nx-gls+1:nx) 399 nzb_local(:,nx+1:nx+gls) = nzb_local(:,0:gls-1) 394 400 395 401 GOTO 12 … … 413 419 !-- case in the user interface. There, the subroutine user_init_grid 414 420 !-- checks which of these two conditions applies. 415 CALL user_init_grid( nzb_local )421 CALL user_init_grid( gls, nzb_local ) 416 422 417 423 END SELECT 424 425 ! 426 !-- Test output of nzb_local -1:ny+1,-1:nx+1 427 WRITE (9,*) '*** nzb_local ***' 428 DO j = ny+1, -1, -1 429 WRITE (9,'(194(1X,I2))') ( nzb_local(j,i), i = -1, nx+1 ) 430 ENDDO 418 431 419 432 ! … … 478 491 !-- nzb_s_outer: 479 492 !-- extend nzb_local east-/westwards first, then north-/southwards 480 nzb_tmp = nzb_local 493 nzb_tmp = nzb_local(-1:ny+1,-1:nx+1) 481 494 DO j = -1, ny + 1 482 495 DO i = 0, nx … … 510 523 !-- nzb_u_inner: 511 524 !-- extend nzb_local rightwards only 512 nzb_tmp = nzb_local 525 nzb_tmp = nzb_local(-1:ny+1,-1:nx+1) 513 526 DO j = -1, ny + 1 514 527 DO i = 0, nx + 1 … … 542 555 !-- nzb_v_inner: 543 556 !-- extend nzb_local northwards only 544 nzb_tmp = nzb_local 557 nzb_tmp = nzb_local(-1:ny+1,-1:nx+1) 545 558 DO i = -1, nx + 1 546 559 DO j = 0, ny + 1 … … 728 741 729 742 ! 743 !-- Calculate wall flag arrays for the multigrid method 744 IF ( psolver == 'multigrid' ) THEN 745 ! 746 !-- Gridpoint increment of the current level 747 inc = 1 748 749 DO l = maximum_grid_level, 1 , -1 750 751 nxl_l = nxl_mg(l) 752 nxr_l = nxr_mg(l) 753 nys_l = nys_mg(l) 754 nyn_l = nyn_mg(l) 755 nzt_l = nzt_mg(l) 756 757 ! 758 !-- Assign the flag level to be calculated 759 SELECT CASE ( l ) 760 CASE ( 1 ) 761 flags => wall_flags_1 762 CASE ( 2 ) 763 flags => wall_flags_2 764 CASE ( 3 ) 765 flags => wall_flags_3 766 CASE ( 4 ) 767 flags => wall_flags_4 768 CASE ( 5 ) 769 flags => wall_flags_5 770 CASE ( 6 ) 771 flags => wall_flags_6 772 CASE ( 7 ) 773 flags => wall_flags_7 774 CASE ( 8 ) 775 flags => wall_flags_8 776 CASE ( 9 ) 777 flags => wall_flags_9 778 CASE ( 10 ) 779 flags => wall_flags_10 780 END SELECT 781 782 ! 783 !-- Depending on the grid level, set the respective bits in case of 784 !-- neighbouring walls 785 !-- Bit 0: wall to the bottom 786 !-- Bit 1: wall to the top (not realized in remaining PALM code so far) 787 !-- Bit 2: wall to the south 788 !-- Bit 3: wall to the north 789 !-- Bit 4: wall to the left 790 !-- Bit 5: wall to the right 791 !-- Bit 6: inside / outside building (1/0) 792 793 flags = 0 794 795 DO i = nxl_l-1, nxr_l+1 796 DO j = nys_l-1, nyn_l+1 797 DO k = nzb, nzt_l+1 798 799 ! 800 !-- Inside/outside building (inside building does not need 801 !-- further tests for walls) 802 IF ( k*inc <= nzb_local(j*inc,i*inc) ) THEN 803 804 flags(k,j,i) = IBSET( flags(k,j,i), 6 ) 805 806 ELSE 807 ! 808 !-- Bottom wall 809 IF ( (k-1)*inc <= nzb_local(j*inc,i*inc) ) THEN 810 flags(k,j,i) = IBSET( flags(k,j,i), 0 ) 811 ENDIF 812 ! 813 !-- South wall 814 IF ( k*inc <= nzb_local((j-1)*inc,i*inc) ) THEN 815 flags(k,j,i) = IBSET( flags(k,j,i), 2 ) 816 ENDIF 817 ! 818 !-- North wall 819 IF ( k*inc <= nzb_local((j+1)*inc,i*inc) ) THEN 820 flags(k,j,i) = IBSET( flags(k,j,i), 3 ) 821 ENDIF 822 ! 823 !-- Left wall 824 IF ( k*inc <= nzb_local(j*inc,(i-1)*inc) ) THEN 825 flags(k,j,i) = IBSET( flags(k,j,i), 4 ) 826 ENDIF 827 ! 828 !-- Right wall 829 IF ( k*inc <= nzb_local(j*inc,(i+1)*inc) ) THEN 830 flags(k,j,i) = IBSET( flags(k,j,i), 5 ) 831 ENDIF 832 833 ENDIF 834 835 ENDDO 836 ENDDO 837 ENDDO 838 839 ! 840 !-- Test output of flag arrays 841 i = nxl_l 842 WRITE (9,*) ' ' 843 WRITE (9,*) '*** mg level ', l, ' ***', mg_switch_to_pe0_level 844 WRITE (9,*) ' inc=', inc, ' i =', nxl_l 845 WRITE (9,*) ' nxl_l',nxl_l,' nxr_l=',nxr_l,' nys_l=',nys_l,' nyn_l=',nyn_l 846 DO k = nzt_l+1, nzb, -1 847 WRITE (9,'(194(1X,I2))') ( flags(k,j,i), j = nys_l-1, nyn_l+1 ) 848 ENDDO 849 850 inc = inc * 2 851 852 ENDDO 853 854 ENDIF 855 856 ! 730 857 !-- In case of topography: limit near-wall mixing length l_wall further: 731 858 !-- Go through all points of the subdomain one by one and look for the closest -
palm/trunk/SOURCE/init_pegrid.f90
r110 r114 4 4 ! Actual revisions: 5 5 ! ----------------- 6 ! Allocation of wall flag arrays for multigrid solver 6 7 ! TEST OUTPUT (TO BE REMOVED) logging mpi2 ierr values 7 8 ! … … 915 916 ENDIF 916 917 918 ! 919 !-- Allocate wall flag arrays used in the multigrid solver 920 IF ( psolver == 'multigrid' ) THEN 921 922 DO i = maximum_grid_level, 1, -1 923 924 SELECT CASE ( i ) 925 926 CASE ( 1 ) 927 ALLOCATE( wall_flags_1(nzb:nzt_mg(i)+1, & 928 nys_mg(i)-1:nyn_mg(i)+1, & 929 nxl_mg(i)-1:nxr_mg(i)+1) ) 930 931 CASE ( 2 ) 932 ALLOCATE( wall_flags_2(nzb:nzt_mg(i)+1, & 933 nys_mg(i)-1:nyn_mg(i)+1, & 934 nxl_mg(i)-1:nxr_mg(i)+1) ) 935 936 CASE ( 3 ) 937 ALLOCATE( wall_flags_3(nzb:nzt_mg(i)+1, & 938 nys_mg(i)-1:nyn_mg(i)+1, & 939 nxl_mg(i)-1:nxr_mg(i)+1) ) 940 941 CASE ( 4 ) 942 ALLOCATE( wall_flags_4(nzb:nzt_mg(i)+1, & 943 nys_mg(i)-1:nyn_mg(i)+1, & 944 nxl_mg(i)-1:nxr_mg(i)+1) ) 945 946 CASE ( 5 ) 947 ALLOCATE( wall_flags_5(nzb:nzt_mg(i)+1, & 948 nys_mg(i)-1:nyn_mg(i)+1, & 949 nxl_mg(i)-1:nxr_mg(i)+1) ) 950 951 CASE ( 6 ) 952 ALLOCATE( wall_flags_6(nzb:nzt_mg(i)+1, & 953 nys_mg(i)-1:nyn_mg(i)+1, & 954 nxl_mg(i)-1:nxr_mg(i)+1) ) 955 956 CASE ( 7 ) 957 ALLOCATE( wall_flags_7(nzb:nzt_mg(i)+1, & 958 nys_mg(i)-1:nyn_mg(i)+1, & 959 nxl_mg(i)-1:nxr_mg(i)+1) ) 960 961 CASE ( 8 ) 962 ALLOCATE( wall_flags_8(nzb:nzt_mg(i)+1, & 963 nys_mg(i)-1:nyn_mg(i)+1, & 964 nxl_mg(i)-1:nxr_mg(i)+1) ) 965 966 CASE ( 9 ) 967 ALLOCATE( wall_flags_9(nzb:nzt_mg(i)+1, & 968 nys_mg(i)-1:nyn_mg(i)+1, & 969 nxl_mg(i)-1:nxr_mg(i)+1) ) 970 971 CASE ( 10 ) 972 ALLOCATE( wall_flags_10(nzb:nzt_mg(i)+1, & 973 nys_mg(i)-1:nyn_mg(i)+1, & 974 nxl_mg(i)-1:nxr_mg(i)+1) ) 975 976 CASE DEFAULT 977 IF ( myid == 0 ) PRINT*, '+++ init_pegrid: more than 10 ', & 978 ' multigrid levels' 979 CALL local_stop 980 981 END SELECT 982 983 ENDDO 984 985 ENDIF 986 917 987 END SUBROUTINE init_pegrid -
palm/trunk/SOURCE/modules.f90
r110 r114 5 5 ! Actual revisions: 6 6 ! ----------------- 7 ! 7 ! +flags, wall_flags_1..10 8 8 ! 9 9 ! Former revisions: … … 589 589 nzb_v_outer, nzb_w_inner, nzb_w_outer, nzb_2d 590 590 591 INTEGER, DIMENSION(:,:,:), POINTER :: flags 592 593 INTEGER, DIMENSION(:,:,:), ALLOCATABLE, TARGET :: wall_flags_1, & 594 wall_flags_2, wall_flags_3, wall_flags_4, wall_flags_5, & 595 wall_flags_6, wall_flags_7, wall_flags_8, wall_flags_9, & 596 wall_flags_10 597 591 598 SAVE 592 599 -
palm/trunk/SOURCE/palm.f90
r110 r114 66 66 INTEGER :: i, run_description_header_i(80) 67 67 68 version = 'PALM 3.4 '68 version = 'PALM 3.4a' 69 69 70 70 #if defined( __parallel ) -
palm/trunk/SOURCE/poismg.f90
r77 r114 8 8 ! Actual revisions: 9 9 ! ----------------- 10 ! 10 ! Boundary conditions at walls are implicitly set using flag arrays. Only 11 ! Neumann BC is allowed. Upper walls are still not realized. 12 ! Bottom and top BCs for array f_mg in restrict removed because boundary 13 ! values are not needed (right hand side of SOR iteration) 11 14 ! 12 15 ! Former revisions: … … 150 153 l = grid_level 151 154 155 ! 156 !-- Choose flag array of this level 157 SELECT CASE ( l ) 158 CASE ( 1 ) 159 flags => wall_flags_1 160 CASE ( 2 ) 161 flags => wall_flags_2 162 CASE ( 3 ) 163 flags => wall_flags_3 164 CASE ( 4 ) 165 flags => wall_flags_4 166 CASE ( 5 ) 167 flags => wall_flags_5 168 CASE ( 6 ) 169 flags => wall_flags_6 170 CASE ( 7 ) 171 flags => wall_flags_7 172 CASE ( 8 ) 173 flags => wall_flags_8 174 CASE ( 9 ) 175 flags => wall_flags_9 176 CASE ( 10 ) 177 flags => wall_flags_10 178 END SELECT 179 152 180 !$OMP PARALLEL PRIVATE (i,j,k) 153 181 !$OMP DO … … 155 183 DO j = nys_mg(l), nyn_mg(l) 156 184 DO k = nzb+1, nzt_mg(l) 157 r(k,j,i) = f_mg(k,j,i) & 158 - ddx2_mg(l) * ( p_mg(k,j,i+1) + p_mg(k,j,i-1) ) & 159 - ddy2_mg(l) * ( p_mg(k,j+1,i) + p_mg(k,j-1,i) ) & 160 - f2_mg(k,l) * p_mg(k+1,j,i) & 161 - f3_mg(k,l) * p_mg(k-1,j,i) & 185 r(k,j,i) = f_mg(k,j,i) & 186 - ddx2_mg(l) * & 187 ( p_mg(k,j,i+1) + IBITS( flags(k,j,i), 5, 1 ) * & 188 ( p_mg(k,j,i) - p_mg(k,j,i+1) ) + & 189 p_mg(k,j,i-1) + IBITS( flags(k,j,i), 4, 1 ) * & 190 ( p_mg(k,j,i) - p_mg(k,j,i-1) ) ) & 191 - ddy2_mg(l) * & 192 ( p_mg(k,j+1,i) + IBITS( flags(k,j,i), 3, 1 ) * & 193 ( p_mg(k,j,i) - p_mg(k,j+1,i) ) + & 194 p_mg(k,j-1,i) + IBITS( flags(k,j,i), 2, 1 ) * & 195 ( p_mg(k,j,i) - p_mg(k,j-1,i) ) ) & 196 - f2_mg(k,l) * p_mg(k+1,j,i) & 197 - f3_mg(k,l) * & 198 ( p_mg(k-1,j,i) + IBITS( flags(k,j,i), 0, 1 ) * & 199 ( p_mg(k,j,i) - p_mg(k-1,j,i) ) ) & 162 200 + f1_mg(k,l) * p_mg(k,j,i) 201 ! 202 !-- Residual within topography should be zero 203 r(k,j,i) = r(k,j,i) * ( 1.0 - IBITS( flags(k,j,i), 6, 1 ) ) 163 204 ENDDO 164 205 ENDDO … … 181 222 182 223 ! 183 !-- Bottom and top boundary conditions 184 IF ( ibc_p_b == 1 ) THEN 185 r(nzb,:,: ) = r(nzb+1,:,:) 186 ELSE 187 r(nzb,:,: ) = 0.0 188 ENDIF 189 224 !-- Top boundary condition 225 !-- A Neumann boundary condition for r is implicitly set in routine restrict 190 226 IF ( ibc_p_t == 1 ) THEN 191 227 r(nzt_mg(l)+1,:,: ) = r(nzt_mg(l),:,:) … … 217 253 INTEGER :: i, ic, j, jc, k, kc, l 218 254 255 REAL :: rkjim, rkjip, rkjmi, rkjmim, rkjmip, rkjpi, rkjpim, rkjpip, & 256 rkmji, rkmjim, rkmjip, rkmjmi, rkmjmim, rkmjmip, rkmjpi, rkmjpim, & 257 rkmjpip 258 219 259 REAL, DIMENSION(nzb:nzt_mg(grid_level)+1, & 220 260 nys_mg(grid_level)-1:nyn_mg(grid_level)+1, & … … 229 269 l = grid_level 230 270 271 ! 272 !-- Choose flag array of the upper level 273 SELECT CASE ( l ) 274 CASE ( 1 ) 275 flags => wall_flags_1 276 CASE ( 2 ) 277 flags => wall_flags_2 278 CASE ( 3 ) 279 flags => wall_flags_3 280 CASE ( 4 ) 281 flags => wall_flags_4 282 CASE ( 5 ) 283 flags => wall_flags_5 284 CASE ( 6 ) 285 flags => wall_flags_6 286 CASE ( 7 ) 287 flags => wall_flags_7 288 CASE ( 8 ) 289 flags => wall_flags_8 290 CASE ( 9 ) 291 flags => wall_flags_9 292 CASE ( 10 ) 293 flags => wall_flags_10 294 END SELECT 295 231 296 !$OMP PARALLEL PRIVATE (i,j,k,ic,jc,kc) 232 297 !$OMP DO … … 237 302 DO kc = nzb+1, nzt_mg(l) 238 303 k = 2*kc-1 304 ! 305 !-- Use implicit Neumann BCs if the respective gridpoint is inside 306 !-- the building 307 rkjim = r(k,j,i-1) + IBITS( flags(k,j,i-1), 6, 1 ) * & 308 ( r(k,j,i) - r(k,j,i-1) ) 309 rkjip = r(k,j,i+1) + IBITS( flags(k,j,i+1), 6, 1 ) * & 310 ( r(k,j,i) - r(k,j,i+1) ) 311 rkjpi = r(k,j+1,i) + IBITS( flags(k,j+1,i), 6, 1 ) * & 312 ( r(k,j,i) - r(k,j+1,i) ) 313 rkjmi = r(k,j-1,i) + IBITS( flags(k,j-1,i), 6, 1 ) * & 314 ( r(k,j,i) - r(k,j-1,i) ) 315 rkjmim = r(k,j-1,i-1) + IBITS( flags(k,j-1,i-1), 6, 1 ) * & 316 ( r(k,j,i) - r(k,j-1,i-1) ) 317 rkjpim = r(k,j+1,i-1) + IBITS( flags(k,j+1,i-1), 6, 1 ) * & 318 ( r(k,j,i) - r(k,j+1,i-1) ) 319 rkjmip = r(k,j-1,i+1) + IBITS( flags(k,j-1,i+1), 6, 1 ) * & 320 ( r(k,j,i) - r(k,j-1,i+1) ) 321 rkjpip = r(k,j+1,i+1) + IBITS( flags(k,j+1,i+1), 6, 1 ) * & 322 ( r(k,j,i) - r(k,j+1,i+1) ) 323 rkmji = r(k-1,j,i) + IBITS( flags(k-1,j,i), 6, 1 ) * & 324 ( r(k,j,i) - r(k-1,j,i) ) 325 rkmjim = r(k-1,j,i-1) + IBITS( flags(k-1,j,i-1), 6, 1 ) * & 326 ( r(k,j,i) - r(k-1,j,i-1) ) 327 rkmjip = r(k-1,j,i+1) + IBITS( flags(k-1,j,i+1), 6, 1 ) * & 328 ( r(k,j,i) - r(k-1,j,i+1) ) 329 rkmjpi = r(k-1,j+1,i) + IBITS( flags(k-1,j+1,i), 6, 1 ) * & 330 ( r(k,j,i) - r(k-1,j+1,i) ) 331 rkmjmi = r(k-1,j-1,i) + IBITS( flags(k-1,j-1,i), 6, 1 ) * & 332 ( r(k,j,i) - r(k-1,j-1,i) ) 333 rkmjmim = r(k-1,j-1,i-1) + IBITS( flags(k-1,j-1,i-1), 6, 1 ) * & 334 ( r(k,j,i) - r(k-1,j-1,i-1) ) 335 rkmjpim = r(k-1,j+1,i-1) + IBITS( flags(k-1,j+1,i-1), 6, 1 ) * & 336 ( r(k,j,i) - r(k-1,j+1,i-1) ) 337 rkmjmip = r(k-1,j-1,i+1) + IBITS( flags(k-1,j-1,i+1), 6, 1 ) * & 338 ( r(k,j,i) - r(k-1,j-1,i+1) ) 339 rkmjpip = r(k-1,j+1,i+1) + IBITS( flags(k-1,j+1,i+1), 6, 1 ) * & 340 ( r(k,j,i) - r(k-1,j+1,i+1) ) 341 239 342 f_mg(kc,jc,ic) = 1.0 / 64.0 * ( & 240 343 8.0 * r(k,j,i) & 241 + 4.0 * ( r (k,j,i-1) + r(k,j,i+1) +&242 r (k,j+1,i) + r(k,j-1,i) )&243 + 2.0 * ( r (k,j-1,i-1) + r(k,j+1,i-1) +&244 r (k,j-1,i+1) + r(k,j+1,i+1) )&245 + 4.0 * r (k-1,j,i)&246 + 2.0 * ( r (k-1,j,i-1) + r(k-1,j,i+1) +&247 r (k-1,j+1,i) + r(k-1,j-1,i) )&248 + ( r (k-1,j-1,i-1) + r(k-1,j+1,i-1) +&249 r (k-1,j-1,i+1) + r(k-1,j+1,i+1) )&344 + 4.0 * ( rkjim + rkjip + & 345 rkjpi + rkjmi ) & 346 + 2.0 * ( rkjmim + rkjpim + & 347 rkjmip + rkjpip ) & 348 + 4.0 * rkmji & 349 + 2.0 * ( rkmjim + rkmjim + & 350 rkmjpi + rkmjmi ) & 351 + ( rkmjmim + rkmjpim + & 352 rkmjmip + rkmjpip ) & 250 353 + 4.0 * r(k+1,j,i) & 251 354 + 2.0 * ( r(k+1,j,i-1) + r(k+1,j,i+1) + & … … 254 357 r(k+1,j-1,i+1) + r(k+1,j+1,i+1) ) & 255 358 ) 359 360 ! f_mg(kc,jc,ic) = 1.0 / 64.0 * ( & 361 ! 8.0 * r(k,j,i) & 362 ! + 4.0 * ( r(k,j,i-1) + r(k,j,i+1) + & 363 ! r(k,j+1,i) + r(k,j-1,i) ) & 364 ! + 2.0 * ( r(k,j-1,i-1) + r(k,j+1,i-1) + & 365 ! r(k,j-1,i+1) + r(k,j+1,i+1) ) & 366 ! + 4.0 * r(k-1,j,i) & 367 ! + 2.0 * ( r(k-1,j,i-1) + r(k-1,j,i+1) + & 368 ! r(k-1,j+1,i) + r(k-1,j-1,i) ) & 369 ! + ( r(k-1,j-1,i-1) + r(k-1,j+1,i-1) + & 370 ! r(k-1,j-1,i+1) + r(k-1,j+1,i+1) ) & 371 ! + 4.0 * r(k+1,j,i) & 372 ! + 2.0 * ( r(k+1,j,i-1) + r(k+1,j,i+1) + & 373 ! r(k+1,j+1,i) + r(k+1,j-1,i) ) & 374 ! + ( r(k+1,j-1,i-1) + r(k+1,j+1,i-1) + & 375 ! r(k+1,j-1,i+1) + r(k+1,j+1,i+1) ) & 376 ! ) 256 377 ENDDO 257 378 ENDDO … … 275 396 ! 276 397 !-- Bottom and top boundary conditions 277 IF ( ibc_p_b == 1 ) THEN278 f_mg(nzb,:,: ) = f_mg(nzb+1,:,:)279 ELSE280 f_mg(nzb,:,: ) = 0.0281 ENDIF282 283 IF ( ibc_p_t == 1 ) THEN284 f_mg(nzt_mg(l)+1,:,: ) = f_mg(nzt_mg(l),:,:)285 ELSE286 f_mg(nzt_mg(l)+1,:,: ) = 0.0287 ENDIF398 ! IF ( ibc_p_b == 1 ) THEN 399 ! f_mg(nzb,:,: ) = f_mg(nzb+1,:,:) 400 ! ELSE 401 ! f_mg(nzb,:,: ) = 0.0 402 ! ENDIF 403 ! 404 ! IF ( ibc_p_t == 1 ) THEN 405 ! f_mg(nzt_mg(l)+1,:,: ) = f_mg(nzt_mg(l),:,:) 406 ! ELSE 407 ! f_mg(nzt_mg(l)+1,:,: ) = 0.0 408 ! ENDIF 288 409 289 410 … … 412 533 LOGICAL :: unroll 413 534 535 REAL :: wall_left, wall_north, wall_right, wall_south, wall_total, wall_top 536 414 537 REAL, DIMENSION(nzb:nzt_mg(grid_level)+1, & 415 538 nys_mg(grid_level)-1:nyn_mg(grid_level)+1, & … … 418 541 419 542 l = grid_level 543 544 ! 545 !-- Choose flag array of this level 546 SELECT CASE ( l ) 547 CASE ( 1 ) 548 flags => wall_flags_1 549 CASE ( 2 ) 550 flags => wall_flags_2 551 CASE ( 3 ) 552 flags => wall_flags_3 553 CASE ( 4 ) 554 flags => wall_flags_4 555 CASE ( 5 ) 556 flags => wall_flags_5 557 CASE ( 6 ) 558 flags => wall_flags_6 559 CASE ( 7 ) 560 flags => wall_flags_7 561 CASE ( 8 ) 562 flags => wall_flags_8 563 CASE ( 9 ) 564 flags => wall_flags_9 565 CASE ( 10 ) 566 flags => wall_flags_10 567 END SELECT 420 568 421 569 unroll = ( MOD( nyn_mg(l)-nys_mg(l)+1, 4 ) == 0 .AND. & … … 434 582 DO j = nys_mg(l) + 2 - colour, nyn_mg(l), 2 435 583 DO k = nzb+1, nzt_mg(l), 2 584 ! p_mg(k,j,i) = 1.0 / f1_mg(k,l) * ( & 585 ! ddx2_mg(l) * ( p_mg(k,j,i+1) + p_mg(k,j,i-1) ) & 586 ! + ddy2_mg(l) * ( p_mg(k,j+1,i) + p_mg(k,j-1,i) ) & 587 ! + f2_mg(k,l) * p_mg(k+1,j,i) & 588 ! + f3_mg(k,l) * p_mg(k-1,j,i) - f_mg(k,j,i) & 589 ! ) 590 436 591 p_mg(k,j,i) = 1.0 / f1_mg(k,l) * ( & 437 ddx2_mg(l) * ( p_mg(k,j,i+1) + p_mg(k,j,i-1) ) & 438 + ddy2_mg(l) * ( p_mg(k,j+1,i) + p_mg(k,j-1,i) ) & 439 + f2_mg(k,l) * p_mg(k+1,j,i) & 440 + f3_mg(k,l) * p_mg(k-1,j,i) - f_mg(k,j,i) & 441 ) 592 ddx2_mg(l) * & 593 ( p_mg(k,j,i+1) + IBITS( flags(k,j,i), 5, 1 ) * & 594 ( p_mg(k,j,i) - p_mg(k,j,i+1) ) + & 595 p_mg(k,j,i-1) + IBITS( flags(k,j,i), 4, 1 ) * & 596 ( p_mg(k,j,i) - p_mg(k,j,i-1) ) ) & 597 + ddy2_mg(l) * & 598 ( p_mg(k,j+1,i) + IBITS( flags(k,j,i), 3, 1 ) * & 599 ( p_mg(k,j,i) - p_mg(k,j+1,i) ) + & 600 p_mg(k,j-1,i) + IBITS( flags(k,j,i), 2, 1 ) * & 601 ( p_mg(k,j,i) - p_mg(k,j-1,i) ) ) & 602 + f2_mg(k,l) * p_mg(k+1,j,i) & 603 + f3_mg(k,l) * & 604 ( p_mg(k-1,j,i) + IBITS( flags(k,j,i), 0, 1 ) * & 605 ( p_mg(k,j,i) - p_mg(k-1,j,i) ) ) & 606 - f_mg(k,j,i) ) 442 607 ENDDO 443 608 ENDDO … … 448 613 DO k = nzb+1, nzt_mg(l), 2 449 614 p_mg(k,j,i) = 1.0 / f1_mg(k,l) * ( & 450 ddx2_mg(l) * ( p_mg(k,j,i+1) + p_mg(k,j,i-1) ) & 451 + ddy2_mg(l) * ( p_mg(k,j+1,i) + p_mg(k,j-1,i) ) & 452 + f2_mg(k,l) * p_mg(k+1,j,i) & 453 + f3_mg(k,l) * p_mg(k-1,j,i) - f_mg(k,j,i) & 454 ) 615 ddx2_mg(l) * & 616 ( p_mg(k,j,i+1) + IBITS( flags(k,j,i), 5, 1 ) * & 617 ( p_mg(k,j,i) - p_mg(k,j,i+1) ) + & 618 p_mg(k,j,i-1) + IBITS( flags(k,j,i), 4, 1 ) * & 619 ( p_mg(k,j,i) - p_mg(k,j,i-1) ) ) & 620 + ddy2_mg(l) * & 621 ( p_mg(k,j+1,i) + IBITS( flags(k,j,i), 3, 1 ) * & 622 ( p_mg(k,j,i) - p_mg(k,j+1,i) ) + & 623 p_mg(k,j-1,i) + IBITS( flags(k,j,i), 2, 1 ) * & 624 ( p_mg(k,j,i) - p_mg(k,j-1,i) ) ) & 625 + f2_mg(k,l) * p_mg(k+1,j,i) & 626 + f3_mg(k,l) * & 627 ( p_mg(k-1,j,i) + IBITS( flags(k,j,i), 0, 1 ) * & 628 ( p_mg(k,j,i) - p_mg(k-1,j,i) ) ) & 629 - f_mg(k,j,i) ) 455 630 ENDDO 456 631 ENDDO … … 461 636 DO k = nzb+2, nzt_mg(l), 2 462 637 p_mg(k,j,i) = 1.0 / f1_mg(k,l) * ( & 463 ddx2_mg(l) * ( p_mg(k,j,i+1) + p_mg(k,j,i-1) ) & 464 + ddy2_mg(l) * ( p_mg(k,j+1,i) + p_mg(k,j-1,i) ) & 465 + f2_mg(k,l) * p_mg(k+1,j,i) & 466 + f3_mg(k,l) * p_mg(k-1,j,i) - f_mg(k,j,i) & 467 ) 638 ddx2_mg(l) * & 639 ( p_mg(k,j,i+1) + IBITS( flags(k,j,i), 5, 1 ) * & 640 ( p_mg(k,j,i) - p_mg(k,j,i+1) ) + & 641 p_mg(k,j,i-1) + IBITS( flags(k,j,i), 4, 1 ) * & 642 ( p_mg(k,j,i) - p_mg(k,j,i-1) ) ) & 643 + ddy2_mg(l) * & 644 ( p_mg(k,j+1,i) + IBITS( flags(k,j,i), 3, 1 ) * & 645 ( p_mg(k,j,i) - p_mg(k,j+1,i) ) + & 646 p_mg(k,j-1,i) + IBITS( flags(k,j,i), 2, 1 ) * & 647 ( p_mg(k,j,i) - p_mg(k,j-1,i) ) ) & 648 + f2_mg(k,l) * p_mg(k+1,j,i) & 649 + f3_mg(k,l) * & 650 ( p_mg(k-1,j,i) + IBITS( flags(k,j,i), 0, 1 ) * & 651 ( p_mg(k,j,i) - p_mg(k-1,j,i) ) ) & 652 - f_mg(k,j,i) ) 468 653 ENDDO 469 654 ENDDO … … 474 659 DO k = nzb+2, nzt_mg(l), 2 475 660 p_mg(k,j,i) = 1.0 / f1_mg(k,l) * ( & 476 ddx2_mg(l) * ( p_mg(k,j,i+1) + p_mg(k,j,i-1) ) & 477 + ddy2_mg(l) * ( p_mg(k,j+1,i) + p_mg(k,j-1,i) ) & 478 + f2_mg(k,l) * p_mg(k+1,j,i) & 479 + f3_mg(k,l) * p_mg(k-1,j,i) - f_mg(k,j,i) & 480 ) 661 ddx2_mg(l) * & 662 ( p_mg(k,j,i+1) + IBITS( flags(k,j,i), 5, 1 ) * & 663 ( p_mg(k,j,i) - p_mg(k,j,i+1) ) + & 664 p_mg(k,j,i-1) + IBITS( flags(k,j,i), 4, 1 ) * & 665 ( p_mg(k,j,i) - p_mg(k,j,i-1) ) ) & 666 + ddy2_mg(l) * & 667 ( p_mg(k,j+1,i) + IBITS( flags(k,j,i), 3, 1 ) * & 668 ( p_mg(k,j,i) - p_mg(k,j+1,i) ) + & 669 p_mg(k,j-1,i) + IBITS( flags(k,j,i), 2, 1 ) * & 670 ( p_mg(k,j,i) - p_mg(k,j-1,i) ) ) & 671 + f2_mg(k,l) * p_mg(k+1,j,i) & 672 + f3_mg(k,l) * & 673 ( p_mg(k-1,j,i) + IBITS( flags(k,j,i), 0, 1 ) * & 674 ( p_mg(k,j,i) - p_mg(k-1,j,i) ) ) & 675 - f_mg(k,j,i) ) 481 676 ENDDO 482 677 ENDDO … … 496 691 j = jj 497 692 p_mg(k,j,i) = 1.0 / f1_mg(k,l) * ( & 498 ddx2_mg(l) * ( p_mg(k,j,i+1) + p_mg(k,j,i-1) ) & 499 + ddy2_mg(l) * ( p_mg(k,j+1,i) + p_mg(k,j-1,i) ) & 500 + f2_mg(k,l) * p_mg(k+1,j,i) & 501 + f3_mg(k,l) * p_mg(k-1,j,i) - f_mg(k,j,i) & 502 ) 693 ddx2_mg(l) * & 694 ( p_mg(k,j,i+1) + IBITS( flags(k,j,i), 5, 1 ) * & 695 ( p_mg(k,j,i) - p_mg(k,j,i+1) ) + & 696 p_mg(k,j,i-1) + IBITS( flags(k,j,i), 4, 1 ) * & 697 ( p_mg(k,j,i) - p_mg(k,j,i-1) ) ) & 698 + ddy2_mg(l) * & 699 ( p_mg(k,j+1,i) + IBITS( flags(k,j,i), 3, 1 ) * & 700 ( p_mg(k,j,i) - p_mg(k,j+1,i) ) + & 701 p_mg(k,j-1,i) + IBITS( flags(k,j,i), 2, 1 ) * & 702 ( p_mg(k,j,i) - p_mg(k,j-1,i) ) ) & 703 + f2_mg(k,l) * p_mg(k+1,j,i) & 704 + f3_mg(k,l) * & 705 ( p_mg(k-1,j,i) + IBITS( flags(k,j,i), 0, 1 ) * & 706 ( p_mg(k,j,i) - p_mg(k-1,j,i) ) ) & 707 - f_mg(k,j,i) ) 503 708 j = jj+2 504 709 p_mg(k,j,i) = 1.0 / f1_mg(k,l) * ( & 505 ddx2_mg(l) * ( p_mg(k,j,i+1) + p_mg(k,j,i-1) ) & 506 + ddy2_mg(l) * ( p_mg(k,j+1,i) + p_mg(k,j-1,i) ) & 507 + f2_mg(k,l) * p_mg(k+1,j,i) & 508 + f3_mg(k,l) * p_mg(k-1,j,i) - f_mg(k,j,i) & 509 ) 510 ! j = jj+4 511 ! p_mg(k,j,i) = 1.0 / f1_mg(k,l) * ( & 512 ! ddx2_mg(l) * ( p_mg(k,j,i+1) + p_mg(k,j,i-1) ) & 513 ! + ddy2_mg(l) * ( p_mg(k,j+1,i) + p_mg(k,j-1,i) ) & 514 ! + f2_mg(k,l) * p_mg(k+1,j,i) & 515 ! + f3_mg(k,l) * p_mg(k-1,j,i) - f_mg(k,j,i) & 516 ! ) 517 ! j = jj+6 518 ! p_mg(k,j,i) = 1.0 / f1_mg(k,l) * ( & 519 ! ddx2_mg(l) * ( p_mg(k,j,i+1) + p_mg(k,j,i-1) ) & 520 ! + ddy2_mg(l) * ( p_mg(k,j+1,i) + p_mg(k,j-1,i) ) & 521 ! + f2_mg(k,l) * p_mg(k+1,j,i) & 522 ! + f3_mg(k,l) * p_mg(k-1,j,i) - f_mg(k,j,i) & 523 ! ) 710 ddx2_mg(l) * & 711 ( p_mg(k,j,i+1) + IBITS( flags(k,j,i), 5, 1 ) * & 712 ( p_mg(k,j,i) - p_mg(k,j,i+1) ) + & 713 p_mg(k,j,i-1) + IBITS( flags(k,j,i), 4, 1 ) * & 714 ( p_mg(k,j,i) - p_mg(k,j,i-1) ) ) & 715 + ddy2_mg(l) * & 716 ( p_mg(k,j+1,i) + IBITS( flags(k,j,i), 3, 1 ) * & 717 ( p_mg(k,j,i) - p_mg(k,j+1,i) ) + & 718 p_mg(k,j-1,i) + IBITS( flags(k,j,i), 2, 1 ) * & 719 ( p_mg(k,j,i) - p_mg(k,j-1,i) ) ) & 720 + f2_mg(k,l) * p_mg(k+1,j,i) & 721 + f3_mg(k,l) * & 722 ( p_mg(k-1,j,i) + IBITS( flags(k,j,i), 0, 1 ) * & 723 ( p_mg(k,j,i) - p_mg(k-1,j,i) ) ) & 724 - f_mg(k,j,i) ) 524 725 ENDDO 525 726 … … 529 730 j =jj 530 731 p_mg(k,j,i) = 1.0 / f1_mg(k,l) * ( & 531 ddx2_mg(l) * ( p_mg(k,j,i+1) + p_mg(k,j,i-1) ) & 532 + ddy2_mg(l) * ( p_mg(k,j+1,i) + p_mg(k,j-1,i) ) & 533 + f2_mg(k,l) * p_mg(k+1,j,i) & 534 + f3_mg(k,l) * p_mg(k-1,j,i) - f_mg(k,j,i) & 535 ) 732 ddx2_mg(l) * & 733 ( p_mg(k,j,i+1) + IBITS( flags(k,j,i), 5, 1 ) * & 734 ( p_mg(k,j,i) - p_mg(k,j,i+1) ) + & 735 p_mg(k,j,i-1) + IBITS( flags(k,j,i), 4, 1 ) * & 736 ( p_mg(k,j,i) - p_mg(k,j,i-1) ) ) & 737 + ddy2_mg(l) * & 738 ( p_mg(k,j+1,i) + IBITS( flags(k,j,i), 3, 1 ) * & 739 ( p_mg(k,j,i) - p_mg(k,j+1,i) ) + & 740 p_mg(k,j-1,i) + IBITS( flags(k,j,i), 2, 1 ) * & 741 ( p_mg(k,j,i) - p_mg(k,j-1,i) ) ) & 742 + f2_mg(k,l) * p_mg(k+1,j,i) & 743 + f3_mg(k,l) * & 744 ( p_mg(k-1,j,i) + IBITS( flags(k,j,i), 0, 1 ) * & 745 ( p_mg(k,j,i) - p_mg(k-1,j,i) ) ) & 746 - f_mg(k,j,i) ) 536 747 j = jj+2 537 748 p_mg(k,j,i) = 1.0 / f1_mg(k,l) * ( & 538 ddx2_mg(l) * ( p_mg(k,j,i+1) + p_mg(k,j,i-1) ) & 539 + ddy2_mg(l) * ( p_mg(k,j+1,i) + p_mg(k,j-1,i) ) & 540 + f2_mg(k,l) * p_mg(k+1,j,i) & 541 + f3_mg(k,l) * p_mg(k-1,j,i) - f_mg(k,j,i) & 542 ) 543 ! j = jj+4 544 ! p_mg(k,j,i) = 1.0 / f1_mg(k,l) * ( & 545 ! ddx2_mg(l) * ( p_mg(k,j,i+1) + p_mg(k,j,i-1) ) & 546 ! + ddy2_mg(l) * ( p_mg(k,j+1,i) + p_mg(k,j-1,i) ) & 547 ! + f2_mg(k,l) * p_mg(k+1,j,i) & 548 ! + f3_mg(k,l) * p_mg(k-1,j,i) - f_mg(k,j,i) & 549 ! ) 550 ! j = jj+6 551 ! p_mg(k,j,i) = 1.0 / f1_mg(k,l) * ( & 552 ! ddx2_mg(l) * ( p_mg(k,j,i+1) + p_mg(k,j,i-1) ) & 553 ! + ddy2_mg(l) * ( p_mg(k,j+1,i) + p_mg(k,j-1,i) ) & 554 ! + f2_mg(k,l) * p_mg(k+1,j,i) & 555 ! + f3_mg(k,l) * p_mg(k-1,j,i) - f_mg(k,j,i) & 556 ! ) 749 ddx2_mg(l) * & 750 ( p_mg(k,j,i+1) + IBITS( flags(k,j,i), 5, 1 ) * & 751 ( p_mg(k,j,i) - p_mg(k,j,i+1) ) + & 752 p_mg(k,j,i-1) + IBITS( flags(k,j,i), 4, 1 ) * & 753 ( p_mg(k,j,i) - p_mg(k,j,i-1) ) ) & 754 + ddy2_mg(l) * & 755 ( p_mg(k,j+1,i) + IBITS( flags(k,j,i), 3, 1 ) * & 756 ( p_mg(k,j,i) - p_mg(k,j+1,i) ) + & 757 p_mg(k,j-1,i) + IBITS( flags(k,j,i), 2, 1 ) * & 758 ( p_mg(k,j,i) - p_mg(k,j-1,i) ) ) & 759 + f2_mg(k,l) * p_mg(k+1,j,i) & 760 + f3_mg(k,l) * & 761 ( p_mg(k-1,j,i) + IBITS( flags(k,j,i), 0, 1 ) * & 762 ( p_mg(k,j,i) - p_mg(k-1,j,i) ) ) & 763 - f_mg(k,j,i) ) 557 764 ENDDO 558 765 … … 562 769 j =jj 563 770 p_mg(k,j,i) = 1.0 / f1_mg(k,l) * ( & 564 ddx2_mg(l) * ( p_mg(k,j,i+1) + p_mg(k,j,i-1) ) & 565 + ddy2_mg(l) * ( p_mg(k,j+1,i) + p_mg(k,j-1,i) ) & 566 + f2_mg(k,l) * p_mg(k+1,j,i) & 567 + f3_mg(k,l) * p_mg(k-1,j,i) - f_mg(k,j,i) & 568 ) 771 ddx2_mg(l) * & 772 ( p_mg(k,j,i+1) + IBITS( flags(k,j,i), 5, 1 ) * & 773 ( p_mg(k,j,i) - p_mg(k,j,i+1) ) + & 774 p_mg(k,j,i-1) + IBITS( flags(k,j,i), 4, 1 ) * & 775 ( p_mg(k,j,i) - p_mg(k,j,i-1) ) ) & 776 + ddy2_mg(l) * & 777 ( p_mg(k,j+1,i) + IBITS( flags(k,j,i), 3, 1 ) * & 778 ( p_mg(k,j,i) - p_mg(k,j+1,i) ) + & 779 p_mg(k,j-1,i) + IBITS( flags(k,j,i), 2, 1 ) * & 780 ( p_mg(k,j,i) - p_mg(k,j-1,i) ) ) & 781 + f2_mg(k,l) * p_mg(k+1,j,i) & 782 + f3_mg(k,l) * & 783 ( p_mg(k-1,j,i) + IBITS( flags(k,j,i), 0, 1 ) * & 784 ( p_mg(k,j,i) - p_mg(k-1,j,i) ) ) & 785 - f_mg(k,j,i) ) 569 786 j = jj+2 570 787 p_mg(k,j,i) = 1.0 / f1_mg(k,l) * ( & 571 ddx2_mg(l) * ( p_mg(k,j,i+1) + p_mg(k,j,i-1) ) & 572 + ddy2_mg(l) * ( p_mg(k,j+1,i) + p_mg(k,j-1,i) ) & 573 + f2_mg(k,l) * p_mg(k+1,j,i) & 574 + f3_mg(k,l) * p_mg(k-1,j,i) - f_mg(k,j,i) & 575 ) 576 ! j = jj+4 577 ! p_mg(k,j,i) = 1.0 / f1_mg(k,l) * ( & 578 ! ddx2_mg(l) * ( p_mg(k,j,i+1) + p_mg(k,j,i-1) ) & 579 ! + ddy2_mg(l) * ( p_mg(k,j+1,i) + p_mg(k,j-1,i) ) & 580 ! + f2_mg(k,l) * p_mg(k+1,j,i) & 581 ! + f3_mg(k,l) * p_mg(k-1,j,i) - f_mg(k,j,i) & 582 ! ) 583 ! j = jj+6 584 ! p_mg(k,j,i) = 1.0 / f1_mg(k,l) * ( & 585 ! ddx2_mg(l) * ( p_mg(k,j,i+1) + p_mg(k,j,i-1) ) & 586 ! + ddy2_mg(l) * ( p_mg(k,j+1,i) + p_mg(k,j-1,i) ) & 587 ! + f2_mg(k,l) * p_mg(k+1,j,i) & 588 ! + f3_mg(k,l) * p_mg(k-1,j,i) - f_mg(k,j,i) & 589 ! ) 788 ddx2_mg(l) * & 789 ( p_mg(k,j,i+1) + IBITS( flags(k,j,i), 5, 1 ) * & 790 ( p_mg(k,j,i) - p_mg(k,j,i+1) ) + & 791 p_mg(k,j,i-1) + IBITS( flags(k,j,i), 4, 1 ) * & 792 ( p_mg(k,j,i) - p_mg(k,j,i-1) ) ) & 793 + ddy2_mg(l) * & 794 ( p_mg(k,j+1,i) + IBITS( flags(k,j,i), 3, 1 ) * & 795 ( p_mg(k,j,i) - p_mg(k,j+1,i) ) + & 796 p_mg(k,j-1,i) + IBITS( flags(k,j,i), 2, 1 ) * & 797 ( p_mg(k,j,i) - p_mg(k,j-1,i) ) ) & 798 + f2_mg(k,l) * p_mg(k+1,j,i) & 799 + f3_mg(k,l) * & 800 ( p_mg(k-1,j,i) + IBITS( flags(k,j,i), 0, 1 ) * & 801 ( p_mg(k,j,i) - p_mg(k-1,j,i) ) ) & 802 - f_mg(k,j,i) ) 590 803 ENDDO 591 804 … … 595 808 j =jj 596 809 p_mg(k,j,i) = 1.0 / f1_mg(k,l) * ( & 597 ddx2_mg(l) * ( p_mg(k,j,i+1) + p_mg(k,j,i-1) ) & 598 + ddy2_mg(l) * ( p_mg(k,j+1,i) + p_mg(k,j-1,i) ) & 599 + f2_mg(k,l) * p_mg(k+1,j,i) & 600 + f3_mg(k,l) * p_mg(k-1,j,i) - f_mg(k,j,i) & 601 ) 810 ddx2_mg(l) * & 811 ( p_mg(k,j,i+1) + IBITS( flags(k,j,i), 5, 1 ) * & 812 ( p_mg(k,j,i) - p_mg(k,j,i+1) ) + & 813 p_mg(k,j,i-1) + IBITS( flags(k,j,i), 4, 1 ) * & 814 ( p_mg(k,j,i) - p_mg(k,j,i-1) ) ) & 815 + ddy2_mg(l) * & 816 ( p_mg(k,j+1,i) + IBITS( flags(k,j,i), 3, 1 ) * & 817 ( p_mg(k,j,i) - p_mg(k,j+1,i) ) + & 818 p_mg(k,j-1,i) + IBITS( flags(k,j,i), 2, 1 ) * & 819 ( p_mg(k,j,i) - p_mg(k,j-1,i) ) ) & 820 + f2_mg(k,l) * p_mg(k+1,j,i) & 821 + f3_mg(k,l) * & 822 ( p_mg(k-1,j,i) + IBITS( flags(k,j,i), 0, 1 ) * & 823 ( p_mg(k,j,i) - p_mg(k-1,j,i) ) ) & 824 - f_mg(k,j,i) ) 602 825 j = jj+2 603 826 p_mg(k,j,i) = 1.0 / f1_mg(k,l) * ( & 604 ddx2_mg(l) * ( p_mg(k,j,i+1) + p_mg(k,j,i-1) ) & 605 + ddy2_mg(l) * ( p_mg(k,j+1,i) + p_mg(k,j-1,i) ) & 606 + f2_mg(k,l) * p_mg(k+1,j,i) & 607 + f3_mg(k,l) * p_mg(k-1,j,i) - f_mg(k,j,i) & 608 ) 609 ! j = jj+4 610 ! p_mg(k,j,i) = 1.0 / f1_mg(k,l) * ( & 611 ! ddx2_mg(l) * ( p_mg(k,j,i+1) + p_mg(k,j,i-1) ) & 612 ! + ddy2_mg(l) * ( p_mg(k,j+1,i) + p_mg(k,j-1,i) ) & 613 ! + f2_mg(k,l) * p_mg(k+1,j,i) & 614 ! + f3_mg(k,l) * p_mg(k-1,j,i) - f_mg(k,j,i) & 615 ! ) 616 ! j = jj+6 617 ! p_mg(k,j,i) = 1.0 / f1_mg(k,l) * ( & 618 ! ddx2_mg(l) * ( p_mg(k,j,i+1) + p_mg(k,j,i-1) ) & 619 ! + ddy2_mg(l) * ( p_mg(k,j+1,i) + p_mg(k,j-1,i) ) & 620 ! + f2_mg(k,l) * p_mg(k+1,j,i) & 621 ! + f3_mg(k,l) * p_mg(k-1,j,i) - f_mg(k,j,i) & 622 ! ) 827 ddx2_mg(l) * & 828 ( p_mg(k,j,i+1) + IBITS( flags(k,j,i), 5, 1 ) * & 829 ( p_mg(k,j,i) - p_mg(k,j,i+1) ) + & 830 p_mg(k,j,i-1) + IBITS( flags(k,j,i), 4, 1 ) * & 831 ( p_mg(k,j,i) - p_mg(k,j,i-1) ) ) & 832 + ddy2_mg(l) * & 833 ( p_mg(k,j+1,i) + IBITS( flags(k,j,i), 3, 1 ) * & 834 ( p_mg(k,j,i) - p_mg(k,j+1,i) ) + & 835 p_mg(k,j-1,i) + IBITS( flags(k,j,i), 2, 1 ) * & 836 ( p_mg(k,j,i) - p_mg(k,j-1,i) ) ) & 837 + f2_mg(k,l) * p_mg(k+1,j,i) & 838 + f3_mg(k,l) * & 839 ( p_mg(k-1,j,i) + IBITS( flags(k,j,i), 0, 1 ) * & 840 ( p_mg(k,j,i) - p_mg(k-1,j,i) ) ) & 841 - f_mg(k,j,i) ) 623 842 ENDDO 624 843 … … 669 888 ENDDO 670 889 890 ! 891 !-- Set pressure within topography and at the topography surfaces 892 !$OMP PARALLEL PRIVATE (i,j,k,wall_left,wall_north,wall_right,wall_south,wall_top,wall_total) 893 !$OMP DO 894 DO i = nxl_mg(l), nxr_mg(l) 895 DO j = nys_mg(l), nyn_mg(l) 896 DO k = nzb, nzt_mg(l) 897 ! 898 !-- First, set pressure inside topography to zero 899 p_mg(k,j,i) = p_mg(k,j,i) * ( 1.0 - IBITS( flags(k,j,i), 6, 1 ) ) 900 ! 901 !-- Second, determine if the gridpoint inside topography is adjacent 902 !-- to a wall and set its value to a value given by the average of 903 !-- those values obtained from Neumann boundary condition 904 wall_left = IBITS( flags(k,j,i-1), 5, 1 ) 905 wall_right = IBITS( flags(k,j,i+1), 4, 1 ) 906 wall_south = IBITS( flags(k,j-1,i), 3, 1 ) 907 wall_north = IBITS( flags(k,j+1,i), 2, 1 ) 908 wall_top = IBITS( flags(k+1,j,i), 0, 1 ) 909 wall_total = wall_left + wall_right + wall_south + wall_north + & 910 wall_top 911 912 IF ( wall_total > 0.0 ) THEN 913 p_mg(k,j,i) = 1.0 / wall_total * & 914 ( wall_left * p_mg(k,j,i-1) + & 915 wall_right * p_mg(k,j,i+1) + & 916 wall_south * p_mg(k,j-1,i) + & 917 wall_north * p_mg(k,j+1,i) + & 918 wall_top * p_mg(k+1,j,i) ) 919 ENDIF 920 ENDDO 921 ENDDO 922 ENDDO 923 !$OMP END PARALLEL 924 925 ! 926 !-- One more time horizontal boundary conditions 927 CALL exchange_horiz( p_mg ) 671 928 672 929 END SUBROUTINE redblack -
palm/trunk/SOURCE/user_interface.f90
r110 r114 4 4 ! Actual revisions: 5 5 ! ----------------- 6 ! 6 ! +argument gls in user_init_grid 7 7 ! 8 8 ! Former revisions: … … 238 238 239 239 240 SUBROUTINE user_init_grid( nzb_local )240 SUBROUTINE user_init_grid( gls, nzb_local ) 241 241 242 242 !------------------------------------------------------------------------------! … … 245 245 ! ------------ 246 246 ! Execution of user-defined grid initializing actions 247 ! First argument gls contains the number of ghost layers, which is > 1 if the 248 ! multigrid method for the pressure solver is used 247 249 !------------------------------------------------------------------------------! 248 250 … … 253 255 IMPLICIT NONE 254 256 255 INTEGER, DIMENSION(-1:ny+1,-1:nx+1) :: nzb_local 257 INTEGER :: gls 258 259 INTEGER, DIMENSION(-gls:ny+gls,-gls:nx+gls) :: nzb_local 256 260 257 261 ! -
palm/trunk/UTIL/combine_plot_fields.f90
- Property svn:keywords set to Id
r108 r114 4 4 ! Actual revisions: 5 5 ! ----------------- 6 ! loop for processing of output by coupled runs, id_string does not contain 7 ! modus any longer 6 ! Bugfix: model_string needed a default value 8 7 ! 9 8 ! Former revisions: 10 9 ! ----------------- 10 ! $Id$ 11 ! 12 ! Aug 07 Loop for processing of output by coupled runs, id_string does not 13 ! contain modus any longer 14 ! 11 15 ! 18/01/06 Output of time-averaged data 12 16 ! … … 106 110 107 111 PRINT*, '' 112 PRINT*, '' 108 113 PRINT*, '*** combine_plot_fields ***' 114 115 ! 116 !-- Find out if a coupled run has been carried out 109 117 INQUIRE( FILE='COUPLING_PORT_OPENED', EXIST=found ) 110 118 IF ( found ) THEN … … 115 123 PRINT*, ' uncoupled run' 116 124 ENDIF 125 126 ! 127 !-- Do everything for each model 117 128 DO model = 1, models 129 ! 130 !-- Set the model string used to identify the filenames 131 model_string = '' 118 132 IF ( models == 2 ) THEN 119 PRINT*, ''120 133 PRINT*, '' 121 134 PRINT*, '*** combine_plot_fields ***' 122 135 IF ( model == 2 ) THEN 123 model_string = "_O"136 model_string = '_O' 124 137 PRINT*, ' now combining ocean data' 125 138 PRINT*, ' ========================' 126 139 ELSE 127 model_string = ""128 140 PRINT*, ' now combining atmosphere data' 129 141 PRINT*, ' ============================='
Note: See TracChangeset
for help on using the changeset viewer.
