Index: /palm/trunk/SOURCE/advec_s_bc.f90
===================================================================
--- /palm/trunk/SOURCE/advec_s_bc.f90	(revision 4487)
+++ /palm/trunk/SOURCE/advec_s_bc.f90	(revision 4488)
@@ -1,20 +1,19 @@
 !> @file advec_s_bc.f90
-!------------------------------------------------------------------------------!
+!--------------------------------------------------------------------------------------------------!
 ! This file is part of the PALM model system.
 !
-! PALM is free software: you can redistribute it and/or modify it under the
-! terms of the GNU General Public License as published by the Free Software
-! Foundation, either version 3 of the License, or (at your option) any later
-! version.
-!
-! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
-! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
-! A PARTICULAR PURPOSE.  See the GNU General Public License for more details.
-!
-! You should have received a copy of the GNU General Public License along with
-! PALM. If not, see <http://www.gnu.org/licenses/>.
+! PALM is free software: you can redistribute it and/or modify it under the terms of the GNU General
+! Public License as published by the Free Software Foundation, either version 3 of the License, or
+! (at your option) any later version.
+!
+! PALM is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the
+! implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General
+! Public License for more details.
+!
+! You should have received a copy of the GNU General Public License along with PALM. If not, see
+! <http://www.gnu.org/licenses/>.
 !
 ! Copyright 1997-2020 Leibniz Universitaet Hannover
-!------------------------------------------------------------------------------!
+! -------------------------------------------------------------------------------------------------!
 !
 ! Current revisions:
@@ -25,12 +24,15 @@
 ! -----------------
 ! $Id$
+! file re-formatted to follow the PALM coding standard
+!
+! 4429 2020-02-27 15:24:30Z raasch
 ! bugfix: cpp-directives added for serial mode
-! 
+!
 ! 4360 2020-01-07 11:25:50Z suehring
 ! Corrected "Former revisions" section
-! 
+!
 ! 3761 2019-02-25 15:31:42Z raasch
 ! unused variables removed
-! 
+!
 ! 3655 2019-01-07 16:51:22Z knoop
 ! nopointer option removed
@@ -45,13 +47,12 @@
 !> Computation in individual steps for each of the three dimensions.
 !> Limiting assumptions:
-!> So far the scheme has been assuming equidistant grid spacing. As this is not
-!> the case in the stretched portion of the z-direction, there dzw(k) is used as
-!> a substitute for a constant grid length. This certainly causes incorrect
-!> results; however, it is hoped that they are not too apparent for weakly 
-!> stretched grids.
+!> So far the scheme has been assuming equidistant grid spacing. As this is not the case in the
+!> stretched portion of the z-direction, there dzw(k) is used as a substitute for a constant grid
+!> length. This certainly causes incorrect results; however, it is hoped that they are not too
+!> apparent for weakly stretched grids.
 !> NOTE: This is a provisional, non-optimised version!
-!------------------------------------------------------------------------------!
-MODULE advec_s_bc_mod
- 
+!--------------------------------------------------------------------------------------------------!
+ MODULE advec_s_bc_mod
+
 
     PRIVATE
@@ -64,29 +65,28 @@
  CONTAINS
 
-!------------------------------------------------------------------------------!
+!--------------------------------------------------------------------------------------------------!
 ! Description:
 ! ------------
 !> @todo Missing subroutine description.
-!------------------------------------------------------------------------------!
+!--------------------------------------------------------------------------------------------------!
     SUBROUTINE advec_s_bc( sk, sk_char )
 
-       USE advection,                                                             &
+       USE advection,                                                                              &
            ONLY:  aex, bex, dex, eex
 
-       USE arrays_3d,                                                             &
+       USE arrays_3d,                                                                              &
            ONLY:  d, ddzw, dzu, dzw, tend, u, v, w
 
-       USE control_parameters,                                                    &
-           ONLY:  dt_3d, bc_pt_t_val, bc_q_t_val, bc_s_t_val, ibc_pt_b, ibc_pt_t, &
-                  ibc_q_t, ibc_s_t, message_string, pt_slope_offset,              &
-                  sloping_surface, u_gtrans, v_gtrans
-
-       USE cpulog,                                                                &
+       USE control_parameters,                                                                     &
+           ONLY:  bc_pt_t_val, bc_q_t_val, bc_s_t_val, dt_3d, ibc_pt_b, ibc_pt_t, ibc_q_t, ibc_s_t,&
+                  message_string, pt_slope_offset, sloping_surface, u_gtrans, v_gtrans
+
+       USE cpulog,                                                                                 &
            ONLY:  cpu_log, log_point_s
 
-       USE grid_variables,                                                        &
+       USE grid_variables,                                                                         &
            ONLY:  ddx, ddy
 
-       USE indices,                                                               &
+       USE indices,                                                                                &
            ONLY:  nx, nxl, nxr, nyn, nys, nzb, nzt
 
@@ -95,5 +95,5 @@
        USE pegrid
 
-       USE statistics,                                                            &
+       USE statistics,                                                                             &
            ONLY:  rmask, statistic_regions, sums_wsts_bc_l
 
@@ -112,37 +112,34 @@
        INTEGER(iwp) ::  type_xz_2 !<
 #endif
-
-       REAL(wp) ::  cim    !<
-       REAL(wp) ::  cimf   !<
-       REAL(wp) ::  cip    !<
-       REAL(wp) ::  cipf   !<
-       REAL(wp) ::  d_new  !<
-       REAL(wp) ::  denomi !< denominator
-       REAL(wp) ::  fminus !<
-       REAL(wp) ::  fplus  !<
-       REAL(wp) ::  f2     !<
-       REAL(wp) ::  f4     !<
-       REAL(wp) ::  f8     !<
-       REAL(wp) ::  f12    !<
-       REAL(wp) ::  f24    !<
-       REAL(wp) ::  f48    !<
-       REAL(wp) ::  f1920  !<
-       REAL(wp) ::  im     !<
-       REAL(wp) ::  ip     !<
-       REAL(wp) ::  m1n    !<
-       REAL(wp) ::  m1z    !<
-       REAL(wp) ::  m2     !<
-       REAL(wp) ::  m3     !<
-       REAL(wp) ::  numera !< numerator
-       REAL(wp) ::  snenn  !<
-       REAL(wp) ::  sterm  !<
-       REAL(wp) ::  tendcy !<
-       REAL(wp) ::  t1     !<
-       REAL(wp) ::  t2     !<
+       REAL(wp) ::  cim       !<
+       REAL(wp) ::  cimf      !<
+       REAL(wp) ::  cip       !<
+       REAL(wp) ::  cipf      !<
+       REAL(wp) ::  d_new     !<
+       REAL(wp) ::  denomi    !< denominator
+       REAL(wp) ::  fminus    !<
+       REAL(wp) ::  fplus     !<
+       REAL(wp) ::  f2        !<
+       REAL(wp) ::  f4        !<
+       REAL(wp) ::  f8        !<
+       REAL(wp) ::  f12       !<
+       REAL(wp) ::  f24       !<
+       REAL(wp) ::  f48       !<
+       REAL(wp) ::  f1920     !<
+       REAL(wp) ::  im        !<
+       REAL(wp) ::  ip        !<
+       REAL(wp) ::  m1n       !<
+       REAL(wp) ::  m1z       !<
+       REAL(wp) ::  m2        !<
+       REAL(wp) ::  m3        !<
+       REAL(wp) ::  numera    !< numerator
+       REAL(wp) ::  snenn     !<
+       REAL(wp) ::  sterm     !<
+       REAL(wp) ::  tendcy    !<
+       REAL(wp) ::  t1        !<
+       REAL(wp) ::  t2        !<
 
        REAL(wp) ::  fmax(2)   !<
        REAL(wp) ::  fmax_l(2) !<
-       
-       REAL(wp), DIMENSION(:,:,:), POINTER ::  sk
 
        REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  a0   !<
@@ -161,6 +158,9 @@
        REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  m1   !<
        REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  sw   !<
-       
+
        REAL(wp), DIMENSION(:,:,:), ALLOCATABLE ::  sk_p !<
+
+       REAL(wp), DIMENSION(:,:,:), POINTER ::  sk !<
+
 
 !
@@ -197,11 +197,11 @@
 !
 !--    Send left boundary, receive right boundary
-       CALL MPI_SENDRECV( sk_p(nzb-2,nys-3,nxl+1), ngp, MPI_REAL, pleft,  0,      &
-                          sk_p(nzb-2,nys-3,nxr+2), ngp, MPI_REAL, pright, 0,      &
+       CALL MPI_SENDRECV( sk_p(nzb-2,nys-3,nxl+1), ngp, MPI_REAL, pleft,  0,                       &
+                          sk_p(nzb-2,nys-3,nxr+2), ngp, MPI_REAL, pright, 0,                       &
                           comm2d, status, ierr )
 !
 !--    Send right boundary, receive left boundary
-       CALL MPI_SENDRECV( sk_p(nzb-2,nys-3,nxr-2), ngp, MPI_REAL, pright, 1,      &
-                          sk_p(nzb-2,nys-3,nxl-3), ngp, MPI_REAL, pleft,  1,      &
+       CALL MPI_SENDRECV( sk_p(nzb-2,nys-3,nxr-2), ngp, MPI_REAL, pright, 1,                       &
+                          sk_p(nzb-2,nys-3,nxl-3), ngp, MPI_REAL, pleft,  1,                       &
                           comm2d, status, ierr )
        CALL cpu_log( log_point_s(11), 'advec_s_bc:sendrecv', 'pause' )
@@ -217,7 +217,7 @@
 
 !
-!--    In case of a sloping surface, the additional gridpoints in x-direction
-!--    of the temperature field at the left and right boundary of the total
-!--    domain must be adjusted by the temperature difference between this distance
+!--    In case of a sloping surface, the additional gridpoints in x-direction of the temperature
+!--    field at the left and right boundary of the total domain must be adjusted by the temperature
+!--    difference between this distance
        IF ( sloping_surface  .AND.  sk_char == 'pt' )  THEN
           IF ( nxl ==  0 )  THEN
@@ -241,5 +241,5 @@
           DO  j = nys, nyn
              DO  k = nzb+1, nzt
-                numera = ABS( sk_p(k,j,i+1) - 2.0_wp * sk_p(k,j,i) + sk_p(k,j,i-1) )
+                numera  = ABS( sk_p(k,j,i+1) - 2.0_wp * sk_p(k,j,i) + sk_p(k,j,i-1) )
                 denomi  = ABS( sk_p(k,j,i+1) - sk_p(k,j,i-1) )
                 fmax_l(1) = MAX( fmax_l(1) , numera )
@@ -259,23 +259,20 @@
 !
 !--    Allocate temporary arrays
-       ALLOCATE( a0(nzb+1:nzt,nxl-1:nxr+1),   a1(nzb+1:nzt,nxl-1:nxr+1),          &
-                 a2(nzb+1:nzt,nxl-1:nxr+1),   a12(nzb+1:nzt,nxl-1:nxr+1),         &
-                 a22(nzb+1:nzt,nxl-1:nxr+1),  immb(nzb+1:nzt,nxl-1:nxr+1),        &
-                 imme(nzb+1:nzt,nxl-1:nxr+1), impb(nzb+1:nzt,nxl-1:nxr+1),        &
-                 impe(nzb+1:nzt,nxl-1:nxr+1), ipmb(nzb+1:nzt,nxl-1:nxr+1),        &
-                 ipme(nzb+1:nzt,nxl-1:nxr+1), ippb(nzb+1:nzt,nxl-1:nxr+1),        &
-                 ippe(nzb+1:nzt,nxl-1:nxr+1), m1(nzb+1:nzt,nxl-2:nxr+2),          &
-                 sw(nzb+1:nzt,nxl-1:nxr+1)                                        &
+       ALLOCATE( a0(nzb+1:nzt,nxl-1:nxr+1),   a1(nzb+1:nzt,nxl-1:nxr+1), a2(nzb+1:nzt,nxl-1:nxr+1),&
+                 a12(nzb+1:nzt,nxl-1:nxr+1),  a22(nzb+1:nzt,nxl-1:nxr+1),                          &
+                 immb(nzb+1:nzt,nxl-1:nxr+1), imme(nzb+1:nzt,nxl-1:nxr+1),                         &
+                 impb(nzb+1:nzt,nxl-1:nxr+1), impe(nzb+1:nzt,nxl-1:nxr+1),                         &
+                 ipmb(nzb+1:nzt,nxl-1:nxr+1), ipme(nzb+1:nzt,nxl-1:nxr+1),                         &
+                 ippb(nzb+1:nzt,nxl-1:nxr+1), ippe(nzb+1:nzt,nxl-1:nxr+1),                         &
+                 m1(nzb+1:nzt,nxl-2:nxr+2),   sw(nzb+1:nzt,nxl-1:nxr+1)                            &
                )
        imme = 0.0_wp; impe = 0.0_wp; ipme = 0.0_wp; ippe = 0.0_wp
 
 !
-!--    Initialise point of time measuring of the exponential portion (this would
-!--    not work if done locally within the loop)
+!--    Initialise point of time measuring of the exponential portion (this would not work if done
+!--    locally within the loop)
        CALL cpu_log( log_point_s(12), 'advec_s_bc:exp', 'start' )
        CALL cpu_log( log_point_s(12), 'advec_s_bc:exp', 'pause' )
 
-!
-!--    Outer loop of all j
        DO  j = nys, nyn
 
@@ -285,15 +282,16 @@
              DO  k = nzb+1, nzt
                 a12(k,i) = 0.5_wp * ( sk_p(k,j,i+1) - sk_p(k,j,i-1) )
-                a22(k,i) = 0.5_wp * ( sk_p(k,j,i+1) - 2.0_wp * sk_p(k,j,i)        &
-                                                    + sk_p(k,j,i-1) )
-                a0(k,i) = ( 9.0_wp * sk_p(k,j,i+2)    - 116.0_wp * sk_p(k,j,i+1)  &
-                            + 2134.0_wp * sk_p(k,j,i) - 116.0_wp * sk_p(k,j,i-1)  &
-                            + 9.0_wp * sk_p(k,j,i-2) ) * f1920
-                a1(k,i) = ( -5.0_wp * sk_p(k,j,i+2)   + 34.0_wp * sk_p(k,j,i+1)   &
-                            - 34.0_wp * sk_p(k,j,i-1) + 5.0_wp * sk_p(k,j,i-2)    &
+                a22(k,i) = 0.5_wp * ( sk_p(k,j,i+1) - 2.0_wp * sk_p(k,j,i) + sk_p(k,j,i-1) )
+                a0(k,i) = ( 9.0_wp * sk_p(k,j,i+2)    - 116.0_wp * sk_p(k,j,i+1)                   &
+                            + 2134.0_wp * sk_p(k,j,i) - 116.0_wp * sk_p(k,j,i-1)                   &
+                            + 9.0_wp * sk_p(k,j,i-2)                                               &
+                          ) * f1920
+                a1(k,i) = ( -5.0_wp * sk_p(k,j,i+2)   + 34.0_wp * sk_p(k,j,i+1)                    &
+                            - 34.0_wp * sk_p(k,j,i-1) + 5.0_wp * sk_p(k,j,i-2)                     &
                           ) * f48
-                a2(k,i) = ( -3.0_wp * sk_p(k,j,i+2) + 36.0_wp * sk_p(k,j,i+1)     &
-                            - 66.0_wp * sk_p(k,j,i) + 36.0_wp * sk_p(k,j,i-1)     &
-                            - 3.0_wp * sk_p(k,j,i-2) ) * f48
+                a2(k,i) = ( -3.0_wp * sk_p(k,j,i+2) + 36.0_wp * sk_p(k,j,i+1)                      &
+                            - 66.0_wp * sk_p(k,j,i) + 36.0_wp * sk_p(k,j,i-1)                      &
+                            - 3.0_wp * sk_p(k,j,i-2)                                               &
+                          ) * f48
              ENDDO
           ENDDO
@@ -301,5 +299,4 @@
 !
 !--       Fluxes using the Bott scheme
-!--       *VOCL LOOP,UNROLL(2)
           DO  i = nxl, nxr
              DO  k = nzb+1, nzt
@@ -308,9 +305,9 @@
                 cipf = 1.0_wp - 2.0_wp * cip
                 cimf = 1.0_wp - 2.0_wp * cim
-                ip   =   a0(k,i)   * f2  * ( 1.0_wp - cipf )                      &
-                       + a1(k,i)   * f8  * ( 1.0_wp - cipf*cipf )                 &
+                ip   =   a0(k,i)   * f2  * ( 1.0_wp - cipf )                                       &
+                       + a1(k,i)   * f8  * ( 1.0_wp - cipf*cipf )                                  &
                        + a2(k,i)   * f24 * ( 1.0_wp - cipf*cipf*cipf )
-                im   =   a0(k,i+1) * f2  * ( 1.0_wp - cimf )                      &
-                       - a1(k,i+1) * f8  * ( 1.0_wp - cimf*cimf )                 &
+                im   =   a0(k,i+1) * f2  * ( 1.0_wp - cimf )                                       &
+                       - a1(k,i+1) * f8  * ( 1.0_wp - cimf*cimf )                                  &
                        + a2(k,i+1) * f24 * ( 1.0_wp - cimf*cimf*cimf )
                 ip   = MAX( ip, 0.0_wp )
@@ -323,9 +320,9 @@
                 cipf = 1.0_wp - 2.0_wp * cip
                 cimf = 1.0_wp - 2.0_wp * cim
-                ip   =   a0(k,i-1) * f2  * ( 1.0_wp - cipf )                      &
-                       + a1(k,i-1) * f8  * ( 1.0_wp - cipf*cipf )                 &
+                ip   =   a0(k,i-1) * f2  * ( 1.0_wp - cipf )                                       &
+                       + a1(k,i-1) * f8  * ( 1.0_wp - cipf*cipf )                                  &
                        + a2(k,i-1) * f24 * ( 1.0_wp - cipf*cipf*cipf )
-                im   =   a0(k,i)   * f2  * ( 1.0_wp - cimf )                      &
-                       - a1(k,i)   * f8  * ( 1.0_wp - cimf*cimf )                 &
+                im   =   a0(k,i)   * f2  * ( 1.0_wp - cimf )                                       &
+                       - a1(k,i)   * f8  * ( 1.0_wp - cimf*cimf )                                  &
                        + a2(k,i)   * f24 * ( 1.0_wp - cimf*cimf*cimf )
                 ip   = MAX( ip, 0.0_wp )
@@ -358,10 +355,8 @@
           DO  i = nxl-1, nxr+1
              DO  k = nzb+1, nzt
-                m2 = 2.0_wp * ABS( a1(k,i) - a12(k,i) ) /                         &
-                     MAX( ABS( a1(k,i) + a12(k,i) ), 1E-35_wp )
+                m2 = 2.0_wp * ABS( a1(k,i) - a12(k,i) ) / MAX( ABS( a1(k,i) + a12(k,i) ), 1E-35_wp )
                 IF ( ABS( a1(k,i) + a12(k,i) ) < fmax(2) )  m2 = 0.0_wp
 
-                m3 = 2.0_wp * ABS( a2(k,i) - a22(k,i) ) /                         &
-                     MAX( ABS( a2(k,i) + a22(k,i) ), 1E-35_wp )
+                m3 = 2.0_wp * ABS( a2(k,i) - a22(k,i) ) / MAX( ABS( a2(k,i) + a22(k,i) ), 1E-35_wp )
                 IF ( ABS( a2(k,i) + a22(k,i) ) < fmax(1) )  m3 = 0.0_wp
 
@@ -369,10 +364,7 @@
                 t2 = 0.35_wp
                 IF ( m1(k,i) == -1.0_wp )  t2 = 0.12_wp
-
-!--             *VOCL STMT,IF(10)
-                IF ( m1(k,i-1) == 1.0_wp .OR. m1(k,i) == 1.0_wp                   &
-                     .OR. m1(k,i+1) == 1.0_wp .OR.  m2 > t2  .OR.  m3 > t2  .OR.  &
-                     ( m1(k,i) > t1  .AND.  m1(k,i-1) /= -1.0_wp  .AND.           &
-                       m1(k,i) /= -1.0_wp  .AND.  m1(k,i+1) /= -1.0_wp )          &
+                IF ( m1(k,i-1) == 1.0_wp  .OR.  m1(k,i) == 1.0_wp  .OR.  m1(k,i+1) == 1.0_wp  .OR. &
+                     m2 > t2  .OR.  m3 > t2  .OR.  ( m1(k,i) > t1 .AND. m1(k,i-1) /= -1.0_wp .AND. &
+                     m1(k,i) /= -1.0_wp .AND. m1(k,i+1) /= -1.0_wp )                               &
                    )  sw(k,i) = 1.0_wp
              ENDDO
@@ -385,5 +377,4 @@
              DO  k = nzb+1, nzt
 
-!--             *VOCL STMT,IF(10)
                 IF ( sw(k,i) == 1.0_wp )  THEN
                    snenn = sk_p(k,j,i+1) - sk_p(k,j,i-1)
@@ -397,9 +388,9 @@
                    cip =  MAX( 0.0_wp, ( u(k,j,i+1) - u_gtrans ) * dt_3d * ddx )
 
-                   ippe(k,i) = sk_p(k,j,i-1) * cip + snenn * (                    &
-                               aex(ix) * cip + bex(ix) / dex(ix) * (              &
-                               eex(ix) -                                          &
-                               EXP( dex(ix)*0.5_wp * ( 1.0_wp - 2.0_wp * cip ) )  &
-                                                                   )              &
+                   ippe(k,i) = sk_p(k,j,i-1) * cip + snenn * (                                     &
+                               aex(ix) * cip + bex(ix) / dex(ix) * (                               &
+                               eex(ix) -                                                           &
+                               EXP( dex(ix)*0.5_wp * ( 1.0_wp - 2.0_wp * cip ) )                   &
+                                                                   )                               &
                                                              )
                    IF ( sterm == 0.0001_wp )  ippe(k,i) = sk_p(k,j,i) * cip
@@ -416,9 +407,9 @@
                    cim = -MIN( 0.0_wp, ( u(k,j,i) - u_gtrans ) * dt_3d * ddx )
 
-                   imme(k,i) = sk_p(k,j,i+1) * cim + snenn * (                    &
-                               aex(ix) * cim + bex(ix) / dex(ix) * (              &
-                               eex(ix) -                                          &
-                               EXP( dex(ix)*0.5_wp * ( 1.0_wp - 2.0_wp * cim ) )  &
-                                                                   )              &
+                   imme(k,i) = sk_p(k,j,i+1) * cim + snenn * (                                     &
+                               aex(ix) * cim + bex(ix) / dex(ix) * (                               &
+                               eex(ix) -                                                           &
+                               EXP( dex(ix)*0.5_wp * ( 1.0_wp - 2.0_wp * cim ) )                   &
+                                                                   )                               &
                                                              )
                    IF ( sterm == 0.0001_wp )  imme(k,i) = sk_p(k,j,i) * cim
@@ -426,5 +417,4 @@
                 ENDIF
 
-!--             *VOCL STMT,IF(10)
                 IF ( sw(k,i+1) == 1.0_wp )  THEN
                    snenn = sk_p(k,j,i) - sk_p(k,j,i+2)
@@ -438,9 +428,9 @@
                    cim = -MIN( 0.0_wp, ( u(k,j,i+1) - u_gtrans ) * dt_3d * ddx )
 
-                   impe(k,i) = sk_p(k,j,i+2) * cim + snenn * (                    &
-                               aex(ix) * cim + bex(ix) / dex(ix) * (              &
-                               eex(ix) -                                          &
-                               EXP( dex(ix)*0.5_wp * ( 1.0_wp - 2.0_wp * cim ) )  &
-                                                                   )              &
+                   impe(k,i) = sk_p(k,j,i+2) * cim + snenn * (                                     &
+                               aex(ix) * cim + bex(ix) / dex(ix) * (                               &
+                               eex(ix) -                                                           &
+                               EXP( dex(ix)*0.5_wp * ( 1.0_wp - 2.0_wp * cim ) )                   &
+                                                                   )                               &
                                                              )
                    IF ( sterm == 0.0001_wp )  impe(k,i) = sk_p(k,j,i+1) * cim
@@ -448,5 +438,4 @@
                 ENDIF
 
-!--             *VOCL STMT,IF(10)
                 IF ( sw(k,i-1) == 1.0_wp )  THEN
                    snenn = sk_p(k,j,i) - sk_p(k,j,i-2)
@@ -460,9 +449,9 @@
                    cip = MAX( 0.0_wp, ( u(k,j,i) - u_gtrans ) * dt_3d * ddx )
 
-                   ipme(k,i) = sk_p(k,j,i-2) * cip + snenn * (                    &
-                               aex(ix) * cip + bex(ix) / dex(ix) * (              &
-                               eex(ix) -                                          &
-                               EXP( dex(ix)*0.5_wp * ( 1.0_wp - 2.0_wp * cip ) )  &
-                                                                   )              &
+                   ipme(k,i) = sk_p(k,j,i-2) * cip + snenn * (                                     &
+                               aex(ix) * cip + bex(ix) / dex(ix) * (                               &
+                               eex(ix) -                                                           &
+                               EXP( dex(ix)*0.5_wp * ( 1.0_wp - 2.0_wp * cip ) )                   &
+                                                                   )                               &
                                                              )
                    IF ( sterm == 0.0001_wp )  ipme(k,i) = sk_p(k,j,i-1) * cip
@@ -478,8 +467,8 @@
           DO  i = nxl, nxr
              DO  k = nzb+1, nzt
-                fplus  = ( 1.0_wp - sw(k,i)   ) * ippb(k,i) + sw(k,i)   * ippe(k,i)  &
-                       - ( 1.0_wp - sw(k,i+1) ) * impb(k,i) - sw(k,i+1) * impe(k,i)
-                fminus = ( 1.0_wp - sw(k,i-1) ) * ipmb(k,i) + sw(k,i-1) * ipme(k,i)  &
-                       - ( 1.0_wp - sw(k,i)   ) * immb(k,i) - sw(k,i)   * imme(k,i)
+                fplus  = ( 1.0_wp - sw(k,i) )     * ippb(k,i) + sw(k,i)   * ippe(k,i)              &
+                         - ( 1.0_wp - sw(k,i+1) ) * impb(k,i) - sw(k,i+1) * impe(k,i)
+                fminus = ( 1.0_wp - sw(k,i-1) )   * ipmb(k,i) + sw(k,i-1) * ipme(k,i)              &
+                         - ( 1.0_wp - sw(k,i) )   * immb(k,i) - sw(k,i)   * imme(k,i)
                 tendcy = fplus - fminus
 !
@@ -490,6 +479,5 @@
 !--             Density correction because of possible remaining divergences
                 d_new = d(k,j,i) - ( u(k,j,i+1) - u(k,j,i) ) * dt_3d * ddx
-                sk_p(k,j,i) = ( ( 1.0_wp + d(k,j,i) ) * sk_p(k,j,i) - tendcy ) /    &
-                              ( 1.0_wp + d_new )
+                sk_p(k,j,i) = ( ( 1.0_wp + d(k,j,i) ) * sk_p(k,j,i) - tendcy ) / ( 1.0_wp + d_new )
                 d(k,j,i)  = d_new
              ENDDO
@@ -500,6 +488,5 @@
 !
 !--    Deallocate temporary arrays
-       DEALLOCATE( a0, a1, a2, a12, a22, immb, imme, impb, impe, ipmb, ipme,      &
-                   ippb, ippe, m1, sw )
+       DEALLOCATE( a0, a1, a2, a12, a22, immb, imme, impb, impe, ipmb, ipme, ippb, ippe, m1, sw )
 
 
@@ -508,5 +495,5 @@
 #if defined( __parallel )
        ngp = ( nzt - nzb + 6 ) * ( nyn - nys + 7 )
-       CALL MPI_TYPE_VECTOR( nxr-nxl+7, 3*(nzt-nzb+6), ngp, MPI_REAL,             &
+       CALL MPI_TYPE_VECTOR( nxr-nxl+7, 3*(nzt-nzb+6), ngp, MPI_REAL,                              &
                              type_xz_2, ierr )
        CALL MPI_TYPE_COMMIT( type_xz_2, ierr )
@@ -514,11 +501,11 @@
 !--    Send front boundary, receive rear boundary
        CALL cpu_log( log_point_s(11), 'advec_s_bc:sendrecv', 'continue' )
-       CALL MPI_SENDRECV( sk_p(nzb-2,nys,nxl-3),   1, type_xz_2, psouth, 0,       &
-                          sk_p(nzb-2,nyn+1,nxl-3), 1, type_xz_2, pnorth, 0,       &
+       CALL MPI_SENDRECV( sk_p(nzb-2,nys,nxl-3),   1, type_xz_2, psouth, 0,                        &
+                          sk_p(nzb-2,nyn+1,nxl-3), 1, type_xz_2, pnorth, 0,                        &
                           comm2d, status, ierr )
 !
 !--    Send rear boundary, receive front boundary
-       CALL MPI_SENDRECV( sk_p(nzb-2,nyn-2,nxl-3), 1, type_xz_2, pnorth, 1,       &
-                          sk_p(nzb-2,nys-3,nxl-3), 1, type_xz_2, psouth, 1,       &
+       CALL MPI_SENDRECV( sk_p(nzb-2,nyn-2,nxl-3), 1, type_xz_2, pnorth, 1,                        &
+                          sk_p(nzb-2,nys-3,nxl-3), 1, type_xz_2, psouth, 1,                        &
                           comm2d, status, ierr )
        CALL MPI_TYPE_FREE( type_xz_2, ierr )
@@ -561,17 +548,14 @@
 !
 !--    Allocate temporary arrays
-       ALLOCATE( a0(nzb+1:nzt,nys-1:nyn+1),   a1(nzb+1:nzt,nys-1:nyn+1),          &
-                 a2(nzb+1:nzt,nys-1:nyn+1),   a12(nzb+1:nzt,nys-1:nyn+1),         &
-                 a22(nzb+1:nzt,nys-1:nyn+1),  immb(nzb+1:nzt,nys-1:nyn+1),        &
-                 imme(nzb+1:nzt,nys-1:nyn+1), impb(nzb+1:nzt,nys-1:nyn+1),        &
-                 impe(nzb+1:nzt,nys-1:nyn+1), ipmb(nzb+1:nzt,nys-1:nyn+1),        &
-                 ipme(nzb+1:nzt,nys-1:nyn+1), ippb(nzb+1:nzt,nys-1:nyn+1),        &
-                 ippe(nzb+1:nzt,nys-1:nyn+1), m1(nzb+1:nzt,nys-2:nyn+2),          &
-                 sw(nzb+1:nzt,nys-1:nyn+1)                                        &
+       ALLOCATE( a0(nzb+1:nzt,nys-1:nyn+1),   a1(nzb+1:nzt,nys-1:nyn+1), a2(nzb+1:nzt,nys-1:nyn+1),&
+                 a12(nzb+1:nzt,nys-1:nyn+1),  a22(nzb+1:nzt,nys-1:nyn+1),                          &
+                 immb(nzb+1:nzt,nys-1:nyn+1), imme(nzb+1:nzt,nys-1:nyn+1),                         &
+                 impb(nzb+1:nzt,nys-1:nyn+1), impe(nzb+1:nzt,nys-1:nyn+1),                         &
+                 ipmb(nzb+1:nzt,nys-1:nyn+1), ipme(nzb+1:nzt,nys-1:nyn+1),                         &
+                 ippb(nzb+1:nzt,nys-1:nyn+1), ippe(nzb+1:nzt,nys-1:nyn+1),                         &
+                 m1(nzb+1:nzt,nys-2:nyn+2),   sw(nzb+1:nzt,nys-1:nyn+1)                            &
                )
        imme = 0.0_wp; impe = 0.0_wp; ipme = 0.0_wp; ippe = 0.0_wp
 
-!
-!--    Outer loop of all i
        DO  i = nxl, nxr
 
@@ -581,15 +565,16 @@
              DO  k = nzb+1, nzt
                 a12(k,j) = 0.5_wp * ( sk_p(k,j+1,i) - sk_p(k,j-1,i) )
-                a22(k,j) = 0.5_wp * ( sk_p(k,j+1,i) - 2.0_wp * sk_p(k,j,i)        &
-                                                    + sk_p(k,j-1,i) )
-                a0(k,j) = ( 9.0_wp * sk_p(k,j+2,i)    - 116.0_wp * sk_p(k,j+1,i)  &
-                            + 2134.0_wp * sk_p(k,j,i) - 116.0_wp * sk_p(k,j-1,i)  &
-                            + 9.0_wp * sk_p(k,j-2,i) ) * f1920
-                a1(k,j) = ( -5.0_wp   * sk_p(k,j+2,i) + 34.0_wp * sk_p(k,j+1,i)   &
-                            - 34.0_wp * sk_p(k,j-1,i) + 5.0_wp  * sk_p(k,j-2,i)   &
+                a22(k,j) = 0.5_wp * ( sk_p(k,j+1,i) - 2.0_wp * sk_p(k,j,i) + sk_p(k,j-1,i) )
+                a0(k,j) = ( 9.0_wp * sk_p(k,j+2,i)    - 116.0_wp * sk_p(k,j+1,i)                   &
+                            + 2134.0_wp * sk_p(k,j,i) - 116.0_wp * sk_p(k,j-1,i)                   &
+                            + 9.0_wp * sk_p(k,j-2,i)                                               &
+                          ) * f1920
+                a1(k,j) = ( -5.0_wp   * sk_p(k,j+2,i) + 34.0_wp * sk_p(k,j+1,i)                    &
+                            - 34.0_wp * sk_p(k,j-1,i) + 5.0_wp  * sk_p(k,j-2,i)                    &
                           ) * f48
-                a2(k,j) = ( -3.0_wp * sk_p(k,j+2,i) + 36.0_wp * sk_p(k,j+1,i)     &
-                            - 66.0_wp * sk_p(k,j,i) + 36.0_wp * sk_p(k,j-1,i)     &
-                            - 3.0_wp * sk_p(k,j-2,i) ) * f48
+                a2(k,j) = ( -3.0_wp * sk_p(k,j+2,i) + 36.0_wp * sk_p(k,j+1,i)                      &
+                            - 66.0_wp * sk_p(k,j,i) + 36.0_wp * sk_p(k,j-1,i)                      &
+                            - 3.0_wp * sk_p(k,j-2,i)                                               &
+                          ) * f48
              ENDDO
           ENDDO
@@ -597,5 +582,4 @@
 !
 !--       Fluxes using the Bott scheme
-!--       *VOCL LOOP,UNROLL(2)
           DO  j = nys, nyn
              DO  k = nzb+1, nzt
@@ -604,9 +588,9 @@
                 cipf = 1.0_wp - 2.0_wp * cip
                 cimf = 1.0_wp - 2.0_wp * cim
-                ip   =   a0(k,j)   * f2  * ( 1.0_wp - cipf )                      &
-                       + a1(k,j)   * f8  * ( 1.0_wp - cipf*cipf )                 &
+                ip   =   a0(k,j)   * f2  * ( 1.0_wp - cipf )                                       &
+                       + a1(k,j)   * f8  * ( 1.0_wp - cipf*cipf )                                  &
                        + a2(k,j)   * f24 * ( 1.0_wp - cipf*cipf*cipf )
-                im   =   a0(k,j+1) * f2  * ( 1.0_wp - cimf )                      &
-                       - a1(k,j+1) * f8  * ( 1.0_wp - cimf*cimf )                 &
+                im   =   a0(k,j+1) * f2  * ( 1.0_wp - cimf )                                       &
+                       - a1(k,j+1) * f8  * ( 1.0_wp - cimf*cimf )                                  &
                        + a2(k,j+1) * f24 * ( 1.0_wp - cimf*cimf*cimf )
                 ip   = MAX( ip, 0.0_wp )
@@ -619,9 +603,9 @@
                 cipf = 1.0_wp - 2.0_wp * cip
                 cimf = 1.0_wp - 2.0_wp * cim
-                ip   =   a0(k,j-1) * f2  * ( 1.0_wp - cipf )                      &
-                       + a1(k,j-1) * f8  * ( 1.0_wp - cipf*cipf )                 &
+                ip   =   a0(k,j-1) * f2  * ( 1.0_wp - cipf )                                       &
+                       + a1(k,j-1) * f8  * ( 1.0_wp - cipf*cipf )                                  &
                        + a2(k,j-1) * f24 * ( 1.0_wp - cipf*cipf*cipf )
-                im   =   a0(k,j)   * f2  * ( 1.0_wp - cimf )                      &
-                       - a1(k,j)   * f8  * ( 1.0_wp - cimf*cimf )                 &
+                im   =   a0(k,j)   * f2  * ( 1.0_wp - cimf )                                       &
+                       - a1(k,j)   * f8  * ( 1.0_wp - cimf*cimf )                                  &
                        + a2(k,j)   * f24 * ( 1.0_wp - cimf*cimf*cimf )
                 ip   = MAX( ip, 0.0_wp )
@@ -654,9 +638,9 @@
           DO  j = nys-1, nyn+1
              DO  k = nzb+1, nzt
-                m2 = 2.0_wp * ABS( a1(k,j) - a12(k,j) ) /                            &
+                m2 = 2.0_wp * ABS( a1(k,j) - a12(k,j) ) /                                          &
                      MAX( ABS( a1(k,j) + a12(k,j) ), 1E-35_wp )
                 IF ( ABS( a1(k,j) + a12(k,j) ) < fmax(2) )  m2 = 0.0_wp
 
-                m3 = 2.0_wp * ABS( a2(k,j) - a22(k,j) ) /                            &
+                m3 = 2.0_wp * ABS( a2(k,j) - a22(k,j) ) /                                          &
                      MAX( ABS( a2(k,j) + a22(k,j) ), 1E-35_wp )
                 IF ( ABS( a2(k,j) + a22(k,j) ) < fmax(1) )  m3 = 0.0_wp
@@ -666,9 +650,7 @@
                 IF ( m1(k,j) == -1.0_wp )  t2 = 0.12_wp
 
-!--             *VOCL STMT,IF(10)
-                IF ( m1(k,j-1) == 1.0_wp .OR. m1(k,j) == 1.0_wp                   &
-                     .OR. m1(k,j+1) == 1.0_wp .OR.  m2 > t2  .OR.  m3 > t2  .OR.  &
-                     ( m1(k,j) > t1  .AND.  m1(k,j-1) /= -1.0_wp  .AND.           &
-                       m1(k,j) /= -1.0_wp  .AND.  m1(k,j+1) /= -1.0_wp )          &
+                IF ( m1(k,j-1) == 1.0_wp  .OR.  m1(k,j) == 1.0_wp  .OR.  m1(k,j+1) == 1.0_wp  .OR. &
+                     m2 > t2  .OR.  m3 > t2  .OR.  ( m1(k,j) > t1 .AND. m1(k,j-1) /= -1.0_wp .AND. &
+                     m1(k,j) /= -1.0_wp .AND. m1(k,j+1) /= -1.0_wp )                               &
                    )  sw(k,j) = 1.0_wp
              ENDDO
@@ -681,5 +663,4 @@
              DO  k = nzb+1, nzt
 
-!--             *VOCL STMT,IF(10)
                 IF ( sw(k,j) == 1.0_wp )  THEN
                    snenn = sk_p(k,j+1,i) - sk_p(k,j-1,i)
@@ -693,9 +674,9 @@
                    cip =  MAX( 0.0_wp, ( v(k,j+1,i) - v_gtrans ) * dt_3d * ddy )
 
-                   ippe(k,j) = sk_p(k,j-1,i) * cip + snenn * (                    &
-                               aex(ix) * cip + bex(ix) / dex(ix) * (              &
-                               eex(ix) -                                          &
-                               EXP( dex(ix)*0.5_wp * ( 1.0_wp - 2.0_wp * cip ) )  &
-                                                                   )              &
+                   ippe(k,j) = sk_p(k,j-1,i) * cip + snenn * (                                     &
+                               aex(ix) * cip + bex(ix) / dex(ix) * (                               &
+                               eex(ix) -                                                           &
+                               EXP( dex(ix)*0.5_wp * ( 1.0_wp - 2.0_wp * cip ) )                   &
+                                                                   )                               &
                                                              )
                    IF ( sterm == 0.0001_wp )  ippe(k,j) = sk_p(k,j,i) * cip
@@ -712,9 +693,9 @@
                    cim = -MIN( 0.0_wp, ( v(k,j,i) - v_gtrans ) * dt_3d * ddy )
 
-                   imme(k,j) = sk_p(k,j+1,i) * cim + snenn * (                    &
-                               aex(ix) * cim + bex(ix) / dex(ix) * (              &
-                               eex(ix) -                                          &
-                               EXP( dex(ix)*0.5_wp * ( 1.0_wp - 2.0_wp * cim ) )  &
-                                                                   )              &
+                   imme(k,j) = sk_p(k,j+1,i) * cim + snenn * (                                     &
+                               aex(ix) * cim + bex(ix) / dex(ix) * (                               &
+                               eex(ix) -                                                           &
+                               EXP( dex(ix)*0.5_wp * ( 1.0_wp - 2.0_wp * cim ) )                   &
+                                                                   )                               &
                                                              )
                    IF ( sterm == 0.0001_wp )  imme(k,j) = sk_p(k,j,i) * cim
@@ -722,5 +703,4 @@
                 ENDIF
 
-!--             *VOCL STMT,IF(10)
                 IF ( sw(k,j+1) == 1.0_wp )  THEN
                    snenn = sk_p(k,j,i) - sk_p(k,j+2,i)
@@ -734,9 +714,9 @@
                    cim = -MIN( 0.0_wp, ( v(k,j+1,i) - v_gtrans ) * dt_3d * ddy )
 
-                   impe(k,j) = sk_p(k,j+2,i) * cim + snenn * (                    &
-                               aex(ix) * cim + bex(ix) / dex(ix) * (              &
-                               eex(ix) -                                          &
-                               EXP( dex(ix)*0.5_wp * ( 1.0_wp - 2.0_wp * cim ) )  &
-                                                                   )              &
+                   impe(k,j) = sk_p(k,j+2,i) * cim + snenn * (                                     &
+                               aex(ix) * cim + bex(ix) / dex(ix) * (                               &
+                               eex(ix) -                                                           &
+                               EXP( dex(ix)*0.5_wp * ( 1.0_wp - 2.0_wp * cim ) )                   &
+                                                                   )                               &
                                                              )
                    IF ( sterm == 0.0001_wp )  impe(k,j) = sk_p(k,j+1,i) * cim
@@ -744,5 +724,4 @@
                 ENDIF
 
-!--             *VOCL STMT,IF(10)
                 IF ( sw(k,j-1) == 1.0_wp )  THEN
                    snenn = sk_p(k,j,i) - sk_p(k,j-2,i)
@@ -756,9 +735,9 @@
                    cip = MAX( 0.0_wp, ( v(k,j,i) - v_gtrans ) * dt_3d * ddy )
 
-                   ipme(k,j) = sk_p(k,j-2,i) * cip + snenn * (                    &
-                               aex(ix) * cip + bex(ix) / dex(ix) * (              &
-                               eex(ix) -                                          &
-                               EXP( dex(ix)*0.5_wp * ( 1.0_wp - 2.0_wp * cip ) )  &
-                                                                   )              &
+                   ipme(k,j) = sk_p(k,j-2,i) * cip + snenn * (                                     &
+                               aex(ix) * cip + bex(ix) / dex(ix) * (                               &
+                               eex(ix) -                                                           &
+                               EXP( dex(ix)*0.5_wp * ( 1.0_wp - 2.0_wp * cip ) )                   &
+                                                                   )                               &
                                                              )
                    IF ( sterm == 0.0001_wp )  ipme(k,j) = sk_p(k,j-1,i) * cip
@@ -774,8 +753,8 @@
           DO  j = nys, nyn
              DO  k = nzb+1, nzt
-                fplus  = ( 1.0_wp - sw(k,j)   ) * ippb(k,j) + sw(k,j)   * ippe(k,j) &
-                       - ( 1.0_wp - sw(k,j+1) ) * impb(k,j) - sw(k,j+1) * impe(k,j)
-                fminus = ( 1.0_wp - sw(k,j-1) ) * ipmb(k,j) + sw(k,j-1) * ipme(k,j) &
-                       - ( 1.0_wp - sw(k,j)   ) * immb(k,j) - sw(k,j)   * imme(k,j)
+                fplus  = ( 1.0_wp - sw(k,j) )     * ippb(k,j) + sw(k,j)   * ippe(k,j)              &
+                         - ( 1.0_wp - sw(k,j+1) ) * impb(k,j) - sw(k,j+1) * impe(k,j)
+                fminus = ( 1.0_wp - sw(k,j-1) )   * ipmb(k,j) + sw(k,j-1) * ipme(k,j)              &
+                         - ( 1.0_wp - sw(k,j) )   * immb(k,j) - sw(k,j)   * imme(k,j)
                 tendcy = fplus - fminus
 !
@@ -786,11 +765,12 @@
 !--             Density correction because of possible remaining divergences
                 d_new = d(k,j,i) - ( v(k,j+1,i) - v(k,j,i) ) * dt_3d * ddy
-                sk_p(k,j,i) = ( ( 1.0_wp + d(k,j,i) ) * sk_p(k,j,i) - tendcy ) /  &
-                              ( 1.0_wp + d_new )
-                d(k,j,i)  = d_new
+                sk_p(k,j,i) = ( ( 1.0_wp + d(k,j,i) ) * sk_p(k,j,i) - tendcy ) / ( 1.0_wp + d_new )
+                d(k,j,i) = d_new
              ENDDO
           ENDDO
 
        ENDDO   ! End of the advection in y-direction
+
+
        CALL cpu_log( log_point_s(11), 'advec_s_bc:sendrecv', 'continue' )
        CALL cpu_log( log_point_s(11), 'advec_s_bc:sendrecv', 'stop' )
@@ -798,6 +778,5 @@
 !
 !--    Deallocate temporary arrays
-       DEALLOCATE( a0, a1, a2, a12, a22, immb, imme, impb, impe, ipmb, ipme,      &
-                   ippb, ippe, m1, sw )
+       DEALLOCATE( a0, a1, a2, a12, a22, immb, imme, impb, impe, ipmb, ipme, ippb, ippe, m1, sw )
 
 
@@ -959,6 +938,6 @@
 
 !
-!--       Cloud water content boundary condition at the bottom boundary:
-!--       Dirichlet (fixed surface rain water content).
+!--       Cloud water content boundary condition at the bottom boundary: Dirichlet (fixed surface
+!--       rain water content).
           DO  i = nxl, nxr
              DO  j = nys, nyn
@@ -981,6 +960,6 @@
 
 !
-!--       Rain water content boundary condition at the bottom boundary:
-!--       Dirichlet (fixed surface rain water content).
+!--       Rain water content boundary condition at the bottom boundary: Dirichlet (fixed surface
+!--       rain water content).
           DO  i = nxl, nxr
              DO  j = nys, nyn
@@ -1003,6 +982,6 @@
 
 !
-!--       Cloud drop concentration boundary condition at the bottom boundary:
-!--       Dirichlet (fixed surface cloud drop concentration).
+!--       Cloud drop concentration boundary condition at the bottom boundary: Dirichlet (fixed
+!--       surface cloud drop concentration).
           DO  i = nxl, nxr
              DO  j = nys, nyn
@@ -1025,6 +1004,6 @@
 
 !
-!--       Rain drop concentration boundary condition at the bottom boundary:
-!--       Dirichlet (fixed surface rain drop concentration).
+!--       Rain drop concentration boundary condition at the bottom boundary: Dirichlet (fixed
+!--       surface rain drop concentration).
           DO  i = nxl, nxr
              DO  j = nys, nyn
@@ -1060,8 +1039,7 @@
        ELSE
 
-          WRITE( message_string, * ) 'no vertical boundary condi',                &
-                                   'tion for variable "', sk_char, '"'
+          WRITE( message_string, * ) 'no vertical boundary condition for variable "', sk_char, '"'
           CALL message( 'advec_s_bc', 'PA0158', 1, 2, 0, 6, 0 )
-         
+
        ENDIF
 
@@ -1090,17 +1068,12 @@
 !
 !--    Allocate temporary arrays
-       ALLOCATE( a0(nzb:nzt+1,nys:nyn),   a1(nzb:nzt+1,nys:nyn),                  &
-                 a2(nzb:nzt+1,nys:nyn),   a12(nzb:nzt+1,nys:nyn),                 &
-                 a22(nzb:nzt+1,nys:nyn),  immb(nzb+1:nzt,nys:nyn),                &
-                 imme(nzb+1:nzt,nys:nyn), impb(nzb+1:nzt,nys:nyn),                &
-                 impe(nzb+1:nzt,nys:nyn), ipmb(nzb+1:nzt,nys:nyn),                &
-                 ipme(nzb+1:nzt,nys:nyn), ippb(nzb+1:nzt,nys:nyn),                &
-                 ippe(nzb+1:nzt,nys:nyn), m1(nzb-1:nzt+2,nys:nyn),                &
-                 sw(nzb:nzt+1,nys:nyn)                                            &
+       ALLOCATE( a0(nzb:nzt+1,nys:nyn),   a1(nzb:nzt+1,nys:nyn),   a2(nzb:nzt+1,nys:nyn),          &
+                 a12(nzb:nzt+1,nys:nyn),  a22(nzb:nzt+1,nys:nyn),  immb(nzb+1:nzt,nys:nyn),        &
+                 imme(nzb+1:nzt,nys:nyn), impb(nzb+1:nzt,nys:nyn), impe(nzb+1:nzt,nys:nyn),        &
+                 ipmb(nzb+1:nzt,nys:nyn), ipme(nzb+1:nzt,nys:nyn), ippb(nzb+1:nzt,nys:nyn),        &
+                 ippe(nzb+1:nzt,nys:nyn), m1(nzb-1:nzt+2,nys:nyn), sw(nzb:nzt+1,nys:nyn)           &
                )
        imme = 0.0_wp; impe = 0.0_wp; ipme = 0.0_wp; ippe = 0.0_wp
 
-!
-!--    Outer loop of all i
        DO  i = nxl, nxr
 
@@ -1110,15 +1083,16 @@
              DO  k = nzb, nzt+1
                 a12(k,j) = 0.5_wp * ( sk_p(k+1,j,i) - sk_p(k-1,j,i) )
-                a22(k,j) = 0.5_wp * ( sk_p(k+1,j,i) - 2.0_wp * sk_p(k,j,i)        &
-                                                    + sk_p(k-1,j,i) )
-                a0(k,j) = ( 9.0_wp * sk_p(k+2,j,i)    - 116.0_wp * sk_p(k+1,j,i)  &
-                            + 2134.0_wp * sk_p(k,j,i) - 116.0_wp * sk_p(k-1,j,i)  &
-                            + 9.0_wp * sk_p(k-2,j,i) ) * f1920
-                a1(k,j) = ( -5.0_wp   * sk_p(k+2,j,i) + 34.0_wp * sk_p(k+1,j,i)   &
-                            - 34.0_wp * sk_p(k-1,j,i) + 5.0_wp  * sk_p(k-2,j,i)   &
+                a22(k,j) = 0.5_wp * ( sk_p(k+1,j,i) - 2.0_wp * sk_p(k,j,i) + sk_p(k-1,j,i) )
+                a0(k,j) = ( 9.0_wp * sk_p(k+2,j,i)    - 116.0_wp * sk_p(k+1,j,i)                   &
+                            + 2134.0_wp * sk_p(k,j,i) - 116.0_wp * sk_p(k-1,j,i)                   &
+                            + 9.0_wp * sk_p(k-2,j,i)                                               &
+                          ) * f1920
+                a1(k,j) = ( -5.0_wp   * sk_p(k+2,j,i) + 34.0_wp * sk_p(k+1,j,i)                    &
+                            - 34.0_wp * sk_p(k-1,j,i) + 5.0_wp  * sk_p(k-2,j,i)                    &
                           ) * f48
-                a2(k,j) = ( -3.0_wp * sk_p(k+2,j,i) + 36.0_wp * sk_p(k+1,j,i)     &
-                            - 66.0_wp * sk_p(k,j,i) + 36.0_wp * sk_p(k-1,j,i)     &
-                            - 3.0_wp * sk_p(k-2,j,i) ) * f48
+                a2(k,j) = ( -3.0_wp * sk_p(k+2,j,i) + 36.0_wp * sk_p(k+1,j,i)                      &
+                            - 66.0_wp * sk_p(k,j,i) + 36.0_wp * sk_p(k-1,j,i)                      &
+                            - 3.0_wp * sk_p(k-2,j,i)                                               &
+                          ) * f48
              ENDDO
           ENDDO
@@ -1126,5 +1100,4 @@
 !
 !--       Fluxes using the Bott scheme
-!--       *VOCL LOOP,UNROLL(2)
           DO  j = nys, nyn
              DO  k = nzb+1, nzt
@@ -1133,9 +1106,9 @@
                 cipf = 1.0_wp - 2.0_wp * cip
                 cimf = 1.0_wp - 2.0_wp * cim
-                ip   =   a0(k,j)   * f2  * ( 1.0_wp - cipf )                      &
-                       + a1(k,j)   * f8  * ( 1.0_wp - cipf*cipf )                 &
+                ip   =   a0(k,j)   * f2  * ( 1.0_wp - cipf )                                       &
+                       + a1(k,j)   * f8  * ( 1.0_wp - cipf*cipf )                                  &
                        + a2(k,j)   * f24 * ( 1.0_wp - cipf*cipf*cipf )
-                im   =   a0(k+1,j) * f2  * ( 1.0_wp - cimf )                      &
-                       - a1(k+1,j) * f8  * ( 1.0_wp - cimf*cimf )                 &
+                im   =   a0(k+1,j) * f2  * ( 1.0_wp - cimf )                                       &
+                       - a1(k+1,j) * f8  * ( 1.0_wp - cimf*cimf )                                  &
                        + a2(k+1,j) * f24 * ( 1.0_wp - cimf*cimf*cimf )
                 ip   = MAX( ip, 0.0_wp )
@@ -1148,9 +1121,9 @@
                 cipf = 1.0_wp - 2.0_wp * cip
                 cimf = 1.0_wp - 2.0_wp * cim
-                ip   =   a0(k-1,j) * f2  * ( 1.0_wp - cipf )                      &
-                       + a1(k-1,j) * f8  * ( 1.0_wp - cipf*cipf )                 &
+                ip   =   a0(k-1,j) * f2  * ( 1.0_wp - cipf )                                       &
+                       + a1(k-1,j) * f8  * ( 1.0_wp - cipf*cipf )                                  &
                        + a2(k-1,j) * f24 * ( 1.0_wp - cipf*cipf*cipf )
-                im   =   a0(k,j)   * f2  * ( 1.0_wp - cimf )                      &
-                       - a1(k,j)   * f8  * ( 1.0_wp - cimf*cimf )                 &
+                im   =   a0(k,j)   * f2  * ( 1.0_wp - cimf )                                       &
+                       - a1(k,j)   * f8  * ( 1.0_wp - cimf*cimf )                                  &
                        + a2(k,j)   * f24 * ( 1.0_wp - cimf*cimf*cimf )
                 ip   = MAX( ip, 0.0_wp )
@@ -1183,9 +1156,9 @@
           DO  j = nys, nyn
              DO  k = nzb, nzt+1
-                m2 = 2.0_wp * ABS( a1(k,j) - a12(k,j) ) /                         &
+                m2 = 2.0_wp * ABS( a1(k,j) - a12(k,j) ) /                                          &
                      MAX( ABS( a1(k,j) + a12(k,j) ), 1E-35_wp )
                 IF ( ABS( a1(k,j) + a12(k,j) ) < fmax(2) )  m2 = 0.0_wp
 
-                m3 = 2.0_wp * ABS( a2(k,j) - a22(k,j) ) /                         &
+                m3 = 2.0_wp * ABS( a2(k,j) - a22(k,j) ) /                                          &
                      MAX( ABS( a2(k,j) + a22(k,j) ), 1E-35_wp )
                 IF ( ABS( a2(k,j) + a22(k,j) ) < fmax(1) )  m3 = 0.0_wp
@@ -1195,9 +1168,7 @@
                 IF ( m1(k,j) == -1.0_wp )  t2 = 0.12_wp
 
-!--             *VOCL STMT,IF(10)
-                IF ( m1(k-1,j) == 1.0_wp .OR. m1(k,j) == 1.0_wp                   &
-                     .OR. m1(k+1,j) == 1.0_wp .OR.  m2 > t2  .OR.  m3 > t2  .OR.  &
-                     ( m1(k,j) > t1  .AND.  m1(k-1,j) /= -1.0_wp  .AND.           &
-                       m1(k,j) /= -1.0_wp  .AND.  m1(k+1,j) /= -1.0_wp )          &
+                IF ( m1(k-1,j) == 1.0_wp  .OR.  m1(k,j) == 1.0_wp  .OR.  m1(k+1,j) == 1.0_wp  .OR. &
+                     m2 > t2  .OR.  m3 > t2  .OR.  ( m1(k,j) > t1 .AND.  m1(k-1,j) /= -1.0_wp .AND.&
+                     m1(k,j) /= -1.0_wp .AND. m1(k+1,j) /= -1.0_wp )                               &
                    )  sw(k,j) = 1.0_wp
              ENDDO
@@ -1210,5 +1181,4 @@
              DO  k = nzb+1, nzt
 
-!--             *VOCL STMT,IF(10)
                 IF ( sw(k,j) == 1.0_wp )  THEN
                    snenn = sk_p(k+1,j,i) - sk_p(k-1,j,i)
@@ -1222,9 +1192,9 @@
                    cip =  MAX( 0.0_wp, w(k,j,i) * dt_3d * ddzw(k) )
 
-                   ippe(k,j) = sk_p(k-1,j,i) * cip + snenn * (                    &
-                               aex(ix) * cip + bex(ix) / dex(ix) * (              &
-                               eex(ix) -                                          &
-                               EXP( dex(ix)*0.5_wp * ( 1.0_wp - 2.0_wp * cip ) )  &
-                                                                   )              &
+                   ippe(k,j) = sk_p(k-1,j,i) * cip + snenn * (                                     &
+                               aex(ix) * cip + bex(ix) / dex(ix) * (                               &
+                               eex(ix) -                                                           &
+                               EXP( dex(ix)*0.5_wp * ( 1.0_wp - 2.0_wp * cip ) )                   &
+                                                                   )                               &
                                                              )
                    IF ( sterm == 0.0001_wp )  ippe(k,j) = sk_p(k,j,i) * cip
@@ -1241,9 +1211,9 @@
                    cim = -MIN( 0.0_wp, w(k-1,j,i) * dt_3d * ddzw(k) )
 
-                   imme(k,j) = sk_p(k+1,j,i) * cim + snenn * (                    &
-                               aex(ix) * cim + bex(ix) / dex(ix) * (              &
-                               eex(ix) -                                          &
-                               EXP( dex(ix)*0.5_wp * ( 1.0_wp - 2.0_wp * cim ) ) &
-                                                                   )              &
+                   imme(k,j) = sk_p(k+1,j,i) * cim + snenn * (                                     &
+                               aex(ix) * cim + bex(ix) / dex(ix) * (                               &
+                               eex(ix) -                                                           &
+                               EXP( dex(ix)*0.5_wp * ( 1.0_wp - 2.0_wp * cim ) )                   &
+                                                                   )                               &
                                                              )
                    IF ( sterm == 0.0001_wp )  imme(k,j) = sk_p(k,j,i) * cim
@@ -1251,5 +1221,4 @@
                 ENDIF
 
-!--             *VOCL STMT,IF(10)
                 IF ( sw(k+1,j) == 1.0_wp )  THEN
                    snenn = sk_p(k,j,i) - sk_p(k+2,j,i)
@@ -1263,9 +1232,9 @@
                    cim = -MIN( 0.0_wp, w(k,j,i) * dt_3d * ddzw(k) )
 
-                   impe(k,j) = sk_p(k+2,j,i) * cim + snenn * (                    &
-                               aex(ix) * cim + bex(ix) / dex(ix) * (              &
-                               eex(ix) -                                          &
-                               EXP( dex(ix)*0.5_wp * ( 1.0_wp - 2.0_wp * cim ) )  &
-                                                                   )              &
+                   impe(k,j) = sk_p(k+2,j,i) * cim + snenn * (                                     &
+                               aex(ix) * cim + bex(ix) / dex(ix) * (                               &
+                               eex(ix) -                                                           &
+                               EXP( dex(ix)*0.5_wp * ( 1.0_wp - 2.0_wp * cim ) )                   &
+                                                                   )                               &
                                                              )
                    IF ( sterm == 0.0001_wp )  impe(k,j) = sk_p(k+1,j,i) * cim
@@ -1273,5 +1242,4 @@
                 ENDIF
 
-!--             *VOCL STMT,IF(10)
                 IF ( sw(k-1,j) == 1.0_wp )  THEN
                    snenn = sk_p(k,j,i) - sk_p(k-2,j,i)
@@ -1285,9 +1253,9 @@
                    cip = MAX( 0.0_wp, w(k-1,j,i) * dt_3d * ddzw(k) )
 
-                   ipme(k,j) = sk_p(k-2,j,i) * cip + snenn * (                    &
-                               aex(ix) * cip + bex(ix) / dex(ix) * (              &
-                               eex(ix) -                                          &
-                               EXP( dex(ix)*0.5_wp * ( 1.0_wp - 2.0_wp * cip ) )  &
-                                                                   )              &
+                   ipme(k,j) = sk_p(k-2,j,i) * cip + snenn * (                                     &
+                               aex(ix) * cip + bex(ix) / dex(ix) * (                               &
+                               eex(ix) -                                                           &
+                               EXP( dex(ix)*0.5_wp * ( 1.0_wp - 2.0_wp * cip ) )                   &
+                                                                   )                               &
                                                              )
                    IF ( sterm == 0.0001_wp )  ipme(k,j) = sk_p(k-1,j,i) * cip
@@ -1303,8 +1271,8 @@
           DO  j = nys, nyn
              DO  k = nzb+1, nzt
-                fplus  = ( 1.0_wp - sw(k,j)   ) * ippb(k,j) + sw(k,j)   * ippe(k,j) &
-                       - ( 1.0_wp - sw(k+1,j) ) * impb(k,j) - sw(k+1,j) * impe(k,j)
-                fminus = ( 1.0_wp - sw(k-1,j) ) * ipmb(k,j) + sw(k-1,j) * ipme(k,j) &
-                       - ( 1.0_wp - sw(k,j)   ) * immb(k,j) - sw(k,j)   * imme(k,j)
+                fplus  = ( 1.0_wp - sw(k,j) )     * ippb(k,j) + sw(k,j)   * ippe(k,j)              &
+                         - ( 1.0_wp - sw(k+1,j) ) * impb(k,j) - sw(k+1,j) * impe(k,j)
+                fminus = ( 1.0_wp - sw(k-1,j) )   * ipmb(k,j) + sw(k-1,j) * ipme(k,j)              &
+                         - ( 1.0_wp - sw(k,j) )   * immb(k,j) - sw(k,j)   * imme(k,j)
                 tendcy = fplus - fminus
 !
@@ -1315,9 +1283,8 @@
 !--             Density correction because of possible remaining divergences
                 d_new = d(k,j,i) - ( w(k,j,i) - w(k-1,j,i) ) * dt_3d * ddzw(k)
-                sk_p(k,j,i) = ( ( 1.0_wp + d(k,j,i) ) * sk_p(k,j,i) - tendcy ) /  &
-                              ( 1.0_wp + d_new )
+                sk_p(k,j,i) = ( ( 1.0_wp + d(k,j,i) ) * sk_p(k,j,i) - tendcy ) / ( 1.0_wp + d_new )
 !
 !--             Store heat flux for subsequent statistics output.
-!--             array m1 is here used as temporary storage
+!--             Array m1 is here used as temporary storage
                 m1(k,j) = fplus / dt_3d * dzw(k)
              ENDDO
@@ -1325,11 +1292,10 @@
 
 !
-!--       Sum up heat flux in order to order to obtain horizontal averages
+!--       Sum up heat flux in order to obtain horizontal averages
           IF ( sk_char == 'pt' )  THEN
              DO  sr = 0, statistic_regions
                 DO  j = nys, nyn
                    DO  k = nzb+1, nzt
-                      sums_wsts_bc_l(k,sr) = sums_wsts_bc_l(k,sr) +               &
-                                             m1(k,j) * rmask(j,i,sr)
+                      sums_wsts_bc_l(k,sr) = sums_wsts_bc_l(k,sr) + m1(k,j) * rmask(j,i,sr)
                    ENDDO
                 ENDDO
@@ -1337,5 +1303,8 @@
           ENDIF
 
-       ENDDO   ! End of the advection in z-direction
+       ENDDO
+!
+!--    End of the advection in z-direction
+
        CALL cpu_log( log_point_s(12), 'advec_s_bc:exp', 'continue' )
        CALL cpu_log( log_point_s(12), 'advec_s_bc:exp', 'stop' )
@@ -1343,6 +1312,5 @@
 !
 !--    Deallocate temporary arrays
-       DEALLOCATE( a0, a1, a2, a12, a22, immb, imme, impb, impe, ipmb, ipme,      &
-                   ippb, ippe, m1, sw )
+       DEALLOCATE( a0, a1, a2, a12, a22, immb, imme, impb, impe, ipmb, ipme, ippb, ippe, m1, sw )
 
 !
Index: /palm/trunk/SOURCE/advec_s_pw.f90
===================================================================
--- /palm/trunk/SOURCE/advec_s_pw.f90	(revision 4487)
+++ /palm/trunk/SOURCE/advec_s_pw.f90	(revision 4488)
@@ -1,20 +1,19 @@
 !> @file advec_s_pw.f90
-!------------------------------------------------------------------------------!
+!--------------------------------------------------------------------------------------------------!
 ! This file is part of the PALM model system.
 !
-! PALM is free software: you can redistribute it and/or modify it under the
-! terms of the GNU General Public License as published by the Free Software
-! Foundation, either version 3 of the License, or (at your option) any later
-! version.
+! PALM is free software: you can redistribute it and/or modify it under the terms of the GNU General
+! Public License as published by the Free Software Foundation, either version 3 of the License, or
+! (at your option) any later version.
 !
-! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
-! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
-! A PARTICULAR PURPOSE.  See the GNU General Public License for more details.
+! PALM is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the
+! implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General
+! Public License for more details.
 !
-! You should have received a copy of the GNU General Public License along with
-! PALM. If not, see <http://www.gnu.org/licenses/>.
+! You should have received a copy of the GNU General Public License along with PALM. If not, see
+! <http://www.gnu.org/licenses/>.
 !
 ! Copyright 1997-2020 Leibniz Universitaet Hannover
-!------------------------------------------------------------------------------!
+!--------------------------------------------------------------------------------------------------!
 !
 ! Current revisions:
@@ -25,12 +24,15 @@
 ! -----------------
 ! $Id$
+! file re-formatted to follow the PALM coding standard
+!
+! 4360 2020-01-07 11:25:50Z suehring
 ! Corrected "Former revisions" section
-! 
+!
 ! 3927 2019-04-23 13:24:29Z raasch
 ! pointer attribute removed from scalar 3d-array for performance reasons
-! 
+!
 ! 3665 2019-01-10 08:28:24Z raasch
 ! unused variables removed
-! 
+!
 ! 3655 2019-01-07 16:51:22Z knoop
 ! nopointer option removed
@@ -42,13 +44,12 @@
 ! Description:
 ! ------------
-!> Advection term for scalar variables using the Piacsek and Williams scheme
-!> (form C3). Contrary to PW itself, for reasons of accuracy their scheme is
-!> slightly modified as follows: the values of those scalars that are used for
-!> the computation of the flux divergence are reduced by the value of the
-!> relevant scalar at the location where the difference is computed (sk(k,j,i)).
+!> Advection term for scalar variables using the Piacsek and Williams scheme (form C3). Contrary to
+!> PW itself, for reasons of accuracy their scheme is slightly modified as follows: the values of
+!> those scalars that are used for the computation of the flux divergence are reduced by the value
+!> of the relevant scalar at the location where the difference is computed (sk(k,j,i)).
 !> NOTE: at the first grid point above the surface computation still takes place!
-!------------------------------------------------------------------------------!
+!--------------------------------------------------------------------------------------------------!
  MODULE advec_s_pw_mod
- 
+
 
     PRIVATE
@@ -59,136 +60,136 @@
        MODULE PROCEDURE advec_s_pw_ij
     END INTERFACE
- 
+
  CONTAINS
 
 
-!------------------------------------------------------------------------------!
+!--------------------------------------------------------------------------------------------------!
 ! Description:
 ! ------------
 !> Call for all grid points
-!------------------------------------------------------------------------------!
-    SUBROUTINE advec_s_pw( sk )
+!--------------------------------------------------------------------------------------------------!
+ SUBROUTINE advec_s_pw( sk )
 
-       USE arrays_3d,                                                          &
-           ONLY:  dd2zu, tend, u, u_stokes_zu, v, v_stokes_zu, w
+    USE arrays_3d,                                                                                 &
+        ONLY:  dd2zu, tend, u, u_stokes_zu, v, v_stokes_zu, w
 
-       USE control_parameters,                                                 &
-           ONLY:  u_gtrans, v_gtrans
+    USE control_parameters,                                                                        &
+        ONLY:  u_gtrans, v_gtrans
 
-       USE grid_variables,                                                     &
-           ONLY:  ddx, ddy
+    USE grid_variables,                                                                            &
+        ONLY:  ddx, ddy
 
-       USE indices,                                                            &
-           ONLY:  nxl, nxlg, nxr, nxrg, nyn, nyng, nys, nysg, nzb, nzt
+    USE indices,                                                                                   &
+        ONLY:  nxl, nxlg, nxr, nxrg, nyn, nyng, nys, nysg, nzb, nzt
 
-       USE kinds
+    USE kinds
 
 
-       IMPLICIT NONE
+    IMPLICIT NONE
 
-       INTEGER(iwp) ::  i !< grid index along x-direction
-       INTEGER(iwp) ::  j !< grid index along y-direction
-       INTEGER(iwp) ::  k !< grid index along z-direction
+    INTEGER(iwp) ::  i !< grid index along x-direction
+    INTEGER(iwp) ::  j !< grid index along y-direction
+    INTEGER(iwp) ::  k !< grid index along z-direction
 
-       REAL(wp)     ::  gu  !< local additional advective velocity
-       REAL(wp)     ::  gv  !< local additional advective velocity
+    REAL(wp)     ::  gu  !< local additional advective velocity
+    REAL(wp)     ::  gv  !< local additional advective velocity
 
-       REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg) ::  sk
- 
+    REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg) ::  sk
 
-       DO  i = nxl, nxr
-          DO  j = nys, nyn
-             DO  k = nzb+1, nzt
+
+    DO  i = nxl, nxr
+       DO  j = nys, nyn
+          DO  k = nzb+1, nzt
 
 !
-!--             Galilean transformation + Stokes drift velocity (ocean only)
-                gu = u_gtrans - u_stokes_zu(k)
-                gv = v_gtrans - v_stokes_zu(k)
+!--          Galilean transformation + Stokes drift velocity (ocean only)
+             gu = u_gtrans - u_stokes_zu(k)
+             gv = v_gtrans - v_stokes_zu(k)
 
-                tend(k,j,i) = tend(k,j,i) +                                    &
-                                     ( -0.5_wp * ( ( u(k,j,i+1) - gu       ) * &
-                                                   ( sk(k,j,i+1) - sk(k,j,i) ) &
-                                                 - ( u(k,j,i)   - gu       ) * &
-                                                   ( sk(k,j,i-1) - sk(k,j,i) ) &
-                                                 ) * ddx                       &
-                                       -0.5_wp * ( ( v(k,j+1,i) - gv       ) * &
-                                                   ( sk(k,j+1,i) - sk(k,j,i) ) &
-                                                 - ( v(k,j,i)   - gv       ) * &
-                                                   ( sk(k,j-1,i) - sk(k,j,i) ) &
-                                                 ) * ddy                       &
-                                       -         (   w(k,j,i)                * &
-                                                   ( sk(k+1,j,i) - sk(k,j,i) ) &
-                                                 -   w(k-1,j,i)              * &
-                                                   ( sk(k-1,j,i) - sk(k,j,i) ) &
-                                                 ) * dd2zu(k)                  &
-                                      )
-             ENDDO
+             tend(k,j,i) = tend(k,j,i) +                                                           &
+                                  ( -0.5_wp * ( ( u(k,j,i+1)  - gu        ) *                      &
+                                                ( sk(k,j,i+1) - sk(k,j,i) )                        &
+                                              - ( u(k,j,i)    - gu        ) *                      &
+                                                ( sk(k,j,i-1) - sk(k,j,i) )                        &
+                                              ) * ddx                                              &
+                                    -0.5_wp * ( ( v(k,j+1,i)  - gv        ) *                      &
+                                                ( sk(k,j+1,i) - sk(k,j,i) )                        &
+                                              - ( v(k,j,i)    - gv        ) *                      &
+                                                ( sk(k,j-1,i) - sk(k,j,i) )                        &
+                                              ) * ddy                                              &
+                                    -         (   w(k,j,i)                  *                      &
+                                                ( sk(k+1,j,i) - sk(k,j,i) )                        &
+                                              -   w(k-1,j,i)                *                      &
+                                                ( sk(k-1,j,i) - sk(k,j,i) )                        &
+                                              ) * dd2zu(k)                                         &
+                                   )
           ENDDO
        ENDDO
+    ENDDO
 
-    END SUBROUTINE advec_s_pw
+ END SUBROUTINE advec_s_pw
 
 
-!------------------------------------------------------------------------------!
+!--------------------------------------------------------------------------------------------------!
 ! Description:
 ! ------------
 !> Call for grid point i,j
-!------------------------------------------------------------------------------!
-    SUBROUTINE advec_s_pw_ij( i, j, sk )
+!--------------------------------------------------------------------------------------------------!
+ SUBROUTINE advec_s_pw_ij( i, j, sk )
 
-       USE arrays_3d,                                                          &
-           ONLY:  dd2zu, tend, u, u_stokes_zu, v, v_stokes_zu, w
+    USE arrays_3d,                                                                                 &
+        ONLY:  dd2zu, tend, u, u_stokes_zu, v, v_stokes_zu, w
 
-       USE control_parameters,                                                 &
-           ONLY:  u_gtrans, v_gtrans
+    USE control_parameters,                                                                        &
+        ONLY:  u_gtrans, v_gtrans
 
-       USE grid_variables,                                                     &
-           ONLY:  ddx, ddy
+    USE grid_variables,                                                                            &
+        ONLY:  ddx, ddy
 
-       USE indices,                                                            &
-           ONLY:  nxlg, nxrg, nyng, nysg, nzb, nzt
+    USE indices,                                                                                   &
+        ONLY:  nxlg, nxrg, nyng, nysg, nzb, nzt
 
-       USE kinds
+    USE kinds
 
 
-       IMPLICIT NONE
+    IMPLICIT NONE
 
-       INTEGER(iwp) ::  i !< grid index along x-direction
-       INTEGER(iwp) ::  j !< grid index along y-direction
-       INTEGER(iwp) ::  k !< grid index along z-direction
+    INTEGER(iwp) ::  i !< grid index along x-direction
+    INTEGER(iwp) ::  j !< grid index along y-direction
+    INTEGER(iwp) ::  k !< grid index along z-direction
 
-       REAL(wp)     ::  gu  !< local additional advective velocity
-       REAL(wp)     ::  gv  !< local additional advective velocity
+    REAL(wp)     ::  gu  !< local additional advective velocity
+    REAL(wp)     ::  gv  !< local additional advective velocity
 
-       REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg) ::  sk
+    REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg) ::  sk
 
 
-       DO  k = nzb+1, nzt
+    DO  k = nzb+1, nzt
 
 !
-!--       Galilean transformation + Stokes drift velocity (ocean only)
-          gu = u_gtrans - u_stokes_zu(k)
-          gv = v_gtrans - v_stokes_zu(k)
+!--    Galilean transformation + Stokes drift velocity (ocean only)
+       gu = u_gtrans - u_stokes_zu(k)
+       gv = v_gtrans - v_stokes_zu(k)
 
-          tend(k,j,i) = tend(k,j,i) +                                          &
-                                    ( -0.5_wp * ( ( u(k,j,i+1) - gu        ) * &
-                                                  ( sk(k,j,i+1) - sk(k,j,i) )  &
-                                                - ( u(k,j,i)   - gu        ) * &
-                                                  ( sk(k,j,i-1) - sk(k,j,i) )  &
-                                                ) * ddx                        &
-                                      -0.5_wp * ( ( v(k,j+1,i) - gv       ) *  &
-                                                  ( sk(k,j+1,i) - sk(k,j,i) )  &
-                                                - ( v(k,j,i)   - gv       ) *  &
-                                                  ( sk(k,j-1,i) - sk(k,j,i) )  &
-                                                ) * ddy                        &
-                                      -         (   w(k,j,i)                *  &
-                                                  ( sk(k+1,j,i) - sk(k,j,i) )  &
-                                                -   w(k-1,j,i)              *  &
-                                                  ( sk(k-1,j,i) - sk(k,j,i) )  &
-                                                ) * dd2zu(k)                   &
-                                    )
-       ENDDO
+       tend(k,j,i) = tend(k,j,i) +                                                                 &
+                                 ( -0.5_wp * ( ( u(k,j,i+1)  - gu        ) *                       &
+                                               ( sk(k,j,i+1) - sk(k,j,i) )                         &
+                                             - ( u(k,j,i)    - gu        ) *                       &
+                                               ( sk(k,j,i-1) - sk(k,j,i) )                         &
+                                             ) * ddx                                               &
+                                   -0.5_wp * ( ( v(k,j+1,i)  - gv        ) *                       &
+                                               ( sk(k,j+1,i) - sk(k,j,i) )                         &
+                                             - ( v(k,j,i)    - gv        ) *                       &
+                                               ( sk(k,j-1,i) - sk(k,j,i) )                         &
+                                             ) * ddy                                               &
+                                   -         (   w(k,j,i)                  *                       &
+                                               ( sk(k+1,j,i) - sk(k,j,i) )                         &
+                                             -   w(k-1,j,i)                *                       &
+                                               ( sk(k-1,j,i) - sk(k,j,i) )                         &
+                                             ) * dd2zu(k)                                          &
+                                 )
+    ENDDO
 
-    END SUBROUTINE advec_s_pw_ij
+ END SUBROUTINE advec_s_pw_ij
 
  END MODULE advec_s_pw_mod
Index: /palm/trunk/SOURCE/advec_s_up.f90
===================================================================
--- /palm/trunk/SOURCE/advec_s_up.f90	(revision 4487)
+++ /palm/trunk/SOURCE/advec_s_up.f90	(revision 4488)
@@ -1,20 +1,19 @@
 !> @file advec_s_up.f90
-!------------------------------------------------------------------------------!
+!--------------------------------------------------------------------------------------------------!
 ! This file is part of the PALM model system.
 !
-! PALM is free software: you can redistribute it and/or modify it under the
-! terms of the GNU General Public License as published by the Free Software
-! Foundation, either version 3 of the License, or (at your option) any later
-! version.
+! PALM is free software: you can redistribute it and/or modify it under the terms of the GNU General
+! Public License as published by the Free Software Foundation, either version 3 of the License, or
+! (at your option) any later version.
 !
-! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
-! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
-! A PARTICULAR PURPOSE.  See the GNU General Public License for more details.
+! PALM is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the
+! implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General
+! Public License for more details.
 !
-! You should have received a copy of the GNU General Public License along with
-! PALM. If not, see <http://www.gnu.org/licenses/>.
+! You should have received a copy of the GNU General Public License along with PALM. If not, see
+! <http://www.gnu.org/licenses/>.
 !
 ! Copyright 1997-2020 Leibniz Universitaet Hannover
-!------------------------------------------------------------------------------!
+!--------------------------------------------------------------------------------------------------!
 !
 ! Current revisions:
@@ -25,12 +24,15 @@
 ! -----------------
 ! $Id$
+! file re-formatted to follow the PALM coding standard
+!
+! 4360 2020-01-07 11:25:50Z suehring
 ! Corrected "Former revisions" section
-! 
+!
 ! 3927 2019-04-23 13:24:29Z raasch
 ! pointer attribute removed from scalar 3d-array for performance reasons
-! 
+!
 ! 3665 2019-01-10 08:28:24Z raasch
 ! unused variables removed
-! 
+!
 ! 3655 2019-01-07 16:51:22Z knoop
 ! nopointer option removed
@@ -45,7 +47,7 @@
 !> NOTE: vertical advection at k=1 still has wrong grid spacing for w>0!
 !>       The same problem occurs for all topography boundaries!
-!------------------------------------------------------------------------------!
+!--------------------------------------------------------------------------------------------------!
  MODULE advec_s_up_mod
- 
+
 
     PRIVATE
@@ -60,150 +62,138 @@
 
 
-!------------------------------------------------------------------------------!
+!--------------------------------------------------------------------------------------------------!
 ! Description:
 ! ------------
 !> Call for all grid points
-!------------------------------------------------------------------------------!
-    SUBROUTINE advec_s_up( sk )
+!--------------------------------------------------------------------------------------------------!
+ SUBROUTINE advec_s_up( sk )
 
-       USE arrays_3d,                                                          &
-           ONLY:  ddzu, tend, u, v, w
+    USE arrays_3d,                                                                                 &
+        ONLY:  ddzu, tend, u, v, w
 
-       USE control_parameters,                                                 &
-           ONLY:  u_gtrans, v_gtrans
+    USE control_parameters,                                                                        &
+        ONLY:  u_gtrans, v_gtrans
 
-       USE grid_variables,                                                     &
-           ONLY:  ddx, ddy
+    USE grid_variables,                                                                            &
+        ONLY:  ddx, ddy
 
-       USE indices,                                                            &
-           ONLY:  nxl, nxlg, nxr, nxrg, nyn, nyng, nys, nysg, nzb, nzt
+    USE indices,                                                                                   &
+        ONLY:  nxl, nxlg, nxr, nxrg, nyn, nyng, nys, nysg, nzb, nzt
 
-       USE kinds
+    USE kinds
 
 
-       IMPLICIT NONE
+    IMPLICIT NONE
 
-       INTEGER(iwp) ::  i !< grid index along x-direction
-       INTEGER(iwp) ::  j !< grid index along y-direction
-       INTEGER(iwp) ::  k !< grid index along z-direction
+    INTEGER(iwp) ::  i !< grid index along x-direction
+    INTEGER(iwp) ::  j !< grid index along y-direction
+    INTEGER(iwp) ::  k !< grid index along z-direction
 
-       REAL(wp) ::  ukomp !< advection velocity along x-direction 
-       REAL(wp) ::  vkomp !< advection velocity along y-direction 
-       REAL(wp) ::  wkomp !< advection velocity along z-direction 
+    REAL(wp) ::  ukomp !< advection velocity along x-direction
+    REAL(wp) ::  vkomp !< advection velocity along y-direction
+    REAL(wp) ::  wkomp !< advection velocity along z-direction
 
-       REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg) ::  sk !< treated scalar
+    REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg) ::  sk !< treated scalar
 
 
-       DO  i = nxl, nxr
-          DO  j = nys, nyn
-             DO  k = nzb+1, nzt
+    DO  i = nxl, nxr
+       DO  j = nys, nyn
+          DO  k = nzb+1, nzt
 !
-!--             x-direction
-                ukomp = 0.5_wp * ( u(k,j,i) + u(k,j,i+1) ) - u_gtrans
-                IF ( ukomp > 0.0_wp )  THEN
-                   tend(k,j,i) = tend(k,j,i) - ukomp *                         &
-                                         ( sk(k,j,i) - sk(k,j,i-1) ) * ddx
-                ELSE
-                   tend(k,j,i) = tend(k,j,i) - ukomp *                         &
-                                         ( sk(k,j,i+1) - sk(k,j,i) ) * ddx
-                ENDIF
+!--          x-direction
+             ukomp = 0.5_wp * ( u(k,j,i) + u(k,j,i+1) ) - u_gtrans
+             IF ( ukomp > 0.0_wp )  THEN
+                tend(k,j,i) = tend(k,j,i) - ukomp * ( sk(k,j,i) - sk(k,j,i-1) ) * ddx
+             ELSE
+                tend(k,j,i) = tend(k,j,i) - ukomp * ( sk(k,j,i+1) - sk(k,j,i) ) * ddx
+             ENDIF
 !
-!--             y-direction
-                vkomp = 0.5_wp * ( v(k,j,i) + v(k,j+1,i) ) - v_gtrans
-                IF ( vkomp > 0.0_wp )  THEN
-                   tend(k,j,i) = tend(k,j,i) - vkomp *                         &
-                                         ( sk(k,j,i) - sk(k,j-1,i) ) * ddy
-                ELSE
-                   tend(k,j,i) = tend(k,j,i) - vkomp *                         &
-                                         ( sk(k,j+1,i) - sk(k,j,i) ) * ddy
-                ENDIF
+!--          y-direction
+             vkomp = 0.5_wp * ( v(k,j,i) + v(k,j+1,i) ) - v_gtrans
+             IF ( vkomp > 0.0_wp )  THEN
+                tend(k,j,i) = tend(k,j,i) - vkomp * ( sk(k,j,i) - sk(k,j-1,i) ) * ddy
+             ELSE
+                tend(k,j,i) = tend(k,j,i) - vkomp * ( sk(k,j+1,i) - sk(k,j,i) ) * ddy
+             ENDIF
 !
-!--             z-direction
-                wkomp = 0.5_wp * ( w(k,j,i) + w(k-1,j,i) )
-                IF ( wkomp > 0.0_wp )  THEN
-                   tend(k,j,i) = tend(k,j,i) - wkomp *                         &
-                                         ( sk(k,j,i) - sk(k-1,j,i) ) * ddzu(k)
-                ELSE
-                   tend(k,j,i) = tend(k,j,i) - wkomp *                         &
-                                         ( sk(k+1,j,i)-sk(k,j,i) ) * ddzu(k+1)
-                ENDIF
+!--          z-direction
+             wkomp = 0.5_wp * ( w(k,j,i) + w(k-1,j,i) )
+             IF ( wkomp > 0.0_wp )  THEN
+                tend(k,j,i) = tend(k,j,i) - wkomp * ( sk(k,j,i) - sk(k-1,j,i) ) * ddzu(k)
+             ELSE
+                tend(k,j,i) = tend(k,j,i) - wkomp * ( sk(k+1,j,i) - sk(k,j,i) ) * ddzu(k+1)
+             ENDIF
 
-             ENDDO
           ENDDO
        ENDDO
+    ENDDO
 
-    END SUBROUTINE advec_s_up
+ END SUBROUTINE advec_s_up
 
 
-!------------------------------------------------------------------------------!
+!--------------------------------------------------------------------------------------------------!
 ! Description:
 ! ------------
 !> Call for grid point i,j
-!------------------------------------------------------------------------------!
-    SUBROUTINE advec_s_up_ij( i, j, sk )
+!--------------------------------------------------------------------------------------------------!
+ SUBROUTINE advec_s_up_ij( i, j, sk )
 
-       USE arrays_3d,                                                          &
-           ONLY:  ddzu, tend, u, v, w
+    USE arrays_3d,                                                                                 &
+        ONLY:  ddzu, tend, u, v, w
 
-       USE control_parameters,                                                 &
-           ONLY:  u_gtrans, v_gtrans
+    USE control_parameters,                                                                        &
+        ONLY:  u_gtrans, v_gtrans
 
-       USE grid_variables,                                                     &
-           ONLY:  ddx, ddy
+    USE grid_variables,                                                                            &
+        ONLY:  ddx, ddy
 
-       USE indices,                                                            &
-           ONLY:  nxlg, nxrg, nyng, nysg, nzb, nzt
+    USE indices,                                                                                   &
+        ONLY:  nxlg, nxrg, nyng, nysg, nzb, nzt
 
-       USE kinds
+    USE kinds
 
 
-       IMPLICIT NONE
+    IMPLICIT NONE
 
-       INTEGER(iwp) ::  i !< grid index along x-direction
-       INTEGER(iwp) ::  j !< grid index along y-direction
-       INTEGER(iwp) ::  k !< grid index along z-direction
+    INTEGER(iwp) ::  i !< grid index along x-direction
+    INTEGER(iwp) ::  j !< grid index along y-direction
+    INTEGER(iwp) ::  k !< grid index along z-direction
 
-       REAL(wp) ::  ukomp !< advection velocity along x-direction 
-       REAL(wp) ::  vkomp !< advection velocity along y-direction 
-       REAL(wp) ::  wkomp !< advection velocity along z-direction 
-       
-       REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg) ::  sk !< treated scalar
+    REAL(wp) ::  ukomp !< advection velocity along x-direction
+    REAL(wp) ::  vkomp !< advection velocity along y-direction
+    REAL(wp) ::  wkomp !< advection velocity along z-direction
+
+    REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg) ::  sk !< treated scalar
 
 
-       DO  k = nzb+1, nzt
+    DO  k = nzb+1, nzt
 !
-!--       x-direction
-          ukomp = 0.5_wp * ( u(k,j,i) + u(k,j,i+1) ) - u_gtrans
-          IF ( ukomp > 0.0_wp )  THEN
-             tend(k,j,i) = tend(k,j,i) - ukomp *                               &
-                                         ( sk(k,j,i) - sk(k,j,i-1) ) * ddx
-          ELSE
-             tend(k,j,i) = tend(k,j,i) - ukomp *                               &
-                                         ( sk(k,j,i+1) - sk(k,j,i) ) * ddx
-          ENDIF
+!--    x-direction
+       ukomp = 0.5_wp * ( u(k,j,i) + u(k,j,i+1) ) - u_gtrans
+       IF ( ukomp > 0.0_wp )  THEN
+          tend(k,j,i) = tend(k,j,i) - ukomp * ( sk(k,j,i) - sk(k,j,i-1) ) * ddx
+       ELSE
+          tend(k,j,i) = tend(k,j,i) - ukomp * ( sk(k,j,i+1) - sk(k,j,i) ) * ddx
+       ENDIF
 !
-!--       y-direction
-          vkomp = 0.5_wp * ( v(k,j,i) + v(k,j+1,i) ) - v_gtrans
-          IF ( vkomp > 0.0_wp )  THEN
-             tend(k,j,i) = tend(k,j,i) - vkomp *                               &
-                                         ( sk(k,j,i) - sk(k,j-1,i) ) * ddy
-          ELSE
-             tend(k,j,i) = tend(k,j,i) - vkomp *                               &
-                                         ( sk(k,j+1,i) - sk(k,j,i) ) * ddy
-          ENDIF
+!--    y-direction
+       vkomp = 0.5_wp * ( v(k,j,i) + v(k,j+1,i) ) - v_gtrans
+       IF ( vkomp > 0.0_wp )  THEN
+          tend(k,j,i) = tend(k,j,i) - vkomp * ( sk(k,j,i) - sk(k,j-1,i) ) * ddy
+       ELSE
+          tend(k,j,i) = tend(k,j,i) - vkomp * ( sk(k,j+1,i) - sk(k,j,i) ) * ddy
+       ENDIF
 !
-!--       z-direction
-          wkomp = 0.5_wp * ( w(k,j,i) + w(k-1,j,i) )
-          IF ( wkomp > 0.0_wp )  THEN
-             tend(k,j,i) = tend(k,j,i) - wkomp *                               &
-                                         ( sk(k,j,i) - sk(k-1,j,i) ) * ddzu(k)
-          ELSE
-             tend(k,j,i) = tend(k,j,i) - wkomp *                               &
-                                         ( sk(k+1,j,i)-sk(k,j,i) ) * ddzu(k+1)
-          ENDIF
+!--    z-direction
+       wkomp = 0.5_wp * ( w(k,j,i) + w(k-1,j,i) )
+       IF ( wkomp > 0.0_wp )  THEN
+          tend(k,j,i) = tend(k,j,i) - wkomp * ( sk(k,j,i) - sk(k-1,j,i) ) * ddzu(k)
+       ELSE
+          tend(k,j,i) = tend(k,j,i) - wkomp * ( sk(k+1,j,i) - sk(k,j,i) ) * ddzu(k+1)
+       ENDIF
 
-       ENDDO
+    ENDDO
 
-    END SUBROUTINE advec_s_up_ij
+ END SUBROUTINE advec_s_up_ij
 
  END MODULE advec_s_up_mod
Index: /palm/trunk/SOURCE/advec_u_pw.f90
===================================================================
--- /palm/trunk/SOURCE/advec_u_pw.f90	(revision 4487)
+++ /palm/trunk/SOURCE/advec_u_pw.f90	(revision 4488)
@@ -1,20 +1,19 @@
 !> @file advec_u_pw.f90
-!------------------------------------------------------------------------------!
+!--------------------------------------------------------------------------------------------------!
 ! This file is part of the PALM model system.
 !
-! PALM is free software: you can redistribute it and/or modify it under the
-! terms of the GNU General Public License as published by the Free Software
-! Foundation, either version 3 of the License, or (at your option) any later
-! version.
+! PALM is free software: you can redistribute it and/or modify it under the terms of the GNU General
+! Public License as published by the Free Software Foundation, either version 3 of the License, or
+! (at your option) any later version.
 !
-! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
-! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
-! A PARTICULAR PURPOSE.  See the GNU General Public License for more details.
+! PALM is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the
+! implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General
+! Public License for more details.
 !
-! You should have received a copy of the GNU General Public License along with
-! PALM. If not, see <http://www.gnu.org/licenses/>.
+! You should have received a copy of the GNU General Public License along with PALM. If not, see
+! <http://www.gnu.org/licenses/>.
 !
 ! Copyright 1997-2020 Leibniz Universitaet Hannover
-!------------------------------------------------------------------------------!
+!--------------------------------------------------------------------------------------------------!
 !
 ! Current revisions:
@@ -25,6 +24,9 @@
 ! -----------------
 ! $Id$
+! file re-formatted to follow the PALM coding standard
+!
+! 4360 2020-01-07 11:25:50Z suehring
 ! Corrected "Former revisions" section
-! 
+!
 ! 3655 2019-01-07 16:51:22Z knoop
 ! variables documented
@@ -37,10 +39,10 @@
 ! ------------
 !> Advection term for u velocity-component using Piacsek and Williams.
-!> Vertical advection at the first grid point above the surface is done with
-!> normal centred differences, because otherwise no information from the surface
-!> would be communicated upwards due to w=0 at K=nzb.
-!------------------------------------------------------------------------------!
+!> Vertical advection at the first grid point above the surface is done with normal centred
+!> differences, because otherwise no information from the surface would be communicated upwards due
+!> to w=0 at K=nzb.
+!--------------------------------------------------------------------------------------------------!
  MODULE advec_u_pw_mod
- 
+
 
     PRIVATE
@@ -51,106 +53,104 @@
        MODULE PROCEDURE advec_u_pw_ij
     END INTERFACE advec_u_pw
- 
+
  CONTAINS
 
 
-!------------------------------------------------------------------------------!
+!--------------------------------------------------------------------------------------------------!
 ! Description:
 ! ------------
 !> Call for all grid points
-!------------------------------------------------------------------------------!
-    SUBROUTINE advec_u_pw
+!--------------------------------------------------------------------------------------------------!
+ SUBROUTINE advec_u_pw
 
-       USE arrays_3d,                                                          &
-           ONLY:  ddzw, tend, u, v, w
+    USE arrays_3d,                                                                                 &
+        ONLY:  ddzw, tend, u, v, w
 
-       USE control_parameters,                                                 &
-           ONLY:  u_gtrans, v_gtrans
+    USE control_parameters,                                                                        &
+        ONLY:  u_gtrans, v_gtrans
 
-       USE grid_variables,                                                     &
-           ONLY:  ddx, ddy
+    USE grid_variables,                                                                            &
+        ONLY:  ddx, ddy
 
-       USE indices,                                                            &
-           ONLY:  nxlu, nxr, nyn, nys, nzb, nzt
+    USE indices,                                                                                   &
+        ONLY:  nxlu, nxr, nyn, nys, nzb, nzt
 
-       USE kinds
+    USE kinds
 
 
-       IMPLICIT NONE
+    IMPLICIT NONE
 
-       INTEGER(iwp) ::  i !< grid index along x-direction
-       INTEGER(iwp) ::  j !< grid index along y-direction
-       INTEGER(iwp) ::  k !< grid index along z-direction
-       
-       REAL(wp)    ::  gu !< Galilei-transformation velocity along x
-       REAL(wp)    ::  gv !< Galilei-transformation velocity along y
- 
-       gu = 2.0_wp * u_gtrans
-       gv = 2.0_wp * v_gtrans
-       DO  i = nxlu, nxr
-          DO  j = nys, nyn
-             DO  k = nzb+1, nzt
-                tend(k,j,i) = tend(k,j,i) - 0.25_wp * (                        &
-                         ( u(k,j,i+1) * ( u(k,j,i+1) + u(k,j,i) - gu )         &
-                         - u(k,j,i-1) * ( u(k,j,i) + u(k,j,i-1) - gu ) ) * ddx &
-                       + ( u(k,j+1,i) * ( v(k,j+1,i) + v(k,j+1,i-1) - gv )     &
-                         - u(k,j-1,i) * ( v(k,j,i) + v(k,j,i-1) - gv ) ) * ddy &
-                       + ( u(k+1,j,i) * ( w(k,j,i) + w(k,j,i-1) )              &
-                         - u(k-1,j,i) * ( w(k-1,j,i) + w(k-1,j,i-1) ) )        &
-                                                                  * ddzw(k)    &
-                                                      )
-             ENDDO
+    INTEGER(iwp) ::  i !< grid index along x-direction
+    INTEGER(iwp) ::  j !< grid index along y-direction
+    INTEGER(iwp) ::  k !< grid index along z-direction
+
+    REAL(wp)     ::  gu !< Galilei-transformation velocity along x
+    REAL(wp)     ::  gv !< Galilei-transformation velocity along y
+
+    gu = 2.0_wp * u_gtrans
+    gv = 2.0_wp * v_gtrans
+    DO  i = nxlu, nxr
+       DO  j = nys, nyn
+          DO  k = nzb+1, nzt
+             tend(k,j,i) = tend(k,j,i) - 0.25_wp * (                                               &
+                           ( u(k,j,i+1)   * ( u(k,j,i+1) + u(k,j,i) - gu )                         &
+                           - u(k,j,i-1)   * ( u(k,j,i)   + u(k,j,i-1) - gu ) ) * ddx               &
+                           + ( u(k,j+1,i) * ( v(k,j+1,i) + v(k,j+1,i-1) - gv )                     &
+                           - u(k,j-1,i)   * ( v(k,j,i)   + v(k,j,i-1) - gv ) ) * ddy               &
+                           + ( u(k+1,j,i) * ( w(k,j,i)   + w(k,j,i-1) )                            &
+                           - u(k-1,j,i)   * ( w(k-1,j,i) + w(k-1,j,i-1) ) )    * ddzw(k)           &
+                                                   )
           ENDDO
        ENDDO
+    ENDDO
 
-    END SUBROUTINE advec_u_pw
+ END SUBROUTINE advec_u_pw
 
 
-!------------------------------------------------------------------------------!
+!--------------------------------------------------------------------------------------------------!
 ! Description:
 ! ------------
 !> Call for grid point i,j
-!------------------------------------------------------------------------------!
-    SUBROUTINE advec_u_pw_ij( i, j )
+!--------------------------------------------------------------------------------------------------!
+ SUBROUTINE advec_u_pw_ij( i, j )
 
-       USE arrays_3d,                                                          &
-           ONLY:  ddzw, tend, u, v, w
+    USE arrays_3d,                                                                                 &
+        ONLY:  ddzw, tend, u, v, w
 
-       USE control_parameters,                                                 &
-           ONLY:  u_gtrans, v_gtrans
+    USE control_parameters,                                                                        &
+        ONLY:  u_gtrans, v_gtrans
 
-       USE grid_variables,                                                     &
-           ONLY:  ddx, ddy
+    USE grid_variables,                                                                            &
+        ONLY:  ddx, ddy
 
-       USE indices,                                                            &
-           ONLY:  nzb, nzt
+    USE indices,                                                                                   &
+        ONLY:  nzb, nzt
 
-       USE kinds
+    USE kinds
 
 
-       IMPLICIT NONE
+    IMPLICIT NONE
 
-       INTEGER(iwp) ::  i !< grid index along x-direction
-       INTEGER(iwp) ::  j !< grid index along y-direction
-       INTEGER(iwp) ::  k !< grid index along z-direction
-       
-       REAL(wp)    ::  gu !< Galilei-transformation velocity along x
-       REAL(wp)    ::  gv !< Galilei-transformation velocity along y
+    INTEGER(iwp) ::  i !< grid index along x-direction
+    INTEGER(iwp) ::  j !< grid index along y-direction
+    INTEGER(iwp) ::  k !< grid index along z-direction
 
-       gu = 2.0_wp * u_gtrans
-       gv = 2.0_wp * v_gtrans
-       DO  k = nzb+1, nzt
-          tend(k,j,i) = tend(k,j,i) - 0.25_wp * (                              &
-                         ( u(k,j,i+1) * ( u(k,j,i+1) + u(k,j,i) - gu )         &
-                         - u(k,j,i-1) * ( u(k,j,i) + u(k,j,i-1) - gu ) ) * ddx &
-                       + ( u(k,j+1,i) * ( v(k,j+1,i) + v(k,j+1,i-1) - gv )     &
-                         - u(k,j-1,i) * ( v(k,j,i) + v(k,j,i-1) - gv ) ) * ddy &
-                       + ( u(k+1,j,i) * ( w(k,j,i) + w(k,j,i-1) )              &
-                         - u(k-1,j,i) * ( w(k-1,j,i) + w(k-1,j,i-1) ) )        &
-                                                                  * ddzw(k)    &
-                                                )
-       ENDDO
+    REAL(wp)     ::  gu !< Galilei-transformation velocity along x
+    REAL(wp)     ::  gv !< Galilei-transformation velocity along y
 
-    END SUBROUTINE advec_u_pw_ij
+    gu = 2.0_wp * u_gtrans
+    gv = 2.0_wp * v_gtrans
+    DO  k = nzb+1, nzt
+       tend(k,j,i) = tend(k,j,i) - 0.25_wp * (                                                     &
+                     ( u(k,j,i+1)   * ( u(k,j,i+1) + u(k,j,i) - gu )                               &
+                     - u(k,j,i-1)   * ( u(k,j,i)   + u(k,j,i-1) - gu ) ) * ddx                     &
+                     + ( u(k,j+1,i) * ( v(k,j+1,i) + v(k,j+1,i-1) - gv )                           &
+                     - u(k,j-1,i)   * ( v(k,j,i)   + v(k,j,i-1) - gv ) ) * ddy                     &
+                     + ( u(k+1,j,i) * ( w(k,j,i)   + w(k,j,i-1) )                                  &
+                     - u(k-1,j,i)   * ( w(k-1,j,i) + w(k-1,j,i-1) ) )    * ddzw(k)                 &
+                                             )
+    ENDDO
+
+ END SUBROUTINE advec_u_pw_ij
 
  END MODULE advec_u_pw_mod
Index: /palm/trunk/SOURCE/advec_u_up.f90
===================================================================
--- /palm/trunk/SOURCE/advec_u_up.f90	(revision 4487)
+++ /palm/trunk/SOURCE/advec_u_up.f90	(revision 4488)
@@ -1,20 +1,19 @@
 !> @file advec_u_up.f90
-!------------------------------------------------------------------------------!
+!--------------------------------------------------------------------------------------------------!
 ! This file is part of the PALM model system.
 !
-! PALM is free software: you can redistribute it and/or modify it under the
-! terms of the GNU General Public License as published by the Free Software
-! Foundation, either version 3 of the License, or (at your option) any later
-! version.
+! PALM is free software: you can redistribute it and/or modify it under the terms of the GNU General
+! Public License as published by the Free Software Foundation, either version 3 of the License, or
+! (at your option) any later version.
 !
-! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
-! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
-! A PARTICULAR PURPOSE.  See the GNU General Public License for more details.
+! PALM is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the
+! implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General
+! Public License for more details.
 !
-! You should have received a copy of the GNU General Public License along with
-! PALM. If not, see <http://www.gnu.org/licenses/>.
+! You should have received a copy of the GNU General Public License along with PALM. If not, see
+! <http://www.gnu.org/licenses/>.
 !
 ! Copyright 1997-2020 Leibniz Universitaet Hannover
-!------------------------------------------------------------------------------!
+!--------------------------------------------------------------------------------------------------!
 !
 ! Current revisions:
@@ -25,6 +24,9 @@
 ! -----------------
 ! $Id$
+! file re-formatted to follow the PALM coding standard
+!
+! 4360 2020-01-07 11:25:50Z suehring
 ! Corrected "Former revisions" section
-! 
+!
 ! 3655 2019-01-07 16:51:22Z knoop
 ! variables documented
@@ -39,7 +41,7 @@
 !> NOTE: vertical advection at k=1 still has wrong grid spacing for w>0!
 !>       The same problem occurs for all topography boundaries!
-!------------------------------------------------------------------------------!
+!--------------------------------------------------------------------------------------------------!
  MODULE advec_u_up_mod
- 
+
 
     PRIVATE
@@ -54,149 +56,134 @@
 
 
-!------------------------------------------------------------------------------!
+!--------------------------------------------------------------------------------------------------!
 ! Description:
 ! ------------
 !> Call for all grid points
-!------------------------------------------------------------------------------!
-    SUBROUTINE advec_u_up
+!--------------------------------------------------------------------------------------------------!
+ SUBROUTINE advec_u_up
 
-       USE arrays_3d,                                                          &
-           ONLY:  ddzu, tend, u, v, w
+    USE arrays_3d,                                                                                 &
+        ONLY:  ddzu, tend, u, v, w
 
-       USE control_parameters,                                                 &
-           ONLY:  u_gtrans, v_gtrans
+    USE control_parameters,                                                                        &
+        ONLY:  u_gtrans, v_gtrans
 
-       USE grid_variables,                                                     &
-           ONLY:  ddx, ddy
+    USE grid_variables,                                                                            &
+        ONLY:  ddx, ddy
 
-       USE indices,                                                            &
-           ONLY:  nxlu, nxr, nyn, nys, nzb, nzt
+    USE indices,                                                                                   &
+        ONLY:  nxlu, nxr, nyn, nys, nzb, nzt
 
-       USE kinds
+    USE kinds
 
 
-       IMPLICIT NONE
+    IMPLICIT NONE
 
-       INTEGER(iwp) ::  i !< grid index along x-direction
-       INTEGER(iwp) ::  j !< grid index along y-direction
-       INTEGER(iwp) ::  k !< grid index along z-direction
+    INTEGER(iwp) ::  i !< grid index along x-direction
+    INTEGER(iwp) ::  j !< grid index along y-direction
+    INTEGER(iwp) ::  k !< grid index along z-direction
 
-       REAL(wp) ::  ukomp !< advection velocity along x-direction 
-       REAL(wp) ::  vkomp !< advection velocity along y-direction 
-       REAL(wp) ::  wkomp !< advection velocity along z-direction 
+    REAL(wp) ::  ukomp !< advection velocity along x-direction
+    REAL(wp) ::  vkomp !< advection velocity along y-direction
+    REAL(wp) ::  wkomp !< advection velocity along z-direction
 
-       
-       DO  i = nxlu, nxr
-          DO  j = nys, nyn
-             DO  k = nzb+1, nzt
+
+    DO  i = nxlu, nxr
+       DO  j = nys, nyn
+          DO  k = nzb+1, nzt
 !
-!--             x-direction
-                ukomp = u(k,j,i) - u_gtrans
-                IF ( ukomp > 0.0_wp )  THEN
-                   tend(k,j,i) = tend(k,j,i) - ukomp *                         &
-                                         ( u(k,j,i) - u(k,j,i-1) ) * ddx
-                ELSE
-                   tend(k,j,i) = tend(k,j,i) - ukomp *                         &
-                                          ( u(k,j,i+1) - u(k,j,i) ) * ddx
-                ENDIF
+!--          x-direction
+             ukomp = u(k,j,i) - u_gtrans
+             IF ( ukomp > 0.0_wp )  THEN
+                tend(k,j,i) = tend(k,j,i) - ukomp * ( u(k,j,i) - u(k,j,i-1) ) * ddx
+             ELSE
+                tend(k,j,i) = tend(k,j,i) - ukomp * ( u(k,j,i+1) - u(k,j,i) ) * ddx
+             ENDIF
 !
-!--             y-direction
-                vkomp = 0.25_wp * ( v(k,j,i)   + v(k,j+1,i) +                  &
-                                 v(k,j,i-1) + v(k,j+1,i-1) ) - v_gtrans
-                IF ( vkomp > 0.0_wp )  THEN
-                   tend(k,j,i) = tend(k,j,i) - vkomp *                         &
-                                         ( u(k,j,i) - u(k,j-1,i) ) * ddy
-                ELSE
-                   tend(k,j,i) = tend(k,j,i) - vkomp *                         &
-                                         ( u(k,j+1,i) - u(k,j,i) ) * ddy
-                ENDIF
+!--          y-direction
+             vkomp = 0.25_wp * ( v(k,j,i) + v(k,j+1,i) + v(k,j,i-1) + v(k,j+1,i-1) ) - v_gtrans
+             IF ( vkomp > 0.0_wp )  THEN
+                tend(k,j,i) = tend(k,j,i) - vkomp * ( u(k,j,i) - u(k,j-1,i) ) * ddy
+             ELSE
+                tend(k,j,i) = tend(k,j,i) - vkomp * ( u(k,j+1,i) - u(k,j,i) ) * ddy
+             ENDIF
 !
-!--             z-direction
-                wkomp = 0.25_wp * ( w(k,j,i)   + w(k-1,j,i) +                  &
-                                 w(k,j,i-1) + w(k-1,j,i-1) )
-                IF ( wkomp > 0.0_wp )  THEN
-                   tend(k,j,i) = tend(k,j,i) - wkomp *                         &
-                                         ( u(k,j,i) - u(k-1,j,i) ) * ddzu(k)
-                ELSE
-                   tend(k,j,i) = tend(k,j,i) - wkomp *                         &
-                                         ( u(k+1,j,i) - u(k,j,i) ) * ddzu(k+1)
-                ENDIF
+!--          z-direction
+             wkomp = 0.25_wp * ( w(k,j,i) + w(k-1,j,i) + w(k,j,i-1) + w(k-1,j,i-1) )
+             IF ( wkomp > 0.0_wp )  THEN
+                tend(k,j,i) = tend(k,j,i) - wkomp * ( u(k,j,i) - u(k-1,j,i) ) * ddzu(k)
+             ELSE
+                tend(k,j,i) = tend(k,j,i) - wkomp * ( u(k+1,j,i) - u(k,j,i) ) * ddzu(k+1)
+             ENDIF
 
-             ENDDO
           ENDDO
        ENDDO
+    ENDDO
 
-    END SUBROUTINE advec_u_up
+ END SUBROUTINE advec_u_up
 
 
-!------------------------------------------------------------------------------!
+!--------------------------------------------------------------------------------------------------!
 ! Description:
 ! ------------
 !> Call for grid point i,j
-!------------------------------------------------------------------------------!
-    SUBROUTINE advec_u_up_ij( i, j )
+!--------------------------------------------------------------------------------------------------!
+ SUBROUTINE advec_u_up_ij( i, j )
 
-       USE arrays_3d,                                                          &
-           ONLY:  ddzu, tend, u, v, w
+    USE arrays_3d,                                                                                 &
+        ONLY:  ddzu, tend, u, v, w
 
-       USE control_parameters,                                                 &
-           ONLY:  u_gtrans, v_gtrans
+    USE control_parameters,                                                                        &
+        ONLY:  u_gtrans, v_gtrans
 
-       USE grid_variables,                                                     &
-           ONLY:  ddx, ddy
+    USE grid_variables,                                                                            &
+        ONLY:  ddx, ddy
 
-       USE indices,                                                            &
-           ONLY:  nzb, nzt
+    USE indices,                                                                                   &
+        ONLY:  nzb, nzt
 
-       USE kinds
+    USE kinds
 
 
-       IMPLICIT NONE
+    IMPLICIT NONE
 
-       INTEGER(iwp) ::  i !< grid index along x-direction
-       INTEGER(iwp) ::  j !< grid index along y-direction
-       INTEGER(iwp) ::  k !< grid index along z-direction
+    INTEGER(iwp) ::  i !< grid index along x-direction
+    INTEGER(iwp) ::  j !< grid index along y-direction
+    INTEGER(iwp) ::  k !< grid index along z-direction
 
-       REAL(wp) ::  ukomp !< advection velocity along x-direction 
-       REAL(wp) ::  vkomp !< advection velocity along y-direction 
-       REAL(wp) ::  wkomp !< advection velocity along z-direction 
+    REAL(wp) ::  ukomp !< advection velocity along x-direction
+    REAL(wp) ::  vkomp !< advection velocity along y-direction
+    REAL(wp) ::  wkomp !< advection velocity along z-direction
 
 
-       DO  k = nzb+1, nzt
+    DO  k = nzb+1, nzt
 !
-!--       x-direction
-          ukomp = u(k,j,i) - u_gtrans
-          IF ( ukomp > 0.0_wp )  THEN
-             tend(k,j,i) = tend(k,j,i) - ukomp *                               &
-                                         ( u(k,j,i) - u(k,j,i-1) ) * ddx
-          ELSE
-             tend(k,j,i) = tend(k,j,i) - ukomp *                               &
-                                         ( u(k,j,i+1) - u(k,j,i) ) * ddx
-          ENDIF
+!--    x-direction
+       ukomp = u(k,j,i) - u_gtrans
+       IF ( ukomp > 0.0_wp )  THEN
+          tend(k,j,i) = tend(k,j,i) - ukomp * ( u(k,j,i) - u(k,j,i-1) ) * ddx
+       ELSE
+          tend(k,j,i) = tend(k,j,i) - ukomp * ( u(k,j,i+1) - u(k,j,i) ) * ddx
+       ENDIF
 !
-!--       y-direction
-          vkomp = 0.25_wp * ( v(k,j,i) + v(k,j+1,i) + v(k,j,i-1) + v(k,j+1,i-1) &
-                         ) - v_gtrans
-          IF ( vkomp > 0.0_wp )  THEN
-             tend(k,j,i) = tend(k,j,i) - vkomp *                               &
-                                         ( u(k,j,i) - u(k,j-1,i) ) * ddy
-          ELSE
-             tend(k,j,i) = tend(k,j,i) - vkomp *                               &
-                                         ( u(k,j+1,i) - u(k,j,i) ) * ddy
-          ENDIF
+!--    y-direction
+       vkomp = 0.25_wp * ( v(k,j,i) + v(k,j+1,i) + v(k,j,i-1) + v(k,j+1,i-1) ) - v_gtrans
+       IF ( vkomp > 0.0_wp )  THEN
+          tend(k,j,i) = tend(k,j,i) - vkomp * ( u(k,j,i) - u(k,j-1,i) ) * ddy
+       ELSE
+          tend(k,j,i) = tend(k,j,i) - vkomp * ( u(k,j+1,i) - u(k,j,i) ) * ddy
+       ENDIF
 !
-!--       z-direction
-          wkomp = 0.25_wp * ( w(k,j,i) + w(k-1,j,i) + w(k,j,i-1) + w(k-1,j,i-1) )
-          IF ( wkomp > 0.0_wp )  THEN
-             tend(k,j,i) = tend(k,j,i) - wkomp *                               &
-                                         ( u(k,j,i) - u(k-1,j,i) ) * ddzu(k)
-          ELSE
-             tend(k,j,i) = tend(k,j,i) - wkomp *                               &
-                                         ( u(k+1,j,i) - u(k,j,i) ) * ddzu(k+1)
-          ENDIF
+!--    z-direction
+       wkomp = 0.25_wp * ( w(k,j,i) + w(k-1,j,i) + w(k,j,i-1) + w(k-1,j,i-1) )
+       IF ( wkomp > 0.0_wp )  THEN
+          tend(k,j,i) = tend(k,j,i) - wkomp * ( u(k,j,i) - u(k-1,j,i) ) * ddzu(k)
+       ELSE
+          tend(k,j,i) = tend(k,j,i) - wkomp * ( u(k+1,j,i) - u(k,j,i) ) * ddzu(k+1)
+       ENDIF
 
-       ENDDO
+    ENDDO
 
-    END SUBROUTINE advec_u_up_ij
+ END SUBROUTINE advec_u_up_ij
 
  END MODULE advec_u_up_mod
Index: /palm/trunk/SOURCE/advec_v_pw.f90
===================================================================
--- /palm/trunk/SOURCE/advec_v_pw.f90	(revision 4487)
+++ /palm/trunk/SOURCE/advec_v_pw.f90	(revision 4488)
@@ -1,21 +1,19 @@
 !> @file advec_v_pw.f90
-!------------------------------------------------------------------------------!
-!------------------------------------------------------------------------------!
+!--------------------------------------------------------------------------------------------------!
 ! This file is part of the PALM model system.
 !
-! PALM is free software: you can redistribute it and/or modify it under the
-! terms of the GNU General Public License as published by the Free Software
-! Foundation, either version 3 of the License, or (at your option) any later
-! version.
+! PALM is free software: you can redistribute it and/or modify it under the terms of the GNU General
+! Public License as published by the Free Software Foundation, either version 3 of the License, or
+! (at your option) any later version.
 !
-! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
-! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
-! A PARTICULAR PURPOSE.  See the GNU General Public License for more details.
+! PALM is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the
+! implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General
+! Public License for more details.
 !
-! You should have received a copy of the GNU General Public License along with
-! PALM. If not, see <http://www.gnu.org/licenses/>.
+! You should have received a copy of the GNU General Public License along with PALM. If not, see
+! <http://www.gnu.org/licenses/>.
 !
 ! Copyright 1997-2020 Leibniz Universitaet Hannover
-!------------------------------------------------------------------------------!
+!--------------------------------------------------------------------------------------------------!
 !
 ! Current revisions:
@@ -26,4 +24,7 @@
 ! -----------------
 ! $Id$
+! file re-formatted to follow the PALM coding standard
+!
+! 4360 2020-01-07 11:25:50Z suehring
 ! Corrected "Former revisions" section
 ! 
@@ -38,8 +39,8 @@
 ! ------------
 !> Advection term for v velocity-component using Piacsek and Williams.
-!> Vertical advection at the first grid point above the surface is done with
-!> normal centred differences, because otherwise no information from the surface
-!> would be communicated upwards due to w=0 at K=nzb.
-!------------------------------------------------------------------------------!
+!> Vertical advection at the first grid point above the surface is done with normal centred
+!> differences, because otherwise no information from the surface would be communicated upwards due
+!> to w=0 at K=nzb.
+!--------------------------------------------------------------------------------------------------!
  MODULE advec_v_pw_mod
  
@@ -56,104 +57,102 @@
 
 
-!------------------------------------------------------------------------------!
+!--------------------------------------------------------------------------------------------------!
 ! Description:
 ! ------------
 !> Call for all grid points
-!------------------------------------------------------------------------------!
-    SUBROUTINE advec_v_pw
+!--------------------------------------------------------------------------------------------------!
+ SUBROUTINE advec_v_pw
 
-       USE arrays_3d,                                                          &
-           ONLY:  ddzw, tend, u, v, w
+    USE arrays_3d,                                                                                 &
+        ONLY:  ddzw, tend, u, v, w
 
-       USE control_parameters,                                                 &
-           ONLY:  u_gtrans, v_gtrans
+    USE control_parameters,                                                                        &
+        ONLY:  u_gtrans, v_gtrans
 
-       USE grid_variables,                                                     &
-           ONLY:  ddx, ddy
+    USE grid_variables,                                                                            &
+        ONLY:  ddx, ddy
 
-       USE indices,                                                            &
-           ONLY:  nxl, nxr, nyn, nysv, nzb, nzt
+    USE indices,                                                                                   &
+        ONLY:  nxl, nxr, nyn, nysv, nzb, nzt
 
-       USE kinds
+    USE kinds
 
 
-       IMPLICIT NONE
+    IMPLICIT NONE
 
-       INTEGER(iwp) ::  i !< grid index along x-direction
-       INTEGER(iwp) ::  j !< grid index along y-direction
-       INTEGER(iwp) ::  k !< grid index along z-direction
-       
-       REAL(wp)    ::  gu !< Galilei-transformation velocity along x
-       REAL(wp)    ::  gv !< Galilei-transformation velocity along y
- 
+    INTEGER(iwp) ::  i !< grid index along x-direction
+    INTEGER(iwp) ::  j !< grid index along y-direction
+    INTEGER(iwp) ::  k !< grid index along z-direction
+    
+    REAL(wp)     ::  gu !< Galilei-transformation velocity along x
+    REAL(wp)     ::  gv !< Galilei-transformation velocity along y
 
-       gu = 2.0_wp * u_gtrans
-       gv = 2.0_wp * v_gtrans
-       DO  i = nxl, nxr
-          DO  j = nysv, nyn
-             DO  k = nzb+1, nzt
-                tend(k,j,i) = tend(k,j,i) - 0.25_wp * (                        &
-                         ( v(k,j,i+1) * ( u(k,j-1,i+1) + u(k,j,i+1) - gu )     &
-                         - v(k,j,i-1) * ( u(k,j-1,i) + u(k,j,i) - gu ) ) * ddx &
-                       + ( v(k,j+1,i) * ( v(k,j+1,i) + v(k,j,i) - gv )         &
-                         - v(k,j-1,i) * ( v(k,j,i) + v(k,j-1,i) - gv ) ) * ddy &
-                       + ( v(k+1,j,i) * ( w(k,j-1,i) + w(k,j,i) )              &
-                         - v(k-1,j,i) * ( w(k-1,j-1,i) + w(k-1,j,i) ) )        &
-                                                                  * ddzw(k)    &
-                                                      )
-             ENDDO
+
+    gu = 2.0_wp * u_gtrans
+    gv = 2.0_wp * v_gtrans
+    DO  i = nxl, nxr
+       DO  j = nysv, nyn
+          DO  k = nzb+1, nzt
+             tend(k,j,i) = tend(k,j,i) - 0.25_wp * (                                               &
+                           ( v(k,j,i+1)   * ( u(k,j-1,i+1) + u(k,j,i+1) - gu )                     &
+                           - v(k,j,i-1)   * ( u(k,j-1,i)   + u(k,j,i)   - gu ) ) * ddx             &
+                           + ( v(k,j+1,i) * ( v(k,j+1,i)   + v(k,j,i)   - gv )                     &
+                           - v(k,j-1,i)   * ( v(k,j,i)     + v(k,j-1,i) - gv ) ) * ddy             &
+                           + ( v(k+1,j,i) * ( w(k,j-1,i)   + w(k,j,i) )                            &
+                           - v(k-1,j,i)   * ( w(k-1,j-1,i) + w(k-1,j,i) ) )      * ddzw(k)         &
+                                                   )
           ENDDO
        ENDDO
+    ENDDO
 
-    END SUBROUTINE advec_v_pw
+ END SUBROUTINE advec_v_pw
 
 
-!------------------------------------------------------------------------------!
+!--------------------------------------------------------------------------------------------------!
 ! Description:
 ! ------------
 !> Call for grid point i,j
-!------------------------------------------------------------------------------!
-    SUBROUTINE advec_v_pw_ij( i, j )
+!--------------------------------------------------------------------------------------------------!
+ SUBROUTINE advec_v_pw_ij( i, j )
 
-       USE arrays_3d,                                                          &
-           ONLY:  ddzw, tend, u, v, w
+    USE arrays_3d,                                                                                 &
+        ONLY:  ddzw, tend, u, v, w
 
-       USE control_parameters,                                                 &
-           ONLY:  u_gtrans, v_gtrans
+    USE control_parameters,                                                                        &
+        ONLY:  u_gtrans, v_gtrans
 
-       USE grid_variables,                                                     &
-           ONLY:  ddx, ddy
+    USE grid_variables,                                                                            &
+        ONLY:  ddx, ddy
 
-       USE indices,                                                            &
-           ONLY:  nzb, nzt
+    USE indices,                                                                                   &
+        ONLY:  nzb, nzt
 
-       USE kinds
+    USE kinds
 
 
-       IMPLICIT NONE
+    IMPLICIT NONE
 
-       INTEGER(iwp) ::  i !< grid index along x-direction
-       INTEGER(iwp) ::  j !< grid index along y-direction
-       INTEGER(iwp) ::  k !< grid index along z-direction
-       
-       REAL(wp)    ::  gu !< Galilei-transformation velocity along x
-       REAL(wp)    ::  gv !< Galilei-transformation velocity along y
+    INTEGER(iwp) ::  i !< grid index along x-direction
+    INTEGER(iwp) ::  j !< grid index along y-direction
+    INTEGER(iwp) ::  k !< grid index along z-direction
+    
+    REAL(wp)     ::  gu !< Galilei-transformation velocity along x
+    REAL(wp)     ::  gv !< Galilei-transformation velocity along y
 
 
-       gu = 2.0_wp * u_gtrans
-       gv = 2.0_wp * v_gtrans
-       DO  k = nzb+1, nzt
-          tend(k,j,i) = tend(k,j,i) - 0.25_wp * (                              &
-                         ( v(k,j,i+1) * ( u(k,j-1,i+1) + u(k,j,i+1) - gu )     &
-                         - v(k,j,i-1) * ( u(k,j-1,i) + u(k,j,i) - gu ) ) * ddx &
-                       + ( v(k,j+1,i) * ( v(k,j+1,i) + v(k,j,i) - gv )         &
-                         - v(k,j-1,i) * ( v(k,j,i) + v(k,j-1,i) - gv ) ) * ddy &
-                       + ( v(k+1,j,i) * ( w(k,j-1,i) + w(k,j,i) )              &
-                         - v(k-1,j,i) * ( w(k-1,j-1,i) + w(k-1,j,i) ) )        &
-                                                                  * ddzw(k)    &
-                                                )
-       ENDDO
- 
-    END SUBROUTINE advec_v_pw_ij
+    gu = 2.0_wp * u_gtrans
+    gv = 2.0_wp * v_gtrans
+    DO  k = nzb+1, nzt
+       tend(k,j,i) = tend(k,j,i) - 0.25_wp * (                                                     &
+                      ( v(k,j,i+1)   * ( u(k,j-1,i+1) + u(k,j,i+1) - gu )                          &
+                      - v(k,j,i-1)   * ( u(k,j-1,i)   + u(k,j,i)   - gu ) ) * ddx                  &
+                      + ( v(k,j+1,i) * ( v(k,j+1,i)   + v(k,j,i)   - gv )                          &
+                      - v(k,j-1,i)   * ( v(k,j,i)     + v(k,j-1,i) - gv ) ) * ddy                  &
+                      + ( v(k+1,j,i) * ( w(k,j-1,i)   + w(k,j,i) )                                 &
+                      - v(k-1,j,i)   * ( w(k-1,j-1,i) + w(k-1,j,i) ) )      * ddzw(k)              &
+                                             )
+    ENDDO
+
+ END SUBROUTINE advec_v_pw_ij
 
  END MODULE advec_v_pw_mod
Index: /palm/trunk/SOURCE/advec_v_up.f90
===================================================================
--- /palm/trunk/SOURCE/advec_v_up.f90	(revision 4487)
+++ /palm/trunk/SOURCE/advec_v_up.f90	(revision 4488)
@@ -1,20 +1,19 @@
 !> @file advec_v_up.f90
-!------------------------------------------------------------------------------!
+!--------------------------------------------------------------------------------------------------!
 ! This file is part of the PALM model system.
 !
-! PALM is free software: you can redistribute it and/or modify it under the
-! terms of the GNU General Public License as published by the Free Software
-! Foundation, either version 3 of the License, or (at your option) any later
-! version.
+! PALM is free software: you can redistribute it and/or modify it under the terms of the GNU General
+! Public License as published by the Free Software Foundation, either version 3 of the License, or
+! (at your option) any later version.
 !
-! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
-! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
-! A PARTICULAR PURPOSE.  See the GNU General Public License for more details.
+! PALM is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the
+! implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General
+! Public License for more details.
 !
-! You should have received a copy of the GNU General Public License along with
-! PALM. If not, see <http://www.gnu.org/licenses/>.
+! You should have received a copy of the GNU General Public License along with PALM. If not, see
+! <http://www.gnu.org/licenses/>.
 !
 ! Copyright 1997-2020 Leibniz Universitaet Hannover
-!------------------------------------------------------------------------------!
+!--------------------------------------------------------------------------------------------------!
 !
 ! Current revisions:
@@ -25,4 +24,7 @@
 ! -----------------
 ! $Id$
+! file re-formatted to follow the PALM coding standard
+!
+! 4360 2020-01-07 11:25:50Z suehring
 ! Corrected "Former revisions" section
 ! 
@@ -39,5 +41,5 @@
 !> NOTE: vertical advection at k=1 still has wrong grid spacing for w>0!
 !>       The same problem occurs for all topography boundaries!
-!------------------------------------------------------------------------------!
+!--------------------------------------------------------------------------------------------------!
  MODULE advec_v_up_mod
  
@@ -54,149 +56,134 @@
 
 
-!------------------------------------------------------------------------------!
+!--------------------------------------------------------------------------------------------------!
 ! Description:
 ! ------------
 !> Call for all grid points
-!------------------------------------------------------------------------------!
-    SUBROUTINE advec_v_up
+!--------------------------------------------------------------------------------------------------!
+ SUBROUTINE advec_v_up
 
-       USE arrays_3d,                                                          &
-           ONLY:  ddzu, tend, u, v, w
+    USE arrays_3d,                                                                                 &
+        ONLY:  ddzu, tend, u, v, w
 
-       USE control_parameters,                                                 &
-           ONLY:  u_gtrans, v_gtrans
+    USE control_parameters,                                                                        &
+        ONLY:  u_gtrans, v_gtrans
 
-       USE grid_variables,                                                     &
-           ONLY:  ddx, ddy
+    USE grid_variables,                                                                            &
+        ONLY:  ddx, ddy
 
-       USE indices,                                                            &
-           ONLY:  nxl, nxr, nyn, nysv, nzb, nzt
+    USE indices,                                                                                   &
+        ONLY:  nxl, nxr, nyn, nysv, nzb, nzt
 
-       USE kinds
+    USE kinds
 
 
-       IMPLICIT NONE
+    IMPLICIT NONE
 
-       INTEGER(iwp) ::  i !< grid index along x-direction
-       INTEGER(iwp) ::  j !< grid index along y-direction
-       INTEGER(iwp) ::  k !< grid index along z-direction
+    INTEGER(iwp) ::  i !< grid index along x-direction
+    INTEGER(iwp) ::  j !< grid index along y-direction
+    INTEGER(iwp) ::  k !< grid index along z-direction
 
-       REAL(wp) ::  ukomp !< advection velocity along x-direction 
-       REAL(wp) ::  vkomp !< advection velocity along y-direction 
-       REAL(wp) ::  wkomp !< advection velocity along z-direction       
+    REAL(wp) ::  ukomp !< advection velocity along x-direction 
+    REAL(wp) ::  vkomp !< advection velocity along y-direction 
+    REAL(wp) ::  wkomp !< advection velocity along z-direction       
 
 
-       DO  i = nxl, nxr
-          DO  j = nysv, nyn
-             DO  k = nzb+1, nzt
+    DO  i = nxl, nxr
+       DO  j = nysv, nyn
+          DO  k = nzb+1, nzt
 !
-!--             x-direction
-                ukomp = 0.25_wp * ( u(k,j,i)   + u(k,j-1,i) +                  &
-                                 u(k,j,i+1) + u(k,j-1,i+1) ) - u_gtrans
-                IF ( ukomp > 0.0_wp )  THEN
-                   tend(k,j,i) = tend(k,j,i) - ukomp *                         &
-                                         ( v(k,j,i) - v(k,j,i-1) ) * ddx
-                ELSE
-                   tend(k,j,i) = tend(k,j,i) - ukomp *                         &
-                                         ( v(k,j,i+1) - v(k,j,i) ) * ddx
-                ENDIF
+!--          x-direction
+             ukomp = 0.25_wp * ( u(k,j,i) + u(k,j-1,i) + u(k,j,i+1) + u(k,j-1,i+1) ) - u_gtrans
+             IF ( ukomp > 0.0_wp )  THEN
+                tend(k,j,i) = tend(k,j,i) - ukomp * ( v(k,j,i) - v(k,j,i-1) ) * ddx
+             ELSE
+                tend(k,j,i) = tend(k,j,i) - ukomp * ( v(k,j,i+1) - v(k,j,i) ) * ddx
+             ENDIF
 !
-!--             y-direction
-                vkomp = v(k,j,i) - v_gtrans
-                IF ( vkomp > 0.0_wp )  THEN
-                   tend(k,j,i) = tend(k,j,i) - vkomp *                         &
-                                         ( v(k,j,i) - v(k,j-1,i) ) * ddy
-                ELSE
-                   tend(k,j,i) = tend(k,j,i) - vkomp *                         &
-                                         ( v(k,j+1,i) - v(k,j,i) ) * ddy
-                ENDIF
+!--          y-direction
+             vkomp = v(k,j,i) - v_gtrans
+             IF ( vkomp > 0.0_wp )  THEN
+                tend(k,j,i) = tend(k,j,i) - vkomp * ( v(k,j,i) - v(k,j-1,i) ) * ddy
+             ELSE
+                tend(k,j,i) = tend(k,j,i) - vkomp * ( v(k,j+1,i) - v(k,j,i) ) * ddy
+             ENDIF
 !
-!--             z-direction
-                wkomp = 0.25_wp * ( w(k,j,i)  + w(k-1,j,i) +                   &
-                                 w(k,j-1,i) + w(k-1,j-1,i) )
-                IF ( wkomp > 0.0_wp )  THEN
-                   tend(k,j,i) = tend(k,j,i) - wkomp *                         &
-                                         ( v(k,j,i) - v(k-1,j,i) ) * ddzu(k)
-                ELSE
-                   tend(k,j,i) = tend(k,j,i) - wkomp *                         &
-                                         ( v(k+1,j,i) - v(k,j,i) ) * ddzu(k+1)
-                ENDIF
+!--          z-direction
+             wkomp = 0.25_wp * ( w(k,j,i) + w(k-1,j,i) + w(k,j-1,i) + w(k-1,j-1,i) )
+             IF ( wkomp > 0.0_wp )  THEN
+                tend(k,j,i) = tend(k,j,i) - wkomp * ( v(k,j,i) - v(k-1,j,i) ) * ddzu(k)
+             ELSE
+                tend(k,j,i) = tend(k,j,i) - wkomp * ( v(k+1,j,i) - v(k,j,i) ) * ddzu(k+1)
+             ENDIF
 
-             ENDDO
           ENDDO
        ENDDO
+    ENDDO
 
-    END SUBROUTINE advec_v_up
+ END SUBROUTINE advec_v_up
 
 
-!------------------------------------------------------------------------------!
+!--------------------------------------------------------------------------------------------------!
 ! Description:
 ! ------------
 !> Call for grid point i,j
-!------------------------------------------------------------------------------!
-    SUBROUTINE advec_v_up_ij( i, j )
+!--------------------------------------------------------------------------------------------------!
+ SUBROUTINE advec_v_up_ij( i, j )
 
-       USE arrays_3d,                                                          &
-           ONLY:  ddzu, tend, u, v, w
+    USE arrays_3d,                                                                                 &
+        ONLY:  ddzu, tend, u, v, w
 
-       USE control_parameters,                                                 &
-           ONLY:  u_gtrans, v_gtrans
+    USE control_parameters,                                                                        &
+        ONLY:  u_gtrans, v_gtrans
 
-       USE grid_variables,                                                     &
-           ONLY:  ddx, ddy
+    USE grid_variables,                                                                            &
+        ONLY:  ddx, ddy
 
-       USE indices,                                                            &
-           ONLY:  nzb, nzt
+    USE indices,                                                                                   &
+        ONLY:  nzb, nzt
 
-       USE kinds
+    USE kinds
 
 
-       IMPLICIT NONE
+    IMPLICIT NONE
 
-       INTEGER(iwp) ::  i !< grid index along x-direction
-       INTEGER(iwp) ::  j !< grid index along y-direction
-       INTEGER(iwp) ::  k !< grid index along z-direction
+    INTEGER(iwp) ::  i !< grid index along x-direction
+    INTEGER(iwp) ::  j !< grid index along y-direction
+    INTEGER(iwp) ::  k !< grid index along z-direction
 
-       REAL(wp) ::  ukomp !< advection velocity along x-direction 
-       REAL(wp) ::  vkomp !< advection velocity along y-direction 
-       REAL(wp) ::  wkomp !< advection velocity along z-direction 
+    REAL(wp) ::  ukomp !< advection velocity along x-direction 
+    REAL(wp) ::  vkomp !< advection velocity along y-direction 
+    REAL(wp) ::  wkomp !< advection velocity along z-direction 
 
 
-       DO  k = nzb+1, nzt
+    DO  k = nzb+1, nzt
 !
-!--       x-direction
-          ukomp = 0.25_wp * ( u(k,j,i) + u(k,j-1,i) + u(k,j,i+1) + u(k,j-1,i+1) &
-                         ) - u_gtrans
-          IF ( ukomp > 0.0_wp )  THEN
-             tend(k,j,i) = tend(k,j,i) - ukomp *                               &
-                                         ( v(k,j,i) - v(k,j,i-1) ) * ddx
-          ELSE
-             tend(k,j,i) = tend(k,j,i) - ukomp *                               &
-                                         ( v(k,j,i+1) - v(k,j,i) ) * ddx
-          ENDIF
+!--    x-direction
+       ukomp = 0.25_wp * ( u(k,j,i) + u(k,j-1,i) + u(k,j,i+1) + u(k,j-1,i+1) ) - u_gtrans
+       IF ( ukomp > 0.0_wp )  THEN
+          tend(k,j,i) = tend(k,j,i) - ukomp * ( v(k,j,i) - v(k,j,i-1) ) * ddx
+       ELSE
+          tend(k,j,i) = tend(k,j,i) - ukomp * ( v(k,j,i+1) - v(k,j,i) ) * ddx
+       ENDIF
 !
-!--       y-direction
-          vkomp = v(k,j,i) - v_gtrans
-          IF ( vkomp > 0.0_wp )  THEN
-             tend(k,j,i) = tend(k,j,i) - vkomp *                               &
-                                         ( v(k,j,i) - v(k,j-1,i) ) * ddy
-          ELSE
-             tend(k,j,i) = tend(k,j,i) - vkomp *                               &
-                                         ( v(k,j+1,i) - v(k,j,i) ) * ddy
-          ENDIF
+!--    y-direction
+       vkomp = v(k,j,i) - v_gtrans
+       IF ( vkomp > 0.0_wp )  THEN
+          tend(k,j,i) = tend(k,j,i) - vkomp * ( v(k,j,i) - v(k,j-1,i) ) * ddy
+       ELSE
+          tend(k,j,i) = tend(k,j,i) - vkomp * ( v(k,j+1,i) - v(k,j,i) ) * ddy
+       ENDIF
 !
-!--       z-direction
-          wkomp = 0.25_wp * ( w(k,j,i) + w(k-1,j,i) + w(k,j-1,i) + w(k-1,j-1,i) )
-          IF ( wkomp > 0.0_wp )  THEN
-             tend(k,j,i) = tend(k,j,i) - wkomp *                               &
-                                         ( v(k,j,i) - v(k-1,j,i) ) * ddzu(k)
-          ELSE
-             tend(k,j,i) = tend(k,j,i) - wkomp *                               &
-                                         ( v(k+1,j,i) - v(k,j,i) ) * ddzu(k+1)
-          ENDIF
+!--    z-direction
+       wkomp = 0.25_wp * ( w(k,j,i) + w(k-1,j,i) + w(k,j-1,i) + w(k-1,j-1,i) )
+       IF ( wkomp > 0.0_wp )  THEN
+          tend(k,j,i) = tend(k,j,i) - wkomp * ( v(k,j,i) - v(k-1,j,i) ) * ddzu(k)
+       ELSE
+          tend(k,j,i) = tend(k,j,i) - wkomp * ( v(k+1,j,i) - v(k,j,i) ) * ddzu(k+1)
+       ENDIF
 
-       ENDDO
+    ENDDO
 
-    END SUBROUTINE advec_v_up_ij
+ END SUBROUTINE advec_v_up_ij
 
  END MODULE advec_v_up_mod
Index: /palm/trunk/SOURCE/advec_w_pw.f90
===================================================================
--- /palm/trunk/SOURCE/advec_w_pw.f90	(revision 4487)
+++ /palm/trunk/SOURCE/advec_w_pw.f90	(revision 4488)
@@ -1,20 +1,19 @@
 !> @file advec_w_pw.f90
-!------------------------------------------------------------------------------!
+!--------------------------------------------------------------------------------------------------!
 ! This file is part of the PALM model system.
 !
-! PALM is free software: you can redistribute it and/or modify it under the
-! terms of the GNU General Public License as published by the Free Software
-! Foundation, either version 3 of the License, or (at your option) any later
-! version.
+! PALM is free software: you can redistribute it and/or modify it under the terms of the GNU General
+! Public License as published by the Free Software Foundation, either version 3 of the License, or
+! (at your option) any later version.
 !
-! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
-! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
-! A PARTICULAR PURPOSE.  See the GNU General Public License for more details.
+! PALM is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the
+! implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General
+! Public License for more details.
 !
-! You should have received a copy of the GNU General Public License along with
-! PALM. If not, see <http://www.gnu.org/licenses/>.
+! You should have received a copy of the GNU General Public License along with PALM. If not, see
+! <http://www.gnu.org/licenses/>.
 !
 ! Copyright 1997-2020 Leibniz Universitaet Hannover
-!------------------------------------------------------------------------------!
+!--------------------------------------------------------------------------------------------------!
 !
 ! Current revisions:
@@ -25,4 +24,7 @@
 ! -----------------
 ! $Id$
+! file re-formatted to follow the PALM coding standard
+!
+! 4360 2020-01-07 11:25:50Z suehring
 ! Corrected "Former revisions" section
 ! 
@@ -37,8 +39,8 @@
 ! ------------
 !> Advection term for w velocity-component using Piacsek and Williams.
-!> Vertical advection at the first grid point above the surface is done with
-!> normal centred differences, because otherwise no information from the surface
-!> would be communicated upwards due to w=0 at k=nzb.
-!------------------------------------------------------------------------------!
+!> Vertical advection at the first grid point above the surface is done with normal centred
+!> differences, because otherwise no information from the surface would be communicated upwards due
+!> to w=0 at k=nzb.
+!--------------------------------------------------------------------------------------------------!
  MODULE advec_w_pw_mod
  
@@ -55,102 +57,100 @@
 
 
-!------------------------------------------------------------------------------!
+!--------------------------------------------------------------------------------------------------!
 ! Description:
 ! ------------
 !> Call for all grid points
-!------------------------------------------------------------------------------!
-    SUBROUTINE advec_w_pw
+!--------------------------------------------------------------------------------------------------!
+ SUBROUTINE advec_w_pw
 
-       USE arrays_3d,                                                          &
-           ONLY:  ddzu, tend, u, v, w
+    USE arrays_3d,                                                                                 &
+        ONLY:  ddzu, tend, u, v, w
 
-       USE control_parameters,                                                 &
-           ONLY:  u_gtrans, v_gtrans
+    USE control_parameters,                                                                        &
+        ONLY:  u_gtrans, v_gtrans
 
-       USE grid_variables,                                                     &
-           ONLY:  ddx, ddy
+    USE grid_variables,                                                                            &
+        ONLY:  ddx, ddy
 
-       USE indices,                                                            &
-           ONLY:  nxl, nxr, nyn, nys, nzb, nzt
+    USE indices,                                                                                   &
+        ONLY:  nxl, nxr, nyn, nys, nzb, nzt
 
-       USE kinds
+    USE kinds
 
 
-       IMPLICIT NONE
+    IMPLICIT NONE
 
-       INTEGER(iwp) ::  i !< grid index along x-direction
-       INTEGER(iwp) ::  j !< grid index along y-direction
-       INTEGER(iwp) ::  k !< grid index along z-direction
-       
-       REAL(wp)    ::  gu !< Galilei-transformation velocity along x
-       REAL(wp)    ::  gv !< Galilei-transformation velocity along y
+    INTEGER(iwp) ::  i !< grid index along x-direction
+    INTEGER(iwp) ::  j !< grid index along y-direction
+    INTEGER(iwp) ::  k !< grid index along z-direction
+    
+    REAL(wp)     ::  gu !< Galilei-transformation velocity along x
+    REAL(wp)     ::  gv !< Galilei-transformation velocity along y
 
- 
-       gu = 2.0_wp * u_gtrans
-       gv = 2.0_wp * v_gtrans
-       DO  i = nxl, nxr
-          DO  j = nys, nyn
-             DO  k = nzb+1, nzt
-                tend(k,j,i) = tend(k,j,i) - 0.25_wp * (                        &
-                         ( w(k,j,i+1) * ( u(k+1,j,i+1) + u(k,j,i+1) - gu )     &
-                         - w(k,j,i-1) * ( u(k+1,j,i) + u(k,j,i) - gu ) ) * ddx &
-                       + ( w(k,j+1,i) * ( v(k+1,j+1,i) + v(k,j+1,i) - gv )     &
-                         - w(k,j-1,i) * ( v(k+1,j,i) + v(k,j,i) - gv ) ) * ddy &
-                       + ( w(k+1,j,i) * ( w(k+1,j,i) + w(k,j,i) )              &
-                         - w(k-1,j,i) * ( w(k,j,i) + w(k-1,j,i) ) )            &
-                                                                  * ddzu(k+1)  &
-                                                      )
-             ENDDO
+
+    gu = 2.0_wp * u_gtrans
+    gv = 2.0_wp * v_gtrans
+    DO  i = nxl, nxr
+       DO  j = nys, nyn
+          DO  k = nzb+1, nzt
+             tend(k,j,i) = tend(k,j,i) - 0.25_wp * (                                               &
+                           ( w(k,j,i+1)   * ( u(k+1,j,i+1) + u(k,j,i+1) - gu )                     &
+                           - w(k,j,i-1)   * ( u(k+1,j,i)   + u(k,j,i) - gu ) ) * ddx               &
+                           + ( w(k,j+1,i) * ( v(k+1,j+1,i) + v(k,j+1,i) - gv )                     &
+                           - w(k,j-1,i)   * ( v(k+1,j,i)   + v(k,j,i) - gv ) ) * ddy               &
+                           + ( w(k+1,j,i) * ( w(k+1,j,i)   + w(k,j,i) )                            &
+                           - w(k-1,j,i)   * ( w(k,j,i)     + w(k-1,j,i) ) )    * ddzu(k+1)         &
+                                                   )
           ENDDO
        ENDDO
+    ENDDO
 
-    END SUBROUTINE advec_w_pw
+ END SUBROUTINE advec_w_pw
 
 
-!------------------------------------------------------------------------------!
+!--------------------------------------------------------------------------------------------------!
 ! Description:
 ! ------------
 !> Call for grid point i,j
-!------------------------------------------------------------------------------!
-    SUBROUTINE advec_w_pw_ij( i, j )
+!--------------------------------------------------------------------------------------------------!
+ SUBROUTINE advec_w_pw_ij( i, j )
 
-       USE arrays_3d,                                                          &
-           ONLY:  ddzu, tend, u, v, w
+    USE arrays_3d,                                                                                 &
+        ONLY:  ddzu, tend, u, v, w
 
-       USE control_parameters,                                                 &
-           ONLY:  u_gtrans, v_gtrans
+    USE control_parameters,                                                                        &
+        ONLY:  u_gtrans, v_gtrans
 
-       USE grid_variables,                                                     &
-           ONLY:  ddx, ddy
+    USE grid_variables,                                                                            &
+        ONLY:  ddx, ddy
 
-       USE indices,                                                            &
-           ONLY:  nzb, nzt
+    USE indices,                                                                                   &
+        ONLY:  nzb, nzt
 
-       USE kinds
+    USE kinds
 
 
-       IMPLICIT NONE
+    IMPLICIT NONE
 
-       INTEGER(iwp) ::  i !< grid index along x-direction
-       INTEGER(iwp) ::  j !< grid index along y-direction
-       INTEGER(iwp) ::  k !< grid index along z-direction
-       
-       REAL(wp)    ::  gu !< Galilei-transformation velocity along x
-       REAL(wp)    ::  gv !< Galilei-transformation velocity along y
+    INTEGER(iwp) ::  i !< grid index along x-direction
+    INTEGER(iwp) ::  j !< grid index along y-direction
+    INTEGER(iwp) ::  k !< grid index along z-direction
+    
+    REAL(wp)     ::  gu !< Galilei-transformation velocity along x
+    REAL(wp)     ::  gv !< Galilei-transformation velocity along y
 
-       gu = 2.0_wp * u_gtrans
-       gv = 2.0_wp * v_gtrans
-       DO  k = nzb+1, nzt
-          tend(k,j,i) = tend(k,j,i) - 0.25_wp * (                              &
-                         ( w(k,j,i+1) * ( u(k+1,j,i+1) + u(k,j,i+1) - gu )     &
-                         - w(k,j,i-1) * ( u(k+1,j,i) + u(k,j,i) - gu ) ) * ddx &
-                       + ( w(k,j+1,i) * ( v(k+1,j+1,i) + v(k,j+1,i) - gv )     &
-                         - w(k,j-1,i) * ( v(k+1,j,i) + v(k,j,i) - gv ) ) * ddy &
-                       + ( w(k+1,j,i) * ( w(k+1,j,i) + w(k,j,i) )              &
-                         - w(k-1,j,i) * ( w(k,j,i) + w(k-1,j,i) ) )            &
-                                                                  * ddzu(k+1)  &
-                                                )
-       ENDDO
-    END SUBROUTINE advec_w_pw_ij
+    gu = 2.0_wp * u_gtrans
+    gv = 2.0_wp * v_gtrans
+    DO  k = nzb+1, nzt
+       tend(k,j,i) = tend(k,j,i) - 0.25_wp * (                                                     &
+                      ( w(k,j,i+1)   * ( u(k+1,j,i+1) + u(k,j,i+1) - gu )                          &
+                      - w(k,j,i-1)   * ( u(k+1,j,i)   + u(k,j,i) - gu ) ) * ddx                    &
+                      + ( w(k,j+1,i) * ( v(k+1,j+1,i) + v(k,j+1,i) - gv )                          &
+                      - w(k,j-1,i)   * ( v(k+1,j,i)   + v(k,j,i) - gv ) ) * ddy                    &
+                      + ( w(k+1,j,i) * ( w(k+1,j,i)   + w(k,j,i) )                                 &
+                      - w(k-1,j,i)   * ( w(k,j,i)     + w(k-1,j,i) ) )    * ddzu(k+1)              &
+                                             )
+    ENDDO
+ END SUBROUTINE advec_w_pw_ij
 
  END MODULE advec_w_pw_mod
