Index: palm/trunk/SOURCE/init_pegrid.f90
===================================================================
--- palm/trunk/SOURCE/init_pegrid.f90	(revision 708)
+++ palm/trunk/SOURCE/init_pegrid.f90	(revision 709)
@@ -4,5 +4,5 @@
 ! Current revisions:
 ! -----------------
-! 
+! formatting adjustments
 ! 
 ! ATTENTION: nnz_x undefined problem still has to be solved!!!!!!!!
@@ -577,17 +577,7 @@
 
           CALL MPI_OPEN_PORT( MPI_INFO_NULL, port_name, ierr )
-!
-!--       TEST OUTPUT (TO BE REMOVED)
-          WRITE(9,*)  TRIM( coupling_mode ),  &
-               ', ierr after MPI_OPEN_PORT: ', ierr
-          CALL LOCAL_FLUSH( 9 )
 
           CALL MPI_PUBLISH_NAME( 'palm_coupler', MPI_INFO_NULL, port_name, &
                                  ierr )
-!
-!--       TEST OUTPUT (TO BE REMOVED)
-          WRITE(9,*)  TRIM( coupling_mode ),  &
-               ', ierr after MPI_PUBLISH_NAME: ', ierr
-          CALL LOCAL_FLUSH( 9 )
 
 !
@@ -614,10 +604,4 @@
 
           CALL MPI_LOOKUP_NAME( 'palm_coupler', MPI_INFO_NULL, port_name, ierr )
-!
-!--       TEST OUTPUT (TO BE REMOVED)
-          WRITE(9,*)  TRIM( coupling_mode ),  &
-               ', ierr after MPI_LOOKUP_NAME: ', ierr
-          CALL LOCAL_FLUSH( 9 )
-
 
        ENDIF
@@ -631,21 +615,12 @@
     IF ( coupling_mode == 'atmosphere_to_ocean' )  THEN
 
-       PRINT*, '... before COMM_ACCEPT'
        CALL MPI_COMM_ACCEPT( port_name, MPI_INFO_NULL, 0, MPI_COMM_WORLD, &
                              comm_inter, ierr )
-       PRINT*, '--- ierr = ', ierr
-       PRINT*, '--- comm_inter atmosphere = ', comm_inter
-
        coupling_mode_remote = 'ocean_to_atmosphere'
 
     ELSEIF ( coupling_mode == 'ocean_to_atmosphere' )  THEN
 
-       IF ( myid == 0 )  PRINT*, '*** read: ', port_name, '  ierr = ', ierr
-       PRINT*, '... before COMM_CONNECT'
        CALL MPI_COMM_CONNECT( port_name, MPI_INFO_NULL, 0, MPI_COMM_WORLD, &
                               comm_inter, ierr )
-       PRINT*, '--- ierr = ', ierr
-       PRINT*, '--- comm_inter ocean      = ', comm_inter
-
        coupling_mode_remote = 'atmosphere_to_ocean'
 
@@ -654,26 +629,20 @@
 
 ! 
-!-- Determine the number of ghost points
-    IF (scalar_advec == 'ws-scheme' .OR. momentum_advec == 'ws-scheme') THEN
+!-- Determine the number of ghost point layers
+    IF ( scalar_advec == 'ws-scheme' .OR. momentum_advec == 'ws-scheme' )  THEN
        nbgp = 3
     ELSE
        nbgp = 1
-    END IF 
-
-!
-!-- In case of coupled runs, create a new MPI derived datatype for the
-!-- exchange of surface (xy) data .
-!-- Gridpoint number for the exchange of ghost points (xy-plane)
-
+    ENDIF 
+
+!
+!-- Create a new MPI derived datatype for the exchange of surface (xy) data,
+!-- which is needed for coupled atmosphere-ocean runs.
+!-- First, calculate number of grid points of an xy-plane.
     ngp_xy  = ( nxr - nxl + 1 + 2 * nbgp ) * ( nyn - nys + 1 + 2 * nbgp )
-
-!
-!-- Define a new MPI derived datatype for the exchange of ghost points in
-!-- y-direction for 2D-arrays (line)
     CALL MPI_TYPE_VECTOR( ngp_xy, 1, nzt-nzb+2, MPI_REAL, type_xy, ierr )
     CALL MPI_TYPE_COMMIT( type_xy, ierr )
 
-
-    IF ( TRIM( coupling_mode ) .NE. 'uncoupled' ) THEN
+    IF ( TRIM( coupling_mode ) /= 'uncoupled' )  THEN
    
 !
@@ -685,22 +654,23 @@
           ny_a = ny
 
-          IF ( myid == 0 ) THEN
-             CALL MPI_SEND( nx_a, 1, MPI_INTEGER, numprocs, 1, &
-                            comm_inter, ierr )
-             CALL MPI_SEND( ny_a, 1, MPI_INTEGER, numprocs, 2, &
-                            comm_inter, ierr )
-             CALL MPI_SEND( pdims, 2, MPI_INTEGER, numprocs, 3, &
-                            comm_inter, ierr )
-             CALL MPI_RECV( nx_o, 1, MPI_INTEGER, numprocs, 4, &
-                            comm_inter, status, ierr )
-             CALL MPI_RECV( ny_o, 1, MPI_INTEGER, numprocs, 5, &
-                            comm_inter, status, ierr )
-             CALL MPI_RECV( pdims_remote, 2, MPI_INTEGER, numprocs, 6, &
+          IF ( myid == 0 )  THEN
+
+             CALL MPI_SEND( nx_a, 1, MPI_INTEGER, numprocs, 1, comm_inter,  &
+                            ierr )
+             CALL MPI_SEND( ny_a, 1, MPI_INTEGER, numprocs, 2, comm_inter,  &
+                            ierr )
+             CALL MPI_SEND( pdims, 2, MPI_INTEGER, numprocs, 3, comm_inter, &
+                            ierr )
+             CALL MPI_RECV( nx_o, 1, MPI_INTEGER, numprocs, 4, comm_inter,  &
+                            status, ierr )
+             CALL MPI_RECV( ny_o, 1, MPI_INTEGER, numprocs, 5, comm_inter,  &
+                            status, ierr )
+             CALL MPI_RECV( pdims_remote, 2, MPI_INTEGER, numprocs, 6,      &
                             comm_inter, status, ierr )
           ENDIF
 
-          CALL MPI_BCAST( nx_o, 1, MPI_INTEGER, 0, comm2d, ierr)
-          CALL MPI_BCAST( ny_o, 1, MPI_INTEGER, 0, comm2d, ierr) 
-          CALL MPI_BCAST( pdims_remote, 2, MPI_INTEGER, 0, comm2d, ierr)
+          CALL MPI_BCAST( nx_o, 1, MPI_INTEGER, 0, comm2d, ierr )
+          CALL MPI_BCAST( ny_o, 1, MPI_INTEGER, 0, comm2d, ierr ) 
+          CALL MPI_BCAST( pdims_remote, 2, MPI_INTEGER, 0, comm2d, ierr )
        
        ELSEIF ( coupling_mode == 'ocean_to_atmosphere' )  THEN
@@ -710,16 +680,14 @@
 
           IF ( myid == 0 ) THEN
-             CALL MPI_RECV( nx_a, 1, MPI_INTEGER, 0, 1, &
-                            comm_inter, status, ierr )
-             CALL MPI_RECV( ny_a, 1, MPI_INTEGER, 0, 2, &
-                            comm_inter, status, ierr )
-             CALL MPI_RECV( pdims_remote, 2, MPI_INTEGER, 0, 3, &
-                            comm_inter, status, ierr )
-             CALL MPI_SEND( nx_o, 1, MPI_INTEGER, 0, 4, &
-                            comm_inter, ierr )
-             CALL MPI_SEND( ny_o, 1, MPI_INTEGER, 0, 5, &
-                            comm_inter, ierr )
-             CALL MPI_SEND( pdims, 2, MPI_INTEGER, 0, 6, &
-                            comm_inter, ierr )
+
+             CALL MPI_RECV( nx_a, 1, MPI_INTEGER, 0, 1, comm_inter, status, &
+                            ierr )
+             CALL MPI_RECV( ny_a, 1, MPI_INTEGER, 0, 2, comm_inter, status, &
+                            ierr )
+             CALL MPI_RECV( pdims_remote, 2, MPI_INTEGER, 0, 3, comm_inter, &
+                            status, ierr )
+             CALL MPI_SEND( nx_o, 1, MPI_INTEGER, 0, 4, comm_inter, ierr )
+             CALL MPI_SEND( ny_o, 1, MPI_INTEGER, 0, 5, comm_inter, ierr )
+             CALL MPI_SEND( pdims, 2, MPI_INTEGER, 0, 6, comm_inter, ierr )
           ENDIF
 
@@ -730,12 +698,11 @@
        ENDIF
   
-       ngp_a = (nx_a+1+2*nbgp)*(ny_a+1+2*nbgp)
-       ngp_o = (nx_o+1+2*nbgp)*(ny_o+1+2*nbgp)
-
-!
-!--    determine if the horizontal grid and the number of PEs 
-!--    in ocean and atmosphere is same or not
-!--    (different number of PEs still not implemented)
-       IF ( nx_o == nx_a .AND. ny_o == ny_a .AND.  &
+       ngp_a = ( nx_a+1 + 2 * nbgp ) * ( ny_a+1 + 2 * nbgp )
+       ngp_o = ( nx_o+1 + 2 * nbgp ) * ( ny_o+1 + 2 * nbgp )
+
+!
+!--    Determine if the horizontal grid and the number of PEs in ocean and
+!--    atmosphere is same or not
+       IF ( nx_o == nx_a  .AND.  ny_o == ny_a  .AND.  &
             pdims(1) == pdims_remote(1) .AND. pdims(2) == pdims_remote(2) ) &
        THEN
@@ -748,6 +715,9 @@
 !--    Determine the target PEs for the exchange between ocean and
 !--    atmosphere (comm2d)
-       IF ( coupling_topology == 0) THEN
-          IF ( TRIM( coupling_mode ) .EQ. 'atmosphere_to_ocean' ) THEN
+       IF ( coupling_topology == 0 )  THEN
+!
+!--       In case of identical topologies, every atmosphere PE has exactly one
+!--       ocean PE counterpart and vice versa
+          IF ( TRIM( coupling_mode ) == 'atmosphere_to_ocean' ) THEN
              target_id = myid + numprocs
           ELSE
@@ -756,18 +726,19 @@
 
        ELSE
-
 !
 !--       In case of nonequivalent topology in ocean and atmosphere only for
 !--       PE0 in ocean and PE0 in atmosphere a target_id is needed, since
-!--       data echxchange between ocean and atmosphere will be done only by
-!--       those PEs.    
-          IF ( myid == 0 ) THEN
-             IF ( TRIM( coupling_mode ) .EQ. 'atmosphere_to_ocean' ) THEN
+!--       data echxchange between ocean and atmosphere will be done only
+!--       between these PEs.    
+          IF ( myid == 0 )  THEN
+
+             IF ( TRIM( coupling_mode ) == 'atmosphere_to_ocean' )  THEN
                 target_id = numprocs 
              ELSE
                 target_id = 0
              ENDIF
- print*, coupling_mode, myid, " -> ", target_id, "numprocs: ", numprocs 
+
           ENDIF
+
        ENDIF
 
@@ -861,7 +832,9 @@
 !
 !--    Find out, if the total domain allows more levels. These additional
-!--    levels are processed on PE0 only.
+!--    levels are identically processed on all PEs.
        IF ( numprocs > 1  .AND.  mg_switch_to_pe0_level /= -1 )  THEN
+
           IF ( mg_levels_z > MIN( mg_levels_x, mg_levels_y ) )  THEN
+
              mg_switch_to_pe0_level_l = maximum_grid_level
 
@@ -889,7 +862,10 @@
                 mg_switch_to_pe0_level_l = 0
              ENDIF
+
           ELSE
+
              mg_switch_to_pe0_level_l = 0
              maximum_grid_level_l = maximum_grid_level
+
           ENDIF
 
@@ -920,4 +896,5 @@
 
              ENDIF
+
           ENDIF
 
@@ -939,5 +916,4 @@
 !--          Save the grid size of the subdomain at the switch level, because
 !--          it is needed in poismg.
-!--          Array bounds of the local subdomain grids are gathered on PE0
              ind(1) = nxl_l; ind(2) = nxr_l
              ind(3) = nys_l; ind(4) = nyn_l
@@ -953,5 +929,5 @@
              DEALLOCATE( ind_all )
 !
-!--          Calculate the grid size of the total domain gathered on PE0
+!--          Calculate the grid size of the total domain
              nxr_l = ( nxr_l-nxl_l+1 ) * pdims(1) - 1
              nxl_l = 0
@@ -1006,9 +982,11 @@
 
 !
-!-- Define a new MPI derived datatype for the exchange of ghost points in
-!-- y-direction for 2D-arrays (line)
-    CALL MPI_TYPE_VECTOR( nxr-nxl+1+2*nbgp, nbgp, ngp_y, MPI_REAL, type_x, ierr )
+!-- Define new MPI derived datatypes for the exchange of ghost points in
+!-- x- and y-direction for 2D-arrays (line)
+    CALL MPI_TYPE_VECTOR( nxr-nxl+1+2*nbgp, nbgp, ngp_y, MPI_REAL, type_x, &
+                          ierr )
     CALL MPI_TYPE_COMMIT( type_x, ierr )
-    CALL MPI_TYPE_VECTOR( nxr-nxl+1+2*nbgp, nbgp, ngp_y, MPI_INTEGER, type_x_int, ierr )
+    CALL MPI_TYPE_VECTOR( nxr-nxl+1+2*nbgp, nbgp, ngp_y, MPI_INTEGER, &
+                          type_x_int, ierr )
     CALL MPI_TYPE_COMMIT( type_x_int, ierr )
 
@@ -1029,37 +1007,43 @@
 
     nxl_l = nxl; nxr_l = nxr; nys_l = nys; nyn_l = nyn; nzb_l = nzb; nzt_l = nzt
+
 !
 !-- Discern between the model grid, which needs nbgp ghost points and
 !-- grid levels for the multigrid scheme. In the latter case only one
 !-- ghost point is necessary.
-!-- First definition of mpi-vectors for exchange of ghost layers on normal 
+!-- First definition of MPI-datatypes for exchange of ghost layers on normal 
 !-- grid. The following loop is needed for data exchange in poismg.f90.
 !
 !-- Determine number of grid points of yz-layer for exchange
     ngp_yz(0) = (nzt - nzb + 2) * (nyn - nys + 1 + 2 * nbgp)
-!
-!-- Define a new mpi datatype for the exchange of left - right boundaries.
-!-- Indeed the data are connected in the physical memory and no mpi-vector
-!-- is necessary, but the data exchange between left and right PE's using
-!-- mpi-vectors is 10% faster than without.
+
+!
+!-- Define an MPI-datatype for the exchange of left/right boundaries.
+!-- Although data are contiguous in physical memory (which does not
+!-- necessarily require an MPI-derived datatype), the data exchange between
+!-- left and right PE's using the MPI-derived type is 10% faster than without.
     CALL MPI_TYPE_VECTOR( nxr-nxl+1+2*nbgp, nbgp*(nzt-nzb+2), ngp_yz(0), &
-                             MPI_REAL, type_xz(0), ierr )
+                          MPI_REAL, type_xz(0), ierr )
     CALL MPI_TYPE_COMMIT( type_xz(0), ierr )
 
-    CALL MPI_TYPE_VECTOR( nbgp, ngp_yz(0), ngp_yz(0), MPI_REAL, type_yz(0), ierr) 
+    CALL MPI_TYPE_VECTOR( nbgp, ngp_yz(0), ngp_yz(0), MPI_REAL, type_yz(0), &
+                          ierr ) 
     CALL MPI_TYPE_COMMIT( type_yz(0), ierr )
-!
-!-- Definition of mpi-vectors for multigrid
+
+!
+!-- Definition of MPI-datatypes for multigrid method (coarser level grids)
     IF ( psolver == 'multigrid' )  THEN
 !    
-!--   The definition of mpi-vectors as aforementioned, but only 1 ghost point is used.
-       DO i = maximum_grid_level, 1 , -1
+!--    Definition of MPI-datatyoe as above, but only 1 ghost level is used
+       DO  i = maximum_grid_level, 1 , -1
+
           ngp_yz(i) = (nzt_l - nzb_l + 2) * (nyn_l - nys_l + 3)
 
           CALL MPI_TYPE_VECTOR( nxr_l-nxl_l+3, nzt_l-nzb_l+2, ngp_yz(i), &
-                             MPI_REAL, type_xz(i), ierr )
+                                MPI_REAL, type_xz(i), ierr )
           CALL MPI_TYPE_COMMIT( type_xz(i), ierr )
 
-          CALL MPI_TYPE_VECTOR( 1, ngp_yz(i), ngp_yz(i), MPI_REAL, type_yz(i), ierr)
+          CALL MPI_TYPE_VECTOR( 1, ngp_yz(i), ngp_yz(i), MPI_REAL, type_yz(i), &
+                                ierr )
           CALL MPI_TYPE_COMMIT( type_yz(i), ierr )
 
@@ -1069,6 +1053,8 @@
           nyn_l = nyn_l / 2
           nzt_l = nzt_l / 2
+
        ENDDO
-    END IF
+
+    ENDIF
 #endif
 
