source: palm/trunk/SOURCE/advec_s_ups.f90 @ 844

 SUBROUTINE advec_s_ups( s, var_char )

!------------------------------------------------------------------------------!
! Current revisions:
! -----------------
! 
!
! Former revisions:
! -----------------
! $Id: advec_s_ups.f90 484 2010-02-05 07:36:54Z gryschka $
!
! 164 2008-05-15 08:46:15Z raasch
! Arguments removed from transpose routines
!
! February 2007
! RCS Log replace by Id keyword, revision history cleaned up
!
! Revision 1.6  2004/04/30 08:02:43  raasch
! Enlarged transposition arrays introduced
!
! Revision 1.1  1999/02/05 08:44:47  raasch
! Initial revision
!
!
! Description:
! ------------
! Upstream-Spline advection of scalar quantities (potential temperature,
! turbulent kinetic energy). The advection process is divided into three 
! subsequent steps, one for each of the dimensions. The result is stored as a 
! tendency in array tend. The computation of the cubic splines and the possible
! execution of the Long-filter require that all grid points of the relevant 
! dimension are available. For model runs on more than one PE therefore both the
! advected and the advecting quantities are transposed accordingly.
!
! Actual arguments:
! s        = scalar quantity to be advected (remains unchanged in this UP)
! var_char = character string specifying the quantity to be advected
!
! Internally used arrays:
! v_ad     = scalar quantity to be advected, initialized = s at the beginning,
!            also being used as temporary storage after each time step
! d        = advecting component (u, v, or w)
!------------------------------------------------------------------------------!

    USE advection
    USE arrays_3d
    USE cpulog
    USE grid_variables
    USE indices
    USE interfaces
    USE control_parameters

    IMPLICIT NONE

    CHARACTER (LEN=*) ::  var_char

    INTEGER ::  i, j, k
    REAL    ::  s(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1)
    REAL, DIMENSION(:,:,:), ALLOCATABLE ::  v_ad

    CALL cpu_log( log_point_s(16), 'advec_s_ups', 'start' )

#if defined( __parallel )

!
!-- Advection of the scalar in x-direction:
!-- Store the scalar in temporary array v_ad (component to be advected,
!-- boundaries are not used because they disturb the transposition)
    ALLOCATE( v_ad(nzb+1:nzta,nys:nyna,nxl:nxra) )
    v_ad = 0.0
    v_ad(nzb+1:nzt,nys:nyn,nxl:nxr) = s(nzb+1:nzt,nys:nyn,nxl:nxr)

!
!-- Enlarge the size of tend, used as a working array for the transpositions
    IF ( nxra > nxr  .OR.  nyna > nyn  .OR.  nza > nz )  THEN
       DEALLOCATE( tend )
       ALLOCATE( tend(1:nza,nys:nyna,nxl:nxra) )
    ENDIF

!
!-- Transpose the component to be advected: z --> x
    CALL transpose_zx( v_ad, tend, v_ad )

#else

!
!-- Advection of the scalar in x-direction:
!-- Store the scalar in temporary array v_ad (component to be advected)
    ALLOCATE( v_ad(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1) )
    v_ad(:,:,:) = s(:,:,:)

#endif

!
!-- Advecting component (u) must be averaged out on the scalar's grid
    DO  i = nxl, nxr
       DO  j = nys, nyn
          DO  k = nzb+1, nzt
              d(k,j,i) = 0.5 * ( u(k,j,i) + u(k,j,i+1) ) - u_gtrans
          ENDDO
       ENDDO
    ENDDO

#if defined( __parallel )

!
!-- Transpose the advecting componnet: z --> x
    CALL transpose_zx( d, tend, d )

#endif

!
!-- Upstream-Spline advection of the scalar in x-direction
    CALL spline_x( v_ad, d, var_char )

!
!-- Advection of the scalar in y-direction:
!-- advecting component (v) must be averaged out on the scalar's grid
    DO  i = nxl, nxr
       DO  j = nys, nyn
          DO  k = nzb+1, nzt
              d(k,j,i) = 0.5 * ( v(k,j,i) + v(k,j+1,i) ) - v_gtrans
          ENDDO
       ENDDO
    ENDDO

#if defined( __parallel )
    
!
!-- Transpose the advecting component: z --> y
    CALL transpose_zx( d, tend, d )
    CALL transpose_xy( d, tend, d )

!
!-- Transpose the component to be advected: x --> y
    CALL transpose_xy( v_ad, tend, v_ad )

#endif

!
!-- Upstream-Spline advection of the scalar in y-direction
    CALL spline_y( v_ad, d, var_char )

!
!-- Advection of the scalar in z-direction:
!-- the advecting component (w) must be averaged out on the scalar's grid
!-- (weighted for non-equidistant grid)
    d = 0.0
    DO  i = nxl, nxr
       DO  j = nys, nyn
          DO  k = nzb+1, nzt
             d(k,j,i) =  ( w(k,j,i) * ( zu(k) - zw(k-1) ) + &
                           w(k-1,j,i) * ( zw(k) - zu(k) ) ) * ddzw(k)
          ENDDO
       ENDDO
    ENDDO

#if defined( __parallel )

!
!-- Transpose the component to be advected: y --> z (= y --> x + x --> z)
    CALL transpose_yx( v_ad, tend, v_ad )
    CALL transpose_xz( v_ad, tend, v_ad )

!
!-- Resize tend to its normal size
    IF ( nxra > nxr  .OR.  nyna > nyn  .OR.  nza > nz )  THEN
       DEALLOCATE( tend )
       ALLOCATE( tend(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1) )
    ENDIF

#endif

!
!-- Upstream-Spline advection of the scalar in z-direction
    CALL spline_z( v_ad, d, dzu, spl_tri_zu, var_char )

!
!-- Compute the tendency term
    DO  i = nxl, nxr
       DO  j = nys, nyn
          DO  k = nzb+1, nzt
             tend(k,j,i) = ( v_ad(k,j,i) - s(k,j,i) ) / dt_3d 
          ENDDO
       ENDDO
    ENDDO

    DEALLOCATE( v_ad )

    CALL cpu_log( log_point_s(16), 'advec_s_ups', 'stop' )

 END SUBROUTINE advec_s_ups