- Timestamp:
- Dec 20, 2017 5:32:50 PM (7 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
palm/trunk/SOURCE/vertical_nesting_mod.f90
r2696 r2712 26 26 ! ----------------- 27 27 ! $Id$ 28 ! Formatting and clean-up (SadiqHuq) 29 ! 30 ! 2696 2017-12-14 17:12:51Z kanani 28 31 ! renamed diffusivities to tcm_diffusivities (TG) 29 32 ! … … 48 51 !> 49 52 !> Definition of parameters and variables for vertical nesting 53 !> The horizontal extent of the parent (Coarse Grid) and the child (Fine Grid) 54 !> have to be identical. The vertical extent of the FG should be smaller than CG. 55 !> Only integer nesting ratio supported. Odd nesting ratio preferred 56 !> The code follows MPI-1 standards. The available processors are split into 57 !> two groups using MPI_COMM_SPLIT. Exchange of data from CG to FG is called 58 !> interpolation. FG initialization by interpolation is done once at the start. 59 !> FG boundary conditions are set by interpolated at every timestep. 60 !> Exchange of data from CG to FG is called anterpolation, the two-way interaction 61 !> occurs at every timestep. 62 !> vnest_start_time set in PARIN allows delayed start of the coupling 63 !> after spin-up of the CG 50 64 !> 51 65 !> @todo Ensure that code can be compiled for serial and parallel mode. Please … … 64 78 IMPLICIT NONE 65 79 66 INTEGER(iwp),DIMENSION(3,2) :: bdims = 0 67 INTEGER(iwp),DIMENSION(3,2) :: bdims_rem = 0 !> Add description. It should not be longer than up to this point | 68 !> If really necessary, a second line can be added like this. 69 INTEGER(iwp) :: cg_nprocs, fg_nprocs 70 INTEGER(iwp),DIMENSION(:,:,:),ALLOCATABLE:: c2f_dims_cg, f2c_dims_cg 71 INTEGER(iwp),DIMENSION(:),ALLOCATABLE :: c2f_dims_fg, f2c_dims_fg 72 INTEGER(iwp) :: TYPE_VNEST_BC, TYPE_VNEST_ANTER 73 74 INTEGER(iwp),DIMENSION(:,:),ALLOCATABLE :: f_rnk_lst, c_rnk_lst 75 INTEGER(iwp),DIMENSION(3) :: cfratio 76 77 INTEGER(iwp) :: nxc, nxf, nyc, nyf, nzc, nzf 78 INTEGER(iwp) :: ngp_c, ngp_f 79 80 INTEGER(iwp) :: target_idex, n_cell_c, n_cell_f 81 INTEGER(iwp),DIMENSION(2) :: pdims_partner 82 INTEGER(iwp),DIMENSION(2) :: offset,map_coord 83 84 REAL(wp) :: dxc, dyc, dxf, dyf, dzc, dzf,dtc,dtf 85 86 REAL(wp),DIMENSION(:,:,:), ALLOCATABLE :: work3d 87 REAL(wp),DIMENSION(:,:), ALLOCATABLE :: work2dshf 88 REAL(wp),DIMENSION(:,:), ALLOCATABLE :: work2dusws 89 REAL(wp),DIMENSION(:,:), ALLOCATABLE :: work2dvsws 90 REAL(wp),DIMENSION(:,:), ALLOCATABLE :: work2dts 91 REAL(wp),DIMENSION(:,:), ALLOCATABLE :: work2dus 92 REAL(wp),DIMENSION(:,:), ALLOCATABLE :: work2dz0 93 94 95 REAL(wp), DIMENSION(:), ALLOCATABLE :: zuc, zuf, zwc, zwf 96 REAL(wp), DIMENSION(:,:,:), POINTER :: interpol3d,anterpol3d 97 ! REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: interpol3d 98 99 LOGICAL :: vnest_init = .FALSE., vnested = .FALSE., & 100 vnest_twi = .FALSE., vnest_couple_rk3 = .FALSE. 101 102 ! PARIN 103 REAL(wp) :: vnest_start_time = 9999999.9 80 LOGICAL :: vnested = .FALSE. !> set to true when 81 !> mrun is called with -N option 82 LOGICAL :: vnest_init = .FALSE. !> set to true when FG is initialized 83 REAL(wp) :: vnest_start_time = 9999999.9 !> simulated time when FG should be 84 !> initialized. Should be 85 !> identical in PARIN & PARIN_N 86 87 88 89 INTEGER(iwp),DIMENSION(3,2) :: bdims = 0 !> sub-domain grid topology of current PE 90 INTEGER(iwp),DIMENSION(3,2) :: bdims_rem = 0 !> sub-domain grid topology of partner PE 91 INTEGER(iwp) :: cg_nprocs !> no. of PE in CG. Set by mrun -N 92 INTEGER(iwp) :: fg_nprocs !> no. of PE in FG. Set by mrun -N 93 INTEGER(iwp) :: TYPE_VNEST_BC !> derived contiguous data type for interpolation 94 INTEGER(iwp) :: TYPE_VNEST_ANTER !> derived contiguous data type for anterpolation 95 INTEGER(iwp),DIMENSION(:,:,:),ALLOCATABLE :: c2f_dims_cg !> One CG PE sends data to multiple FG PEs 96 !> list of grid-topology of partners 97 INTEGER(iwp),DIMENSION(:,:,:),ALLOCATABLE :: f2c_dims_cg !> One CG PE receives data from multiple FG PEs 98 !> list of grid-topology of partners 99 INTEGER(iwp),DIMENSION(:),ALLOCATABLE :: c2f_dims_fg !> One FG PE sends data to multiple CG PE 100 !> list of grid-topology of partner 101 INTEGER(iwp),DIMENSION(:),ALLOCATABLE :: f2c_dims_fg !> One FG PE sends data to only one CG PE 102 !> list of grid-topology of partner 103 104 INTEGER(iwp),DIMENSION(:,:),ALLOCATABLE :: f_rnk_lst !> list storing rank of FG PE denoted by pdims 105 INTEGER(iwp),DIMENSION(:,:),ALLOCATABLE :: c_rnk_lst !> list storing rank of CG PE denoted by pdims 106 INTEGER(iwp),DIMENSION(3) :: cfratio !> Nesting ratio in x,y and z-directions 107 108 INTEGER(iwp) :: nxc !> no. of CG grid points in x-direction 109 INTEGER(iwp) :: nxf !> no. of FG grid points in x-direction 110 INTEGER(iwp) :: nyc !> no. of CG grid points in y-direction 111 INTEGER(iwp) :: nyf !> no. of FG grid points in y-direction 112 INTEGER(iwp) :: nzc !> no. of CG grid points in z-direction 113 INTEGER(iwp) :: nzf !> no. of FG grid points in z-direction 114 INTEGER(iwp) :: ngp_c !> no. of CG grid points in one vertical level 115 INTEGER(iwp) :: ngp_f !> no. of FG grid points in one vertical level 116 117 INTEGER(iwp) :: n_cell_c !> total no. of CG grid points in a PE 118 INTEGER(iwp) :: n_cell_f !> total no. of FG grid points in a PE 119 INTEGER(iwp),DIMENSION(2) :: pdims_partner !> processor topology of partner PE 120 INTEGER(iwp) :: target_idex !> temporary variable 121 INTEGER(iwp),DIMENSION(2) :: offset !> temporary variable 122 INTEGER(iwp),DIMENSION(2) :: map_coord !> temporary variable 123 124 REAL(wp) :: dxc !> CG grid pacing in x-direction 125 REAL(wp) :: dyc !> FG grid pacing in x-direction 126 REAL(wp) :: dxf !> CG grid pacing in y-direction 127 REAL(wp) :: dyf !> FG grid pacing in y-direction 128 REAL(wp) :: dzc !> CG grid pacing in z-direction 129 REAL(wp) :: dzf !> FG grid pacing in z-direction 130 REAL(wp) :: dtc !> dt calculated for CG 131 REAL(wp) :: dtf !> dt calculated for FG 132 133 REAL(wp), DIMENSION(:), ALLOCATABLE :: zuc !> CG vertical u-levels 134 REAL(wp), DIMENSION(:), ALLOCATABLE :: zuf !> FG vertical u-levels 135 REAL(wp), DIMENSION(:), ALLOCATABLE :: zwc !> CG vertical w-levels 136 REAL(wp), DIMENSION(:), ALLOCATABLE :: zwf !> FG vertical w-levels 137 REAL(wp), DIMENSION(:,:,:), POINTER :: interpol3d !> pointers to simplify function calls 138 REAL(wp), DIMENSION(:,:,:), POINTER :: anterpol3d !> pointers to simplify function calls 139 140 141 REAL(wp),DIMENSION(:,:,:), ALLOCATABLE :: work3d !> temporary array for exchange of 3D data 142 REAL(wp),DIMENSION(:,:), ALLOCATABLE :: work2dusws !> temporary array for exchange of 2D data 143 REAL(wp),DIMENSION(:,:), ALLOCATABLE :: work2dvsws !> temporary array for exchange of 2D data 144 REAL(wp),DIMENSION(:,:), ALLOCATABLE :: work2dts !> temporary array for exchange of 2D data 145 REAL(wp),DIMENSION(:,:), ALLOCATABLE :: work2dus !> temporary array for exchange of 2D data 104 146 105 147 SAVE … … 112 154 113 155 !-- Public constants and variables 114 PUBLIC vnested, vnest_init, vnest_twi, vnest_couple_rk3, & 115 vnest_start_time 156 PUBLIC vnested, vnest_init, vnest_start_time 116 157 117 158 PRIVATE bdims, bdims_rem, & 118 work3d, work2dshf, work2dusws, work2dvsws, & 119 work2dts, work2dus, work2dz0, & 159 work3d, work2dusws, work2dvsws, work2dts, work2dus, & 120 160 dxc, dyc, dxf, dyf, dzc, dzf, dtc, dtf, & 121 161 zuc, zuf, zwc, zwf, interpol3d, anterpol3d, & … … 176 216 177 217 SUBROUTINE vnest_init_fine 218 #if defined( __parallel ) 178 219 179 220 !--------------------------------------------------------------------------------! … … 191 232 USE interfaces 192 233 USE pegrid 193 USE surface_mod, 234 USE surface_mod, & 194 235 ONLY : surf_def_h, surf_def_v 195 USE turbulence_closure_mod, 236 USE turbulence_closure_mod, & 196 237 ONLY : tcm_diffusivities 238 197 239 198 240 IMPLICIT NONE 199 241 200 REAL(wp) :: time_since_reference_point_rem 201 202 INTEGER(iwp) :: i, j, k,im,jn,ko 203 INTEGER(iwp) :: if, jf, kf 204 205 #if defined( __parallel ) 206 207 if (myid ==0 )print *, ' TIME TO INIT FINE from COARSE', simulated_time 242 REAL(wp) :: time_since_reference_point_rem 243 INTEGER(iwp) :: i 244 INTEGER(iwp) :: j 245 INTEGER(iwp) :: k 246 INTEGER(iwp) :: im 247 INTEGER(iwp) :: jn 248 INTEGER(iwp) :: ko 249 INTEGER(iwp) :: iif 250 INTEGER(iwp) :: jjf 251 INTEGER(iwp) :: kkf 252 253 254 if (myid ==0 ) print *, ' TIME TO INIT FINE from COARSE', simulated_time 208 255 209 256 ! … … 344 391 (bdims(2,2)-bdims(2,1)+3) 345 392 346 !-- WARNING 347 !-- shf,z0 not interpolated 348 !-- line commented in interpolate_to_fine_flux 349 !-- FG needs to read it's own data file 350 !MERGE-WIP CALL MPI_SEND(surf_def_h(0)%shf ( bdims(2,1)-1:bdims(2,2)+1, & 351 !MERGE-WIP bdims(1,1)-1:bdims(1,2)+1),& 352 !MERGE-WIP n_cell_c, MPI_REAL, target_idex, & 353 !MERGE-WIP 109, comm_inter, ierr ) 354 !MERGE-WIP 355 !MERGE-WIP CALL MPI_SEND(surf_def_h(0)%usws( bdims(2,1)-1:bdims(2,2)+1, & 356 !MERGE-WIP bdims(1,1)-1:bdims(1,2)+1),& 357 !MERGE-WIP n_cell_c, MPI_REAL, target_idex, & 358 !MERGE-WIP 110, comm_inter, ierr ) 359 !MERGE-WIP 360 !MERGE-WIP CALL MPI_SEND(surf_def_h(0)%vsws( bdims(2,1)-1:bdims(2,2)+1, & 361 !MERGE-WIP bdims(1,1)-1:bdims(1,2)+1),& 362 !MERGE-WIP n_cell_c, MPI_REAL, target_idex, & 363 !MERGE-WIP 111, comm_inter, ierr ) 364 !MERGE-WIP 365 !MERGE CALL MPI_SEND(ts ( bdims(2,1)-1:bdims(2,2)+1, & 366 !MERGE bdims(1,1)-1:bdims(1,2)+1),& 367 !MERGE n_cell_c, MPI_REAL, target_idex, & 368 !MERGE 112, comm_inter, ierr ) 369 !MERGE 370 !MERGE CALL MPI_SEND(us ( bdims(2,1)-1:bdims(2,2)+1, & 371 !MERGE bdims(1,1)-1:bdims(1,2)+1),& 372 !MERGE n_cell_c, MPI_REAL, target_idex, & 373 !MERGE 113, comm_inter, ierr ) 374 !MERGE 375 !MERGE CALL MPI_SEND(z0 ( bdims(2,1)-1:bdims(2,2)+1, & 376 !MERGE bdims(1,1)-1:bdims(1,2)+1),& 377 !MERGE n_cell_c, MPI_REAL, target_idex, & 378 !MERGE 114, comm_inter, ierr ) 393 ! 394 !-- shf and z0 for CG / FG need to initialized in input file or user_code 395 !-- TODO 396 !-- initialization of usws, vsws, ts and us not vital to vnest FG 397 !-- variables are not compatible with the new surface layer module 398 ! 399 ! CALL MPI_SEND(surf_def_h(0)%usws( bdims(2,1)-1:bdims(2,2)+1, & 400 ! bdims(1,1)-1:bdims(1,2)+1),& 401 ! n_cell_c, MPI_REAL, target_idex, & 402 ! 110, comm_inter, ierr ) 403 ! 404 ! CALL MPI_SEND(surf_def_h(0)%vsws( bdims(2,1)-1:bdims(2,2)+1, & 405 ! bdims(1,1)-1:bdims(1,2)+1),& 406 ! n_cell_c, MPI_REAL, target_idex, & 407 ! 111, comm_inter, ierr ) 408 ! 409 ! CALL MPI_SEND(ts ( bdims(2,1)-1:bdims(2,2)+1, & 410 ! bdims(1,1)-1:bdims(1,2)+1),& 411 ! n_cell_c, MPI_REAL, target_idex, & 412 ! 112, comm_inter, ierr ) 413 ! 414 ! CALL MPI_SEND(us ( bdims(2,1)-1:bdims(2,2)+1, & 415 ! bdims(1,1)-1:bdims(1,2)+1),& 416 ! n_cell_c, MPI_REAL, target_idex, & 417 ! 113, comm_inter, ierr ) 418 ! 379 419 ENDIF 380 420 … … 465 505 (bdims_rem(2,2)-bdims_rem(2,1)+3) 466 506 467 ALLOCATE( work2dshf ( bdims_rem(2,1)-1:bdims_rem(2,2)+1, &468 bdims_rem(1,1)-1:bdims_rem(1,2)+1) )469 507 ALLOCATE( work2dusws ( bdims_rem(2,1)-1:bdims_rem(2,2)+1, & 470 508 bdims_rem(1,1)-1:bdims_rem(1,2)+1) ) … … 475 513 ALLOCATE( work2dus ( bdims_rem(2,1)-1:bdims_rem(2,2)+1, & 476 514 bdims_rem(1,1)-1:bdims_rem(1,2)+1) ) 477 ALLOCATE( work2dz0 ( bdims_rem(2,1)-1:bdims_rem(2,2)+1, & 478 bdims_rem(1,1)-1:bdims_rem(1,2)+1) ) 479 480 !MERGE-WIP CALL MPI_RECV( work2dshf ,n_cell_c, MPI_REAL, target_idex, 109, & 481 !MERGE-WIP comm_inter,status, ierr ) 482 !MERGE-WIP 483 !MERGE-WIP CALL MPI_RECV( work2dusws,n_cell_c, MPI_REAL, target_idex, 110, & 484 !MERGE-WIP comm_inter,status, ierr ) 485 !MERGE-WIP 486 !MERGE-WIP CALL MPI_RECV( work2dvsws,n_cell_c, MPI_REAL, target_idex, 111, & 487 !MERGE-WIP comm_inter,status, ierr ) 488 !MERGE-WIP 489 !MERGE CALL MPI_RECV( work2dts ,n_cell_c, MPI_REAL, target_idex, 112, & 490 !MERGE comm_inter,status, ierr ) 491 !MERGE 492 !MERGE CALL MPI_RECV( work2dus ,n_cell_c, MPI_REAL, target_idex, 113, & 493 !MERGE comm_inter,status, ierr ) 494 !MERGE 495 !MERGE CALL MPI_RECV( work2dz0 ,n_cell_c, MPI_REAL, target_idex, 114, & 496 !MERGE comm_inter,status, ierr ) 497 !MERGE 498 !MERGE-WIP CALL interpolate_to_fine_flux ( 108 ) 499 500 DEALLOCATE( work2dshf ) 515 516 ! 517 !-- shf and z0 for CG / FG need to initialized in input file or user_code 518 !-- TODO 519 !-- initialization of usws, vsws, ts and us not vital to vnest FG 520 !-- variables are not compatible with the new surface layer module 521 ! 522 ! CALL MPI_RECV( work2dusws,n_cell_c, MPI_REAL, target_idex, 110, & 523 ! comm_inter,status, ierr ) 524 ! 525 ! CALL MPI_RECV( work2dvsws,n_cell_c, MPI_REAL, target_idex, 111, & 526 ! comm_inter,status, ierr ) 527 ! 528 ! CALL MPI_RECV( work2dts ,n_cell_c, MPI_REAL, target_idex, 112, & 529 ! comm_inter,status, ierr ) 530 ! 531 ! CALL MPI_RECV( work2dus ,n_cell_c, MPI_REAL, target_idex, 113, & 532 ! comm_inter,status, ierr ) 533 ! 534 ! CALL interpolate_to_fine_flux ( 108 ) 535 501 536 DEALLOCATE( work2dusws ) 502 537 DEALLOCATE( work2dvsws ) 503 538 DEALLOCATE( work2dts ) 504 539 DEALLOCATE( work2dus ) 505 DEALLOCATE( work2dz0 )506 540 ENDIF 507 541 508 542 IF ( .NOT. constant_diffusion ) THEN 509 DO k f = bdims(3,1)+1,bdims(3,2)+1510 DO j f = bdims(2,1),bdims(2,2)511 DO i f = bdims(1,1),bdims(1,2)512 513 IF ( e(k f,jf,if) < 0.0 ) THEN514 e(k f,jf,if) = 1E-15_wp543 DO kkf = bdims(3,1)+1,bdims(3,2)+1 544 DO jjf = bdims(2,1),bdims(2,2) 545 DO iif = bdims(1,1),bdims(1,2) 546 547 IF ( e(kkf,jjf,iif) < 0.0 ) THEN 548 e(kkf,jjf,iif) = 1E-15_wp 515 549 END IF 516 550 … … 586 620 587 621 if (myid==0) print *, '** Fine Initalized ** simulated_time:', simulated_time 588 #endif 622 589 623 CONTAINS 590 624 591 625 SUBROUTINE interpolate_to_fine_w( tag ) 592 593 #if defined( __parallel )594 626 595 627 USE arrays_3d … … 602 634 IMPLICIT NONE 603 635 604 INTEGER(iwp), intent(in) :: tag 605 INTEGER(iwp) :: i, j, k 606 INTEGER(iwp) :: if, jf, kf 607 INTEGER(iwp) :: bottomx, topx 608 INTEGER(iwp) :: bottomy, topy 609 INTEGER(iwp) :: bottomz, topz 610 REAL(wp) :: eps, alpha, eminus, edot, eplus 611 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: wprs, wprf 612 613 INTEGER(iwp) :: nzbottom, nztop 636 INTEGER(iwp), intent(in) :: tag 637 INTEGER(iwp) :: i 638 INTEGER(iwp) :: j 639 INTEGER(iwp) :: k 640 INTEGER(iwp) :: iif 641 INTEGER(iwp) :: jjf 642 INTEGER(iwp) :: kkf 643 INTEGER(iwp) :: nzbottom 644 INTEGER(iwp) :: nztop 645 INTEGER(iwp) :: bottomx 646 INTEGER(iwp) :: bottomy 647 INTEGER(iwp) :: bottomz 648 INTEGER(iwp) :: topx 649 INTEGER(iwp) :: topy 650 INTEGER(iwp) :: topz 651 REAL(wp) :: eps 652 REAL(wp) :: alpha 653 REAL(wp) :: eminus 654 REAL(wp) :: edot 655 REAL(wp) :: eplus 656 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: wprs 657 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: wprf 614 658 615 659 … … 632 676 topx = (nxf+1)/(nxc+1) * (i+1) - 1 633 677 634 DO i f = bottomx, topx635 636 eps = ( i f * dxf + 0.5 * dxf - i * dxc - 0.5 * dxc ) / dxc678 DO iif = bottomx, topx 679 680 eps = ( iif * dxf + 0.5 * dxf - i * dxc - 0.5 * dxc ) / dxc 637 681 alpha = ( ( dxf / dxc )**2.0 - 1.0 ) / 24.0 638 682 eminus = eps * ( eps - 1.0 ) / 2.0 + alpha … … 640 684 eplus = eps * ( eps + 1.0 ) / 2.0 + alpha 641 685 642 wprf(k,j,i f) = eminus * work3d(k,j,i-1) &686 wprf(k,j,iif) = eminus * work3d(k,j,i-1) & 643 687 + edot * work3d(k,j,i) & 644 688 + eplus * work3d(k,j,i+1) … … 657 701 topy = (nyf+1)/(nyc+1) * (j+1) - 1 658 702 659 DO i f = nxl, nxr660 DO j f = bottomy, topy703 DO iif = nxl, nxr 704 DO jjf = bottomy, topy 661 705 662 eps = ( j f * dyf + 0.5 * dyf - j * dyc - 0.5 * dyc ) / dyc706 eps = ( jjf * dyf + 0.5 * dyf - j * dyc - 0.5 * dyc ) / dyc 663 707 alpha = ( ( dyf / dyc )**2.0 - 1.0 ) / 24.0 664 708 eminus = eps * ( eps - 1.0 ) / 2.0 + alpha … … 666 710 eplus = eps * ( eps + 1.0 ) / 2.0 + alpha 667 711 668 wprs(k,j f,if) = eminus * wprf(k,j-1,if) &669 + edot * wprf(k,j,i f) &670 + eplus * wprf(k,j+1,i f)712 wprs(k,jjf,iif) = eminus * wprf(k,j-1,iif) & 713 + edot * wprf(k,j,iif) & 714 + eplus * wprf(k,j+1,iif) 671 715 672 716 END DO … … 684 728 topz = (dzc/dzf) * (k+1) - 1 685 729 686 DO j f = nys, nyn687 DO i f = nxl, nxr688 DO k f = bottomz, topz689 690 w(k f,jf,if) = wprs(k,jf,if) + ( zwf(kf) - zwc(k) ) &691 * ( wprs(k+1,j f,if) - wprs(k,jf,if) ) / dzc730 DO jjf = nys, nyn 731 DO iif = nxl, nxr 732 DO kkf = bottomz, topz 733 734 w(kkf,jjf,iif) = wprs(k,jjf,iif) + ( zwf(kkf) - zwc(k) ) & 735 * ( wprs(k+1,jjf,iif) - wprs(k,jjf,iif) ) / dzc 692 736 693 737 END DO … … 697 741 END DO 698 742 699 DO j f = nys, nyn700 DO i f = nxl, nxr701 702 w(nzt,j f,if) = wprs(nztop,jf,if)743 DO jjf = nys, nyn 744 DO iif = nxl, nxr 745 746 w(nzt,jjf,iif) = wprs(nztop,jjf,iif) 703 747 704 748 END DO … … 709 753 DEALLOCATE( wprf, wprs ) 710 754 711 #endif712 755 END SUBROUTINE interpolate_to_fine_w 713 756 714 757 SUBROUTINE interpolate_to_fine_u( tag ) 715 758 716 #if defined( __parallel )717 759 718 760 USE arrays_3d … … 725 767 IMPLICIT NONE 726 768 727 INTEGER(iwp), intent(in) :: tag 728 INTEGER(iwp) :: i, j, k 729 INTEGER(iwp) :: if, jf, kf 730 INTEGER(iwp) :: bottomx, topx 731 INTEGER(iwp) :: bottomy, topy 732 INTEGER(iwp) :: bottomz, topz 733 REAL(wp) :: eps, alpha, eminus, edot, eplus 734 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: uprs, uprf 735 736 INTEGER(iwp) :: nzbottom, nztop 769 INTEGER(iwp), intent(in) :: tag 770 INTEGER(iwp) :: i 771 INTEGER(iwp) :: j 772 INTEGER(iwp) :: k 773 INTEGER(iwp) :: iif 774 INTEGER(iwp) :: jjf 775 INTEGER(iwp) :: kkf 776 INTEGER(iwp) :: nzbottom 777 INTEGER(iwp) :: nztop 778 INTEGER(iwp) :: bottomx 779 INTEGER(iwp) :: bottomy 780 INTEGER(iwp) :: bottomz 781 INTEGER(iwp) :: topx 782 INTEGER(iwp) :: topy 783 INTEGER(iwp) :: topz 784 REAL(wp) :: eps 785 REAL(wp) :: alpha 786 REAL(wp) :: eminus 787 REAL(wp) :: edot 788 REAL(wp) :: eplus 789 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: uprf 790 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: uprs 791 737 792 738 793 … … 756 811 757 812 DO i = bdims_rem(1,1)-1, bdims_rem(1,2)+1 758 DO j f = bottomy, topy759 760 eps = ( j f * dyf + 0.5 * dyf - j * dyc - 0.5 * dyc ) / dyc813 DO jjf = bottomy, topy 814 815 eps = ( jjf * dyf + 0.5 * dyf - j * dyc - 0.5 * dyc ) / dyc 761 816 alpha = ( ( dyf / dyc )**2.0 - 1.0 ) / 24.0 762 817 eminus = eps * ( eps - 1.0 ) / 2.0 + alpha … … 764 819 eplus = eps * ( eps + 1.0 ) / 2.0 + alpha 765 820 766 uprf(k,j f,i) = eminus * work3d(k,j-1,i) &821 uprf(k,jjf,i) = eminus * work3d(k,j-1,i) & 767 822 + edot * work3d(k,j,i) & 768 823 + eplus * work3d(k,j+1,i) … … 782 837 topz = (dzc/dzf) * k 783 838 784 DO j f = nys, nyn839 DO jjf = nys, nyn 785 840 DO i = bdims_rem(1,1)-1, bdims_rem(1,2)+1 786 DO k f = bottomz, topz841 DO kkf = bottomz, topz 787 842 788 eps = ( zuf(k f) - zuc(k) ) / dzc843 eps = ( zuf(kkf) - zuc(k) ) / dzc 789 844 alpha = ( ( dzf / dzc )**2.0 - 1.0 ) / 24.0 790 845 eminus = eps * ( eps - 1.0 ) / 2.0 + alpha … … 792 847 eplus = eps * ( eps + 1.0 ) / 2.0 + alpha 793 848 794 uprs(k f,jf,i) = eminus * uprf(k-1,jf,i) &795 + edot * uprf(k,j f,i) &796 + eplus * uprf(k+1,j f,i)849 uprs(kkf,jjf,i) = eminus * uprf(k-1,jjf,i) & 850 + edot * uprf(k,jjf,i) & 851 + eplus * uprf(k+1,jjf,i) 797 852 798 853 END DO … … 802 857 END DO 803 858 804 DO j f = nys, nyn859 DO jjf = nys, nyn 805 860 DO i = bdims_rem(1,1)-1, bdims_rem(1,2)+1 806 861 … … 811 866 eplus = eps * ( eps + 1.0 ) / 2.0 + alpha 812 867 813 uprs(nzt+1,j f,i) = eminus * uprf(nztop,jf,i) &814 + edot * uprf(nztop+1,j f,i) &815 + eplus * uprf(nztop+2,j f,i)868 uprs(nzt+1,jjf,i) = eminus * uprf(nztop,jjf,i) & 869 + edot * uprf(nztop+1,jjf,i) & 870 + eplus * uprf(nztop+2,jjf,i) 816 871 817 872 END DO … … 821 876 !-- Interpolation in x-direction (linear) 822 877 823 DO k f = nzb+1, nzt+1824 DO j f = nys, nyn878 DO kkf = nzb+1, nzt+1 879 DO jjf = nys, nyn 825 880 DO i = bdims_rem(1,1), bdims_rem(1,2) 826 881 … … 828 883 topx = (nxf+1)/(nxc+1) * (i+1) - 1 829 884 830 DO i f = bottomx, topx831 u(k f,jf,if) = uprs(kf,jf,i) + (if * dxf - i * dxc ) &832 * ( uprs(k f,jf,i+1) - uprs(kf,jf,i) ) / dxc885 DO iif = bottomx, topx 886 u(kkf,jjf,iif) = uprs(kkf,jjf,i) + ( iif * dxf - i * dxc ) & 887 * ( uprs(kkf,jjf,i+1) - uprs(kkf,jjf,i) ) / dxc 833 888 END DO 834 889 … … 843 898 DEALLOCATE( uprf, uprs ) 844 899 845 #endif846 900 END SUBROUTINE interpolate_to_fine_u 847 901 … … 849 903 SUBROUTINE interpolate_to_fine_v( tag ) 850 904 851 #if defined( __parallel )852 905 853 906 USE arrays_3d … … 859 912 860 913 IMPLICIT NONE 914 915 INTEGER(iwp), intent(in) :: tag 916 INTEGER(iwp) :: i 917 INTEGER(iwp) :: j 918 INTEGER(iwp) :: k 919 INTEGER(iwp) :: iif 920 INTEGER(iwp) :: jjf 921 INTEGER(iwp) :: kkf 922 INTEGER(iwp) :: nzbottom 923 INTEGER(iwp) :: nztop 924 INTEGER(iwp) :: bottomx 925 INTEGER(iwp) :: bottomy 926 INTEGER(iwp) :: bottomz 927 INTEGER(iwp) :: topx 928 INTEGER(iwp) :: topy 929 INTEGER(iwp) :: topz 930 REAL(wp) :: eps 931 REAL(wp) :: alpha 932 REAL(wp) :: eminus 933 REAL(wp) :: edot 934 REAL(wp) :: eplus 935 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: vprs 936 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: vprf 861 937 862 INTEGER(iwp), intent(in) :: tag863 INTEGER(iwp) :: i, j, k864 INTEGER(iwp) :: if, jf, kf865 INTEGER(iwp) :: bottomx, topx866 INTEGER(iwp) :: bottomy, topy867 INTEGER(iwp) :: bottomz, topz868 REAL(wp) :: eps, alpha, eminus, edot, eplus869 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: vprs, vprf870 871 INTEGER(iwp) :: nzbottom, nztop872 873 938 874 939 nzbottom = bdims_rem (3,1) … … 890 955 topx = (nxf+1)/(nxc+1) * (i+1) - 1 891 956 892 DO i f = bottomx, topx893 894 eps = ( i f * dxf + 0.5 * dxf - i * dxc - 0.5 * dxc ) / dxc957 DO iif = bottomx, topx 958 959 eps = ( iif * dxf + 0.5 * dxf - i * dxc - 0.5 * dxc ) / dxc 895 960 alpha = ( ( dxf / dxc )**2.0 - 1.0 ) / 24.0 896 961 eminus = eps * ( eps - 1.0 ) / 2.0 + alpha … … 898 963 eplus = eps * ( eps + 1.0 ) / 2.0 + alpha 899 964 900 vprf(k,j,i f) = eminus * work3d(k,j,i-1) &965 vprf(k,j,iif) = eminus * work3d(k,j,i-1) & 901 966 + edot * work3d(k,j,i) & 902 967 + eplus * work3d(k,j,i+1) … … 917 982 918 983 DO j = bdims_rem(2,1)-1, bdims_rem(2,2)+1 919 DO i f = nxl, nxr920 DO k f = bottomz, topz921 922 eps = ( zuf(k f) - zuc(k) ) / dzc984 DO iif = nxl, nxr 985 DO kkf = bottomz, topz 986 987 eps = ( zuf(kkf) - zuc(k) ) / dzc 923 988 alpha = ( ( dzf / dzc )**2.0 - 1.0 ) / 24.0 924 989 eminus = eps * ( eps - 1.0 ) / 2.0 + alpha … … 926 991 eplus = eps * ( eps + 1.0 ) / 2.0 + alpha 927 992 928 vprs(k f,j,if) = eminus * vprf(k-1,j,if) &929 + edot * vprf(k,j,i f) &930 + eplus * vprf(k+1,j,i f)993 vprs(kkf,j,iif) = eminus * vprf(k-1,j,iif) & 994 + edot * vprf(k,j,iif) & 995 + eplus * vprf(k+1,j,iif) 931 996 932 997 END DO … … 937 1002 938 1003 DO j = bdims_rem(2,1)-1, bdims_rem(2,2)+1 939 DO i f = nxl, nxr1004 DO iif = nxl, nxr 940 1005 941 1006 eps = ( zuf(nzt+1) - zuc(nztop+1) ) / dzc … … 945 1010 eplus = eps * ( eps + 1.0 ) / 2.0 + alpha 946 1011 947 vprs(nzt+1,j,i f) = eminus * vprf(nztop,j,if) &948 + edot * vprf(nztop+1,j,i f) &949 + eplus * vprf(nztop+2,j,i f)1012 vprs(nzt+1,j,iif) = eminus * vprf(nztop,j,iif) & 1013 + edot * vprf(nztop+1,j,iif) & 1014 + eplus * vprf(nztop+2,j,iif) 950 1015 951 1016 END DO … … 955 1020 !-- Interpolation in y-direction (linear) 956 1021 957 DO k f = nzb+1, nzt+11022 DO kkf = nzb+1, nzt+1 958 1023 DO j = bdims_rem(2,1), bdims_rem(2,2) 959 1024 … … 961 1026 topy = (nyf+1)/(nyc+1) * (j+1) - 1 962 1027 963 DO i f = nxl, nxr964 DO j f = bottomy, topy965 v (k f,jf,if) = vprs(kf,j,if) + (jf * dyf - j * dyc ) &966 * ( vprs(k f,j+1,if) - vprs(kf,j,if) ) / dyc1028 DO iif = nxl, nxr 1029 DO jjf = bottomy, topy 1030 v (kkf,jjf,iif) = vprs(kkf,j,iif) + ( jjf * dyf - j * dyc ) & 1031 * ( vprs(kkf,j+1,iif) - vprs(kkf,j,iif) ) / dyc 967 1032 END DO 968 1033 END DO … … 977 1042 DEALLOCATE( vprf, vprs ) 978 1043 979 #endif980 1044 END SUBROUTINE interpolate_to_fine_v 981 1045 … … 983 1047 SUBROUTINE interpolate_to_fine_s( tag ) 984 1048 985 #if defined( __parallel )986 1049 987 1050 USE arrays_3d … … 993 1056 994 1057 IMPLICIT NONE 995 996 997 INTEGER(iwp), intent(in) :: tag 998 INTEGER(iwp) :: i, j, k 999 INTEGER(iwp) :: if, jf, kf 1000 INTEGER(iwp) :: bottomx, topx 1001 INTEGER(iwp) :: bottomy, topy 1002 INTEGER(iwp) :: bottomz, topz 1003 REAL(wp) :: eps, alpha, eminus, edot, eplus 1004 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: ptprs, ptprf 1005 1006 INTEGER(iwp) :: nzbottom, nztop 1058 1059 INTEGER(iwp), intent(in) :: tag 1060 INTEGER(iwp) :: i 1061 INTEGER(iwp) :: j 1062 INTEGER(iwp) :: k 1063 INTEGER(iwp) :: iif 1064 INTEGER(iwp) :: jjf 1065 INTEGER(iwp) :: kkf 1066 INTEGER(iwp) :: nzbottom 1067 INTEGER(iwp) :: nztop 1068 INTEGER(iwp) :: bottomx 1069 INTEGER(iwp) :: bottomy 1070 INTEGER(iwp) :: bottomz 1071 INTEGER(iwp) :: topx 1072 INTEGER(iwp) :: topy 1073 INTEGER(iwp) :: topz 1074 REAL(wp) :: eps 1075 REAL(wp) :: alpha 1076 REAL(wp) :: eminus 1077 REAL(wp) :: edot 1078 REAL(wp) :: eplus 1079 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: ptprs 1080 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: ptprf 1007 1081 1008 1082 … … 1026 1100 topx = (nxf+1)/(nxc+1) * (i+1) - 1 1027 1101 1028 DO i f = bottomx, topx1029 1030 eps = ( i f * dxf + 0.5 * dxf - i * dxc - 0.5 * dxc ) / dxc1102 DO iif = bottomx, topx 1103 1104 eps = ( iif * dxf + 0.5 * dxf - i * dxc - 0.5 * dxc ) / dxc 1031 1105 alpha = ( ( dxf / dxc )**2.0 - 1.0 ) / 24.0 1032 1106 eminus = eps * ( eps - 1.0 ) / 2.0 + alpha … … 1034 1108 eplus = eps * ( eps + 1.0 ) / 2.0 + alpha 1035 1109 1036 ptprf(k,j,i f) = eminus * work3d(k,j,i-1) &1110 ptprf(k,j,iif) = eminus * work3d(k,j,i-1) & 1037 1111 + edot * work3d(k,j,i) & 1038 1112 + eplus * work3d(k,j,i+1) … … 1052 1126 topy = (nyf+1)/(nyc+1) * (j+1) - 1 1053 1127 1054 DO i f = nxl, nxr1055 DO j f = bottomy, topy1056 1057 eps = ( j f * dyf + 0.5 * dyf - j * dyc - 0.5 * dyc ) / dyc1128 DO iif = nxl, nxr 1129 DO jjf = bottomy, topy 1130 1131 eps = ( jjf * dyf + 0.5 * dyf - j * dyc - 0.5 * dyc ) / dyc 1058 1132 alpha = ( ( dyf / dyc )**2.0 - 1.0 ) / 24.0 1059 1133 eminus = eps * ( eps - 1.0 ) / 2.0 + alpha … … 1061 1135 eplus = eps * ( eps + 1.0 ) / 2.0 + alpha 1062 1136 1063 ptprs(k,j f,if) = eminus * ptprf(k,j-1,if) &1064 + edot * ptprf(k,j,i f) &1065 + eplus * ptprf(k,j+1,i f)1137 ptprs(k,jjf,iif) = eminus * ptprf(k,j-1,iif) & 1138 + edot * ptprf(k,j,iif) & 1139 + eplus * ptprf(k,j+1,iif) 1066 1140 1067 1141 END DO … … 1079 1153 topz = (dzc/dzf) * k 1080 1154 1081 DO j f = nys, nyn1082 DO i f = nxl, nxr1083 DO k f = bottomz, topz1155 DO jjf = nys, nyn 1156 DO iif = nxl, nxr 1157 DO kkf = bottomz, topz 1084 1158 1085 eps = ( zuf(k f) - zuc(k) ) / dzc1159 eps = ( zuf(kkf) - zuc(k) ) / dzc 1086 1160 alpha = ( ( dzf / dzc )**2.0 - 1.0 ) / 24.0 1087 1161 eminus = eps * ( eps - 1.0 ) / 2.0 + alpha … … 1089 1163 eplus = eps * ( eps + 1.0 ) / 2.0 + alpha 1090 1164 1091 interpol3d(k f,jf,if) = eminus * ptprs(k-1,jf,if) &1092 + edot * ptprs(k,j f,if) &1093 + eplus * ptprs(k+1,j f,if)1165 interpol3d(kkf,jjf,iif) = eminus * ptprs(k-1,jjf,iif) & 1166 + edot * ptprs(k,jjf,iif) & 1167 + eplus * ptprs(k+1,jjf,iif) 1094 1168 1095 1169 END DO … … 1099 1173 END DO 1100 1174 1101 DO j f = nys, nyn1102 DO i f = nxl, nxr1175 DO jjf = nys, nyn 1176 DO iif = nxl, nxr 1103 1177 1104 1178 eps = ( zuf(nzt+1) - zuc(nztop+1) ) / dzc … … 1108 1182 eplus = eps * ( eps + 1.0 ) / 2.0 + alpha 1109 1183 1110 interpol3d(nzt+1,j f,if) = eminus * ptprs(nztop,jf,if) &1111 + edot * ptprs(nztop+1,j f,if) &1112 + eplus * ptprs(nztop+2,j f,if)1184 interpol3d(nzt+1,jjf,iif) = eminus * ptprs(nztop,jjf,iif) & 1185 + edot * ptprs(nztop+1,jjf,iif) & 1186 + eplus * ptprs(nztop+2,jjf,iif) 1113 1187 1114 1188 END DO … … 1118 1192 DEALLOCATE( ptprf, ptprs ) 1119 1193 1120 #endif1121 1194 END SUBROUTINE interpolate_to_fine_s 1122 1195 … … 1124 1197 SUBROUTINE interpolate_to_fine_kh( tag ) 1125 1198 1126 #if defined( __parallel )1127 1199 1128 1200 USE arrays_3d … … 1134 1206 1135 1207 IMPLICIT NONE 1136 1137 1138 INTEGER(iwp), intent(in) :: tag 1139 INTEGER(iwp) :: i, j, k 1140 INTEGER(iwp) :: if, jf, kf 1141 INTEGER(iwp) :: bottomx, topx 1142 INTEGER(iwp) :: bottomy, topy 1143 INTEGER(iwp) :: bottomz, topz 1144 REAL(wp) :: eps, alpha, eminus, edot, eplus 1145 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: uprs, uprf 1146 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: vprs, vprf 1147 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: wprs, wprf 1148 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: ptprs, ptprf 1149 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: eprs, eprf 1150 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: kmprs, kmprf 1151 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: khprs, khprf 1152 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: shfpr, uswspr, vswspr 1153 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: tspr, uspr, z0pr 1154 1155 INTEGER(iwp) :: nzbottom, nztop 1208 1209 INTEGER(iwp), intent(in) :: tag 1210 INTEGER(iwp) :: i 1211 INTEGER(iwp) :: j 1212 INTEGER(iwp) :: k 1213 INTEGER(iwp) :: iif 1214 INTEGER(iwp) :: jjf 1215 INTEGER(iwp) :: kkf 1216 INTEGER(iwp) :: nzbottom 1217 INTEGER(iwp) :: nztop 1218 INTEGER(iwp) :: bottomx 1219 INTEGER(iwp) :: bottomy 1220 INTEGER(iwp) :: bottomz 1221 INTEGER(iwp) :: topx 1222 INTEGER(iwp) :: topy 1223 INTEGER(iwp) :: topz 1224 REAL(wp) :: eps 1225 REAL(wp) :: alpha 1226 REAL(wp) :: eminus 1227 REAL(wp) :: edot 1228 REAL(wp) :: eplus 1229 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: uprs 1230 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: vprs 1231 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: wprs 1232 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: ptprs 1233 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: eprs 1234 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: kmprs 1235 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: khprs 1236 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: tspr 1237 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: uprf 1238 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: vprf 1239 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: wprf 1240 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: ptprf 1241 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: eprf 1242 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: kmprf 1243 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: khprf 1244 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: uswspr 1245 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: vswspr 1246 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: uspr 1156 1247 1157 1248 … … 1178 1269 topx = (nxf+1)/(nxc+1) * (i+1) - 1 1179 1270 1180 DO i f = bottomx, topx1181 1182 eps = ( i f * dxf + 0.5 * dxf - i * dxc - 0.5 * dxc ) / dxc1271 DO iif = bottomx, topx 1272 1273 eps = ( iif * dxf + 0.5 * dxf - i * dxc - 0.5 * dxc ) / dxc 1183 1274 alpha = ( ( dxf / dxc )**2.0 - 1.0 ) / 24.0 1184 1275 eminus = eps * ( eps - 1.0 ) / 2.0 + alpha … … 1186 1277 eplus = eps * ( eps + 1.0 ) / 2.0 + alpha 1187 1278 1188 ptprf(k,j,i f) = eminus * work3d(k,j,i-1) &1279 ptprf(k,j,iif) = eminus * work3d(k,j,i-1) & 1189 1280 + edot * work3d(k,j,i) & 1190 1281 + eplus * work3d(k,j,i+1) … … 1204 1295 topy = (nyf+1)/(nyc+1) * (j+1) - 1 1205 1296 1206 DO i f = nxl, nxr1207 DO j f = bottomy, topy1208 1209 eps = ( j f * dyf + 0.5 * dyf - j * dyc - 0.5 * dyc ) / dyc1297 DO iif = nxl, nxr 1298 DO jjf = bottomy, topy 1299 1300 eps = ( jjf * dyf + 0.5 * dyf - j * dyc - 0.5 * dyc ) / dyc 1210 1301 alpha = ( ( dyf / dyc )**2.0 - 1.0 ) / 24.0 1211 1302 eminus = eps * ( eps - 1.0 ) / 2.0 + alpha … … 1213 1304 eplus = eps * ( eps + 1.0 ) / 2.0 + alpha 1214 1305 1215 ptprs(k,j f,if) = eminus * ptprf(k,j-1,if) &1216 + edot * ptprf(k,j,i f) &1217 + eplus * ptprf(k,j+1,i f)1306 ptprs(k,jjf,iif) = eminus * ptprf(k,j-1,iif) & 1307 + edot * ptprf(k,j,iif) & 1308 + eplus * ptprf(k,j+1,iif) 1218 1309 1219 1310 END DO … … 1231 1322 topz = (dzc/dzf) * k 1232 1323 1233 DO j f = nys, nyn1234 DO i f = nxl, nxr1235 DO k f = bottomz, topz1324 DO jjf = nys, nyn 1325 DO iif = nxl, nxr 1326 DO kkf = bottomz, topz 1236 1327 1237 eps = ( zuf(k f) - zuc(k) ) / dzc1328 eps = ( zuf(kkf) - zuc(k) ) / dzc 1238 1329 alpha = ( ( dzf / dzc )**2.0 - 1.0 ) / 24.0 1239 1330 eminus = eps * ( eps - 1.0 ) / 2.0 + alpha … … 1241 1332 eplus = eps * ( eps + 1.0 ) / 2.0 + alpha 1242 1333 1243 kh(k f,jf,if) = eminus * ptprs(k-1,jf,if) &1244 + edot * ptprs(k,j f,if) &1245 + eplus * ptprs(k+1,j f,if)1334 kh(kkf,jjf,iif) = eminus * ptprs(k-1,jjf,iif) & 1335 + edot * ptprs(k,jjf,iif) & 1336 + eplus * ptprs(k+1,jjf,iif) 1246 1337 1247 1338 END DO … … 1251 1342 END DO 1252 1343 1253 DO j f = nys, nyn1254 DO i f = nxl, nxr1344 DO jjf = nys, nyn 1345 DO iif = nxl, nxr 1255 1346 1256 1347 eps = ( zuf(nzt+1) - zuc(nztop+1) ) / dzc … … 1260 1351 eplus = eps * ( eps + 1.0 ) / 2.0 + alpha 1261 1352 1262 kh(nzt+1,j f,if) = eminus * ptprs(nztop,jf,if) &1263 + edot * ptprs(nztop+1,j f,if) &1264 + eplus * ptprs(nztop+2,j f,if)1353 kh(nzt+1,jjf,iif) = eminus * ptprs(nztop,jjf,iif) & 1354 + edot * ptprs(nztop+1,jjf,iif) & 1355 + eplus * ptprs(nztop+2,jjf,iif) 1265 1356 1266 1357 END DO … … 1270 1361 DEALLOCATE( ptprf, ptprs ) 1271 1362 1272 #endif1273 1363 END SUBROUTINE interpolate_to_fine_kh 1274 1364 1275 1365 SUBROUTINE interpolate_to_fine_km( tag ) 1276 1366 1277 #if defined( __parallel )1278 1367 1279 1368 USE arrays_3d … … 1285 1374 1286 1375 IMPLICIT NONE 1287 1288 1289 INTEGER(iwp), intent(in) :: tag 1290 INTEGER(iwp) :: i, j, k 1291 INTEGER(iwp) :: if, jf, kf 1292 INTEGER(iwp) :: bottomx, topx 1293 INTEGER(iwp) :: bottomy, topy 1294 INTEGER(iwp) :: bottomz, topz 1295 REAL(wp) :: eps, alpha, eminus, edot, eplus 1296 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: uprs, uprf 1297 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: vprs, vprf 1298 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: wprs, wprf 1299 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: ptprs, ptprf 1300 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: eprs, eprf 1301 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: kmprs, kmprf 1302 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: khprs, khprf 1303 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: shfpr, uswspr, vswspr 1304 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: tspr, uspr, z0pr 1305 1306 INTEGER(iwp) :: nzbottom, nztop 1376 1377 INTEGER(iwp), intent(in) :: tag 1378 INTEGER(iwp) :: i 1379 INTEGER(iwp) :: j 1380 INTEGER(iwp) :: k 1381 INTEGER(iwp) :: iif 1382 INTEGER(iwp) :: jjf 1383 INTEGER(iwp) :: kkf 1384 INTEGER(iwp) :: nzbottom 1385 INTEGER(iwp) :: nztop 1386 INTEGER(iwp) :: bottomx 1387 INTEGER(iwp) :: bottomy 1388 INTEGER(iwp) :: bottomz 1389 INTEGER(iwp) :: topx 1390 INTEGER(iwp) :: topy 1391 INTEGER(iwp) :: topz 1392 REAL(wp) :: eps 1393 REAL(wp) :: alpha 1394 REAL(wp) :: eminus 1395 REAL(wp) :: edot 1396 REAL(wp) :: eplus 1397 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: uprs 1398 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: vprs 1399 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: wprs 1400 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: ptprs 1401 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: eprs 1402 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: kmprs 1403 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: khprs 1404 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: vprf 1405 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: wprf 1406 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: ptprf 1407 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: eprf 1408 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: kmprf 1409 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: khprf 1410 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: uswspr 1411 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: vswspr 1412 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: tspr 1413 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: uspr 1307 1414 1308 1415 … … 1329 1436 topx = (nxf+1)/(nxc+1) * (i+1) - 1 1330 1437 1331 DO i f = bottomx, topx1332 1333 eps = ( i f * dxf + 0.5 * dxf - i * dxc - 0.5 * dxc ) / dxc1438 DO iif = bottomx, topx 1439 1440 eps = ( iif * dxf + 0.5 * dxf - i * dxc - 0.5 * dxc ) / dxc 1334 1441 alpha = ( ( dxf / dxc )**2.0 - 1.0 ) / 24.0 1335 1442 eminus = eps * ( eps - 1.0 ) / 2.0 + alpha … … 1337 1444 eplus = eps * ( eps + 1.0 ) / 2.0 + alpha 1338 1445 1339 ptprf(k,j,i f) = eminus * work3d(k,j,i-1) &1446 ptprf(k,j,iif) = eminus * work3d(k,j,i-1) & 1340 1447 + edot * work3d(k,j,i) & 1341 1448 + eplus * work3d(k,j,i+1) … … 1355 1462 topy = (nyf+1)/(nyc+1) * (j+1) - 1 1356 1463 1357 DO i f = nxl, nxr1358 DO j f = bottomy, topy1359 1360 eps = ( j f * dyf + 0.5 * dyf - j * dyc - 0.5 * dyc ) / dyc1464 DO iif = nxl, nxr 1465 DO jjf = bottomy, topy 1466 1467 eps = ( jjf * dyf + 0.5 * dyf - j * dyc - 0.5 * dyc ) / dyc 1361 1468 alpha = ( ( dyf / dyc )**2.0 - 1.0 ) / 24.0 1362 1469 eminus = eps * ( eps - 1.0 ) / 2.0 + alpha … … 1364 1471 eplus = eps * ( eps + 1.0 ) / 2.0 + alpha 1365 1472 1366 ptprs(k,j f,if) = eminus * ptprf(k,j-1,if) &1367 + edot * ptprf(k,j,i f) &1368 + eplus * ptprf(k,j+1,i f)1473 ptprs(k,jjf,iif) = eminus * ptprf(k,j-1,iif) & 1474 + edot * ptprf(k,j,iif) & 1475 + eplus * ptprf(k,j+1,iif) 1369 1476 1370 1477 END DO … … 1382 1489 topz = (dzc/dzf) * k 1383 1490 1384 DO j f = nys, nyn1385 DO i f = nxl, nxr1386 DO k f = bottomz, topz1491 DO jjf = nys, nyn 1492 DO iif = nxl, nxr 1493 DO kkf = bottomz, topz 1387 1494 1388 eps = ( zuf(k f) - zuc(k) ) / dzc1495 eps = ( zuf(kkf) - zuc(k) ) / dzc 1389 1496 alpha = ( ( dzf / dzc )**2.0 - 1.0 ) / 24.0 1390 1497 eminus = eps * ( eps - 1.0 ) / 2.0 + alpha … … 1392 1499 eplus = eps * ( eps + 1.0 ) / 2.0 + alpha 1393 1500 1394 km(k f,jf,if) = eminus * ptprs(k-1,jf,if) &1395 + edot * ptprs(k,j f,if) &1396 + eplus * ptprs(k+1,j f,if)1501 km(kkf,jjf,iif) = eminus * ptprs(k-1,jjf,iif) & 1502 + edot * ptprs(k,jjf,iif) & 1503 + eplus * ptprs(k+1,jjf,iif) 1397 1504 1398 1505 END DO … … 1402 1509 END DO 1403 1510 1404 DO j f = nys, nyn1405 DO i f = nxl, nxr1511 DO jjf = nys, nyn 1512 DO iif = nxl, nxr 1406 1513 1407 1514 eps = ( zuf(nzt+1) - zuc(nztop+1) ) / dzc … … 1411 1518 eplus = eps * ( eps + 1.0 ) / 2.0 + alpha 1412 1519 1413 km(nzt+1,j f,if) = eminus * ptprs(nztop,jf,if) &1414 + edot * ptprs(nztop+1,j f,if) &1415 + eplus * ptprs(nztop+2,j f,if)1520 km(nzt+1,jjf,iif) = eminus * ptprs(nztop,jjf,iif) & 1521 + edot * ptprs(nztop+1,jjf,iif) & 1522 + eplus * ptprs(nztop+2,jjf,iif) 1416 1523 1417 1524 END DO … … 1421 1528 DEALLOCATE( ptprf, ptprs ) 1422 1529 1423 #endif1424 1530 END SUBROUTINE interpolate_to_fine_km 1425 1531 … … 1429 1535 SUBROUTINE interpolate_to_fine_flux( tag ) 1430 1536 1431 #if defined( __parallel )1432 1537 1433 1538 USE arrays_3d … … 1439 1544 1440 1545 IMPLICIT NONE 1441 1442 1443 INTEGER(iwp), intent(in) :: tag 1444 INTEGER(iwp) :: i, j, k 1445 INTEGER(iwp) :: if, jf, kf 1446 INTEGER(iwp) :: bottomx, topx 1447 INTEGER(iwp) :: bottomy, topy 1448 INTEGER(iwp) :: bottomz, topz 1449 REAL(wp) :: eps, alpha, eminus, edot, eplus 1450 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: shfpr, uswspr, vswspr 1451 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: tspr, uspr, z0pr 1452 1453 1454 INTEGER(iwp) :: nzbottom, nztop 1455 1456 ALLOCATE( shfpr (bdims_rem(2,1)-1:bdims_rem(2,2)+1,nxl:nxr) ) 1546 1547 INTEGER(iwp), intent(in) :: tag 1548 INTEGER(iwp) :: i 1549 INTEGER(iwp) :: j 1550 INTEGER(iwp) :: k 1551 INTEGER(iwp) :: iif 1552 INTEGER(iwp) :: jjf 1553 INTEGER(iwp) :: kkf 1554 INTEGER(iwp) :: nzbottom 1555 INTEGER(iwp) :: nztop 1556 INTEGER(iwp) :: bottomx 1557 INTEGER(iwp) :: bottomy 1558 INTEGER(iwp) :: bottomz 1559 INTEGER(iwp) :: topx 1560 INTEGER(iwp) :: topy 1561 INTEGER(iwp) :: topz 1562 REAL(wp) :: eps 1563 REAL(wp) :: alpha 1564 REAL(wp) :: eminus 1565 REAL(wp) :: edot 1566 REAL(wp) :: eplus 1567 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: uswspr 1568 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: vswspr 1569 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: tspr 1570 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: uspr 1571 1457 1572 ALLOCATE( uswspr(bdims_rem(2,1)-1:bdims_rem(2,2)+1,nxl:nxr) ) 1458 1573 ALLOCATE( vswspr(bdims_rem(2,1)-1:bdims_rem(2,2)+1,nxl:nxr) ) 1459 1574 ALLOCATE( tspr (bdims_rem(2,1)-1:bdims_rem(2,2)+1,nxl:nxr) ) 1460 1575 ALLOCATE( uspr (bdims_rem(2,1)-1:bdims_rem(2,2)+1,nxl:nxr) ) 1461 ALLOCATE( z0pr (bdims_rem(2,1)-1:bdims_rem(2,2)+1,nxl:nxr) )1462 1576 1463 1577 ! … … 1473 1587 topx = (nxf+1)/(nxc+1) * (i+1) - 1 1474 1588 1475 DO i f = bottomx, topx1476 1477 eps = ( i f * dxf + 0.5 * dxf - i * dxc - 0.5 * dxc ) / dxc1589 DO iif = bottomx, topx 1590 1591 eps = ( iif * dxf + 0.5 * dxf - i * dxc - 0.5 * dxc ) / dxc 1478 1592 alpha = ( ( dxf / dxc )**2.0 - 1.0 ) / 24.0 1479 1593 eminus = eps * ( eps - 1.0 ) / 2.0 + alpha … … 1481 1595 eplus = eps * ( eps + 1.0 ) / 2.0 + alpha 1482 1596 1483 shfpr(j,if) = eminus * work2dshf(j,i-1) & 1484 + edot * work2dshf(j,i) & 1485 + eplus * work2dshf(j,i+1) 1486 1487 uswspr(j,if) = eminus * work2dusws(j,i-1) & 1597 uswspr(j,iif) = eminus * work2dusws(j,i-1) & 1488 1598 + edot * work2dusws(j,i) & 1489 1599 + eplus * work2dusws(j,i+1) 1490 1600 1491 vswspr(j,i f) = eminus * work2dvsws(j,i-1) &1601 vswspr(j,iif) = eminus * work2dvsws(j,i-1) & 1492 1602 + edot * work2dvsws(j,i) & 1493 1603 + eplus * work2dvsws(j,i+1) 1494 1604 1495 tspr(j,i f) = eminus * work2dts(j,i-1) &1605 tspr(j,iif) = eminus * work2dts(j,i-1) & 1496 1606 + edot * work2dts(j,i) & 1497 1607 + eplus * work2dts(j,i+1) 1498 1608 1499 uspr(j,i f) = eminus * work2dus(j,i-1) &1609 uspr(j,iif) = eminus * work2dus(j,i-1) & 1500 1610 + edot * work2dus(j,i) & 1501 1611 + eplus * work2dus(j,i+1) 1502 1503 z0pr(j,if) = eminus * work2dz0(j,i-1) &1504 + edot * work2dz0(j,i) &1505 + eplus * work2dz0(j,i+1)1506 1612 1507 1613 END DO … … 1518 1624 topy = (nyf+1)/(nyc+1) * (j+1) - 1 1519 1625 1520 DO i f = nxl, nxr1521 DO j f = bottomy, topy1522 1523 eps = ( j f * dyf + 0.5 * dyf - j * dyc - 0.5 * dyc ) / dyc1626 DO iif = nxl, nxr 1627 DO jjf = bottomy, topy 1628 1629 eps = ( jjf * dyf + 0.5 * dyf - j * dyc - 0.5 * dyc ) / dyc 1524 1630 alpha = ( ( dyf / dyc )**2.0 - 1.0 ) / 24.0 1525 1631 eminus = eps * ( eps - 1.0 ) / 2.0 + alpha 1526 1632 edot = ( 1.0 - eps**2.0 ) - 2.0 * alpha 1527 1633 eplus = eps * ( eps + 1.0 ) / 2.0 + alpha 1528 1529 !-- WARNING 1530 !-- shf,z0 not interpolated 1531 !-- line commented in interpolate_to_fine_flux 1532 !-- FG needs to read it's own data file 1533 !MERGE-WIP! surf_def_h(0)%shf(jf,if) = eminus * shfpr(j-1,if) & 1534 !MERGE-WIP! + edot * shfpr(j,if) & 1535 !MERGE-WIP! + eplus * shfpr(j+1,if) 1536 !MERGE-WIP 1537 !MERGE-WIP surf_def_h(0)%usws(jf,if) = eminus * uswspr(j-1,if) & 1538 !MERGE-WIP + edot * uswspr(j,if) & 1539 !MERGE-WIP + eplus * uswspr(j+1,if) 1540 !MERGE-WIP 1541 !MERGE-WIP surf_def_h(0)%vsws(jf,if) = eminus * vswspr(j-1,if) & 1542 !MERGE-WIP + edot * vswspr(j,if) & 1543 !MERGE-WIP + eplus * vswspr(j+1,if) 1544 !MERGE-WIP 1545 !MERGE ts(jf,if) = eminus * tspr(j-1,if) & 1546 !MERGE + edot * tspr(j,if) & 1547 !MERGE + eplus * tspr(j+1,if) 1548 !MERGE 1549 !MERGE us(jf,if) = eminus * uspr(j-1,if) & 1550 !MERGE + edot * uspr(j,if) & 1551 !MERGE + eplus * uspr(j+1,if) 1552 !MERGE 1553 !MERGE! z0(jf,if) = eminus * z0pr(j-1,if) & 1554 !MERGE! + edot * z0pr(j,if) & 1555 !MERGE! + eplus * z0pr(j+1,if) 1634 1635 !! 1636 !!-- TODO 1637 !-- variables are not compatible with the new surface layer module 1638 ! 1639 ! surf_def_h(0)%usws(jjf,iif) = eminus * uswspr(j-1,if) & 1640 ! + edot * uswspr(j,iif) & 1641 ! + eplus * uswspr(j+1,iif) 1642 ! 1643 ! surf_def_h(0)%vsws(jjf,iif) = eminus * vswspr(j-1,if) & 1644 ! + edot * vswspr(j,iif) & 1645 ! + eplus * vswspr(j+1,iif) 1646 ! 1647 ! ts(jjf,iif) = eminus * tspr(j-1,if) & 1648 ! + edot * tspr(j,iif) & 1649 ! + eplus * tspr(j+1,iif) 1650 ! 1651 ! us(jjf,iif) = eminus * uspr(j-1,if) & 1652 ! + edot * uspr(j,iif) & 1653 ! + eplus * uspr(j+1,iif) 1556 1654 1557 1655 END DO … … 1561 1659 1562 1660 1563 DEALLOCATE( shfpr, uswspr, vswspr ) 1564 DEALLOCATE( tspr, uspr, z0pr ) 1565 1566 1567 #endif 1661 DEALLOCATE( uswspr, vswspr ) 1662 DEALLOCATE( tspr, uspr ) 1663 1664 1568 1665 END SUBROUTINE interpolate_to_fine_flux 1569 1666 1570 1667 1668 #endif 1571 1669 END SUBROUTINE vnest_init_fine 1572 1670 1573 1671 SUBROUTINE vnest_boundary_conds 1672 #if defined( __parallel ) 1574 1673 !------------------------------------------------------------------------------! 1575 1674 ! Description: … … 1591 1690 IMPLICIT NONE 1592 1691 1593 INTEGER(iwp) :: i, j, k 1594 INTEGER(iwp) :: if, jf 1595 1596 REAL(wp) :: c_max, denom 1597 1598 #if defined( __parallel ) 1692 INTEGER(iwp) :: i 1693 INTEGER(iwp) :: j 1694 INTEGER(iwp) :: iif 1695 INTEGER(iwp) :: jjf 1696 REAL(wp) :: c_max 1697 REAL(wp) :: denom 1698 1599 1699 1600 1700 ! … … 1707 1807 1708 1808 !-- TKE Neumann BC for FG top 1709 DO j f = nys, nyn1710 DO i f = nxl, nxr1711 e(nzt+1,j f,if) = e(nzt,jf,if)1809 DO jjf = nys, nyn 1810 DO iif = nxl, nxr 1811 e(nzt+1,jjf,iif) = e(nzt,jjf,iif) 1712 1812 END DO 1713 1813 END DO … … 1731 1831 1732 1832 1733 #endif1734 1833 CONTAINS 1735 1834 1736 1835 SUBROUTINE vnest_set_topbc_w 1737 1836 1738 #if defined( __parallel )1739 1837 1740 1838 USE arrays_3d … … 1746 1844 1747 1845 IMPLICIT NONE 1748 1749 INTEGER(iwp) :: i, j, k 1750 INTEGER(iwp) :: if, jf 1751 INTEGER(iwp) :: bottomx, topx 1752 INTEGER(iwp) :: bottomy, topy 1753 REAL(wp) :: eps, alpha, eminus, edot, eplus 1754 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: wprf 1755 1846 1847 INTEGER(iwp) :: i 1848 INTEGER(iwp) :: j 1849 INTEGER(iwp) :: k 1850 INTEGER(iwp) :: iif 1851 INTEGER(iwp) :: jjf 1852 INTEGER(iwp) :: kkf 1853 INTEGER(iwp) :: bottomx 1854 INTEGER(iwp) :: bottomy 1855 INTEGER(iwp) :: topx 1856 INTEGER(iwp) :: topy 1857 REAL(wp) :: eps 1858 REAL(wp) :: alpha 1859 REAL(wp) :: eminus 1860 REAL(wp) :: edot 1861 REAL(wp) :: eplus 1862 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: wprf 1863 1756 1864 1757 1865 ALLOCATE( wprf(bdims_rem(2,1)-1:bdims_rem(2,2)+1,nxl:nxr) ) … … 1773 1881 topx = (nxf+1)/(nxc+1) * (i+1) - 1 1774 1882 1775 DO i f = bottomx, topx1776 1777 eps = (i f * dxf + 0.5 * dxf - i * dxc - 0.5 * dxc) / dxc1883 DO iif = bottomx, topx 1884 1885 eps = (iif * dxf + 0.5 * dxf - i * dxc - 0.5 * dxc) / dxc 1778 1886 alpha = ( (dxf/dxc)**2.0 - 1.0) / 24.0 1779 1887 eminus = eps * ( eps - 1.0 ) / 2.0 + alpha 1780 1888 edot = ( 1.0 - eps**2.0 ) - 2.0 * alpha 1781 1889 eplus = eps * ( eps + 1.0 ) / 2.0 + alpha 1782 wprf(j,i f) = eminus * work3d(bdims_rem(3,1),j,i-1) &1890 wprf(j,iif) = eminus * work3d(bdims_rem(3,1),j,i-1) & 1783 1891 + edot * work3d(bdims_rem(3,1),j,i) & 1784 1892 + eplus * work3d(bdims_rem(3,1),j,i+1) … … 1798 1906 topy = (nyf+1)/(nyc+1) * (j+1) - 1 1799 1907 1800 DO i f = nxl, nxr1801 1802 DO j f = bottomy, topy1803 1804 eps = (j f * dyf + 0.5 * dyf - j * dyc - 0.5 * dyc) / dyc1908 DO iif = nxl, nxr 1909 1910 DO jjf = bottomy, topy 1911 1912 eps = (jjf * dyf + 0.5 * dyf - j * dyc - 0.5 * dyc) / dyc 1805 1913 1806 1914 alpha = ( (dyf/dyc)**2.0 - 1.0) / 24.0 … … 1812 1920 eplus = eps * ( eps + 1.0 ) / 2.0 + alpha 1813 1921 1814 w(nzt,j f,if) = eminus * wprf(j-1,if) &1815 + edot * wprf(j,i f) &1816 + eplus * wprf(j+1,i f)1922 w(nzt,jjf,iif) = eminus * wprf(j-1,iif) & 1923 + edot * wprf(j,iif) & 1924 + eplus * wprf(j+1,iif) 1817 1925 1818 1926 END DO … … 1823 1931 1824 1932 DEALLOCATE( wprf ) 1825 #endif1826 1933 1827 1934 END SUBROUTINE vnest_set_topbc_w … … 1830 1937 SUBROUTINE vnest_set_topbc_u 1831 1938 1832 #if defined( __parallel )1833 1939 1834 1940 USE arrays_3d … … 1840 1946 1841 1947 IMPLICIT NONE 1842 1843 INTEGER(iwp) :: i, j, k 1844 INTEGER(iwp) :: if, jf 1845 INTEGER(iwp) :: bottomx, topx 1846 INTEGER(iwp) :: bottomy, topy 1847 REAL(wp) :: eps, alpha, eminus, edot, eplus 1848 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: uprf 1849 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: uprs 1850 1851 1948 1949 INTEGER(iwp) :: i 1950 INTEGER(iwp) :: j 1951 INTEGER(iwp) :: k 1952 INTEGER(iwp) :: iif 1953 INTEGER(iwp) :: jjf 1954 INTEGER(iwp) :: kkf 1955 INTEGER(iwp) :: bottomx 1956 INTEGER(iwp) :: bottomy 1957 INTEGER(iwp) :: topx 1958 INTEGER(iwp) :: topy 1959 REAL(wp) :: eps 1960 REAL(wp) :: alpha 1961 REAL(wp) :: eminus 1962 REAL(wp) :: edot 1963 REAL(wp) :: eplus 1964 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: uprf 1965 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: uprs 1852 1966 1853 1967 ALLOCATE( uprf(bdims_rem(3,1):bdims_rem(3,2),nys:nyn,bdims_rem(1,1)-1:bdims_rem(1,2)+1) ) … … 1865 1979 1866 1980 DO i = bdims_rem(1,1)-1, bdims_rem(1,2)+1 1867 DO j f = bottomy, topy1868 1869 eps = (j f * dyf + 0.5 * dyf - j * dyc - 0.5 * dyc) / dyc1981 DO jjf = bottomy, topy 1982 1983 eps = (jjf * dyf + 0.5 * dyf - j * dyc - 0.5 * dyc) / dyc 1870 1984 alpha = ( (dyf/dyc)**2.0 - 1.0) / 24.0 1871 1985 eminus = eps * ( eps - 1.0 ) / 2.0 + alpha … … 1873 1987 eplus = eps * ( eps + 1.0 ) / 2.0 + alpha 1874 1988 1875 uprf(k,j f,i) = eminus * work3d(k,j-1,i) &1989 uprf(k,jjf,i) = eminus * work3d(k,j-1,i) & 1876 1990 + edot * work3d(k,j,i) & 1877 1991 + eplus * work3d(k,j+1,i) … … 1885 1999 !-- Interpolation in z-direction 1886 2000 1887 DO j f = nys, nyn2001 DO jjf = nys, nyn 1888 2002 DO i = bdims_rem(1,1)-1, bdims_rem(1,2)+1 1889 2003 eps = ( zuf(nzt+1) - zuc(bdims_rem(3,1)+1) ) / dzc … … 1892 2006 edot = ( 1.0 - eps**2.0 ) - 2.0 * alpha 1893 2007 eplus = eps * ( eps + 1.0 ) / 2.0 + alpha 1894 uprs(j f,i) = eminus * uprf(bdims_rem(3,1),jf,i) &1895 + edot * uprf(bdims_rem(3,1)+1,j f,i) &1896 + eplus * uprf(bdims_rem(3,1)+2,j f,i)2008 uprs(jjf,i) = eminus * uprf(bdims_rem(3,1),jjf,i) & 2009 + edot * uprf(bdims_rem(3,1)+1,jjf,i) & 2010 + eplus * uprf(bdims_rem(3,1)+2,jjf,i) 1897 2011 END DO 1898 2012 END DO … … 1901 2015 !-- Interpolation in x-direction 1902 2016 1903 DO j f = nys, nyn2017 DO jjf = nys, nyn 1904 2018 DO i = bdims_rem(1,1), bdims_rem(1,2) 1905 2019 … … 1907 2021 topx = (nxf+1)/(nxc+1) * (i+1) - 1 1908 2022 1909 DO i f = bottomx, topx1910 u(nzt+1,j f,if) = uprs(jf,i) + ( if * dxf - i * dxc ) * ( uprs(jf,i+1) - uprs(jf,i) ) / dxc2023 DO iif = bottomx, topx 2024 u(nzt+1,jjf,iif) = uprs(jjf,i) + ( iif * dxf - i * dxc ) * ( uprs(jjf,i+1) - uprs(jjf,i) ) / dxc 1911 2025 END DO 1912 2026 … … 1917 2031 1918 2032 DEALLOCATE ( uprf, uprs ) 1919 #endif1920 2033 1921 2034 END SUBROUTINE vnest_set_topbc_u … … 1924 2037 SUBROUTINE vnest_set_topbc_v 1925 2038 1926 #if defined( __parallel )1927 2039 1928 2040 USE arrays_3d … … 1934 2046 1935 2047 IMPLICIT NONE 1936 1937 INTEGER(iwp) :: i, j, k 1938 INTEGER(iwp) :: if, jf 1939 INTEGER(iwp) :: bottomx, topx 1940 INTEGER(iwp) :: bottomy, topy 1941 REAL(wp) :: eps, alpha, eminus, edot, eplus 1942 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: vprf 1943 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: vprs 2048 2049 INTEGER(iwp) :: i 2050 INTEGER(iwp) :: j 2051 INTEGER(iwp) :: k 2052 INTEGER(iwp) :: iif 2053 INTEGER(iwp) :: jjf 2054 INTEGER(iwp) :: kkf 2055 INTEGER(iwp) :: bottomx 2056 INTEGER(iwp) :: bottomy 2057 INTEGER(iwp) :: topx 2058 INTEGER(iwp) :: topy 2059 REAL(wp) :: eps 2060 REAL(wp) :: alpha 2061 REAL(wp) :: eminus 2062 REAL(wp) :: edot 2063 REAL(wp) :: eplus 2064 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: vprf 2065 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: vprs 1944 2066 1945 2067 … … 1963 2085 topx = (nxf+1)/(nxc+1) * (i+1) - 1 1964 2086 1965 DO i f = bottomx, topx1966 1967 eps = (i f * dxf + 0.5 * dxf - i * dxc - 0.5 * dxc) / dxc2087 DO iif = bottomx, topx 2088 2089 eps = (iif * dxf + 0.5 * dxf - i * dxc - 0.5 * dxc) / dxc 1968 2090 alpha = ( (dxf/dxc)**2.0 - 1.0) / 24.0 1969 2091 eminus = eps * ( eps - 1.0 ) / 2.0 + alpha 1970 2092 edot = ( 1.0 - eps**2.0 ) - 2.0 * alpha 1971 2093 eplus = eps * ( eps + 1.0 ) / 2.0 + alpha 1972 vprf(k,j,i f) = eminus * work3d(k,j,i-1) &2094 vprf(k,j,iif) = eminus * work3d(k,j,i-1) & 1973 2095 + edot * work3d(k,j,i) & 1974 2096 + eplus * work3d(k,j,i+1) … … 1983 2105 1984 2106 DO j = bdims_rem(2,1)-1, bdims_rem(2,2)+1 1985 DO i f = nxl, nxr2107 DO iif = nxl, nxr 1986 2108 1987 2109 eps = ( zuf(nzt+1) - zuc(bdims_rem(3,1)+1) ) / dzc … … 1990 2112 edot = ( 1.0 - eps**2.0 ) - 2.0 * alpha 1991 2113 eplus = eps * ( eps + 1.0 ) / 2.0 + alpha 1992 vprs(j,i f) = eminus * vprf(bdims_rem(3,1),j,if) &1993 + edot * vprf(bdims_rem(3,1)+1,j,i f) &1994 + eplus * vprf(bdims_rem(3,1)+2,j,i f)2114 vprs(j,iif) = eminus * vprf(bdims_rem(3,1),j,iif) & 2115 + edot * vprf(bdims_rem(3,1)+1,j,iif) & 2116 + eplus * vprf(bdims_rem(3,1)+2,j,iif) 1995 2117 1996 2118 END DO … … 2001 2123 2002 2124 DO j = bdims_rem(2,1), bdims_rem(2,2) 2003 DO i f = nxl, nxr2125 DO iif = nxl, nxr 2004 2126 2005 2127 bottomy = (nyf+1)/(nyc+1) * j 2006 2128 topy = (nyf+1)/(nyc+1) * (j+1) - 1 2007 2129 2008 DO j f = bottomy, topy2009 2010 v(nzt+1,j f,if) = vprs(j,if) + ( jf * dyf - j * dyc ) * ( vprs(j+1,if) - vprs(j,if) ) / dyc2130 DO jjf = bottomy, topy 2131 2132 v(nzt+1,jjf,iif) = vprs(j,iif) + ( jjf * dyf - j * dyc ) * ( vprs(j+1,iif) - vprs(j,iif) ) / dyc 2011 2133 2012 2134 END DO … … 2018 2140 2019 2141 2020 #endif2021 2142 2022 2143 END SUBROUTINE vnest_set_topbc_v … … 2025 2146 SUBROUTINE vnest_set_topbc_s 2026 2147 2027 #if defined( __parallel )2028 2148 2029 2149 USE arrays_3d … … 2035 2155 2036 2156 IMPLICIT NONE 2037 2038 INTEGER(iwp) :: i, j, k 2039 INTEGER(iwp) :: if, jf 2040 INTEGER(iwp) :: bottomx, topx 2041 INTEGER(iwp) :: bottomy, topy 2042 REAL(wp) :: eps, alpha, eminus, edot, eplus 2043 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: ptprf, ptprs 2044 2157 2158 INTEGER(iwp) :: i 2159 INTEGER(iwp) :: j 2160 INTEGER(iwp) :: k 2161 INTEGER(iwp) :: iif 2162 INTEGER(iwp) :: jjf 2163 INTEGER(iwp) :: kkf 2164 INTEGER(iwp) :: bottomx 2165 INTEGER(iwp) :: bottomy 2166 INTEGER(iwp) :: topx 2167 INTEGER(iwp) :: topy 2168 REAL(wp) :: eps 2169 REAL(wp) :: alpha 2170 REAL(wp) :: eminus 2171 REAL(wp) :: edot 2172 REAL(wp) :: eplus 2173 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: ptprf 2174 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: ptprs 2175 2045 2176 2046 2177 … … 2064 2195 topx = (nxf+1)/(nxc+1) *(i+1) - 1 2065 2196 2066 DO i f = bottomx, topx2067 2068 eps = (i f * dxf + 0.5 * dxf - i * dxc - 0.5 * dxc) / dxc2197 DO iif = bottomx, topx 2198 2199 eps = (iif * dxf + 0.5 * dxf - i * dxc - 0.5 * dxc) / dxc 2069 2200 2070 2201 alpha = ( (dxf/dxc)**2.0 - 1.0) / 24.0 … … 2076 2207 eplus = eps * ( eps + 1.0 ) / 2.0 + alpha 2077 2208 2078 ptprf(k,j,i f) = eminus * work3d(k,j,i-1) &2209 ptprf(k,j,iif) = eminus * work3d(k,j,i-1) & 2079 2210 + edot * work3d(k,j,i) & 2080 2211 + eplus * work3d(k,j,i+1) … … 2097 2228 topy = (nyf+1)/(nyc+1) * (j+1) - 1 2098 2229 2099 DO i f = nxl, nxr2100 2101 DO j f = bottomy, topy2102 2103 eps = (j f * dyf + 0.5 * dyf - j * dyc - 0.5 * dyc) / dyc2230 DO iif = nxl, nxr 2231 2232 DO jjf = bottomy, topy 2233 2234 eps = (jjf * dyf + 0.5 * dyf - j * dyc - 0.5 * dyc) / dyc 2104 2235 2105 2236 alpha = ( (dyf/dyc)**2.0 - 1.0) / 24.0 … … 2111 2242 eplus = eps * ( eps + 1.0 ) / 2.0 + alpha 2112 2243 2113 ptprs(k,j f,if) = eminus * ptprf(k,j-1,if) &2114 + edot * ptprf(k,j,i f) &2115 + eplus * ptprf(k,j+1,i f)2244 ptprs(k,jjf,iif) = eminus * ptprf(k,j-1,iif) & 2245 + edot * ptprf(k,j,iif) & 2246 + eplus * ptprf(k,j+1,iif) 2116 2247 END DO 2117 2248 … … 2125 2256 !-- Interpolation in z-direction 2126 2257 2127 DO j f = nys, nyn2128 DO i f = nxl, nxr2258 DO jjf = nys, nyn 2259 DO iif = nxl, nxr 2129 2260 2130 2261 eps = ( zuf(nzt+1) - zuc(bdims_rem(3,1)+1) ) / dzc … … 2138 2269 eplus = eps * ( eps + 1.0 ) / 2.0 + alpha 2139 2270 2140 interpol3d (nzt+1,j f,if) = eminus * ptprs(bdims_rem(3,1),jf,if) &2141 + edot * ptprs(bdims_rem(3,1)+1,j f,if) &2142 + eplus * ptprs(bdims_rem(3,1)+2,j f,if)2271 interpol3d (nzt+1,jjf,iif) = eminus * ptprs(bdims_rem(3,1),jjf,iif) & 2272 + edot * ptprs(bdims_rem(3,1)+1,jjf,iif) & 2273 + eplus * ptprs(bdims_rem(3,1)+2,jjf,iif) 2143 2274 2144 2275 END DO … … 2148 2279 2149 2280 2150 #endif2151 2281 2152 2282 END SUBROUTINE vnest_set_topbc_s 2283 #endif 2153 2284 END SUBROUTINE vnest_boundary_conds 2154 2285 2155 2286 2156 2287 SUBROUTINE vnest_boundary_conds_khkm 2288 #if defined( __parallel ) 2157 2289 2158 2290 !--------------------------------------------------------------------------------! … … 2175 2307 IMPLICIT NONE 2176 2308 2177 INTEGER(iwp) :: i, j, k 2178 INTEGER(iwp) :: if, jf 2179 2180 REAL(wp) :: c_max, denom 2181 2182 #if defined( __parallel ) 2309 INTEGER(iwp) :: i 2310 INTEGER(iwp) :: j 2311 INTEGER(iwp) :: iif 2312 INTEGER(iwp) :: jjf 2313 REAL(wp) :: c_max 2314 REAL(wp) :: denom 2315 2183 2316 2184 2317 IF ( coupling_mode == 'vnested_crse' ) THEN … … 2265 2398 2266 2399 ! Neumann BC for FG kh 2267 DO j f = nys, nyn2268 DO i f = nxl, nxr2269 kh(nzt+1,j f,if) = kh(nzt,jf,if)2400 DO jjf = nys, nyn 2401 DO iif = nxl, nxr 2402 kh(nzt+1,jjf,iif) = kh(nzt,jjf,iif) 2270 2403 END DO 2271 2404 END DO … … 2275 2408 2276 2409 ! Neumann BC for FG kh 2277 DO j f = nys, nyn2278 DO i f = nxl, nxr2279 km(nzt+1,j f,if) = km(nzt,jf,if)2410 DO jjf = nys, nyn 2411 DO iif = nxl, nxr 2412 km(nzt+1,jjf,iif) = km(nzt,jjf,iif) 2280 2413 END DO 2281 2414 END DO … … 2284 2417 ! 2285 2418 !-- The following evaluation can only be performed, if the fine grid is situated below the inversion 2286 !! DO j f = nys-1, nyn+12287 !! DO i f = nxl-1, nxr+12419 !! DO jjf = nys-1, nyn+1 2420 !! DO iif = nxl-1, nxr+1 2288 2421 !! 2289 !! km(nzt+1,j f,if) = 0.1 * l_grid(nzt+1) * SQRT( e(nzt+1,jf,if) )2290 !! kh(nzt+1,j f,if) = 3.0 * km(nzt+1,jf,if)2422 !! km(nzt+1,jjf,iif) = 0.1 * l_grid(nzt+1) * SQRT( e(nzt+1,jjf,iif) ) 2423 !! kh(nzt+1,jjf,iif) = 3.0 * km(nzt+1,jjf,iif) 2291 2424 !! 2292 2425 !! END DO … … 2300 2433 ENDIF 2301 2434 2302 #endif2303 2435 2304 2436 CONTAINS … … 2306 2438 SUBROUTINE vnest_set_topbc_kh 2307 2439 2308 #if defined( __parallel )2309 2440 2310 2441 USE arrays_3d … … 2316 2447 2317 2448 IMPLICIT NONE 2318 2319 INTEGER(iwp) :: i, j, k 2320 INTEGER(iwp) :: if, jf 2321 INTEGER(iwp) :: bottomx, topx 2322 INTEGER(iwp) :: bottomy, topy 2323 REAL(wp) :: eps, alpha, eminus, edot, eplus 2324 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: ptprf, ptprs 2325 2449 2450 INTEGER(iwp) :: i 2451 INTEGER(iwp) :: j 2452 INTEGER(iwp) :: k 2453 INTEGER(iwp) :: iif 2454 INTEGER(iwp) :: jjf 2455 INTEGER(iwp) :: kkf 2456 INTEGER(iwp) :: bottomx 2457 INTEGER(iwp) :: bottomy 2458 INTEGER(iwp) :: topx 2459 INTEGER(iwp) :: topy 2460 REAL(wp) :: eps 2461 REAL(wp) :: alpha 2462 REAL(wp) :: eminus 2463 REAL(wp) :: edot 2464 REAL(wp) :: eplus 2465 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: ptprf 2466 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: ptprs 2467 2326 2468 2327 2469 … … 2345 2487 topx = (nxf+1)/(nxc+1) *(i+1) - 1 2346 2488 2347 DO i f = bottomx, topx2348 2349 eps = (i f * dxf + 0.5 * dxf - i * dxc - 0.5 * dxc) / dxc2489 DO iif = bottomx, topx 2490 2491 eps = (iif * dxf + 0.5 * dxf - i * dxc - 0.5 * dxc) / dxc 2350 2492 2351 2493 alpha = ( (dxf/dxc)**2.0 - 1.0) / 24.0 … … 2357 2499 eplus = eps * ( eps + 1.0 ) / 2.0 + alpha 2358 2500 2359 ptprf(k,j,i f) = eminus * work3d(k,j,i-1) &2501 ptprf(k,j,iif) = eminus * work3d(k,j,i-1) & 2360 2502 + edot * work3d(k,j,i) & 2361 2503 + eplus * work3d(k,j,i+1) … … 2378 2520 topy = (nyf+1)/(nyc+1) * (j+1) - 1 2379 2521 2380 DO i f = nxl, nxr2381 2382 DO j f = bottomy, topy2383 2384 eps = (j f * dyf + 0.5 * dyf - j * dyc - 0.5 * dyc) / dyc2522 DO iif = nxl, nxr 2523 2524 DO jjf = bottomy, topy 2525 2526 eps = (jjf * dyf + 0.5 * dyf - j * dyc - 0.5 * dyc) / dyc 2385 2527 2386 2528 alpha = ( (dyf/dyc)**2.0 - 1.0) / 24.0 … … 2392 2534 eplus = eps * ( eps + 1.0 ) / 2.0 + alpha 2393 2535 2394 ptprs(k,j f,if) = eminus * ptprf(k,j-1,if) &2395 + edot * ptprf(k,j,i f) &2396 + eplus * ptprf(k,j+1,i f)2536 ptprs(k,jjf,iif) = eminus * ptprf(k,j-1,iif) & 2537 + edot * ptprf(k,j,iif) & 2538 + eplus * ptprf(k,j+1,iif) 2397 2539 END DO 2398 2540 … … 2406 2548 !-- Interpolation in z-direction 2407 2549 2408 DO j f = nys, nyn2409 DO i f = nxl, nxr2550 DO jjf = nys, nyn 2551 DO iif = nxl, nxr 2410 2552 2411 2553 eps = ( zuf(nzt+1) - zuc(bdims_rem(3,1)+1) ) / dzc … … 2419 2561 eplus = eps * ( eps + 1.0 ) / 2.0 + alpha 2420 2562 2421 kh (nzt+1,j f,if) = eminus * ptprs(bdims_rem(3,1),jf,if) &2422 + edot * ptprs(bdims_rem(3,1)+1,j f,if) &2423 + eplus * ptprs(bdims_rem(3,1)+2,j f,if)2563 kh (nzt+1,jjf,iif) = eminus * ptprs(bdims_rem(3,1),jjf,iif) & 2564 + edot * ptprs(bdims_rem(3,1)+1,jjf,iif) & 2565 + eplus * ptprs(bdims_rem(3,1)+2,jjf,iif) 2424 2566 2425 2567 END DO … … 2429 2571 2430 2572 2431 #endif2432 2573 2433 2574 END SUBROUTINE vnest_set_topbc_kh … … 2435 2576 SUBROUTINE vnest_set_topbc_km 2436 2577 2437 #if defined( __parallel )2438 2578 2439 2579 USE arrays_3d … … 2445 2585 2446 2586 IMPLICIT NONE 2447 2448 INTEGER(iwp) :: i, j, k 2449 INTEGER(iwp) :: if, jf 2450 INTEGER(iwp) :: bottomx, topx 2451 INTEGER(iwp) :: bottomy, topy 2452 REAL(wp) :: eps, alpha, eminus, edot, eplus 2453 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: ptprf, ptprs 2454 2587 2588 INTEGER(iwp) :: i 2589 INTEGER(iwp) :: j 2590 INTEGER(iwp) :: k 2591 INTEGER(iwp) :: iif 2592 INTEGER(iwp) :: jjf 2593 INTEGER(iwp) :: bottomx 2594 INTEGER(iwp) :: bottomy 2595 INTEGER(iwp) :: bottomz 2596 INTEGER(iwp) :: topx 2597 INTEGER(iwp) :: topy 2598 INTEGER(iwp) :: topz 2599 REAL(wp) :: eps 2600 REAL(wp) :: alpha 2601 REAL(wp) :: eminus 2602 REAL(wp) :: edot 2603 REAL(wp) :: eplus 2604 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: ptprf 2605 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: ptprs 2606 2455 2607 2456 2608 … … 2474 2626 topx = (nxf+1)/(nxc+1) *(i+1) - 1 2475 2627 2476 DO i f = bottomx, topx2477 2478 eps = (i f * dxf + 0.5 * dxf - i * dxc - 0.5 * dxc) / dxc2628 DO iif = bottomx, topx 2629 2630 eps = (iif * dxf + 0.5 * dxf - i * dxc - 0.5 * dxc) / dxc 2479 2631 2480 2632 alpha = ( (dxf/dxc)**2.0 - 1.0) / 24.0 … … 2486 2638 eplus = eps * ( eps + 1.0 ) / 2.0 + alpha 2487 2639 2488 ptprf(k,j,i f) = eminus * work3d(k,j,i-1) &2640 ptprf(k,j,iif) = eminus * work3d(k,j,i-1) & 2489 2641 + edot * work3d(k,j,i) & 2490 2642 + eplus * work3d(k,j,i+1) … … 2507 2659 topy = (nyf+1)/(nyc+1) * (j+1) - 1 2508 2660 2509 DO i f = nxl, nxr2510 2511 DO j f = bottomy, topy2512 2513 eps = (j f * dyf + 0.5 * dyf - j * dyc - 0.5 * dyc) / dyc2661 DO iif = nxl, nxr 2662 2663 DO jjf = bottomy, topy 2664 2665 eps = (jjf * dyf + 0.5 * dyf - j * dyc - 0.5 * dyc) / dyc 2514 2666 2515 2667 alpha = ( (dyf/dyc)**2.0 - 1.0) / 24.0 … … 2521 2673 eplus = eps * ( eps + 1.0 ) / 2.0 + alpha 2522 2674 2523 ptprs(k,j f,if) = eminus * ptprf(k,j-1,if) &2524 + edot * ptprf(k,j,i f) &2525 + eplus * ptprf(k,j+1,i f)2675 ptprs(k,jjf,iif) = eminus * ptprf(k,j-1,iif) & 2676 + edot * ptprf(k,j,iif) & 2677 + eplus * ptprf(k,j+1,iif) 2526 2678 END DO 2527 2679 … … 2535 2687 !-- Interpolation in z-direction 2536 2688 2537 DO j f = nys, nyn2538 DO i f = nxl, nxr2689 DO jjf = nys, nyn 2690 DO iif = nxl, nxr 2539 2691 2540 2692 eps = ( zuf(nzt+1) - zuc(bdims_rem(3,1)+1) ) / dzc … … 2548 2700 eplus = eps * ( eps + 1.0 ) / 2.0 + alpha 2549 2701 2550 km (nzt+1,j f,if) = eminus * ptprs(bdims_rem(3,1),jf,if) &2551 + edot * ptprs(bdims_rem(3,1)+1,j f,if) &2552 + eplus * ptprs(bdims_rem(3,1)+2,j f,if)2702 km (nzt+1,jjf,iif) = eminus * ptprs(bdims_rem(3,1),jjf,iif) & 2703 + edot * ptprs(bdims_rem(3,1)+1,jjf,iif) & 2704 + eplus * ptprs(bdims_rem(3,1)+2,jjf,iif) 2553 2705 2554 2706 END DO … … 2558 2710 2559 2711 2560 #endif2561 2712 2562 2713 END SUBROUTINE vnest_set_topbc_km 2563 2714 2564 2715 2716 #endif 2565 2717 END SUBROUTINE vnest_boundary_conds_khkm 2566 2718 … … 2568 2720 2569 2721 SUBROUTINE vnest_anterpolate 2722 2723 #if defined( __parallel ) 2570 2724 2571 2725 !--------------------------------------------------------------------------------! … … 2587 2741 IMPLICIT NONE 2588 2742 2589 REAL(wp) :: time_since_reference_point_rem 2590 INTEGER(iwp) :: i, j, k, im, jn, ko 2591 2592 !--- INTEGER(iwp) :: j !< grid index y direction 2593 !-- INTEGER(iwp) :: k !< grid index z direction 2594 INTEGER(iwp) :: kb !< variable to set respective boundary value, depends on facing. 2595 INTEGER(iwp) :: l !< running index boundary type, for up- and downward-facing walls 2596 INTEGER(iwp) :: m !< running index surface elements 2597 2598 #if defined( __parallel ) 2743 REAL(wp) :: time_since_reference_point_rem 2744 INTEGER(iwp) :: i 2745 INTEGER(iwp) :: j 2746 INTEGER(iwp) :: k 2747 INTEGER(iwp) :: im 2748 INTEGER(iwp) :: jn 2749 INTEGER(iwp) :: ko 2750 INTEGER(iwp) :: kb !< variable to set respective boundary value, depends on facing. 2751 INTEGER(iwp) :: l !< running index boundary type, for up- and downward-facing walls 2752 INTEGER(iwp) :: m !< running index surface elements 2753 2599 2754 2600 2755 … … 2835 2990 2836 2991 2837 #endif2838 2992 2839 2993 CONTAINS 2840 2994 SUBROUTINE anterpolate_to_crse_u( tag ) 2841 2995 2842 #if defined( __parallel )2843 2996 2844 2997 USE arrays_3d … … 2850 3003 2851 3004 IMPLICIT NONE 2852 2853 INTEGER(iwp) :: i, j, k 2854 INTEGER(iwp) :: if, jf, kf 2855 INTEGER(iwp) :: bottomx, topx 2856 INTEGER(iwp) :: bottomy, topy 2857 INTEGER(iwp) :: bottomz, topz 2858 REAL(wp) :: aweight 2859 INTEGER(iwp), intent(in) :: tag 3005 3006 INTEGER(iwp), intent(in) :: tag 3007 INTEGER(iwp) :: i 3008 INTEGER(iwp) :: j 3009 INTEGER(iwp) :: k 3010 INTEGER(iwp) :: iif 3011 INTEGER(iwp) :: jjf 3012 INTEGER(iwp) :: kkf 3013 INTEGER(iwp) :: bottomx 3014 INTEGER(iwp) :: bottomy 3015 INTEGER(iwp) :: bottomz 3016 INTEGER(iwp) :: topx 3017 INTEGER(iwp) :: topy 3018 INTEGER(iwp) :: topz 3019 REAL(wp) :: aweight 2860 3020 2861 3021 ! … … 2876 3036 DO i = bdims_rem(1,1),bdims_rem(1,2) 2877 3037 2878 i f = (nxf+1) / (nxc+1) * i3038 iif = (nxf+1) / (nxc+1) * i 2879 3039 2880 3040 aweight = 0.0 2881 3041 2882 DO k f = bottomz, topz2883 DO j f = bottomy, topy2884 2885 aweight = aweight + anterpol3d(k f,jf,if) * &3042 DO kkf = bottomz, topz 3043 DO jjf = bottomy, topy 3044 3045 aweight = aweight + anterpol3d(kkf,jjf,iif) * & 2886 3046 (dzf/dzc) * (dyf/dyc) 2887 3047 … … 2898 3058 2899 3059 2900 #endif2901 3060 2902 3061 END SUBROUTINE anterpolate_to_crse_u … … 2905 3064 SUBROUTINE anterpolate_to_crse_v( tag ) 2906 3065 2907 #if defined( __parallel )2908 3066 2909 3067 USE arrays_3d … … 2915 3073 2916 3074 IMPLICIT NONE 2917 2918 INTEGER(iwp) :: i, j, k 2919 INTEGER(iwp) :: if, jf, kf 2920 INTEGER(iwp) :: bottomx, topx 2921 INTEGER(iwp) :: bottomy, topy 2922 INTEGER(iwp) :: bottomz, topz 2923 REAL(wp) :: aweight 2924 INTEGER(iwp), intent(in) :: tag 3075 3076 INTEGER(iwp), intent(in) :: tag 3077 INTEGER(iwp) :: i 3078 INTEGER(iwp) :: j 3079 INTEGER(iwp) :: k 3080 INTEGER(iwp) :: iif 3081 INTEGER(iwp) :: jjf 3082 INTEGER(iwp) :: kkf 3083 INTEGER(iwp) :: bottomx 3084 INTEGER(iwp) :: bottomy 3085 INTEGER(iwp) :: bottomz 3086 INTEGER(iwp) :: topx 3087 INTEGER(iwp) :: topy 3088 INTEGER(iwp) :: topz 3089 REAL(wp) :: aweight 3090 2925 3091 ! 2926 3092 !-- Anterpolation of the velocity components v … … 2935 3101 DO j = bdims_rem(2,1), bdims_rem(2,2) 2936 3102 2937 j f = (nyf+1) / (nyc+1) * j3103 jjf = (nyf+1) / (nyc+1) * j 2938 3104 2939 3105 DO i = bdims_rem(1,1), bdims_rem(1,2) … … 2944 3110 aweight = 0.0 2945 3111 2946 DO k f = bottomz, topz2947 DO i f = bottomx, topx2948 2949 aweight = aweight + anterpol3d(k f,jf,if) * &3112 DO kkf = bottomz, topz 3113 DO iif = bottomx, topx 3114 3115 aweight = aweight + anterpol3d(kkf,jjf,iif) * & 2950 3116 (dzf/dzc) * (dxf/dxc) 2951 3117 … … 2961 3127 2962 3128 2963 #endif2964 3129 2965 3130 END SUBROUTINE anterpolate_to_crse_v … … 2968 3133 SUBROUTINE anterpolate_to_crse_w( tag ) 2969 3134 2970 #if defined( __parallel )2971 3135 2972 3136 USE arrays_3d … … 2978 3142 2979 3143 IMPLICIT NONE 2980 2981 INTEGER(iwp) :: i, j, k 2982 INTEGER(iwp) :: if, jf, kf 2983 INTEGER(iwp) :: bottomx, topx 2984 INTEGER(iwp) :: bottomy, topy 2985 INTEGER(iwp) :: bottomz, topz 2986 REAL(wp) :: aweight 2987 INTEGER(iwp), intent(in) :: tag 3144 3145 INTEGER(iwp), intent(in) :: tag 3146 INTEGER(iwp) :: i 3147 INTEGER(iwp) :: j 3148 INTEGER(iwp) :: k 3149 INTEGER(iwp) :: iif 3150 INTEGER(iwp) :: jjf 3151 INTEGER(iwp) :: kkf 3152 INTEGER(iwp) :: bottomx 3153 INTEGER(iwp) :: bottomy 3154 INTEGER(iwp) :: bottomz 3155 INTEGER(iwp) :: topx 3156 INTEGER(iwp) :: topy 3157 INTEGER(iwp) :: topz 3158 REAL(wp) :: aweight 3159 2988 3160 ! 2989 3161 !-- Anterpolation of the velocity components w … … 2993 3165 DO k = bdims_rem(3,1), bdims_rem(3,2)-1 2994 3166 2995 k f = cfratio(3) * k3167 kkf = cfratio(3) * k 2996 3168 2997 3169 DO j = bdims_rem(2,1), bdims_rem(2,2) … … 3007 3179 aweight = 0.0 3008 3180 3009 DO j f = bottomy, topy3010 DO i f = bottomx, topx3011 3012 aweight = aweight + anterpol3d (k f,jf,if) * &3181 DO jjf = bottomy, topy 3182 DO iif = bottomx, topx 3183 3184 aweight = aweight + anterpol3d (kkf,jjf,iif) * & 3013 3185 (dxf/dxc) * (dyf/dyc) 3014 3186 … … 3024 3196 END DO 3025 3197 3026 #endif3027 3198 3028 3199 END SUBROUTINE anterpolate_to_crse_w … … 3031 3202 SUBROUTINE anterpolate_to_crse_s( tag ) 3032 3203 3033 #if defined( __parallel )3034 3204 3035 3205 USE arrays_3d … … 3041 3211 3042 3212 IMPLICIT NONE 3043 3044 INTEGER(iwp) :: i, j, k 3045 INTEGER(iwp) :: if, jf, kf 3046 INTEGER(iwp) :: bottomx, topx 3047 INTEGER(iwp) :: bottomy, topy 3048 INTEGER(iwp) :: bottomz, topz 3049 REAL(wp) :: aweight 3050 INTEGER(iwp), intent(in) :: tag 3213 3214 INTEGER(iwp), intent(in) :: tag 3215 INTEGER(iwp) :: i 3216 INTEGER(iwp) :: j 3217 INTEGER(iwp) :: k 3218 INTEGER(iwp) :: iif 3219 INTEGER(iwp) :: jjf 3220 INTEGER(iwp) :: kkf 3221 INTEGER(iwp) :: bottomx 3222 INTEGER(iwp) :: bottomy 3223 INTEGER(iwp) :: bottomz 3224 INTEGER(iwp) :: topx 3225 INTEGER(iwp) :: topy 3226 INTEGER(iwp) :: topz 3227 REAL(wp) :: aweight 3051 3228 3052 3229 ! … … 3071 3248 aweight = 0.0 3072 3249 3073 DO k f = bottomz, topz3074 DO j f = bottomy, topy3075 DO i f = bottomx, topx3076 3077 aweight = aweight + anterpol3d(k f,jf,if) * &3250 DO kkf = bottomz, topz 3251 DO jjf = bottomy, topy 3252 DO iif = bottomx, topx 3253 3254 aweight = aweight + anterpol3d(kkf,jjf,iif) * & 3078 3255 (dzf/dzc) * (dyf/dyc) * (dxf/dxc) 3079 3256 … … 3090 3267 END DO 3091 3268 3092 #endif3093 3269 3094 3270 END SUBROUTINE anterpolate_to_crse_s 3271 #endif 3095 3272 END SUBROUTINE vnest_anterpolate 3096 3273 … … 3098 3275 3099 3276 SUBROUTINE vnest_anterpolate_e 3277 #if defined( __parallel ) 3100 3278 3101 3279 !--------------------------------------------------------------------------------! … … 3115 3293 IMPLICIT NONE 3116 3294 3117 REAL(wp) :: time_since_reference_point_rem 3118 INTEGER(iwp) :: i, j, k, im, jn, ko 3119 3120 #if defined( __parallel ) 3295 REAL(wp) :: time_since_reference_point_rem 3296 INTEGER(iwp) :: i 3297 INTEGER(iwp) :: j 3298 INTEGER(iwp) :: k 3299 INTEGER(iwp) :: im 3300 INTEGER(iwp) :: jn 3301 INTEGER(iwp) :: ko 3302 3121 3303 3122 3304 ! … … 3243 3425 ENDIF 3244 3426 3245 #endif3246 3427 3247 3428 CONTAINS … … 3253 3434 SUBROUTINE anterpolate_to_crse_e( tag ) 3254 3435 3255 #if defined( __parallel )3256 3436 3257 3437 USE arrays_3d … … 3263 3443 3264 3444 IMPLICIT NONE 3265 3266 INTEGER(iwp) :: i, j, k 3267 INTEGER(iwp) :: if, jf, kf 3268 INTEGER(iwp) :: bottomx, topx 3269 INTEGER(iwp) :: bottomy, topy 3270 INTEGER(iwp) :: bottomz, topz 3271 REAL(wp) :: aweight_a, aweight_b, aweight_c, aweight_d, aweight_e 3272 REAL(wp) :: energ 3273 INTEGER(iwp), intent(in) :: tag 3274 3445 3446 INTEGER(iwp), intent(in) :: tag 3447 INTEGER(iwp) :: i 3448 INTEGER(iwp) :: j 3449 INTEGER(iwp) :: k 3450 INTEGER(iwp) :: iif 3451 INTEGER(iwp) :: jjf 3452 INTEGER(iwp) :: kkf 3453 INTEGER(iwp) :: bottomx 3454 INTEGER(iwp) :: bottomy 3455 INTEGER(iwp) :: bottomz 3456 INTEGER(iwp) :: topx 3457 INTEGER(iwp) :: topy 3458 INTEGER(iwp) :: topz 3459 REAL(wp) :: aweight_a 3460 REAL(wp) :: aweight_b 3461 REAL(wp) :: aweight_c 3462 REAL(wp) :: aweight_d 3463 REAL(wp) :: aweight_e 3464 REAL(wp) :: energ 3275 3465 3276 3466 DO k = bdims_rem(3,1)+1, bdims_rem(3,2) … … 3295 3485 aweight_e = 0.0 3296 3486 3297 DO k f = bottomz, topz3298 DO j f = bottomy, topy3299 DO i f = bottomx, topx3300 3301 aweight_a = aweight_a + anterpol3d(k f,jf,if) * &3487 DO kkf = bottomz, topz 3488 DO jjf = bottomy, topy 3489 DO iif = bottomx, topx 3490 3491 aweight_a = aweight_a + anterpol3d(kkf,jjf,iif) * & 3302 3492 (dzf/dzc) * (dyf/dyc) * (dxf/dxc) 3303 3493 3304 3494 3305 energ = ( 0.5 * ( u(k f,jf,if) + u(kf,jf,if+1) ) )**2.0 + &3306 ( 0.5 * ( v(k f,jf,if) + v(kf,jf+1,if) ) )**2.0 + &3307 ( 0.5 * ( w(k f-1,jf,if) + w(kf,jf,if) ) )**2.03495 energ = ( 0.5 * ( u(kkf,jjf,iif) + u(kkf,jjf,iif+1) ) )**2.0 + & 3496 ( 0.5 * ( v(kkf,jjf,iif) + v(kkf,jjf+1,iif) ) )**2.0 + & 3497 ( 0.5 * ( w(kkf-1,jjf,iif) + w(kkf,jjf,iif) ) )**2.0 3308 3498 3309 3499 aweight_b = aweight_b + energ * & 3310 3500 (dzf/dzc) * (dyf/dyc) * (dxf/dxc) 3311 3501 3312 aweight_c = aweight_c + 0.5 * ( u(k f,jf,if) + u(kf,jf,if+1) ) * &3502 aweight_c = aweight_c + 0.5 * ( u(kkf,jjf,iif) + u(kkf,jjf,iif+1) ) * & 3313 3503 (dzf/dzc) * (dyf/dyc) * (dxf/dxc) 3314 3504 3315 aweight_d = aweight_d + 0.5 * ( v(k f,jf,if) + v(kf,jf+1,if) ) * &3505 aweight_d = aweight_d + 0.5 * ( v(kkf,jjf,iif) + v(kkf,jjf+1,iif) ) * & 3316 3506 (dzf/dzc) * (dyf/dyc) * (dxf/dxc) 3317 3507 3318 aweight_e = aweight_e + 0.5 * ( w(k f-1,jf,if) + w(kf,jf,if) ) * &3508 aweight_e = aweight_e + 0.5 * ( w(kkf-1,jjf,iif) + w(kkf,jjf,iif) ) * & 3319 3509 (dzf/dzc) * (dyf/dyc) * (dxf/dxc) 3320 3510 … … 3336 3526 3337 3527 3338 #endif3339 3528 3340 3529 END SUBROUTINE anterpolate_to_crse_e 3530 #endif 3341 3531 END SUBROUTINE vnest_anterpolate_e 3342 3532 3343 3533 SUBROUTINE vnest_init_pegrid_rank 3534 #if defined( __parallel ) 3344 3535 ! Domain decomposition and exchange of grid variables between coarse and fine 3345 3536 ! Given processor coordinates as index f_rnk_lst(pcoord(1), pcoord(2)) … … 3369 3560 IMPLICIT NONE 3370 3561 3371 INTEGER(iwp) :: dest_rnk3372 INTEGER(iwp) :: i !<3562 INTEGER(iwp) :: dest_rnk 3563 INTEGER(iwp) :: i !< 3373 3564 3374 3565 IF (myid == 0) THEN … … 3429 3620 3430 3621 3622 #endif 3431 3623 3432 3624 END SUBROUTINE vnest_init_pegrid_rank … … 3434 3626 3435 3627 SUBROUTINE vnest_init_pegrid_domain 3628 #if defined( __parallel ) 3436 3629 3437 3630 USE control_parameters, & … … 3451 3644 IMPLICIT NONE 3452 3645 3453 INTEGER(iwp) :: dest_rnk 3454 INTEGER(iwp) :: i, j !< 3455 INTEGER(iwp) :: tempx, tempy !< 3456 INTEGER(iwp) :: TYPE_INT_YZ, SIZEOFREAL 3457 INTEGER(iwp) :: MTV_X,MTV_Y,MTV_Z,MTV_RX,MTV_RY,MTV_RZ 3646 INTEGER(iwp) :: dest_rnk 3647 INTEGER(iwp) :: i !< 3648 INTEGER(iwp) :: j !< 3649 INTEGER(iwp) :: tempx 3650 INTEGER(iwp) :: tempy 3651 INTEGER(iwp) :: TYPE_INT_YZ 3652 INTEGER(iwp) :: SIZEOFREAL 3653 INTEGER(iwp) :: MTV_X 3654 INTEGER(iwp) :: MTV_Y 3655 INTEGER(iwp) :: MTV_Z 3656 INTEGER(iwp) :: MTV_RX 3657 INTEGER(iwp) :: MTV_RY 3658 INTEGER(iwp) :: MTV_RZ 3458 3659 3459 3660 ! … … 3692 3893 3693 3894 ENDIF 3694 3895 #endif 3695 3896 END SUBROUTINE vnest_init_pegrid_domain 3696 3897 … … 3698 3899 SUBROUTINE vnest_init_grid 3699 3900 3901 #if defined( __parallel ) 3700 3902 USE arrays_3d, & 3701 3903 ONLY: zu, zw … … 3761 3963 ENDIF 3762 3964 3965 #endif 3763 3966 END SUBROUTINE vnest_init_grid 3764 3967 3765 3968 3766 3969 SUBROUTINE vnest_check_parameters 3767 3768 USE arrays_3d, & 3769 ONLY: zu, zw 3770 3771 USE control_parameters, & 3772 ONLY: coupling_mode 3773 3774 USE indices, & 3775 ONLY: nzt 3776 3777 USE kinds 3778 3779 USE pegrid 3970 #if defined( __parallel ) 3971 3972 USE pegrid, & 3973 ONLY: myid 3780 3974 3781 3975 IMPLICIT NONE 3782 3976 3783 3784 3977 IF (myid==0) PRINT*, '*** vnest: check parameters not implemented yet ***' 3785 3978 3786 3979 #endif 3787 3980 END SUBROUTINE vnest_check_parameters 3788 3981 … … 3790 3983 SUBROUTINE vnest_timestep_sync 3791 3984 3985 #if defined( __parallel ) 3792 3986 USE control_parameters, & 3793 ONLY: coupling_mode, dt_3d, dt_coupling3987 ONLY: coupling_mode, dt_3d, old_dt 3794 3988 3795 3989 USE interfaces … … 3825 4019 !-- Identical timestep for coarse and fine grids 3826 4020 dt_3d = MIN( dtc, dtf ) 3827 !-- Nest coupling at every timestep 3828 dt_coupling = dt_3d 3829 4021 old_dt = dt_3d 4022 #endif 3830 4023 END SUBROUTINE vnest_timestep_sync 3831 4024 3832 4025 SUBROUTINE vnest_deallocate 4026 #if defined( __parallel ) 3833 4027 USE control_parameters, & 3834 4028 ONLY: coupling_mode … … 3846 4040 IF ( ALLOCATED (f2c_dims_fg) ) DEALLOCATE (f2c_dims_fg) 3847 4041 ENDIF 3848 4042 #endif 3849 4043 END SUBROUTINE vnest_deallocate 3850 4044
Note: See TracChangeset
for help on using the changeset viewer.