Changeset 1797 for palm/trunk/SOURCE


Ignore:
Timestamp:
Mar 21, 2016 4:50:28 PM (9 years ago)
Author:
raasch
Message:

Introduction of different data transfer modes; restart mechanism adjusted for nested runs; parameter consistency checks for nested runs; further formatting cleanup

Location:
palm/trunk/SOURCE
Files:
8 edited

Legend:

Unmodified
Added
Removed
  • palm/trunk/SOURCE/Makefile

    r1792 r1797  
    2020# Current revisions:
    2121# ------------------
    22 #
     22# dependcy for check_for_restart updated
    2323#
    2424# Former revisions:
     
    342342calc_precipitation.o: modules.o mod_kinds.o
    343343calc_radiation.o: modules.o mod_kinds.o
    344 check_for_restart.o: modules.o mod_kinds.o
     344check_for_restart.o: modules.o mod_kinds.o pmc_interface.o
    345345check_open.o: modules.o mod_kinds.o mod_particle_attributes.o netcdf_interface.o
    346346check_parameters.o: modules.o mod_kinds.o land_surface_model.o \
  • palm/trunk/SOURCE/check_for_restart.f90

    r1683 r1797  
    1919! Current revisions:
    2020! -----------------
    21 !
     21! check now accounts for nesting mode
    2222!
    2323! Former revisions:
     
    7070               termination_time_needed, time_restart,                          &
    7171               time_since_reference_point, write_binary
     72
    7273    USE kinds
     74
    7375    USE pegrid
    7476
     77    USE pmc_interface,                                                         &
     78        ONLY:  comm_world_nesting, cpl_id, nested_run
     79
    7580    IMPLICIT NONE
    7681
    77 
    78     LOGICAL :: terminate_run_l           !<
    79     LOGICAL :: do_stop_now = .FALSE.     !<
    80     LOGICAL :: do_restart_now = .FALSE.  !<
     82    INTEGER ::  global_communicator       !< global communicator to be used here
     83
     84    LOGICAL ::  terminate_run_l           !<
     85    LOGICAL ::  do_stop_now = .FALSE.     !<
     86    LOGICAL ::  do_restart_now = .FALSE.  !<
    8187
    8288    REAL(wp) ::  remaining_time !<
     
    96102    ENDIF
    97103
     104!
     105!-- Set the global communicator to be used (depends on the mode in which PALM is
     106!-- running)
     107    IF ( nested_run )  THEN
     108       global_communicator = comm_world_nesting
     109    ELSE
     110       global_communicator = comm2d
     111    ENDIF
     112
    98113#if defined( __parallel )
    99114!
    100115!-- Make a logical OR for all processes. Stop the model run if at least
    101 !-- one processor has reached the time limit.
    102     IF ( collective_wait )  CALL MPI_BARRIER( comm2d, ierr )
    103     CALL MPI_ALLREDUCE( terminate_run_l, terminate_run, 1, MPI_LOGICAL,        &
    104                         MPI_LOR, comm2d, ierr )
     116!-- one process has reached the time limit.
     117    IF ( collective_wait )  CALL MPI_BARRIER( global_communicator, ierr )
     118    CALL MPI_ALLREDUCE( terminate_run_l, terminate_run, 1, MPI_LOGICAL,     &
     119                        MPI_LOR, global_communicator, ierr )
    105120#else
    106121    terminate_run = terminate_run_l
     
    164179#if defined( __parallel )
    165180!
    166 !-- Make a logical OR for all processes. Stop the model run if at least
    167 !-- one processor has reached the time limit.
    168     IF ( collective_wait )  CALL MPI_BARRIER( comm2d, ierr )
     181!-- Make a logical OR for all processes. Stop the model run if a flag file has
     182!-- been detected above.
     183    IF ( collective_wait )  CALL MPI_BARRIER( global_communicator, ierr )
    169184    CALL MPI_ALLREDUCE( terminate_run_l, terminate_run, 1, MPI_LOGICAL,        &
    170                         MPI_LOR, comm2d, ierr )
     185                        MPI_LOR, global_communicator, ierr )
    171186#else
    172187    terminate_run = terminate_run_l
     
    256271!-- the start of a continuation run, except if the user forced to stop the
    257272!-- run without restart
    258     IF ( terminate_run  .AND.  myid == 0 .AND. .NOT. do_stop_now)  THEN
     273    IF ( terminate_run  .AND.  myid == 0  .AND.  cpl_id == 1  .AND.            &
     274         .NOT. do_stop_now)  THEN
    259275
    260276       OPEN ( 90, FILE='CONTINUE_RUN', FORM='FORMATTED' )
  • palm/trunk/SOURCE/header.f90

    r1792 r1797  
    1919! Current revisions:
    2020! -----------------
    21 !
     21! output of nesting datatransfer mode
    2222!
    2323! Former revisions:
     
    291291
    292292    USE pmc_interface,                                                         &
    293         ONLY:  nested_run, nesting_mode
     293        ONLY:  nested_run, nesting_datatransfer_mode, nesting_mode
    294294
    295295    USE radiation_model_mod,                                                   &
     
    483483    IF ( nested_run )  THEN
    484484
    485        WRITE ( io, 600 )  TRIM( nesting_mode )
     485       WRITE ( io, 600 )  TRIM( nesting_mode ),                                &
     486                          TRIM( nesting_datatransfer_mode )
    486487       CALL pmc_get_model_info( ncpl = ncpl, cpl_id = my_cpl_id )
    487488
     
    24762477600 FORMAT (/' Nesting informations:'/ &
    24772478            ' --------------------'/ &
    2478             ' Nesting mode:                     ',A// &
     2479            ' Nesting mode:                     ',A/ &
     2480            ' Nesting-datatransfer mode:        ',A// &
    24792481            ' Nest id  parent  number   lower left coordinates   name'/ &
    24802482            ' (*=me)     id    of PEs      x (m)     y (m)' )
  • palm/trunk/SOURCE/pmc_client.f90

    r1792 r1797  
    2020! Current revisions:
    2121! ------------------
    22 !
     22! introduction of different datatransfer modes
    2323!
    2424! Former revisions:
     
    401401       do i=1,me%inter_npes
    402402          aPE => me%PEs(i)
    403           ar  => aPE%array_list(next_array_in_list)    !actual array is last array in list
     403          ar  => aPE%array_list(next_array_in_list)
    404404          ar%NrDims    = NrDims
    405405          ar%A_dim     = dims
     
    575575      INTEGER(kind=MPI_ADDRESS_KIND)    ::  target_disp
    576576
    577       t1 = PMC_Time()
    578       CALL MPI_Barrier(me%intra_comm, ierr)                         ! Wait for server to fill buffer
    579       t2 = PMC_Time()-t1
    580       if(present(WaitTime)) WaitTime = t2
     577!
     578!--   Synchronization of the model is done in pmci_client_synchronize and pmci_server_synchronize
     579!--   Therefor the RMA window can be filled without sychronization at this point and a barrier
     580!--   is not necessary
     581!--   Please note that WaitTime has to be set in PMC_S_FillBuffer AND PMC_C_GetBuffer
     582      if(present(WaitTime))  then
     583         t1 = PMC_Time()
     584         CALL MPI_Barrier(me%intra_comm, ierr)
     585         t2 = PMC_Time()
     586         WaitTime = t2-t1
     587      end if
    581588
    582589      CALL MPI_Barrier(me%intra_comm, ierr)                         ! Wait for buffer is filled
  • palm/trunk/SOURCE/pmc_handle_communicator.f90

    r1792 r1797  
    2020! Current revisions:
    2121! ------------------
    22 !
     22! introduction of different datatransfer modes,
     23! export of comm_world_nesting
    2324!
    2425! Former revisions:
     
    9091   ! Coupler Setup
    9192
     93   INTEGER                                    :: m_world_comm !global nesting communicator
    9294   INTEGER                                    :: m_my_CPL_id  !Coupler id of this model
    9395   INTEGER                                    :: m_Parent_id  !Coupler id of parent of this model
     
    125127 CONTAINS
    126128
    127    SUBROUTINE pmc_init_model( comm, nesting_mode, pmc_status )
     129   SUBROUTINE pmc_init_model( comm, nesting_datatransfer_mode, nesting_mode,   &
     130                              pmc_status )
    128131
    129132      USE control_parameters,                                                  &
     
    136139
    137140      CHARACTER(LEN=7), INTENT(OUT) ::  nesting_mode
     141      CHARACTER(LEN=7), INTENT(OUT) ::  nesting_datatransfer_mode
    138142
    139143      INTEGER, INTENT(OUT)                ::  comm
     
    148152      pmc_status   = pmc_status_ok
    149153      comm         = -1
     154      m_world_comm = MPI_COMM_WORLD
    150155      m_my_cpl_id  = -1
    151156      clientcount  =  0
     
    159164      IF ( m_world_rank == 0 )  THEN
    160165
    161          CALL read_coupling_layout( nesting_mode, pmc_status )
     166         CALL read_coupling_layout( nesting_datatransfer_mode, nesting_mode,   &
     167                                    pmc_status )
    162168
    163169         IF ( pmc_status /= pmc_no_namelist_found  .AND.                       &
     
    221227         CALL MPI_BCAST( m_couplers(i)%lower_left_x, 1, MPI_REAL,    0, MPI_COMM_WORLD, istat )
    222228         CALL MPI_BCAST( m_couplers(i)%lower_left_y, 1, MPI_REAL,    0, MPI_COMM_WORLD, istat )
     229!--      TO_DO: the next two calls can and should probably moved outside this loop
    223230         CALL MPI_BCAST( nesting_mode, LEN( nesting_mode ), MPI_CHARACTER, 0, MPI_COMM_WORLD, istat )
     231         CALL MPI_BCAST( nesting_datatransfer_mode, LEN(nesting_datatransfer_mode), MPI_CHARACTER, 0, MPI_COMM_WORLD, istat )
    224232      ENDDO
    225233
     
    329337!
    330338!-- Provide module private variables of the pmc for PALM
    331     SUBROUTINE pmc_get_model_info( cpl_id, cpl_name, cpl_parent_id,            &
    332                                    lower_left_x, lower_left_y, ncpl, npe_total,&
    333                                    request_for_cpl_id )
     339    SUBROUTINE pmc_get_model_info( comm_world_nesting, cpl_id, cpl_name,       &
     340                                   cpl_parent_id, lower_left_x, lower_left_y,  &
     341                                   ncpl, npe_total, request_for_cpl_id )
    334342
    335343      USE kinds
     
    341349      INTEGER, INTENT(IN), OPTIONAL ::  request_for_cpl_id
    342350
     351      INTEGER, INTENT(OUT), OPTIONAL ::  comm_world_nesting
    343352      INTEGER, INTENT(OUT), OPTIONAL ::  cpl_id
    344353      INTEGER, INTENT(OUT), OPTIONAL ::  cpl_parent_id
     
    364373!
    365374!--   Return the requested information
     375      IF ( PRESENT( comm_world_nesting )  )  THEN
     376         comm_world_nesting = m_world_comm
     377      ENDIF
    366378      IF ( PRESENT( cpl_id )        )  THEN
    367          cpl_id        = requested_cpl_id
     379         cpl_id = requested_cpl_id
    368380      ENDIF
    369381      IF ( PRESENT( cpl_parent_id ) )  THEN
     
    371383      ENDIF
    372384      IF ( PRESENT( cpl_name )      )  THEN
    373          cpl_name      = m_couplers(requested_cpl_id)%name
     385         cpl_name = m_couplers(requested_cpl_id)%name
    374386      ENDIF
    375387      IF ( PRESENT( ncpl )          )  THEN
    376          ncpl          = m_ncpl
     388         ncpl = m_ncpl
    377389      ENDIF
    378390      IF ( PRESENT( npe_total )     )  THEN
    379          npe_total     = m_couplers(requested_cpl_id)%npe_total
     391         npe_total = m_couplers(requested_cpl_id)%npe_total
    380392      ENDIF
    381393      IF ( PRESENT( lower_left_x )  )  THEN
    382          lower_left_x  = m_couplers(requested_cpl_id)%lower_left_x
     394         lower_left_x = m_couplers(requested_cpl_id)%lower_left_x
    383395      ENDIF
    384396      IF ( PRESENT( lower_left_y )  )  THEN
    385          lower_left_y  = m_couplers(requested_cpl_id)%lower_left_y
     397         lower_left_y = m_couplers(requested_cpl_id)%lower_left_y
    386398      ENDIF
    387399
     
    400412
    401413
    402  SUBROUTINE read_coupling_layout( nesting_mode, pmc_status )
     414 SUBROUTINE read_coupling_layout( nesting_datatransfer_mode, nesting_mode,     &
     415                                  pmc_status )
    403416
    404417    IMPLICIT NONE
    405418
    406     CHARACTER(LEN=7) ::  nesting_mode
     419    CHARACTER(LEN=7), INTENT(INOUT) ::  nesting_mode
     420    CHARACTER(LEN=7), INTENT(INOUT) ::  nesting_datatransfer_mode
    407421
    408422    INTEGER, INTENT(INOUT) ::  pmc_status
     
    412426
    413427
    414     NAMELIST /nestpar/  domain_layouts, nesting_mode
     428    NAMELIST /nestpar/  domain_layouts, nesting_datatransfer_mode, nesting_mode
    415429
    416430!
  • palm/trunk/SOURCE/pmc_interface.f90

    r1792 r1797  
    2020! Current revisions:
    2121! ------------------
    22 !
     22! introduction of different datatransfer modes,
     23! further formatting cleanup, parameter checks added (including mismatches
     24! between root and client model settings),
     25! +routine pmci_check_setting_mismatches
     26! comm_world_nesting introduced
    2327!
    2428! Former revisions:
     
    145149!
    146150!-- Coupler setup
     151    INTEGER(iwp), SAVE      ::  comm_world_nesting     !:
    147152    INTEGER(iwp), SAVE      ::  cpl_id  = 1            !:
    148153    CHARACTER(LEN=32), SAVE ::  cpl_name               !:
     
    152157!
    153158!-- Control parameters, will be made input parameters later
     159    CHARACTER(LEN=7), SAVE ::  nesting_datatransfer_mode = 'mixed'  !: steering
     160                                                         !: parameter for data-
     161                                                         !: transfer mode
    154162    CHARACTER(LEN=7), SAVE ::  nesting_mode = 'two-way'  !: steering parameter
    155163                                                         !: for 1- or 2-way nesting
     
    308316
    309317
    310     INTERFACE pmci_client_datatrans
    311        MODULE PROCEDURE pmci_client_datatrans
     318    INTERFACE pmci_check_setting_mismatches
     319       MODULE PROCEDURE pmci_check_setting_mismatches
    312320    END INTERFACE
    313321
     
    320328    END INTERFACE
    321329
     330    INTERFACE pmci_datatrans
     331       MODULE PROCEDURE pmci_datatrans
     332    END INTERFACE pmci_datatrans
     333
    322334    INTERFACE pmci_ensure_nest_mass_conservation
    323335       MODULE PROCEDURE pmci_ensure_nest_mass_conservation
     
    346358    PUBLIC anterp_relax_length_l, anterp_relax_length_r,                       &
    347359           anterp_relax_length_s, anterp_relax_length_n,                       &
    348            anterp_relax_length_t, client_to_server, cpl_id, nested_run,        &
    349            nesting_mode, server_to_client
    350     PUBLIC pmci_client_datatrans
     360           anterp_relax_length_t, client_to_server, comm_world_nesting,        &
     361           cpl_id, nested_run, nesting_datatransfer_mode, nesting_mode,        &
     362           server_to_client
    351363    PUBLIC pmci_client_initialize
    352364    PUBLIC pmci_client_synchronize
     365    PUBLIC pmci_datatrans
    353366    PUBLIC pmci_ensure_nest_mass_conservation
    354367    PUBLIC pmci_init
    355368    PUBLIC pmci_modelconfiguration
    356     PUBLIC pmci_server_datatrans
    357369    PUBLIC pmci_server_initialize
    358370    PUBLIC pmci_server_synchronize
     
    365377 SUBROUTINE pmci_init( world_comm )
    366378
     379    USE control_parameters,                                                  &
     380        ONLY:  message_string
     381
    367382    IMPLICIT NONE
    368383
     
    376391
    377392
    378     CALL pmc_init_model( world_comm, nesting_mode, pmc_status )
     393    CALL pmc_init_model( world_comm, nesting_datatransfer_mode, nesting_mode,  &
     394                         pmc_status )
    379395
    380396    IF ( pmc_status == pmc_no_namelist_found )  THEN
     
    390406
    391407!
     408!-- Check steering parameter values
     409    IF ( TRIM( nesting_mode ) /= 'one-way'  .AND.                              &
     410         TRIM( nesting_mode ) /= 'two-way' )                                   &
     411    THEN
     412       message_string = 'illegal nesting mode: ' // TRIM( nesting_mode )
     413       CALL message( 'pmci_init', 'PA0417', 3, 2, 0, 6, 0 )
     414    ENDIF
     415
     416    IF ( TRIM( nesting_datatransfer_mode ) /= 'cascade'  .AND.                 &
     417         TRIM( nesting_datatransfer_mode ) /= 'mixed'    .AND.                 &
     418         TRIM( nesting_datatransfer_mode ) /= 'overlap' )                      &
     419    THEN
     420       message_string = 'illegal nesting datatransfer mode: '                  &
     421                        // TRIM( nesting_datatransfer_mode )
     422       CALL message( 'pmci_init', 'PA0418', 3, 2, 0, 6, 0 )
     423    ENDIF
     424
     425!
    392426!-- Set the general steering switch which tells PALM that its a nested run
    393427    nested_run = .TRUE.
     
    396430!-- Get some variables required by the pmc-interface (and in some cases in the
    397431!-- PALM code out of the pmci) out of the pmc-core
    398     CALL pmc_get_model_info( cpl_id = cpl_id, cpl_parent_id = cpl_parent_id,   &
     432    CALL pmc_get_model_info( comm_world_nesting = comm_world_nesting,          &
     433                             cpl_id = cpl_id, cpl_parent_id = cpl_parent_id,   &
    399434                             cpl_name = cpl_name, npe_total = cpl_npe_total,   &
    400435                             lower_left_x = lower_left_coord_x,                &
     
    445480!-- Initialize PMC Server
    446481    CALL pmci_setup_server
     482!
     483!-- Check for mismatches between settings of master and client variables
     484!-- (e.g., all clients have to follow the end_time settings of the root master)
     485    CALL pmci_check_setting_mismatches
    447486
    448487    CALL location_message( 'finished', .TRUE. )
     
    701740       ELSE
    702741!
    703 !--       TO_DO: Klaus: comment why thie dummy allocation is required
     742!--       TO_DO: Klaus: comment why this dummy allocation is required
    704743          ALLOCATE( index_list(6,1) )
    705744          CALL pmc_s_set_2d_index_list( client_id, index_list )
     
    789828                                     SIZE(define_coarse_grid_real), 0, 21, ierr )
    790829          CALL pmc_recv_from_server( define_coarse_grid_int,  3, 0, 22, ierr )
    791 
    792 !
    793 !--       Receive also the dz-,zu- and zw-arrays here.
    794 !--       TO_DO: what is the meaning of above comment
    795830!
    796831!--        Debug-printouts - keep them
     
    852887       
    853888!
    854 !--    TO_DO: give comments what is happening here
     889!--    Find the index bounds for the nest domain in the coarse-grid index space
    855890       CALL pmci_map_fine_to_coarse_grid
     891!
     892!--    TO_DO: Klaus give a comment what is happening here
    856893       CALL pmc_c_get_2d_index_list
    857894
     
    879916
    880917!
    881 !--    Two-way coupling
    882 !--    TO_DO: comment what is happening here
     918!--    Two-way coupling.
     919!--    Precompute the index arrays and relaxation functions for the
     920!--    anterpolation
    883921       IF ( nesting_mode == 'two-way' )  THEN
    884922          CALL pmci_init_anterp_tophat
     
    895933
    896934    SUBROUTINE pmci_map_fine_to_coarse_grid
    897 
     935!
     936!--    Determine index bounds of interpolation/anterpolation area in the coarse
     937!--    grid index space
    898938       IMPLICIT NONE
    899939
     
    907947
    908948!
    909 !--    Determine indices of interpolation/anterpolation area in the coarse grid
    910 !--    If the fine- and coarse grid nodes do not match.
     949!--    If the fine- and coarse grid nodes do not match:
    911950       loffset = MOD( coord_x(nxl), cg%dx )
    912951       xexl    = cg%dx + loffset
     
    9911030!--    Precomputation of the interpolation coefficients and client-array indices
    9921031!--    to be used by the interpolation routines interp_tril_lr, interp_tril_ns
    993 !--    and interp_tril_t. Constant dz is still assumed.
     1032!--    and interp_tril_t.
    9941033
    9951034       IMPLICIT NONE
     
    11991238!--    Left boundary
    12001239       IF ( nest_bound_l )  THEN
     1240
    12011241          ALLOCATE( logc_u_l(nzb:nzt_topo_nestbc_l,nys:nyn,1:2) )
    12021242          ALLOCATE( logc_v_l(nzb:nzt_topo_nestbc_l,nys:nyn,1:2) )
     
    12351275
    12361276          ENDDO
     1277
    12371278       ENDIF
    12381279
     
    12401281!--    Right boundary
    12411282       IF ( nest_bound_r )  THEN
     1283
    12421284          ALLOCATE( logc_u_r(nzb:nzt_topo_nestbc_r,nys:nyn,1:2) )
    12431285          ALLOCATE( logc_v_r(nzb:nzt_topo_nestbc_r,nys:nyn,1:2) )
     
    12521294          DO  j = nys, nyn
    12531295!
    1254 !--          Right boundary for u.
     1296!--          Right boundary for u
    12551297             i   = nxr + 1
    12561298             kb  = nzb_u_inner(j,i)
     
    12621304             logc_ratio_u_r(k,j,1,0:ncorr-1) = lcr(0:ncorr-1)
    12631305             lcr(0:ncorr-1) = 1.0_wp
    1264 
    1265 !
    1266 !--          Right boundary for v.
     1306!
     1307!--          Right boundary for v
    12671308             i   = nxr + 1
    12681309             kb  = nzb_v_inner(j,i)
     
    12761317
    12771318          ENDDO
     1319
    12781320       ENDIF
    12791321
     
    12811323!--    South boundary
    12821324       IF ( nest_bound_s )  THEN
     1325
    12831326          ALLOCATE( logc_u_s(nzb:nzt_topo_nestbc_s,nxl:nxr,1:2) )
    12841327          ALLOCATE( logc_v_s(nzb:nzt_topo_nestbc_s,nxl:nxr,1:2) )
     
    12911334          direction      = 1
    12921335          inc            = 1
     1336
    12931337          DO  i = nxl, nxr
    12941338!
    1295 !--          South boundary for u.
     1339!--          South boundary for u
    12961340             j   = -1
    12971341             kb  =  nzb_u_inner(j,i)
     
    13031347             logc_ratio_u_s(k,i,1,0:ncorr-1) = lcr(0:ncorr-1)
    13041348             lcr(0:ncorr-1) = 1.0_wp
    1305 
    13061349!
    13071350!--          South boundary for v
     
    13151358             logc_ratio_v_s(k,i,1,0:ncorr-1) = lcr(0:ncorr-1)
    13161359             lcr(0:ncorr-1) = 1.0_wp
    1317           ENDDO
     1360
     1361          ENDDO
     1362
    13181363       ENDIF
    13191364
     
    13211366!--    North boundary
    13221367       IF ( nest_bound_n )  THEN
     1368
    13231369          ALLOCATE( logc_u_n(nzb:nzt_topo_nestbc_n,nxl:nxr,1:2) )
    13241370          ALLOCATE( logc_v_n(nzb:nzt_topo_nestbc_n,nxl:nxr,1:2) )
     
    13311377          direction      = 1
    13321378          inc            = 1
     1379
    13331380          DO  i = nxl, nxr
    13341381!
    1335 !--          North boundary for u.
     1382!--          North boundary for u
    13361383             j   = nyn + 1
    13371384             kb  = nzb_u_inner(j,i)
     
    13431390             logc_ratio_u_n(k,i,1,0:ncorr-1) = lcr(0:ncorr-1)
    13441391             lcr(0:ncorr-1) = 1.0_wp
    1345 
    1346 !
    1347 !--          North boundary for v.
     1392!
     1393!--          North boundary for v
    13481394             j   = nyn + 1
    13491395             kb  = nzb_v_inner(j,i)
     
    13571403
    13581404          ENDDO
     1405
    13591406       ENDIF
    13601407
    13611408!       
    1362 !--    Then vertical walls and corners if necessary.
     1409!--    Then vertical walls and corners if necessary
    13631410       IF ( topography /= 'flat' )  THEN
     1411
    13641412          kb = 0       ! kb is not used when direction > 1
    13651413!       
    13661414!--       Left boundary
    13671415          IF ( nest_bound_l )  THEN
     1416
    13681417             ALLOCATE( logc_w_l(nzb:nzt_topo_nestbc_l,nys:nyn,1:2) )
    13691418             ALLOCATE( logc_ratio_w_l(nzb:nzt_topo_nestbc_l,nys:nyn,1:2,       &
     
    13741423             DO  j = nys, nyn
    13751424                DO  k = nzb, nzt_topo_nestbc_l
    1376 
    13771425!
    13781426!--                Wall for u on the south side, but not on the north side
     
    13921440                      lcr(0:ncorr-1) = 1.0_wp
    13931441                   ENDIF
     1442
    13941443!
    13951444!--                Wall for u on the north side, but not on the south side
    13961445                   i  = 0
    1397 !--                TO_DO: routine must be indentet by 1 space from here on,
    1398 !--                       and long lines must be wrapped
    1399                   IF ( ( nzb_u_outer(j,i) > nzb_u_outer(j-1,i) ) .AND.     &
    1400                        ( nzb_u_outer(j,i) == nzb_u_outer(j+1,i) ) )  THEN
    1401                      inc        = -1                 
    1402                      wall_index =  j + 1
    1403                      CALL pmci_define_loglaw_correction_parameters( lc, lcr, k, j, inc, wall_index, z0(j,i), kb, direction, ncorr )
    1404 
    1405 !
    1406 !--                  The direction of the wall-normal index is stored as the sign of the logc-element.
    1407                      logc_u_l(k,j,2) = inc * lc
    1408                      logc_ratio_u_l(k,j,2,0:ncorr-1) = lcr(0:ncorr-1)
    1409                      lcr(0:ncorr-1) = 1.0_wp
    1410                   ENDIF
    1411 
    1412 !
    1413 !--               Wall for w on the south side, but not on the north side.
    1414                   i  = -1
    1415                   IF ( ( nzb_w_outer(j,i) > nzb_w_outer(j+1,i) )  .AND.    &
    1416                        ( nzb_w_outer(j,i) == nzb_w_outer(j-1,i) ) )  THEN
    1417                      inc        =  1
    1418                      wall_index =  j
    1419                      CALL pmci_define_loglaw_correction_parameters( lc, lcr, k, j, inc, wall_index, z0(j,i), kb, direction, ncorr )
    1420 
    1421 !
    1422 !--                  The direction of the wall-normal index is stored as the sign of the logc-element.
    1423                      logc_w_l(k,j,2) = inc * lc
    1424                      logc_ratio_w_l(k,j,2,0:ncorr-1) = lcr(0:ncorr-1)
    1425                      lcr(0:ncorr-1) = 1.0_wp
    1426                   ENDIF
    1427 
    1428 !
    1429 !--               Wall for w on the north side, but not on the south side.
    1430                   i  = -1
    1431                   IF ( ( nzb_w_outer(j,i) > nzb_w_outer(j-1,i) ) .AND.     &
    1432                        ( nzb_w_outer(j,i) == nzb_w_outer(j+1,i) ) )  THEN
    1433                      inc        = -1
    1434                      wall_index =  j+1
    1435                      CALL pmci_define_loglaw_correction_parameters( lc, lcr, k, j, inc, wall_index, z0(j,i), kb, direction, ncorr )
    1436 
    1437 !
    1438 !--                  The direction of the wall-normal index is stored as the sign of the logc-element.
    1439                      logc_w_l(k,j,2) = inc * lc
    1440                      logc_ratio_w_l(k,j,2,0:ncorr-1) = lcr(0:ncorr-1)
    1441                      lcr(0:ncorr-1) = 1.0_wp
    1442                   ENDIF
    1443                ENDDO
    1444             ENDDO
    1445          ENDIF   !  IF ( nest_bound_l )
     1446                   IF ( ( nzb_u_outer(j,i) > nzb_u_outer(j-1,i) ) .AND.        &
     1447                        ( nzb_u_outer(j,i) == nzb_u_outer(j+1,i) ) )  THEN
     1448                      inc        = -1
     1449                      wall_index =  j + 1
     1450                      CALL pmci_define_loglaw_correction_parameters( lc, lcr,  &
     1451                          k, j, inc, wall_index, z0(j,i), kb, direction, ncorr )
     1452!
     1453!--                   The direction of the wall-normal index is stored as the
     1454!--                   sign of the logc-element.
     1455                      logc_u_l(k,j,2) = inc * lc
     1456                      logc_ratio_u_l(k,j,2,0:ncorr-1) = lcr(0:ncorr-1)
     1457                      lcr(0:ncorr-1) = 1.0_wp
     1458                   ENDIF
     1459
     1460!
     1461!--                Wall for w on the south side, but not on the north side.
     1462                   i  = -1
     1463                   IF ( ( nzb_w_outer(j,i) > nzb_w_outer(j+1,i) )  .AND.       &
     1464                        ( nzb_w_outer(j,i) == nzb_w_outer(j-1,i) ) )  THEN
     1465                      inc        =  1
     1466                      wall_index =  j
     1467                      CALL pmci_define_loglaw_correction_parameters( lc, lcr,  &
     1468                          k, j, inc, wall_index, z0(j,i), kb, direction, ncorr )
     1469!
     1470!--                   The direction of the wall-normal index is stored as the
     1471!--                   sign of the logc-element.
     1472                      logc_w_l(k,j,2) = inc * lc
     1473                      logc_ratio_w_l(k,j,2,0:ncorr-1) = lcr(0:ncorr-1)
     1474                      lcr(0:ncorr-1) = 1.0_wp
     1475                   ENDIF
     1476
     1477!
     1478!--                Wall for w on the north side, but not on the south side.
     1479                   i  = -1
     1480                   IF ( ( nzb_w_outer(j,i) > nzb_w_outer(j-1,i) )  .AND.       &
     1481                        ( nzb_w_outer(j,i) == nzb_w_outer(j+1,i) ) )  THEN
     1482                      inc        = -1
     1483                      wall_index =  j+1
     1484                      CALL pmci_define_loglaw_correction_parameters( lc, lcr,  &
     1485                          k, j, inc, wall_index, z0(j,i), kb, direction, ncorr )
     1486!
     1487!--                   The direction of the wall-normal index is stored as the
     1488!--                   sign of the logc-element.
     1489                      logc_w_l(k,j,2) = inc * lc
     1490                      logc_ratio_w_l(k,j,2,0:ncorr-1) = lcr(0:ncorr-1)
     1491                      lcr(0:ncorr-1) = 1.0_wp
     1492                   ENDIF
     1493
     1494                ENDDO
     1495             ENDDO
     1496
     1497          ENDIF   !  IF ( nest_bound_l )
    14461498
    14471499!       
    1448 !--      Right boundary.
    1449          IF ( nest_bound_r )  THEN
    1450             ALLOCATE( logc_w_r(nzb:nzt_topo_nestbc_r,nys:nyn,1:2) )
    1451             ALLOCATE( logc_ratio_w_r(nzb:nzt_topo_nestbc_r,nys:nyn,1:2,0:ncorr-1) )
    1452             logc_w_r       = 0
    1453             logc_ratio_w_r = 1.0_wp
    1454             direction      = 2
    1455             i  = nxr + 1       
    1456             DO  j = nys, nyn
    1457                DO  k = nzb, nzt_topo_nestbc_r
    1458 
    1459 !
    1460 !--               Wall for u on the south side, but not on the north side.
    1461                   IF ( ( nzb_u_outer(j,i) > nzb_u_outer(j+1,i) )  .AND.    &
    1462                        ( nzb_u_outer(j,i) == nzb_u_outer(j-1,i) ) )  THEN
    1463                      inc        = 1
    1464                      wall_index = j
    1465                      CALL pmci_define_loglaw_correction_parameters( lc, lcr, k, j, inc, wall_index, z0(j,i), kb, direction, ncorr )
    1466 
    1467 !
    1468 !--                  The direction of the wall-normal index is stored as the sign of the logc-element.
    1469                      logc_u_r(k,j,2) = inc * lc
    1470                      logc_ratio_u_r(k,j,2,0:ncorr-1) = lcr(0:ncorr-1)
    1471                      lcr(0:ncorr-1) = 1.0_wp
    1472                   ENDIF
    1473 
    1474 !
    1475 !--               Wall for u on the north side, but not on the south side.
    1476                   IF ( ( nzb_u_outer(j,i) > nzb_u_outer(j-1,i) )  .AND.    &
    1477                        ( nzb_u_outer(j,i) == nzb_u_outer(j+1,i) ) )  THEN
    1478                      inc        = -1                 
    1479                      wall_index =  j+1
    1480                      CALL pmci_define_loglaw_correction_parameters( lc, lcr, k, j, inc, wall_index, z0(j,i), kb, direction, ncorr )
    1481 
    1482 !
    1483 !--                  The direction of the wall-normal index is stored as the sign of the logc-element.
    1484                      logc_u_r(k,j,2) = inc * lc
    1485                      logc_ratio_u_r(k,j,2,0:ncorr-1) = lcr(0:ncorr-1)
    1486                      lcr(0:ncorr-1) = 1.0_wp
    1487                   ENDIF
    1488 
    1489 !
    1490 !--               Wall for w on the south side, but not on the north side.
    1491                   IF ( ( nzb_w_outer(j,i) > nzb_w_outer(j+1,i) )  .AND.    &
    1492                        ( nzb_w_outer(j,i) == nzb_w_outer(j-1,i) ) )  THEN
    1493                      inc        =  1
    1494                      wall_index =  j
    1495                      CALL pmci_define_loglaw_correction_parameters( lc, lcr, k, j, inc, wall_index, z0(j,i), kb, direction, ncorr )
    1496 
    1497 !
    1498 !--                  The direction of the wall-normal index is stored as the sign of the logc-element.
    1499                      logc_w_r(k,j,2) = inc * lc
    1500                      logc_ratio_w_r(k,j,2,0:ncorr-1) = lcr(0:ncorr-1)
    1501                      lcr(0:ncorr-1) = 1.0_wp
    1502                   ENDIF
    1503 !
    1504 !--               Wall for w on the north side, but not on the south side.
    1505                   IF ( ( nzb_w_outer(j,i) > nzb_w_outer(j-1,i) ) .AND.     &
    1506                        ( nzb_w_outer(j,i) == nzb_w_outer(j+1,i) ) )  THEN
    1507                      inc        = -1
    1508                      wall_index =  j+1
    1509                      CALL pmci_define_loglaw_correction_parameters( lc, lcr, k, j, inc, wall_index, z0(j,i), kb, direction, ncorr )
    1510 
    1511 !
    1512 !--                  The direction of the wall-normal index is stored as the sign of the logc-element.
    1513                      logc_w_r(k,j,2) = inc * lc
    1514                      logc_ratio_w_r(k,j,2,0:ncorr-1) = lcr(0:ncorr-1)
    1515                      lcr(0:ncorr-1) = 1.0_wp
    1516                   ENDIF
    1517                ENDDO
    1518             ENDDO
    1519          ENDIF   !  IF ( nest_bound_r )
     1500!--       Right boundary
     1501          IF ( nest_bound_r )  THEN
     1502
     1503             ALLOCATE( logc_w_r(nzb:nzt_topo_nestbc_r,nys:nyn,1:2) )
     1504             ALLOCATE( logc_ratio_w_r(nzb:nzt_topo_nestbc_r,nys:nyn,1:2,       &
     1505                                      0:ncorr-1) )
     1506             logc_w_r       = 0
     1507             logc_ratio_w_r = 1.0_wp
     1508             direction      = 2
     1509             i  = nxr + 1
     1510
     1511             DO  j = nys, nyn
     1512                DO  k = nzb, nzt_topo_nestbc_r
     1513!
     1514!--                Wall for u on the south side, but not on the north side
     1515                   IF ( ( nzb_u_outer(j,i) > nzb_u_outer(j+1,i) )  .AND.       &
     1516                        ( nzb_u_outer(j,i) == nzb_u_outer(j-1,i) ) )  THEN
     1517                      inc        = 1
     1518                      wall_index = j
     1519                      CALL pmci_define_loglaw_correction_parameters( lc, lcr,  &
     1520                          k, j, inc, wall_index, z0(j,i), kb, direction, ncorr )
     1521!
     1522!--                   The direction of the wall-normal index is stored as the
     1523!--                   sign of the logc-element.
     1524                      logc_u_r(k,j,2) = inc * lc
     1525                      logc_ratio_u_r(k,j,2,0:ncorr-1) = lcr(0:ncorr-1)
     1526                      lcr(0:ncorr-1) = 1.0_wp
     1527                   ENDIF
     1528
     1529!
     1530!--                Wall for u on the north side, but not on the south side
     1531                   IF ( ( nzb_u_outer(j,i) > nzb_u_outer(j-1,i) )  .AND.       &
     1532                        ( nzb_u_outer(j,i) == nzb_u_outer(j+1,i) ) )  THEN
     1533                      inc        = -1
     1534                      wall_index =  j+1
     1535                      CALL pmci_define_loglaw_correction_parameters( lc, lcr,  &
     1536                          k, j, inc, wall_index, z0(j,i), kb, direction, ncorr )
     1537!
     1538!--                   The direction of the wall-normal index is stored as the
     1539!--                   sign of the logc-element.
     1540                      logc_u_r(k,j,2) = inc * lc
     1541                      logc_ratio_u_r(k,j,2,0:ncorr-1) = lcr(0:ncorr-1)
     1542                      lcr(0:ncorr-1) = 1.0_wp
     1543                   ENDIF
     1544
     1545!
     1546!--                Wall for w on the south side, but not on the north side
     1547                   IF ( ( nzb_w_outer(j,i) > nzb_w_outer(j+1,i) )  .AND.       &
     1548                        ( nzb_w_outer(j,i) == nzb_w_outer(j-1,i) ) )  THEN
     1549                      inc        =  1
     1550                      wall_index =  j
     1551                      CALL pmci_define_loglaw_correction_parameters( lc, lcr,  &
     1552                          k, j, inc, wall_index, z0(j,i), kb, direction, ncorr )
     1553!
     1554!--                   The direction of the wall-normal index is stored as the
     1555!--                   sign of the logc-element.
     1556                      logc_w_r(k,j,2) = inc * lc
     1557                      logc_ratio_w_r(k,j,2,0:ncorr-1) = lcr(0:ncorr-1)
     1558                      lcr(0:ncorr-1) = 1.0_wp
     1559                   ENDIF
     1560
     1561!
     1562!--                Wall for w on the north side, but not on the south side
     1563                   IF ( ( nzb_w_outer(j,i) > nzb_w_outer(j-1,i) )  .AND.       &
     1564                        ( nzb_w_outer(j,i) == nzb_w_outer(j+1,i) ) )  THEN
     1565                      inc        = -1
     1566                      wall_index =  j+1
     1567                      CALL pmci_define_loglaw_correction_parameters( lc, lcr,  &
     1568                          k, j, inc, wall_index, z0(j,i), kb, direction, ncorr )
     1569
     1570!
     1571!--                   The direction of the wall-normal index is stored as the
     1572!--                   sign of the logc-element.
     1573                      logc_w_r(k,j,2) = inc * lc
     1574                      logc_ratio_w_r(k,j,2,0:ncorr-1) = lcr(0:ncorr-1)
     1575                      lcr(0:ncorr-1) = 1.0_wp
     1576                   ENDIF
     1577
     1578                ENDDO
     1579             ENDDO
     1580
     1581          ENDIF   !  IF ( nest_bound_r )
    15201582
    15211583!       
    1522 !--      South boundary.
    1523          IF ( nest_bound_s )  THEN
    1524             ALLOCATE( logc_w_s(nzb:nzt_topo_nestbc_s, nxl:nxr, 1:2) )
    1525             ALLOCATE( logc_ratio_w_s(nzb:nzt_topo_nestbc_s, nxl:nxr, 1:2, 0:ncorr-1) )
    1526             logc_w_s       = 0
    1527             logc_ratio_w_s = 1.0_wp 
    1528             direction      = 3
    1529             DO  i = nxl, nxr
    1530                DO  k = nzb, nzt_topo_nestbc_s
    1531 
    1532 !
    1533 !--               Wall for v on the left side, but not on the right side.
    1534                   j  = 0
    1535                   IF ( ( nzb_v_outer(j,i) > nzb_v_outer(j,i+1) )  .AND.    &
    1536                        ( nzb_v_outer(j,i) == nzb_v_outer(j,i-1) ) )  THEN
    1537                      inc        =  1
    1538                      wall_index =  i
    1539                      CALL pmci_define_loglaw_correction_parameters( lc, lcr, k, i, inc, wall_index, z0(j,i), kb, direction, ncorr )
    1540 
    1541 !
    1542 !--                  The direction of the wall-normal index is stored as the sign of the logc-element.
    1543                      logc_v_s(k,i,2) = inc * lc
    1544                      logc_ratio_v_s(k,i,2,0:ncorr-1) = lcr(0:ncorr-1)
    1545                      lcr(0:ncorr-1) = 1.0_wp
    1546                   ENDIF
    1547 !
    1548 !--               Wall for v on the right side, but not on the left side.
    1549                   j  = 0
    1550                   IF ( ( nzb_v_outer(j,i) > nzb_v_outer(j,i-1) )  .AND.    &
    1551                        ( nzb_v_outer(j,i) == nzb_v_outer(j,i+1) ) )  THEN
    1552                      inc        = -1
    1553                      wall_index =  i+1
    1554                      CALL pmci_define_loglaw_correction_parameters( lc, lcr, k, i, inc, wall_index, z0(j,i), kb, direction, ncorr )
    1555 
    1556 !
    1557 !--                  The direction of the wall-normal index is stored as the sign of the logc-element.
    1558                      logc_v_s(k,i,2) = inc * lc
    1559                      logc_ratio_v_s(k,i,2,0:ncorr-1) = lcr(0:ncorr-1)
    1560                      lcr(0:ncorr-1) = 1.0_wp
    1561                   ENDIF
    1562 
    1563 !
    1564 !--               Wall for w on the left side, but not on the right side.
    1565                   j  = -1
    1566                   IF ( ( nzb_w_outer(j,i) > nzb_w_outer(j,i+1) )  .AND.    &
    1567                        ( nzb_w_outer(j,i) == nzb_w_outer(j,i-1) ) )  THEN
    1568                      inc        =  1
    1569                      wall_index =  i
    1570                      CALL pmci_define_loglaw_correction_parameters( lc, lcr, k, i, inc, wall_index, z0(j,i), kb, direction, ncorr )
    1571 
    1572 !
    1573 !--                  The direction of the wall-normal index is stored as the sign of the logc-element.
    1574                      logc_w_s(k,i,2) = inc * lc
    1575                      logc_ratio_w_s(k,i,2,0:ncorr - 1) = lcr(0:ncorr-1)
    1576                      lcr(0:ncorr-1) = 1.0_wp
    1577                   ENDIF
    1578 
    1579 !
    1580 !--               Wall for w on the right side, but not on the left side.
    1581                   j  = -1
    1582                   IF ( ( nzb_w_outer(j,i) > nzb_w_outer(j,i-1) )  .AND.    &
    1583                        ( nzb_w_outer(j,i) == nzb_w_outer(j,i+1) ) )  THEN
    1584                      inc        = -1
    1585                      wall_index =  i+1
    1586                      CALL pmci_define_loglaw_correction_parameters( lc, lcr, k, i, inc, wall_index, z0(j,i), kb, direction, ncorr )
    1587 
    1588 !
    1589 !--                  The direction of the wall-normal index is stored as the sign of the logc-element.
    1590                      logc_w_s(k,i,2) = inc * lc
    1591                      logc_ratio_w_s(k,i,2,0:ncorr-1) = lcr(0:ncorr-1)
    1592                      lcr(0:ncorr-1) = 1.0_wp
    1593                   ENDIF
    1594                ENDDO
    1595             ENDDO
    1596          ENDIF   !  IF (nest_bound_s )
     1584!--       South boundary
     1585          IF ( nest_bound_s )  THEN
     1586
     1587             ALLOCATE( logc_w_s(nzb:nzt_topo_nestbc_s, nxl:nxr, 1:2) )
     1588             ALLOCATE( logc_ratio_w_s(nzb:nzt_topo_nestbc_s,nxl:nxr,1:2,       &
     1589                                      0:ncorr-1) )
     1590             logc_w_s       = 0
     1591             logc_ratio_w_s = 1.0_wp
     1592             direction      = 3
     1593
     1594             DO  i = nxl, nxr
     1595                DO  k = nzb, nzt_topo_nestbc_s
     1596!
     1597!--                Wall for v on the left side, but not on the right side
     1598                   j  = 0
     1599                   IF ( ( nzb_v_outer(j,i) > nzb_v_outer(j,i+1) )  .AND.       &
     1600                        ( nzb_v_outer(j,i) == nzb_v_outer(j,i-1) ) )  THEN
     1601                      inc        =  1
     1602                      wall_index =  i
     1603                      CALL pmci_define_loglaw_correction_parameters( lc, lcr,  &
     1604                          k, i, inc, wall_index, z0(j,i), kb, direction, ncorr )
     1605!
     1606!--                   The direction of the wall-normal index is stored as the
     1607!--                   sign of the logc-element.
     1608                      logc_v_s(k,i,2) = inc * lc
     1609                      logc_ratio_v_s(k,i,2,0:ncorr-1) = lcr(0:ncorr-1)
     1610                      lcr(0:ncorr-1) = 1.0_wp
     1611                   ENDIF
     1612
     1613!
     1614!--                Wall for v on the right side, but not on the left side
     1615                   j  = 0
     1616                   IF ( ( nzb_v_outer(j,i) > nzb_v_outer(j,i-1) )  .AND.       &
     1617                        ( nzb_v_outer(j,i) == nzb_v_outer(j,i+1) ) )  THEN
     1618                      inc        = -1
     1619                      wall_index =  i+1
     1620                      CALL pmci_define_loglaw_correction_parameters( lc, lcr,  &
     1621                          k, i, inc, wall_index, z0(j,i), kb, direction, ncorr )
     1622!
     1623!--                   The direction of the wall-normal index is stored as the
     1624!--                   sign of the logc-element.
     1625                      logc_v_s(k,i,2) = inc * lc
     1626                      logc_ratio_v_s(k,i,2,0:ncorr-1) = lcr(0:ncorr-1)
     1627                      lcr(0:ncorr-1) = 1.0_wp
     1628                   ENDIF
     1629
     1630!
     1631!--                Wall for w on the left side, but not on the right side
     1632                   j  = -1
     1633                   IF ( ( nzb_w_outer(j,i) > nzb_w_outer(j,i+1) )  .AND.       &
     1634                        ( nzb_w_outer(j,i) == nzb_w_outer(j,i-1) ) )  THEN
     1635                      inc        =  1
     1636                      wall_index =  i
     1637                      CALL pmci_define_loglaw_correction_parameters( lc, lcr,  &
     1638                          k, i, inc, wall_index, z0(j,i), kb, direction, ncorr )
     1639!
     1640!--                   The direction of the wall-normal index is stored as the
     1641!--                   sign of the logc-element.
     1642                      logc_w_s(k,i,2) = inc * lc
     1643                      logc_ratio_w_s(k,i,2,0:ncorr - 1) = lcr(0:ncorr-1)
     1644                      lcr(0:ncorr-1) = 1.0_wp
     1645                   ENDIF
     1646
     1647!
     1648!--                Wall for w on the right side, but not on the left side
     1649                   j  = -1
     1650                   IF ( ( nzb_w_outer(j,i) > nzb_w_outer(j,i-1) )  .AND.       &
     1651                        ( nzb_w_outer(j,i) == nzb_w_outer(j,i+1) ) )  THEN
     1652                      inc        = -1
     1653                      wall_index =  i+1
     1654                      CALL pmci_define_loglaw_correction_parameters( lc, lcr,  &
     1655                          k, i, inc, wall_index, z0(j,i), kb, direction, ncorr )
     1656!
     1657!--                   The direction of the wall-normal index is stored as the
     1658!--                   sign of the logc-element.
     1659                      logc_w_s(k,i,2) = inc * lc
     1660                      logc_ratio_w_s(k,i,2,0:ncorr-1) = lcr(0:ncorr-1)
     1661                      lcr(0:ncorr-1) = 1.0_wp
     1662                   ENDIF
     1663
     1664                ENDDO
     1665             ENDDO
     1666
     1667          ENDIF   !  IF (nest_bound_s )
    15971668
    15981669!       
    1599 !--      North boundary.
    1600          IF ( nest_bound_n )  THEN
    1601             ALLOCATE( logc_w_n(nzb:nzt_topo_nestbc_n, nxl:nxr, 1:2) )
    1602             ALLOCATE( logc_ratio_w_n(nzb:nzt_topo_nestbc_n, nxl:nxr, 1:2, 0:ncorr-1) )
    1603             logc_w_n       = 0
    1604             logc_ratio_w_n = 1.0_wp
    1605             direction      = 3
    1606             j  = nyn + 1
    1607             DO  i = nxl, nxr
    1608                DO  k = nzb, nzt_topo_nestbc_n
    1609 
    1610 !
    1611 !--               Wall for v on the left side, but not on the right side.
    1612                   IF ( ( nzb_v_outer(j,i) > nzb_v_outer(j,i+1) )  .AND.    &
    1613                        ( nzb_v_outer(j,i) == nzb_v_outer(j,i-1) ) )  THEN
    1614                      inc        = 1
    1615                      wall_index = i
    1616                      CALL pmci_define_loglaw_correction_parameters( lc, lcr, k, i, inc, wall_index, z0(j,i), kb, direction, ncorr )
    1617 
    1618 !
    1619 !--                  The direction of the wall-normal index is stored as the sign of the logc-element.
    1620                      logc_v_n(k,i,2) = inc * lc
    1621                      logc_ratio_v_n(k,i,2,0:ncorr-1) = lcr(0:ncorr-1)
    1622                      lcr(0:ncorr-1) = 1.0_wp
    1623                   ENDIF
    1624 
    1625 !
    1626 !--               Wall for v on the right side, but not on the left side.
    1627                   IF ( ( nzb_v_outer(j,i) > nzb_v_outer(j,i-1) )  .AND.    &
    1628                        ( nzb_v_outer(j,i) == nzb_v_outer(j,i+1) ) )  THEN
    1629                      inc        = -1                 
    1630                      wall_index =  i + 1
    1631                      CALL pmci_define_loglaw_correction_parameters( lc, lcr, k, i, inc, wall_index, z0(j,i), kb, direction, ncorr )
    1632 
    1633 !
    1634 !--                  The direction of the wall-normal index is stored as the sign of the logc-element.
    1635                      logc_v_n(k,i,2) = inc * lc
    1636                      logc_ratio_v_n(k,i,2,0:ncorr-1) = lcr(0:ncorr-1)
    1637                      lcr(0:ncorr-1) = 1.0_wp
    1638                   ENDIF
    1639 
    1640 !
    1641 !--               Wall for w on the left side, but not on the right side.
    1642                   IF ( ( nzb_w_outer(j,i) > nzb_w_outer(j,i+1) )  .AND.    &
    1643                        ( nzb_w_outer(j,i) == nzb_w_outer(j,i-1) ) )  THEN
    1644                      inc        = 1
    1645                      wall_index = i
    1646                      CALL pmci_define_loglaw_correction_parameters( lc, lcr, k, i, inc, wall_index, z0(j,i), kb, direction, ncorr )
    1647 
    1648 !
    1649 !--                  The direction of the wall-normal index is stored as the sign of the logc-element.
    1650                      logc_w_n(k,i,2) = inc * lc
    1651                      logc_ratio_w_n(k,i,2,0:ncorr-1) = lcr(0:ncorr-1)
    1652                      lcr(0:ncorr-1) = 1.0_wp
    1653                   ENDIF
    1654 
    1655 !
    1656 !--               Wall for w on the right side, but not on the left side.
    1657                   IF ( ( nzb_w_outer(j,i) > nzb_w_outer(j,i-1) ) .AND.     &
    1658                        ( nzb_w_outer(j,i) == nzb_w_outer(j,i+1) ) )  THEN
    1659                      inc        = -1
    1660                      wall_index =  i+1
    1661                      CALL pmci_define_loglaw_correction_parameters( lc, lcr, k, i, inc, wall_index, z0(j,i), kb, direction, ncorr )
    1662 
    1663 !
    1664 !--                  The direction of the wall-normal index is stored as the sign of the logc-element.
    1665                      logc_w_n(k,i,2) = inc * lc
    1666                      logc_ratio_w_n(k,i,2,0:ncorr-1) = lcr(0:ncorr-1)
    1667                      lcr(0:ncorr-1) = 1.0_wp
    1668                   ENDIF
    1669                ENDDO
    1670             ENDDO
    1671          ENDIF   !  IF ( nest_bound_n )
    1672       ENDIF   !  IF ( topography /= 'flat' )
    1673 
    1674    END SUBROUTINE pmci_init_loglaw_correction
     1670!--       North boundary
     1671          IF ( nest_bound_n )  THEN
     1672
     1673             ALLOCATE( logc_w_n(nzb:nzt_topo_nestbc_n, nxl:nxr, 1:2) )
     1674             ALLOCATE( logc_ratio_w_n(nzb:nzt_topo_nestbc_n,nxl:nxr,1:2,       &
     1675                                      0:ncorr-1) )
     1676             logc_w_n       = 0
     1677             logc_ratio_w_n = 1.0_wp
     1678             direction      = 3
     1679             j  = nyn + 1
     1680
     1681             DO  i = nxl, nxr
     1682                DO  k = nzb, nzt_topo_nestbc_n
     1683!
     1684!--                Wall for v on the left side, but not on the right side
     1685                   IF ( ( nzb_v_outer(j,i) > nzb_v_outer(j,i+1) )  .AND.       &
     1686                        ( nzb_v_outer(j,i) == nzb_v_outer(j,i-1) ) )  THEN
     1687                      inc        = 1
     1688                      wall_index = i
     1689                      CALL pmci_define_loglaw_correction_parameters( lc, lcr,  &
     1690                          k, i, inc, wall_index, z0(j,i), kb, direction, ncorr )
     1691!
     1692!--                   The direction of the wall-normal index is stored as the
     1693!--                   sign of the logc-element.
     1694                      logc_v_n(k,i,2) = inc * lc
     1695                      logc_ratio_v_n(k,i,2,0:ncorr-1) = lcr(0:ncorr-1)
     1696                      lcr(0:ncorr-1) = 1.0_wp
     1697                   ENDIF
     1698
     1699!
     1700!--                Wall for v on the right side, but not on the left side
     1701                   IF ( ( nzb_v_outer(j,i) > nzb_v_outer(j,i-1) )  .AND.       &
     1702                        ( nzb_v_outer(j,i) == nzb_v_outer(j,i+1) ) )  THEN
     1703                      inc        = -1
     1704                      wall_index =  i + 1
     1705                      CALL pmci_define_loglaw_correction_parameters( lc, lcr,  &
     1706                          k, i, inc, wall_index, z0(j,i), kb, direction, ncorr )
     1707!
     1708!--                   The direction of the wall-normal index is stored as the
     1709!--                   sign of the logc-element.
     1710                      logc_v_n(k,i,2) = inc * lc
     1711                      logc_ratio_v_n(k,i,2,0:ncorr-1) = lcr(0:ncorr-1)
     1712                      lcr(0:ncorr-1) = 1.0_wp
     1713                   ENDIF
     1714
     1715!
     1716!--                Wall for w on the left side, but not on the right side
     1717                   IF ( ( nzb_w_outer(j,i) > nzb_w_outer(j,i+1) )  .AND.       &
     1718                        ( nzb_w_outer(j,i) == nzb_w_outer(j,i-1) ) )  THEN
     1719                      inc        = 1
     1720                      wall_index = i
     1721                      CALL pmci_define_loglaw_correction_parameters( lc, lcr,  &
     1722                          k, i, inc, wall_index, z0(j,i), kb, direction, ncorr )
     1723!
     1724!--                   The direction of the wall-normal index is stored as the
     1725!--                   sign of the logc-element.
     1726                      logc_w_n(k,i,2) = inc * lc
     1727                      logc_ratio_w_n(k,i,2,0:ncorr-1) = lcr(0:ncorr-1)
     1728                      lcr(0:ncorr-1) = 1.0_wp
     1729                   ENDIF
     1730
     1731!
     1732!--                Wall for w on the right side, but not on the left side
     1733                   IF ( ( nzb_w_outer(j,i) > nzb_w_outer(j,i-1) )  .AND.       &
     1734                        ( nzb_w_outer(j,i) == nzb_w_outer(j,i+1) ) )  THEN
     1735                      inc        = -1
     1736                      wall_index =  i+1
     1737                      CALL pmci_define_loglaw_correction_parameters( lc, lcr,  &
     1738                          k, i, inc, wall_index, z0(j,i), kb, direction, ncorr )
     1739!
     1740!--                   The direction of the wall-normal index is stored as the
     1741!--                   sign of the logc-element.
     1742                      logc_w_n(k,i,2) = inc * lc
     1743                      logc_ratio_w_n(k,i,2,0:ncorr-1) = lcr(0:ncorr-1)
     1744                      lcr(0:ncorr-1) = 1.0_wp
     1745                   ENDIF
     1746
     1747                ENDDO
     1748             ENDDO
     1749
     1750          ENDIF   !  IF ( nest_bound_n )
     1751
     1752       ENDIF   !  IF ( topography /= 'flat' )
     1753
     1754    END SUBROUTINE pmci_init_loglaw_correction
    16751755
    16761756
     
    19011981
    19021982
    1903 !-- TO_DO:  indentation and wrap long lines from here on to the end of the file
    1904    SUBROUTINE pmci_init_anterp_tophat
    1905 !
    1906 !--   Precomputation of the client-array indices for
    1907 !--   corresponding coarse-grid array index and the
    1908 !--   Under-relaxation coefficients to be used by anterp_tophat.
    1909 
    1910       IMPLICIT NONE
    1911 
    1912       INTEGER(iwp) ::  i        !: Fine-grid index
    1913       INTEGER(iwp) ::  ii       !: Coarse-grid index
    1914       INTEGER(iwp) ::  istart   !:
    1915       INTEGER(iwp) ::  j        !: Fine-grid index
    1916       INTEGER(iwp) ::  jj       !: Coarse-grid index
    1917       INTEGER(iwp) ::  jstart   !:
    1918       INTEGER(iwp) ::  k        !: Fine-grid index
    1919       INTEGER(iwp) ::  kk       !: Coarse-grid index
    1920       INTEGER(iwp) ::  kstart   !:
    1921       REAL(wp)     ::  xi       !:
    1922       REAL(wp)     ::  eta      !:
    1923       REAL(wp)     ::  zeta     !:
     1983
     1984
     1985    SUBROUTINE pmci_init_anterp_tophat
     1986!
     1987!--    Precomputation of the client-array indices for
     1988!--    corresponding coarse-grid array index and the
     1989!--    Under-relaxation coefficients to be used by anterp_tophat.
     1990
     1991       IMPLICIT NONE
     1992
     1993       INTEGER(iwp) ::  i        !: Fine-grid index
     1994       INTEGER(iwp) ::  ii       !: Coarse-grid index
     1995       INTEGER(iwp) ::  istart   !:
     1996       INTEGER(iwp) ::  j        !: Fine-grid index
     1997       INTEGER(iwp) ::  jj       !: Coarse-grid index
     1998       INTEGER(iwp) ::  jstart   !:
     1999       INTEGER(iwp) ::  k        !: Fine-grid index
     2000       INTEGER(iwp) ::  kk       !: Coarse-grid index
     2001       INTEGER(iwp) ::  kstart   !:
     2002       REAL(wp)     ::  xi       !:
     2003       REAL(wp)     ::  eta      !:
     2004       REAL(wp)     ::  zeta     !:
    19242005     
    1925 
    1926 !
    1927 !--   Default values:
    1928       IF ( anterp_relax_length_l < 0.0_wp )  THEN
    1929          anterp_relax_length_l = 0.1_wp * ( nx + 1 ) * dx
    1930       ENDIF
    1931       IF ( anterp_relax_length_r < 0.0_wp )  THEN
    1932          anterp_relax_length_r = 0.1_wp * ( nx + 1 ) * dx
    1933       ENDIF
    1934       IF ( anterp_relax_length_s < 0.0_wp )  THEN
    1935          anterp_relax_length_s = 0.1_wp * ( ny + 1 ) * dy
    1936       ENDIF
    1937       IF ( anterp_relax_length_n < 0.0_wp )  THEN
    1938          anterp_relax_length_n = 0.1_wp * ( ny + 1 ) * dy
    1939       ENDIF
    1940       IF ( anterp_relax_length_t < 0.0_wp )  THEN
    1941          anterp_relax_length_t = 0.1_wp * zu(nzt)
    1942       ENDIF
    1943 
    1944 !
    1945 !--   First determine kctu and kctw that are the coarse-grid upper bounds for index k.
    1946       kk = 0
    1947       DO  WHILE ( cg%zu(kk) < zu(nzt) )
    1948          kk = kk + 1
    1949       ENDDO
    1950       kctu = kk - 1
    1951 
    1952       kk = 0
    1953       DO  WHILE ( cg%zw(kk) < zw(nzt-1) )
    1954          kk = kk + 1
    1955       ENDDO
    1956       kctw = kk - 1
    1957 
    1958       ALLOCATE( iflu(icl:icr) )
    1959       ALLOCATE( iflo(icl:icr) )
    1960       ALLOCATE( ifuu(icl:icr) )
    1961       ALLOCATE( ifuo(icl:icr) )
    1962       ALLOCATE( jflv(jcs:jcn) )
    1963       ALLOCATE( jflo(jcs:jcn) )
    1964       ALLOCATE( jfuv(jcs:jcn) )
    1965       ALLOCATE( jfuo(jcs:jcn) )
    1966       ALLOCATE( kflw(0:kctw) )
    1967       ALLOCATE( kflo(0:kctu) )
    1968       ALLOCATE( kfuw(0:kctw) )
    1969       ALLOCATE( kfuo(0:kctu) )
    1970 
    1971 !
    1972 !--   i-indices of u for each l-index value.   
    1973       istart = nxlg
    1974       DO  ii = icl, icr 
    1975          i = istart
    1976          DO  WHILE ( ( coord_x(i)  <  cg%coord_x(ii) - 0.5_wp * cg%dx )  .AND.  ( i < nxrg ) )
    1977             i = i + 1
    1978          ENDDO
    1979          iflu(ii) = MIN( MAX( i, nxlg ), nxrg )
    1980          DO  WHILE ( ( coord_x(i)  <  cg%coord_x(ii) + 0.5_wp * cg%dx )  .AND.  ( i < nxrg ) )
    1981             i = i + 1
    1982          ENDDO
    1983          ifuu(ii) = MIN( MAX( i, nxlg ), nxrg )
    1984          istart = iflu(ii)
    1985       ENDDO
    1986 
    1987 !
    1988 !--   i-indices of others for each l-index value.
    1989       istart = nxlg
    1990       DO  ii = icl, icr   
    1991          i = istart
    1992          DO  WHILE ( ( coord_x(i) + 0.5_wp * dx  <  cg%coord_x(ii) )  .AND.  ( i < nxrg ) )
    1993             i = i + 1
    1994          ENDDO
    1995          iflo(ii) = MIN( MAX( i, nxlg ), nxrg )
    1996          DO  WHILE ( ( coord_x(i) + 0.5_wp * dx  <  cg%coord_x(ii) + cg%dx )  .AND.  ( i < nxrg ) )
    1997             i = i + 1
    1998          ENDDO
    1999          ifuo(ii) = MIN(MAX(i,nxlg),nxrg)
    2000          istart = iflo(ii)
    2001       ENDDO
    2002 
    2003 !
    2004 !--   j-indices of v for each m-index value.
    2005       jstart = nysg
    2006       DO  jj = jcs, jcn
    2007          j = jstart
    2008          DO  WHILE ( ( coord_y(j)  <  cg%coord_y(jj) - 0.5_wp * cg%dy )  .AND.  ( j < nyng ) )
    2009             j = j + 1
    2010          ENDDO
    2011          jflv(jj) = MIN( MAX( j, nysg ), nyng )
    2012          DO  WHILE ( ( coord_y(j)  <  cg%coord_y(jj) + 0.5_wp * cg%dy )  .AND.  ( j < nyng ) )
    2013             j = j + 1
    2014          ENDDO
    2015          jfuv(jj) = MIN( MAX( j, nysg ), nyng )
    2016          jstart = jflv(jj)
    2017       ENDDO
    2018 
    2019 !
    2020 !--   j-indices of others for each m-index value.
    2021       jstart = nysg
    2022       DO  jj = jcs, jcn
    2023          j = jstart
    2024          DO  WHILE ( ( coord_y(j) + 0.5_wp * dy  <  cg%coord_y(jj) )  .AND.  ( j < nyng ) ) 
    2025             j = j + 1
    2026          ENDDO
    2027          jflo(jj) = MIN( MAX( j, nysg ), nyng )
    2028          DO  WHILE ( ( coord_y(j) + 0.5_wp * dy  <  cg%coord_y(jj) + cg%dy )  .AND.  ( j < nyng ) )
    2029             j = j + 1
    2030          ENDDO
    2031          jfuo(jj) = MIN( MAX( j, nysg ), nyng )
    2032          jstart = jflv(jj)
    2033       ENDDO
    2034 
    2035 !
    2036 !--   k-indices of w for each n-index value.
    2037       kstart  = 0
    2038       kflw(0) = 0
    2039       kfuw(0) = 0
    2040       DO  kk = 1, kctw
    2041          k = kstart
    2042          DO  WHILE ( ( zw(k) < cg%zw(kk) - 0.5_wp * cg%dzw(kk) )  .AND.  ( k < nzt ) )
    2043             k = k + 1
    2044          ENDDO
    2045          kflw(kk) = MIN( MAX( k, 1 ), nzt + 1 )
    2046          DO  WHILE ( ( zw(k) < cg%zw(kk) + 0.5_wp * cg%dzw(kk+1) )  .AND.  ( k < nzt ) )
    2047             k = k + 1
    2048          ENDDO
    2049          kfuw(kk) = MIN( MAX( k, 1 ), nzt + 1 )
    2050          kstart = kflw(kk)
    2051       ENDDO
    2052 
    2053 !
    2054 !--   k-indices of others for each n-index value.
    2055       kstart  = 0
    2056       kflo(0) = 0
    2057       kfuo(0) = 0
    2058       DO  kk = 1, kctu
    2059          k = kstart
    2060          DO  WHILE ( ( zu(k) < cg%zu(kk) - 0.5_wp * cg%dzu(kk) )  .AND.  ( k < nzt ) )
    2061             k = k + 1
    2062          ENDDO
    2063          kflo(kk) = MIN( MAX( k, 1 ), nzt + 1 )
    2064          DO  WHILE ( ( zu(k) < cg%zu(kk) + 0.5_wp * cg%dzu(kk+1) )  .AND.  ( k < nzt ) ) 
    2065             k = k + 1
    2066          ENDDO
    2067          kfuo(kk) = MIN( MAX( k-1, 1 ), nzt + 1 )
    2068          kstart = kflo(kk)
    2069       ENDDO
     2006!
     2007!--    Default values:
     2008       IF ( anterp_relax_length_l < 0.0_wp )  THEN
     2009          anterp_relax_length_l = 0.1_wp * ( nx + 1 ) * dx
     2010       ENDIF
     2011       IF ( anterp_relax_length_r < 0.0_wp )  THEN
     2012          anterp_relax_length_r = 0.1_wp * ( nx + 1 ) * dx
     2013       ENDIF
     2014       IF ( anterp_relax_length_s < 0.0_wp )  THEN
     2015          anterp_relax_length_s = 0.1_wp * ( ny + 1 ) * dy
     2016       ENDIF
     2017       IF ( anterp_relax_length_n < 0.0_wp )  THEN
     2018          anterp_relax_length_n = 0.1_wp * ( ny + 1 ) * dy
     2019       ENDIF
     2020       IF ( anterp_relax_length_t < 0.0_wp )  THEN
     2021          anterp_relax_length_t = 0.1_wp * zu(nzt)
     2022       ENDIF
     2023
     2024!
     2025!--    First determine kctu and kctw that are the coarse-grid upper bounds for
     2026!--    index k
     2027       kk = 0
     2028       DO  WHILE ( cg%zu(kk) < zu(nzt) )
     2029          kk = kk + 1
     2030       ENDDO
     2031       kctu = kk - 1
     2032
     2033       kk = 0
     2034       DO  WHILE ( cg%zw(kk) < zw(nzt-1) )
     2035          kk = kk + 1
     2036       ENDDO
     2037       kctw = kk - 1
     2038
     2039       ALLOCATE( iflu(icl:icr) )
     2040       ALLOCATE( iflo(icl:icr) )
     2041       ALLOCATE( ifuu(icl:icr) )
     2042       ALLOCATE( ifuo(icl:icr) )
     2043       ALLOCATE( jflv(jcs:jcn) )
     2044       ALLOCATE( jflo(jcs:jcn) )
     2045       ALLOCATE( jfuv(jcs:jcn) )
     2046       ALLOCATE( jfuo(jcs:jcn) )
     2047       ALLOCATE( kflw(0:kctw) )
     2048       ALLOCATE( kflo(0:kctu) )
     2049       ALLOCATE( kfuw(0:kctw) )
     2050       ALLOCATE( kfuo(0:kctu) )
     2051
     2052!
     2053!--    i-indices of u for each ii-index value
     2054       istart = nxlg
     2055       DO  ii = icl, icr
     2056          i = istart
     2057          DO  WHILE ( ( coord_x(i) < cg%coord_x(ii) - 0.5_wp * cg%dx )  .AND.  &
     2058                      ( i < nxrg ) )
     2059             i = i + 1
     2060          ENDDO
     2061          iflu(ii) = MIN( MAX( i, nxlg ), nxrg )
     2062          DO  WHILE ( ( coord_x(i) < cg%coord_x(ii) + 0.5_wp * cg%dx )  .AND.  &
     2063                      ( i < nxrg ) )
     2064             i = i + 1
     2065          ENDDO
     2066          ifuu(ii) = MIN( MAX( i, nxlg ), nxrg )
     2067          istart = iflu(ii)
     2068       ENDDO
     2069
     2070!
     2071!--    i-indices of others for each ii-index value
     2072       istart = nxlg
     2073       DO  ii = icl, icr
     2074          i = istart
     2075          DO  WHILE ( ( coord_x(i) + 0.5_wp * dx < cg%coord_x(ii) )  .AND.     &
     2076                      ( i < nxrg ) )
     2077             i = i + 1
     2078          ENDDO
     2079          iflo(ii) = MIN( MAX( i, nxlg ), nxrg )
     2080          DO  WHILE ( ( coord_x(i) + 0.5_wp * dx < cg%coord_x(ii) + cg%dx )    &
     2081                      .AND.  ( i < nxrg ) )
     2082             i = i + 1
     2083          ENDDO
     2084          ifuo(ii) = MIN(MAX(i,nxlg),nxrg)
     2085          istart = iflo(ii)
     2086       ENDDO
     2087
     2088!
     2089!--    j-indices of v for each jj-index value
     2090       jstart = nysg
     2091       DO  jj = jcs, jcn
     2092          j = jstart
     2093          DO  WHILE ( ( coord_y(j) < cg%coord_y(jj) - 0.5_wp * cg%dy )  .AND.  &
     2094                      ( j < nyng ) )
     2095             j = j + 1
     2096          ENDDO
     2097          jflv(jj) = MIN( MAX( j, nysg ), nyng )
     2098          DO  WHILE ( ( coord_y(j) < cg%coord_y(jj) + 0.5_wp * cg%dy )  .AND.  &
     2099                      ( j < nyng ) )
     2100             j = j + 1
     2101          ENDDO
     2102          jfuv(jj) = MIN( MAX( j, nysg ), nyng )
     2103          jstart = jflv(jj)
     2104       ENDDO
     2105
     2106!
     2107!--    j-indices of others for each jj-index value
     2108       jstart = nysg
     2109       DO  jj = jcs, jcn
     2110          j = jstart
     2111          DO  WHILE ( ( coord_y(j) + 0.5_wp * dy < cg%coord_y(jj) )  .AND.     &
     2112                      ( j < nyng ) )
     2113             j = j + 1
     2114          ENDDO
     2115          jflo(jj) = MIN( MAX( j, nysg ), nyng )
     2116          DO  WHILE ( ( coord_y(j) + 0.5_wp * dy < cg%coord_y(jj) + cg%dy )    &
     2117                      .AND.  ( j < nyng ) )
     2118             j = j + 1
     2119          ENDDO
     2120          jfuo(jj) = MIN( MAX( j, nysg ), nyng )
     2121          jstart = jflv(jj)
     2122       ENDDO
     2123
     2124!
     2125!--    k-indices of w for each kk-index value
     2126       kstart  = 0
     2127       kflw(0) = 0
     2128       kfuw(0) = 0
     2129       DO  kk = 1, kctw
     2130          k = kstart
     2131          DO  WHILE ( ( zw(k) < cg%zw(kk) - 0.5_wp * cg%dzw(kk) )  .AND.       &
     2132                      ( k < nzt ) )
     2133             k = k + 1
     2134          ENDDO
     2135          kflw(kk) = MIN( MAX( k, 1 ), nzt + 1 )
     2136          DO  WHILE ( ( zw(k) < cg%zw(kk) + 0.5_wp * cg%dzw(kk+1) )  .AND.     &
     2137                      ( k < nzt ) )
     2138             k = k + 1
     2139          ENDDO
     2140          kfuw(kk) = MIN( MAX( k, 1 ), nzt + 1 )
     2141          kstart = kflw(kk)
     2142       ENDDO
     2143
     2144!
     2145!--    k-indices of others for each kk-index value
     2146       kstart  = 0
     2147       kflo(0) = 0
     2148       kfuo(0) = 0
     2149       DO  kk = 1, kctu
     2150          k = kstart
     2151          DO  WHILE ( ( zu(k) < cg%zu(kk) - 0.5_wp * cg%dzu(kk) )  .AND.       &
     2152                      ( k < nzt ) )
     2153             k = k + 1
     2154          ENDDO
     2155          kflo(kk) = MIN( MAX( k, 1 ), nzt + 1 )
     2156          DO  WHILE ( ( zu(k) < cg%zu(kk) + 0.5_wp * cg%dzu(kk+1) )  .AND.     &
     2157                      ( k < nzt ) )
     2158             k = k + 1
     2159          ENDDO
     2160          kfuo(kk) = MIN( MAX( k-1, 1 ), nzt + 1 )
     2161          kstart = kflo(kk)
     2162       ENDDO
    20702163     
    20712164!
    2072 !--   Spatial under-relaxation coefficients
    2073       ALLOCATE( frax(icl:icr) )
    2074 
    2075       DO  ii = icl, icr
    2076          IF ( nest_bound_l ) THEN
    2077             xi   = ( MAX( 0.0_wp, ( cg%coord_x(ii) - lower_left_coord_x ) ) / anterp_relax_length_l )**4
    2078          ELSEIF ( nest_bound_r ) THEN
    2079             xi   = ( MAX( 0.0_wp, ( lower_left_coord_x + ( nx + 1 ) * dx - cg%coord_x(ii) ) ) / anterp_relax_length_r )**4
    2080          ELSE
    2081             xi   = 999999.9_wp
    2082          ENDIF
    2083          frax(ii)  = xi / ( 1.0_wp + xi )
    2084       ENDDO
    2085 
    2086       ALLOCATE( fray(jcs:jcn) )
    2087 
    2088       DO  jj = jcs, jcn
    2089          IF ( nest_bound_s ) THEN           
    2090             eta  = ( MAX( 0.0_wp, ( cg%coord_y(jj) - lower_left_coord_y ) ) / anterp_relax_length_s )**4
    2091          ELSEIF ( nest_bound_n ) THEN
    2092             eta  = ( MAX( 0.0_wp, ( lower_left_coord_y + ( ny + 1 ) * dy - cg%coord_y(jj)) ) / anterp_relax_length_n )**4
    2093          ELSE
    2094             eta  = 999999.9_wp
    2095          ENDIF
    2096          fray(jj)  = eta / ( 1.0_wp + eta )
    2097       ENDDO
     2165!--    Spatial under-relaxation coefficients
     2166       ALLOCATE( frax(icl:icr) )
     2167
     2168       DO  ii = icl, icr
     2169          IF ( nest_bound_l )  THEN
     2170             xi = ( MAX( 0.0_wp, ( cg%coord_x(ii) - lower_left_coord_x ) ) /   &
     2171                    anterp_relax_length_l )**4
     2172          ELSEIF ( nest_bound_r )  THEN
     2173             xi = ( MAX( 0.0_wp, ( lower_left_coord_x + ( nx + 1 ) * dx -      &
     2174                                   cg%coord_x(ii) ) ) /                        &
     2175                    anterp_relax_length_r )**4
     2176          ELSE
     2177             xi = 999999.9_wp
     2178          ENDIF
     2179          frax(ii) = xi / ( 1.0_wp + xi )
     2180       ENDDO
     2181
     2182       ALLOCATE( fray(jcs:jcn) )
     2183
     2184       DO  jj = jcs, jcn
     2185          IF ( nest_bound_s )  THEN
     2186             eta = ( MAX( 0.0_wp, ( cg%coord_y(jj) - lower_left_coord_y ) ) /  &
     2187                     anterp_relax_length_s )**4
     2188          ELSEIF ( nest_bound_n )  THEN
     2189             eta = ( MAX( 0.0_wp, ( lower_left_coord_y + ( ny + 1 ) * dy -     &
     2190                                    cg%coord_y(jj)) ) /                        &
     2191                     anterp_relax_length_n )**4
     2192          ELSE
     2193             eta = 999999.9_wp
     2194          ENDIF
     2195          fray(jj) = eta / ( 1.0_wp + eta )
     2196       ENDDO
    20982197     
    2099       ALLOCATE( fraz(0:kctu) )
    2100       DO  kk = 0, kctu       
    2101          zeta = ( ( zu(nzt) - cg%zu(kk) ) / anterp_relax_length_t )**4       
    2102          fraz(kk) = zeta / ( 1.0_wp + zeta )
    2103       ENDDO
    2104 
    2105    END SUBROUTINE pmci_init_anterp_tophat
    2106 
    2107 
    2108 
    2109    SUBROUTINE pmci_init_tkefactor
    2110 
    2111 !
    2112 !--   Computes the scaling factor for the SGS TKE from coarse grid to be used
    2113 !--   as BC for the fine grid. Based on the Kolmogorov energy spectrum
    2114 !--   for the inertial subrange and assumption of sharp cut-off of the resolved
    2115 !--   energy spectrum. Near the surface, the reduction of TKE is made
    2116 !--   smaller than further away from the surface.
    2117 !
    2118 !                      Antti Hellsten 4.3.2015
    2119 !
    2120 !--   Extended for non-flat topography and variable dz.
    2121 !
    2122 !                      Antti Hellsten 26.3.2015
    2123 !
    2124 !--   The current near-wall adaptation can be replaced by a new one which
    2125 !--   uses a step function [0,1] based on the logc-arrays. AH 30.12.2015
    2126       IMPLICIT NONE
    2127       REAL(wp), PARAMETER  ::  cfw = 0.2_wp          !:
    2128       REAL(wp), PARAMETER  ::  c_tkef = 0.6_wp       !:
    2129       REAL(wp)             ::  fw                    !:
    2130       REAL(wp), PARAMETER  ::  fw0 = 0.9_wp          !:
    2131       REAL(wp)             ::  glsf                  !:
    2132       REAL(wp)             ::  glsc                  !:
    2133       REAL(wp)             ::  height                !:
    2134       REAL(wp), PARAMETER  ::  p13 = 1.0_wp/3.0_wp   !:
    2135       REAL(wp), PARAMETER  ::  p23 = 2.0_wp/3.0_wp   !:
    2136       INTEGER(iwp)         ::  k                     !:
    2137       INTEGER(iwp)         ::  kc                    !:
     2198       ALLOCATE( fraz(0:kctu) )
     2199       DO  kk = 0, kctu
     2200          zeta = ( ( zu(nzt) - cg%zu(kk) ) / anterp_relax_length_t )**4
     2201          fraz(kk) = zeta / ( 1.0_wp + zeta )
     2202       ENDDO
     2203
     2204    END SUBROUTINE pmci_init_anterp_tophat
     2205
     2206
     2207
     2208    SUBROUTINE pmci_init_tkefactor
     2209
     2210!
     2211!--    Computes the scaling factor for the SGS TKE from coarse grid to be used
     2212!--    as BC for the fine grid. Based on the Kolmogorov energy spectrum
     2213!--    for the inertial subrange and assumption of sharp cut-off of the resolved
     2214!--    energy spectrum. Near the surface, the reduction of TKE is made
     2215!--    smaller than further away from the surface.
     2216
     2217       IMPLICIT NONE
     2218       REAL(wp), PARAMETER ::  cfw = 0.2_wp          !:
     2219       REAL(wp), PARAMETER ::  c_tkef = 0.6_wp       !:
     2220       REAL(wp)            ::  fw                    !:
     2221       REAL(wp), PARAMETER ::  fw0 = 0.9_wp          !:
     2222       REAL(wp)            ::  glsf                  !:
     2223       REAL(wp)            ::  glsc                  !:
     2224       REAL(wp)            ::  height                !:
     2225       REAL(wp), PARAMETER ::  p13 = 1.0_wp/3.0_wp   !:
     2226       REAL(wp), PARAMETER ::  p23 = 2.0_wp/3.0_wp   !:
     2227       INTEGER(iwp)        ::  k                     !:
     2228       INTEGER(iwp)        ::  kc                    !:
    21382229       
    21392230
    2140       IF ( nest_bound_l )  THEN
    2141          ALLOCATE( tkefactor_l(nzb:nzt+1,nysg:nyng) )   
    2142          tkefactor_l = 0.0_wp
    2143          i = nxl - 1
    2144          DO  j = nysg, nyng
    2145             DO  k = nzb_s_inner(j,i) + 1, nzt
    2146                kc     = kco(k+1)
    2147                glsf   = ( dx * dy * dzu(k) )**p13
    2148                glsc   = ( cg%dx * cg%dy *cg%dzu(kc) )**p13
    2149                height = zu(k) - zu(nzb_s_inner(j,i))
    2150                fw     = EXP( -cfw * height / glsf )
    2151                tkefactor_l(k,j) = c_tkef * ( fw0 * fw + ( 1.0_wp - fw ) * ( glsf / glsc )**p23 )
    2152             ENDDO
    2153             tkefactor_l(nzb_s_inner(j,i),j) = c_tkef * fw0
    2154          ENDDO
    2155       ENDIF
    2156 
    2157       IF ( nest_bound_r )  THEN
    2158          ALLOCATE( tkefactor_r(nzb:nzt+1,nysg:nyng) )   
    2159          tkefactor_r = 0.0_wp
    2160          i = nxr + 1
    2161          DO  j = nysg, nyng
    2162             DO  k = nzb_s_inner(j,i) + 1, nzt
    2163                kc     = kco(k+1)
    2164                glsf   = ( dx * dy * dzu(k) )**p13
    2165                glsc   = ( cg%dx * cg%dy * cg%dzu(kc) )**p13
    2166                height = zu(k) - zu(nzb_s_inner(j,i))
    2167                fw     = EXP( -cfw * height / glsf )
    2168                tkefactor_r(k,j) = c_tkef * (fw0 * fw + ( 1.0_wp - fw ) * ( glsf / glsc )**p23 )
    2169             ENDDO
    2170             tkefactor_r(nzb_s_inner(j,i),j) = c_tkef * fw0
    2171          ENDDO
    2172       ENDIF
     2231       IF ( nest_bound_l )  THEN
     2232          ALLOCATE( tkefactor_l(nzb:nzt+1,nysg:nyng) )
     2233          tkefactor_l = 0.0_wp
     2234          i = nxl - 1
     2235          DO  j = nysg, nyng
     2236             DO  k = nzb_s_inner(j,i) + 1, nzt
     2237                kc     = kco(k+1)
     2238                glsf   = ( dx * dy * dzu(k) )**p13
     2239                glsc   = ( cg%dx * cg%dy *cg%dzu(kc) )**p13
     2240                height = zu(k) - zu(nzb_s_inner(j,i))
     2241                fw     = EXP( -cfw * height / glsf )
     2242                tkefactor_l(k,j) = c_tkef * ( fw0 * fw + ( 1.0_wp - fw ) *     &
     2243                                              ( glsf / glsc )**p23 )
     2244             ENDDO
     2245             tkefactor_l(nzb_s_inner(j,i),j) = c_tkef * fw0
     2246          ENDDO
     2247       ENDIF
     2248
     2249       IF ( nest_bound_r )  THEN
     2250          ALLOCATE( tkefactor_r(nzb:nzt+1,nysg:nyng) )
     2251          tkefactor_r = 0.0_wp
     2252          i = nxr + 1
     2253          DO  j = nysg, nyng
     2254             DO  k = nzb_s_inner(j,i) + 1, nzt
     2255                kc     = kco(k+1)
     2256                glsf   = ( dx * dy * dzu(k) )**p13
     2257                glsc   = ( cg%dx * cg%dy * cg%dzu(kc) )**p13
     2258                height = zu(k) - zu(nzb_s_inner(j,i))
     2259                fw     = EXP( -cfw * height / glsf )
     2260                tkefactor_r(k,j) = c_tkef * ( fw0 * fw + ( 1.0_wp - fw ) *     &
     2261                                              ( glsf / glsc )**p23 )
     2262             ENDDO
     2263             tkefactor_r(nzb_s_inner(j,i),j) = c_tkef * fw0
     2264          ENDDO
     2265       ENDIF
    21732266
    21742267      IF ( nest_bound_s )  THEN
    2175          ALLOCATE( tkefactor_s(nzb:nzt+1,nxlg:nxrg) )
    2176          tkefactor_s = 0.0_wp
    2177          j = nys - 1
    2178          DO  i = nxlg, nxrg
    2179             DO  k = nzb_s_inner(j,i) + 1, nzt
    2180                kc     = kco(k+1)
    2181                glsf   = ( dx * dy * dzu(k) )**p13
    2182                glsc   = ( cg%dx * cg%dy * cg%dzu(kc) ) ** p13
    2183                height = zu(k) - zu(nzb_s_inner(j,i))
    2184                fw     = EXP( -cfw*height / glsf )
    2185                tkefactor_s(k,i) = c_tkef * ( fw0 * fw + ( 1.0_wp - fw ) * ( glsf / glsc )**p23 )
    2186             ENDDO
    2187             tkefactor_s(nzb_s_inner(j,i),i) = c_tkef * fw0
    2188          ENDDO
    2189       ENDIF
    2190 
    2191       IF ( nest_bound_n )  THEN
    2192          ALLOCATE( tkefactor_n(nzb:nzt+1,nxlg:nxrg) )
    2193          tkefactor_n = 0.0_wp   
    2194          j = nyn + 1
    2195          DO  i = nxlg, nxrg
    2196             DO  k = nzb_s_inner(j,i)+1, nzt
    2197                kc     = kco(k+1)
    2198                glsf   = ( dx * dy * dzu(k) )**p13
    2199                glsc   = ( cg%dx * cg%dy * cg%dzu(kc) )**p13
    2200                height = zu(k) - zu(nzb_s_inner(j,i))
    2201                fw     = EXP( -cfw * height / glsf )
    2202                tkefactor_n(k,i) = c_tkef * ( fw0 * fw + ( 1.0_wp - fw ) * ( glsf / glsc )**p23 )
    2203             ENDDO
    2204             tkefactor_n(nzb_s_inner(j,i),i) = c_tkef * fw0
    2205          ENDDO
    2206       ENDIF
    2207 
    2208       ALLOCATE( tkefactor_t(nysg:nyng,nxlg:nxrg) )
    2209       k = nzt
    2210       DO  i = nxlg, nxrg
    2211          DO  j = nysg, nyng
    2212             kc     = kco(k+1)
    2213             glsf   = ( dx * dy * dzu(k) )**p13
    2214             glsc   = ( cg%dx * cg%dy * cg%dzu(kc) )**p13
    2215             height = zu(k) - zu(nzb_s_inner(j,i))
    2216             fw     = EXP( -cfw * height / glsf )
    2217             tkefactor_t(j,i) = c_tkef * ( fw0 * fw + ( 1.0_wp - fw ) * ( glsf / glsc )**p23 )
    2218          ENDDO
    2219       ENDDO
     2268          ALLOCATE( tkefactor_s(nzb:nzt+1,nxlg:nxrg) )
     2269          tkefactor_s = 0.0_wp
     2270          j = nys - 1
     2271          DO  i = nxlg, nxrg
     2272             DO  k = nzb_s_inner(j,i) + 1, nzt
     2273                kc     = kco(k+1)
     2274                glsf   = ( dx * dy * dzu(k) )**p13
     2275                glsc   = ( cg%dx * cg%dy * cg%dzu(kc) ) ** p13
     2276                height = zu(k) - zu(nzb_s_inner(j,i))
     2277                fw     = EXP( -cfw*height / glsf )
     2278                tkefactor_s(k,i) = c_tkef * ( fw0 * fw + ( 1.0_wp - fw ) *     &
     2279                                              ( glsf / glsc )**p23 )
     2280             ENDDO
     2281             tkefactor_s(nzb_s_inner(j,i),i) = c_tkef * fw0
     2282          ENDDO
     2283       ENDIF
     2284
     2285       IF ( nest_bound_n )  THEN
     2286          ALLOCATE( tkefactor_n(nzb:nzt+1,nxlg:nxrg) )
     2287          tkefactor_n = 0.0_wp
     2288          j = nyn + 1
     2289          DO  i = nxlg, nxrg
     2290             DO  k = nzb_s_inner(j,i)+1, nzt
     2291                kc     = kco(k+1)
     2292                glsf   = ( dx * dy * dzu(k) )**p13
     2293                glsc   = ( cg%dx * cg%dy * cg%dzu(kc) )**p13
     2294                height = zu(k) - zu(nzb_s_inner(j,i))
     2295                fw     = EXP( -cfw * height / glsf )
     2296                tkefactor_n(k,i) = c_tkef * ( fw0 * fw + ( 1.0_wp - fw ) *     &
     2297                                              ( glsf / glsc )**p23 )
     2298             ENDDO
     2299             tkefactor_n(nzb_s_inner(j,i),i) = c_tkef * fw0
     2300          ENDDO
     2301       ENDIF
     2302
     2303       ALLOCATE( tkefactor_t(nysg:nyng,nxlg:nxrg) )
     2304       k = nzt
     2305       DO  i = nxlg, nxrg
     2306          DO  j = nysg, nyng
     2307             kc     = kco(k+1)
     2308             glsf   = ( dx * dy * dzu(k) )**p13
     2309             glsc   = ( cg%dx * cg%dy * cg%dzu(kc) )**p13
     2310             height = zu(k) - zu(nzb_s_inner(j,i))
     2311             fw     = EXP( -cfw * height / glsf )
     2312             tkefactor_t(j,i) = c_tkef * ( fw0 * fw + ( 1.0_wp - fw ) *        &
     2313                                           ( glsf / glsc )**p23 )
     2314          ENDDO
     2315       ENDDO
    22202316     
    22212317    END SUBROUTINE pmci_init_tkefactor
     
    23662462
    23672463!
    2368 !-- List of array names, which can be coupled.
    2369 !-- TO_DO: Antti: what is the meaning of the next line?
    2370 !-- AH: Note that the k-range of the *c arrays is changed from 1:nz to 0:nz+1.
     2464!-- List of array names, which can be coupled
    23712465    IF ( TRIM( name ) == "u" )  THEN
    23722466       IF ( .NOT. ALLOCATED( uc ) )  ALLOCATE( uc(0:nzc+1, js:je, is:ie) )
     
    24062500       ELSE
    24072501!
    2408 !--       Avoid others to continue
     2502!--       Prevent others from continuing
    24092503          CALL MPI_BARRIER( comm2d, ierr )
    24102504       ENDIF
     
    26502744
    26512745
     2746 SUBROUTINE pmci_check_setting_mismatches
     2747!
     2748!-- Check for mismatches between settings of master and client variables
     2749!-- (e.g., all clients have to follow the end_time settings of the root model).
     2750!-- The root model overwrites variables in the other models, so these variables
     2751!-- only need to be set once in file PARIN.
     2752
     2753#if defined( __parallel )
     2754
     2755    USE control_parameters,                                                    &
     2756        ONLY:  dt_restart, end_time, message_string, restart_time, time_restart
     2757
     2758    IMPLICIT NONE
     2759
     2760    INTEGER ::  ierr
     2761
     2762    REAL(wp) ::  dt_restart_root
     2763    REAL(wp) ::  end_time_root
     2764    REAL(wp) ::  restart_time_root
     2765    REAL(wp) ::  time_restart_root
     2766
     2767!
     2768!-- Check the time to be simulated.
     2769!-- Here, and in the following, the root process communicates the respective
     2770!-- variable to all others, and its value will then be compared with the local
     2771!-- values.
     2772    IF ( pmc_is_rootmodel() )  end_time_root = end_time
     2773    CALL MPI_BCAST( end_time_root, 1, MPI_REAL, 0, comm_world_nesting, ierr )
     2774
     2775    IF ( .NOT. pmc_is_rootmodel() )  THEN
     2776       IF ( end_time /= end_time_root )  THEN
     2777          WRITE( message_string, * )  'mismatch between root model and ',      &
     2778               'client settings &   end_time(root) = ', end_time_root,         &
     2779               ' &   end_time(client) = ', end_time, ' & client value is set', &
     2780               ' to root value'
     2781          CALL message( 'pmci_check_setting_mismatches', 'PA0419', 0, 1, 0, 6, &
     2782                        0 )
     2783          end_time = end_time_root
     2784       ENDIF
     2785    ENDIF
     2786
     2787!
     2788!-- Same for restart time
     2789    IF ( pmc_is_rootmodel() )  restart_time_root = restart_time
     2790    CALL MPI_BCAST( restart_time_root, 1, MPI_REAL, 0, comm_world_nesting, ierr )
     2791
     2792    IF ( .NOT. pmc_is_rootmodel() )  THEN
     2793       IF ( restart_time /= restart_time_root )  THEN
     2794          WRITE( message_string, * )  'mismatch between root model and ',      &
     2795               'client settings &   restart_time(root) = ', restart_time_root, &
     2796               ' &   restart_time(client) = ', restart_time, ' & client ',     &
     2797               'value is set to root value'
     2798          CALL message( 'pmci_check_setting_mismatches', 'PA0419', 0, 1, 0, 6, &
     2799                        0 )
     2800          restart_time = restart_time_root
     2801       ENDIF
     2802    ENDIF
     2803
     2804!
     2805!-- Same for dt_restart
     2806    IF ( pmc_is_rootmodel() )  dt_restart_root = dt_restart
     2807    CALL MPI_BCAST( dt_restart_root, 1, MPI_REAL, 0, comm_world_nesting, ierr )
     2808
     2809    IF ( .NOT. pmc_is_rootmodel() )  THEN
     2810       IF ( dt_restart /= dt_restart_root )  THEN
     2811          WRITE( message_string, * )  'mismatch between root model and ',      &
     2812               'client settings &   dt_restart(root) = ', dt_restart_root,     &
     2813               ' &   dt_restart(client) = ', dt_restart, ' & client ',         &
     2814               'value is set to root value'
     2815          CALL message( 'pmci_check_setting_mismatches', 'PA0419', 0, 1, 0, 6, &
     2816                        0 )
     2817          dt_restart = dt_restart_root
     2818       ENDIF
     2819    ENDIF
     2820
     2821!
     2822!-- Same for time_restart
     2823    IF ( pmc_is_rootmodel() )  time_restart_root = time_restart
     2824    CALL MPI_BCAST( time_restart_root, 1, MPI_REAL, 0, comm_world_nesting, ierr )
     2825
     2826    IF ( .NOT. pmc_is_rootmodel() )  THEN
     2827       IF ( time_restart /= time_restart_root )  THEN
     2828          WRITE( message_string, * )  'mismatch between root model and ',      &
     2829               'client settings &   time_restart(root) = ', time_restart_root, &
     2830               ' &   time_restart(client) = ', time_restart, ' & client ',     &
     2831               'value is set to root value'
     2832          CALL message( 'pmci_check_setting_mismatches', 'PA0419', 0, 1, 0, 6, &
     2833                        0 )
     2834          time_restart = time_restart_root
     2835       ENDIF
     2836    ENDIF
     2837
     2838#endif
     2839
     2840 END SUBROUTINE pmci_check_setting_mismatches
     2841
     2842
     2843
    26522844 SUBROUTINE pmci_ensure_nest_mass_conservation
    26532845
     
    27872979    REAL(wp), DIMENSION(1) ::  dtl         !:
    27882980
     2981
     2982    CALL cpu_log( log_point_s(70), 'pmc sync', 'start' )
     2983
    27892984!
    27902985!-- First find the smallest native time step of all the clients of the current
     
    28133008    ENDDO
    28143009
     3010    CALL cpu_log( log_point_s(70), 'pmc sync', 'stop' )
     3011
    28153012#endif
    28163013 END SUBROUTINE pmci_server_synchronize
     
    28363033    dtl(1) = dt_3d
    28373034    IF ( cpl_id > 1 )  THEN
     3035
     3036       CALL cpu_log( log_point_s(70), 'pmc sync', 'start' )
     3037
    28383038       IF ( myid==0 )  THEN
    28393039          CALL pmc_send_to_server( dtl, SIZE( dtl ), 0, 101, ierr )
     
    28453045!--    Broadcast the unified time step to all server processes
    28463046       CALL MPI_BCAST( dt_3d, 1, MPI_REAL, 0, comm2d, ierr )
     3047
     3048       CALL cpu_log( log_point_s(70), 'pmc sync', 'stop' )
     3049
    28473050    ENDIF
    28483051
     
    28533056
    28543057 SUBROUTINE pmci_set_swaplevel( swaplevel )
     3058!
     3059!-- After each Runge-Kutta sub-timestep, alternately set buffer one or buffer
     3060!-- two active
    28553061
    28563062    IMPLICIT NONE
     
    28623068    INTEGER(iwp)            ::  m          !:
    28633069
    2864 !
    2865 !-- After each timestep, alternately set buffer one or buffer two active
    28663070    DO  m = 1, SIZE( pmc_server_for_client )-1
    28673071       client_id = pmc_server_for_client(m)
     
    28703074
    28713075 END SUBROUTINE pmci_set_swaplevel
     3076
     3077
     3078
     3079 SUBROUTINE pmci_datatrans( local_nesting_mode )
     3080!
     3081!-- Althoug nesting_mode is a variable of this model, pass it as an argument to
     3082!-- allow for example to force one-way initialization phase
     3083
     3084    IMPLICIT NONE
     3085
     3086    INTEGER(iwp)           ::  ierr   !:
     3087    INTEGER(iwp)           ::  istat  !:
     3088
     3089    CHARACTER(LEN=*),INTENT(IN) ::  local_nesting_mode
     3090
     3091    IF ( local_nesting_mode == 'one-way' )  THEN
     3092
     3093       CALL pmci_client_datatrans( server_to_client )
     3094       CALL pmci_server_datatrans( server_to_client )
     3095
     3096    ELSE
     3097
     3098       IF( nesting_datatransfer_mode == 'cascade' )  THEN
     3099
     3100          CALL pmci_client_datatrans( server_to_client )
     3101          CALL pmci_server_datatrans( server_to_client )
     3102
     3103          CALL pmci_server_datatrans( client_to_server )
     3104          CALL pmci_client_datatrans( client_to_server )
     3105
     3106       ELSEIF( nesting_datatransfer_mode == 'overlap')  THEN
     3107
     3108          CALL pmci_server_datatrans( server_to_client )
     3109          CALL pmci_client_datatrans( server_to_client )
     3110
     3111          CALL pmci_client_datatrans( client_to_server )
     3112          CALL pmci_server_datatrans( client_to_server )
     3113
     3114       ELSEIF( TRIM( nesting_datatransfer_mode ) == 'mixed' )  THEN
     3115
     3116          CALL pmci_server_datatrans( server_to_client )
     3117          CALL pmci_client_datatrans( server_to_client )
     3118
     3119          CALL pmci_server_datatrans( client_to_server )
     3120          CALL pmci_client_datatrans( client_to_server )
     3121
     3122       ENDIF
     3123
     3124    ENDIF
     3125
     3126 END SUBROUTINE pmci_datatrans
     3127
    28723128
    28733129
     
    28903146    REAL(wp), DIMENSION(1) ::  dtl         !:
    28913147
    2892 !
    2893 !-- First find the smallest native time step of all the clients of the current
    2894 !-- server.
    2895     dtl(1) = 999999.9_wp
    2896     DO  m = 1, SIZE( pmc_server_for_client )-1
    2897        client_id = pmc_server_for_client(m)
    2898        IF ( myid==0 )  THEN
    2899           CALL pmc_recv_from_client( client_id, dtc, SIZE( dtc ), 0, 101, ierr )
    2900           dtl(1) = MIN( dtl(1), dtc(1) )
    2901           dt_3d  = dtl(1)
    2902        ENDIF
    2903     ENDDO
    2904 
    2905 !
    2906 !-- Broadcast the unified time step to all server processes
    2907     CALL MPI_BCAST( dt_3d, 1, MPI_REAL, 0, comm2d, ierr )
    29083148
    29093149    DO  m = 1, SIZE( PMC_Server_for_Client )-1
    29103150       client_id = PMC_Server_for_Client(m)
    2911        CALL cpu_log( log_point_s(70), 'pmc sync', 'start' )
    2912 
    2913 !
    2914 !--    Send the new time step to all the clients of the current server
    2915        IF ( myid == 0 )  THEN
    2916           CALL pmc_send_to_client( client_id, dtl, SIZE( dtl ), 0, 102, ierr )
    2917        ENDIF
    2918        CALL cpu_log( log_point_s(70), 'pmc sync', 'stop' )
    29193151       
    29203152       IF ( direction == server_to_client )  THEN
    29213153          CALL cpu_log( log_point_s(71), 'pmc server send', 'start' )
    2922           CALL pmc_s_fillbuffer( client_id, waittime=waittime )
     3154          CALL pmc_s_fillbuffer( client_id )
    29233155          CALL cpu_log( log_point_s(71), 'pmc server send', 'stop' )
    29243156       ELSE
     
    29363168!
    29373169!--          Inside buildings/topography reset velocities and TKE back to zero.
    2938 !--          TO_DO: at least temperature should be included here immediately
    29393170!--          Other scalars (pt, q, s, km, kh, p, sa, ...) are ignored at
    29403171!--          present, maybe revise later.
    29413172             DO   i = nxlg, nxrg
    29423173                DO   j = nysg, nyng
    2943                    u(nzb:nzb_u_inner(j,i),j,i) = 0.0_wp
    2944                    v(nzb:nzb_v_inner(j,i),j,i) = 0.0_wp
    2945                    w(nzb:nzb_w_inner(j,i),j,i) = 0.0_wp
    2946                    e(nzb:nzb_s_inner(j,i),j,i) = 0.0_wp
     3174                   u(nzb:nzb_u_inner(j,i),j,i)  = 0.0_wp
     3175                   v(nzb:nzb_v_inner(j,i),j,i)  = 0.0_wp
     3176                   w(nzb:nzb_w_inner(j,i),j,i)  = 0.0_wp
     3177                   e(nzb:nzb_s_inner(j,i),j,i)  = 0.0_wp
     3178                   pt(nzb:nzb_s_inner(j,i),j,i) = 0.0_wp
     3179                   IF ( humidity  .OR.  passive_scalar )  THEN
     3180                      q(nzb:nzb_s_inner(j,i),j,i) = 0.0_wp
     3181                   ENDIF
    29473182                ENDDO
    29483183             ENDDO
     
    29693204    INTEGER(iwp) ::  jcn         !:
    29703205   
    2971     REAL(wp) ::  waittime    !:
    2972 
    29733206    REAL(wp), DIMENSION(1) ::  dtl         !:
    29743207    REAL(wp), DIMENSION(1) ::  dts         !:
     
    29773210    dtl = dt_3d
    29783211    IF ( cpl_id > 1 )  THEN
    2979        CALL cpu_log( log_point_s(70), 'pmc sync', 'start' )
    2980        IF ( myid==0 )  THEN
    2981           CALL pmc_send_to_server( dtl, SIZE( dtl ), 0, 101, ierr )
    2982           CALL pmc_recv_from_server( dts, SIZE( dts ), 0, 102, ierr )
    2983           dt_3d = dts(1)
    2984        ENDIF
    2985 
    2986 !
    2987 !--    Broadcast the unified time step to all server processes.
    2988        CALL MPI_BCAST( dt_3d, 1, MPI_REAL, 0, comm2d, ierr )
    2989        CALL cpu_log( log_point_s(70), 'pmc sync', 'stop' )
    2990 
    29913212!
    29923213!--    Client domain boundaries in the server indice space.
     
    29993220
    30003221          CALL cpu_log( log_point_s(73), 'pmc client recv', 'start' )
    3001           CALL pmc_c_getbuffer( WaitTime = WaitTime )
     3222          CALL pmc_c_getbuffer( )
    30023223          CALL cpu_log( log_point_s(73), 'pmc client recv', 'stop' )
    30033224
     3225          CALL cpu_log( log_point_s(75), 'pmc interpolation', 'start' )
    30043226          CALL pmci_interpolation
     3227          CALL cpu_log( log_point_s(75), 'pmc interpolation', 'stop' )
    30053228
    30063229       ELSE
    30073230!
    3008 !--       direction == server_to_client
     3231!--       direction == client_to_server
     3232          CALL cpu_log( log_point_s(76), 'pmc anterpolation', 'start' )
    30093233          CALL pmci_anterpolation
     3234          CALL cpu_log( log_point_s(76), 'pmc anterpolation', 'stop' )
    30103235
    30113236          CALL cpu_log( log_point_s(74), 'pmc client send', 'start' )
    3012           CALL pmc_c_putbuffer( WaitTime = WaitTime )
     3237          CALL pmc_c_putbuffer( )
    30133238          CALL cpu_log( log_point_s(74), 'pmc client send', 'stop' )
    30143239
     
    32343459
    32353460      CALL pmci_anterp_tophat( u,  uc,  kctu, iflu, ifuu, jflo, jfuo, kflo,    &
    3236                                kfuo, nzb_u_inner, 'u' )
     3461                               kfuo, 'u' )
    32373462      CALL pmci_anterp_tophat( v,  vc,  kctu, iflo, ifuo, jflv, jfuv, kflo,    &
    3238                                kfuo, nzb_v_inner, 'v' )
     3463                               kfuo, 'v' )
    32393464      CALL pmci_anterp_tophat( w,  wc,  kctw, iflo, ifuo, jflo, jfuo, kflw,    &
    3240                                kfuw, nzb_w_inner, 'w' )
     3465                               kfuw, 'w' )
    32413466      CALL pmci_anterp_tophat( pt, ptc, kctu, iflo, ifuo, jflo, jfuo, kflo,    &
    3242                                kfuo, nzb_s_inner, 's' )
     3467                               kfuo, 's' )
    32433468      IF ( humidity  .OR.  passive_scalar )  THEN
    32443469         CALL pmci_anterp_tophat( q, qc, kctu, iflo, ifuo, jflo, jfuo, kflo,   &
    3245                                   kfuo, nzb_s_inner, 's' )
     3470                                  kfuo, 's' )
    32463471      ENDIF
    32473472
     
    32563481!--   Interpolation of ghost-node values used as the client-domain boundary
    32573482!--   conditions. This subroutine handles the left and right boundaries. It is
    3258 !--   based on trilinear interpolation. Constant dz is still assumed.
    3259 !--   TO_DO:  constant dz is an important restriction and should be checked
    3260 !--           somewhere in order to let users not run into this trap
     3483!--   based on trilinear interpolation.
     3484
    32613485      IMPLICIT NONE
    32623486
    3263 !--   TO_DO: wrap long lines in this and the remaining interp_tril routines
    3264       REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg), INTENT(INOUT)          ::  f            !:
    3265       REAL(wp), DIMENSION(0:cg%nz+1,jcs:jcn,icl:icr), INTENT(IN)                 ::  fc           !:
    3266       REAL(wp), DIMENSION(nzb:nzt_topo_nestbc,nys:nyn,1:2,0:ncorr-1), INTENT(IN) ::  logc_ratio   !:
    3267       REAL(wp), DIMENSION(nxlg:nxrg), INTENT(IN)                                 ::  r1x          !:
    3268       REAL(wp), DIMENSION(nxlg:nxrg), INTENT(IN)                                 ::  r2x          !:
    3269       REAL(wp), DIMENSION(nysg:nyng), INTENT(IN)                                 ::  r1y          !:
    3270       REAL(wp), DIMENSION(nysg:nyng), INTENT(IN)                                 ::  r2y          !:
    3271       REAL(wp), DIMENSION(nzb:nzt+1), INTENT(IN)                                 ::  r1z          !:
    3272       REAL(wp), DIMENSION(nzb:nzt+1), INTENT(IN)                                 ::  r2z          !:
     3487      REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg),                      &
     3488                                      INTENT(INOUT) ::  f       !:
     3489      REAL(wp), DIMENSION(0:cg%nz+1,jcs:jcn,icl:icr),                          &
     3490                                      INTENT(IN)    ::  fc      !:
     3491      REAL(wp), DIMENSION(nzb:nzt_topo_nestbc,nys:nyn,1:2,0:ncorr-1),          &
     3492                                      INTENT(IN)    ::  logc_ratio   !:
     3493      REAL(wp), DIMENSION(nxlg:nxrg), INTENT(IN)    ::  r1x     !:
     3494      REAL(wp), DIMENSION(nxlg:nxrg), INTENT(IN)    ::  r2x     !:
     3495      REAL(wp), DIMENSION(nysg:nyng), INTENT(IN)    ::  r1y     !:
     3496      REAL(wp), DIMENSION(nysg:nyng), INTENT(IN)    ::  r2y     !:
     3497      REAL(wp), DIMENSION(nzb:nzt+1), INTENT(IN)    ::  r1z     !:
     3498      REAL(wp), DIMENSION(nzb:nzt+1), INTENT(IN)    ::  r2z     !:
    32733499     
    3274       INTEGER(iwp), DIMENSION(nxlg:nxrg), INTENT(IN)                             ::  ic           !:
    3275       INTEGER(iwp), DIMENSION(nysg:nyng), INTENT(IN)                             ::  jc           !:
    3276       INTEGER(iwp), DIMENSION(nysg:nyng,nxlg:nxrg), INTENT(IN)                   ::  kb           !:
    3277       INTEGER(iwp), DIMENSION(nzb:nzt+1), INTENT(IN)                             ::  kc           !:
    3278       INTEGER(iwp), DIMENSION(nzb:nzt_topo_nestbc,nys:nyn,1:2), INTENT(IN)       ::  logc         !:
    3279       INTEGER(iwp)                                                               ::  nzt_topo_nestbc   !:
     3500      INTEGER(iwp), DIMENSION(nxlg:nxrg), INTENT(IN)           ::  ic     !:
     3501      INTEGER(iwp), DIMENSION(nysg:nyng), INTENT(IN)           ::  jc     !:
     3502      INTEGER(iwp), DIMENSION(nysg:nyng,nxlg:nxrg), INTENT(IN) ::  kb     !:
     3503      INTEGER(iwp), DIMENSION(nzb:nzt+1), INTENT(IN)           ::  kc     !:
     3504      INTEGER(iwp), DIMENSION(nzb:nzt_topo_nestbc,nys:nyn,1:2),                &
     3505                                          INTENT(IN)           ::  logc   !:
     3506      INTEGER(iwp) ::  nzt_topo_nestbc   !:
    32803507
    32813508      CHARACTER(LEN=1),INTENT(IN) ::  edge   !:
     
    33123539     
    33133540!
    3314 !--   Check which edge is to be handled: left or right. Note the assumption that the same PE never
    3315 !--   holds both left and right nest boundaries. Should this be changed?
     3541!--   Check which edge is to be handled
    33163542      IF ( edge == 'l' )  THEN
    3317          IF ( var == 'u' )  THEN   !  For u, nxl is a ghost node, but not for the other variables.
     3543!
     3544!--      For u, nxl is a ghost node, but not for the other variables
     3545         IF ( var == 'u' )  THEN
    33183546            i  = nxl
    33193547            ib = nxl - 1
     
    33273555      ENDIF
    33283556     
    3329       DO  j = nys, nyn + 1
    3330          DO  k = kb(j,i), nzt + 1
     3557      DO  j = nys, nyn+1
     3558         DO  k = kb(j,i), nzt+1
    33313559            l = ic(i)
    33323560            m = jc(j)
     
    33453573!--   Generalized log-law-correction algorithm.
    33463574!--   Doubly two-dimensional index arrays logc(:,:,1:2) and log-ratio arrays
    3347 !--   logc_ratio(:,:,1:2,0:ncorr-1) have been precomputed in subroutine pmci_init_loglaw_correction.
     3575!--   logc_ratio(:,:,1:2,0:ncorr-1) have been precomputed in subroutine
     3576!--   pmci_init_loglaw_correction.
    33483577!
    33493578!--   Solid surface below the node
    33503579      IF ( var == 'u' .OR. var == 'v' )  THEN           
    33513580         DO  j = nys, nyn
    3352             k = kb(j,i) + 1
     3581            k = kb(j,i)+1
    33533582            IF ( ( logc(k,j,1) /= 0 )  .AND.  ( logc(k,j,2) == 0 ) )  THEN
    33543583               k1 = logc(k,j,1)
     
    33623591
    33633592!
    3364 !--   In case of non-flat topography, also vertical walls and corners need to be treated.
    3365 !--   Only single and double wall nodes are corrected. Triple and higher-multiple wall nodes
    3366 !--   are not corrected as the log law would not be valid anyway in such locations.
     3593!--   In case of non-flat topography, also vertical walls and corners need to be
     3594!--   treated. Only single and double wall nodes are corrected. Triple and
     3595!--   higher-multiple wall nodes are not corrected as the log law would not be
     3596!--   valid anyway in such locations.
    33673597      IF ( topography /= 'flat' )  THEN
    33683598         IF ( var == 'u' .OR. var == 'w' )  THEN                 
     
    33713601!--         Solid surface only on south/north side of the node                   
    33723602            DO  j = nys, nyn
    3373                DO  k = kb(j,i) + 1, nzt_topo_nestbc
     3603               DO  k = kb(j,i)+1, nzt_topo_nestbc
    33743604                  IF ( ( logc(k,j,2) /= 0 )  .AND.  ( logc(k,j,1) == 0 ) )  THEN
    33753605
    33763606!
    3377 !--                  Direction of the wall-normal index is carried in as the sign of logc.           
     3607!--                  Direction of the wall-normal index is carried in as the
     3608!--                  sign of logc
    33783609                     jinc = SIGN( 1, logc(k,j,2) )
    33793610                     j1   = ABS( logc(k,j,2) )
    3380                      DO  jcorr=0, ncorr - 1
     3611                     DO  jcorr = 0, ncorr-1
    33813612                        jco = j + jinc * jcorr
    33823613                        f(k,jco,i) = logc_ratio(k,j,2,jcorr) * f(k,j1,i)
     
    33963627                  jinc = SIGN( 1, logc(k,j,2) )
    33973628                  j1   = ABS( logc(k,j,2) )                 
    3398                   DO  jcorr = 0, ncorr - 1
     3629                  DO  jcorr = 0, ncorr-1
    33993630                     jco = j + jinc * jcorr
    3400                      DO  kcorr = 0, ncorr - 1
     3631                     DO  kcorr = 0, ncorr-1
    34013632                        kco = k + kcorr
    3402                         f(kco,jco,i) = 0.5_wp * ( logc_ratio(k,j,1,kcorr) * f(k1,j,i)  &
    3403                                      + logc_ratio(k,j,2,jcorr) * f(k,j1,i) )
     3633                        f(kco,jco,i) = 0.5_wp * ( logc_ratio(k,j,1,kcorr) *    &
     3634                                                  f(k1,j,i)                    &
     3635                                                + logc_ratio(k,j,2,jcorr) *    &
     3636                                                  f(k,j1,i) )
    34043637                     ENDDO
    34053638                  ENDDO
     
    34203653            ENDDO
    34213654         ELSEIF ( edge == 'r' )  THEN           
    3422             DO  j = nys, nyn + 1
    3423                DO  k = kb(j,i), nzt + 1
     3655            DO  j = nys, nyn+1
     3656               DO  k = kb(j,i), nzt+1
    34243657                  f(k,j,i) = tkefactor_r(k,j) * f(k,j,i)
    34253658               ENDDO
     
    34313664!--   Store the boundary values also into the other redundant ghost node layers
    34323665      IF ( edge == 'l' )  THEN
    3433          DO ibgp = -nbgp, ib
     3666         DO  ibgp = -nbgp, ib
    34343667            f(0:nzt+1,nysg:nyng,ibgp) = f(0:nzt+1,nysg:nyng,i)
    34353668         ENDDO
    3436       ELSE IF ( edge == 'r' )  THEN
    3437          DO ibgp = ib, nx + nbgp
     3669      ELSEIF ( edge == 'r' )  THEN
     3670         DO  ibgp = ib, nx+nbgp
    34383671            f(0:nzt+1,nysg:nyng,ibgp) = f(0:nzt+1,nysg:nyng,i)
    34393672         ENDDO
     
    34443677
    34453678
    3446    SUBROUTINE pmci_interp_tril_sn( f, fc, ic, jc, kc, r1x, r2x, r1y, r2y, r1z, r2z, kb, logc, logc_ratio, &
    3447                                   nzt_topo_nestbc, edge, var )
     3679   SUBROUTINE pmci_interp_tril_sn( f, fc, ic, jc, kc, r1x, r2x, r1y, r2y, r1z, &
     3680                                   r2z, kb, logc, logc_ratio,                  &
     3681                                   nzt_topo_nestbc, edge, var )
    34483682
    34493683!
     
    34513685!--   conditions. This subroutine handles the south and north boundaries.
    34523686!--   This subroutine is based on trilinear interpolation.
    3453 !--   Constant dz is still assumed.
    3454 !
    3455 !--                                      Antti Hellsten 22.2.2015.
    3456 !
    3457 !--   Rewritten so that all the coefficients and client-array indices are
    3458 !--   precomputed in the initialization phase by pmci_init_interp_tril.
    3459 !
    3460 !--                                      Antti Hellsten 3.3.2015.
    3461 !
    3462 !--   Constant dz no more assumed.
    3463 !--                                      Antti Hellsten 23.3.2015.
    3464 !
    3465 !--   Adapted for non-flat topography. However, the near-wall velocities
    3466 !--   are log-corrected only over horifontal surfaces, not yet near vertical
    3467 !--   walls.
    3468 !--                                      Antti Hellsten 26.3.2015.
    3469 !
    3470 !--   Indexing in the principal direction (j) is changed. Now, the nest-boundary
    3471 !--   values are interpolated only into the first ghost-node layers on each later
    3472 !--   boundary. These values are then simply copied to the second ghost-node layer.
    3473 !
    3474 !--                                      Antti Hellsten 6.10.2015.
     3687
    34753688      IMPLICIT NONE
    3476       REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg), INTENT(INOUT)          ::  f                 !:
    3477       REAL(wp), DIMENSION(0:cg%nz+1,jcs:jcn,icl:icr), INTENT(IN)                 ::  fc                !:
    3478       REAL(wp), DIMENSION(nzb:nzt_topo_nestbc,nxl:nxr,1:2,0:ncorr-1), INTENT(IN) ::  logc_ratio        !:
    3479       REAL(wp), DIMENSION(nxlg:nxrg), INTENT(IN)                                 ::  r1x               !:
    3480       REAL(wp), DIMENSION(nxlg:nxrg), INTENT(IN)                                 ::  r2x               !:
    3481       REAL(wp), DIMENSION(nysg:nyng), INTENT(IN)                                 ::  r1y               !:
    3482       REAL(wp), DIMENSION(nysg:nyng), INTENT(IN)                                 ::  r2y               !:
    3483       REAL(wp), DIMENSION(nzb:nzt+1), INTENT(IN)                                 ::  r1z               !:
    3484       REAL(wp), DIMENSION(nzb:nzt+1), INTENT(IN)                                 ::  r2z               !:
     3689
     3690      REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg),                      &
     3691                                      INTENT(INOUT) ::  f             !:
     3692      REAL(wp), DIMENSION(0:cg%nz+1,jcs:jcn,icl:icr),                          &
     3693                                      INTENT(IN)    ::  fc            !:
     3694      REAL(wp), DIMENSION(nzb:nzt_topo_nestbc,nxl:nxr,1:2,0:ncorr-1),          &
     3695                                      INTENT(IN)    ::  logc_ratio    !:
     3696      REAL(wp), DIMENSION(nxlg:nxrg), INTENT(IN)    ::  r1x           !:
     3697      REAL(wp), DIMENSION(nxlg:nxrg), INTENT(IN)    ::  r2x           !:
     3698      REAL(wp), DIMENSION(nysg:nyng), INTENT(IN)    ::  r1y           !:
     3699      REAL(wp), DIMENSION(nysg:nyng), INTENT(IN)    ::  r2y           !:
     3700      REAL(wp), DIMENSION(nzb:nzt+1), INTENT(IN)    ::  r1z           !:
     3701      REAL(wp), DIMENSION(nzb:nzt+1), INTENT(IN)    ::  r2z           !:
    34853702     
    3486       INTEGER(iwp), DIMENSION(nxlg:nxrg), INTENT(IN)                             ::  ic                !:
    3487       INTEGER(iwp), DIMENSION(nysg:nyng), INTENT(IN)                             ::  jc                !:
    3488       INTEGER(iwp), DIMENSION(nysg:nyng,nxlg:nxrg), INTENT(IN)                   ::  kb                !:
    3489       INTEGER(iwp), DIMENSION(nzb:nzt+1), INTENT(IN)                             ::  kc                !:
    3490       INTEGER(iwp), DIMENSION(nzb:nzt_topo_nestbc,nxl:nxr,1:2), INTENT(IN)       ::  logc              !:
    3491       INTEGER(iwp)                                                               ::  nzt_topo_nestbc   !:
     3703      INTEGER(iwp), DIMENSION(nxlg:nxrg), INTENT(IN)           ::  ic    !:
     3704      INTEGER(iwp), DIMENSION(nysg:nyng), INTENT(IN)           ::  jc    !:
     3705      INTEGER(iwp), DIMENSION(nysg:nyng,nxlg:nxrg), INTENT(IN) ::  kb    !:
     3706      INTEGER(iwp), DIMENSION(nzb:nzt+1), INTENT(IN)           ::  kc    !:
     3707      INTEGER(iwp), DIMENSION(nzb:nzt_topo_nestbc,nxl:nxr,1:2),                &
     3708                                          INTENT(IN)           ::  logc  !:
     3709      INTEGER(iwp) ::  nzt_topo_nestbc   !:
    34923710
    34933711      CHARACTER(LEN=1), INTENT(IN) ::  edge   !:
     
    35213739     
    35223740!
    3523 !--   Check which edge is to be handled: south or north. Note the assumption that the same PE never
    3524 !--   holds both south and north nest boundaries. Should this be changed?
     3741!--   Check which edge is to be handled: south or north
    35253742      IF ( edge == 's' )  THEN
    3526          IF ( var == 'v' )  THEN   !  For v, nys is a ghost node, but not for the other variables.
     3743!
     3744!--      For v, nys is a ghost node, but not for the other variables
     3745         IF ( var == 'v' )  THEN
    35273746            j  = nys
    35283747            jb = nys - 1
     
    35363755      ENDIF
    35373756
    3538       DO  i = nxl, nxr + 1
    3539          DO  k = kb(j,i), nzt + 1
     3757      DO  i = nxl, nxr+1
     3758         DO  k = kb(j,i), nzt+1
    35403759            l = ic(i)
    35413760            m = jc(j)
     
    35543773!--   Generalized log-law-correction algorithm.
    35553774!--   Multiply two-dimensional index arrays logc(:,:,1:2) and log-ratio arrays
    3556 !--   logc_ratio(:,:,1:2,0:ncorr-1) have been precomputed in subroutine pmci_init_loglaw_correction.
     3775!--   logc_ratio(:,:,1:2,0:ncorr-1) have been precomputed in subroutine
     3776!--   pmci_init_loglaw_correction.
    35573777!
    35583778!--   Solid surface below the node
     
    35623782            IF ( ( logc(k,i,1) /= 0 )  .AND.  ( logc(k,i,2) == 0 ) )  THEN
    35633783               k1 = logc(k,i,1)
    3564                DO  kcorr = 0, ncorr - 1
     3784               DO  kcorr = 0, ncorr-1
    35653785                  kco = k + kcorr
    35663786                  f(kco,j,i) = logc_ratio(k,i,1,kcorr) * f(k1,j,i)
     
    35713791
    35723792!
    3573 !--   In case of non-flat topography, also vertical walls and corners need to be treated.
    3574 !--   Only single and double wall nodes are corrected.
    3575 !--   Triple and higher-multiple wall nodes are not corrected as it would be extremely complicated
    3576 !--   and the log law would not be valid anyway in such locations.
     3793!--   In case of non-flat topography, also vertical walls and corners need to be
     3794!--   treated. Only single and double wall nodes are corrected.
     3795!--   Triple and higher-multiple wall nodes are not corrected as it would be
     3796!--   extremely complicated and the log law would not be valid anyway in such
     3797!--   locations.
    35773798      IF ( topography /= 'flat' )  THEN
    35783799         IF ( var == 'v' .OR. var == 'w' )  THEN
     
    35853806
    35863807!
    3587 !--                  Direction of the wall-normal index is carried in as the sign of logc.           
     3808!--                  Direction of the wall-normal index is carried in as the
     3809!--                  sign of logc
    35883810                     iinc = SIGN( 1, logc(k,i,2) )
    35893811                     i1  = ABS( logc(k,i,2) )
    3590                      DO  icorr = 0, ncorr - 1
     3812                     DO  icorr = 0, ncorr-1
    35913813                        ico = i + iinc * icorr
    35923814                        f(k,j,ico) = logc_ratio(k,i,2,icorr) * f(k,j,i1)
     
    36063828                  iinc = SIGN( 1, logc(k,i,2) )
    36073829                  i1   = ABS( logc(k,i,2) )
    3608                   DO  icorr = 0, ncorr - 1
     3830                  DO  icorr = 0, ncorr-1
    36093831                     ico = i + iinc * icorr
    3610                      DO  kcorr = 0, ncorr - 1
     3832                     DO  kcorr = 0, ncorr-1
    36113833                        kco = k + kcorr
    3612                         f(kco,i,ico) = 0.5_wp * ( logc_ratio(k,i,1,kcorr) * f(k1,j,i)  &
    3613                                      + logc_ratio(k,i,2,icorr) * f(k,j,i1) )
     3834                        f(kco,i,ico) = 0.5_wp * ( logc_ratio(k,i,1,kcorr) *    &
     3835                                                  f(k1,j,i)  &
     3836                                                + logc_ratio(k,i,2,icorr) *    &
     3837                                                  f(k,j,i1) )
    36143838                     ENDDO
    36153839                  ENDDO
     
    36253849         IF ( edge == 's' )  THEN
    36263850            DO  i = nxl, nxr + 1
    3627                DO  k = kb(j,i), nzt + 1
     3851               DO  k = kb(j,i), nzt+1
    36283852                  f(k,j,i) = tkefactor_s(k,i) * f(k,j,i)
    36293853               ENDDO
     
    36313855         ELSEIF ( edge == 'n' )  THEN
    36323856            DO  i = nxl, nxr + 1
    3633                DO  k = kb(j,i), nzt + 1
     3857               DO  k = kb(j,i), nzt+1
    36343858                  f(k,j,i) = tkefactor_n(k,i) * f(k,j,i)
    36353859               ENDDO
     
    36413865!--   Store the boundary values also into the other redundant ghost node layers
    36423866      IF ( edge == 's' )  THEN
    3643          DO jbgp = -nbgp, jb
     3867         DO  jbgp = -nbgp, jb
    36443868            f(0:nzt+1,jbgp,nxlg:nxrg) = f(0:nzt+1,j,nxlg:nxrg)
    36453869         ENDDO
    3646       ELSE IF ( edge == 'n' )  THEN
    3647          DO jbgp = jb, ny + nbgp
     3870      ELSEIF ( edge == 'n' )  THEN
     3871         DO  jbgp = jb, ny+nbgp
    36483872            f(0:nzt+1,jbgp,nxlg:nxrg) = f(0:nzt+1,j,nxlg:nxrg)
    36493873         ENDDO
     
    36543878 
    36553879
    3656    SUBROUTINE pmci_interp_tril_t( f, fc, ic, jc, kc, r1x, r2x, r1y, r2y, r1z, r2z, var )
     3880   SUBROUTINE pmci_interp_tril_t( f, fc, ic, jc, kc, r1x, r2x, r1y, r2y, r1z,  &
     3881                                  r2z, var )
    36573882
    36583883!
     
    36603885!--   conditions. This subroutine handles the top boundary.
    36613886!--   This subroutine is based on trilinear interpolation.
    3662 !--   Constant dz is still assumed.
     3887
    36633888      IMPLICIT NONE
    3664       REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg), INTENT(INOUT) ::  f     !:
    3665       REAL(wp), DIMENSION(0:cg%nz+1,jcs:jcn,icl:icr), INTENT(IN)        ::  fc    !:
    3666       REAL(wp), DIMENSION(nxlg:nxrg), INTENT(IN)                        ::  r1x   !:
    3667       REAL(wp), DIMENSION(nxlg:nxrg), INTENT(IN)                        ::  r2x   !:
    3668       REAL(wp), DIMENSION(nysg:nyng), INTENT(IN)                        ::  r1y   !:
    3669       REAL(wp), DIMENSION(nysg:nyng), INTENT(IN)                        ::  r2y   !:
    3670       REAL(wp), DIMENSION(nzb:nzt+1), INTENT(IN)                        ::  r1z   !:
    3671       REAL(wp), DIMENSION(nzb:nzt+1), INTENT(IN)                        ::  r2z   !:
     3889
     3890      REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg),                      &
     3891                                      INTENT(INOUT) ::  f     !:
     3892      REAL(wp), DIMENSION(0:cg%nz+1,jcs:jcn,icl:icr),                          &
     3893                                      INTENT(IN)    ::  fc    !:
     3894      REAL(wp), DIMENSION(nxlg:nxrg), INTENT(IN)    ::  r1x   !:
     3895      REAL(wp), DIMENSION(nxlg:nxrg), INTENT(IN)    ::  r2x   !:
     3896      REAL(wp), DIMENSION(nysg:nyng), INTENT(IN)    ::  r1y   !:
     3897      REAL(wp), DIMENSION(nysg:nyng), INTENT(IN)    ::  r2y   !:
     3898      REAL(wp), DIMENSION(nzb:nzt+1), INTENT(IN)    ::  r1z   !:
     3899      REAL(wp), DIMENSION(nzb:nzt+1), INTENT(IN)    ::  r2z   !:
    36723900     
    3673       INTEGER(iwp), DIMENSION(nxlg:nxrg), INTENT(IN)                    ::  ic    !:
    3674       INTEGER(iwp), DIMENSION(nysg:nyng), INTENT(IN)                    ::  jc    !:
    3675       INTEGER(iwp), DIMENSION(nzb:nzt+1), INTENT(IN)                    ::  kc    !:
     3901      INTEGER(iwp), DIMENSION(nxlg:nxrg), INTENT(IN) ::  ic    !:
     3902      INTEGER(iwp), DIMENSION(nysg:nyng), INTENT(IN) ::  jc    !:
     3903      INTEGER(iwp), DIMENSION(nzb:nzt+1), INTENT(IN) ::  kc    !:
    36763904     
    36773905      CHARACTER(LEN=1), INTENT(IN) :: var   !:
     
    37013929      ENDIF
    37023930     
    3703       DO  i = nxl - 1, nxr + 1
    3704          DO  j = nys - 1, nyn + 1
     3931      DO  i = nxl-1, nxr+1
     3932         DO  j = nys-1, nyn+1
    37053933            l = ic(i)
    37063934            m = jc(j)
     
    37243952!
    37253953!--   Rescale if f is the TKE.
    3726 !--   It is assumed that the bottom surface never reaches the top
    3727 !---  boundary of a nest domain.
    3728       IF ( var == 'e')  THEN
     3954!--   It is assumed that the bottom surface never reaches the top boundary of a
     3955!--   nest domain.
     3956      IF ( var == 'e' )  THEN
    37293957         DO  i = nxl, nxr
    37303958            DO  j = nys, nyn
     
    37463974!--    subroutine handles the left and right boundaries. However, this operation
    37473975!--    is only needed in case of one-way coupling.
     3976
    37483977       IMPLICIT NONE
    37493978
     
    37854014       ENDIF
    37864015
    3787        DO  j = nys, nyn + 1
    3788           DO  k = kb(j,i), nzt +1
     4016       DO  j = nys, nyn+1
     4017          DO  k = kb(j,i), nzt+1
    37894018             vdotnor = outnor * u(k,j,ied)
    37904019!
     
    38064035          ENDDO
    38074036       ELSEIF ( edge == 'r' )  THEN
    3808           DO ibgp = ib, nx + nbgp
     4037          DO ibgp = ib, nx+nbgp
    38094038             f(0:nzt+1,nysg:nyng,ibgp) = f(0:nzt+1,nysg:nyng,i)
    38104039          ENDDO
     
    38224051!--    interpolated values by values extrapolated from the domain. This
    38234052!--    subroutine handles the south and north boundaries.
     4053
    38244054       IMPLICIT NONE
    38254055
     
    38614091       ENDIF
    38624092
    3863        DO  i = nxl, nxr + 1
    3864           DO  k = kb(j,i), nzt + 1
     4093       DO  i = nxl, nxr+1
     4094          DO  k = kb(j,i), nzt+1
    38654095             vdotnor = outnor * v(k,jed,i)
    38664096!
     
    38784108!--    Store the boundary values also into the redundant ghost node layers.
    38794109       IF ( edge == 's' )  THEN
    3880           DO jbgp = -nbgp, jb
     4110          DO  jbgp = -nbgp, jb
    38814111             f(0:nzt+1,jbgp,nxlg:nxrg) = f(0:nzt+1,j,nxlg:nxrg)
    38824112          ENDDO
    38834113       ELSEIF ( edge == 'n' )  THEN
    3884           DO jbgp = jb, ny + nbgp
     4114          DO  jbgp = jb, ny+nbgp
    38854115             f(0:nzt+1,jbgp,nxlg:nxrg) = f(0:nzt+1,j,nxlg:nxrg)
    38864116          ENDDO
     
    38964126!--    conditions. This subroutine handles the top boundary. It is based on
    38974127!--    trilinear interpolation.
     4128
    38984129       IMPLICIT NONE
    38994130
     
    39414172
    39424173    SUBROUTINE pmci_anterp_tophat( f, fc, kct, ifl, ifu, jfl, jfu, kfl, kfu,   &
    3943                                    kb, var )
     4174                                   var )
    39444175!
    39454176!--    Anterpolation of internal-node values to be used as the server-domain
     
    39474178!--    integration of the fine-grid values contained within the coarse-grid
    39484179!--    cell.
     4180
    39494181       IMPLICIT NONE
    39504182
     
    39684200       INTEGER(iwp), INTENT(IN) ::  kct   !:
    39694201
    3970        INTEGER(iwp), DIMENSION(icl:icr), INTENT(IN)             ::  ifl   !:
    3971        INTEGER(iwp), DIMENSION(icl:icr), INTENT(IN)             ::  ifu   !:
    3972        INTEGER(iwp), DIMENSION(jcs:jcn), INTENT(IN)             ::  jfl   !:
    3973        INTEGER(iwp), DIMENSION(jcs:jcn), INTENT(IN)             ::  jfu   !:
    3974 !--    TO_DO: is the next line really unnecessary?
    3975        INTEGER(iwp), DIMENSION(nysg:nyng,nxlg:nxrg), INTENT(IN) ::  kb    !: may be unnecessary
    3976        INTEGER(iwp), DIMENSION(0:kct), INTENT(IN)               ::  kfl   !:
    3977        INTEGER(iwp), DIMENSION(0:kct), INTENT(IN)               ::  kfu   !:
     4202       INTEGER(iwp), DIMENSION(icl:icr), INTENT(IN) ::  ifl   !:
     4203       INTEGER(iwp), DIMENSION(icl:icr), INTENT(IN) ::  ifu   !:
     4204       INTEGER(iwp), DIMENSION(jcs:jcn), INTENT(IN) ::  jfl   !:
     4205       INTEGER(iwp), DIMENSION(jcs:jcn), INTENT(IN) ::  jfu   !:
     4206       INTEGER(iwp), DIMENSION(0:kct), INTENT(IN)   ::  kfl   !:
     4207       INTEGER(iwp), DIMENSION(0:kct), INTENT(IN)   ::  kfu   !:
    39784208
    39794209
     
    40494279                   ENDIF
    40504280                ENDIF
    4051 
    4052                 fc(kk,jj,ii) = ( 1.0_wp - fra ) * fc(kk,jj,ii) +                 &
     4281!
     4282!--             TO DO: introduce 3-d coarse grid array for precomputed
     4283!--             1/REAL(nfc) values
     4284                fc(kk,jj,ii) = ( 1.0_wp - fra ) * fc(kk,jj,ii) +               &
    40534285                               fra * cellsum / REAL( nfc, KIND = wp )
    40544286
  • palm/trunk/SOURCE/pmc_server.f90

    r1792 r1797  
    2020! Current revisions:
    2121! ------------------
    22 !
     22! introduction of different datatransfer modes
    2323!
    2424! Former revisions:
     
    498498      REAL(wp), POINTER, DIMENSION(:,:,:) ::  data_3d
    499499
    500       t1 = PMC_Time()
    501       CALL MPI_Barrier(Clients(ClientId)%intra_comm, ierr)              ! Wait for buffer empty
    502       t2 = PMC_Time()
    503       if(present(WaitTime)) WaitTime = t2-t1
     500!--   Synchronization of the model is done in pmci_client_synchronize and pmci_server_synchronize
     501!--   Therefor the RMA window cann be filled without sychronization at this point and the barrier
     502!--   is not necessary
     503!--   Please note that WaitTime has to be set in PMC_S_FillBuffer AND PMC_C_GetBuffer
     504
     505      if(present(WaitTime)) then
     506        t1 = PMC_Time()
     507        CALL MPI_Barrier(Clients(ClientId)%intra_comm, ierr)
     508        t2 = PMC_Time()
     509        WaitTime = t2-t1
     510      end if
    504511
    505512      do ip=1,Clients(ClientId)%inter_npes
  • palm/trunk/SOURCE/time_integration.f90

    r1792 r1797  
    1919! Current revisions:
    2020! ------------------
    21 !
     21! introduction of different datatransfer modes
    2222!
    2323! Former revisions:
     
    258258    USE pmc_interface,                                                         &
    259259        ONLY:  client_to_server, nested_run, nesting_mode,                     &
    260                pmci_ensure_nest_mass_conservation, pmci_client_datatrans,      &
    261                pmci_client_synchronize, pmci_server_datatrans,                 &
     260               pmci_client_synchronize, pmci_datatrans,                        &
     261               pmci_ensure_nest_mass_conservation,                             &
    262262               pmci_server_synchronize, server_to_client
    263263
     
    673673
    674674          IF ( nested_run )  THEN
    675 !
    676 !--          TO_DO: try to give more meaningful comments here
    677 !--          Domain nesting
    678 !--          Note that the nesting operations are omitted intentionally on the
    679 !--          first two RK-substeps.
     675
    680676             CALL cpu_log( log_point(60), 'nesting', 'start' )
    681677!
    682 !--          From server to client commmunication ( direction=SERVER_TO_CLIENT )
    683              CALL pmci_server_datatrans( server_to_client )
    684              CALL pmci_client_datatrans( server_to_client )
     678!--          Domain nesting. The data transfer subroutines pmci_server_datatrans
     679!--          and pmci_client_datatatrans are called inside the wrapper
     680!--          subroutine pmci_datatrans according to the control parameters
     681!--          nesting_mode and nesting_datatransfer_mode.
     682!--          TO_DO: why is nesting_mode given as a parameter here?
     683             CALL pmci_datatrans( nesting_mode )
    685684
    686685             IF ( nesting_mode == 'two-way' )  THEN
    687 !
    688 !--             From client to server commmunication ( direction=CLIENT_TO_SERVER )
    689                 CALL pmci_server_datatrans( client_to_server )
    690                 CALL pmci_client_datatrans( client_to_server )
    691686!
    692687!--             Exchange_horiz is needed for all server-domains after the
Note: See TracChangeset for help on using the changeset viewer.