Changeset 1791
- Timestamp:
- Mar 11, 2016 10:41:25 AM (9 years ago)
- Location:
- palm/trunk/SOURCE
- Files:
-
- 7 edited
-
Makefile (modified) (2 diffs)
-
header.f90 (modified) (6 diffs)
-
pmc_client.f90 (modified) (2 diffs)
-
pmc_handle_communicator.f90 (modified) (15 diffs)
-
pmc_interface.f90 (modified) (60 diffs)
-
pmc_server.f90 (modified) (3 diffs)
-
time_integration.f90 (modified) (3 diffs)
Legend:
- Unmodified
- Added
- Removed
-
palm/trunk/SOURCE/Makefile
r1787 r1791 20 20 # Current revisions: 21 21 # ------------------ 22 # 22 # dependencies of header changed 23 23 # 24 24 # Former revisions: … … 381 381 global_min_max.o: modules.o mod_kinds.o 382 382 header.o: modules.o cpulog.o mod_kinds.o netcdf_interface.o land_surface_model.o\ 383 plant_canopy_model.o pmc_ interface.o radiation_model.o spectrum.o \384 subsidence.o383 plant_canopy_model.o pmc_handle_communicator.o pmc_interface.o \ 384 radiation_model.o spectrum.o subsidence.o 385 385 impact_of_latent_heat.o: modules.o mod_kinds.o 386 386 inflow_turbulence.o: modules.o cpulog.o mod_kinds.o -
palm/trunk/SOURCE/header.f90
r1789 r1791 19 19 ! Current revisions: 20 20 ! ----------------- 21 ! 22 ! 21 ! output of nesting informations of all domains 22 ! 23 23 ! Former revisions: 24 24 ! ----------------- … … 284 284 plant_canopy 285 285 286 USE pmc_handle_communicator, & 287 ONLY: pmc_get_model_info 288 286 289 USE pmc_interface, & 287 ONLY: cpl_id, cpl_parent_id, cpl_name, lower_left_coord_x, & 288 lower_left_coord_y, nested_run, nesting_mode 290 ONLY: nested_run, nesting_mode 289 291 290 292 USE radiation_model_mod, & … … 311 313 312 314 CHARACTER (LEN=26) :: ver_rev !< 315 316 CHARACTER (LEN=32) :: cpl_name !< 313 317 314 318 CHARACTER (LEN=40) :: output_format !< … … 338 342 CHARACTER (LEN=1), DIMENSION(1:3) :: dir = (/ 'x', 'y', 'z' /) !< 339 343 340 INTEGER(iwp) :: av !< 341 INTEGER(iwp) :: bh !< 342 INTEGER(iwp) :: blx !< 343 INTEGER(iwp) :: bly !< 344 INTEGER(iwp) :: bxl !< 345 INTEGER(iwp) :: bxr !< 346 INTEGER(iwp) :: byn !< 347 INTEGER(iwp) :: bys !< 348 INTEGER(iwp) :: ch !< 349 INTEGER(iwp) :: count !< 350 INTEGER(iwp) :: cwx !< 351 INTEGER(iwp) :: cwy !< 352 INTEGER(iwp) :: cxl !< 353 INTEGER(iwp) :: cxr !< 354 INTEGER(iwp) :: cyn !< 355 INTEGER(iwp) :: cys !< 356 INTEGER(iwp) :: dim !< 357 INTEGER(iwp) :: i !< 358 INTEGER(iwp) :: io !< 359 INTEGER(iwp) :: j !< 360 INTEGER(iwp) :: k !< 361 INTEGER(iwp) :: l !< 362 INTEGER(iwp) :: ll !< 363 INTEGER(iwp) :: mpi_type !< 344 INTEGER(iwp) :: av !< 345 INTEGER(iwp) :: bh !< 346 INTEGER(iwp) :: blx !< 347 INTEGER(iwp) :: bly !< 348 INTEGER(iwp) :: bxl !< 349 INTEGER(iwp) :: bxr !< 350 INTEGER(iwp) :: byn !< 351 INTEGER(iwp) :: bys !< 352 INTEGER(iwp) :: ch !< 353 INTEGER(iwp) :: count !< 354 INTEGER(iwp) :: cpl_parent_id !< 355 INTEGER(iwp) :: cwx !< 356 INTEGER(iwp) :: cwy !< 357 INTEGER(iwp) :: cxl !< 358 INTEGER(iwp) :: cxr !< 359 INTEGER(iwp) :: cyn !< 360 INTEGER(iwp) :: cys !< 361 INTEGER(iwp) :: dim !< 362 INTEGER(iwp) :: i !< 363 INTEGER(iwp) :: io !< 364 INTEGER(iwp) :: j !< 365 INTEGER(iwp) :: k !< 366 INTEGER(iwp) :: l !< 367 INTEGER(iwp) :: ll !< 368 INTEGER(iwp) :: mpi_type !< 369 INTEGER(iwp) :: my_cpl_id !< 370 INTEGER(iwp) :: n !< 371 INTEGER(iwp) :: ncpl !< 372 INTEGER(iwp) :: npe_total !< 364 373 365 374 REAL(wp) :: canopy_height !< canopy height (in m) 366 375 REAL(wp) :: cpuseconds_per_simulated_second !< 376 REAL(wp) :: lower_left_coord_x !< x-coordinate of nest domain 377 REAL(wp) :: lower_left_coord_y !< y-coordinate of nest domain 367 378 368 379 ! … … 468 479 !-- Nesting informations 469 480 IF ( nested_run ) THEN 470 WRITE ( io, 600 ) cpl_id, TRIM( cpl_name ), cpl_parent_id, & 471 nesting_mode, lower_left_coord_x, lower_left_coord_y 481 482 WRITE ( io, 600 ) TRIM( nesting_mode ) 483 CALL pmc_get_model_info( ncpl = ncpl, cpl_id = my_cpl_id ) 484 485 DO n = 1, ncpl 486 CALL pmc_get_model_info( request_for_cpl_id = n, cpl_name = cpl_name,& 487 cpl_parent_id = cpl_parent_id, & 488 lower_left_x = lower_left_coord_x, & 489 lower_left_y = lower_left_coord_y, & 490 npe_total = npe_total ) 491 IF ( n == my_cpl_id ) THEN 492 char1 = '*' 493 ELSE 494 char1 = ' ' 495 ENDIF 496 WRITE ( io, 601 ) TRIM( char1 ), n, cpl_parent_id, npe_total, & 497 lower_left_coord_x, lower_left_coord_y, & 498 TRIM( cpl_name ) 499 ENDDO 472 500 ENDIF 473 501 WRITE ( io, 99 ) … … 2443 2471 513 FORMAT (' --> Scalar advection via Wicker-Skamarock-Scheme 5th order ' // & 2444 2472 '+ monotonic adjustment') 2445 600 FORMAT (/' Nesting informations:'/ & 2446 ' Nest id / name: ',I2.2,' / ',A,' (parent id: ',I2.2,')'/ & 2447 ' Nesting mode: ',A/ & 2448 ' Lower left corner coordinates: ','x = ',F8.2,' m, y = ',F8.2,' m'/) 2473 600 FORMAT (/' Nesting informations:'/ & 2474 ' --------------------'/ & 2475 ' Nesting mode: ',A// & 2476 ' Nest id parent number lower left coordinates name'/ & 2477 ' (*=me) id of PEs x (m) y (m)' ) 2478 601 FORMAT (2X,A1,1X,I2.2,6X,I2.2,5X,I5,5X,F8.2,2X,F8.2,5X,A) 2449 2479 2450 2480 END SUBROUTINE header -
palm/trunk/SOURCE/pmc_client.f90
r1787 r1791 20 20 ! Current revisions: 21 21 ! ------------------ 22 ! 22 ! Debug write-statement commented out 23 23 ! 24 24 ! Former revisions: … … 150 150 end do 151 151 152 if(me%model_rank == 0) write(0,'(a,5i6)') 'PMC_ClientInit ',me%model_rank,me%model_npes,me%inter_npes,me%intra_rank152 ! if(me%model_rank == 0) write(0,'(a,5i6)') 'PMC_ClientInit ',me%model_rank,me%model_npes,me%inter_npes,me%intra_rank 153 153 154 154 return -
palm/trunk/SOURCE/pmc_handle_communicator.f90
r1787 r1791 20 20 ! Current revisions: 21 21 ! ------------------ 22 ! 22 ! m_nrofcpl renamed m_ncpl, 23 ! pmc_get_local_model_info renamed pmc_get_model_info, some keywords also 24 ! renamed and some added, 25 ! debug write-statements commented out 23 26 ! 24 27 ! Former revisions: … … 86 89 INTEGER :: m_my_CPL_id !Coupler id of this model 87 90 INTEGER :: m_Parent_id !Coupler id of parent of this model 88 INTEGER :: m_NrOfCpl !Number of Coupler in layout file 89 TYPE(PMC_layout),DIMENSION(PMC_MAX_MODELL) :: m_couplers !Information of all coupler 91 INTEGER :: m_ncpl !Number of Couplers in layout file 92 93 TYPE(PMC_layout),DIMENSION(PMC_MAX_MODELL) :: m_couplers !Information of all couplers 90 94 91 95 ! MPI settings … … 110 114 END INTERFACE pmc_is_rootmodel 111 115 112 INTERFACE PMC_get_local_model_info113 MODULE PROCEDURE PMC_get_local_model_info114 END INTERFACE PMC_get_local_model_info115 116 PUBLIC pmc_get_ local_model_info, pmc_init_model, pmc_is_rootmodel116 INTERFACE pmc_get_model_info 117 MODULE PROCEDURE pmc_get_model_info 118 END INTERFACE pmc_get_model_info 119 120 PUBLIC pmc_get_model_info, pmc_init_model, pmc_is_rootmodel 117 121 118 122 CONTAINS … … 160 164 !-- Calculate start PE of every model 161 165 start_pe(1) = 0 162 DO i = 2, m_n rofcpl+1166 DO i = 2, m_ncpl+1 163 167 start_pe(i) = start_pe(i-1) + m_couplers(i-1)%npe_total 164 168 ENDDO … … 167 171 !-- The number of cores provided with the run must be the same as the 168 172 !-- total sum of cores required by all nest domains 169 IF ( start_pe(m_n rofcpl+1) /= m_world_npes ) THEN173 IF ( start_pe(m_ncpl+1) /= m_world_npes ) THEN 170 174 WRITE ( message_string, '(A,I6,A,I6,A)' ) & 171 175 'nesting-setup requires more MPI procs (', & 172 start_pe(m_n rofcpl+1), ') than provided (',&176 start_pe(m_ncpl+1), ') than provided (', & 173 177 m_world_npes,')' 174 178 CALL message( 'pmc_init_model', 'PA0229', 3, 2, 0, 6, 0 ) … … 202 206 ENDIF 203 207 204 CALL MPI_BCAST( m_n rofcpl, 1,MPI_INTEGER, 0, MPI_COMM_WORLD, istat)205 CALL MPI_BCAST( start_pe, m_n rofcpl+1, MPI_INTEGER, 0, MPI_COMM_WORLD, istat)208 CALL MPI_BCAST( m_ncpl, 1, MPI_INTEGER, 0, MPI_COMM_WORLD, istat) 209 CALL MPI_BCAST( start_pe, m_ncpl+1, MPI_INTEGER, 0, MPI_COMM_WORLD, istat) 206 210 207 211 ! 208 212 !-- Broadcast coupling layout 209 DO i = 1, m_n rofcpl213 DO i = 1, m_ncpl 210 214 CALL MPI_BCAST( m_couplers(i)%name, LEN( m_couplers(i)%name ), MPI_CHARACTER, 0, MPI_COMM_WORLD, istat ) 211 215 CALL MPI_BCAST( m_couplers(i)%id, 1, MPI_INTEGER, 0, MPI_COMM_WORLD, istat ) … … 219 223 ! 220 224 !-- Assign global MPI processes to individual models by setting the couple id 221 DO i = 1, m_n rofcpl225 DO i = 1, m_ncpl 222 226 IF ( m_world_rank >= start_pe(i) .AND. m_world_rank < start_pe(i+1) ) & 223 227 THEN … … 241 245 ! 242 246 !-- Broadcast (from PE 0) the parent id and id of every model 243 DO i = 1, m_n rofcpl247 DO i = 1, m_ncpl 244 248 CALL MPI_BCAST( m_couplers(i)%parent_id, 1, MPI_INTEGER, 0, & 245 249 MPI_COMM_WORLD, istat ) … … 257 261 !-- different colors. 258 262 !-- The grouping was done above with MPI_COMM_SPLIT 259 DO i = 2, m_n rofcpl263 DO i = 2, m_ncpl 260 264 261 265 IF ( m_couplers(i)%parent_id == m_my_cpl_id ) THEN … … 292 296 293 297 clientcount = 0 294 DO i = 2, m_n rofcpl298 DO i = 2, m_ncpl 295 299 IF ( activeserver(i) == 1 ) THEN 296 300 clientcount = clientcount + 1 … … 321 325 322 326 ! 323 !-- Make module private variables available to palm 324 SUBROUTINE pmc_get_local_model_info( my_cpl_id, my_cpl_parent_id, cpl_name, & 325 npe_total, lower_left_x, lower_left_y ) 327 !-- Provide module private variables of the pmc for PALM 328 SUBROUTINE pmc_get_model_info( cpl_id, cpl_name, cpl_parent_id, & 329 lower_left_x, lower_left_y, ncpl, npe_total,& 330 request_for_cpl_id ) 326 331 327 332 USE kinds … … 330 335 331 336 CHARACTER(LEN=*), INTENT(OUT), OPTIONAL :: cpl_name 332 INTEGER, INTENT(OUT), OPTIONAL :: my_cpl_id 333 INTEGER, INTENT(OUT), OPTIONAL :: my_cpl_parent_id 334 INTEGER, INTENT(OUT), OPTIONAL :: npe_total 335 REAL(wp), INTENT(OUT), OPTIONAL :: lower_left_x 336 REAL(wp), INTENT(OUT), OPTIONAL :: lower_left_y 337 338 IF ( PRESENT( my_cpl_id ) ) my_cpl_id = m_my_cpl_id 339 IF ( PRESENT( my_cpl_parent_id ) ) my_cpl_parent_id = m_couplers(my_cpl_id)%parent_id 340 IF ( PRESENT( cpl_name ) ) cpl_name = m_couplers(my_cpl_id)%name 341 IF ( PRESENT( npe_total ) ) npe_total = m_couplers(my_cpl_id)%npe_total 342 IF ( PRESENT( lower_left_x ) ) lower_left_x = m_couplers(my_cpl_id)%lower_left_x 343 IF ( PRESENT( lower_left_y ) ) lower_left_y = m_couplers(my_cpl_id)%lower_left_y 344 345 END SUBROUTINE pmc_get_local_model_info 337 338 INTEGER, INTENT(IN), OPTIONAL :: request_for_cpl_id 339 340 INTEGER, INTENT(OUT), OPTIONAL :: cpl_id 341 INTEGER, INTENT(OUT), OPTIONAL :: cpl_parent_id 342 INTEGER, INTENT(OUT), OPTIONAL :: ncpl 343 INTEGER, INTENT(OUT), OPTIONAL :: npe_total 344 345 INTEGER :: requested_cpl_id 346 347 REAL(wp), INTENT(OUT), OPTIONAL :: lower_left_x 348 REAL(wp), INTENT(OUT), OPTIONAL :: lower_left_y 349 350 ! 351 !-- Set the requested coupler id 352 IF ( PRESENT( request_for_cpl_id ) ) THEN 353 requested_cpl_id = request_for_cpl_id 354 ! 355 !-- Check for allowed range of values 356 IF ( requested_cpl_id < 1 .OR. requested_cpl_id > m_ncpl ) RETURN 357 ELSE 358 requested_cpl_id = m_my_cpl_id 359 ENDIF 360 361 ! 362 !-- Return the requested information 363 IF ( PRESENT( cpl_id ) ) THEN 364 cpl_id = requested_cpl_id 365 ENDIF 366 IF ( PRESENT( cpl_parent_id ) ) THEN 367 cpl_parent_id = m_couplers(requested_cpl_id)%parent_id 368 ENDIF 369 IF ( PRESENT( cpl_name ) ) THEN 370 cpl_name = m_couplers(requested_cpl_id)%name 371 ENDIF 372 IF ( PRESENT( ncpl ) ) THEN 373 ncpl = m_ncpl 374 ENDIF 375 IF ( PRESENT( npe_total ) ) THEN 376 npe_total = m_couplers(requested_cpl_id)%npe_total 377 ENDIF 378 IF ( PRESENT( lower_left_x ) ) THEN 379 lower_left_x = m_couplers(requested_cpl_id)%lower_left_x 380 ENDIF 381 IF ( PRESENT( lower_left_y ) ) THEN 382 lower_left_y = m_couplers(requested_cpl_id)%lower_left_y 383 ENDIF 384 385 END SUBROUTINE pmc_get_model_info 346 386 347 387 … … 364 404 365 405 INTEGER, INTENT(INOUT) :: pmc_status 366 INTEGER :: i, istat , iunit406 INTEGER :: i, istat 367 407 368 408 TYPE(pmc_layout), DIMENSION(pmc_max_modell) :: domain_layouts … … 374 414 !-- Initialize some coupling variables 375 415 domain_layouts(1:pmc_max_modell)%id = -1 376 m_nrofcpl = 0 377 iunit = 345 416 m_ncpl = 0 378 417 379 418 pmc_status = pmc_status_ok … … 412 451 !-- Get the number of nested models given in the nestpar-NAMELIST 413 452 DO i = 1, pmc_max_modell 414 415 IF ( m_couplers(i)%id /= -1 .AND. i <= pmc_max_modell ) THEN 416 WRITE ( 0, '(A,A,1X,3I7,1X,2F10.2)' ) 'Set up Model ', & 417 TRIM( m_couplers(i)%name ), m_couplers(i)%id, & 418 m_couplers(i)%Parent_id, m_couplers(i)%npe_total, & 419 m_couplers(i)%lower_left_x, & 420 m_couplers(i)%lower_left_y 421 ELSE 422 ! 423 !-- When id=-1 is found for the first time, the list of domains is 424 !-- finished (or latest after pmc_max_modell entries 425 m_nrofcpl = i - 1 426 EXIT 453 ! 454 !-- When id=-1 is found for the first time, the list of domains is finished 455 IF ( m_couplers(i)%id == -1 .OR. i == pmc_max_modell ) THEN 456 IF ( m_couplers(i)%id == -1 ) THEN 457 m_ncpl = i - 1 458 EXIT 459 ELSE 460 m_ncpl = pmc_max_modell 461 ENDIF 427 462 ENDIF 428 463 -
palm/trunk/SOURCE/pmc_interface.f90
r1784 r1791 20 20 ! Current revisions: 21 21 ! ------------------ 22 ! 22 ! routine pmci_update_new removed, 23 ! pmc_get_local_model_info renamed pmc_get_model_info, some keywords also 24 ! renamed, 25 ! filling up redundant ghost points introduced, 26 ! some index bound variables renamed, 27 ! further formatting cleanup 23 28 ! 24 29 ! Former revisions: … … 111 116 112 117 USE pmc_handle_communicator, & 113 ONLY: pmc_get_ local_model_info, pmc_init_model, pmc_is_rootmodel,&118 ONLY: pmc_get_model_info, pmc_init_model, pmc_is_rootmodel, & 114 119 pmc_no_namelist_found, pmc_server_for_client 115 120 … … 128 133 IMPLICIT NONE 129 134 130 !-- TO_DO: a lot of lines (including comments) in this file exceed the 80 char 131 !-- limit. Try to reduce as much as possible 132 133 !-- TO_DO: are all of these variables following now really PUBLIC? 134 !-- Klaus and I guess they are not 135 PRIVATE !: Note that the default publicity is here set to private. 135 PRIVATE 136 136 137 137 ! 138 138 !-- Constants 139 INTEGER(iwp), PARAMETER , PUBLIC:: client_to_server = 2 !:140 INTEGER(iwp), PARAMETER , PUBLIC:: server_to_client = 1 !:139 INTEGER(iwp), PARAMETER :: client_to_server = 2 !: 140 INTEGER(iwp), PARAMETER :: server_to_client = 1 !: 141 141 142 142 ! 143 143 !-- Coupler setup 144 INTEGER(iwp), PUBLIC,SAVE :: cpl_id = 1 !:145 CHARACTER(LEN=32), PUBLIC,SAVE :: cpl_name !:146 INTEGER(iwp), PUBLIC,SAVE :: cpl_npe_total !:147 INTEGER(iwp), PUBLIC,SAVE :: cpl_parent_id !:144 INTEGER(iwp), SAVE :: cpl_id = 1 !: 145 CHARACTER(LEN=32), SAVE :: cpl_name !: 146 INTEGER(iwp), SAVE :: cpl_npe_total !: 147 INTEGER(iwp), SAVE :: cpl_parent_id !: 148 148 149 149 ! 150 150 !-- Control parameters, will be made input parameters later 151 CHARACTER(LEN=7), PUBLIC, SAVE :: nesting_mode = 'two-way' !: steering parameter for one- or two-way nesting 152 153 LOGICAL, PUBLIC, SAVE :: nested_run = .FALSE. !: general switch if nested run or not 154 155 REAL(wp), PUBLIC, SAVE :: anterp_relax_length_l = -1.0_wp !: 156 REAL(wp), PUBLIC, SAVE :: anterp_relax_length_r = -1.0_wp !: 157 REAL(wp), PUBLIC, SAVE :: anterp_relax_length_s = -1.0_wp !: 158 REAL(wp), PUBLIC, SAVE :: anterp_relax_length_n = -1.0_wp !: 159 REAL(wp), PUBLIC, SAVE :: anterp_relax_length_t = -1.0_wp !: 151 CHARACTER(LEN=7), SAVE :: nesting_mode = 'two-way' !: steering parameter 152 !: for 1- or 2-way nesting 153 154 LOGICAL, SAVE :: nested_run = .FALSE. !: general switch 155 156 REAL(wp), SAVE :: anterp_relax_length_l = -1.0_wp !: 157 REAL(wp), SAVE :: anterp_relax_length_r = -1.0_wp !: 158 REAL(wp), SAVE :: anterp_relax_length_s = -1.0_wp !: 159 REAL(wp), SAVE :: anterp_relax_length_n = -1.0_wp !: 160 REAL(wp), SAVE :: anterp_relax_length_t = -1.0_wp !: 160 161 161 162 ! 162 163 !-- Geometry 163 REAL(wp), PUBLIC, SAVE:: area_t !:164 REAL(wp), SAVE :: area_t !: 164 165 REAL(wp), SAVE, DIMENSION(:), ALLOCATABLE :: coord_x !: 165 166 REAL(wp), SAVE, DIMENSION(:), ALLOCATABLE :: coord_y !: 166 REAL(wp), PUBLIC, SAVE:: lower_left_coord_x !:167 REAL(wp), PUBLIC, SAVE:: lower_left_coord_y !:167 REAL(wp), SAVE :: lower_left_coord_x !: 168 REAL(wp), SAVE :: lower_left_coord_y !: 168 169 169 170 ! 170 171 !-- Client coarse data arrays 171 REAL(wp), SAVE, DIMENSION(:,:,:), ALLOCATABLE, TARGET, PUBLIC :: ec !: 172 REAL(wp), SAVE, DIMENSION(:,:,:), ALLOCATABLE, TARGET, PUBLIC :: ptc !: 173 REAL(wp), SAVE, DIMENSION(:,:,:), ALLOCATABLE, TARGET, PUBLIC :: uc !: 174 REAL(wp), SAVE, DIMENSION(:,:,:), ALLOCATABLE, TARGET, PUBLIC :: vc !: 175 REAL(wp), SAVE, DIMENSION(:,:,:), ALLOCATABLE, TARGET, PUBLIC :: wc !: 176 REAL(wp), SAVE, DIMENSION(:,:,:), ALLOCATABLE, TARGET, PUBLIC :: qc !: 177 178 INTEGER(iwp), DIMENSION(5) :: coarse_bound !: 179 REAL(wp), PUBLIC, SAVE :: xexl !: 180 REAL(wp), PUBLIC, SAVE :: xexr !: 181 REAL(wp), PUBLIC, SAVE :: yexs !: 182 REAL(wp), PUBLIC, SAVE :: yexn !: 183 REAL(wp), PUBLIC, SAVE, DIMENSION(:,:), ALLOCATABLE :: tkefactor_l !: 184 REAL(wp), PUBLIC, SAVE, DIMENSION(:,:), ALLOCATABLE :: tkefactor_n !: 185 REAL(wp), PUBLIC, SAVE, DIMENSION(:,:), ALLOCATABLE :: tkefactor_r !: 186 REAL(wp), PUBLIC, SAVE, DIMENSION(:,:), ALLOCATABLE :: tkefactor_s !: 187 REAL(wp), PUBLIC, SAVE, DIMENSION(:,:), ALLOCATABLE :: tkefactor_t !: 188 189 ! 190 !-- Client interpolation coefficients and client-array indices to be precomputed and stored. 191 INTEGER(iwp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:) :: ico !: 192 INTEGER(iwp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:) :: icu !: 193 INTEGER(iwp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:) :: jco !: 194 INTEGER(iwp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:) :: jcv !: 195 INTEGER(iwp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:) :: kco !: 196 INTEGER(iwp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:) :: kcw !: 197 REAL(wp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:) :: r1xo !: 198 REAL(wp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:) :: r2xo !: 199 REAL(wp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:) :: r1xu !: 200 REAL(wp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:) :: r2xu !: 201 REAL(wp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:) :: r1yo !: 202 REAL(wp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:) :: r2yo !: 203 REAL(wp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:) :: r1yv !: 204 REAL(wp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:) :: r2yv !: 205 REAL(wp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:) :: r1zo !: 206 REAL(wp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:) :: r2zo !: 207 REAL(wp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:) :: r1zw !: 208 REAL(wp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:) :: r2zw !: 209 210 ! 211 !-- Client index arrays and log-ratio arrays for the log-law near-wall corrections. 212 !-- These are not truly 3-D arrays but multiply 2-D arrays. 213 INTEGER(iwp), PUBLIC, SAVE :: ncorr !: ncorr is the 4th dimension of the log_ratio-arrays 214 INTEGER(iwp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:,:,:) :: logc_u_l !: 215 INTEGER(iwp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:,:,:) :: logc_u_n !: 216 INTEGER(iwp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:,:,:) :: logc_u_r !: 217 INTEGER(iwp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:,:,:) :: logc_u_s !: 218 INTEGER(iwp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:,:,:) :: logc_v_l !: 219 INTEGER(iwp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:,:,:) :: logc_v_n !: 220 INTEGER(iwp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:,:,:) :: logc_v_r !: 221 INTEGER(iwp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:,:,:) :: logc_v_s !: 222 INTEGER(iwp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:,:,:) :: logc_w_l !: 223 INTEGER(iwp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:,:,:) :: logc_w_n !: 224 INTEGER(iwp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:,:,:) :: logc_w_r !: 225 INTEGER(iwp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:,:,:) :: logc_w_s !: 226 REAL(wp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:,:,:,:) :: logc_ratio_u_l !: 227 REAL(wp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:,:,:,:) :: logc_ratio_u_n !: 228 REAL(wp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:,:,:,:) :: logc_ratio_u_r !: 229 REAL(wp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:,:,:,:) :: logc_ratio_u_s !: 230 REAL(wp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:,:,:,:) :: logc_ratio_v_l !: 231 REAL(wp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:,:,:,:) :: logc_ratio_v_n !: 232 REAL(wp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:,:,:,:) :: logc_ratio_v_r !: 233 REAL(wp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:,:,:,:) :: logc_ratio_v_s !: 234 REAL(wp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:,:,:,:) :: logc_ratio_w_l !: 235 REAL(wp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:,:,:,:) :: logc_ratio_w_n !: 236 REAL(wp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:,:,:,:) :: logc_ratio_w_r !: 237 REAL(wp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:,:,:,:) :: logc_ratio_w_s !: 172 INTEGER(iwp), DIMENSION(5) :: coarse_bound !: 173 174 REAL(wp), SAVE :: xexl !: 175 REAL(wp), SAVE :: xexr !: 176 REAL(wp), SAVE :: yexs !: 177 REAL(wp), SAVE :: yexn !: 178 REAL(wp), SAVE, DIMENSION(:,:), ALLOCATABLE :: tkefactor_l !: 179 REAL(wp), SAVE, DIMENSION(:,:), ALLOCATABLE :: tkefactor_n !: 180 REAL(wp), SAVE, DIMENSION(:,:), ALLOCATABLE :: tkefactor_r !: 181 REAL(wp), SAVE, DIMENSION(:,:), ALLOCATABLE :: tkefactor_s !: 182 REAL(wp), SAVE, DIMENSION(:,:), ALLOCATABLE :: tkefactor_t !: 183 184 REAL(wp), SAVE, DIMENSION(:,:,:), ALLOCATABLE, TARGET :: ec !: 185 REAL(wp), SAVE, DIMENSION(:,:,:), ALLOCATABLE, TARGET :: ptc !: 186 REAL(wp), SAVE, DIMENSION(:,:,:), ALLOCATABLE, TARGET :: uc !: 187 REAL(wp), SAVE, DIMENSION(:,:,:), ALLOCATABLE, TARGET :: vc !: 188 REAL(wp), SAVE, DIMENSION(:,:,:), ALLOCATABLE, TARGET :: wc !: 189 REAL(wp), SAVE, DIMENSION(:,:,:), ALLOCATABLE, TARGET :: qc !: 190 191 ! 192 !-- Client interpolation coefficients and client-array indices to be precomputed 193 !-- and stored. 194 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: ico !: 195 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: icu !: 196 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: jco !: 197 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: jcv !: 198 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: kco !: 199 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: kcw !: 200 REAL(wp), SAVE, ALLOCATABLE, DIMENSION(:) :: r1xo !: 201 REAL(wp), SAVE, ALLOCATABLE, DIMENSION(:) :: r2xo !: 202 REAL(wp), SAVE, ALLOCATABLE, DIMENSION(:) :: r1xu !: 203 REAL(wp), SAVE, ALLOCATABLE, DIMENSION(:) :: r2xu !: 204 REAL(wp), SAVE, ALLOCATABLE, DIMENSION(:) :: r1yo !: 205 REAL(wp), SAVE, ALLOCATABLE, DIMENSION(:) :: r2yo !: 206 REAL(wp), SAVE, ALLOCATABLE, DIMENSION(:) :: r1yv !: 207 REAL(wp), SAVE, ALLOCATABLE, DIMENSION(:) :: r2yv !: 208 REAL(wp), SAVE, ALLOCATABLE, DIMENSION(:) :: r1zo !: 209 REAL(wp), SAVE, ALLOCATABLE, DIMENSION(:) :: r2zo !: 210 REAL(wp), SAVE, ALLOCATABLE, DIMENSION(:) :: r1zw !: 211 REAL(wp), SAVE, ALLOCATABLE, DIMENSION(:) :: r2zw !: 212 213 ! 214 !-- Client index arrays and log-ratio arrays for the log-law near-wall 215 !-- corrections. These are not truly 3-D arrays but multiple 2-D arrays. 216 INTEGER(iwp), SAVE :: ncorr !: 4th dimension of the log_ratio-arrays 217 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:,:,:) :: logc_u_l !: 218 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:,:,:) :: logc_u_n !: 219 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:,:,:) :: logc_u_r !: 220 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:,:,:) :: logc_u_s !: 221 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:,:,:) :: logc_v_l !: 222 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:,:,:) :: logc_v_n !: 223 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:,:,:) :: logc_v_r !: 224 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:,:,:) :: logc_v_s !: 225 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:,:,:) :: logc_w_l !: 226 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:,:,:) :: logc_w_n !: 227 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:,:,:) :: logc_w_r !: 228 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:,:,:) :: logc_w_s !: 229 REAL(wp), SAVE, ALLOCATABLE, DIMENSION(:,:,:,:) :: logc_ratio_u_l !: 230 REAL(wp), SAVE, ALLOCATABLE, DIMENSION(:,:,:,:) :: logc_ratio_u_n !: 231 REAL(wp), SAVE, ALLOCATABLE, DIMENSION(:,:,:,:) :: logc_ratio_u_r !: 232 REAL(wp), SAVE, ALLOCATABLE, DIMENSION(:,:,:,:) :: logc_ratio_u_s !: 233 REAL(wp), SAVE, ALLOCATABLE, DIMENSION(:,:,:,:) :: logc_ratio_v_l !: 234 REAL(wp), SAVE, ALLOCATABLE, DIMENSION(:,:,:,:) :: logc_ratio_v_n !: 235 REAL(wp), SAVE, ALLOCATABLE, DIMENSION(:,:,:,:) :: logc_ratio_v_r !: 236 REAL(wp), SAVE, ALLOCATABLE, DIMENSION(:,:,:,:) :: logc_ratio_v_s !: 237 REAL(wp), SAVE, ALLOCATABLE, DIMENSION(:,:,:,:) :: logc_ratio_w_l !: 238 REAL(wp), SAVE, ALLOCATABLE, DIMENSION(:,:,:,:) :: logc_ratio_w_n !: 239 REAL(wp), SAVE, ALLOCATABLE, DIMENSION(:,:,:,:) :: logc_ratio_w_r !: 240 REAL(wp), SAVE, ALLOCATABLE, DIMENSION(:,:,:,:) :: logc_ratio_w_s !: 238 241 239 242 ! 240 243 !-- Upper bounds for k in anterpolation. 241 INTEGER(iwp), PUBLIC, SAVE :: kceu !:242 INTEGER(iwp), PUBLIC, SAVE :: kcew !:244 INTEGER(iwp), SAVE :: kctu !: 245 INTEGER(iwp), SAVE :: kctw !: 243 246 244 247 ! 245 248 !-- Upper bound for k in log-law correction in interpolation. 246 INTEGER(iwp), PUBLIC,SAVE :: nzt_topo_nestbc_l !:247 INTEGER(iwp), PUBLIC,SAVE :: nzt_topo_nestbc_n !:248 INTEGER(iwp), PUBLIC,SAVE :: nzt_topo_nestbc_r !:249 INTEGER(iwp), PUBLIC,SAVE :: nzt_topo_nestbc_s !:249 INTEGER(iwp), SAVE :: nzt_topo_nestbc_l !: 250 INTEGER(iwp), SAVE :: nzt_topo_nestbc_n !: 251 INTEGER(iwp), SAVE :: nzt_topo_nestbc_r !: 252 INTEGER(iwp), SAVE :: nzt_topo_nestbc_s !: 250 253 251 254 ! 252 255 !-- Number of ghost nodes in coarse-grid arrays for i and j in anterpolation. 253 INTEGER(iwp), PUBLIC,SAVE :: nhll !:254 INTEGER(iwp), PUBLIC,SAVE :: nhlr !:255 INTEGER(iwp), PUBLIC,SAVE :: nhls !:256 INTEGER(iwp), PUBLIC,SAVE :: nhln !:256 INTEGER(iwp), SAVE :: nhll !: 257 INTEGER(iwp), SAVE :: nhlr !: 258 INTEGER(iwp), SAVE :: nhls !: 259 INTEGER(iwp), SAVE :: nhln !: 257 260 258 261 ! 259 262 !-- Spatial under-relaxation coefficients for anterpolation. 260 REAL(wp), PUBLIC,SAVE, ALLOCATABLE, DIMENSION(:) :: frax !:261 REAL(wp), PUBLIC,SAVE, ALLOCATABLE, DIMENSION(:) :: fray !:262 REAL(wp), PUBLIC,SAVE, ALLOCATABLE, DIMENSION(:) :: fraz !:263 REAL(wp), SAVE, ALLOCATABLE, DIMENSION(:) :: frax !: 264 REAL(wp), SAVE, ALLOCATABLE, DIMENSION(:) :: fray !: 265 REAL(wp), SAVE, ALLOCATABLE, DIMENSION(:) :: fraz !: 263 266 264 267 ! 265 268 !-- Client-array indices to be precomputed and stored for anterpolation. 266 INTEGER(iwp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:) :: iflu !: 267 INTEGER(iwp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:) :: ifuu !: 268 INTEGER(iwp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:) :: iflo !: 269 INTEGER(iwp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:) :: ifuo !: 270 INTEGER(iwp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:) :: jflv !: 271 INTEGER(iwp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:) :: jfuv !: 272 INTEGER(iwp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:) :: jflo !: 273 INTEGER(iwp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:) :: jfuo !: 274 INTEGER(iwp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:) :: kflw !: 275 INTEGER(iwp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:) :: kfuw !: 276 INTEGER(iwp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:) :: kflo !: 277 INTEGER(iwp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:) :: kfuo !: 278 279 ! 280 !-- Module private variables. 269 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: iflu !: 270 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: ifuu !: 271 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: iflo !: 272 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: ifuo !: 273 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: jflv !: 274 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: jfuv !: 275 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: jflo !: 276 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: jfuo !: 277 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: kflw !: 278 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: kfuw !: 279 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: kflo !: 280 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: kfuo !: 281 281 282 INTEGER(iwp), DIMENSION(3) :: define_coarse_grid_int !: 282 283 REAL(wp), DIMENSION(7) :: define_coarse_grid_real !: … … 301 302 END TYPE coarsegrid_def 302 303 303 TYPE(coarsegrid_def), SAVE 304 TYPE(coarsegrid_def), SAVE :: cg !: 304 305 305 306 … … 340 341 END INTERFACE pmci_set_swaplevel 341 342 342 INTERFACE pmci_update_new343 MODULE PROCEDURE pmci_update_new344 END INTERFACE345 343 PUBLIC anterp_relax_length_l, anterp_relax_length_r, & 344 anterp_relax_length_s, anterp_relax_length_n, & 345 anterp_relax_length_t, client_to_server, cpl_id, nested_run, & 346 nesting_mode, server_to_client 346 347 PUBLIC pmci_client_datatrans 347 348 PUBLIC pmci_client_initialize … … 354 355 PUBLIC pmci_server_synchronize 355 356 PUBLIC pmci_set_swaplevel 356 PUBLIC pmci_update_new357 357 358 358 … … 393 393 !-- Get some variables required by the pmc-interface (and in some cases in the 394 394 !-- PALM code out of the pmci) out of the pmc-core 395 CALL pmc_get_local_model_info( my_cpl_id = cpl_id, & 396 my_cpl_parent_id = cpl_parent_id, & 397 cpl_name = cpl_name, & 398 npe_total = cpl_npe_total, & 399 lower_left_x = lower_left_coord_x, & 400 lower_left_y = lower_left_coord_y ) 395 CALL pmc_get_model_info( cpl_id = cpl_id, cpl_parent_id = cpl_parent_id, & 396 cpl_name = cpl_name, npe_total = cpl_npe_total, & 397 lower_left_x = lower_left_coord_x, & 398 lower_left_y = lower_left_coord_y ) 401 399 ! 402 400 !-- Set the steering switch which tells the models that they are nested (of 403 !-- course the root domain (cpl_id = 1 401 !-- course the root domain (cpl_id = 1) is not nested) 404 402 IF ( cpl_id >= 2 ) THEN 405 403 nest_domain = .TRUE. … … 435 433 436 434 CALL location_message( 'setup the nested model configuration', .FALSE. ) 437 CALL pmci_setup_coordinates !: Compute absolute coordinates valid for all models 438 CALL pmci_setup_client !: Initialize PMC Client (Must be called before pmc_setup_server) 439 CALL pmci_setup_server !: Initialize PMC Server 435 ! 436 !-- Compute absolute coordinates for all models 437 CALL pmci_setup_coordinates 438 ! 439 !-- Initialize the client (must be called before pmc_setup_server) 440 CALL pmci_setup_client 441 ! 442 !-- Initialize PMC Server 443 CALL pmci_setup_server 444 440 445 CALL location_message( 'finished', .TRUE. ) 441 446 … … 449 454 IMPLICIT NONE 450 455 451 INTEGER(iwp) :: client_id !: 452 INTEGER(iwp) :: ierr !: 453 INTEGER(iwp) :: i !: 454 INTEGER(iwp) :: j !: 455 INTEGER(iwp) :: k !: 456 INTEGER(iwp) :: m !: 457 INTEGER(iwp) :: nomatch !: 458 INTEGER(iwp) :: nx_cl !: 459 INTEGER(iwp) :: ny_cl !: 460 INTEGER(iwp) :: nz_cl !: 461 INTEGER(iwp), DIMENSION(5) :: val !: 462 REAL(wp), DIMENSION(1) :: fval !: 463 REAL(wp) :: dx_cl !: 464 REAL(wp) :: dy_cl !: 465 REAL(wp) :: xez !: 466 REAL(wp) :: yez !: 467 CHARACTER(len=32) :: mychannel 468 CHARACTER(len=32) :: myname 456 CHARACTER(LEN=32) :: myname 457 458 INTEGER(iwp) :: client_id !: 459 INTEGER(iwp) :: ierr !: 460 INTEGER(iwp) :: i !: 461 INTEGER(iwp) :: j !: 462 INTEGER(iwp) :: k !: 463 INTEGER(iwp) :: m !: 464 INTEGER(iwp) :: nomatch !: 465 INTEGER(iwp) :: nx_cl !: 466 INTEGER(iwp) :: ny_cl !: 467 INTEGER(iwp) :: nz_cl !: 468 469 INTEGER(iwp), DIMENSION(5) :: val !: 470 471 REAL(wp) :: dx_cl !: 472 REAL(wp) :: dy_cl !: 473 REAL(wp) :: xez !: 474 REAL(wp) :: yez !: 475 476 REAL(wp), DIMENSION(1) :: fval !: 477 469 478 REAL(wp), DIMENSION(:), ALLOCATABLE :: cl_coord_x !: 470 479 REAL(wp), DIMENSION(:), ALLOCATABLE :: cl_coord_y !: … … 472 481 473 482 ! 474 ! Initialize the PMCserver483 ! Initialize the pmc server 475 484 CALL pmc_serverinit 476 485 … … 478 487 !-- Get coordinates from all clients 479 488 DO m = 1, SIZE( pmc_server_for_client ) - 1 489 480 490 client_id = pmc_server_for_client(m) 481 491 IF ( myid == 0 ) THEN … … 510 520 CALL pmc_recv_from_client( client_id, cl_coord_y, SIZE( cl_coord_y ),& 511 521 0, 12, ierr ) 512 WRITE ( 0, * ) 'receive from pmc Client ', client_id, nx_cl, ny_cl522 ! WRITE ( 0, * ) 'receive from pmc Client ', client_id, nx_cl, ny_cl 513 523 514 524 define_coarse_grid_real(1) = lower_left_coord_x … … 529 539 xez = ( nbgp + 1 ) * dx 530 540 yez = ( nbgp + 1 ) * dy 531 IF ( cl_coord_x(0) < define_coarse_grid_real(1) + xez ) nomatch = 1 532 IF ( cl_coord_x(nx_cl + 1) > define_coarse_grid_real(5) - xez ) nomatch = 1 533 IF ( cl_coord_y(0) < define_coarse_grid_real(2) + yez ) nomatch = 1 534 IF ( cl_coord_y(ny_cl + 1) > define_coarse_grid_real(6) - yez ) nomatch = 1 541 IF ( ( cl_coord_x(0) < define_coarse_grid_real(1) + xez ) .OR. & 542 ( cl_coord_x(nx_cl+1) > define_coarse_grid_real(5) - xez ) .OR. & 543 ( cl_coord_y(0) < define_coarse_grid_real(2) + yez ) .OR. & 544 ( cl_coord_y(ny_cl+1) > define_coarse_grid_real(6) - yez ) ) & 545 THEN 546 nomatch = 1 547 ENDIF 535 548 536 549 DEALLOCATE( cl_coord_x ) … … 539 552 ! 540 553 !-- Send coarse grid information to client 541 CALL pmc_send_to_client( client_id, Define_coarse_grid_real, &542 SIZE( define_coarse_grid_real), 0,&543 21,ierr )544 CALL pmc_send_to_client( client_id, Define_coarse_grid_int, 3, 0, &554 CALL pmc_send_to_client( client_id, define_coarse_grid_real, & 555 SIZE( define_coarse_grid_real ), 0, 21, & 556 ierr ) 557 CALL pmc_send_to_client( client_id, define_coarse_grid_int, 3, 0, & 545 558 22, ierr ) 546 559 547 560 ! 548 561 !-- Send local grid to client 549 CALL pmc_send_to_client( client_id, coord_x, nx+1+2*nbgp, 0, 24, ierr ) 550 CALL pmc_send_to_client( client_id, coord_y, ny+1+2*nbgp, 0, 25, ierr ) 562 CALL pmc_send_to_client( client_id, coord_x, nx+1+2*nbgp, 0, 24, & 563 ierr ) 564 CALL pmc_send_to_client( client_id, coord_y, ny+1+2*nbgp, 0, 25, & 565 ierr ) 551 566 552 567 ! 553 568 !-- Also send the dzu-, dzw-, zu- and zw-arrays here 554 CALL pmc_send_to_client( client_id, dzu, nz_cl +1, 0, 26, ierr )555 CALL pmc_send_to_client( client_id, dzw, nz_cl +1, 0, 27, ierr )556 CALL pmc_send_to_client( client_id, zu, nz_cl +2, 0, 28, ierr )557 CALL pmc_send_to_client( client_id, zw, nz_cl +2, 0, 29, ierr )569 CALL pmc_send_to_client( client_id, dzu, nz_cl+1, 0, 26, ierr ) 570 CALL pmc_send_to_client( client_id, dzw, nz_cl+1, 0, 27, ierr ) 571 CALL pmc_send_to_client( client_id, zu, nz_cl+2, 0, 28, ierr ) 572 CALL pmc_send_to_client( client_id, zw, nz_cl+2, 0, 29, ierr ) 558 573 559 574 ENDIF … … 567 582 568 583 CALL MPI_BCAST( nz_cl, 1, MPI_INTEGER, 0, comm2d, ierr ) 569 584 585 ! 586 !-- TO_DO: Klaus: please give a comment what is done here 570 587 CALL pmci_create_index_list 571 588 572 589 ! 573 590 !-- Include couple arrays into server content 591 !-- TO_DO: Klaus: please give a more meaningful comment 574 592 CALL pmc_s_clear_next_array_list 575 593 DO WHILE ( pmc_s_getnextarray( client_id, myname ) ) … … 587 605 IMPLICIT NONE 588 606 607 INTEGER(iwp) :: i !: 608 INTEGER(iwp) :: ic !: 609 INTEGER(iwp) :: ierr !: 610 INTEGER(iwp) :: j !: 611 INTEGER(iwp) :: k !: 612 INTEGER(iwp) :: m !: 613 INTEGER(iwp) :: n !: 614 INTEGER(iwp) :: npx !: 615 INTEGER(iwp) :: npy !: 616 INTEGER(iwp) :: nrx !: 617 INTEGER(iwp) :: nry !: 618 INTEGER(iwp) :: px !: 619 INTEGER(iwp) :: py !: 620 INTEGER(iwp) :: server_pe !: 621 622 INTEGER(iwp), DIMENSION(2) :: scoord !: 623 INTEGER(iwp), DIMENSION(2) :: size_of_array !: 624 589 625 INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE :: coarse_bound_all !: 590 INTEGER(iwp) :: i !:591 INTEGER(iwp) :: ic !:592 INTEGER(iwp) :: ierr !:593 626 INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE :: index_list !: 594 INTEGER(iwp) :: j !:595 INTEGER(iwp) :: k !:596 INTEGER(iwp) :: m !:597 INTEGER(iwp) :: n !:598 INTEGER(iwp) :: npx !:599 INTEGER(iwp) :: npy !:600 INTEGER(iwp) :: nrx !:601 INTEGER(iwp) :: nry !:602 INTEGER(iwp) :: px !:603 INTEGER(iwp) :: py !:604 INTEGER(iwp), DIMENSION(2) :: scoord !:605 INTEGER(iwp) :: server_pe !:606 INTEGER(iwp), DIMENSION(2) :: size_of_array !:607 608 627 609 628 IF ( myid == 0 ) THEN 629 !-- TO_DO: Klaus: give more specific comment what size_of_array stands for 610 630 CALL pmc_recv_from_client( client_id, size_of_array, 2, 0, 40, ierr ) 611 631 ALLOCATE( coarse_bound_all(size_of_array(1),size_of_array(2)) ) … … 628 648 CALL MPI_COMM_SIZE( comm1dx, npx, ierr ) 629 649 CALL MPI_COMM_SIZE( comm1dy, npy, ierr ) 630 631 nrx = nxr - nxl + 1 ! +1 in index because FORTRAN indexing starts with 1, palm with 0 650 ! 651 !-- The +1 in index is because PALM starts with nx=0 652 nrx = nxr - nxl + 1 632 653 nry = nyn - nys + 1 633 654 ic = 0 634 DO k = 1, size_of_array(2) ! loop over all client PEs 635 DO j = coarse_bound_all(3,k), coarse_bound_all(4,k) ! area in y required by actual client PE 636 DO i = coarse_bound_all(1,k), coarse_bound_all(2,k) ! area in x required by actual client PE 655 ! 656 !-- Loop over all client PEs 657 DO k = 1, size_of_array(2) 658 ! 659 !-- Area along y required by actual client PE 660 DO j = coarse_bound_all(3,k), coarse_bound_all(4,k) 661 ! 662 !-- Area along x required by actual client PE 663 DO i = coarse_bound_all(1,k), coarse_bound_all(2,k) 664 637 665 px = i / nrx 638 666 py = j / nry … … 642 670 643 671 ic = ic + 1 644 index_list(1,ic) = i - ( px * nrx ) + 1 + nbgp ! First index in Server Array 645 index_list(2,ic) = j - ( py * nry ) + 1 + nbgp ! Second index in Server Array 646 index_list(3,ic) = i - coarse_bound_all(1,k) + 1 ! x Index client coarse grid 647 index_list(4,ic) = j - coarse_bound_all(3,k) + 1 ! y Index client coarse grid 648 index_list(5,ic) = k - 1 ! PE Number client 649 index_list(6,ic) = server_pe ! PE Number server 672 ! 673 !-- First index in server array 674 index_list(1,ic) = i - ( px * nrx ) + 1 + nbgp 675 ! 676 !-- Second index in server array 677 index_list(2,ic) = j - ( py * nry ) + 1 + nbgp 678 ! 679 !-- x index of client coarse grid 680 index_list(3,ic) = i - coarse_bound_all(1,k) + 1 681 ! 682 !-- y index of client coarse grid 683 index_list(4,ic) = j - coarse_bound_all(3,k) + 1 684 ! 685 !-- PE number of client 686 index_list(5,ic) = k - 1 687 ! 688 !-- PE number of server 689 index_list(6,ic) = server_pe 690 650 691 ENDDO 651 692 ENDDO 652 693 ENDDO 694 ! 695 !-- TO_DO: Klaus: comment what is done here 653 696 CALL pmc_s_set_2d_index_list( client_id, index_list(:,1:ic) ) 697 654 698 ELSE 655 ALLOCATE( index_list(6,1) ) ! Dummy allocate 699 ! 700 !-- TO_DO: Klaus: comment why thie dummy allocation is required 701 ALLOCATE( index_list(6,1) ) 656 702 CALL pmc_s_set_2d_index_list( client_id, index_list ) 657 703 ENDIF … … 671 717 IMPLICIT NONE 672 718 673 CHARACTER(LEN=DA_Namelen) :: myname !: 674 675 INTEGER(iwp) :: i !: 676 INTEGER(iwp) :: ierr !: 677 INTEGER(iwp) :: icl !: 678 INTEGER(iwp) :: icr !: 679 INTEGER(iwp) :: j !: 680 INTEGER(iwp) :: jcn !: 681 INTEGER(iwp) :: jcs !: 719 CHARACTER(LEN=da_namelen) :: myname !: 720 721 INTEGER(iwp) :: i !: 722 INTEGER(iwp) :: ierr !: 723 INTEGER(iwp) :: icl !: 724 INTEGER(iwp) :: icr !: 725 INTEGER(iwp) :: j !: 726 INTEGER(iwp) :: jcn !: 727 INTEGER(iwp) :: jcs !: 728 682 729 INTEGER(iwp), DIMENSION(5) :: val !: 683 730 684 REAL(wp), DIMENSION(1) :: fval !: 685 REAL(wp) :: xcs !: 686 REAL(wp) :: xce !: 687 REAL(wp) :: ycs !: 688 REAL(wp) :: yce !: 731 REAL(wp) :: xcs !: 732 REAL(wp) :: xce !: 733 REAL(wp) :: ycs !: 734 REAL(wp) :: yce !: 735 736 REAL(wp), DIMENSION(1) :: fval !: 689 737 690 738 ! 691 739 !-- TO_DO: describe what is happening in this if-clause 692 740 !-- Root Model does not have Server and is not a client 693 741 IF ( .NOT. pmc_is_rootmodel() ) THEN 742 694 743 CALL pmc_clientinit 695 744 ! … … 725 774 726 775 IF ( myid == 0 ) THEN 776 727 777 CALL pmc_send_to_server( val, SIZE( val ), 0, 123, ierr ) 728 778 CALL pmc_send_to_server( fval, SIZE( fval ), 0, 124, ierr ) … … 732 782 ! 733 783 !-- Receive Coarse grid information. 784 !-- TO_DO: find shorter and more meaningful name for define_coarse_grid_real 734 785 CALL pmc_recv_from_server( define_coarse_grid_real, & 735 786 SIZE(define_coarse_grid_real), 0, 21, ierr ) … … 738 789 ! 739 790 !-- Receive also the dz-,zu- and zw-arrays here. 740 !-- TO_DO: what is the meaning of above comment + remove write statements 741 !-- and give this informations in header 742 WRITE(0,*) 'Coarse grid from Server ' 743 WRITE(0,*) 'startx_tot = ',define_coarse_grid_real(1) 744 WRITE(0,*) 'starty_tot = ',define_coarse_grid_real(2) 745 WRITE(0,*) 'endx_tot = ',define_coarse_grid_real(5) 746 WRITE(0,*) 'endy_tot = ',define_coarse_grid_real(6) 747 WRITE(0,*) 'dx = ',define_coarse_grid_real(3) 748 WRITE(0,*) 'dy = ',define_coarse_grid_real(4) 749 WRITE(0,*) 'dz = ',define_coarse_grid_real(7) 750 WRITE(0,*) 'nx_coarse = ',define_coarse_grid_int(1) 751 WRITE(0,*) 'ny_coarse = ',define_coarse_grid_int(2) 752 WRITE(0,*) 'nz_coarse = ',define_coarse_grid_int(3) 791 !-- TO_DO: what is the meaning of above comment 792 ! 793 !-- Debug-printouts - keep them 794 ! WRITE(0,*) 'Coarse grid from Server ' 795 ! WRITE(0,*) 'startx_tot = ',define_coarse_grid_real(1) 796 ! WRITE(0,*) 'starty_tot = ',define_coarse_grid_real(2) 797 ! WRITE(0,*) 'endx_tot = ',define_coarse_grid_real(5) 798 ! WRITE(0,*) 'endy_tot = ',define_coarse_grid_real(6) 799 ! WRITE(0,*) 'dx = ',define_coarse_grid_real(3) 800 ! WRITE(0,*) 'dy = ',define_coarse_grid_real(4) 801 ! WRITE(0,*) 'dz = ',define_coarse_grid_real(7) 802 ! WRITE(0,*) 'nx_coarse = ',define_coarse_grid_int(1) 803 ! WRITE(0,*) 'ny_coarse = ',define_coarse_grid_int(2) 804 ! WRITE(0,*) 'nz_coarse = ',define_coarse_grid_int(3) 753 805 ENDIF 754 806 … … 765 817 766 818 ! 767 !-- Get Server coordinates on coarse grid819 !-- Get server coordinates on coarse grid 768 820 ALLOCATE( cg%coord_x(-nbgp:cg%nx+nbgp) ) 769 821 ALLOCATE( cg%coord_y(-nbgp:cg%ny+nbgp) ) … … 776 828 ! 777 829 !-- Get coarse grid coordinates and vales of the z-direction from server 778 IF ( myid == 0) THEN 779 CALL pmc_recv_from_server( cg%coord_x, cg%nx + 1 + 2 * nbgp, 0, 24, & 780 ierr ) 781 CALL pmc_recv_from_server( cg%coord_y, cg%ny + 1 + 2 * nbgp, 0, 25, & 782 ierr ) 830 IF ( myid == 0) THEN 831 832 CALL pmc_recv_from_server( cg%coord_x, cg%nx+1+2*nbgp, 0, 24, ierr ) 833 CALL pmc_recv_from_server( cg%coord_y, cg%ny+1+2*nbgp, 0, 25, ierr ) 783 834 CALL pmc_recv_from_server( cg%dzu, cg%nz + 1, 0, 26, ierr ) 784 835 CALL pmc_recv_from_server( cg%dzw, cg%nz + 1, 0, 27, ierr ) 785 836 CALL pmc_recv_from_server( cg%zu, cg%nz + 2, 0, 28, ierr ) 786 837 CALL pmc_recv_from_server( cg%zw, cg%nz + 2, 0, 29, ierr ) 787 ENDIF 788 789 ! 790 !-- and broadcast this information 791 CALL MPI_BCAST( cg%coord_x, cg%nx + 1 + 2 * nbgp, MPI_REAL, 0, comm2d, & 792 ierr ) 793 CALL MPI_BCAST( cg%coord_y, cg%ny + 1 + 2 * nbgp, MPI_REAL, 0, comm2d, & 794 ierr ) 795 CALL MPI_BCAST( cg%dzu, cg%nz + 1, MPI_REAL, 0, comm2d, ierr ) 796 CALL MPI_BCAST( cg%dzw, cg%nz + 1, MPI_REAL, 0, comm2d, ierr ) 797 CALL MPI_BCAST( cg%zu, cg%nz + 2, MPI_REAL, 0, comm2d, ierr ) 798 CALL MPI_BCAST( cg%zw, cg%nz + 2, MPI_REAL, 0, comm2d, ierr ) 838 839 ENDIF 840 841 ! 842 !-- Broadcast this information 843 CALL MPI_BCAST( cg%coord_x, cg%nx+1+2*nbgp, MPI_REAL, 0, comm2d, ierr ) 844 CALL MPI_BCAST( cg%coord_y, cg%ny+1+2*nbgp, MPI_REAL, 0, comm2d, ierr ) 845 CALL MPI_BCAST( cg%dzu, cg%nz+1, MPI_REAL, 0, comm2d, ierr ) 846 CALL MPI_BCAST( cg%dzw, cg%nz+1, MPI_REAL, 0, comm2d, ierr ) 847 CALL MPI_BCAST( cg%zu, cg%nz+2, MPI_REAL, 0, comm2d, ierr ) 848 CALL MPI_BCAST( cg%zw, cg%nz+2, MPI_REAL, 0, comm2d, ierr ) 799 849 850 ! 851 !-- TO_DO: give comments what is happening here 800 852 CALL pmci_map_fine_to_coarse_grid 801 853 CALL pmc_c_get_2d_index_list 802 854 803 855 ! 804 !-- Include couple arrays into client content. 856 !-- Include couple arrays into client content 857 !-- TO_DO: Klaus: better explain the above comment (what is client content?) 805 858 CALL pmc_c_clear_next_array_list 806 859 DO WHILE ( pmc_c_getnextarray( myname ) ) … … 823 876 824 877 ! 825 !-- Two-way coupling 878 !-- Two-way coupling 879 !-- TO_DO: comment what is happening here 826 880 IF ( nesting_mode == 'two-way' ) THEN 827 881 CALL pmci_init_anterp_tophat … … 837 891 CONTAINS 838 892 839 840 893 SUBROUTINE pmci_map_fine_to_coarse_grid 841 894 842 843 844 INTEGER(iwp), DIMENSION(5,numprocs):: coarse_bound_all !:845 INTEGER(iwp), DIMENSION(2):: size_of_array !:895 IMPLICIT NONE 896 897 INTEGER(iwp), DIMENSION(5,numprocs) :: coarse_bound_all !: 898 INTEGER(iwp), DIMENSION(2) :: size_of_array !: 846 899 847 REAL(wp) :: coarse_dx!:848 REAL(wp) :: coarse_dy!:849 REAL(wp) :: loffset !:850 REAL(wp) :: noffset !:851 REAL(wp) :: roffset !: 852 REAL(wp) :: soffset !: 853 854 ! 855 !-- Determine indices of interpolation/anterpolation area in coarse grid. 856 coarse_dx = cg%dx857 coarse_dy = cg%dy 858 859 loffset = MOD( coord_x(nxl), coarse_dx ) ! If the fine- and coarse grid nodes do not match.860 xexl = coarse_dx + loffset861 nhll = CEILING( xexl / coarse_dx ) ! This is needed in the anterpolation phase.862 xcs = coord_x(nxl) - xexl863 DO i = 0, cg%nx864 IF ( cg%coord_x(i) > xcs ) THEN865 icl = MAX( -1, i-1 )866 EXIT867 ENDIF 868 ENDDO 869 870 roffset = MOD( coord_x(nxr + 1), coarse_dx ) ! If the fine- and coarse grid nodes do not match.871 xexr = coarse_dx + roffset 872 nhlr = CEILING( xexr / coarse_dx ) ! This is needed in the anterpolation phase. 873 xce = coord_x(nxr) + xexr874 DO i = cg%nx, 0 , -1875 IF ( cg%coord_x(i) < xce ) THEN876 icr = MIN( cg%nx + 1, i + 1 )877 EXIT878 ENDIF879 ENDDO880 881 soffset = MOD( coord_y(nys), coarse_dy ) ! If the fine- and coarse grid nodes do not match 882 yexs = coarse_dy + soffset 883 nhls = CEILING( yexs / coarse_dy ) ! This is needed in the anterpolation phase.884 ycs = coord_y(nys) - yexs885 DO j = 0, cg%ny 886 IF ( cg%coord_y(j) > ycs ) THEN 887 jcs = MAX( -nbgp, j - 1 )888 EXIT889 ENDIF890 ENDDO891 892 noffset = MOD( coord_y(nyn + 1), coarse_dy ) ! If the fine- and coarse grid nodes do not match893 yexn = coarse_dy + noffset894 nhln = CEILING( yexn / coarse_dy ) ! This is needed in the anterpolation phase.895 yce = coord_y(nyn) + yexn 896 DO j = cg%ny, 0, -1 897 IF ( cg%coord_y(j) < yce ) THEN898 jcn = MIN( cg%ny + nbgp, j + 1 )899 EXIT 900 ENDIF 901 ENDDO902 903 coarse_bound(1) = icl904 coarse_bound(2) = icr905 coarse_bound(3) = jcs906 coarse_bound(4) = jcn907 coarse_bound(5) = myid908 ! 909 !-- Note that MPI_Gather receives data from all processes in the rank order. 910 CALL MPI_GATHER( coarse_bound, 5, MPI_INTEGER, coarse_bound_all, 5, &911 MPI_INTEGER, 0, comm2d, ierr )912 913 IF ( myid == 0 ) THEN914 size_of_array(1) = SIZE( coarse_bound_all, 1 )915 size_of_array(2) = SIZE( coarse_bound_all, 2 ) 916 CALL pmc_send_to_server( size_of_array, 2, 0, 40, ierr ) 917 CALL pmc_send_to_server( coarse_bound_all, SIZE( coarse_bound_all ), 0, 41, ierr ) 918 ENDIF919 920 END SUBROUTINE pmci_map_fine_to_coarse_grid 921 922 923 924 SUBROUTINE pmci_init_interp_tril925 926 ! 927 !-- Precomputation of the interpolation coefficients and client-array indices 928 !-- to be used by the interpolation routines interp_tril_lr, interp_tril_ns and 929 !-- interp_tril_t. Constant dz is still assumed. 930 ! 931 !-- Antti Hellsten 3.3.2015. 932 ! 933 !-- Modified for variable dz, but not yet tested. 934 ! -- Antti Hellsten 27.3.2015.935 ! 936 IMPLICIT NONE 937 938 INTEGER(iwp) :: i !: 939 INTEGER(iwp) :: i1 !:940 INTEGER(iwp) :: j !: 941 INTEGER(iwp) :: j1!:942 INTEGER(iwp) :: k!:943 INTEGER(iwp) :: kc!:944 945 REAL(wp) :: coarse_dx!:946 REAL(wp) :: coarse_dy!:947 REAL(wp) :: coarse_dz !: 948 REAL(wp) :: xb !:949 REAL(wp) :: xcsu !:950 REAL(wp) :: xfso !:951 REAL(wp) :: xcso !:952 REAL(wp) :: xfsu !:953 REAL(wp) :: yb !:954 REAL(wp) :: ycso !:955 REAL(wp) :: ycsv !:956 REAL(wp) :: yfso !:957 REAL(wp) :: yfsv !:958 REAL(wp) :: zcso !:959 REAL(wp) :: zcsw !:960 REAL(wp) :: zfso !:961 REAL(wp) :: zfsw !:900 REAL(wp) :: loffset !: 901 REAL(wp) :: noffset !: 902 REAL(wp) :: roffset !: 903 REAL(wp) :: soffset !: 904 905 ! 906 !-- Determine indices of interpolation/anterpolation area in the coarse grid 907 !-- If the fine- and coarse grid nodes do not match. 908 loffset = MOD( coord_x(nxl), cg%dx ) 909 xexl = cg%dx + loffset 910 ! 911 !-- This is needed in the anterpolation phase 912 nhll = CEILING( xexl / cg%dx ) 913 xcs = coord_x(nxl) - xexl 914 DO i = 0, cg%nx 915 IF ( cg%coord_x(i) > xcs ) THEN 916 icl = MAX( -1, i-1 ) 917 EXIT 918 ENDIF 919 ENDDO 920 ! 921 !-- If the fine- and coarse grid nodes do not match 922 roffset = MOD( coord_x(nxr+1), cg%dx ) 923 xexr = cg%dx + roffset 924 ! 925 !-- This is needed in the anterpolation phase 926 nhlr = CEILING( xexr / cg%dx ) 927 xce = coord_x(nxr) + xexr 928 DO i = cg%nx, 0 , -1 929 IF ( cg%coord_x(i) < xce ) THEN 930 icr = MIN( cg%nx+1, i+1 ) 931 EXIT 932 ENDIF 933 ENDDO 934 ! 935 !-- If the fine- and coarse grid nodes do not match 936 soffset = MOD( coord_y(nys), cg%dy ) 937 yexs = cg%dy + soffset 938 ! 939 !-- This is needed in the anterpolation phase 940 nhls = CEILING( yexs / cg%dy ) 941 ycs = coord_y(nys) - yexs 942 DO j = 0, cg%ny 943 IF ( cg%coord_y(j) > ycs ) THEN 944 jcs = MAX( -nbgp, j-1 ) 945 EXIT 946 ENDIF 947 ENDDO 948 ! 949 !-- If the fine- and coarse grid nodes do not match 950 noffset = MOD( coord_y(nyn+1), cg%dy ) 951 yexn = cg%dy + noffset 952 ! 953 !-- This is needed in the anterpolation phase 954 nhln = CEILING( yexn / cg%dy ) 955 yce = coord_y(nyn) + yexn 956 DO j = cg%ny, 0, -1 957 IF ( cg%coord_y(j) < yce ) THEN 958 jcn = MIN( cg%ny + nbgp, j+1 ) 959 EXIT 960 ENDIF 961 ENDDO 962 963 coarse_bound(1) = icl 964 coarse_bound(2) = icr 965 coarse_bound(3) = jcs 966 coarse_bound(4) = jcn 967 coarse_bound(5) = myid 968 ! 969 !-- Note that MPI_Gather receives data from all processes in the rank order 970 !-- TO_DO: refer to the line where this fact becomes important 971 CALL MPI_GATHER( coarse_bound, 5, MPI_INTEGER, coarse_bound_all, 5, & 972 MPI_INTEGER, 0, comm2d, ierr ) 973 974 IF ( myid == 0 ) THEN 975 size_of_array(1) = SIZE( coarse_bound_all, 1 ) 976 size_of_array(2) = SIZE( coarse_bound_all, 2 ) 977 CALL pmc_send_to_server( size_of_array, 2, 0, 40, ierr ) 978 CALL pmc_send_to_server( coarse_bound_all, SIZE( coarse_bound_all ), & 979 0, 41, ierr ) 980 ENDIF 981 982 END SUBROUTINE pmci_map_fine_to_coarse_grid 983 984 985 986 SUBROUTINE pmci_init_interp_tril 987 ! 988 !-- Precomputation of the interpolation coefficients and client-array indices 989 !-- to be used by the interpolation routines interp_tril_lr, interp_tril_ns 990 !-- and interp_tril_t. Constant dz is still assumed. 991 992 IMPLICIT NONE 993 994 INTEGER(iwp) :: i !: 995 INTEGER(iwp) :: i1 !: 996 INTEGER(iwp) :: j !: 997 INTEGER(iwp) :: j1 !: 998 INTEGER(iwp) :: k !: 999 INTEGER(iwp) :: kc !: 1000 1001 REAL(wp) :: xb !: 1002 REAL(wp) :: xcsu !: 1003 REAL(wp) :: xfso !: 1004 REAL(wp) :: xcso !: 1005 REAL(wp) :: xfsu !: 1006 REAL(wp) :: yb !: 1007 REAL(wp) :: ycso !: 1008 REAL(wp) :: ycsv !: 1009 REAL(wp) :: yfso !: 1010 REAL(wp) :: yfsv !: 1011 REAL(wp) :: zcso !: 1012 REAL(wp) :: zcsw !: 1013 REAL(wp) :: zfso !: 1014 REAL(wp) :: zfsw !: 962 1015 963 1016 964 coarse_dx = cg%dx 965 coarse_dy = cg%dy 966 coarse_dz = cg%dz 967 xb = nxl * dx 968 yb = nys * dy 1017 xb = nxl * dx 1018 yb = nys * dy 969 1019 970 ALLOCATE( icu(nxlg:nxrg) )971 ALLOCATE( ico(nxlg:nxrg) )972 ALLOCATE( jcv(nysg:nyng) )973 ALLOCATE( jco(nysg:nyng) )974 ALLOCATE( kcw(nzb:nzt+1) )975 ALLOCATE( kco(nzb:nzt+1) )976 ALLOCATE( r1xu(nxlg:nxrg) )977 ALLOCATE( r2xu(nxlg:nxrg) )978 ALLOCATE( r1xo(nxlg:nxrg) )979 ALLOCATE( r2xo(nxlg:nxrg) )980 ALLOCATE( r1yv(nysg:nyng) )981 ALLOCATE( r2yv(nysg:nyng) )982 ALLOCATE( r1yo(nysg:nyng) )983 ALLOCATE( r2yo(nysg:nyng) )984 ALLOCATE( r1zw(nzb:nzt+1) )985 ALLOCATE( r2zw(nzb:nzt+1) )986 ALLOCATE( r1zo(nzb:nzt+1) )987 ALLOCATE( r2zo(nzb:nzt+1) )988 989 ! 990 !-- Note that the node coordinates xfs... and xcs... are relative991 !-- to the lower-left-bottom corner of the fc-array, not the actual992 !-- client domain corner.993 DO i = nxlg, nxrg994 xfsu = coord_x(i) - ( lower_left_coord_x + xb - xexl )995 xfso = coord_x(i) + 0.5_wp * dx - ( lower_left_coord_x + xb - xexl )996 icu(i) = icl + FLOOR( xfsu / coarse_dx )997 ico(i) = icl + FLOOR( ( xfso - 0.5_wp * coarse_dx ) / coarse_dx )998 xcsu = ( icu(i) - icl ) * coarse_dx999 xcso = ( ico(i) - icl ) * coarse_dx + 0.5_wp * coarse_dx1000 r2xu(i) = ( xfsu - xcsu ) / coarse_dx1001 r2xo(i) = ( xfso - xcso ) / coarse_dx1002 r1xu(i) = 1.0_wp - r2xu(i)1003 r1xo(i) = 1.0_wp - r2xo(i)1004 ENDDO1005 1006 DO j = nysg, nyng1007 yfsv = coord_y(j) - ( lower_left_coord_y + yb - yexs )1008 yfso = coord_y(j) + 0.5_wp * dy - ( lower_left_coord_y + yb - yexs )1009 jcv(j) = jcs + FLOOR( yfsv/coarse_dy )1010 jco(j) = jcs + FLOOR( ( yfso -0.5_wp * coarse_dy ) / coarse_dy )1011 ycsv = ( jcv(j) - jcs ) * coarse_dy1012 ycso = ( jco(j) - jcs ) * coarse_dy + 0.5_wp * coarse_dy1013 r2yv(j) = ( yfsv - ycsv ) / coarse_dy1014 r2yo(j) = ( yfso - ycso ) / coarse_dy1015 r1yv(j) = 1.0_wp - r2yv(j)1016 r1yo(j) = 1.0_wp - r2yo(j)1017 ENDDO1018 1019 DO k = nzb, nzt + 11020 zfsw = zw(k)1021 zfso = zu(k)1022 1023 kc = 01024 DO WHILE ( cg%zw(kc) <= zfsw )1025 kc = kc + 11026 ENDDO1027 kcw(k) = kc - 11020 ALLOCATE( icu(nxlg:nxrg) ) 1021 ALLOCATE( ico(nxlg:nxrg) ) 1022 ALLOCATE( jcv(nysg:nyng) ) 1023 ALLOCATE( jco(nysg:nyng) ) 1024 ALLOCATE( kcw(nzb:nzt+1) ) 1025 ALLOCATE( kco(nzb:nzt+1) ) 1026 ALLOCATE( r1xu(nxlg:nxrg) ) 1027 ALLOCATE( r2xu(nxlg:nxrg) ) 1028 ALLOCATE( r1xo(nxlg:nxrg) ) 1029 ALLOCATE( r2xo(nxlg:nxrg) ) 1030 ALLOCATE( r1yv(nysg:nyng) ) 1031 ALLOCATE( r2yv(nysg:nyng) ) 1032 ALLOCATE( r1yo(nysg:nyng) ) 1033 ALLOCATE( r2yo(nysg:nyng) ) 1034 ALLOCATE( r1zw(nzb:nzt+1) ) 1035 ALLOCATE( r2zw(nzb:nzt+1) ) 1036 ALLOCATE( r1zo(nzb:nzt+1) ) 1037 ALLOCATE( r2zo(nzb:nzt+1) ) 1038 1039 ! 1040 !-- Note that the node coordinates xfs... and xcs... are relative to the 1041 !-- lower-left-bottom corner of the fc-array, not the actual client domain 1042 !-- corner 1043 DO i = nxlg, nxrg 1044 xfsu = coord_x(i) - ( lower_left_coord_x + xb - xexl ) 1045 xfso = coord_x(i) + 0.5_wp * dx - ( lower_left_coord_x + xb - xexl ) 1046 icu(i) = icl + FLOOR( xfsu / cg%dx ) 1047 ico(i) = icl + FLOOR( ( xfso - 0.5_wp * cg%dx ) / cg%dx ) 1048 xcsu = ( icu(i) - icl ) * cg%dx 1049 xcso = ( ico(i) - icl ) * cg%dx + 0.5_wp * cg%dx 1050 r2xu(i) = ( xfsu - xcsu ) / cg%dx 1051 r2xo(i) = ( xfso - xcso ) / cg%dx 1052 r1xu(i) = 1.0_wp - r2xu(i) 1053 r1xo(i) = 1.0_wp - r2xo(i) 1054 ENDDO 1055 1056 DO j = nysg, nyng 1057 yfsv = coord_y(j) - ( lower_left_coord_y + yb - yexs ) 1058 yfso = coord_y(j) + 0.5_wp * dy - ( lower_left_coord_y + yb - yexs ) 1059 jcv(j) = jcs + FLOOR( yfsv / cg%dy ) 1060 jco(j) = jcs + FLOOR( ( yfso -0.5_wp * cg%dy ) / cg%dy ) 1061 ycsv = ( jcv(j) - jcs ) * cg%dy 1062 ycso = ( jco(j) - jcs ) * cg%dy + 0.5_wp * cg%dy 1063 r2yv(j) = ( yfsv - ycsv ) / cg%dy 1064 r2yo(j) = ( yfso - ycso ) / cg%dy 1065 r1yv(j) = 1.0_wp - r2yv(j) 1066 r1yo(j) = 1.0_wp - r2yo(j) 1067 ENDDO 1068 1069 DO k = nzb, nzt + 1 1070 zfsw = zw(k) 1071 zfso = zu(k) 1072 1073 kc = 0 1074 DO WHILE ( cg%zw(kc) <= zfsw ) 1075 kc = kc + 1 1076 ENDDO 1077 kcw(k) = kc - 1 1028 1078 1029 kc = 01030 DO WHILE ( cg%zu(kc) <= zfso )1031 kc = kc + 11032 ENDDO1033 kco(k) = kc - 11034 1035 zcsw = cg%zw(kcw(k))1036 zcso = cg%zu(kco(k))1037 r2zw(k) = ( zfsw - zcsw ) / cg%dzw(kcw(k) + 1)1038 r2zo(k) = ( zfso - zcso ) / cg%dzu(kco(k) + 1)1039 r1zw(k) = 1.0_wp - r2zw(k)1040 r1zo(k) = 1.0_wp - r2zo(k)1041 ENDDO1079 kc = 0 1080 DO WHILE ( cg%zu(kc) <= zfso ) 1081 kc = kc + 1 1082 ENDDO 1083 kco(k) = kc - 1 1084 1085 zcsw = cg%zw(kcw(k)) 1086 zcso = cg%zu(kco(k)) 1087 r2zw(k) = ( zfsw - zcsw ) / cg%dzw(kcw(k)+1) 1088 r2zo(k) = ( zfso - zcso ) / cg%dzu(kco(k)+1) 1089 r1zw(k) = 1.0_wp - r2zw(k) 1090 r1zo(k) = 1.0_wp - r2zo(k) 1091 ENDDO 1042 1092 1043 END SUBROUTINE pmci_init_interp_tril 1044 1045 1046 1047 SUBROUTINE pmci_init_loglaw_correction 1048 1049 ! 1050 !-- Precomputation of the index and log-ratio arrays for the log-law corrections 1051 !-- for near-wall nodes after the nest-BC interpolation. 1052 !-- These are used by the interpolation routines interp_tril_lr and interp_tril_ns. 1053 ! 1054 !-- Antti Hellsten 30.12.2015. 1055 ! 1056 IMPLICIT NONE 1057 INTEGER(iwp) :: direction !: Wall normal index: 1=k, 2=j, 3=i. 1058 INTEGER(iwp) :: i !: 1059 INTEGER(iwp) :: icorr !: 1060 INTEGER(iwp) :: inc !: Wall outward-normal index increment -1 or 1, for direction=1, inc=1 always. 1061 INTEGER(iwp) :: iw !: 1062 INTEGER(iwp) :: j !: 1063 INTEGER(iwp) :: jcorr !: 1064 INTEGER(iwp) :: jw !: 1065 INTEGER(iwp) :: k !: 1066 INTEGER(iwp) :: kb !: 1067 INTEGER(iwp) :: kcorr !: 1068 INTEGER(iwp) :: lc !: 1069 INTEGER(iwp) :: ni !: 1070 INTEGER(iwp) :: nj !: 1071 INTEGER(iwp) :: nk !: 1072 INTEGER(iwp) :: nzt_topo_max !: 1073 INTEGER(iwp) :: wall_index !: Index of the wall-node coordinate 1074 1075 REAL(wp), ALLOCATABLE, DIMENSION(:) :: lcr !: 1076 1077 ! 1078 !-- First determine the maximum k-index needed for the near-wall corrections. 1079 !-- This maximum is individual for each boundary to minimize the storage requirements 1080 !-- and to minimize the corresponding loop k-range in the interpolation routines. 1081 nzt_topo_nestbc_l = nzb 1082 IF ( nest_bound_l ) THEN 1083 DO i = nxl - 1, nxl 1084 DO j = nys, nyn 1085 nzt_topo_nestbc_l = MAX( nzt_topo_nestbc_l, nzb_u_inner(j,i), & 1086 nzb_v_inner(j,i), nzb_w_inner(j,i) ) 1087 ENDDO 1088 ENDDO 1089 nzt_topo_nestbc_l = nzt_topo_nestbc_l + 1 1090 ENDIF 1093 END SUBROUTINE pmci_init_interp_tril 1094 1095 1096 1097 SUBROUTINE pmci_init_loglaw_correction 1098 ! 1099 !-- Precomputation of the index and log-ratio arrays for the log-law 1100 !-- corrections for near-wall nodes after the nest-BC interpolation. 1101 !-- These are used by the interpolation routines interp_tril_lr and 1102 !-- interp_tril_ns. 1103 1104 IMPLICIT NONE 1105 1106 INTEGER(iwp) :: direction !: Wall normal index: 1=k, 2=j, 3=i. 1107 INTEGER(iwp) :: i !: 1108 INTEGER(iwp) :: icorr !: 1109 INTEGER(iwp) :: inc !: Wall outward-normal index increment -1 1110 !: or 1, for direction=1, inc=1 always 1111 INTEGER(iwp) :: iw !: 1112 INTEGER(iwp) :: j !: 1113 INTEGER(iwp) :: jcorr !: 1114 INTEGER(iwp) :: jw !: 1115 INTEGER(iwp) :: k !: 1116 INTEGER(iwp) :: kb !: 1117 INTEGER(iwp) :: kcorr !: 1118 INTEGER(iwp) :: lc !: 1119 INTEGER(iwp) :: ni !: 1120 INTEGER(iwp) :: nj !: 1121 INTEGER(iwp) :: nk !: 1122 INTEGER(iwp) :: nzt_topo_max !: 1123 INTEGER(iwp) :: wall_index !: Index of the wall-node coordinate 1124 1125 REAL(wp), ALLOCATABLE, DIMENSION(:) :: lcr !: 1126 1127 ! 1128 !-- First determine the maximum k-index needed for the near-wall corrections. 1129 !-- This maximum is individual for each boundary to minimize the storage 1130 !-- requirements and to minimize the corresponding loop k-range in the 1131 !-- interpolation routines. 1132 nzt_topo_nestbc_l = nzb 1133 IF ( nest_bound_l ) THEN 1134 DO i = nxl-1, nxl 1135 DO j = nys, nyn 1136 nzt_topo_nestbc_l = MAX( nzt_topo_nestbc_l, nzb_u_inner(j,i), & 1137 nzb_v_inner(j,i), nzb_w_inner(j,i) ) 1138 ENDDO 1139 ENDDO 1140 nzt_topo_nestbc_l = nzt_topo_nestbc_l + 1 1141 ENDIF 1091 1142 1092 nzt_topo_nestbc_r = nzb 1093 IF ( nest_bound_r ) THEN 1094 i = nxr + 1 1095 DO j = nys, nyn 1096 nzt_topo_nestbc_r = MAX( nzt_topo_nestbc_r, nzb_u_inner(j,i), & 1097 nzb_v_inner(j,i), nzb_w_inner(j,i) ) 1098 ENDDO 1099 nzt_topo_nestbc_r = nzt_topo_nestbc_r + 1 1100 ENDIF 1101 1102 nzt_topo_nestbc_s = nzb 1103 IF ( nest_bound_s ) THEN 1104 DO j = nys - 1, nys 1105 DO i = nxl, nxr 1106 nzt_topo_nestbc_s = MAX( nzt_topo_nestbc_s, nzb_u_inner(j,i), & 1107 nzb_v_inner(j,i), nzb_w_inner(j,i) ) 1108 ENDDO 1109 ENDDO 1110 nzt_topo_nestbc_s = nzt_topo_nestbc_s + 1 1111 ENDIF 1112 1113 nzt_topo_nestbc_n = nzb 1114 IF ( nest_bound_n ) THEN 1115 j = nyn + 1 1116 DO i = nxl, nxr 1117 nzt_topo_nestbc_n = MAX( nzt_topo_nestbc_n, nzb_u_inner(j,i), & 1118 nzb_v_inner(j,i), nzb_w_inner(j,i) ) 1119 ENDDO 1120 nzt_topo_nestbc_n = nzt_topo_nestbc_n + 1 1121 ENDIF 1122 1123 ! 1124 !-- Then determine the maximum number of near-wall nodes per wall point based 1125 !-- on the grid-spacing ratios. 1126 nzt_topo_max = MAX( nzt_topo_nestbc_l, nzt_topo_nestbc_r, & 1127 nzt_topo_nestbc_s, nzt_topo_nestbc_n ) 1128 1129 ! 1130 !-- Note that the outer division must be integer division. 1131 ni = CEILING( cg%dx / dx ) / 2 1132 nj = CEILING( cg%dy / dy ) / 2 1133 nk = 1 1134 DO k = 1, nzt_topo_max 1135 nk = MAX( nk, CEILING( cg%dzu(k) / dzu(k) ) ) 1136 ENDDO 1137 nk = nk / 2 ! Note that this must be integer division. 1138 ncorr = MAX( ni, nj, nk ) 1139 1140 ALLOCATE( lcr(0:ncorr - 1) ) 1141 lcr = 1.0_wp 1142 1143 ! 1144 !-- First horizontal walls 1145 !-- Left boundary 1146 IF ( nest_bound_l ) THEN 1147 ALLOCATE( logc_u_l(nzb:nzt_topo_nestbc_l, nys:nyn, 1:2) ) 1148 ALLOCATE( logc_v_l(nzb:nzt_topo_nestbc_l, nys:nyn, 1:2) ) 1149 ALLOCATE( logc_ratio_u_l(nzb:nzt_topo_nestbc_l, nys:nyn, 1:2, 0:ncorr-1) ) 1150 ALLOCATE( logc_ratio_v_l(nzb:nzt_topo_nestbc_l, nys:nyn, 1:2, 0:ncorr-1) ) 1151 logc_u_l = 0 1152 logc_v_l = 0 1153 logc_ratio_u_l = 1.0_wp 1154 logc_ratio_v_l = 1.0_wp 1155 direction = 1 1156 inc = 1 1157 1158 DO j = nys, nyn 1159 ! 1160 !-- Left boundary for u 1161 i = 0 1162 kb = nzb_u_inner(j,i) 1163 k = kb + 1 1164 wall_index = kb 1165 CALL pmci_define_loglaw_correction_parameters( lc, lcr, k, j, inc, & 1166 wall_index, z0(j,i), kb, direction, ncorr ) 1167 logc_u_l(k,j,1) = lc 1168 logc_ratio_u_l(k,j,1,0:ncorr-1) = lcr(0:ncorr-1) 1169 lcr(0:ncorr-1) = 1.0_wp 1170 1171 ! 1172 !-- Left boundary for v 1173 i = -1 1174 kb = nzb_v_inner(j,i) 1175 k = kb + 1 1176 wall_index = kb 1177 CALL pmci_define_loglaw_correction_parameters( lc, lcr, k, j, inc, & 1178 wall_index, z0(j,i), kb, direction, ncorr ) 1179 logc_v_l(k,j,1) = lc 1180 logc_ratio_v_l(k,j,1,0:ncorr-1) = lcr(0:ncorr-1) 1181 lcr(0:ncorr-1) = 1.0_wp 1182 1183 ENDDO 1184 ENDIF 1185 1186 ! 1187 !-- Right boundary 1188 IF ( nest_bound_r ) THEN 1189 ALLOCATE( logc_u_r(nzb:nzt_topo_nestbc_r,nys:nyn,1:2) ) 1190 ALLOCATE( logc_v_r(nzb:nzt_topo_nestbc_r,nys:nyn,1:2) ) 1191 ALLOCATE( logc_ratio_u_r(nzb:nzt_topo_nestbc_r,nys:nyn,1:2,0:ncorr-1) ) 1192 ALLOCATE( logc_ratio_v_r(nzb:nzt_topo_nestbc_r,nys:nyn,1:2,0:ncorr-1) ) 1193 logc_u_r = 0 1194 logc_v_r = 0 1195 logc_ratio_u_r = 1.0_wp 1196 logc_ratio_v_r = 1.0_wp 1197 direction = 1 1198 inc = 1 1199 DO j = nys, nyn 1200 ! 1201 !-- Right boundary for u. 1202 i = nxr + 1 1203 kb = nzb_u_inner(j,i) 1204 k = kb + 1 1205 wall_index = kb 1206 CALL pmci_define_loglaw_correction_parameters( lc, lcr, k, j, inc, & 1207 wall_index, z0(j,i), kb, direction, ncorr ) 1208 logc_u_r(k,j,1) = lc 1209 logc_ratio_u_r(k,j,1,0:ncorr-1) = lcr(0:ncorr-1) 1210 lcr(0:ncorr-1) = 1.0_wp 1211 1212 ! 1213 !-- Right boundary for v. 1214 i = nxr + 1 1215 kb = nzb_v_inner(j,i) 1216 k = kb + 1 1217 wall_index = kb 1218 CALL pmci_define_loglaw_correction_parameters( lc, lcr, k, j, inc, & 1219 wall_index, z0(j,i), kb, direction, ncorr ) 1220 logc_v_r(k,j,1) = lc 1221 logc_ratio_v_r(k,j,1,0:ncorr-1) = lcr(0:ncorr-1) 1222 lcr(0:ncorr-1) = 1.0_wp 1223 1224 ENDDO 1225 ENDIF 1226 1227 ! 1228 !-- South boundary 1229 IF ( nest_bound_s ) THEN 1230 ALLOCATE( logc_u_s(nzb:nzt_topo_nestbc_s,nxl:nxr,1:2) ) 1231 ALLOCATE( logc_v_s(nzb:nzt_topo_nestbc_s,nxl:nxr,1:2) ) 1232 ALLOCATE( logc_ratio_u_s(nzb:nzt_topo_nestbc_s,nxl:nxr,1:2,0:ncorr-1) ) 1233 ALLOCATE( logc_ratio_v_s(nzb:nzt_topo_nestbc_s,nxl:nxr,1:2,0:ncorr-1) ) 1234 logc_u_s = 0 1235 logc_v_s = 0 1236 logc_ratio_u_s = 1.0_wp 1237 logc_ratio_v_s = 1.0_wp 1238 direction = 1 1239 inc = 1 1240 DO i = nxl, nxr 1241 ! 1242 !-- South boundary for u. 1243 j = -1 1244 kb = nzb_u_inner(j,i) 1245 k = kb + 1 1246 wall_index = kb 1247 CALL pmci_define_loglaw_correction_parameters( lc, lcr, k, j, inc, & 1248 wall_index, z0(j,i), kb, direction, ncorr ) 1249 logc_u_s(k,i,1) = lc 1250 logc_ratio_u_s(k,i,1,0:ncorr-1) = lcr(0:ncorr-1) 1251 lcr(0:ncorr-1) = 1.0_wp 1252 1253 ! 1254 !-- South boundary for v 1255 j = 0 1256 kb = nzb_v_inner(j,i) 1257 k = kb + 1 1258 wall_index = kb 1259 CALL pmci_define_loglaw_correction_parameters( lc, lcr, k, j, inc, & 1260 wall_index, z0(j,i), kb, direction, ncorr ) 1261 logc_v_s(k,i,1) = lc 1262 logc_ratio_v_s(k,i,1,0:ncorr-1) = lcr(0:ncorr-1) 1263 lcr(0:ncorr-1) = 1.0_wp 1264 ENDDO 1265 ENDIF 1266 1267 ! 1268 !-- North boundary 1269 IF ( nest_bound_n ) THEN 1270 ALLOCATE( logc_u_n(nzb:nzt_topo_nestbc_n,nxl:nxr,1:2) ) 1271 ALLOCATE( logc_v_n(nzb:nzt_topo_nestbc_n,nxl:nxr,1:2) ) 1272 ALLOCATE( logc_ratio_u_n(nzb:nzt_topo_nestbc_n,nxl:nxr,1:2,0:ncorr-1) ) 1273 ALLOCATE( logc_ratio_v_n(nzb:nzt_topo_nestbc_n,nxl:nxr,1:2,0:ncorr-1) ) 1274 logc_u_n = 0 1275 logc_v_n = 0 1276 logc_ratio_u_n = 1.0_wp 1277 logc_ratio_v_n = 1.0_wp 1278 direction = 1 1279 inc = 1 1280 DO i = nxl, nxr 1281 ! 1282 !-- North boundary for u. 1283 j = nyn + 1 1284 kb = nzb_u_inner(j,i) 1285 k = kb + 1 1286 wall_index = kb 1287 CALL pmci_define_loglaw_correction_parameters( lc, lcr, k, j, inc, & 1288 wall_index, z0(j,i), kb, direction, ncorr ) 1289 logc_u_n(k,i,1) = lc 1290 logc_ratio_u_n(k,i,1,0:ncorr-1) = lcr(0:ncorr-1) 1291 lcr(0:ncorr-1) = 1.0_wp 1292 1293 ! 1294 !-- North boundary for v. 1295 j = nyn + 1 1296 kb = nzb_v_inner(j,i) 1297 k = kb + 1 1298 wall_index = kb 1299 CALL pmci_define_loglaw_correction_parameters( lc, lcr, k, j, inc, & 1300 wall_index, z0(j,i), kb, direction, ncorr ) 1301 logc_v_n(k,i,1) = lc 1302 logc_ratio_v_n(k,i,1,0:ncorr-1) = lcr(0:ncorr-1) 1303 lcr(0:ncorr-1) = 1.0_wp 1304 1305 ENDDO 1306 ENDIF 1143 nzt_topo_nestbc_r = nzb 1144 IF ( nest_bound_r ) THEN 1145 i = nxr + 1 1146 DO j = nys, nyn 1147 nzt_topo_nestbc_r = MAX( nzt_topo_nestbc_r, nzb_u_inner(j,i), & 1148 nzb_v_inner(j,i), nzb_w_inner(j,i) ) 1149 ENDDO 1150 nzt_topo_nestbc_r = nzt_topo_nestbc_r + 1 1151 ENDIF 1152 1153 nzt_topo_nestbc_s = nzb 1154 IF ( nest_bound_s ) THEN 1155 DO j = nys-1, nys 1156 DO i = nxl, nxr 1157 nzt_topo_nestbc_s = MAX( nzt_topo_nestbc_s, nzb_u_inner(j,i), & 1158 nzb_v_inner(j,i), nzb_w_inner(j,i) ) 1159 ENDDO 1160 ENDDO 1161 nzt_topo_nestbc_s = nzt_topo_nestbc_s + 1 1162 ENDIF 1163 1164 nzt_topo_nestbc_n = nzb 1165 IF ( nest_bound_n ) THEN 1166 j = nyn + 1 1167 DO i = nxl, nxr 1168 nzt_topo_nestbc_n = MAX( nzt_topo_nestbc_n, nzb_u_inner(j,i), & 1169 nzb_v_inner(j,i), nzb_w_inner(j,i) ) 1170 ENDDO 1171 nzt_topo_nestbc_n = nzt_topo_nestbc_n + 1 1172 ENDIF 1173 1174 ! 1175 !-- Then determine the maximum number of near-wall nodes per wall point based 1176 !-- on the grid-spacing ratios. 1177 nzt_topo_max = MAX( nzt_topo_nestbc_l, nzt_topo_nestbc_r, & 1178 nzt_topo_nestbc_s, nzt_topo_nestbc_n ) 1179 1180 ! 1181 !-- Note that the outer division must be integer division. 1182 ni = CEILING( cg%dx / dx ) / 2 1183 nj = CEILING( cg%dy / dy ) / 2 1184 nk = 1 1185 DO k = 1, nzt_topo_max 1186 nk = MAX( nk, CEILING( cg%dzu(k) / dzu(k) ) ) 1187 ENDDO 1188 nk = nk / 2 ! Note that this must be integer division. 1189 ncorr = MAX( ni, nj, nk ) 1190 1191 ALLOCATE( lcr(0:ncorr-1) ) 1192 lcr = 1.0_wp 1193 1194 ! 1195 !-- First horizontal walls 1196 !-- Left boundary 1197 IF ( nest_bound_l ) THEN 1198 ALLOCATE( logc_u_l(nzb:nzt_topo_nestbc_l,nys:nyn,1:2) ) 1199 ALLOCATE( logc_v_l(nzb:nzt_topo_nestbc_l,nys:nyn,1:2) ) 1200 ALLOCATE( logc_ratio_u_l(nzb:nzt_topo_nestbc_l,nys:nyn,1:2,0:ncorr-1) ) 1201 ALLOCATE( logc_ratio_v_l(nzb:nzt_topo_nestbc_l,nys:nyn,1:2,0:ncorr-1) ) 1202 logc_u_l = 0 1203 logc_v_l = 0 1204 logc_ratio_u_l = 1.0_wp 1205 logc_ratio_v_l = 1.0_wp 1206 direction = 1 1207 inc = 1 1208 1209 DO j = nys, nyn 1210 ! 1211 !-- Left boundary for u 1212 i = 0 1213 kb = nzb_u_inner(j,i) 1214 k = kb + 1 1215 wall_index = kb 1216 CALL pmci_define_loglaw_correction_parameters( lc, lcr, k, j, & 1217 inc, wall_index, z0(j,i), kb, direction, ncorr ) 1218 logc_u_l(k,j,1) = lc 1219 logc_ratio_u_l(k,j,1,0:ncorr-1) = lcr(0:ncorr-1) 1220 lcr(0:ncorr-1) = 1.0_wp 1221 ! 1222 !-- Left boundary for v 1223 i = -1 1224 kb = nzb_v_inner(j,i) 1225 k = kb + 1 1226 wall_index = kb 1227 CALL pmci_define_loglaw_correction_parameters( lc, lcr, k, j, & 1228 inc, wall_index, z0(j,i), kb, direction, ncorr ) 1229 logc_v_l(k,j,1) = lc 1230 logc_ratio_v_l(k,j,1,0:ncorr-1) = lcr(0:ncorr-1) 1231 lcr(0:ncorr-1) = 1.0_wp 1232 1233 ENDDO 1234 ENDIF 1235 1236 ! 1237 !-- Right boundary 1238 IF ( nest_bound_r ) THEN 1239 ALLOCATE( logc_u_r(nzb:nzt_topo_nestbc_r,nys:nyn,1:2) ) 1240 ALLOCATE( logc_v_r(nzb:nzt_topo_nestbc_r,nys:nyn,1:2) ) 1241 ALLOCATE( logc_ratio_u_r(nzb:nzt_topo_nestbc_r,nys:nyn,1:2,0:ncorr-1) ) 1242 ALLOCATE( logc_ratio_v_r(nzb:nzt_topo_nestbc_r,nys:nyn,1:2,0:ncorr-1) ) 1243 logc_u_r = 0 1244 logc_v_r = 0 1245 logc_ratio_u_r = 1.0_wp 1246 logc_ratio_v_r = 1.0_wp 1247 direction = 1 1248 inc = 1 1249 DO j = nys, nyn 1250 ! 1251 !-- Right boundary for u. 1252 i = nxr + 1 1253 kb = nzb_u_inner(j,i) 1254 k = kb + 1 1255 wall_index = kb 1256 CALL pmci_define_loglaw_correction_parameters( lc, lcr, k, j, & 1257 inc, wall_index, z0(j,i), kb, direction, ncorr ) 1258 logc_u_r(k,j,1) = lc 1259 logc_ratio_u_r(k,j,1,0:ncorr-1) = lcr(0:ncorr-1) 1260 lcr(0:ncorr-1) = 1.0_wp 1261 1262 ! 1263 !-- Right boundary for v. 1264 i = nxr + 1 1265 kb = nzb_v_inner(j,i) 1266 k = kb + 1 1267 wall_index = kb 1268 CALL pmci_define_loglaw_correction_parameters( lc, lcr, k, j, & 1269 inc, wall_index, z0(j,i), kb, direction, ncorr ) 1270 logc_v_r(k,j,1) = lc 1271 logc_ratio_v_r(k,j,1,0:ncorr-1) = lcr(0:ncorr-1) 1272 lcr(0:ncorr-1) = 1.0_wp 1273 1274 ENDDO 1275 ENDIF 1276 1277 ! 1278 !-- South boundary 1279 IF ( nest_bound_s ) THEN 1280 ALLOCATE( logc_u_s(nzb:nzt_topo_nestbc_s,nxl:nxr,1:2) ) 1281 ALLOCATE( logc_v_s(nzb:nzt_topo_nestbc_s,nxl:nxr,1:2) ) 1282 ALLOCATE( logc_ratio_u_s(nzb:nzt_topo_nestbc_s,nxl:nxr,1:2,0:ncorr-1) ) 1283 ALLOCATE( logc_ratio_v_s(nzb:nzt_topo_nestbc_s,nxl:nxr,1:2,0:ncorr-1) ) 1284 logc_u_s = 0 1285 logc_v_s = 0 1286 logc_ratio_u_s = 1.0_wp 1287 logc_ratio_v_s = 1.0_wp 1288 direction = 1 1289 inc = 1 1290 DO i = nxl, nxr 1291 ! 1292 !-- South boundary for u. 1293 j = -1 1294 kb = nzb_u_inner(j,i) 1295 k = kb + 1 1296 wall_index = kb 1297 CALL pmci_define_loglaw_correction_parameters( lc, lcr, k, j, & 1298 inc, wall_index, z0(j,i), kb, direction, ncorr ) 1299 logc_u_s(k,i,1) = lc 1300 logc_ratio_u_s(k,i,1,0:ncorr-1) = lcr(0:ncorr-1) 1301 lcr(0:ncorr-1) = 1.0_wp 1302 1303 ! 1304 !-- South boundary for v 1305 j = 0 1306 kb = nzb_v_inner(j,i) 1307 k = kb + 1 1308 wall_index = kb 1309 CALL pmci_define_loglaw_correction_parameters( lc, lcr, k, j, & 1310 inc, wall_index, z0(j,i), kb, direction, ncorr ) 1311 logc_v_s(k,i,1) = lc 1312 logc_ratio_v_s(k,i,1,0:ncorr-1) = lcr(0:ncorr-1) 1313 lcr(0:ncorr-1) = 1.0_wp 1314 ENDDO 1315 ENDIF 1316 1317 ! 1318 !-- North boundary 1319 IF ( nest_bound_n ) THEN 1320 ALLOCATE( logc_u_n(nzb:nzt_topo_nestbc_n,nxl:nxr,1:2) ) 1321 ALLOCATE( logc_v_n(nzb:nzt_topo_nestbc_n,nxl:nxr,1:2) ) 1322 ALLOCATE( logc_ratio_u_n(nzb:nzt_topo_nestbc_n,nxl:nxr,1:2,0:ncorr-1) ) 1323 ALLOCATE( logc_ratio_v_n(nzb:nzt_topo_nestbc_n,nxl:nxr,1:2,0:ncorr-1) ) 1324 logc_u_n = 0 1325 logc_v_n = 0 1326 logc_ratio_u_n = 1.0_wp 1327 logc_ratio_v_n = 1.0_wp 1328 direction = 1 1329 inc = 1 1330 DO i = nxl, nxr 1331 ! 1332 !-- North boundary for u. 1333 j = nyn + 1 1334 kb = nzb_u_inner(j,i) 1335 k = kb + 1 1336 wall_index = kb 1337 CALL pmci_define_loglaw_correction_parameters( lc, lcr, k, j, & 1338 inc, wall_index, z0(j,i), kb, direction, ncorr ) 1339 logc_u_n(k,i,1) = lc 1340 logc_ratio_u_n(k,i,1,0:ncorr-1) = lcr(0:ncorr-1) 1341 lcr(0:ncorr-1) = 1.0_wp 1342 1343 ! 1344 !-- North boundary for v. 1345 j = nyn + 1 1346 kb = nzb_v_inner(j,i) 1347 k = kb + 1 1348 wall_index = kb 1349 CALL pmci_define_loglaw_correction_parameters( lc, lcr, k, j, & 1350 inc, wall_index, z0(j,i), kb, direction, ncorr ) 1351 logc_v_n(k,i,1) = lc 1352 logc_ratio_v_n(k,i,1,0:ncorr-1) = lcr(0:ncorr-1) 1353 lcr(0:ncorr-1) = 1.0_wp 1354 1355 ENDDO 1356 ENDIF 1307 1357 1308 1358 ! 1309 !-- Then vertical walls and corners if necessary.1310 IF ( topography /= 'flat' ) THEN1311 kb = 0 ! kb is not used when direction > 11359 !-- Then vertical walls and corners if necessary. 1360 IF ( topography /= 'flat' ) THEN 1361 kb = 0 ! kb is not used when direction > 1 1312 1362 ! 1313 !-- Left boundary 1314 IF ( nest_bound_l ) THEN 1315 ALLOCATE( logc_w_l(nzb:nzt_topo_nestbc_l,nys:nyn,1:2) ) 1316 ALLOCATE( logc_ratio_w_l(nzb:nzt_topo_nestbc_l,nys:nyn,1:2,0:ncorr-1) ) 1317 logc_w_l = 0 1318 logc_ratio_w_l = 1.0_wp 1319 direction = 2 1320 DO j = nys, nyn 1321 DO k = nzb, nzt_topo_nestbc_l 1322 1323 ! 1324 !-- Wall for u on the south side, but not on the north side. 1325 i = 0 1326 IF ( ( nzb_u_outer(j,i) > nzb_u_outer(j+1,i) ) .AND. & 1327 ( nzb_u_outer(j,i) == nzb_u_outer(j-1,i) ) ) THEN 1328 inc = 1 1329 wall_index = j 1330 CALL pmci_define_loglaw_correction_parameters( lc, lcr, k, j, inc, wall_index, z0(j,i), kb, direction, ncorr ) 1331 1332 ! 1333 !-- The direction of the wall-normal index is stored as the sign of the logc-element. 1334 logc_u_l(k,j,2) = inc * lc 1335 logc_ratio_u_l(k,j,2,0:ncorr-1) = lcr(0:ncorr-1) 1336 lcr(0:ncorr-1) = 1.0_wp 1337 ENDIF 1338 1339 ! 1340 !-- Wall for u on the north side, but not on the south side. 1341 i = 0 1363 !-- Left boundary 1364 IF ( nest_bound_l ) THEN 1365 ALLOCATE( logc_w_l(nzb:nzt_topo_nestbc_l,nys:nyn,1:2) ) 1366 ALLOCATE( logc_ratio_w_l(nzb:nzt_topo_nestbc_l,nys:nyn,1:2, & 1367 0:ncorr-1) ) 1368 logc_w_l = 0 1369 logc_ratio_w_l = 1.0_wp 1370 direction = 2 1371 DO j = nys, nyn 1372 DO k = nzb, nzt_topo_nestbc_l 1373 1374 ! 1375 !-- Wall for u on the south side, but not on the north side 1376 i = 0 1377 IF ( ( nzb_u_outer(j,i) > nzb_u_outer(j+1,i) ) .AND. & 1378 ( nzb_u_outer(j,i) == nzb_u_outer(j-1,i) ) ) & 1379 THEN 1380 inc = 1 1381 wall_index = j 1382 CALL pmci_define_loglaw_correction_parameters( lc, lcr, & 1383 k, j, inc, wall_index, z0(j,i), kb, direction, ncorr ) 1384 ! 1385 !-- The direction of the wall-normal index is stored as the 1386 !-- sign of the logc-element. 1387 logc_u_l(k,j,2) = inc * lc 1388 logc_ratio_u_l(k,j,2,0:ncorr-1) = lcr(0:ncorr-1) 1389 lcr(0:ncorr-1) = 1.0_wp 1390 ENDIF 1391 ! 1392 !-- Wall for u on the north side, but not on the south side 1393 i = 0 1394 !-- TO_DO: routine must be indentet by 1 space from here on, 1395 !-- and long lines must be wrapped 1342 1396 IF ( ( nzb_u_outer(j,i) > nzb_u_outer(j-1,i) ) .AND. & 1343 1397 ( nzb_u_outer(j,i) == nzb_u_outer(j+1,i) ) ) THEN … … 1619 1673 1620 1674 1621 SUBROUTINE pmci_define_loglaw_correction_parameters( lc, lcr, k, ij, inc, wall_index, z0_l, kb, direction, ncorr ) 1622 IMPLICIT NONE 1623 INTEGER(iwp), INTENT(IN) :: direction !: 1624 INTEGER(iwp), INTENT(IN) :: ij !: 1625 INTEGER(iwp), INTENT(IN) :: inc !: 1626 INTEGER(iwp), INTENT(IN) :: k !: 1627 INTEGER(iwp), INTENT(IN) :: kb !: 1628 INTEGER(iwp), INTENT(OUT) :: lc !: 1629 INTEGER(iwp), INTENT(IN) :: ncorr !: 1630 INTEGER(iwp), INTENT(IN) :: wall_index !: 1631 REAL(wp), DIMENSION(0:ncorr-1), INTENT(OUT) :: lcr !: 1632 REAL(wp), INTENT(IN) :: z0_l !: 1675 SUBROUTINE pmci_define_loglaw_correction_parameters( lc, lcr, k, ij, inc, & 1676 wall_index, z0_l, kb, direction, ncorr ) 1677 1678 IMPLICIT NONE 1679 1680 INTEGER(iwp), INTENT(IN) :: direction !: 1681 INTEGER(iwp), INTENT(IN) :: ij !: 1682 INTEGER(iwp), INTENT(IN) :: inc !: 1683 INTEGER(iwp), INTENT(IN) :: k !: 1684 INTEGER(iwp), INTENT(IN) :: kb !: 1685 INTEGER(iwp), INTENT(OUT) :: lc !: 1686 INTEGER(iwp), INTENT(IN) :: ncorr !: 1687 INTEGER(iwp), INTENT(IN) :: wall_index !: 1688 1689 INTEGER(iwp) :: alcorr !: 1690 INTEGER(iwp) :: corr_index !: 1691 INTEGER(iwp) :: lcorr !: 1692 1693 LOGICAL :: more !: 1694 1695 REAL(wp), DIMENSION(0:ncorr-1), INTENT(OUT) :: lcr !: 1696 REAL(wp), INTENT(IN) :: z0_l !: 1633 1697 1634 INTEGER(iwp) :: alcorr !: 1635 INTEGER(iwp) :: corr_index !: 1636 INTEGER(iwp) :: lcorr !: 1637 REAL(wp) :: logvelc1 !: 1638 LOGICAL :: more !: 1698 REAL(wp) :: logvelc1 !: 1639 1699 1640 1700 1641 SELECT CASE ( direction )1642 1643 CASE (1) ! k1644 more = .TRUE.1645 lcorr = 01646 DO WHILE ( more .AND. lcorr <= ncorr -1 )1647 corr_index = k + lcorr1648 IF ( lcorr == 0 ) THEN1649 CALL pmci_find_logc_pivot_k( lc, logvelc1, z0_l, kb )1650 ENDIF1701 SELECT CASE ( direction ) 1702 1703 CASE (1) ! k 1704 more = .TRUE. 1705 lcorr = 0 1706 DO WHILE ( more .AND. lcorr <= ncorr-1 ) 1707 corr_index = k + lcorr 1708 IF ( lcorr == 0 ) THEN 1709 CALL pmci_find_logc_pivot_k( lc, logvelc1, z0_l, kb ) 1710 ENDIF 1651 1711 1652 IF ( corr_index < lc ) THEN 1653 lcr(lcorr) = LOG( ( zu(k) - zw(kb) ) / z0_l ) / logvelc1 1654 more = .TRUE. 1655 ELSE 1656 lcr(lcorr) = 1.0 1657 more = .FALSE. 1658 ENDIF 1659 lcorr = lcorr + 1 1660 ENDDO 1661 1662 CASE (2) ! j 1663 more = .TRUE. 1664 lcorr = 0 1665 alcorr = 0 1666 DO WHILE ( more .AND. alcorr <= ncorr - 1 ) 1667 corr_index = ij + lcorr ! In this case (direction = 2) ij is j 1668 IF ( lcorr == 0 ) THEN 1669 CALL pmci_find_logc_pivot_j( lc, logvelc1, ij, wall_index, z0_l, inc ) 1670 ENDIF 1671 1672 ! 1673 !-- The role of inc here is to make the comparison operation "<" valid in both directions. 1674 IF ( inc * corr_index < inc * lc ) THEN 1675 lcr(alcorr) = LOG( ABS( coord_y(corr_index) + 0.5_wp * dy & 1676 - coord_y(wall_index) ) / z0_l ) / logvelc1 1677 more = .TRUE. 1678 ELSE 1679 lcr(alcorr) = 1.0_wp 1680 more = .FALSE. 1681 ENDIF 1682 lcorr = lcorr + inc 1683 alcorr = ABS( lcorr ) 1684 ENDDO 1685 1686 CASE (3) ! i 1687 more = .TRUE. 1688 lcorr = 0 1689 alcorr = 0 1690 DO WHILE ( more .AND. alcorr <= ncorr - 1 ) 1691 corr_index = ij + lcorr ! In this case (direction = 3) ij is i 1692 IF ( lcorr == 0 ) THEN 1693 CALL pmci_find_logc_pivot_i( lc, logvelc1, ij, wall_index, z0_l, inc ) 1694 ENDIF 1695 1696 ! 1697 !-- The role of inc here is to make the comparison operation "<" valid in both directions. 1698 IF ( inc * corr_index < inc * lc ) THEN 1699 lcr(alcorr) = LOG( ABS( coord_x(corr_index) + 0.5_wp * dx & 1700 - coord_x(wall_index) ) / z0_l ) / logvelc1 1701 more = .TRUE. 1702 ELSE 1703 lcr(alcorr) = 1.0_wp 1704 more = .FALSE. 1705 ENDIF 1706 lcorr = lcorr + inc 1707 alcorr = ABS( lcorr ) 1708 ENDDO 1709 1710 END SELECT 1711 1712 END SUBROUTINE pmci_define_loglaw_correction_parameters 1713 1714 1715 1716 SUBROUTINE pmci_find_logc_pivot_k( lc, logzc1, z0_l, kb ) 1717 1718 ! 1719 !-- Finds the pivot node and te log-law factor for near-wall nodes for 1720 !-- which the wall-parallel velocity components will be log-law corrected 1721 !-- after interpolation. This subroutine is only for horizontal walls. 1722 ! 1723 !-- Antti Hellsten 30.12.2015 1724 IMPLICIT NONE 1725 REAL(wp),INTENT(OUT) :: logzc1 !: 1726 REAL(wp), INTENT(IN) :: z0_l !: 1727 INTEGER(iwp), INTENT(IN) :: kb !: 1728 INTEGER(iwp), INTENT(OUT) :: lc !: 1712 IF ( corr_index < lc ) THEN 1713 lcr(lcorr) = LOG( ( zu(k) - zw(kb) ) / z0_l ) / logvelc1 1714 more = .TRUE. 1715 ELSE 1716 lcr(lcorr) = 1.0 1717 more = .FALSE. 1718 ENDIF 1719 lcorr = lcorr + 1 1720 ENDDO 1721 1722 CASE (2) ! j 1723 more = .TRUE. 1724 lcorr = 0 1725 alcorr = 0 1726 DO WHILE ( more .AND. alcorr <= ncorr-1 ) 1727 corr_index = ij + lcorr ! In this case (direction = 2) ij is j 1728 IF ( lcorr == 0 ) THEN 1729 CALL pmci_find_logc_pivot_j( lc, logvelc1, ij, wall_index, & 1730 z0_l, inc ) 1731 ENDIF 1732 1733 ! 1734 !-- The role of inc here is to make the comparison operation "<" 1735 !-- valid in both directions 1736 IF ( inc * corr_index < inc * lc ) THEN 1737 lcr(alcorr) = LOG( ABS( coord_y(corr_index) + 0.5_wp * dy & 1738 - coord_y(wall_index) ) / z0_l ) & 1739 / logvelc1 1740 more = .TRUE. 1741 ELSE 1742 lcr(alcorr) = 1.0_wp 1743 more = .FALSE. 1744 ENDIF 1745 lcorr = lcorr + inc 1746 alcorr = ABS( lcorr ) 1747 ENDDO 1748 1749 CASE (3) ! i 1750 more = .TRUE. 1751 lcorr = 0 1752 alcorr = 0 1753 DO WHILE ( more .AND. alcorr <= ncorr-1 ) 1754 corr_index = ij + lcorr ! In this case (direction = 3) ij is i 1755 IF ( lcorr == 0 ) THEN 1756 CALL pmci_find_logc_pivot_i( lc, logvelc1, ij, wall_index, & 1757 z0_l, inc ) 1758 ENDIF 1759 ! 1760 !-- The role of inc here is to make the comparison operation "<" 1761 !-- valid in both directions 1762 IF ( inc * corr_index < inc * lc ) THEN 1763 lcr(alcorr) = LOG( ABS( coord_x(corr_index) + 0.5_wp * dx & 1764 - coord_x(wall_index) ) / z0_l ) & 1765 / logvelc1 1766 more = .TRUE. 1767 ELSE 1768 lcr(alcorr) = 1.0_wp 1769 more = .FALSE. 1770 ENDIF 1771 lcorr = lcorr + inc 1772 alcorr = ABS( lcorr ) 1773 ENDDO 1774 1775 END SELECT 1776 1777 END SUBROUTINE pmci_define_loglaw_correction_parameters 1778 1779 1780 1781 SUBROUTINE pmci_find_logc_pivot_k( lc, logzc1, z0_l, kb ) 1782 ! 1783 !-- Finds the pivot node and te log-law factor for near-wall nodes for 1784 !-- which the wall-parallel velocity components will be log-law corrected 1785 !-- after interpolation. This subroutine is only for horizontal walls. 1786 1787 IMPLICIT NONE 1788 1789 INTEGER(iwp), INTENT(IN) :: kb !: 1790 INTEGER(iwp), INTENT(OUT) :: lc !: 1791 1792 INTEGER(iwp) :: kbc !: 1793 INTEGER(iwp) :: k1 !: 1794 1795 REAL(wp),INTENT(OUT) :: logzc1 !: 1796 REAL(wp), INTENT(IN) :: z0_l !: 1797 1798 REAL(wp) :: zuc1 !: 1799 1800 1801 kbc = nzb + 1 1802 ! 1803 !-- kbc is the first coarse-grid point above the surface 1804 DO WHILE ( cg%zu(kbc) < zu(kb) ) 1805 kbc = kbc + 1 1806 ENDDO 1807 zuc1 = cg%zu(kbc) 1808 k1 = kb + 1 1809 DO WHILE ( zu(k1) < zuc1 ) ! Important: must be <, not <= 1810 k1 = k1 + 1 1811 ENDDO 1812 logzc1 = LOG( (zu(k1) - zw(kb) ) / z0_l ) 1813 lc = k1 1814 1815 END SUBROUTINE pmci_find_logc_pivot_k 1816 1817 1818 1819 SUBROUTINE pmci_find_logc_pivot_j( lc, logyc1, j, jw, z0_l, inc ) 1820 ! 1821 !-- Finds the pivot node and te log-law factor for near-wall nodes for 1822 !-- which the wall-parallel velocity components will be log-law corrected 1823 !-- after interpolation. This subroutine is only for vertical walls on 1824 !-- south/north sides of the node. 1825 1826 IMPLICIT NONE 1827 1828 INTEGER(iwp), INTENT(IN) :: inc !: increment must be 1 or -1. 1829 INTEGER(iwp), INTENT(IN) :: j !: 1830 INTEGER(iwp), INTENT(IN) :: jw !: 1831 INTEGER(iwp), INTENT(OUT) :: lc !: 1832 1833 INTEGER(iwp) :: j1 !: 1834 1835 REAL(wp), INTENT(IN) :: z0_l !: 1836 1837 REAL(wp) :: logyc1 !: 1838 REAL(wp) :: yc1 !: 1839 1840 ! 1841 !-- yc1 is the y-coordinate of the first coarse-grid u- and w-nodes out from 1842 !-- the wall 1843 yc1 = coord_y(jw) + 0.5_wp * inc * cg%dy 1844 ! 1845 !-- j1 is the first fine-grid index further away from the wall than yc1 1846 j1 = j 1847 ! 1848 !-- Important: must be <, not <= 1849 DO WHILE ( inc * ( coord_y(j1) + 0.5_wp * dy ) < inc * yc1 ) 1850 j1 = j1 + inc 1851 ENDDO 1852 1853 logyc1 = LOG( ABS( coord_y(j1) + 0.5_wp * dy - coord_y(jw) ) / z0_l ) 1854 lc = j1 1855 1856 END SUBROUTINE pmci_find_logc_pivot_j 1857 1858 1859 1860 SUBROUTINE pmci_find_logc_pivot_i( lc, logxc1, i, iw, z0_l, inc ) 1861 ! 1862 !-- Finds the pivot node and the log-law factor for near-wall nodes for 1863 !-- which the wall-parallel velocity components will be log-law corrected 1864 !-- after interpolation. This subroutine is only for vertical walls on 1865 !-- south/north sides of the node. 1866 1867 IMPLICIT NONE 1868 1869 INTEGER(iwp), INTENT(IN) :: i !: 1870 INTEGER(iwp), INTENT(IN) :: inc !: increment must be 1 or -1. 1871 INTEGER(iwp), INTENT(IN) :: iw !: 1872 INTEGER(iwp), INTENT(OUT) :: lc !: 1873 1874 INTEGER(iwp) :: i1 !: 1875 1876 REAL(wp), INTENT(IN) :: z0_l !: 1877 1878 REAL(wp) :: logxc1 !: 1879 REAL(wp) :: xc1 !: 1880 1881 ! 1882 !-- xc1 is the x-coordinate of the first coarse-grid v- and w-nodes out from 1883 !-- the wall 1884 xc1 = coord_x(iw) + 0.5_wp * inc * cg%dx 1885 ! 1886 !-- i1 is the first fine-grid index futher away from the wall than xc1. 1887 i1 = i 1888 ! 1889 !-- Important: must be <, not <= 1890 DO WHILE ( inc * ( coord_x(i1) + 0.5_wp *dx ) < inc * xc1 ) 1891 i1 = i1 + inc 1892 ENDDO 1729 1893 1730 REAL(wp) :: zuc1 !: 1731 INTEGER(iwp) :: kbc !: 1732 INTEGER(iwp) :: k1 !: 1733 1734 1735 kbc = nzb + 1 1736 DO WHILE ( cg%zu(kbc) < zu(kb) ) ! kbc is the first coarse-grid point above the surface. 1737 kbc = kbc + 1 1738 ENDDO 1739 zuc1 = cg%zu(kbc) 1740 k1 = kb + 1 1741 DO WHILE ( zu(k1) < zuc1 ) ! Important: must be <, not <= 1742 k1 = k1 + 1 1743 ENDDO 1744 logzc1 = LOG( (zu(k1) - zw(kb) ) / z0_l ) 1745 lc = k1 1746 1747 END SUBROUTINE pmci_find_logc_pivot_k 1748 1749 1750 1751 SUBROUTINE pmci_find_logc_pivot_j( lc, logyc1, j, jw, z0_l, inc ) 1752 1753 ! 1754 !-- Finds the pivot node and te log-law factor for near-wall nodes for 1755 !-- which the wall-parallel velocity components will be log-law corrected 1756 !-- after interpolation. This subroutine is only for vertical walls on 1757 !-- south/north sides of the node. 1758 ! 1759 !-- Antti Hellsten 5.1.2016 1760 IMPLICIT NONE 1761 REAL(wp), INTENT(IN) :: z0_l !: 1762 INTEGER(iwp), INTENT(IN) :: inc !: increment must be 1 or -1. 1763 INTEGER(iwp), INTENT(IN) :: j !: 1764 INTEGER(iwp), INTENT(IN) :: jw !: 1765 INTEGER(iwp), INTENT(OUT) :: lc !: 1766 1767 REAL(wp) :: logyc1 !: 1768 REAL(wp) :: yc1 !: 1769 INTEGER(iwp) :: j1 !: 1770 1771 ! 1772 !-- yc1 is the y-coordinate of the first coarse-grid u- and w-nodes out from the wall. 1773 yc1 = coord_y(jw) + 0.5_wp * inc * cg%dy 1774 1775 ! 1776 !-- j1 is the first fine-grid index futher away from the wall than yc1. 1777 j1 = j 1778 DO WHILE ( inc * ( coord_y(j1) + 0.5_wp * dy ) < inc * yc1 ) ! Important: must be <, not <= 1779 j1 = j1 + inc 1780 ENDDO 1781 1782 logyc1 = LOG( ABS( coord_y(j1) + 0.5_wp * dy - coord_y(jw) ) / z0_l ) 1783 lc = j1 1784 1785 END SUBROUTINE pmci_find_logc_pivot_j 1786 1787 1788 1789 SUBROUTINE pmci_find_logc_pivot_i( lc, logxc1, i, iw, z0_l, inc ) 1790 1791 ! 1792 !-- Finds the pivot node and the log-law factor for near-wall nodes for 1793 !-- which the wall-parallel velocity components will be log-law corrected 1794 !-- after interpolation. This subroutine is only for vertical walls on 1795 !-- south/north sides of the node. 1796 ! 1797 !-- Antti Hellsten 8.1.2016 1798 1799 IMPLICIT NONE 1800 REAL(wp), INTENT(IN) :: z0_l !: 1801 INTEGER(iwp), INTENT(IN) :: i !: 1802 INTEGER(iwp), INTENT(IN) :: inc !: increment must be 1 or -1. 1803 INTEGER(iwp), INTENT(IN) :: iw !: 1804 INTEGER(iwp), INTENT(OUT) :: lc !: 1805 1806 REAL(wp) :: logxc1 !: 1807 REAL(wp) :: xc1 !: 1808 INTEGER(iwp) :: i1 !: 1809 1810 ! 1811 !-- xc1 is the x-coordinate of the first coarse-grid v- and w-nodes out from the wall. 1812 xc1 = coord_x(iw) + 0.5_wp *inc * cg%dx 1813 1814 ! 1815 !-- i1 is the first fine-grid index futher away from the wall than xc1. 1816 i1 = i 1817 DO WHILE ( inc * ( coord_x(i1) + 0.5_wp *dx ) < inc *xc1 ) ! Important: must be <, not <= 1818 i1 = i1 + inc 1819 ENDDO 1820 1821 logxc1 = LOG( ABS( coord_x(i1) + 0.5_wp*dx - coord_x(iw) ) / z0_l ) 1822 lc = i1 1823 1824 END SUBROUTINE pmci_find_logc_pivot_i 1825 1826 1827 1894 logxc1 = LOG( ABS( coord_x(i1) + 0.5_wp*dx - coord_x(iw) ) / z0_l ) 1895 lc = i1 1896 1897 END SUBROUTINE pmci_find_logc_pivot_i 1898 1899 1900 !-- TO_DO: indentation and wrap long lines from here on to the end of the file 1828 1901 SUBROUTINE pmci_init_anterp_tophat 1829 1902 ! … … 1831 1904 !-- corresponding coarse-grid array index and the 1832 1905 !-- Under-relaxation coefficients to be used by anterp_tophat. 1833 ! 1834 !-- Antti Hellsten 9.10.2015. 1906 1835 1907 IMPLICIT NONE 1836 INTEGER(iwp) :: i !: 1908 1909 INTEGER(iwp) :: i !: Fine-grid index 1910 INTEGER(iwp) :: ii !: Coarse-grid index 1837 1911 INTEGER(iwp) :: istart !: 1838 INTEGER(iwp) :: j !: 1912 INTEGER(iwp) :: j !: Fine-grid index 1913 INTEGER(iwp) :: jj !: Coarse-grid index 1839 1914 INTEGER(iwp) :: jstart !: 1840 INTEGER(iwp) :: k !: 1915 INTEGER(iwp) :: k !: Fine-grid index 1916 INTEGER(iwp) :: kk !: Coarse-grid index 1841 1917 INTEGER(iwp) :: kstart !: 1842 INTEGER(iwp) :: l !:1843 INTEGER(iwp) :: m !:1844 INTEGER(iwp) :: n !:1845 1918 REAL(wp) :: xi !: 1846 1919 REAL(wp) :: eta !: 1847 1920 REAL(wp) :: zeta !: 1848 1921 1922 1849 1923 ! 1850 1924 !-- Default values: … … 1866 1940 1867 1941 ! 1868 !-- First determine kc eu and kcew that are the coarse-grid upper bounds for index k.1869 n= 01870 DO WHILE ( cg%zu( n) < zu(nzt) )1871 n = n+ 11942 !-- First determine kctu and kctw that are the coarse-grid upper bounds for index k. 1943 kk = 0 1944 DO WHILE ( cg%zu(kk) < zu(nzt) ) 1945 kk = kk + 1 1872 1946 ENDDO 1873 kc eu = n- 11874 1875 n= 01876 DO WHILE ( cg%zw( n) < zw(nzt-1) )1877 n = n+ 11947 kctu = kk - 1 1948 1949 kk = 0 1950 DO WHILE ( cg%zw(kk) < zw(nzt-1) ) 1951 kk = kk + 1 1878 1952 ENDDO 1879 kc ew = n- 11953 kctw = kk - 1 1880 1954 1881 1955 ALLOCATE( iflu(icl:icr) ) … … 1887 1961 ALLOCATE( jfuv(jcs:jcn) ) 1888 1962 ALLOCATE( jfuo(jcs:jcn) ) 1889 ALLOCATE( kflw(0:kc ew) )1890 ALLOCATE( kflo(0:kc eu) )1891 ALLOCATE( kfuw(0:kc ew) )1892 ALLOCATE( kfuo(0:kc eu) )1963 ALLOCATE( kflw(0:kctw) ) 1964 ALLOCATE( kflo(0:kctu) ) 1965 ALLOCATE( kfuw(0:kctw) ) 1966 ALLOCATE( kfuo(0:kctu) ) 1893 1967 1894 1968 ! 1895 1969 !-- i-indices of u for each l-index value. 1896 1970 istart = nxlg 1897 DO l= icl, icr1971 DO ii = icl, icr 1898 1972 i = istart 1899 DO WHILE ( ( coord_x(i) < cg%coord_x( l) - 0.5_wp * cg%dx ) .AND. ( i < nxrg ) )1973 DO WHILE ( ( coord_x(i) < cg%coord_x(ii) - 0.5_wp * cg%dx ) .AND. ( i < nxrg ) ) 1900 1974 i = i + 1 1901 1975 ENDDO 1902 iflu( l) = MIN( MAX( i, nxlg ), nxrg )1903 DO WHILE ( ( coord_x(i) < cg%coord_x( l) + 0.5_wp * cg%dx ) .AND. ( i < nxrg ) )1976 iflu(ii) = MIN( MAX( i, nxlg ), nxrg ) 1977 DO WHILE ( ( coord_x(i) < cg%coord_x(ii) + 0.5_wp * cg%dx ) .AND. ( i < nxrg ) ) 1904 1978 i = i + 1 1905 1979 ENDDO 1906 ifuu( l) = MIN( MAX( i, nxlg ), nxrg )1907 istart = iflu( l)1980 ifuu(ii) = MIN( MAX( i, nxlg ), nxrg ) 1981 istart = iflu(ii) 1908 1982 ENDDO 1909 1983 … … 1911 1985 !-- i-indices of others for each l-index value. 1912 1986 istart = nxlg 1913 DO l= icl, icr1987 DO ii = icl, icr 1914 1988 i = istart 1915 DO WHILE ( ( coord_x(i) + 0.5_wp * dx < cg%coord_x( l) ) .AND. ( i < nxrg ) )1989 DO WHILE ( ( coord_x(i) + 0.5_wp * dx < cg%coord_x(ii) ) .AND. ( i < nxrg ) ) 1916 1990 i = i + 1 1917 1991 ENDDO 1918 iflo( l) = MIN( MAX( i, nxlg ), nxrg )1919 DO WHILE ( ( coord_x(i) + 0.5_wp * dx < cg%coord_x( l) + cg%dx ) .AND. ( i < nxrg ) )1992 iflo(ii) = MIN( MAX( i, nxlg ), nxrg ) 1993 DO WHILE ( ( coord_x(i) + 0.5_wp * dx < cg%coord_x(ii) + cg%dx ) .AND. ( i < nxrg ) ) 1920 1994 i = i + 1 1921 1995 ENDDO 1922 ifuo( l) = MIN(MAX(i,nxlg),nxrg)1923 istart = iflo( l)1996 ifuo(ii) = MIN(MAX(i,nxlg),nxrg) 1997 istart = iflo(ii) 1924 1998 ENDDO 1925 1999 … … 1927 2001 !-- j-indices of v for each m-index value. 1928 2002 jstart = nysg 1929 DO m= jcs, jcn2003 DO jj = jcs, jcn 1930 2004 j = jstart 1931 DO WHILE ( ( coord_y(j) < cg%coord_y( m) - 0.5_wp * cg%dy ) .AND. ( j < nyng ) )2005 DO WHILE ( ( coord_y(j) < cg%coord_y(jj) - 0.5_wp * cg%dy ) .AND. ( j < nyng ) ) 1932 2006 j = j + 1 1933 2007 ENDDO 1934 jflv( m) = MIN( MAX( j, nysg ), nyng )1935 DO WHILE ( ( coord_y(j) < cg%coord_y( m) + 0.5_wp * cg%dy ) .AND. ( j < nyng ) )2008 jflv(jj) = MIN( MAX( j, nysg ), nyng ) 2009 DO WHILE ( ( coord_y(j) < cg%coord_y(jj) + 0.5_wp * cg%dy ) .AND. ( j < nyng ) ) 1936 2010 j = j + 1 1937 2011 ENDDO 1938 jfuv( m) = MIN( MAX( j, nysg ), nyng )1939 jstart = jflv( m)2012 jfuv(jj) = MIN( MAX( j, nysg ), nyng ) 2013 jstart = jflv(jj) 1940 2014 ENDDO 1941 2015 … … 1943 2017 !-- j-indices of others for each m-index value. 1944 2018 jstart = nysg 1945 DO m= jcs, jcn2019 DO jj = jcs, jcn 1946 2020 j = jstart 1947 DO WHILE ( ( coord_y(j) + 0.5_wp * dy < cg%coord_y( m) ) .AND. ( j < nyng ) )2021 DO WHILE ( ( coord_y(j) + 0.5_wp * dy < cg%coord_y(jj) ) .AND. ( j < nyng ) ) 1948 2022 j = j + 1 1949 2023 ENDDO 1950 jflo( m) = MIN( MAX( j, nysg ), nyng )1951 DO WHILE ( ( coord_y(j) + 0.5_wp * dy < cg%coord_y( m) + cg%dy ) .AND. ( j < nyng ) )2024 jflo(jj) = MIN( MAX( j, nysg ), nyng ) 2025 DO WHILE ( ( coord_y(j) + 0.5_wp * dy < cg%coord_y(jj) + cg%dy ) .AND. ( j < nyng ) ) 1952 2026 j = j + 1 1953 2027 ENDDO 1954 jfuo( m) = MIN( MAX( j, nysg ), nyng )1955 jstart = jflv( m)2028 jfuo(jj) = MIN( MAX( j, nysg ), nyng ) 2029 jstart = jflv(jj) 1956 2030 ENDDO 1957 2031 … … 1961 2035 kflw(0) = 0 1962 2036 kfuw(0) = 0 1963 DO n = 1, kcew2037 DO kk = 1, kctw 1964 2038 k = kstart 1965 DO WHILE ( ( zw(k) < cg%zw( n) - 0.5_wp * cg%dzw(n) ) .AND. ( k < nzt ) )2039 DO WHILE ( ( zw(k) < cg%zw(kk) - 0.5_wp * cg%dzw(kk) ) .AND. ( k < nzt ) ) 1966 2040 k = k + 1 1967 2041 ENDDO 1968 kflw( n) = MIN( MAX( k, 1 ), nzt + 1 )1969 DO WHILE ( ( zw(k) < cg%zw( n) + 0.5_wp * cg%dzw(n+1) ) .AND. ( k < nzt ) )2042 kflw(kk) = MIN( MAX( k, 1 ), nzt + 1 ) 2043 DO WHILE ( ( zw(k) < cg%zw(kk) + 0.5_wp * cg%dzw(kk+1) ) .AND. ( k < nzt ) ) 1970 2044 k = k + 1 1971 2045 ENDDO 1972 kfuw( n) = MIN( MAX( k, 1 ), nzt + 1 )1973 kstart = kflw( n)2046 kfuw(kk) = MIN( MAX( k, 1 ), nzt + 1 ) 2047 kstart = kflw(kk) 1974 2048 ENDDO 1975 2049 … … 1979 2053 kflo(0) = 0 1980 2054 kfuo(0) = 0 1981 DO n = 1, kceu2055 DO kk = 1, kctu 1982 2056 k = kstart 1983 DO WHILE ( ( zu(k) < cg%zu( n) - 0.5_wp * cg%dzu(n) ) .AND. ( k < nzt ) )2057 DO WHILE ( ( zu(k) < cg%zu(kk) - 0.5_wp * cg%dzu(kk) ) .AND. ( k < nzt ) ) 1984 2058 k = k + 1 1985 2059 ENDDO 1986 kflo( n) = MIN( MAX( k, 1 ), nzt + 1 )1987 DO WHILE ( ( zu(k) < cg%zu( n) + 0.5_wp * cg%dzu(n+1) ) .AND. ( k < nzt ) )2060 kflo(kk) = MIN( MAX( k, 1 ), nzt + 1 ) 2061 DO WHILE ( ( zu(k) < cg%zu(kk) + 0.5_wp * cg%dzu(kk+1) ) .AND. ( k < nzt ) ) 1988 2062 k = k + 1 1989 2063 ENDDO 1990 kfuo( n) = MIN( MAX( k-1, 1 ), nzt + 1 )1991 kstart = kflo( n)2064 kfuo(kk) = MIN( MAX( k-1, 1 ), nzt + 1 ) 2065 kstart = kflo(kk) 1992 2066 ENDDO 1993 2067 … … 1995 2069 !-- Spatial under-relaxation coefficients 1996 2070 ALLOCATE( frax(icl:icr) ) 1997 DO l = icl, icr 2071 2072 DO ii = icl, icr 1998 2073 IF ( nest_bound_l ) THEN 1999 xi = ( ( cg%coord_x(l) - lower_left_coord_x) / anterp_relax_length_l )**42074 xi = ( MAX( 0.0_wp, ( cg%coord_x(ii) - lower_left_coord_x ) ) / anterp_relax_length_l )**4 2000 2075 ELSEIF ( nest_bound_r ) THEN 2001 xi = ( ( lower_left_coord_x + ( nx + 1 ) * dx - cg%coord_x(l) ) / anterp_relax_length_r )**42076 xi = ( MAX( 0.0_wp, ( lower_left_coord_x + ( nx + 1 ) * dx - cg%coord_x(ii) ) ) / anterp_relax_length_r )**4 2002 2077 ELSE 2003 2078 xi = 999999.9_wp 2004 2079 ENDIF 2005 frax( l) = xi / ( 1.0_wp + xi )2080 frax(ii) = xi / ( 1.0_wp + xi ) 2006 2081 ENDDO 2007 2082 2008 2083 ALLOCATE( fray(jcs:jcn) ) 2009 DO m = jcs, jcn 2084 2085 DO jj = jcs, jcn 2010 2086 IF ( nest_bound_s ) THEN 2011 eta = ( ( cg%coord_y(m) - lower_left_coord_y) / anterp_relax_length_s )**42087 eta = ( MAX( 0.0_wp, ( cg%coord_y(jj) - lower_left_coord_y ) ) / anterp_relax_length_s )**4 2012 2088 ELSEIF ( nest_bound_n ) THEN 2013 eta = ( (lower_left_coord_y + ( ny + 1 ) * dy - cg%coord_y(m)) / anterp_relax_length_n )**42089 eta = ( MAX( 0.0_wp, ( lower_left_coord_y + ( ny + 1 ) * dy - cg%coord_y(jj)) ) / anterp_relax_length_n )**4 2014 2090 ELSE 2015 2091 eta = 999999.9_wp 2016 2092 ENDIF 2017 fray( m) = eta / ( 1.0_wp + eta )2093 fray(jj) = eta / ( 1.0_wp + eta ) 2018 2094 ENDDO 2019 2095 2020 ALLOCATE( fraz(0:kc eu) )2021 DO n = 0, kceu2022 zeta = ( ( zu(nzt) - cg%zu( n) ) / anterp_relax_length_t )**42023 fraz( n) = zeta / ( 1.0_wp + zeta )2096 ALLOCATE( fraz(0:kctu) ) 2097 DO kk = 0, kctu 2098 zeta = ( ( zu(nzt) - cg%zu(kk) ) / anterp_relax_length_t )**4 2099 fraz(kk) = zeta / ( 1.0_wp + zeta ) 2024 2100 ENDDO 2025 2101 … … 2173 2249 2174 2250 2251 2252 SUBROUTINE pmci_set_array_pointer( name, client_id, nz_cl ) 2253 2254 IMPLICIT NONE 2255 2256 INTEGER, INTENT(IN) :: client_id !: 2257 INTEGER, INTENT(IN) :: nz_cl !: 2258 CHARACTER(LEN=*), INTENT(IN) :: name !: 2259 2260 #if defined( __parallel ) 2261 INTEGER(iwp) :: ierr !: 2262 INTEGER(iwp) :: istat !: 2263 2264 REAL(wp), POINTER, DIMENSION(:,:) :: p_2d !: 2265 REAL(wp), POINTER, DIMENSION(:,:) :: p_2d_sec !: 2266 REAL(wp), POINTER, DIMENSION(:,:,:) :: p_3d !: 2267 REAL(wp), POINTER, DIMENSION(:,:,:) :: p_3d_sec !: 2268 2269 2270 NULLIFY( p_3d ) 2271 NULLIFY( p_2d ) 2272 2273 ! 2274 !-- List of array names, which can be coupled. 2275 !-- In case of 3D please change also the second array for the pointer version 2276 IF ( TRIM(name) == "u" ) p_3d => u 2277 IF ( TRIM(name) == "v" ) p_3d => v 2278 IF ( TRIM(name) == "w" ) p_3d => w 2279 IF ( TRIM(name) == "e" ) p_3d => e 2280 IF ( TRIM(name) == "pt" ) p_3d => pt 2281 IF ( TRIM(name) == "q" ) p_3d => q 2282 ! 2283 !-- Next line is just an example for a 2D array (not active for coupling!) 2284 !-- Please note, that z0 has to be declared as TARGET array in modules.f90 2285 ! IF ( TRIM(name) == "z0" ) p_2d => z0 2286 2287 #if defined( __nopointer ) 2288 IF ( ASSOCIATED( p_3d ) ) THEN 2289 CALL pmc_s_set_dataarray( client_id, p_3d, nz_cl, nz ) 2290 ELSEIF ( ASSOCIATED( p_2d ) ) THEN 2291 CALL pmc_s_set_dataarray( client_id, p_2d ) 2292 ELSE 2293 ! 2294 !-- Give only one message for the root domain 2295 IF ( myid == 0 .AND. cpl_id == 1 ) THEN 2296 2297 message_string = 'pointer for array "' // TRIM( name ) // & 2298 '" can''t be associated' 2299 CALL message( 'pmci_set_array_pointer', 'PA0117', 3, 2, 0, 6, 0 ) 2300 ELSE 2301 ! 2302 !-- Avoid others to continue 2303 CALL MPI_BARRIER( comm2d, ierr ) 2304 ENDIF 2305 ENDIF 2306 #else 2307 IF ( TRIM(name) == "u" ) p_3d_sec => u_2 2308 IF ( TRIM(name) == "v" ) p_3d_sec => v_2 2309 IF ( TRIM(name) == "w" ) p_3d_sec => w_2 2310 IF ( TRIM(name) == "e" ) p_3d_sec => e_2 2311 IF ( TRIM(name) == "pt" ) p_3d_sec => pt_2 2312 IF ( TRIM(name) == "q" ) p_3d_sec => q_2 2313 2314 IF ( ASSOCIATED( p_3d ) ) THEN 2315 CALL pmc_s_set_dataarray( client_id, p_3d, nz_cl, nz, & 2316 array_2 = p_3d_sec ) 2317 ELSEIF ( ASSOCIATED( p_2d ) ) THEN 2318 CALL pmc_s_set_dataarray( client_id, p_2d ) 2319 ELSE 2320 ! 2321 !-- Give only one message for the root domain 2322 IF ( myid == 0 .AND. cpl_id == 1 ) THEN 2323 2324 message_string = 'pointer for array "' // TRIM( name ) // & 2325 '" can''t be associated' 2326 CALL message( 'pmci_set_array_pointer', 'PA0117', 3, 2, 0, 6, 0 ) 2327 ELSE 2328 ! 2329 !-- Avoid others to continue 2330 CALL MPI_BARRIER( comm2d, ierr ) 2331 ENDIF 2332 2333 ENDIF 2334 #endif 2335 2336 #endif 2337 END SUBROUTINE pmci_set_array_pointer 2338 2339 2340 2341 SUBROUTINE pmci_create_client_arrays( name, is, ie, js, je, nzc ) 2342 2343 IMPLICIT NONE 2344 2345 CHARACTER(LEN=*), INTENT(IN) :: name !: 2346 2347 INTEGER(iwp), INTENT(IN) :: ie !: 2348 INTEGER(iwp), INTENT(IN) :: is !: 2349 INTEGER(iwp), INTENT(IN) :: je !: 2350 INTEGER(iwp), INTENT(IN) :: js !: 2351 INTEGER(iwp), INTENT(IN) :: nzc !: Note that nzc is cg%nz 2352 2353 #if defined( __parallel ) 2354 INTEGER(iwp) :: ierr !: 2355 INTEGER(iwp) :: istat !: 2356 2357 REAL(wp), POINTER,DIMENSION(:,:) :: p_2d !: 2358 REAL(wp), POINTER,DIMENSION(:,:,:) :: p_3d !: 2359 2360 2361 NULLIFY( p_3d ) 2362 NULLIFY( p_2d ) 2363 2364 ! 2365 !-- List of array names, which can be coupled. 2366 !-- TO_DO: Antti: what is the meaning of the next line? 2367 !-- AH: Note that the k-range of the *c arrays is changed from 1:nz to 0:nz+1. 2368 IF ( TRIM( name ) == "u" ) THEN 2369 IF ( .NOT. ALLOCATED( uc ) ) ALLOCATE( uc(0:nzc+1, js:je, is:ie) ) 2370 p_3d => uc 2371 ELSEIF ( TRIM( name ) == "v" ) THEN 2372 IF ( .NOT. ALLOCATED( vc ) ) ALLOCATE( vc(0:nzc+1, js:je, is:ie) ) 2373 p_3d => vc 2374 ELSEIF ( TRIM( name ) == "w" ) THEN 2375 IF ( .NOT. ALLOCATED( wc ) ) ALLOCATE( wc(0:nzc+1, js:je, is:ie) ) 2376 p_3d => wc 2377 ELSEIF ( TRIM( name ) == "e" ) THEN 2378 IF ( .NOT. ALLOCATED( ec ) ) ALLOCATE( ec(0:nzc+1, js:je, is:ie) ) 2379 p_3d => ec 2380 ELSEIF ( TRIM( name ) == "pt") THEN 2381 IF ( .NOT. ALLOCATED( ptc ) ) ALLOCATE( ptc(0:nzc+1, js:je, is:ie) ) 2382 p_3d => ptc 2383 ELSEIF ( TRIM( name ) == "q") THEN 2384 IF ( .NOT. ALLOCATED( qc ) ) ALLOCATE( qc(0:nzc+1, js:je, is:ie) ) 2385 p_3d => qc 2386 !ELSEIF (trim(name) == "z0") then 2387 !IF (.not.allocated(z0c)) allocate(z0c(js:je, is:ie)) 2388 !p_2d => z0c 2389 ENDIF 2390 2391 IF ( ASSOCIATED( p_3d ) ) THEN 2392 CALL pmc_c_set_dataarray( p_3d ) 2393 ELSEIF ( ASSOCIATED( p_2d ) ) THEN 2394 CALL pmc_c_set_dataarray( p_2d ) 2395 ELSE 2396 ! 2397 !-- Give only one message for the first client domain 2398 IF ( myid == 0 .AND. cpl_id == 2 ) THEN 2399 2400 message_string = 'pointer for array "' // TRIM( name ) // & 2401 '" can''t be associated' 2402 CALL message( 'pmci_create_client_arrays', 'PA0170', 3, 2, 0, 6, 0 ) 2403 ELSE 2404 ! 2405 !-- Avoid others to continue 2406 CALL MPI_BARRIER( comm2d, ierr ) 2407 ENDIF 2408 ENDIF 2409 2410 #endif 2411 END SUBROUTINE pmci_create_client_arrays 2412 2413 2414 2415 SUBROUTINE pmci_server_initialize 2416 !-- TO_DO: add general explanations about what this subroutine does 2417 #if defined( __parallel ) 2418 IMPLICIT NONE 2419 2420 INTEGER(iwp) :: client_id !: 2421 INTEGER(iwp) :: m !: 2422 2423 REAL(wp) :: waittime !: 2424 2425 2426 DO m = 1, SIZE( pmc_server_for_client ) - 1 2427 client_id = pmc_server_for_client(m) 2428 CALL pmc_s_fillbuffer( client_id, waittime=waittime ) 2429 ENDDO 2430 2431 #endif 2432 END SUBROUTINE pmci_server_initialize 2433 2434 2435 2436 SUBROUTINE pmci_client_initialize 2437 !-- TO_DO: add general explanations about what this subroutine does 2438 #if defined( __parallel ) 2439 IMPLICIT NONE 2440 2441 INTEGER(iwp) :: i !: 2442 INTEGER(iwp) :: icl !: 2443 INTEGER(iwp) :: icr !: 2444 INTEGER(iwp) :: j !: 2445 INTEGER(iwp) :: jcn !: 2446 INTEGER(iwp) :: jcs !: 2447 2448 REAL(wp) :: waittime !: 2449 2450 ! 2451 !-- Root id is never a client 2452 IF ( cpl_id > 1 ) THEN 2453 2454 ! 2455 !-- Client domain boundaries in the server index space 2456 icl = coarse_bound(1) 2457 icr = coarse_bound(2) 2458 jcs = coarse_bound(3) 2459 jcn = coarse_bound(4) 2460 2461 ! 2462 !-- Get data from server 2463 CALL pmc_c_getbuffer( waittime = waittime ) 2464 2465 ! 2466 !-- The interpolation. 2467 CALL pmci_interp_tril_all ( u, uc, icu, jco, kco, r1xu, r2xu, r1yo, & 2468 r2yo, r1zo, r2zo, nzb_u_inner, 'u' ) 2469 CALL pmci_interp_tril_all ( v, vc, ico, jcv, kco, r1xo, r2xo, r1yv, & 2470 r2yv, r1zo, r2zo, nzb_v_inner, 'v' ) 2471 CALL pmci_interp_tril_all ( w, wc, ico, jco, kcw, r1xo, r2xo, r1yo, & 2472 r2yo, r1zw, r2zw, nzb_w_inner, 'w' ) 2473 CALL pmci_interp_tril_all ( e, ec, ico, jco, kco, r1xo, r2xo, r1yo, & 2474 r2yo, r1zo, r2zo, nzb_s_inner, 'e' ) 2475 CALL pmci_interp_tril_all ( pt, ptc, ico, jco, kco, r1xo, r2xo, r1yo, & 2476 r2yo, r1zo, r2zo, nzb_s_inner, 's' ) 2477 IF ( humidity .OR. passive_scalar ) THEN 2478 CALL pmci_interp_tril_all ( q, qc, ico, jco, kco, r1xo, r2xo, r1yo, & 2479 r2yo, r1zo, r2zo, nzb_s_inner, 's' ) 2480 ENDIF 2481 2482 IF ( topography /= 'flat' ) THEN 2483 ! 2484 !-- Inside buildings set velocities and TKE back to zero. 2485 !-- Other scalars (pt, q, s, km, kh, p, sa, ...) are ignored at present, 2486 !-- maybe revise later. 2487 DO i = nxlg, nxrg 2488 DO j = nysg, nyng 2489 u(nzb:nzb_u_inner(j,i),j,i) = 0.0_wp 2490 v(nzb:nzb_v_inner(j,i),j,i) = 0.0_wp 2491 w(nzb:nzb_w_inner(j,i),j,i) = 0.0_wp 2492 e(nzb:nzb_s_inner(j,i),j,i) = 0.0_wp 2493 u_p(nzb:nzb_u_inner(j,i),j,i) = 0.0_wp 2494 v_p(nzb:nzb_v_inner(j,i),j,i) = 0.0_wp 2495 w_p(nzb:nzb_w_inner(j,i),j,i) = 0.0_wp 2496 e_p(nzb:nzb_s_inner(j,i),j,i) = 0.0_wp 2497 ENDDO 2498 ENDDO 2499 ENDIF 2500 ENDIF 2501 2502 2503 CONTAINS 2504 2505 2506 SUBROUTINE pmci_interp_tril_all( f, fc, ic, jc, kc, r1x, r2x, r1y, r2y, & 2507 r1z, r2z, kb, var ) 2508 ! 2509 !-- Interpolation of the internal values for the client-domain initialization 2510 !-- This subroutine is based on trilinear interpolation. 2511 !-- Coding based on interp_tril_lr/sn/t 2512 IMPLICIT NONE 2513 2514 CHARACTER(LEN=1), INTENT(IN) :: var !: 2515 2516 INTEGER(iwp), DIMENSION(nxlg:nxrg), INTENT(IN) :: ic !: 2517 INTEGER(iwp), DIMENSION(nysg:nyng), INTENT(IN) :: jc !: 2518 INTEGER(iwp), DIMENSION(nzb:nzt+1), INTENT(IN) :: kc !: 2519 INTEGER(iwp), DIMENSION(nysg:nyng,nxlg:nxrg), INTENT(IN) :: kb !: 2520 2521 INTEGER(iwp) :: i !: 2522 INTEGER(iwp) :: ib !: 2523 INTEGER(iwp) :: ie !: 2524 INTEGER(iwp) :: j !: 2525 INTEGER(iwp) :: jb !: 2526 INTEGER(iwp) :: je !: 2527 INTEGER(iwp) :: k !: 2528 INTEGER(iwp) :: k1 !: 2529 INTEGER(iwp) :: kbc !: 2530 INTEGER(iwp) :: l !: 2531 INTEGER(iwp) :: m !: 2532 INTEGER(iwp) :: n !: 2533 2534 REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg), INTENT(INOUT) :: f !: 2535 REAL(wp), DIMENSION(0:cg%nz+1,jcs:jcn,icl:icr), INTENT(IN) :: fc !: 2536 REAL(wp), DIMENSION(nxlg:nxrg), INTENT(IN) :: r1x !: 2537 REAL(wp), DIMENSION(nxlg:nxrg), INTENT(IN) :: r2x !: 2538 REAL(wp), DIMENSION(nysg:nyng), INTENT(IN) :: r1y !: 2539 REAL(wp), DIMENSION(nysg:nyng), INTENT(IN) :: r2y !: 2540 REAL(wp), DIMENSION(nzb:nzt+1), INTENT(IN) :: r1z !: 2541 REAL(wp), DIMENSION(nzb:nzt+1), INTENT(IN) :: r2z !: 2542 2543 REAL(wp) :: fk !: 2544 REAL(wp) :: fkj !: 2545 REAL(wp) :: fkjp !: 2546 REAL(wp) :: fkp !: 2547 REAL(wp) :: fkpj !: 2548 REAL(wp) :: fkpjp !: 2549 REAL(wp) :: logratio !: 2550 REAL(wp) :: logzuc1 !: 2551 REAL(wp) :: zuc1 !: 2552 2553 2554 ib = nxl 2555 ie = nxr 2556 jb = nys 2557 je = nyn 2558 IF ( nest_bound_l ) THEN 2559 ib = nxl - 1 2560 ! 2561 !-- For u, nxl is a ghost node, but not for the other variables 2562 IF ( var == 'u' ) THEN 2563 ib = nxl 2564 ENDIF 2565 ENDIF 2566 IF ( nest_bound_s ) THEN 2567 jb = nys - 1 2568 ! 2569 !-- For v, nys is a ghost node, but not for the other variables 2570 IF ( var == 'v' ) THEN 2571 jb = nys 2572 ENDIF 2573 ENDIF 2574 IF ( nest_bound_r ) THEN 2575 ie = nxr + 1 2576 ENDIF 2577 IF ( nest_bound_n ) THEN 2578 je = nyn + 1 2579 ENDIF 2580 2581 ! 2582 !-- Trilinear interpolation. 2583 DO i = ib, ie 2584 DO j = jb, je 2585 DO k = kb(j,i), nzt + 1 2586 l = ic(i) 2587 m = jc(j) 2588 n = kc(k) 2589 fkj = r1x(i) * fc(n,m,l) + r2x(i) * fc(n,m,l+1) 2590 fkjp = r1x(i) * fc(n,m+1,l) + r2x(i) * fc(n,m+1,l+1) 2591 fkpj = r1x(i) * fc(n+1,m,l) + r2x(i) * fc(n+1,m,l+1) 2592 fkpjp = r1x(i) * fc(n+1,m+1,l) + r2x(i) * fc(n+1,m+1,l+1) 2593 fk = r1y(j) * fkj + r2y(j) * fkjp 2594 fkp = r1y(j) * fkpj + r2y(j) * fkpjp 2595 f(k,j,i) = r1z(k) * fk + r2z(k) * fkp 2596 ENDDO 2597 ENDDO 2598 ENDDO 2599 2600 ! 2601 !-- Correct the interpolated values of u and v in near-wall nodes, i.e. in 2602 !-- the nodes below the coarse-grid nodes with k=1. The corrction is only 2603 !-- made over horizontal wall surfaces in this phase. For the nest boundary 2604 !-- conditions, a corresponding correction is made for all vertical walls, 2605 !-- too. 2606 IF ( var == 'u' .OR. var == 'v' ) THEN 2607 DO i = ib, nxr 2608 DO j = jb, nyn 2609 kbc = 1 2610 ! 2611 !-- kbc is the first coarse-grid point above the surface 2612 DO WHILE ( cg%zu(kbc) < zu(kb(j,i)) ) 2613 kbc = kbc + 1 2614 ENDDO 2615 zuc1 = cg%zu(kbc) 2616 k1 = kb(j,i) + 1 2617 DO WHILE ( zu(k1) < zuc1 ) 2618 k1 = k1 + 1 2619 ENDDO 2620 logzuc1 = LOG( ( zu(k1) - zu(kb(j,i)) ) / z0(j,i) ) 2621 2622 k = kb(j,i) + 1 2623 DO WHILE ( zu(k) < zuc1 ) 2624 logratio = ( LOG( ( zu(k) - zu(kb(j,i)) ) / z0(j,i)) ) / logzuc1 2625 f(k,j,i) = logratio * f(k1,j,i) 2626 k = k + 1 2627 ENDDO 2628 f(kb(j,i),j,i) = 0.0_wp 2629 ENDDO 2630 ENDDO 2631 2632 ELSEIF ( var == 'w' ) THEN 2633 2634 DO i = ib, nxr 2635 DO j = jb, nyn 2636 f(kb(j,i),j,i) = 0.0_wp 2637 ENDDO 2638 ENDDO 2639 2640 ENDIF 2641 2642 END SUBROUTINE pmci_interp_tril_all 2643 2644 #endif 2645 END SUBROUTINE pmci_client_initialize 2646 2647 2648 2649 SUBROUTINE pmci_ensure_nest_mass_conservation 2650 2651 #if defined( __parallel ) 2652 ! 2653 !-- Adjust the volume-flow rate through the top boundary so that the net volume 2654 !-- flow through all boundaries of the current nest domain becomes zero. 2655 IMPLICIT NONE 2656 2657 INTEGER(iwp) :: i !: 2658 INTEGER(iwp) :: ierr !: 2659 INTEGER(iwp) :: j !: 2660 INTEGER(iwp) :: k !: 2661 2662 REAL(wp) :: dxdy !: 2663 REAL(wp) :: innor !: 2664 REAL(wp) :: w_lt !: 2665 REAL(wp), DIMENSION(1:3) :: volume_flow_l !: 2666 2667 ! 2668 !-- Sum up the volume flow through the left/right boundaries 2669 volume_flow(1) = 0.0_wp 2670 volume_flow_l(1) = 0.0_wp 2671 2672 IF ( nest_bound_l ) THEN 2673 i = 0 2674 innor = dy 2675 DO j = nys, nyn 2676 DO k = nzb_u_inner(j,i)+1, nzt 2677 volume_flow_l(1) = volume_flow_l(1) + innor * u(k,j,i) * dzw(k) 2678 ENDDO 2679 ENDDO 2680 ENDIF 2681 2682 IF ( nest_bound_r ) THEN 2683 i = nx + 1 2684 innor = -dy 2685 DO j = nys, nyn 2686 DO k = nzb_u_inner(j,i)+1, nzt 2687 volume_flow_l(1) = volume_flow_l(1) + innor * u(k,j,i) * dzw(k) 2688 ENDDO 2689 ENDDO 2690 ENDIF 2691 2692 #if defined( __parallel ) 2693 IF ( collective_wait ) CALL MPI_BARRIER( comm2d, ierr ) 2694 CALL MPI_ALLREDUCE( volume_flow_l(1), volume_flow(1), 1, MPI_REAL, & 2695 MPI_SUM, comm2d, ierr ) 2696 #else 2697 volume_flow(1) = volume_flow_l(1) 2698 #endif 2699 2700 ! 2701 !-- Sum up the volume flow through the south/north boundaries 2702 volume_flow(2) = 0.0_wp 2703 volume_flow_l(2) = 0.0_wp 2704 2705 IF ( nest_bound_s ) THEN 2706 j = 0 2707 innor = dx 2708 DO i = nxl, nxr 2709 DO k = nzb_v_inner(j,i)+1, nzt 2710 volume_flow_l(2) = volume_flow_l(2) + innor * v(k,j,i) * dzw(k) 2711 ENDDO 2712 ENDDO 2713 ENDIF 2714 2715 IF ( nest_bound_n ) THEN 2716 j = ny + 1 2717 innor = -dx 2718 DO i = nxl, nxr 2719 DO k = nzb_v_inner(j,i)+1, nzt 2720 volume_flow_l(2) = volume_flow_l(2) + innor * v(k,j,i) * dzw(k) 2721 ENDDO 2722 ENDDO 2723 ENDIF 2724 2725 #if defined( __parallel ) 2726 IF ( collective_wait ) CALL MPI_BARRIER( comm2d, ierr ) 2727 CALL MPI_ALLREDUCE( volume_flow_l(2), volume_flow(2), 1, MPI_REAL, & 2728 MPI_SUM, comm2d, ierr ) 2729 #else 2730 volume_flow(2) = volume_flow_l(2) 2731 #endif 2732 2733 ! 2734 !-- Sum up the volume flow through the top boundary 2735 volume_flow(3) = 0.0_wp 2736 volume_flow_l(3) = 0.0_wp 2737 dxdy = dx * dy 2738 k = nzt 2739 DO i = nxl, nxr 2740 DO j = nys, nyn 2741 volume_flow_l(3) = volume_flow_l(3) - w(k,j,i) * dxdy 2742 ENDDO 2743 ENDDO 2744 2745 #if defined( __parallel ) 2746 IF ( collective_wait ) CALL MPI_BARRIER( comm2d, ierr ) 2747 CALL MPI_ALLREDUCE( volume_flow_l(3), volume_flow(3), 1, MPI_REAL, & 2748 MPI_SUM, comm2d, ierr ) 2749 #else 2750 volume_flow(3) = volume_flow_l(3) 2751 #endif 2752 2753 ! 2754 !-- Correct the top-boundary value of w 2755 w_lt = (volume_flow(1) + volume_flow(2) + volume_flow(3)) / area_t 2756 DO i = nxl, nxr 2757 DO j = nys, nyn 2758 DO k = nzt, nzt + 1 2759 w(k,j,i) = w(k,j,i) + w_lt 2760 ENDDO 2761 ENDDO 2762 ENDDO 2763 2764 #endif 2765 END SUBROUTINE pmci_ensure_nest_mass_conservation 2766 2767 2768 !-- TO_DO: the timestep sycnchronization could be done easier using 2769 !-- an MPI_ALLREDUCE with MIN over MPI_COMM_WORLD 2175 2770 SUBROUTINE pmci_server_synchronize 2176 2771 2177 2772 #if defined( __parallel ) 2178 2773 ! 2179 !-- Unify the time steps for each model and synchronize. 2180 !-- This is based on the assumption that the native time step 2181 !-- (original dt_3d) of any server is always larger than the smallest 2182 !-- native time step of it s clients. 2774 !-- Unify the time steps for each model and synchronize. This is based on the 2775 !-- assumption that the native time step (original dt_3d) of any server is 2776 !-- always larger than the smallest native time step of it s clients. 2183 2777 IMPLICIT NONE 2184 INTEGER(iwp) :: client_id !: 2778 2779 INTEGER(iwp) :: client_id !: 2780 INTEGER(iwp) :: ierr !: 2781 INTEGER(iwp) :: m !: 2782 2185 2783 REAL(wp), DIMENSION(1) :: dtc !: 2186 2784 REAL(wp), DIMENSION(1) :: dtl !: 2187 INTEGER(iwp) :: ierr !: 2188 INTEGER(iwp) :: m !: 2189 2190 ! 2191 !-- First find the smallest native time step of all the clients of the current server. 2785 2786 ! 2787 !-- First find the smallest native time step of all the clients of the current 2788 !-- server. 2192 2789 dtl(1) = 999999.9_wp 2193 DO m = 1, SIZE( PMC_Server_for_Client ) -12790 DO m = 1, SIZE( PMC_Server_for_Client )-1 2194 2791 client_id = PMC_Server_for_Client(m) 2195 2792 IF ( myid == 0 ) THEN … … 2201 2798 2202 2799 ! 2203 !-- Broadcast the unified time step to all server processes .2800 !-- Broadcast the unified time step to all server processes 2204 2801 CALL MPI_BCAST( dt_3d, 1, MPI_REAL, 0, comm2d, ierr ) 2205 2802 2206 2803 ! 2207 !-- Send the new time step to all the clients of the current server .2804 !-- Send the new time step to all the clients of the current server 2208 2805 DO m = 1, SIZE( PMC_Server_for_Client ) - 1 2209 2806 client_id = PMC_Server_for_Client(m) … … 2222 2819 #if defined( __parallel ) 2223 2820 ! 2224 !-- Unify the time steps for each model and synchronize. 2225 !-- This is based on the assumption that the native time step 2226 !-- (original dt_3d) of any server is always larger than the smallest 2227 !-- native time step of it s clients. 2821 !-- Unify the time steps for each model and synchronize. This is based on the 2822 !-- assumption that the native time step (original dt_3d) of any server is 2823 !-- always larger than the smallest native time step of it s clients. 2228 2824 2229 2825 IMPLICIT NONE 2826 2827 INTEGER(iwp) :: ierr !: 2828 2230 2829 REAL(wp), DIMENSION(1) :: dtl !: 2231 2830 REAL(wp), DIMENSION(1) :: dts !: 2232 INTEGER(iwp) :: ierr !:2233 2831 2234 2832 2235 2833 dtl(1) = dt_3d 2236 IF ( cpl_id > 1 ) THEN ! Root id is never a client2834 IF ( cpl_id > 1 ) THEN 2237 2835 IF ( myid==0 ) THEN 2238 2836 CALL pmc_send_to_server( dtl, SIZE( dtl ), 0, 101, ierr ) … … 2242 2840 2243 2841 ! 2244 !-- Broadcast the unified time step to all server processes .2842 !-- Broadcast the unified time step to all server processes 2245 2843 CALL MPI_BCAST( dt_3d, 1, MPI_REAL, 0, comm2d, ierr ) 2246 2844 ENDIF … … 2255 2853 IMPLICIT NONE 2256 2854 2257 INTEGER(iwp),INTENT(IN) :: swaplevel !: swaplevel (1 or 2) of PALM's timestep 2855 INTEGER(iwp),INTENT(IN) :: swaplevel !: swaplevel (1 or 2) of PALM's 2856 !: timestep 2258 2857 2259 2858 INTEGER(iwp) :: client_id !: … … 2278 2877 2279 2878 #if defined( __parallel ) 2280 INTEGER(iwp) :: client_id !: 2281 INTEGER(iwp) :: i !: 2282 INTEGER(iwp) :: j !: 2283 INTEGER(iwp) :: ierr !: 2284 INTEGER(iwp) :: m !: 2285 REAL(wp) :: waittime !: 2286 REAL(wp), DIMENSION(1) :: dtc !: 2287 REAL(wp), DIMENSION(1) :: dtl !: 2288 2289 ! 2290 !-- First find the smallest native time step of all the clients of the current server. 2879 INTEGER(iwp) :: client_id !: 2880 INTEGER(iwp) :: i !: 2881 INTEGER(iwp) :: j !: 2882 INTEGER(iwp) :: ierr !: 2883 INTEGER(iwp) :: m !: 2884 2885 REAL(wp) :: waittime !: 2886 REAL(wp), DIMENSION(1) :: dtc !: 2887 REAL(wp), DIMENSION(1) :: dtl !: 2888 2889 ! 2890 !-- First find the smallest native time step of all the clients of the current 2891 !-- server. 2291 2892 dtl(1) = 999999.9_wp 2292 DO m = 1, SIZE( PMC_Server_for_Client ) -12293 client_id = PMC_Server_for_Client(m)2893 DO m = 1, SIZE( pmc_server_for_client )-1 2894 client_id = pmc_server_for_client(m) 2294 2895 IF ( myid==0 ) THEN 2295 2896 CALL pmc_recv_from_client( client_id, dtc, SIZE( dtc ), 0, 101, ierr ) … … 2300 2901 2301 2902 ! 2302 !-- Broadcast the unified time step to all server processes .2903 !-- Broadcast the unified time step to all server processes 2303 2904 CALL MPI_BCAST( dt_3d, 1, MPI_REAL, 0, comm2d, ierr ) 2304 2905 2305 DO m = 1, SIZE( PMC_Server_for_Client ) -12906 DO m = 1, SIZE( PMC_Server_for_Client )-1 2306 2907 client_id = PMC_Server_for_Client(m) 2307 CALL cpu_log( log_point_s(70), ' PMC modelsync', 'start' )2308 2309 ! 2310 !-- Send the new time step to all the clients of the current server .2311 IF ( myid == 0 ) THEN2908 CALL cpu_log( log_point_s(70), 'pmc sync', 'start' ) 2909 2910 ! 2911 !-- Send the new time step to all the clients of the current server 2912 IF ( myid == 0 ) THEN 2312 2913 CALL pmc_send_to_client( client_id, dtl, SIZE( dtl ), 0, 102, ierr ) 2313 2914 ENDIF 2314 CALL cpu_log( log_point_s(70), ' PMC modelsync', 'stop' )2915 CALL cpu_log( log_point_s(70), 'pmc sync', 'stop' ) 2315 2916 2316 2917 IF ( direction == server_to_client ) THEN 2317 CALL cpu_log( log_point_s(71), ' PMC Server Send', 'start' )2918 CALL cpu_log( log_point_s(71), 'pmc server send', 'start' ) 2318 2919 CALL pmc_s_fillbuffer( client_id, waittime=waittime ) 2319 CALL cpu_log( log_point_s(71), 'PMC Server Send', 'stop' ) 2320 ELSE ! Communication from client to server. 2321 CALL cpu_log( log_point_s(72), 'PMC Server Recv', 'start' ) 2920 CALL cpu_log( log_point_s(71), 'pmc server send', 'stop' ) 2921 ELSE 2922 ! 2923 !-- Communication from client to server 2924 CALL cpu_log( log_point_s(72), 'pmc server recv', 'start' ) 2322 2925 client_id = pmc_server_for_client(m) 2323 2926 CALL pmc_s_getdata_from_buffer( client_id ) 2324 CALL cpu_log( log_point_s(72), ' PMC Server Recv', 'stop' )2325 2326 ! 2327 !-- The anterpolated data is now available in u etc .2927 CALL cpu_log( log_point_s(72), 'pmc server recv', 'stop' ) 2928 2929 ! 2930 !-- The anterpolated data is now available in u etc 2328 2931 IF ( topography /= 'flat' ) THEN 2329 2932 2330 2933 ! 2331 2934 !-- Inside buildings/topography reset velocities and TKE back to zero. 2332 !-- Other scalars (pt, q, s, km, kh, p, sa, ...) are ignored at present, 2333 !-- maybe revise later. 2935 !-- TO_DO: at least temperature should be included here immediately 2936 !-- Other scalars (pt, q, s, km, kh, p, sa, ...) are ignored at 2937 !-- present, maybe revise later. 2334 2938 DO i = nxlg, nxrg 2335 2939 DO j = nysg, nyng … … 2337 2941 v(nzb:nzb_v_inner(j,i),j,i) = 0.0_wp 2338 2942 w(nzb:nzb_w_inner(j,i),j,i) = 0.0_wp 2339 e(nzb:nzb_ w_inner(j,i),j,i) = 0.0_wp2943 e(nzb:nzb_s_inner(j,i),j,i) = 0.0_wp 2340 2944 ENDDO 2341 2945 ENDDO … … 2356 2960 2357 2961 #if defined( __parallel ) 2358 INTEGER(iwp) 2359 INTEGER(iwp) 2360 INTEGER(iwp) 2361 INTEGER(iwp) 2362 INTEGER(iwp) 2962 INTEGER(iwp) :: ierr !: 2963 INTEGER(iwp) :: icl !: 2964 INTEGER(iwp) :: icr !: 2965 INTEGER(iwp) :: jcs !: 2966 INTEGER(iwp) :: jcn !: 2363 2967 2364 REAL(wp), DIMENSION(1) :: dtl !: 2365 REAL(wp), DIMENSION(1) :: dts !: 2366 REAL(wp) :: waittime !: 2968 REAL(wp) :: waittime !: 2969 2970 REAL(wp), DIMENSION(1) :: dtl !: 2971 REAL(wp), DIMENSION(1) :: dts !: 2367 2972 2368 2973 2369 2974 dtl = dt_3d 2370 IF ( cpl_id > 1 ) THEN ! Root id is never a client2371 CALL cpu_log( log_point_s(70), ' PMC modelsync', 'start' )2975 IF ( cpl_id > 1 ) THEN 2976 CALL cpu_log( log_point_s(70), 'pmc sync', 'start' ) 2372 2977 IF ( myid==0 ) THEN 2373 2978 CALL pmc_send_to_server( dtl, SIZE( dtl ), 0, 101, ierr ) … … 2379 2984 !-- Broadcast the unified time step to all server processes. 2380 2985 CALL MPI_BCAST( dt_3d, 1, MPI_REAL, 0, comm2d, ierr ) 2381 CALL cpu_log( log_point_s(70), ' PMC modelsync', 'stop' )2986 CALL cpu_log( log_point_s(70), 'pmc sync', 'stop' ) 2382 2987 2383 2988 ! … … 2389 2994 2390 2995 IF ( direction == server_to_client ) THEN 2391 CALL cpu_log( log_point_s(73), 'PMC Client Recv', 'start' ) 2996 2997 CALL cpu_log( log_point_s(73), 'pmc client recv', 'start' ) 2392 2998 CALL pmc_c_getbuffer( WaitTime = WaitTime ) 2393 CALL cpu_log( log_point_s(73), ' PMC Client Recv', 'stop' )2999 CALL cpu_log( log_point_s(73), 'pmc client recv', 'stop' ) 2394 3000 2395 3001 CALL pmci_interpolation 2396 3002 2397 ELSE ! IF ( direction == server_to_client ) 2398 3003 ELSE 3004 ! 3005 !-- direction == server_to_client 2399 3006 CALL pmci_anterpolation 2400 3007 2401 CALL cpu_log( log_point_s(74), ' PMC Client Send', 'start' )3008 CALL cpu_log( log_point_s(74), 'pmc client send', 'start' ) 2402 3009 CALL pmc_c_putbuffer( WaitTime = WaitTime ) 2403 CALL cpu_log( log_point_s(74), 'PMC Client Send', 'stop' ) 3010 CALL cpu_log( log_point_s(74), 'pmc client send', 'stop' ) 3011 2404 3012 ENDIF 2405 3013 ENDIF … … 2407 3015 CONTAINS 2408 3016 2409 2410 SUBROUTINE pmci_interpolation 2411 2412 ! 2413 !-- A wrapper routine for all interpolation and extrapolation actions. 2414 IMPLICIT NONE 2415 2416 ! 2417 !-- Add IF-condition here: IF not vertical nesting 2418 IF ( nest_bound_l ) THEN ! Left border pe 2419 CALL pmci_interp_tril_lr( u, uc, icu, jco, kco, r1xu, r2xu, r1yo, r2yo, r1zo, r2zo, & 2420 nzb_u_inner, logc_u_l, logc_ratio_u_l, nzt_topo_nestbc_l, 'l', 'u' ) 2421 CALL pmci_interp_tril_lr( v, vc, ico, jcv, kco, r1xo, r2xo, r1yv, r2yv, r1zo, r2zo, & 2422 nzb_v_inner, logc_v_l, logc_ratio_v_l, nzt_topo_nestbc_l, 'l', 'v' ) 2423 CALL pmci_interp_tril_lr( w, wc, ico, jco, kcw, r1xo, r2xo, r1yo, r2yo, r1zw, r2zw, & 2424 nzb_w_inner, logc_w_l, logc_ratio_w_l, nzt_topo_nestbc_l, 'l', 'w' ) 2425 CALL pmci_interp_tril_lr( e, ec, ico, jco, kco, r1xo, r2xo, r1yo, r2yo, r1zo, r2zo, & 2426 nzb_s_inner, logc_u_l, logc_ratio_u_l, nzt_topo_nestbc_l, 'l', 'e' ) 2427 CALL pmci_interp_tril_lr( pt, ptc, ico, jco, kco, r1xo, r2xo, r1yo, r2yo, r1zo, r2zo, & 2428 nzb_s_inner, logc_u_l, logc_ratio_u_l, nzt_topo_nestbc_l, 'l', 's' ) 2429 IF ( humidity .OR. passive_scalar ) THEN 2430 CALL pmci_interp_tril_lr( q, qc, ico, jco, kco, r1xo, r2xo, r1yo, r2yo, r1zo, r2zo, & 2431 nzb_s_inner, logc_u_l, logc_ratio_u_l, nzt_topo_nestbc_l, 'l', 's' ) 2432 ENDIF 2433 IF ( nesting_mode == 'one-way' ) THEN 2434 CALL pmci_extrap_ifoutflow_lr( u, nzb_u_inner, 'l', 'u' ) 2435 CALL pmci_extrap_ifoutflow_lr( v, nzb_v_inner, 'l', 'v' ) 2436 CALL pmci_extrap_ifoutflow_lr( w, nzb_w_inner, 'l', 'w' ) 2437 CALL pmci_extrap_ifoutflow_lr( e, nzb_s_inner, 'l', 'e' ) 2438 CALL pmci_extrap_ifoutflow_lr( pt,nzb_s_inner, 'l', 's' ) 2439 IF ( humidity .OR. passive_scalar ) THEN 2440 CALL pmci_extrap_ifoutflow_lr( q, nzb_s_inner, 'l', 's' ) 2441 ENDIF 2442 ENDIF 2443 ENDIF 2444 IF ( nest_bound_r ) THEN ! Right border pe 2445 CALL pmci_interp_tril_lr( u, uc, icu, jco, kco, r1xu, r2xu, r1yo, r2yo, r1zo, r2zo, & 2446 nzb_u_inner, logc_u_r, logc_ratio_u_r, nzt_topo_nestbc_r, 'r', 'u' ) 2447 CALL pmci_interp_tril_lr( v, vc, ico, jcv, kco, r1xo, r2xo, r1yv, r2yv, r1zo, r2zo, & 2448 nzb_v_inner, logc_v_r, logc_ratio_v_r, nzt_topo_nestbc_r, 'r', 'v' ) 2449 CALL pmci_interp_tril_lr( w, wc, ico, jco, kcw, r1xo, r2xo, r1yo, r2yo, r1zw, r2zw, & 2450 nzb_w_inner, logc_w_r, logc_ratio_w_r, nzt_topo_nestbc_r, 'r', 'w' ) 2451 CALL pmci_interp_tril_lr( e, ec, ico, jco, kco, r1xo, r2xo, r1yo, r2yo, r1zo, r2zo, & 2452 nzb_s_inner, logc_u_r, logc_ratio_u_r, nzt_topo_nestbc_r, 'r', 'e' ) 2453 CALL pmci_interp_tril_lr( pt, ptc, ico, jco, kco, r1xo, r2xo, r1yo, r2yo, r1zo, r2zo, & 2454 nzb_s_inner, logc_u_r, logc_ratio_u_r, nzt_topo_nestbc_r, 'r', 's' ) 2455 IF ( humidity .OR. passive_scalar ) THEN 2456 CALL pmci_interp_tril_lr( q, qc, ico, jco, kco, r1xo, r2xo, r1yo, r2yo, r1zo, r2zo, & 2457 nzb_s_inner, logc_u_r, logc_ratio_u_r, nzt_topo_nestbc_r, 'r', 's' ) 2458 ENDIF 2459 IF ( nesting_mode == 'one-way' ) THEN 2460 CALL pmci_extrap_ifoutflow_lr( u, nzb_u_inner, 'r', 'u' ) 3017 SUBROUTINE pmci_interpolation 3018 3019 ! 3020 !-- A wrapper routine for all interpolation and extrapolation actions 3021 IMPLICIT NONE 3022 3023 ! 3024 !-- Add IF-condition here: IF not vertical nesting 3025 !-- Left border pe: 3026 IF ( nest_bound_l ) THEN 3027 CALL pmci_interp_tril_lr( u, uc, icu, jco, kco, r1xu, r2xu, r1yo, & 3028 r2yo, r1zo, r2zo, nzb_u_inner, logc_u_l, & 3029 logc_ratio_u_l, nzt_topo_nestbc_l, 'l', & 3030 'u' ) 3031 CALL pmci_interp_tril_lr( v, vc, ico, jcv, kco, r1xo, r2xo, r1yv, & 3032 r2yv, r1zo, r2zo, nzb_v_inner, logc_v_l, & 3033 logc_ratio_v_l, nzt_topo_nestbc_l, 'l', & 3034 'v' ) 3035 CALL pmci_interp_tril_lr( w, wc, ico, jco, kcw, r1xo, r2xo, r1yo, & 3036 r2yo, r1zw, r2zw, nzb_w_inner, logc_w_l, & 3037 logc_ratio_w_l, nzt_topo_nestbc_l, 'l', & 3038 'w' ) 3039 CALL pmci_interp_tril_lr( e, ec, ico, jco, kco, r1xo, r2xo, r1yo, & 3040 r2yo, r1zo, r2zo, nzb_s_inner, logc_u_l, & 3041 logc_ratio_u_l, nzt_topo_nestbc_l, 'l', & 3042 'e' ) 3043 CALL pmci_interp_tril_lr( pt, ptc, ico, jco, kco, r1xo, r2xo, r1yo, & 3044 r2yo, r1zo, r2zo, nzb_s_inner, logc_u_l, & 3045 logc_ratio_u_l, nzt_topo_nestbc_l, 'l', & 3046 's' ) 3047 IF ( humidity .OR. passive_scalar ) THEN 3048 CALL pmci_interp_tril_lr( q, qc, ico, jco, kco, r1xo, r2xo, r1yo, & 3049 r2yo, r1zo, r2zo, nzb_s_inner, & 3050 logc_u_l, logc_ratio_u_l, & 3051 nzt_topo_nestbc_l, 'l', 's' ) 3052 ENDIF 3053 3054 IF ( nesting_mode == 'one-way' ) THEN 3055 CALL pmci_extrap_ifoutflow_lr( u, nzb_u_inner, 'l', 'u' ) 3056 CALL pmci_extrap_ifoutflow_lr( v, nzb_v_inner, 'l', 'v' ) 3057 CALL pmci_extrap_ifoutflow_lr( w, nzb_w_inner, 'l', 'w' ) 3058 CALL pmci_extrap_ifoutflow_lr( e, nzb_s_inner, 'l', 'e' ) 3059 CALL pmci_extrap_ifoutflow_lr( pt,nzb_s_inner, 'l', 's' ) 3060 IF ( humidity .OR. passive_scalar ) THEN 3061 CALL pmci_extrap_ifoutflow_lr( q, nzb_s_inner, 'l', 's' ) 3062 ENDIF 3063 ENDIF 3064 3065 ENDIF 3066 ! 3067 !-- Right border pe 3068 IF ( nest_bound_r ) THEN 3069 CALL pmci_interp_tril_lr( u, uc, icu, jco, kco, r1xu, r2xu, r1yo, & 3070 r2yo, r1zo, r2zo, nzb_u_inner, logc_u_r, & 3071 logc_ratio_u_r, nzt_topo_nestbc_r, 'r', & 3072 'u' ) 3073 CALL pmci_interp_tril_lr( v, vc, ico, jcv, kco, r1xo, r2xo, r1yv, & 3074 r2yv, r1zo, r2zo, nzb_v_inner, logc_v_r, & 3075 logc_ratio_v_r, nzt_topo_nestbc_r, 'r', & 3076 'v' ) 3077 CALL pmci_interp_tril_lr( w, wc, ico, jco, kcw, r1xo, r2xo, r1yo, & 3078 r2yo, r1zw, r2zw, nzb_w_inner, logc_w_r, & 3079 logc_ratio_w_r, nzt_topo_nestbc_r, 'r', & 3080 'w' ) 3081 CALL pmci_interp_tril_lr( e, ec, ico, jco, kco, r1xo, r2xo, r1yo, & 3082 r2yo, r1zo, r2zo, nzb_s_inner, logc_u_r, & 3083 logc_ratio_u_r, nzt_topo_nestbc_r, 'r', & 3084 'e' ) 3085 CALL pmci_interp_tril_lr( pt, ptc, ico, jco, kco, r1xo, r2xo, r1yo, & 3086 r2yo, r1zo, r2zo, nzb_s_inner, logc_u_r, & 3087 logc_ratio_u_r, nzt_topo_nestbc_r, 'r', & 3088 's' ) 3089 IF ( humidity .OR. passive_scalar ) THEN 3090 CALL pmci_interp_tril_lr( q, qc, ico, jco, kco, r1xo, r2xo, r1yo, & 3091 r2yo, r1zo, r2zo, nzb_s_inner, & 3092 logc_u_r, logc_ratio_u_r, & 3093 nzt_topo_nestbc_r, 'r', 's' ) 3094 ENDIF 3095 3096 IF ( nesting_mode == 'one-way' ) THEN 3097 CALL pmci_extrap_ifoutflow_lr( u, nzb_u_inner, 'r', 'u' ) 2461 3098 CALL pmci_extrap_ifoutflow_lr( v, nzb_v_inner, 'r', 'v' ) 2462 3099 CALL pmci_extrap_ifoutflow_lr( w, nzb_w_inner, 'r', 'w' ) … … 2467 3104 ENDIF 2468 3105 ENDIF 2469 ENDIF 2470 IF ( nest_bound_s ) THEN ! South border pe 2471 CALL pmci_interp_tril_sn( u, uc, icu, jco, kco, r1xu, r2xu, r1yo, r2yo, r1zo, r2zo, & 2472 nzb_u_inner, logc_u_s, logc_ratio_u_s, nzt_topo_nestbc_s, 's', 'u' ) 2473 CALL pmci_interp_tril_sn( v, vc, ico, jcv, kco, r1xo, r2xo, r1yv, r2yv, r1zo, r2zo, & 2474 nzb_v_inner, logc_v_s, logc_ratio_v_s, nzt_topo_nestbc_s, 's', 'v' ) 2475 CALL pmci_interp_tril_sn( w, wc, ico, jco, kcw, r1xo, r2xo, r1yo, r2yo, r1zw, r2zw, & 2476 nzb_w_inner, logc_w_s, logc_ratio_w_s, nzt_topo_nestbc_s, 's', 'w' ) 2477 CALL pmci_interp_tril_sn( e, ec, ico, jco, kco, r1xo, r2xo, r1yo, r2yo, r1zo, r2zo, & 2478 nzb_s_inner, logc_u_s, logc_ratio_u_s, nzt_topo_nestbc_s, 's', 'e' ) 2479 CALL pmci_interp_tril_sn( pt, ptc, ico, jco, kco, r1xo, r2xo, r1yo, r2yo, r1zo, r2zo, & 2480 nzb_s_inner, logc_u_s, logc_ratio_u_s, nzt_topo_nestbc_s, 's', 's' ) 3106 3107 ENDIF 3108 ! 3109 !-- South border pe 3110 IF ( nest_bound_s ) THEN 3111 CALL pmci_interp_tril_sn( u, uc, icu, jco, kco, r1xu, r2xu, r1yo, & 3112 r2yo, r1zo, r2zo, nzb_u_inner, logc_u_s, & 3113 logc_ratio_u_s, nzt_topo_nestbc_s, 's', & 3114 'u' ) 3115 CALL pmci_interp_tril_sn( v, vc, ico, jcv, kco, r1xo, r2xo, r1yv, & 3116 r2yv, r1zo, r2zo, nzb_v_inner, logc_v_s, & 3117 logc_ratio_v_s, nzt_topo_nestbc_s, 's', & 3118 'v' ) 3119 CALL pmci_interp_tril_sn( w, wc, ico, jco, kcw, r1xo, r2xo, r1yo, & 3120 r2yo, r1zw, r2zw, nzb_w_inner, logc_w_s, & 3121 logc_ratio_w_s, nzt_topo_nestbc_s, 's', & 3122 'w' ) 3123 CALL pmci_interp_tril_sn( e, ec, ico, jco, kco, r1xo, r2xo, r1yo, & 3124 r2yo, r1zo, r2zo, nzb_s_inner, logc_u_s, & 3125 logc_ratio_u_s, nzt_topo_nestbc_s, 's', & 3126 'e' ) 3127 CALL pmci_interp_tril_sn( pt, ptc, ico, jco, kco, r1xo, r2xo, r1yo, & 3128 r2yo, r1zo, r2zo, nzb_s_inner, logc_u_s, & 3129 logc_ratio_u_s, nzt_topo_nestbc_s, 's', & 3130 's' ) 2481 3131 IF ( humidity .OR. passive_scalar ) THEN 2482 CALL pmci_interp_tril_sn( q, qc, ico, jco, kco, r1xo, r2xo, r1yo, r2yo, r1zo, r2zo, & 2483 nzb_s_inner, logc_u_s, logc_ratio_u_s, nzt_topo_nestbc_s, 's', 's' ) 3132 CALL pmci_interp_tril_sn( q, qc, ico, jco, kco, r1xo, r2xo, r1yo, & 3133 r2yo, r1zo, r2zo, nzb_s_inner, & 3134 logc_u_s, logc_ratio_u_s, & 3135 nzt_topo_nestbc_s, 's', 's' ) 2484 3136 ENDIF 3137 2485 3138 IF ( nesting_mode == 'one-way' ) THEN 2486 3139 CALL pmci_extrap_ifoutflow_sn( u, nzb_u_inner, 's', 'u' ) … … 2493 3146 ENDIF 2494 3147 ENDIF 2495 ENDIF 2496 IF ( nest_bound_n ) THEN ! North border pe 2497 CALL pmci_interp_tril_sn( u, uc, icu, jco, kco, r1xu, r2xu, r1yo, r2yo, r1zo, r2zo, & 2498 nzb_u_inner, logc_u_n, logc_ratio_u_n, nzt_topo_nestbc_n, 'n', 'u' ) 2499 CALL pmci_interp_tril_sn( v, vc, ico, jcv, kco, r1xo, r2xo, r1yv, r2yv, r1zo, r2zo, & 2500 nzb_v_inner, logc_v_n, logc_ratio_v_n, nzt_topo_nestbc_n, 'n', 'v' ) 2501 CALL pmci_interp_tril_sn( w, wc, ico, jco, kcw, r1xo, r2xo, r1yo, r2yo, r1zw, r2zw, & 2502 nzb_w_inner, logc_w_n, logc_ratio_w_n, nzt_topo_nestbc_n, 'n', 'w' ) 2503 CALL pmci_interp_tril_sn( e, ec, ico,jco,kco,r1xo,r2xo,r1yo,r2yo,r1zo,r2zo, & 2504 nzb_s_inner, logc_u_n, logc_ratio_u_n, nzt_topo_nestbc_n, 'n', 'e' ) 2505 CALL pmci_interp_tril_sn( pt, ptc, ico, jco, kco, r1xo, r2xo, r1yo, r2yo, r1zo, r2zo, & 2506 nzb_s_inner, logc_u_n, logc_ratio_u_n, nzt_topo_nestbc_n, 'n', 's' ) 3148 3149 ENDIF 3150 ! 3151 !-- North border pe 3152 IF ( nest_bound_n ) THEN 3153 CALL pmci_interp_tril_sn( u, uc, icu, jco, kco, r1xu, r2xu, r1yo, & 3154 r2yo, r1zo, r2zo, nzb_u_inner, logc_u_n, & 3155 logc_ratio_u_n, nzt_topo_nestbc_n, 'n', & 3156 'u' ) 3157 CALL pmci_interp_tril_sn( v, vc, ico, jcv, kco, r1xo, r2xo, r1yv, & 3158 r2yv, r1zo, r2zo, nzb_v_inner, logc_v_n, & 3159 logc_ratio_v_n, nzt_topo_nestbc_n, 'n', & 3160 'v' ) 3161 CALL pmci_interp_tril_sn( w, wc, ico, jco, kcw, r1xo, r2xo, r1yo, & 3162 r2yo, r1zw, r2zw, nzb_w_inner, logc_w_n, & 3163 logc_ratio_w_n, nzt_topo_nestbc_n, 'n', & 3164 'w' ) 3165 CALL pmci_interp_tril_sn( e, ec, ico, jco, kco, r1xo, r2xo, r1yo, & 3166 r2yo, r1zo, r2zo, nzb_s_inner, logc_u_n, & 3167 logc_ratio_u_n, nzt_topo_nestbc_n, 'n', & 3168 'e' ) 3169 CALL pmci_interp_tril_sn( pt, ptc, ico, jco, kco, r1xo, r2xo, r1yo, & 3170 r2yo, r1zo, r2zo, nzb_s_inner, logc_u_n, & 3171 logc_ratio_u_n, nzt_topo_nestbc_n, 'n', & 3172 's' ) 2507 3173 IF ( humidity .OR. passive_scalar ) THEN 2508 CALL pmci_interp_tril_sn( q, qc, ico,jco,kco,r1xo,r2xo,r1yo,r2yo,r1zo,r2zo, & 2509 nzb_s_inner, logc_u_n, logc_ratio_u_n, nzt_topo_nestbc_n, 'n', 's' ) 3174 CALL pmci_interp_tril_sn( q, qc, ico, jco, kco, r1xo, r2xo, r1yo, & 3175 r2yo, r1zo, r2zo, nzb_s_inner, & 3176 logc_u_n, logc_ratio_u_n, & 3177 nzt_topo_nestbc_n, 'n', 's' ) 2510 3178 ENDIF 3179 2511 3180 IF ( nesting_mode == 'one-way' ) THEN 2512 3181 CALL pmci_extrap_ifoutflow_sn( u, nzb_u_inner, 'n', 'u' ) … … 2518 3187 CALL pmci_extrap_ifoutflow_sn( q, nzb_s_inner, 'n', 's' ) 2519 3188 ENDIF 3189 2520 3190 ENDIF 2521 ENDIF 2522 2523 ! 2524 !-- All PEs are top-border PEs 2525 CALL pmci_interp_tril_t( u, uc, icu, jco, kco, r1xu, r2xu, r1yo, r2yo, r1zo, r2zo, 'u' ) 2526 CALL pmci_interp_tril_t( v, vc, ico, jcv, kco, r1xo, r2xo, r1yv, r2yv, r1zo, r2zo, 'v' ) 2527 CALL pmci_interp_tril_t( w, wc, ico, jco, kcw, r1xo, r2xo, r1yo, r2yo, r1zw, r2zw, 'w' ) 2528 CALL pmci_interp_tril_t( e, ec, ico, jco, kco, r1xo, r2xo, r1yo, r2yo, r1zo, r2zo, 'e' ) 2529 CALL pmci_interp_tril_t( pt, ptc, ico, jco, kco, r1xo, r2xo, r1yo, r2yo, r1zo, r2zo, 's' ) 2530 IF ( humidity .OR. passive_scalar ) THEN 2531 CALL pmci_interp_tril_t( q, qc, ico, jco, kco, r1xo, r2xo, r1yo, r2yo, r1zo, r2zo, 's' ) 2532 ENDIF 3191 3192 ENDIF 3193 3194 ! 3195 !-- All PEs are top-border PEs 3196 CALL pmci_interp_tril_t( u, uc, icu, jco, kco, r1xu, r2xu, r1yo, & 3197 r2yo, r1zo, r2zo, 'u' ) 3198 CALL pmci_interp_tril_t( v, vc, ico, jcv, kco, r1xo, r2xo, r1yv, & 3199 r2yv, r1zo, r2zo, 'v' ) 3200 CALL pmci_interp_tril_t( w, wc, ico, jco, kcw, r1xo, r2xo, r1yo, & 3201 r2yo, r1zw, r2zw, 'w' ) 3202 CALL pmci_interp_tril_t( e, ec, ico, jco, kco, r1xo, r2xo, r1yo, & 3203 r2yo, r1zo, r2zo, 'e' ) 3204 CALL pmci_interp_tril_t( pt, ptc, ico, jco, kco, r1xo, r2xo, r1yo, & 3205 r2yo, r1zo, r2zo, 's' ) 3206 IF ( humidity .OR. passive_scalar ) THEN 3207 CALL pmci_interp_tril_t( q, qc, ico, jco, kco, r1xo, r2xo, r1yo, & 3208 r2yo, r1zo, r2zo, 's' ) 3209 ENDIF 3210 2533 3211 IF ( nesting_mode == 'one-way' ) THEN 2534 3212 CALL pmci_extrap_ifoutflow_t( u, 'u' ) … … 2541 3219 ENDIF 2542 3220 ENDIF 3221 2543 3222 END SUBROUTINE pmci_interpolation 2544 3223 … … 2551 3230 IMPLICIT NONE 2552 3231 2553 CALL pmci_anterp_tophat( u, uc, kceu, iflu, ifuu, jflo, jfuo, kflo, kfuo, nzb_u_inner, 'u' ) 2554 CALL pmci_anterp_tophat( v, vc, kceu, iflo, ifuo, jflv, jfuv, kflo, kfuo, nzb_v_inner, 'v' ) 2555 CALL pmci_anterp_tophat( w, wc, kcew, iflo, ifuo, jflo, jfuo, kflw, kfuw, nzb_w_inner, 'w' ) 2556 CALL pmci_anterp_tophat( pt, ptc, kceu, iflo, ifuo, jflo, jfuo, kflo, kfuo, nzb_s_inner, 's' ) 3232 CALL pmci_anterp_tophat( u, uc, kctu, iflu, ifuu, jflo, jfuo, kflo, & 3233 kfuo, nzb_u_inner, 'u' ) 3234 CALL pmci_anterp_tophat( v, vc, kctu, iflo, ifuo, jflv, jfuv, kflo, & 3235 kfuo, nzb_v_inner, 'v' ) 3236 CALL pmci_anterp_tophat( w, wc, kctw, iflo, ifuo, jflo, jfuo, kflw, & 3237 kfuw, nzb_w_inner, 'w' ) 3238 CALL pmci_anterp_tophat( pt, ptc, kctu, iflo, ifuo, jflo, jfuo, kflo, & 3239 kfuo, nzb_s_inner, 's' ) 2557 3240 IF ( humidity .OR. passive_scalar ) THEN 2558 CALL pmci_anterp_tophat( q, qc, kceu, iflo, ifuo, jflo, jfuo, kflo, kfuo, nzb_s_inner, 's' ) 3241 CALL pmci_anterp_tophat( q, qc, kctu, iflo, ifuo, jflo, jfuo, kflo, & 3242 kfuo, nzb_s_inner, 's' ) 2559 3243 ENDIF 2560 3244 … … 2563 3247 2564 3248 2565 SUBROUTINE pmci_interp_tril_lr( f, fc, ic, jc, kc, r1x, r2x, r1y, r2y, r1z, r2z, kb, logc, logc_ratio,&2566 nzt_topo_nestbc, edge, var )2567 3249 SUBROUTINE pmci_interp_tril_lr( f, fc, ic, jc, kc, r1x, r2x, r1y, r2y, r1z, & 3250 r2z, kb, logc, logc_ratio, nzt_topo_nestbc, & 3251 edge, var ) 2568 3252 ! 2569 3253 !-- Interpolation of ghost-node values used as the client-domain boundary 2570 !-- conditions. This subroutine handles the left and right boundaries. 2571 !-- This subroutine is based on trilinear interpolation. 2572 !-- Constant dz is still assumed. 2573 ! 2574 !-- Antti Hellsten 22.2.2015. 2575 ! 2576 !-- Rewritten so that all the coefficients and client-array indices are 2577 !-- precomputed in the initialization phase by pmci_init_interp_tril. 2578 ! 2579 ! Antti Hellsten 3.3.2015. 2580 ! 2581 !-- Constant dz no more assumed. 2582 ! Antti Hellsten 23.3.2015. 2583 ! 2584 !-- Adapted for non-flat topography. However, the near-wall velocities 2585 !-- are log-corrected only over horizontal surfaces, not yet near vertical 2586 !-- walls. 2587 !-- Antti Hellsten 26.3.2015. 2588 ! 2589 !-- Indexing in the principal direction (i) is changed. Now, the nest-boundary 2590 !-- values are interpolated only into the first ghost-node layers on each later 2591 !-- boundary. These values are then simply copied to the second ghost-node layer. 2592 ! 2593 !-- Antti Hellsten 6.10.2015. 3254 !-- conditions. This subroutine handles the left and right boundaries. It is 3255 !-- based on trilinear interpolation. Constant dz is still assumed. 3256 !-- TO_DO: constant dz is an important restriction and should be checked 3257 !-- somewhere in order to let users not run into this trap 2594 3258 IMPLICIT NONE 3259 3260 !-- TO_DO: wrap long lines in this and the remaining interp_tril routines 2595 3261 REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg), INTENT(INOUT) :: f !: 2596 3262 REAL(wp), DIMENSION(0:cg%nz+1,jcs:jcn,icl:icr), INTENT(IN) :: fc !: … … 2615 3281 INTEGER(iwp) :: i !: 2616 3282 INTEGER(iwp) :: ib !: 3283 INTEGER(iwp) :: ibgp !: 2617 3284 INTEGER(iwp) :: iw !: 2618 3285 INTEGER(iwp) :: j !: … … 2759 3426 2760 3427 ! 2761 !-- Store the boundary values also into the second ghost node layer. 2762 f(0:nzt+1,nys:nyn+1,ib) = f(0:nzt+1,nys:nyn+1,i) 3428 !-- Store the boundary values also into the other redundant ghost node layers 3429 IF ( edge == 'l' ) THEN 3430 DO ibgp = -nbgp, ib 3431 f(0:nzt+1,nysg:nyng,ibgp) = f(0:nzt+1,nysg:nyng,i) 3432 ENDDO 3433 ELSE IF ( edge == 'r' ) THEN 3434 DO ibgp = ib, nx + nbgp 3435 f(0:nzt+1,nysg:nyng,ibgp) = f(0:nzt+1,nysg:nyng,i) 3436 ENDDO 3437 ENDIF 2763 3438 2764 3439 END SUBROUTINE pmci_interp_tril_lr … … 2823 3498 INTEGER(iwp) :: j !: 2824 3499 INTEGER(iwp) :: jb !: 3500 INTEGER(iwp) :: jbgp !: 2825 3501 INTEGER(iwp) :: k !: 2826 3502 INTEGER(iwp) :: kcorr !: … … 2960 3636 2961 3637 ! 2962 !-- Store the boundary values also into the second ghost node layer. 2963 f(0:nzt+1,jb,nxl:nxr+1) = f(0:nzt+1,j,nxl:nxr+1) 3638 !-- Store the boundary values also into the other redundant ghost node layers 3639 IF ( edge == 's' ) THEN 3640 DO jbgp = -nbgp, jb 3641 f(0:nzt+1,jbgp,nxlg:nxrg) = f(0:nzt+1,j,nxlg:nxrg) 3642 ENDDO 3643 ELSE IF ( edge == 'n' ) THEN 3644 DO jbgp = jb, ny + nbgp 3645 f(0:nzt+1,jbgp,nxlg:nxrg) = f(0:nzt+1,j,nxlg:nxrg) 3646 ENDDO 3647 ENDIF 2964 3648 2965 3649 END SUBROUTINE pmci_interp_tril_sn … … 2974 3658 !-- This subroutine is based on trilinear interpolation. 2975 3659 !-- Constant dz is still assumed. 2976 !2977 !-- Antti Hellsten 23.2.2015.2978 !2979 !2980 !-- Rewritten so that all the coefficients and client-array indices are2981 !-- precomputed in the initialization phase by pmci_init_interp_tril.2982 !2983 !-- Antti Hellsten 3.3.2015.2984 !2985 !-- Constant dz no more assumed.2986 !-- Antti Hellsten 23.3.2015.2987 !2988 !-- Indexing in the principal direction (k) is changed. Now, the nest-boundary2989 !-- values are interpolated only into the first ghost-node layer. Actually there is2990 !-- only one ghost-node layer in the k-direction.2991 !2992 !-- Antti Hellsten 6.10.2015.2993 !2994 3660 IMPLICIT NONE 2995 3661 REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg), INTENT(INOUT) :: f !: … … 3069 3735 3070 3736 3071 SUBROUTINE pmci_extrap_ifoutflow_lr( f, kb, edge, var ) 3072 3073 ! 3074 !-- After the interpolation of ghost-node values for the client-domain boundary 3075 !-- conditions, this subroutine checks if there is a local outflow through the 3076 !-- boundary. In that case this subroutine overwrites the interpolated values 3077 !-- by values extrapolated from the domain. This subroutine handles the left and 3078 !-- right boundaries. 3079 !-- However, this operation is only needed in case of one-way coupling. 3080 ! 3081 !-- Antti Hellsten 9.3.2015. 3082 ! 3083 !-- Indexing in the principal direction (i) is changed. Now, the nest-boundary 3084 !-- values are interpolated only into the first ghost-node layers on each later 3085 !-- boundary. These values are then simply copied to the second ghost-node layer. 3086 ! 3087 !-- Antti Hellsten 6.10.2015. 3088 ! 3089 IMPLICIT NONE 3090 REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg), INTENT(INOUT) :: f !: 3091 INTEGER(iwp), DIMENSION(nysg:nyng,nxlg:nxrg), INTENT(IN) :: kb !: 3092 3093 CHARACTER(LEN=1),INTENT(IN) :: edge !: 3094 CHARACTER(LEN=1),INTENT(IN) :: var !: 3095 3096 INTEGER(iwp) :: i !: 3097 INTEGER(iwp) :: ib !: 3098 INTEGER(iwp) :: ied !: 3099 INTEGER(iwp) :: j !: 3100 INTEGER(iwp) :: k !: 3737 SUBROUTINE pmci_extrap_ifoutflow_lr( f, kb, edge, var ) 3738 ! 3739 !-- After the interpolation of ghost-node values for the client-domain 3740 !-- boundary conditions, this subroutine checks if there is a local outflow 3741 !-- through the boundary. In that case this subroutine overwrites the 3742 !-- interpolated values by values extrapolated from the domain. This 3743 !-- subroutine handles the left and right boundaries. However, this operation 3744 !-- is only needed in case of one-way coupling. 3745 IMPLICIT NONE 3746 3747 CHARACTER(LEN=1),INTENT(IN) :: edge !: 3748 CHARACTER(LEN=1),INTENT(IN) :: var !: 3749 3750 INTEGER(iwp) :: i !: 3751 INTEGER(iwp) :: ib !: 3752 INTEGER(iwp) :: ibgp !: 3753 INTEGER(iwp) :: ied !: 3754 INTEGER(iwp) :: j !: 3755 INTEGER(iwp) :: k !: 3101 3756 3102 REAL(wp) :: outnor !: 3103 REAL(wp) :: vdotnor !: 3104 3105 ! 3106 !-- Check which edge is to be handled: left or right. 3107 IF ( edge == 'l' ) THEN 3108 IF ( var == 'u' ) THEN 3109 i = nxl 3110 ib = nxl - 1 3111 ied = nxl + 1 3112 ELSE 3113 i = nxl - 1 3114 ib = nxl - 2 3115 ied = nxl 3116 ENDIF 3117 outnor = -1.0_wp 3118 ELSEIF ( edge == 'r' ) THEN 3119 i = nxr + 1 3120 ib = nxr + 2 3121 ied = nxr 3122 outnor = 1.0_wp 3123 ENDIF 3124 3125 DO j = nys, nyn + 1 3126 DO k = kb(j,i), nzt +1 3127 vdotnor = outnor * u(k,j,ied) 3128 IF ( vdotnor > 0.0_wp ) THEN ! Local outflow. 3129 f(k,j,i) = f(k,j,ied) 3130 ENDIF 3131 ENDDO 3132 IF ( (var == 'u' ) .OR. (var == 'v' ) .OR. (var == 'w') ) THEN 3133 f(kb(j,i),j,i) = 0.0_wp 3134 ENDIF 3135 ENDDO 3136 3137 ! 3138 !-- Store the updated boundary values also into the second ghost node layer. 3139 f(0:nzt,nys:nyn+1,ib) = f(0:nzt,nys:nyn+1,i) 3140 3141 END SUBROUTINE pmci_extrap_ifoutflow_lr 3142 3143 3144 3145 SUBROUTINE pmci_extrap_ifoutflow_sn( f, kb, edge, var ) 3146 ! 3147 !-- After the interpolation of ghost-node values for the client-domain boundary 3148 !-- conditions, this subroutine checks if there is a local outflow through the 3149 !-- boundary. In that case this subroutine overwrites the interpolated values 3150 !-- by values extrapolated from the domain. This subroutine handles the south and 3151 !-- north boundaries. 3152 ! 3153 !-- Antti Hellsten 9.3.2015. 3154 ! 3155 !-- Indexing in the principal direction (j) is changed. Now, the nest-boundary 3156 !-- values are interpolated only into the first ghost-node layers on each later 3157 !-- boundary. These values are then simply copied to the second ghost-node layer. 3158 ! 3159 !-- Antti Hellsten 6.10.2015. 3160 3161 IMPLICIT NONE 3162 REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg), INTENT(INOUT) :: f !: 3163 INTEGER(iwp), DIMENSION(nysg:nyng,nxlg:nxrg), INTENT(IN) :: kb !: 3164 CHARACTER(LEN=1), INTENT(IN) :: edge !: 3165 CHARACTER(LEN=1), INTENT(IN) :: var !: 3757 INTEGER(iwp), DIMENSION(nysg:nyng,nxlg:nxrg), INTENT(IN) :: kb !: 3758 3759 REAL(wp) :: outnor !: 3760 REAL(wp) :: vdotnor !: 3761 3762 REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg), INTENT(INOUT) :: f !: 3763 3764 ! 3765 !-- Check which edge is to be handled: left or right 3766 IF ( edge == 'l' ) THEN 3767 IF ( var == 'u' ) THEN 3768 i = nxl 3769 ib = nxl - 1 3770 ied = nxl + 1 3771 ELSE 3772 i = nxl - 1 3773 ib = nxl - 2 3774 ied = nxl 3775 ENDIF 3776 outnor = -1.0_wp 3777 ELSEIF ( edge == 'r' ) THEN 3778 i = nxr + 1 3779 ib = nxr + 2 3780 ied = nxr 3781 outnor = 1.0_wp 3782 ENDIF 3783 3784 DO j = nys, nyn + 1 3785 DO k = kb(j,i), nzt +1 3786 vdotnor = outnor * u(k,j,ied) 3787 ! 3788 !-- Local outflow 3789 IF ( vdotnor > 0.0_wp ) THEN 3790 f(k,j,i) = f(k,j,ied) 3791 ENDIF 3792 ENDDO 3793 IF ( (var == 'u' ) .OR. (var == 'v' ) .OR. (var == 'w') ) THEN 3794 f(kb(j,i),j,i) = 0.0_wp 3795 ENDIF 3796 ENDDO 3797 3798 ! 3799 !-- Store the boundary values also into the redundant ghost node layers. 3800 IF ( edge == 'l' ) THEN 3801 DO ibgp = -nbgp, ib 3802 f(0:nzt+1,nysg:nyng,ibgp) = f(0:nzt+1,nysg:nyng,i) 3803 ENDDO 3804 ELSEIF ( edge == 'r' ) THEN 3805 DO ibgp = ib, nx + nbgp 3806 f(0:nzt+1,nysg:nyng,ibgp) = f(0:nzt+1,nysg:nyng,i) 3807 ENDDO 3808 ENDIF 3809 3810 END SUBROUTINE pmci_extrap_ifoutflow_lr 3811 3812 3813 3814 SUBROUTINE pmci_extrap_ifoutflow_sn( f, kb, edge, var ) 3815 ! 3816 !-- After the interpolation of ghost-node values for the client-domain 3817 !-- boundary conditions, this subroutine checks if there is a local outflow 3818 !-- through the boundary. In that case this subroutine overwrites the 3819 !-- interpolated values by values extrapolated from the domain. This 3820 !-- subroutine handles the south and north boundaries. 3821 IMPLICIT NONE 3822 3823 CHARACTER(LEN=1), INTENT(IN) :: edge !: 3824 CHARACTER(LEN=1), INTENT(IN) :: var !: 3166 3825 3167 INTEGER(iwp) :: i !: 3168 INTEGER(iwp) :: j !: 3169 INTEGER(iwp) :: jb !: 3170 INTEGER(iwp) :: jed !: 3171 INTEGER(iwp) :: k !: 3172 REAL(wp) :: outnor !: 3173 REAL(wp) :: vdotnor !: 3174 3175 ! 3176 !-- Check which edge is to be handled: left or right. 3177 IF ( edge == 's' ) THEN 3178 IF ( var == 'v' ) THEN 3179 j = nys 3180 jb = nys - 1 3181 jed = nys + 1 3182 ELSE 3183 j = nys - 1 3184 jb = nys - 2 3185 jed = nys 3186 ENDIF 3187 outnor = -1.0_wp 3188 ELSEIF ( edge == 'n' ) THEN 3189 j = nyn + 1 3190 jb = nyn + 2 3191 jed = nyn 3192 outnor = 1.0_wp 3193 ENDIF 3194 3195 DO i = nxl, nxr + 1 3196 DO k = kb(j,i), nzt + 1 3197 vdotnor = outnor * v(k,jed,i) 3198 IF ( vdotnor > 0.0_wp ) THEN ! Local outflow. 3199 f(k,j,i) = f(k,jed,i) 3200 ENDIF 3201 ENDDO 3202 IF ( (var == 'u' ) .OR. (var == 'v' ) .OR. (var == 'w') ) THEN 3203 f(kb(j,i),j,i) = 0.0_wp 3204 ENDIF 3205 ENDDO 3206 3207 ! 3208 !-- Store the updated boundary values also into the second ghost node layer. 3209 f(0:nzt,jb,nxl:nxr+1) = f(0:nzt,j,nxl:nxr+1) 3210 3211 END SUBROUTINE pmci_extrap_ifoutflow_sn 3826 INTEGER(iwp) :: i !: 3827 INTEGER(iwp) :: j !: 3828 INTEGER(iwp) :: jb !: 3829 INTEGER(iwp) :: jbgp !: 3830 INTEGER(iwp) :: jed !: 3831 INTEGER(iwp) :: k !: 3832 3833 INTEGER(iwp), DIMENSION(nysg:nyng,nxlg:nxrg), INTENT(IN) :: kb !: 3834 3835 REAL(wp) :: outnor !: 3836 REAL(wp) :: vdotnor !: 3837 3838 REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg), INTENT(INOUT) :: f !: 3839 3840 ! 3841 !-- Check which edge is to be handled: left or right 3842 IF ( edge == 's' ) THEN 3843 IF ( var == 'v' ) THEN 3844 j = nys 3845 jb = nys - 1 3846 jed = nys + 1 3847 ELSE 3848 j = nys - 1 3849 jb = nys - 2 3850 jed = nys 3851 ENDIF 3852 outnor = -1.0_wp 3853 ELSEIF ( edge == 'n' ) THEN 3854 j = nyn + 1 3855 jb = nyn + 2 3856 jed = nyn 3857 outnor = 1.0_wp 3858 ENDIF 3859 3860 DO i = nxl, nxr + 1 3861 DO k = kb(j,i), nzt + 1 3862 vdotnor = outnor * v(k,jed,i) 3863 ! 3864 !-- Local outflow 3865 IF ( vdotnor > 0.0_wp ) THEN 3866 f(k,j,i) = f(k,jed,i) 3867 ENDIF 3868 ENDDO 3869 IF ( (var == 'u' ) .OR. (var == 'v' ) .OR. (var == 'w') ) THEN 3870 f(kb(j,i),j,i) = 0.0_wp 3871 ENDIF 3872 ENDDO 3873 3874 ! 3875 !-- Store the boundary values also into the redundant ghost node layers. 3876 IF ( edge == 's' ) THEN 3877 DO jbgp = -nbgp, jb 3878 f(0:nzt+1,jbgp,nxlg:nxrg) = f(0:nzt+1,j,nxlg:nxrg) 3879 ENDDO 3880 ELSEIF ( edge == 'n' ) THEN 3881 DO jbgp = jb, ny + nbgp 3882 f(0:nzt+1,jbgp,nxlg:nxrg) = f(0:nzt+1,j,nxlg:nxrg) 3883 ENDDO 3884 ENDIF 3885 3886 END SUBROUTINE pmci_extrap_ifoutflow_sn 3212 3887 3213 3888 3214 3889 3215 SUBROUTINE pmci_extrap_ifoutflow_t( f, var ) 3216 3217 ! 3218 !-- Interpolation of ghost-node values used as the client-domain boundary 3219 !-- conditions. This subroutine handles the top boundary. 3220 !-- This subroutine is based on trilinear interpolation. 3221 ! 3222 !-- Antti Hellsten 23.2.2015. 3223 ! 3224 ! 3225 !-- Rewritten so that all the coefficients and client-array indices are 3226 !-- precomputed in the initialization phase by init_interp_tril. 3227 ! 3228 !-- Antti Hellsten 3.3.2015. 3229 ! 3230 !-- Indexing in the principal direction (k) is changed. Now, the nest-boundary 3231 !-- values are extrapolated only into the first ghost-node layer. Actually there is 3232 !-- only one ghost-node layer in the k-direction. 3233 ! 3234 !-- Antti Hellsten 6.10.2015. 3235 IMPLICIT NONE 3236 REAL(wp), DIMENSION(nzb:nzt+1,nys-nbgp:nyn+nbgp,nxl-nbgp:nxr+nbgp), INTENT(INOUT) :: f !: 3237 CHARACTER(LEN=1), INTENT(IN) :: var !: 3890 SUBROUTINE pmci_extrap_ifoutflow_t( f, var ) 3891 ! 3892 !-- Interpolation of ghost-node values used as the client-domain boundary 3893 !-- conditions. This subroutine handles the top boundary. It is based on 3894 !-- trilinear interpolation. 3895 IMPLICIT NONE 3896 3897 CHARACTER(LEN=1), INTENT(IN) :: var !: 3238 3898 3239 INTEGER(iwp) :: i !: 3240 INTEGER(iwp) :: j !: 3241 INTEGER(iwp) :: k !: 3242 INTEGER(iwp) :: ked !: 3243 REAL(wp) :: vdotnor !: 3899 INTEGER(iwp) :: i !: 3900 INTEGER(iwp) :: j !: 3901 INTEGER(iwp) :: k !: 3902 INTEGER(iwp) :: ked !: 3903 3904 REAL(wp) :: vdotnor !: 3905 3906 REAL(wp), DIMENSION(nzb:nzt+1,nys-nbgp:nyn+nbgp,nxl-nbgp:nxr+nbgp), & 3907 INTENT(INOUT) :: f !: 3244 3908 3245 3909 3246 IF ( var == 'w' ) THEN 3247 k = nzt 3248 ked = nzt - 1 3249 ELSE 3250 k = nzt + 1 3251 ked = nzt 3252 ENDIF 3253 3254 DO i = nxl, nxr 3255 DO j = nys, nyn 3256 vdotnor = w(ked,j,i) 3257 IF ( vdotnor > 0.0_wp ) THEN !: Local outflow. 3258 f(k,j,i) = f(ked,j,i) 3259 ENDIF 3260 ENDDO 3261 ENDDO 3262 3263 ! 3264 !-- Just fill up the second ghost-node layer for w. 3265 IF ( var == 'w' ) THEN 3266 f(nzt+1,:,:) = f(nzt,:,:) 3267 ENDIF 3268 3269 END SUBROUTINE pmci_extrap_ifoutflow_t 3270 3271 3272 3273 SUBROUTINE pmci_anterp_tophat( f, fc, kce, ifl, ifu, jfl, jfu, kfl, kfu, kb, var ) 3274 ! 3275 !-- Anterpolation of internal-node values to be used as the server-domain 3276 !-- values. This subroutine is based on the first-order numerical 3277 !-- integration of the fine-grid values contained within the coarse-grid 3278 !-- cell. 3279 ! 3280 !-- Antti Hellsten 16.9.2015. 3281 ! 3282 IMPLICIT NONE 3283 REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg), INTENT(IN) :: f !: 3284 REAL(wp), DIMENSION(0:cg%nz+1,jcs:jcn,icl:icr), INTENT(INOUT) :: fc !: 3285 INTEGER(iwp), DIMENSION(icl:icr), INTENT(IN) :: ifl !: 3286 INTEGER(iwp), DIMENSION(icl:icr), INTENT(IN) :: ifu !: 3287 INTEGER(iwp), DIMENSION(jcs:jcn), INTENT(IN) :: jfl !: 3288 INTEGER(iwp), DIMENSION(jcs:jcn), INTENT(IN) :: jfu !: 3289 INTEGER(iwp), DIMENSION(nysg:nyng,nxlg:nxrg), INTENT(IN) :: kb !: may be unnecessary 3290 INTEGER(iwp), INTENT(IN) :: kce !: 3291 INTEGER(iwp), DIMENSION(0:kce), INTENT(IN) :: kfl !: 3292 INTEGER(iwp), DIMENSION(0:kce), INTENT(IN) :: kfu !: 3293 CHARACTER(LEN=1), INTENT(IN) :: var !: 3294 3295 INTEGER(iwp) :: i !: 3296 INTEGER(iwp) :: icb !: 3297 INTEGER(iwp) :: ice !: 3298 INTEGER(iwp) :: ifc !: 3299 INTEGER(iwp) :: ijfc !: 3300 INTEGER(iwp) :: j !: 3301 INTEGER(iwp) :: jcb !: 3302 INTEGER(iwp) :: jce !: 3303 INTEGER(iwp) :: k !: 3304 INTEGER(iwp) :: kcb !: 3305 INTEGER(iwp) :: l !: 3306 INTEGER(iwp) :: m !: 3307 INTEGER(iwp) :: n !: 3308 INTEGER(iwp) :: nfc !: 3309 REAL(wp) :: cellsum !: 3310 REAL(wp) :: f1f !: 3311 REAL(wp) :: fra !: 3910 IF ( var == 'w' ) THEN 3911 k = nzt 3912 ked = nzt - 1 3913 ELSE 3914 k = nzt + 1 3915 ked = nzt 3916 ENDIF 3917 3918 DO i = nxl, nxr 3919 DO j = nys, nyn 3920 vdotnor = w(ked,j,i) 3921 ! 3922 !-- Local outflow 3923 IF ( vdotnor > 0.0_wp ) THEN 3924 f(k,j,i) = f(ked,j,i) 3925 ENDIF 3926 ENDDO 3927 ENDDO 3928 3929 ! 3930 !-- Just fill up the second ghost-node layer for w 3931 IF ( var == 'w' ) THEN 3932 f(nzt+1,:,:) = f(nzt,:,:) 3933 ENDIF 3934 3935 END SUBROUTINE pmci_extrap_ifoutflow_t 3936 3937 3938 3939 SUBROUTINE pmci_anterp_tophat( f, fc, kct, ifl, ifu, jfl, jfu, kfl, kfu, & 3940 kb, var ) 3941 ! 3942 !-- Anterpolation of internal-node values to be used as the server-domain 3943 !-- values. This subroutine is based on the first-order numerical 3944 !-- integration of the fine-grid values contained within the coarse-grid 3945 !-- cell. 3946 IMPLICIT NONE 3947 3948 CHARACTER(LEN=1), INTENT(IN) :: var !: 3949 3950 INTEGER(iwp) :: i !: Fine-grid index 3951 INTEGER(iwp) :: ii !: Coarse-grid index 3952 INTEGER(iwp) :: iclp !: 3953 INTEGER(iwp) :: icrm !: 3954 INTEGER(iwp) :: ifc !: 3955 INTEGER(iwp) :: ijfc !: 3956 INTEGER(iwp) :: j !: Fine-grid index 3957 INTEGER(iwp) :: jj !: Coarse-grid index 3958 INTEGER(iwp) :: jcnm !: 3959 INTEGER(iwp) :: jcsp !: 3960 INTEGER(iwp) :: k !: Fine-grid index 3961 INTEGER(iwp) :: kk !: Coarse-grid index 3962 INTEGER(iwp) :: kcb !: 3963 INTEGER(iwp) :: nfc !: 3964 3965 INTEGER(iwp), INTENT(IN) :: kct !: 3966 3967 INTEGER(iwp), DIMENSION(icl:icr), INTENT(IN) :: ifl !: 3968 INTEGER(iwp), DIMENSION(icl:icr), INTENT(IN) :: ifu !: 3969 INTEGER(iwp), DIMENSION(jcs:jcn), INTENT(IN) :: jfl !: 3970 INTEGER(iwp), DIMENSION(jcs:jcn), INTENT(IN) :: jfu !: 3971 !-- TO_DO: is the next line really unnecessary? 3972 INTEGER(iwp), DIMENSION(nysg:nyng,nxlg:nxrg), INTENT(IN) :: kb !: may be unnecessary 3973 INTEGER(iwp), DIMENSION(0:kct), INTENT(IN) :: kfl !: 3974 INTEGER(iwp), DIMENSION(0:kct), INTENT(IN) :: kfu !: 3975 3976 3977 REAL(wp) :: cellsum !: 3978 REAL(wp) :: f1f !: 3979 REAL(wp) :: fra !: 3980 3981 REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg), INTENT(IN) :: f !: 3982 REAL(wp), DIMENSION(0:cg%nz+1,jcs:jcn,icl:icr), INTENT(INOUT) :: fc !: 3312 3983 3313 3984 3314 icb = icl 3315 ice = icr 3316 jcb = jcs 3317 jce = jcn 3318 3319 ! 3320 !-- Define the index bounds icb, ice, jcb and jce. 3321 !-- Note that kcb is simply zero and kce enters here as a parameter and it is 3322 !-- determined in init_anterp_tophat 3323 IF ( nest_bound_l ) THEN 3324 IF ( var == 'u' ) THEN 3325 icb = icl + nhll + 1 3326 ELSE 3327 icb = icl + nhll 3328 ENDIF 3329 ENDIF 3330 IF ( nest_bound_r ) THEN 3331 ice = icr - nhlr 3332 ENDIF 3333 3334 IF ( nest_bound_s ) THEN 3335 IF ( var == 'v' ) THEN 3336 jcb = jcs + nhls + 1 3337 ELSE 3338 jcb = jcs + nhls 3339 ENDIF 3340 ENDIF 3341 IF ( nest_bound_n ) THEN 3342 jce = jcn - nhln 3343 ENDIF 3344 kcb = 0 3345 3346 ! 3347 !-- Note that l,m, and n are coarse-grid indices and i,j, and k are fine-grid indices. 3348 DO l = icb, ice 3349 ifc = ifu(l) - ifl(l) + 1 3350 DO m = jcb, jce 3351 ijfc = ifc * ( jfu(m) - jfl(m) +1 ) 3352 3353 ! 3354 !-- How to deal with the lower bound of k in case of non-flat topography? 3355 !kcb = MIN( kb(jfl(m),ifl(l)), kb(jfu(m),ifl(l)), kb(jfl(m),ifu(l)), kb(jfu(m),ifu(l)) ) ! Something wrong with this. 3356 DO n = kcb, kce 3357 nfc = ijfc * ( kfu(n) - kfl(n) + 1 ) 3358 cellsum = 0.0 3359 DO i = ifl(l), ifu(l) 3360 DO j = jfl(m), jfu(m) 3361 DO k = kfl(n), kfu(n) 3362 cellsum = cellsum + f(k,j,i) 3363 ENDDO 3364 ENDDO 3365 ENDDO 3366 3367 ! 3368 !-- Spatial under-relaxation. 3369 fra = frax(l) * fray(m) * fraz(n) 3370 !-- TO_DO: why not KIND=wp ? 3371 fc(n,m,l) = ( 1.0_wp - fra ) * fc(n,m,l) + fra * cellsum / REAL( nfc, KIND=KIND(cellsum) ) 3372 ENDDO 3373 ENDDO 3374 ENDDO 3375 3376 END SUBROUTINE pmci_anterp_tophat 3985 ! 3986 !-- Initialize the index bounds for anterpolation 3987 iclp = icl 3988 icrm = icr 3989 jcsp = jcs 3990 jcnm = jcn 3991 3992 ! 3993 !-- Define the index bounds iclp, icrm, jcsp and jcnm. 3994 !-- Note that kcb is simply zero and kct enters here as a parameter and it is 3995 !-- determined in pmci_init_anterp_tophat 3996 IF ( nest_bound_l ) THEN 3997 IF ( var == 'u' ) THEN 3998 iclp = icl + nhll + 1 3999 ELSE 4000 iclp = icl + nhll 4001 ENDIF 4002 ENDIF 4003 IF ( nest_bound_r ) THEN 4004 icrm = icr - nhlr 4005 ENDIF 4006 4007 IF ( nest_bound_s ) THEN 4008 IF ( var == 'v' ) THEN 4009 jcsp = jcs + nhls + 1 4010 ELSE 4011 jcsp = jcs + nhls 4012 ENDIF 4013 ENDIF 4014 IF ( nest_bound_n ) THEN 4015 jcnm = jcn - nhln 4016 ENDIF 4017 kcb = 0 4018 4019 ! 4020 !-- Note that l,m, and n are coarse-grid indices and i,j, and k are fine-grid 4021 !-- indices. 4022 DO ii = iclp, icrm 4023 ifc = ifu(ii) - ifl(ii) + 1 4024 DO jj = jcsp, jcnm 4025 ijfc = ifc * ( jfu(jj) - jfl(jj) + 1 ) 4026 ! 4027 !-- For simplicity anterpolate within buildings too 4028 DO kk = kcb, kct 4029 nfc = ijfc * ( kfu(kk) - kfl(kk) + 1 ) 4030 cellsum = 0.0_wp 4031 DO i = ifl(ii), ifu(ii) 4032 DO j = jfl(jj), jfu(jj) 4033 DO k = kfl(kk), kfu(kk) 4034 cellsum = cellsum + f(k,j,i) 4035 ENDDO 4036 ENDDO 4037 ENDDO 4038 ! 4039 !-- Spatial under-relaxation. 4040 fra = frax(ii) * fray(jj) * fraz(kk) 4041 ! 4042 !-- Block out the fine-grid corner patches from the anterpolation 4043 IF ( ( ifl(ii) < nxl ) .OR. ( ifu(ii) > nxr ) ) THEN 4044 IF ( ( jfl(jj) < nys ) .OR. ( jfu(jj) > nyn ) ) THEN 4045 fra = 0.0_wp 4046 ENDIF 4047 ENDIF 4048 4049 fc(kk,jj,ii) = ( 1.0_wp - fra ) * fc(kk,jj,ii) + & 4050 fra * cellsum / REAL( nfc, KIND = wp ) 4051 4052 ENDDO 4053 ENDDO 4054 ENDDO 4055 4056 END SUBROUTINE pmci_anterp_tophat 3377 4057 3378 4058 #endif 3379 4059 END SUBROUTINE pmci_client_datatrans 3380 4060 3381 3382 3383 SUBROUTINE pmci_update_new3384 3385 #if defined( __parallel )3386 !3387 !-- Copy the interpolated/anterpolated boundary values to the _p3388 !-- arrays, too, to make sure the interpolated/anterpolated boundary3389 !-- values are carried over from one RK inner step to another.3390 !-- So far works only with the cpp-switch __nopointer.3391 !3392 !-- Antti Hellsten 8.3.20153393 !3394 3395 !-- Just debugging3396 w(nzt+1,:,:) = w(nzt,:,:)3397 3398 u_p = u3399 v_p = v3400 w_p = w3401 e_p = e3402 pt_p = pt3403 IF ( humidity .OR. passive_scalar ) THEN3404 q_p = q3405 ENDIF3406 3407 !3408 !-- TO_DO: Find out later if nesting would work without __nopointer.3409 #endif3410 3411 END SUBROUTINE pmci_update_new3412 3413 3414 3415 SUBROUTINE pmci_set_array_pointer( name, client_id, nz_cl )3416 3417 IMPLICIT NONE3418 3419 INTEGER, INTENT(IN) :: client_id !:3420 INTEGER, INTENT(IN) :: nz_cl !:3421 CHARACTER(LEN=*), INTENT(IN) :: name !:3422 3423 #if defined( __parallel )3424 REAL(wp), POINTER, DIMENSION(:,:) :: p_2d !:3425 REAL(wp), POINTER, DIMENSION(:,:) :: p_2d_sec !:3426 REAL(wp), POINTER, DIMENSION(:,:,:) :: p_3d !:3427 REAL(wp), POINTER, DIMENSION(:,:,:) :: p_3d_sec !:3428 INTEGER(iwp) :: ierr !:3429 INTEGER(iwp) :: istat !:3430 3431 3432 NULLIFY( p_3d )3433 NULLIFY( p_2d )3434 3435 !3436 !-- List of array names, which can be coupled3437 !-- In case of 3D please change also the second array for the pointer version3438 IF ( TRIM(name) == "u" ) p_3d => u3439 IF ( TRIM(name) == "v" ) p_3d => v3440 IF ( TRIM(name) == "w" ) p_3d => w3441 IF ( TRIM(name) == "e" ) p_3d => e3442 IF ( TRIM(name) == "pt" ) p_3d => pt3443 IF ( TRIM(name) == "q" ) p_3d => q3444 !3445 !-- This is just an example for a 2D array, not active for coupling3446 !-- Please note, that z0 has to be declared as TARGET array in modules.f903447 ! IF ( TRIM(name) == "z0" ) p_2d => z03448 3449 #if defined( __nopointer )3450 IF ( ASSOCIATED( p_3d ) ) THEN3451 CALL pmc_s_set_dataarray( client_id, p_3d, nz_cl, nz )3452 ELSEIF ( ASSOCIATED( p_2d ) ) THEN3453 CALL pmc_s_set_dataarray( client_id, p_2d )3454 ELSE3455 !3456 !-- Give only one message for the root domain3457 IF ( myid == 0 .AND. cpl_id == 1 ) THEN3458 3459 message_string = 'pointer for array "' // TRIM( name ) // &3460 '" can''t be associated'3461 CALL message( 'pmci_set_array_pointer', 'PA0117', 3, 2, 0, 6, 0 )3462 ELSE3463 !3464 !-- Avoid others to continue3465 CALL MPI_BARRIER( comm2d, ierr )3466 ENDIF3467 ENDIF3468 #else3469 !-- TO_DO: Why aren't the other pointers (p_3d) not set to u_1, v_1, etc.??3470 IF ( TRIM(name) == "u" ) p_3d_sec => u_23471 IF ( TRIM(name) == "v" ) p_3d_sec => v_23472 IF ( TRIM(name) == "w" ) p_3d_sec => w_23473 IF ( TRIM(name) == "e" ) p_3d_sec => e_23474 IF ( TRIM(name) == "pt" ) p_3d_sec => pt_23475 IF ( TRIM(name) == "q" ) p_3d_sec => q_23476 3477 IF ( ASSOCIATED( p_3d ) ) THEN3478 CALL pmc_s_set_dataarray( client_id, p_3d, nz_cl, nz, &3479 array_2 = p_3d_sec )3480 ELSEIF ( ASSOCIATED( p_2d ) ) THEN3481 CALL pmc_s_set_dataarray( client_id, p_2d )3482 ELSE3483 !3484 !-- Give only one message for the root domain3485 IF ( myid == 0 .AND. cpl_id == 1 ) THEN3486 3487 message_string = 'pointer for array "' // TRIM( name ) // &3488 '" can''t be associated'3489 CALL message( 'pmci_set_array_pointer', 'PA0117', 3, 2, 0, 6, 0 )3490 ELSE3491 !3492 !-- Avoid others to continue3493 CALL MPI_BARRIER( comm2d, ierr )3494 ENDIF3495 3496 ENDIF3497 #endif3498 3499 #endif3500 END SUBROUTINE pmci_set_array_pointer3501 3502 3503 3504 SUBROUTINE pmci_create_client_arrays( name, is, ie, js, je, nzc )3505 3506 IMPLICIT NONE3507 3508 INTEGER(iwp), INTENT(IN) :: ie !:3509 INTEGER(iwp), INTENT(IN) :: is !:3510 INTEGER(iwp), INTENT(IN) :: je !:3511 INTEGER(iwp), INTENT(IN) :: js !:3512 INTEGER(iwp), INTENT(IN) :: nzc !: Note that nzc is cg%nz3513 CHARACTER(LEN=*), INTENT(IN) :: name !:3514 3515 #if defined( __parallel )3516 REAL(wp), POINTER,DIMENSION(:,:) :: p_2d !:3517 REAL(wp), POINTER,DIMENSION(:,:,:) :: p_3d !:3518 INTEGER(iwp) :: ierr !:3519 INTEGER(iwp) :: istat !:3520 3521 3522 NULLIFY( p_3d )3523 NULLIFY( p_2d )3524 3525 !3526 !-- List of array names, which can be coupled.3527 !-- AH: Note that the k-range of the *c arrays is changed from 1:nz to 0:nz+1.3528 IF ( TRIM( name ) == "u" ) THEN3529 IF ( .NOT. ALLOCATED( uc ) ) ALLOCATE( uc(0:nzc+1, js:je, is:ie) )3530 p_3d => uc3531 ELSEIF ( TRIM( name ) == "v" ) THEN3532 IF ( .NOT. ALLOCATED( vc ) ) ALLOCATE( vc(0:nzc+1, js:je, is:ie) )3533 p_3d => vc3534 ELSEIF ( TRIM( name ) == "w" ) THEN3535 IF ( .NOT. ALLOCATED( wc ) ) ALLOCATE( wc(0:nzc+1, js:je, is:ie) )3536 p_3d => wc3537 ELSEIF ( TRIM( name ) == "e" ) THEN3538 IF ( .NOT. ALLOCATED( ec ) ) ALLOCATE( ec(0:nzc+1, js:je, is:ie) )3539 p_3d => ec3540 ELSEIF ( TRIM( name ) == "pt") THEN3541 IF ( .NOT. ALLOCATED( ptc ) ) ALLOCATE( ptc(0:nzc+1, js:je, is:ie) )3542 p_3d => ptc3543 ELSEIF ( TRIM( name ) == "q") THEN3544 IF ( .NOT. ALLOCATED( qc ) ) ALLOCATE( qc(0:nzc+1, js:je, is:ie) )3545 p_3d => qc3546 !ELSEIF (trim(name) == "z0") then3547 !IF (.not.allocated(z0c)) allocate(z0c(js:je, is:ie))3548 !p_2d => z0c3549 ENDIF3550 3551 IF ( ASSOCIATED( p_3d ) ) THEN3552 CALL pmc_c_set_dataarray( p_3d )3553 ELSEIF ( ASSOCIATED( p_2d ) ) THEN3554 CALL pmc_c_set_dataarray( p_2d )3555 ELSE3556 !3557 !-- Give only one message for the first client domain3558 IF ( myid == 0 .AND. cpl_id == 2 ) THEN3559 3560 message_string = 'pointer for array "' // TRIM( name ) // &3561 '" can''t be associated'3562 CALL message( 'pmci_create_client_arrays', 'PA0170', 3, 2, 0, 6, 0 )3563 ELSE3564 !3565 !-- Avoid others to continue3566 CALL MPI_BARRIER( comm2d, ierr )3567 ENDIF3568 ENDIF3569 3570 #endif3571 END SUBROUTINE pmci_create_client_arrays3572 3573 3574 3575 SUBROUTINE pmci_server_initialize3576 3577 #if defined( __parallel )3578 IMPLICIT NONE3579 3580 INTEGER(iwp) :: client_id !:3581 INTEGER(iwp) :: m !:3582 REAL(wp) :: waittime !:3583 3584 3585 DO m = 1, SIZE( pmc_server_for_client ) - 13586 client_id = pmc_server_for_client(m)3587 CALL pmc_s_fillbuffer( client_id, waittime=waittime )3588 ENDDO3589 3590 #endif3591 END SUBROUTINE pmci_server_initialize3592 3593 3594 3595 SUBROUTINE pmci_client_initialize3596 3597 #if defined( __parallel )3598 IMPLICIT NONE3599 3600 INTEGER(iwp) :: i !:3601 INTEGER(iwp) :: icl !:3602 INTEGER(iwp) :: icr !:3603 INTEGER(iwp) :: j !:3604 INTEGER(iwp) :: jcn !:3605 INTEGER(iwp) :: jcs !:3606 REAL(wp) :: waittime !:3607 3608 3609 IF ( cpl_id > 1 ) THEN ! Root id is never a client3610 3611 !3612 !-- Client domain boundaries in the server indice space.3613 icl = coarse_bound(1)3614 icr = coarse_bound(2)3615 jcs = coarse_bound(3)3616 jcn = coarse_bound(4)3617 3618 !3619 !-- Get data from server3620 CALL pmc_c_getbuffer( waittime = waittime )3621 3622 !3623 !-- The interpolation.3624 CALL pmci_interp_tril_all ( u, uc, icu, jco, kco, r1xu, r2xu, r1yo, r2yo, r1zo, r2zo, nzb_u_inner, 'u' )3625 CALL pmci_interp_tril_all ( v, vc, ico, jcv, kco, r1xo, r2xo, r1yv, r2yv, r1zo, r2zo, nzb_v_inner, 'v' )3626 CALL pmci_interp_tril_all ( w, wc, ico, jco, kcw, r1xo, r2xo, r1yo, r2yo, r1zw, r2zw, nzb_w_inner, 'w' )3627 CALL pmci_interp_tril_all ( e, ec, ico, jco, kco, r1xo, r2xo, r1yo, r2yo, r1zo, r2zo, nzb_s_inner, 'e' )3628 CALL pmci_interp_tril_all ( pt, ptc, ico, jco, kco, r1xo, r2xo, r1yo, r2yo, r1zo, r2zo, nzb_s_inner, 's' )3629 IF ( humidity .OR. passive_scalar ) THEN3630 CALL pmci_interp_tril_all ( q, qc, ico, jco, kco, r1xo, r2xo, r1yo, r2yo, r1zo, r2zo, nzb_s_inner, 's' )3631 ENDIF3632 3633 IF ( topography /= 'flat' ) THEN3634 3635 !3636 !-- Inside buildings set velocities and TKE back to zero.3637 !-- Other scalars (pt, q, s, km, kh, p, sa, ...) are ignored at present,3638 !-- maybe revise later.3639 DO i = nxlg, nxrg3640 DO j = nysg, nyng3641 u(nzb:nzb_u_inner(j,i),j,i) = 0.0_wp3642 v(nzb:nzb_v_inner(j,i),j,i) = 0.0_wp3643 w(nzb:nzb_w_inner(j,i),j,i) = 0.0_wp3644 e(nzb:nzb_s_inner(j,i),j,i) = 0.0_wp3645 u_p(nzb:nzb_u_inner(j,i),j,i) = 0.0_wp3646 v_p(nzb:nzb_v_inner(j,i),j,i) = 0.0_wp3647 w_p(nzb:nzb_w_inner(j,i),j,i) = 0.0_wp3648 e_p(nzb:nzb_s_inner(j,i),j,i) = 0.0_wp3649 ENDDO3650 ENDDO3651 ENDIF3652 ENDIF3653 3654 3655 CONTAINS3656 3657 3658 SUBROUTINE pmci_interp_tril_all( f, fc, ic, jc, kc, r1x, r2x, r1y, r2y, r1z, r2z, kb, var )3659 3660 !3661 !-- Interpolation of the internal values for the client-domain initialization.3662 !-- This subroutine is based on trilinear interpolation.3663 !-- Coding based on interp_tril_lr/sn/t3664 !3665 !-- Antti Hellsten 20.10.2015.3666 IMPLICIT NONE3667 REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg), INTENT(INOUT) :: f !:3668 REAL(wp), DIMENSION(0:cg%nz+1,jcs:jcn,icl:icr), INTENT(IN) :: fc !:3669 REAL(wp), DIMENSION(nxlg:nxrg), INTENT(IN) :: r1x !:3670 REAL(wp), DIMENSION(nxlg:nxrg), INTENT(IN) :: r2x !:3671 REAL(wp), DIMENSION(nysg:nyng), INTENT(IN) :: r1y !:3672 REAL(wp), DIMENSION(nysg:nyng), INTENT(IN) :: r2y !:3673 REAL(wp), DIMENSION(nzb:nzt+1), INTENT(IN) :: r1z !:3674 REAL(wp), DIMENSION(nzb:nzt+1), INTENT(IN) :: r2z !:3675 INTEGER(iwp), DIMENSION(nxlg:nxrg), INTENT(IN) :: ic !:3676 INTEGER(iwp), DIMENSION(nysg:nyng), INTENT(IN) :: jc !:3677 INTEGER(iwp), DIMENSION(nzb:nzt+1), INTENT(IN) :: kc !:3678 INTEGER(iwp), DIMENSION(nysg:nyng,nxlg:nxrg), INTENT(IN) :: kb !:3679 CHARACTER(LEN=1), INTENT(IN) :: var !:3680 3681 INTEGER(iwp) :: i !:3682 INTEGER(iwp) :: ib !:3683 INTEGER(iwp) :: ie !:3684 INTEGER(iwp) :: j !:3685 INTEGER(iwp) :: jb !:3686 INTEGER(iwp) :: je !:3687 INTEGER(iwp) :: k !:3688 INTEGER(iwp) :: k1 !:3689 INTEGER(iwp) :: kbc !:3690 INTEGER(iwp) :: l !:3691 INTEGER(iwp) :: m !:3692 INTEGER(iwp) :: n !:3693 REAL(wp) :: fk !:3694 REAL(wp) :: fkj !:3695 REAL(wp) :: fkjp !:3696 REAL(wp) :: fkp !:3697 REAL(wp) :: fkpj !:3698 REAL(wp) :: fkpjp !:3699 REAL(wp) :: logratio !:3700 REAL(wp) :: logzuc1 !:3701 REAL(wp) :: zuc1 !:3702 3703 3704 ib = nxl3705 ie = nxr3706 jb = nys3707 je = nyn3708 IF ( nest_bound_l ) THEN3709 ib = nxl - 13710 IF ( var == 'u' ) THEN ! For u, nxl is a ghost node, but not for the other variables.3711 ib = nxl3712 ENDIF3713 ENDIF3714 IF ( nest_bound_s ) THEN3715 jb = nys - 13716 IF ( var == 'v' ) THEN ! For v, nys is a ghost node, but not for the other variables.3717 jb = nys3718 ENDIF3719 ENDIF3720 IF ( nest_bound_r ) THEN3721 ie = nxr + 13722 ENDIF3723 IF ( nest_bound_n ) THEN3724 je = nyn + 13725 ENDIF3726 3727 !3728 !-- Trilinear interpolation.3729 DO i = ib, ie3730 DO j = jb, je3731 DO k = kb(j,i), nzt + 13732 l = ic(i)3733 m = jc(j)3734 n = kc(k)3735 fkj = r1x(i) * fc(n,m,l) + r2x(i) * fc(n,m,l+1)3736 fkjp = r1x(i) * fc(n,m+1,l) + r2x(i) * fc(n,m+1,l+1)3737 fkpj = r1x(i) * fc(n+1,m,l) + r2x(i) * fc(n+1,m,l+1)3738 fkpjp = r1x(i) * fc(n+1,m+1,l) + r2x(i) * fc(n+1,m+1,l+1)3739 fk = r1y(j) * fkj + r2y(j) * fkjp3740 fkp = r1y(j) * fkpj + r2y(j) * fkpjp3741 f(k,j,i) = r1z(k) * fk + r2z(k) * fkp3742 ENDDO3743 ENDDO3744 ENDDO3745 3746 !3747 !-- Correct the interpolated values of u and v in near-wall nodes, i.e. in3748 !-- the nodes below the coarse-grid nodes with k=1. The corrction is only made3749 !-- over horizontal wall surfaces in this phase. For the nest boundary conditions,3750 !-- a corresponding corrections is made for all vertical walls, too.3751 IF ( var == 'u' .OR. var == 'v' ) THEN3752 DO i = ib, nxr3753 DO j = jb, nyn3754 kbc = 13755 DO WHILE ( cg%zu(kbc) < zu(kb(j,i)) ) ! kbc is the first coarse-grid point above the surface.3756 kbc = kbc + 13757 ENDDO3758 zuc1 = cg%zu(kbc)3759 k1 = kb(j,i) + 13760 DO WHILE ( zu(k1) < zuc1 )3761 k1 = k1 + 13762 ENDDO3763 logzuc1 = LOG( ( zu(k1) - zu(kb(j,i)) ) / z0(j,i) )3764 3765 k = kb(j,i) + 13766 DO WHILE ( zu(k) < zuc1 )3767 logratio = ( LOG( ( zu(k) - zu(kb(j,i)) ) / z0(j,i)) ) / logzuc13768 f(k,j,i) = logratio * f(k1,j,i)3769 k = k + 13770 ENDDO3771 f(kb(j,i),j,i) = 0.0_wp3772 ENDDO3773 ENDDO3774 ELSEIF ( var == 'w' ) THEN3775 DO i = ib, nxr3776 DO j = jb, nyn3777 f(kb(j,i),j,i) = 0.0_wp3778 ENDDO3779 ENDDO3780 ENDIF3781 3782 END SUBROUTINE pmci_interp_tril_all3783 3784 #endif3785 END SUBROUTINE pmci_client_initialize3786 3787 3788 3789 SUBROUTINE pmci_ensure_nest_mass_conservation3790 3791 #if defined( __parallel )3792 !3793 !-- Adjust the volume-flow rate through the top boundary3794 !-- so that the net volume flow through all boundaries3795 !-- of the current nest domain becomes zero.3796 IMPLICIT NONE3797 3798 INTEGER(iwp) :: i !:3799 INTEGER(iwp) :: ierr !:3800 INTEGER(iwp) :: j !:3801 INTEGER(iwp) :: k !:3802 REAL(wp) :: dxdy !:3803 REAL(wp) :: innor !:3804 REAL(wp), DIMENSION(1:3) :: volume_flow_l !:3805 REAL(wp) :: w_lt !:3806 3807 !3808 !-- Sum up the volume flow through the left/right boundaries.3809 volume_flow(1) = 0.0_wp3810 volume_flow_l(1) = 0.0_wp3811 3812 IF ( nest_bound_l ) THEN3813 i = 03814 innor = dy3815 DO j = nys, nyn3816 DO k = nzb_u_inner(j,i) + 1, nzt3817 volume_flow_l(1) = volume_flow_l(1) + innor * u(k,j,i) * dzw(k)3818 ENDDO3819 ENDDO3820 ENDIF3821 3822 IF ( nest_bound_r ) THEN3823 i = nx + 13824 innor = -dy3825 DO j = nys, nyn3826 DO k = nzb_u_inner(j,i) + 1, nzt3827 volume_flow_l(1) = volume_flow_l(1) + innor * u(k,j,i) * dzw(k)3828 ENDDO3829 ENDDO3830 ENDIF3831 3832 #if defined( __parallel )3833 IF ( collective_wait ) CALL MPI_BARRIER( comm2d, ierr )3834 CALL MPI_ALLREDUCE( volume_flow_l(1), volume_flow(1), 1, MPI_REAL, &3835 MPI_SUM, comm2d, ierr )3836 #else3837 volume_flow(1) = volume_flow_l(1)3838 #endif3839 3840 !3841 !-- Sum up the volume flow through the south/north boundaries.3842 volume_flow(2) = 0.0_wp3843 volume_flow_l(2) = 0.0_wp3844 3845 IF ( nest_bound_s ) THEN3846 j = 03847 innor = dx3848 DO i = nxl, nxr3849 DO k = nzb_v_inner(j,i) + 1, nzt3850 volume_flow_l(2) = volume_flow_l(2) + innor * v(k,j,i) * dzw(k)3851 ENDDO3852 ENDDO3853 ENDIF3854 3855 IF ( nest_bound_n ) THEN3856 j = ny + 13857 innor = -dx3858 DO i = nxl, nxr3859 DO k = nzb_v_inner(j,i)+1, nzt3860 volume_flow_l(2) = volume_flow_l(2) + innor * v(k,j,i) * dzw(k)3861 ENDDO3862 ENDDO3863 ENDIF3864 3865 #if defined( __parallel )3866 IF ( collective_wait ) CALL MPI_BARRIER( comm2d, ierr )3867 CALL MPI_ALLREDUCE( volume_flow_l(2), volume_flow(2), 1, MPI_REAL, &3868 MPI_SUM, comm2d, ierr )3869 #else3870 volume_flow(2) = volume_flow_l(2)3871 #endif3872 3873 !3874 !-- Sum up the volume flow through the top boundary.3875 volume_flow(3) = 0.0_wp3876 volume_flow_l(3) = 0.0_wp3877 dxdy = dx * dy3878 k = nzt3879 DO i = nxl, nxr3880 DO j = nys, nyn3881 volume_flow_l(3) = volume_flow_l(3) - w(k,j,i) * dxdy3882 ENDDO3883 ENDDO3884 3885 #if defined( __parallel )3886 IF ( collective_wait ) CALL MPI_BARRIER( comm2d, ierr )3887 CALL MPI_ALLREDUCE( volume_flow_l(3), volume_flow(3), 1, MPI_REAL, &3888 MPI_SUM, comm2d, ierr )3889 #else3890 volume_flow(3) = volume_flow_l(3)3891 #endif3892 3893 !3894 !-- Correct the top-boundary value of w.3895 w_lt = (volume_flow(1) + volume_flow(2) + volume_flow(3)) / area_t3896 DO i = nxl, nxr3897 DO j = nys, nyn3898 DO k = nzt, nzt + 13899 w(k,j,i) = w(k,j,i) + w_lt3900 ENDDO3901 ENDDO3902 ENDDO3903 3904 #endif3905 END SUBROUTINE pmci_ensure_nest_mass_conservation3906 3907 4061 END MODULE pmc_interface -
palm/trunk/SOURCE/pmc_server.f90
r1787 r1791 20 20 ! Current revisions: 21 21 ! ------------------ 22 ! 22 ! Debug write-statements commented out 23 23 ! 24 24 ! Former revisions: … … 136 136 137 137 do i=1,size(PMC_Server_for_Client)-1 138 if(m_model_comm == 0) write(0,*) 'PMC_Server: Initialize client Id',PMC_Server_for_Client(i)138 ! if(m_model_comm == 0) write(0,*) 'PMC_Server: Initialize client Id',PMC_Server_for_Client(i) 139 139 140 140 ClientId = PMC_Server_for_Client(i) … … 152 152 CALL MPI_Comm_rank (Clients(ClientId)%intra_comm, Clients(ClientId)%intra_rank, istat); 153 153 154 write(9,*) 'ClientId ',i,ClientId,m_world_rank, Clients(ClientId)%inter_npes154 ! write(9,*) 'ClientId ',i,ClientId,m_world_rank, Clients(ClientId)%inter_npes 155 155 156 156 ALLOCATE (Clients(ClientId)%PEs(Clients(ClientId)%inter_npes)) -
palm/trunk/SOURCE/time_integration.f90
r1787 r1791 19 19 ! Current revisions: 20 20 ! ------------------ 21 ! 21 ! call of pmci_update_new removed 22 22 ! 23 23 ! Former revisions: … … 257 257 pmci_ensure_nest_mass_conservation, pmci_client_datatrans, & 258 258 pmci_client_synchronize, pmci_server_datatrans, & 259 pmci_server_synchronize, pmci_update_new,server_to_client259 pmci_server_synchronize, server_to_client 260 260 261 261 USE production_e_mod, & … … 703 703 IF ( nest_domain ) THEN 704 704 CALL pmci_ensure_nest_mass_conservation 705 !706 !-- pmc_update_new is not necessary if nesting is made at each707 !-- substep708 CALL pmci_update_new709 705 ENDIF 710 706
Note: See TracChangeset
for help on using the changeset viewer.
