Changeset 63 for palm/trunk/SOURCE
- Timestamp:
- Mar 13, 2007 3:52:49 AM (18 years ago)
- Location:
- palm/trunk/SOURCE
- Files:
-
- 12 edited
-
CURRENT_MODIFICATIONS (modified) (1 diff)
-
advec_particles.f90 (modified) (8 diffs)
-
advec_s_bc.f90 (modified) (43 diffs)
-
check_parameters.f90 (modified) (3 diffs)
-
header.f90 (modified) (3 diffs)
-
init_3d_model.f90 (modified) (1 diff)
-
modules.f90 (modified) (2 diffs)
-
parin.f90 (modified) (2 diffs)
-
prognostic_equations.f90 (modified) (9 diffs)
-
read_var_list.f90 (modified) (2 diffs)
-
time_integration.f90 (modified) (4 diffs)
-
write_var_list.f90 (modified) (2 diffs)
Legend:
- Unmodified
- Added
- Removed
-
palm/trunk/SOURCE/CURRENT_MODIFICATIONS
r61 r63 9 9 10 10 new d3par-parameter dt_max to define the maximum value for the allowed timestep 11 12 new inipar-parameter loop_optimization to control the loop optimization method 11 13 12 14 new inipar-parameter pt_refrence. If given, this value is used as the reference that in buoyancy terms (otherwise, the instantaneous horizontally averaged temperature is used). -
palm/trunk/SOURCE/advec_particles.f90
r60 r63 2315 2315 ENDIF 2316 2316 2317 CALL MPI_SENDRECV( trlp , trlp_count, mpi_particle_type, pleft, 1,&2317 CALL MPI_SENDRECV( trlp(1), trlp_count, mpi_particle_type, pleft, 1, & 2318 2318 particles(number_of_particles+1), trrp_count_recv,& 2319 2319 mpi_particle_type, pright, 1, & … … 2343 2343 ENDIF 2344 2344 2345 CALL MPI_SENDRECV( trlpt, trlpt_count*tlength, MPI_REAL, pleft, 1,& 2345 CALL MPI_SENDRECV( trlpt(1,1,1), trlpt_count*tlength, MPI_REAL, & 2346 pleft, 1, & 2346 2347 particle_tail_coordinates(1,1,number_of_tails+1), & 2347 2348 trrpt_count_recv*tlength, MPI_REAL, pright, 1, & … … 2391 2392 ENDIF 2392 2393 2393 CALL MPI_SENDRECV( trrp , trrp_count, mpi_particle_type, pright, 1,&2394 CALL MPI_SENDRECV( trrp(1), trrp_count, mpi_particle_type, pright, 1,& 2394 2395 particles(number_of_particles+1), trlp_count_recv,& 2395 2396 mpi_particle_type, pleft, 1, & … … 2419 2420 ENDIF 2420 2421 2421 CALL MPI_SENDRECV( trrpt, trrpt_count*tlength, MPI_REAL, pright, & 2422 1, particle_tail_coordinates(1,1,number_of_tails+1), & 2422 CALL MPI_SENDRECV( trrpt(1,1,1), trrpt_count*tlength, MPI_REAL, & 2423 pright, 1, & 2424 particle_tail_coordinates(1,1,number_of_tails+1), & 2423 2425 trlpt_count_recv*tlength, MPI_REAL, pleft, 1, & 2424 2426 comm2d, status, ierr ) … … 2754 2756 ENDIF 2755 2757 2756 CALL MPI_SENDRECV( trsp , trsp_count, mpi_particle_type, psouth, 1,&2758 CALL MPI_SENDRECV( trsp(1), trsp_count, mpi_particle_type, psouth, 1,& 2757 2759 particles(number_of_particles+1), trnp_count_recv,& 2758 2760 mpi_particle_type, pnorth, 1, & … … 2782 2784 ENDIF 2783 2785 2784 CALL MPI_SENDRECV( trspt, trspt_count*tlength, MPI_REAL, psouth, & 2785 1, particle_tail_coordinates(1,1,number_of_tails+1), & 2786 CALL MPI_SENDRECV( trspt(1,1,1), trspt_count*tlength, MPI_REAL, & 2787 psouth, 1, & 2788 particle_tail_coordinates(1,1,number_of_tails+1), & 2786 2789 trnpt_count_recv*tlength, MPI_REAL, pnorth, 1, & 2787 2790 comm2d, status, ierr ) … … 2830 2833 ENDIF 2831 2834 2832 CALL MPI_SENDRECV( trnp , trnp_count, mpi_particle_type, pnorth, 1,&2835 CALL MPI_SENDRECV( trnp(1), trnp_count, mpi_particle_type, pnorth, 1,& 2833 2836 particles(number_of_particles+1), trsp_count_recv,& 2834 2837 mpi_particle_type, psouth, 1, & … … 2858 2861 ENDIF 2859 2862 2860 CALL MPI_SENDRECV( trnpt, trnpt_count*tlength, MPI_REAL, pnorth, & 2861 1, particle_tail_coordinates(1,1,number_of_tails+1), & 2863 CALL MPI_SENDRECV( trnpt(1,1,1), trnpt_count*tlength, MPI_REAL, & 2864 pnorth, 1, & 2865 particle_tail_coordinates(1,1,number_of_tails+1), & 2862 2866 trspt_count_recv*tlength, MPI_REAL, psouth, 1, & 2863 2867 comm2d, status, ierr ) -
palm/trunk/SOURCE/advec_s_bc.f90
r4 r63 1 1 SUBROUTINE advec_s_bc( sk, sk_char ) 2 2 3 !------------------------------------------------------------------------------ -!3 !------------------------------------------------------------------------------! 4 4 ! Actual revisions: 5 5 ! ----------------- 6 ! 6 ! Calculation extended for gridpoint nzt 7 7 ! 8 8 ! Former revisions: … … 29 29 ! stretched grids. 30 30 ! NOTE: This is a provisional, non-optimised version! 31 !------------------------------------------------------------------------------ -!31 !------------------------------------------------------------------------------! 32 32 33 33 USE advection 34 34 USE arrays_3d 35 USE control_parameters 35 36 USE cpulog 36 37 USE grid_variables … … 39 40 USE pegrid 40 41 USE statistics 41 USE control_parameters42 42 43 43 IMPLICIT NONE … … 57 57 REAL, DIMENSION(:,:,:), ALLOCATABLE :: sk_p 58 58 59 #if defined( __ t3eb ) || defined( __t3eh ) || defined( __t3ej2 ) || defined( __t3ej5)59 #if defined( __nec ) 60 60 REAL (kind=4) :: m1n, m1z !Wichtig: Division 61 61 REAL (kind=4), DIMENSION(:,:), ALLOCATABLE :: m1, sw … … 68 68 ! 69 69 !-- Array sk_p requires 2 extra elements for each dimension 70 ALLOCATE( sk_p(nzb-2:nzt+ 2,nys-3:nyn+3,nxl-3:nxr+3) )70 ALLOCATE( sk_p(nzb-2:nzt+3,nys-3:nyn+3,nxl-3:nxr+3) ) 71 71 sk_p = 0.0 72 72 … … 89 89 DO i = nxl-1, nxr+1 90 90 DO j = nys, nyn 91 DO k = nzb +1, nzt-191 DO k = nzb, nzt+1 92 92 sk_p(k,j,i) = sk(k,j,i) 93 93 ENDDO … … 95 95 ENDDO 96 96 #if defined( __parallel ) 97 ngp = 2 * ( nzt - nzb + 5) * ( nyn - nys + 7 )97 ngp = 2 * ( nzt - nzb + 6 ) * ( nyn - nys + 7 ) 98 98 CALL cpu_log( log_point_s(11), 'advec_s_bc:sendrecv', 'start' ) 99 99 ! … … 142 142 DO i = nxl, nxr 143 143 DO j = nys, nyn 144 DO k = nzb+1, nzt -1144 DO k = nzb+1, nzt 145 145 zaehler = ABS( sk_p(k,j,i+1) - 2.0 * sk_p(k,j,i) + sk_p(k,j,i-1) ) 146 146 nenner = ABS( sk_p(k,j,i+1) - sk_p(k,j,i-1) ) … … 160 160 ! 161 161 !-- Allocate temporary arrays 162 ALLOCATE( a0(nzb+1:nzt -1,nxl-1:nxr+1), a1(nzb+1:nzt-1,nxl-1:nxr+1), &163 a2(nzb+1:nzt -1,nxl-1:nxr+1), a12(nzb+1:nzt-1,nxl-1:nxr+1), &164 a22(nzb+1:nzt -1,nxl-1:nxr+1), immb(nzb+1:nzt-1,nxl-1:nxr+1), &165 imme(nzb+1:nzt -1,nxl-1:nxr+1), impb(nzb+1:nzt-1,nxl-1:nxr+1), &166 impe(nzb+1:nzt -1,nxl-1:nxr+1), ipmb(nzb+1:nzt-1,nxl-1:nxr+1), &167 ipme(nzb+1:nzt -1,nxl-1:nxr+1), ippb(nzb+1:nzt-1,nxl-1:nxr+1), &168 ippe(nzb+1:nzt -1,nxl-1:nxr+1), m1(nzb+1:nzt-1,nxl-2:nxr+2), &169 sw(nzb+1:nzt -1,nxl-1:nxr+1)&162 ALLOCATE( a0(nzb+1:nzt,nxl-1:nxr+1), a1(nzb+1:nzt,nxl-1:nxr+1), & 163 a2(nzb+1:nzt,nxl-1:nxr+1), a12(nzb+1:nzt,nxl-1:nxr+1), & 164 a22(nzb+1:nzt,nxl-1:nxr+1), immb(nzb+1:nzt,nxl-1:nxr+1), & 165 imme(nzb+1:nzt,nxl-1:nxr+1), impb(nzb+1:nzt,nxl-1:nxr+1), & 166 impe(nzb+1:nzt,nxl-1:nxr+1), ipmb(nzb+1:nzt,nxl-1:nxr+1), & 167 ipme(nzb+1:nzt,nxl-1:nxr+1), ippb(nzb+1:nzt,nxl-1:nxr+1), & 168 ippe(nzb+1:nzt,nxl-1:nxr+1), m1(nzb+1:nzt,nxl-2:nxr+2), & 169 sw(nzb+1:nzt,nxl-1:nxr+1) & 170 170 ) 171 171 imme = 0.0; impe = 0.0; ipme = 0.0; ippe = 0.0 … … 184 184 !-- Compute polynomial coefficients 185 185 DO i = nxl-1, nxr+1 186 DO k = nzb+1, nzt -1186 DO k = nzb+1, nzt 187 187 a12(k,i) = 0.5 * ( sk_p(k,j,i+1) - sk_p(k,j,i-1) ) 188 188 a22(k,i) = 0.5 * ( sk_p(k,j,i+1) - 2.0 * sk_p(k,j,i) & … … 204 204 !-- *VOCL LOOP,UNROLL(2) 205 205 DO i = nxl, nxr 206 DO k = nzb+1, nzt -1206 DO k = nzb+1, nzt 207 207 cip = MAX( 0.0, ( u(k,j,i+1) - u_gtrans ) * dt_3d * ddx ) 208 208 cim = -MIN( 0.0, ( u(k,j,i+1) - u_gtrans ) * dt_3d * ddx ) … … 240 240 !-- Compute monitor function m1 241 241 DO i = nxl-2, nxr+2 242 DO k = nzb+1, nzt -1242 DO k = nzb+1, nzt 243 243 m1z = ABS( sk_p(k,j,i+1) - 2.0 * sk_p(k,j,i) + sk_p(k,j,i-1) ) 244 244 m1n = ABS( sk_p(k,j,i+1) - sk_p(k,j,i-1) ) … … 258 258 sw = 0.0 259 259 DO i = nxl-1, nxr+1 260 DO k = nzb+1, nzt -1260 DO k = nzb+1, nzt 261 261 m2 = 2.0 * ABS( a1(k,i) - a12(k,i) ) / & 262 262 MAX( ABS( a1(k,i) + a12(k,i) ), 1E-35 ) … … 284 284 CALL cpu_log( log_point_s(12), 'advec_s_bc:exp', 'continue' ) 285 285 DO i = nxl, nxr 286 DO k = nzb+1, nzt -1286 DO k = nzb+1, nzt 287 287 288 288 !-- *VOCL STMT,IF(10) … … 374 374 !-- Prognostic equation 375 375 DO i = nxl, nxr 376 DO k = nzb+1, nzt -1376 DO k = nzb+1, nzt 377 377 fplus = ( 1.0 - sw(k,i) ) * ippb(k,i) + sw(k,i) * ippe(k,i) & 378 378 - ( 1.0 - sw(k,i+1) ) * impb(k,i) - sw(k,i+1) * impe(k,i) … … 404 404 !-- Enlarge boundary of local array cyclically in y-direction 405 405 #if defined( __parallel ) 406 ngp = ( nzt - nzb + 5) * ( nyn - nys + 7 )407 CALL MPI_TYPE_VECTOR( nxr-nxl+7, 3*(nzt-nzb+ 5), ngp, MPI_REAL, &406 ngp = ( nzt - nzb + 6 ) * ( nyn - nys + 7 ) 407 CALL MPI_TYPE_VECTOR( nxr-nxl+7, 3*(nzt-nzb+6), ngp, MPI_REAL, & 408 408 type_xz_2, ierr ) 409 409 CALL MPI_TYPE_COMMIT( type_xz_2, ierr ) … … 423 423 #else 424 424 DO i = nxl, nxr 425 DO k = nzb+1, nzt -1425 DO k = nzb+1, nzt 426 426 sk_p(k,nys-1,i) = sk_p(k,nyn,i) 427 427 sk_p(k,nys-2,i) = sk_p(k,nyn-1,i) … … 439 439 DO i = nxl, nxr 440 440 DO j = nys, nyn 441 DO k = nzb+1, nzt -1441 DO k = nzb+1, nzt 442 442 zaehler = ABS( sk_p(k,j+1,i) - 2.0 * sk_p(k,j,i) + sk_p(k,j-1,i) ) 443 443 nenner = ABS( sk_p(k,j+1,i) - sk_p(k,j-1,i) ) … … 457 457 ! 458 458 !-- Allocate temporary arrays 459 ALLOCATE( a0(nzb+1:nzt -1,nys-1:nyn+1), a1(nzb+1:nzt-1,nys-1:nyn+1), &460 a2(nzb+1:nzt -1,nys-1:nyn+1), a12(nzb+1:nzt-1,nys-1:nyn+1), &461 a22(nzb+1:nzt -1,nys-1:nyn+1), immb(nzb+1:nzt-1,nys-1:nyn+1), &462 imme(nzb+1:nzt -1,nys-1:nyn+1), impb(nzb+1:nzt-1,nys-1:nyn+1), &463 impe(nzb+1:nzt -1,nys-1:nyn+1), ipmb(nzb+1:nzt-1,nys-1:nyn+1), &464 ipme(nzb+1:nzt -1,nys-1:nyn+1), ippb(nzb+1:nzt-1,nys-1:nyn+1), &465 ippe(nzb+1:nzt -1,nys-1:nyn+1), m1(nzb+1:nzt-1,nys-2:nyn+2), &466 sw(nzb+1:nzt -1,nys-1:nyn+1)&459 ALLOCATE( a0(nzb+1:nzt,nys-1:nyn+1), a1(nzb+1:nzt,nys-1:nyn+1), & 460 a2(nzb+1:nzt,nys-1:nyn+1), a12(nzb+1:nzt,nys-1:nyn+1), & 461 a22(nzb+1:nzt,nys-1:nyn+1), immb(nzb+1:nzt,nys-1:nyn+1), & 462 imme(nzb+1:nzt,nys-1:nyn+1), impb(nzb+1:nzt,nys-1:nyn+1), & 463 impe(nzb+1:nzt,nys-1:nyn+1), ipmb(nzb+1:nzt,nys-1:nyn+1), & 464 ipme(nzb+1:nzt,nys-1:nyn+1), ippb(nzb+1:nzt,nys-1:nyn+1), & 465 ippe(nzb+1:nzt,nys-1:nyn+1), m1(nzb+1:nzt,nys-2:nyn+2), & 466 sw(nzb+1:nzt,nys-1:nyn+1) & 467 467 ) 468 468 imme = 0.0; impe = 0.0; ipme = 0.0; ippe = 0.0 … … 475 475 !-- Compute polynomial coefficients 476 476 DO j = nys-1, nyn+1 477 DO k = nzb+1, nzt -1477 DO k = nzb+1, nzt 478 478 a12(k,j) = 0.5 * ( sk_p(k,j+1,i) - sk_p(k,j-1,i) ) 479 479 a22(k,j) = 0.5 * ( sk_p(k,j+1,i) - 2.0 * sk_p(k,j,i) & … … 495 495 !-- *VOCL LOOP,UNROLL(2) 496 496 DO j = nys, nyn 497 DO k = nzb+1, nzt -1497 DO k = nzb+1, nzt 498 498 cip = MAX( 0.0, ( v(k,j+1,i) - v_gtrans ) * dt_3d * ddy ) 499 499 cim = -MIN( 0.0, ( v(k,j+1,i) - v_gtrans ) * dt_3d * ddy ) … … 531 531 !-- Compute monitor function m1 532 532 DO j = nys-2, nyn+2 533 DO k = nzb+1, nzt -1533 DO k = nzb+1, nzt 534 534 m1z = ABS( sk_p(k,j+1,i) - 2.0 * sk_p(k,j,i) + sk_p(k,j-1,i) ) 535 535 m1n = ABS( sk_p(k,j+1,i) - sk_p(k,j-1,i) ) … … 549 549 sw = 0.0 550 550 DO j = nys-1, nyn+1 551 DO k = nzb+1, nzt -1551 DO k = nzb+1, nzt 552 552 m2 = 2.0 * ABS( a1(k,j) - a12(k,j) ) / & 553 553 MAX( ABS( a1(k,j) + a12(k,j) ), 1E-35 ) … … 575 575 CALL cpu_log( log_point_s(12), 'advec_s_bc:exp', 'continue' ) 576 576 DO j = nys, nyn 577 DO k = nzb+1, nzt -1577 DO k = nzb+1, nzt 578 578 579 579 !-- *VOCL STMT,IF(10) … … 665 665 !-- Prognostic equation 666 666 DO j = nys, nyn 667 DO k = nzb+1, nzt -1667 DO k = nzb+1, nzt 668 668 fplus = ( 1.0 - sw(k,j) ) * ippb(k,j) + sw(k,j) * ippe(k,j) & 669 669 - ( 1.0 - sw(k,j+1) ) * impb(k,j) - sw(k,j+1) * impe(k,j) … … 720 720 DO i = nxl, nxr 721 721 DO j = nys, nyn 722 sk_p(nzb,j,i) = sk_p(nzb+1,j,i) 723 sk_p(nzb-1,j,i) = sk_p(nzb+1,j,i) 724 sk_p(nzb-2,j,i) = sk_p(nzb+1,j,i) 722 sk_p(nzb-1,j,i) = sk_p(nzb,j,i) 723 sk_p(nzb-2,j,i) = sk_p(nzb,j,i) 725 724 ENDDO 726 725 ENDDO … … 730 729 ! 731 730 !-- Temperature boundary condition at the top boundary 732 IF ( ibc_pt_t == 0 ) THEN733 ! 734 !-- Dirichlet 731 IF ( ibc_pt_t == 0 .OR. ibc_pt_t == 1 ) THEN 732 ! 733 !-- Dirichlet or Neumann (zero gradient) 735 734 DO i = nxl, nxr 736 735 DO j = nys, nyn 737 sk_p(nzt+ 1,j,i) = sk_p(nzt,j,i)738 sk_p(nzt+ 2,j,i) = sk_p(nzt,j,i)736 sk_p(nzt+2,j,i) = sk_p(nzt+1,j,i) 737 sk_p(nzt+3,j,i) = sk_p(nzt+1,j,i) 739 738 ENDDO 740 739 ENDDO 741 740 742 ELSE 743 ! 744 !-- Neumann: dzu(nzt+2 ) is not defined, thus instead dzu(nzt+1) is used741 ELSEIF ( ibc_pt_t == 2 ) THEN 742 ! 743 !-- Neumann: dzu(nzt+2:3) are not defined, dzu(nzt+1) is used instead 745 744 DO i = nxl, nxr 746 745 DO j = nys, nyn 747 sk_p(nzt,j,i) = sk_p(nzt-1,j,i) + bc_pt_t_val * dzu(nzt)748 sk_p(nzt+1,j,i) = sk_p(nzt,j,i) + bc_pt_t_val * dzu(nzt+1)749 746 sk_p(nzt+2,j,i) = sk_p(nzt+1,j,i) + bc_pt_t_val * dzu(nzt+1) 747 sk_p(nzt+3,j,i) = sk_p(nzt+2,j,i) + bc_pt_t_val * dzu(nzt+1) 750 748 ENDDO 751 749 ENDDO … … 756 754 757 755 ! 758 !-- Temperatureboundary condition at the bottom boundary756 !-- Specific humidity boundary condition at the bottom boundary 759 757 IF ( ibc_q_b == 0 ) THEN 760 758 ! 761 !-- Dirichlet (fixed surface temperature)759 !-- Dirichlet (fixed surface humidity) 762 760 DO i = nxl, nxr 763 761 DO j = nys, nyn … … 772 770 DO i = nxl, nxr 773 771 DO j = nys, nyn 774 sk_p(nzb,j,i) = sk_p(nzb+1,j,i) 775 sk_p(nzb-1,j,i) = sk_p(nzb+1,j,i) 776 sk_p(nzb-2,j,i) = sk_p(nzb+1,j,i) 772 sk_p(nzb-1,j,i) = sk_p(nzb,j,i) 773 sk_p(nzb-2,j,i) = sk_p(nzb,j,i) 777 774 ENDDO 778 775 ENDDO … … 781 778 782 779 ! 783 !-- Temperatureboundary condition at the top boundary780 !-- Specific humidity boundary condition at the top boundary 784 781 IF ( ibc_q_t == 0 ) THEN 785 782 ! … … 787 784 DO i = nxl, nxr 788 785 DO j = nys, nyn 789 sk_p(nzt+ 1,j,i) = sk_p(nzt,j,i)790 sk_p(nzt+ 2,j,i) = sk_p(nzt,j,i)786 sk_p(nzt+2,j,i) = sk_p(nzt+1,j,i) 787 sk_p(nzt+3,j,i) = sk_p(nzt+1,j,i) 791 788 ENDDO 792 789 ENDDO … … 794 791 ELSE 795 792 ! 796 !-- Neumann: dzu(nzt+2 ) is not defined, thus instead dzu(nzt+1) is used793 !-- Neumann: dzu(nzt+2:3) are not defined, dzu(nzt+1) is used instead 797 794 DO i = nxl, nxr 798 795 DO j = nys, nyn 799 sk_p(nzt,j,i) = sk_p(nzt-1,j,i) + bc_q_t_val * dzu(nzt)800 sk_p(nzt+1,j,i) = sk_p(nzt,j,i) + bc_q_t_val * dzu(nzt+1)801 796 sk_p(nzt+2,j,i) = sk_p(nzt+1,j,i) + bc_q_t_val * dzu(nzt+1) 797 sk_p(nzt+3,j,i) = sk_p(nzt+2,j,i) + bc_q_t_val * dzu(nzt+1) 802 798 ENDDO 803 799 ENDDO … … 809 805 ! 810 806 !-- TKE boundary condition at bottom and top boundary (generally Neumann) 811 sk_p(nzb,j,i) = sk_p(nzb+1,j,i) 812 sk_p(nzb-1,j,i) = sk_p(nzb+1,j,i) 813 sk_p(nzb-2,j,i) = sk_p(nzb+1,j,i) 814 sk_p(nzt,j,i) = sk_p(nzt-1,j,i) 815 sk_p(nzt+1,j,i) = sk_p(nzt-1,j,i) 816 sk_p(nzt+2,j,i) = sk_p(nzt-1,j,i) 807 sk_p(nzb-1,j,i) = sk_p(nzb,j,i) 808 sk_p(nzb-2,j,i) = sk_p(nzb,j,i) 809 sk_p(nzt+2,j,i) = sk_p(nzt+1,j,i) 810 sk_p(nzt+3,j,i) = sk_p(nzt+1,j,i) 817 811 818 812 ELSE … … 829 823 DO i = nxl, nxr 830 824 DO j = nys, nyn 831 DO k = nzb, nzt 825 DO k = nzb, nzt+1 832 826 zaehler = ABS( sk_p(k+1,j,i) - 2.0 * sk_p(k,j,i) + sk_p(k-1,j,i) ) 833 827 nenner = ABS( sk_p(k+1,j,i+1) - sk_p(k-1,j,i) ) … … 847 841 ! 848 842 !-- Allocate temporary arrays 849 ALLOCATE( a0(nzb:nzt ,nys:nyn), a1(nzb:nzt,nys:nyn),&850 a2(nzb:nzt ,nys:nyn), a12(nzb:nzt,nys:nyn),&851 a22(nzb:nzt ,nys:nyn), immb(nzb+1:nzt-1,nys:nyn), &852 imme(nzb+1:nzt -1,nys:nyn), impb(nzb+1:nzt-1,nys:nyn), &853 impe(nzb+1:nzt -1,nys:nyn), ipmb(nzb+1:nzt-1,nys:nyn), &854 ipme(nzb+1:nzt -1,nys:nyn), ippb(nzb+1:nzt-1,nys:nyn), &855 ippe(nzb+1:nzt -1,nys:nyn), m1(nzb-1:nzt+1,nys:nyn),&856 sw(nzb:nzt ,nys:nyn)&843 ALLOCATE( a0(nzb:nzt+1,nys:nyn), a1(nzb:nzt+1,nys:nyn), & 844 a2(nzb:nzt+1,nys:nyn), a12(nzb:nzt+1,nys:nyn), & 845 a22(nzb:nzt+1,nys:nyn), immb(nzb+1:nzt,nys:nyn), & 846 imme(nzb+1:nzt,nys:nyn), impb(nzb+1:nzt,nys:nyn), & 847 impe(nzb+1:nzt,nys:nyn), ipmb(nzb+1:nzt,nys:nyn), & 848 ipme(nzb+1:nzt,nys:nyn), ippb(nzb+1:nzt,nys:nyn), & 849 ippe(nzb+1:nzt,nys:nyn), m1(nzb-1:nzt+2,nys:nyn), & 850 sw(nzb:nzt+1,nys:nyn) & 857 851 ) 858 852 imme = 0.0; impe = 0.0; ipme = 0.0; ippe = 0.0 … … 865 859 !-- Compute polynomial coefficients 866 860 DO j = nys, nyn 867 DO k = nzb, nzt 861 DO k = nzb, nzt+1 868 862 a12(k,j) = 0.5 * ( sk_p(k+1,j,i) - sk_p(k-1,j,i) ) 869 863 a22(k,j) = 0.5 * ( sk_p(k+1,j,i) - 2.0 * sk_p(k,j,i) & … … 885 879 !-- *VOCL LOOP,UNROLL(2) 886 880 DO j = nys, nyn 887 DO k = nzb+1, nzt -1881 DO k = nzb+1, nzt 888 882 cip = MAX( 0.0, w(k,j,i) * dt_3d * ddzw(k) ) 889 883 cim = -MIN( 0.0, w(k,j,i) * dt_3d * ddzw(k) ) … … 921 915 !-- Compute monitor function m1 922 916 DO j = nys, nyn 923 DO k = nzb-1, nzt+ 1917 DO k = nzb-1, nzt+2 924 918 m1z = ABS( sk_p(k+1,j,i) - 2.0 * sk_p(k,j,i) + sk_p(k-1,j,i) ) 925 919 m1n = ABS( sk_p(k+1,j,i) - sk_p(k-1,j,i) ) … … 939 933 sw = 0.0 940 934 DO j = nys, nyn 941 DO k = nzb, nzt 935 DO k = nzb, nzt+1 942 936 m2 = 2.0 * ABS( a1(k,j) - a12(k,j) ) / & 943 937 MAX( ABS( a1(k,j) + a12(k,j) ), 1E-35 ) … … 965 959 CALL cpu_log( log_point_s(12), 'advec_s_bc:exp', 'continue' ) 966 960 DO j = nys, nyn 967 DO k = nzb+1, nzt -1961 DO k = nzb+1, nzt 968 962 969 963 !-- *VOCL STMT,IF(10) … … 1055 1049 !-- Prognostic equation 1056 1050 DO j = nys, nyn 1057 DO k = nzb+1, nzt -11051 DO k = nzb+1, nzt 1058 1052 fplus = ( 1.0 - sw(k,j) ) * ippb(k,j) + sw(k,j) * ippe(k,j) & 1059 1053 - ( 1.0 - sw(k+1,j) ) * impb(k,j) - sw(k+1,j) * impe(k,j) … … 1082 1076 DO sr = 0, statistic_regions 1083 1077 DO j = nys, nyn 1084 DO k = nzb+1, nzt -11078 DO k = nzb+1, nzt 1085 1079 sums_wsts_bc_l(k,sr) = sums_wsts_bc_l(k,sr) + & 1086 1080 m1(k,j) * rmask(j,i,sr) … … 1103 1097 DO i = nxl, nxr 1104 1098 DO j = nys, nyn 1105 DO k = nzb+1, nzt -11099 DO k = nzb+1, nzt 1106 1100 tend(k,j,i) = tend(k,j,i) + ( sk_p(k,j,i) - sk(k,j,i) ) / dt_3d 1107 1101 ENDDO -
palm/trunk/SOURCE/check_parameters.f90
r57 r63 5 5 ! ----------------- 6 6 ! "by_user" allowed as initializing action, -data_output_ts, 7 ! leapfrog with non-flat topography not allowed any more, pt_reference is8 ! checked7 ! leapfrog with non-flat topography not allowed any more, loop_optimization 8 ! and pt_reference are checked 9 9 ! 10 10 ! Former revisions: … … 82 82 83 83 ! 84 !-- Check the general loop optimization method 85 IF ( loop_optimization == 'default' ) THEN 86 IF ( host(1:3) == 'nec' ) THEN 87 loop_optimization = 'vector' 88 ELSE 89 loop_optimization = 'cache' 90 ENDIF 91 ENDIF 92 IF ( loop_optimization /= 'noopt' .AND. loop_optimization /= 'cache' & 93 .AND. loop_optimization /= 'vector' ) THEN 94 IF ( myid == 0 ) THEN 95 PRINT*, '+++ check_parameters:' 96 PRINT*, ' illegal value given for loop_optimization: ', & 97 TRIM( loop_optimization ) 98 ENDIF 99 CALL local_stop 100 ENDIF 101 102 ! 84 103 !-- Check topography setting (check for illegal parameter combinations) 85 104 IF ( topography /= 'flat' ) THEN … … 257 276 END SELECT 258 277 259 IF ( scalar_advec /= 'pw-scheme' .AND. timestep_scheme(1:5) == 'runge' )&278 IF ( scalar_advec == 'ups-scheme' .AND. timestep_scheme(1:5) == 'runge' )& 260 279 THEN 261 280 IF ( myid == 0 ) THEN -
palm/trunk/SOURCE/header.f90
r60 r63 4 4 ! Actual revisions: 5 5 ! ----------------- 6 ! Output of netcdf_64bit_3d, particles-package is now part of the default code. 6 ! Output of netcdf_64bit_3d, particles-package is now part of the default code, 7 ! output of the loop optimization method. 7 8 ! 8 9 ! Former revisions: … … 189 190 WRITE ( io, 118 ) 190 191 ENDIF 191 IF ( momentum_advec /= 'ups-scheme' .AND. scalar_advec /= 'ups-scheme' & 192 .AND. scalar_advec /= 'bc-scheme' .AND. host(1:3) /= 'nec' ) THEN 193 WRITE ( io, 139 ) 194 ENDIF 192 193 WRITE ( io, 139 ) TRIM( loop_optimization ) 194 195 195 IF ( galilei_transformation ) THEN 196 196 IF ( use_ug_for_galilei_tr ) THEN … … 1182 1182 I3,')') 1183 1183 138 FORMAT (' Using hybrid version for 1d-domain-decomposition') 1184 139 FORMAT (' --> Prognostic equations are solved in one single loop (fast', & 1185 ' method)') 1184 139 FORMAT (' --> Loop optimization method: ',A) 1186 1185 140 FORMAT (' maximum residual allowed: ',E10.3) 1187 1186 141 FORMAT (' fixed number of multigrid cycles: ',I4) -
palm/trunk/SOURCE/init_3d_model.f90
r51 r63 756 756 ! 757 757 !-- If required, initialize particles 758 CALL init_particles758 IF ( particle_advection ) CALL init_particles 759 759 760 760 ! -
palm/trunk/SOURCE/modules.f90
r57 r63 6 6 ! ----------------- 7 7 ! +array rif_wall 8 ! +netcdf_64bit_3d, zu_s_inner, zw_w_inner, id_var_zusi_*, id_var_zwwi_*, 9 ! ts_value, u_nzb_p1_for_vfc, v_nzb_p1_for_vfc, pt_reference, use_pt_reference, 8 ! +loop_optimization, netcdf_64bit_3d, zu_s_inner, zw_w_inner, id_var_zusi_*, 9 ! id_var_zwwi_*, ts_value, u_nzb_p1_for_vfc, v_nzb_p1_for_vfc, pt_reference, 10 ! use_pt_reference, 10 11 ! +age_m in particle_type 11 12 ! -data_output_ts, dots_n … … 211 212 CHARACTER (LEN=9) :: simulated_time_chr 212 213 CHARACTER (LEN=12) :: version = 'PALM 3.1c' 213 CHARACTER (LEN=16) :: momentum_advec = 'pw-scheme', & 214 CHARACTER (LEN=16) :: loop_optimization = 'default', & 215 momentum_advec = 'pw-scheme', & 214 216 psolver = 'poisfft', & 215 217 scalar_advec = 'pw-scheme' -
palm/trunk/SOURCE/parin.f90
r61 r63 4 4 ! Actual revisions: 5 5 ! ----------------- 6 ! +dt_max, netcdf_64bit_3d in d3par, +pt_reference in inipar, -data_output_ts 6 ! +dt_max, netcdf_64bit_3d in d3par, +loop_optimization, pt_reference in 7 ! inipar, -data_output_ts 7 8 ! 8 9 ! Former revisions: … … 60 61 inflow_disturbance_end, initializing_actions, & 61 62 km_constant, km_damp_max, long_filter_factor, & 62 mixing_length_1d, moisture, momentum_advec, & 63 netcdf_precision, npex, npey, nsor_ini, nx, ny, nz, & 64 omega, outflow_damping_width, overshoot_limit_e, & 63 loop_optimization, mixing_length_1d, moisture, & 64 momentum_advec, netcdf_precision, npex, npey, nsor_ini, & 65 nx, ny, nz, omega, outflow_damping_width, & 66 overshoot_limit_e, & 65 67 overshoot_limit_pt, overshoot_limit_u, & 66 68 overshoot_limit_v, overshoot_limit_w, passive_scalar, & -
palm/trunk/SOURCE/prognostic_equations.f90
r57 r63 4 4 ! Actual revisions: 5 5 ! ----------------- 6 ! z0 removed from arguments in calls of diffusion_u/v/w 6 ! z0 removed from arguments in calls of diffusion_u/v/w, 7 ! subroutine names changed to .._noopt, .._cache, and .._vector 7 8 ! 8 9 ! Former revisions: … … 63 64 64 65 PRIVATE 65 PUBLIC prognostic_equations , prognostic_equations_fast, &66 prognostic_equations_vec 67 68 INTERFACE prognostic_equations 69 MODULE PROCEDURE prognostic_equations 70 END INTERFACE prognostic_equations 71 72 INTERFACE prognostic_equations_ fast73 MODULE PROCEDURE prognostic_equations_ fast74 END INTERFACE prognostic_equations_ fast75 76 INTERFACE prognostic_equations_vec 77 MODULE PROCEDURE prognostic_equations_vec 78 END INTERFACE prognostic_equations_vec 66 PUBLIC prognostic_equations_noopt, prognostic_equations_cache, & 67 prognostic_equations_vector 68 69 INTERFACE prognostic_equations_noopt 70 MODULE PROCEDURE prognostic_equations_noopt 71 END INTERFACE prognostic_equations_noopt 72 73 INTERFACE prognostic_equations_cache 74 MODULE PROCEDURE prognostic_equations_cache 75 END INTERFACE prognostic_equations_cache 76 77 INTERFACE prognostic_equations_vector 78 MODULE PROCEDURE prognostic_equations_vector 79 END INTERFACE prognostic_equations_vector 79 80 80 81 … … 82 83 83 84 84 SUBROUTINE prognostic_equations 85 SUBROUTINE prognostic_equations_noopt 85 86 86 87 !------------------------------------------------------------------------------! … … 615 616 616 617 617 END SUBROUTINE prognostic_equations 618 619 620 SUBROUTINE prognostic_equations_ fast618 END SUBROUTINE prognostic_equations_noopt 619 620 621 SUBROUTINE prognostic_equations_cache 621 622 622 623 !------------------------------------------------------------------------------! … … 996 997 997 998 998 END SUBROUTINE prognostic_equations_ fast999 1000 1001 SUBROUTINE prognostic_equations_vec 999 END SUBROUTINE prognostic_equations_cache 1000 1001 1002 SUBROUTINE prognostic_equations_vector 1002 1003 1003 1004 !------------------------------------------------------------------------------! … … 1242 1243 ! 1243 1244 !-- pt-tendency terms with communication 1245 sat = tsc(1) 1246 sbt = tsc(2) 1244 1247 IF ( scalar_advec == 'bc-scheme' ) THEN 1245 ! 1246 !-- Bott-Chlond scheme always uses Euler time step. Thus: 1247 sat = 1.0 1248 sbt = 1.0 1248 1249 IF ( timestep_scheme(1:5) /= 'runge' ) THEN 1250 ! 1251 !-- Bott-Chlond scheme always uses Euler time step when leapfrog is 1252 !-- switched on. Thus: 1253 sat = 1.0 1254 sbt = 1.0 1255 ENDIF 1249 1256 tend = 0.0 1250 1257 CALL advec_s_bc( pt, 'pt' ) 1251 1258 ELSE 1252 sat = tsc(1)1253 sbt = tsc(2)1254 1259 IF ( tsc(2) /= 2.0 .AND. scalar_advec == 'ups-scheme' ) THEN 1255 1260 tend = 0.0 … … 1340 1345 ! 1341 1346 !-- Scalar/q-tendency terms with communication 1347 sat = tsc(1) 1348 sbt = tsc(2) 1342 1349 IF ( scalar_advec == 'bc-scheme' ) THEN 1343 ! 1344 !-- Bott-Chlond scheme always uses Euler time step. Thus: 1345 sat = 1.0 1346 sbt = 1.0 1350 1351 IF ( timestep_scheme(1:5) /= 'runge' ) THEN 1352 ! 1353 !-- Bott-Chlond scheme always uses Euler time step when leapfrog is 1354 !-- switched on. Thus: 1355 sat = 1.0 1356 sbt = 1.0 1357 ENDIF 1347 1358 tend = 0.0 1348 1359 CALL advec_s_bc( q, 'q' ) 1349 1360 ELSE 1350 sat = tsc(1)1351 sbt = tsc(2)1352 1361 IF ( tsc(2) /= 2.0 ) THEN 1353 1362 IF ( scalar_advec == 'ups-scheme' ) THEN … … 1438 1447 !-- TKE-tendency terms with communication 1439 1448 CALL production_e_init 1449 1450 sat = tsc(1) 1451 sbt = tsc(2) 1440 1452 IF ( .NOT. use_upstream_for_tke ) THEN 1441 1453 IF ( scalar_advec == 'bc-scheme' ) THEN 1442 ! 1443 !-- Bott-Chlond scheme always uses Euler time step. Thus: 1444 sat = 1.0 1445 sbt = 1.0 1454 1455 IF ( timestep_scheme(1:5) /= 'runge' ) THEN 1456 ! 1457 !-- Bott-Chlond scheme always uses Euler time step when leapfrog is 1458 !-- switched on. Thus: 1459 sat = 1.0 1460 sbt = 1.0 1461 ENDIF 1446 1462 tend = 0.0 1447 1463 CALL advec_s_bc( e, 'e' ) 1448 1464 ELSE 1449 sat = tsc(1)1450 sbt = tsc(2)1451 1465 IF ( tsc(2) /= 2.0 ) THEN 1452 1466 IF ( scalar_advec == 'ups-scheme' ) THEN … … 1550 1564 1551 1565 1552 END SUBROUTINE prognostic_equations_vec 1566 END SUBROUTINE prognostic_equations_vector 1553 1567 1554 1568 -
palm/trunk/SOURCE/read_var_list.f90
r57 r63 4 4 ! Actual revisions: 5 5 ! ----------------- 6 ! + pt_reference6 ! +loop_optimization, pt_reference 7 7 ! 8 8 ! Former revisions: … … 208 208 CASE ( 'long_filter_factor' ) 209 209 READ ( 13 ) long_filter_factor 210 CASE ( 'loop_optimization' ) 211 READ ( 13 ) loop_optimization 210 212 CASE ( 'mixing_length_1d' ) 211 213 READ ( 13 ) mixing_length_1d -
palm/trunk/SOURCE/time_integration.f90
r48 r63 5 5 ! ----------------- 6 6 ! Move call of user_actions( 'after_integration' ) below increment of times 7 ! and counters 7 ! and counters, 8 ! calls of prognostic_equations_.. changed to .._noopt, .._cache, and 9 ! .._vector, these calls are now controlled by switch loop_optimization 8 10 ! 9 11 ! Former revisions: … … 98 100 !-- in the other versions a good vectorization is prohibited due to 99 101 !-- inlining problems. 100 IF ( host(1:3) == 'nec' ) THEN101 CALL prognostic_equations_vec 102 IF ( loop_optimization == 'vector' ) THEN 103 CALL prognostic_equations_vector 102 104 ELSE 103 105 IF ( momentum_advec == 'ups-scheme' .OR. & … … 105 107 scalar_advec == 'bc-scheme' ) & 106 108 THEN 107 CALL prognostic_equations 109 CALL prognostic_equations_noopt 108 110 ELSE 109 CALL prognostic_equations_ fast111 CALL prognostic_equations_cache 110 112 ENDIF 111 113 ENDIF … … 114 116 !-- Particle advection (only once during intermediate steps, because 115 117 !-- it uses an Euler-step) 116 IF ( simulated_time >= particle_advection_start .AND. & 118 IF ( particle_advection .AND. & 119 simulated_time >= particle_advection_start .AND. & 117 120 intermediate_timestep_count == 1 ) THEN 118 121 CALL advec_particles -
palm/trunk/SOURCE/write_var_list.f90
r57 r63 4 4 ! Actual revisions: 5 5 ! ----------------- 6 ! + pt_refrence6 ! +loop_optimization, pt_refrence 7 7 ! 8 8 ! Former revisions: … … 174 174 WRITE ( 14 ) 'long_filter_factor ' 175 175 WRITE ( 14 ) long_filter_factor 176 WRITE ( 14 ) 'loop_optimization ' 177 WRITE ( 14 ) loop_optimization 176 178 WRITE ( 14 ) 'mixing_length_1d ' 177 179 WRITE ( 14 ) mixing_length_1d
Note: See TracChangeset
for help on using the changeset viewer.
