scale_atmos_sub_diagnostic.F90

Go to the documentation of this file.

!-------------------------------------------------------------------------------
!-------------------------------------------------------------------------------
#include "inc_openmp.h"
module scale_atmos_diagnostic
  !-----------------------------------------------------------------------------
  !
  !++ used modules
  !
  use scale_precision
  use scale_stdio
  use scale_prof
  use scale_grid_index
  use scale_tracer
  !-----------------------------------------------------------------------------
  implicit none
  private
  !-----------------------------------------------------------------------------
  !
  !++ Public procedure
  !
  public :: atmos_diagnostic_get

  !-----------------------------------------------------------------------------
  !
  !++ Public parameters & variables
  !
  !-----------------------------------------------------------------------------
  !
  !++ Private procedure
  !
  !-----------------------------------------------------------------------------
  !
  !++ Private parameters & variables
  !
  !-----------------------------------------------------------------------------
contains
  !-----------------------------------------------------------------------------
  subroutine atmos_diagnostic_get( &
       POTT, &
       TEMP, &
       PRES, &
       PHYD, &
       W,    &
       U,    &
       V,    &
       N2,   &
       DENS, &
       MOMZ, &
       MOMX, &
       MOMY, &
       RHOT, &
       QTRC )
    use scale_comm, only: &
       comm_vars8, &
       comm_wait
    use scale_const, only: &
       grav => const_grav
    use scale_atmos_thermodyn, only: &
       thermodyn_qd        => atmos_thermodyn_qd, &
       thermodyn_r         => atmos_thermodyn_r,  &
       thermodyn_temp_pres => atmos_thermodyn_temp_pres
    use scale_grid_real, only: &
       cz => real_cz, &
       fz => real_fz
    use scale_gridtrans, only: &
       gsqrt => gtrans_gsqrt, &
       j13g  => gtrans_j13g,  &
       j23g  => gtrans_j23g,  &
       i_xyw, &
       i_xyz
    implicit none

    real(RP), intent(out) :: pott(ka,ia,ja)
    real(RP), intent(out) :: temp(ka,ia,ja)
    real(RP), intent(out) :: pres(ka,ia,ja)
    real(RP), intent(out) :: phyd(ka,ia,ja)
    real(RP), intent(out) :: w   (ka,ia,ja)
    real(RP), intent(out) :: u   (ka,ia,ja)
    real(RP), intent(out) :: v   (ka,ia,ja)
    real(RP), intent(out) :: n2  (ka,ia,ja)
    real(RP), intent(in)  :: dens(ka,ia,ja)
    real(RP), intent(in)  :: momz(ka,ia,ja)
    real(RP), intent(in)  :: momx(ka,ia,ja)
    real(RP), intent(in)  :: momy(ka,ia,ja)
    real(RP), intent(in)  :: rhot(ka,ia,ja)
    real(RP), intent(in)  :: qtrc(ka,ia,ja,qa)

    real(RP) :: ph(ka)
    real(RP) :: rpt(ka)
    real(RP) :: q(qa)
    real(RP) :: qdry, rtot
    real(RP) :: momws

    integer  :: k, i, j
    integer  :: iq
    !---------------------------------------------------------------------------

    call thermodyn_temp_pres( temp(:,:,:),   & ! [OUT]
                              pres(:,:,:),   & ! [OUT]
                              dens(:,:,:),   & ! [IN]
                              rhot(:,:,:),   & ! [IN]
                              qtrc(:,:,:,:), & ! [IN]
                              tracer_cv(:),       & ! [IN]
                              tracer_r(:),       & ! [IN]
                              tracer_mass(:)        ) ! [IN]

    !$omp parallel do private(i,j,k,ph) OMP_SCHEDULE_ collapse(2)
    do j = 1, ja
    do i = 1, ia
       ph(ke) = pres(ke,i,j) - dens(ke,i,j) * grav * ( fz(ke,i,j) - cz(ke,i,j) )
       do k = ke, ks, -1
          ph(k-1) = ph(k) + dens(k,i,j) * grav * ( fz(k,i,j) - fz(k-1,i,j) )
          phyd(k,i,j) = ( ph(k) + ph(k-1) ) * 0.5_rp
       end do
    enddo
    enddo


!OCL XFILL
    !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
    do j = 1, ja
    do i = 1, ia
    do k = ks+1, ke-1
       w(k,i,j) = 0.5_rp * ( momz(k-1,i,j)+momz(k,i,j) ) / dens(k,i,j)
    enddo
    enddo
    enddo
!OCL XFILL
    !$omp parallel do private(i,j) OMP_SCHEDULE_ collapse(2)
    do j = 2, ja
    do i = 2, ia
!       W(KS,i,j) = 0.5_RP * ( MOMZ(KS,i,j) ) / DENS(KS,i,j)

       ! at KS+1/2
       momws = momz(ks,i,j) &
             + ( j13g(ks,i,j,i_xyw) * ( momx(ks,i,j) + momx(ks,i-1,j) + momx(ks+1,i,j) + momx(ks+1,i-1,j) ) &
               + j23g(ks,i,j,i_xyw) * ( momy(ks,i,j) + momy(ks,i,j-1) + momy(ks+1,i,j) + momy(ks+1,i,j-1) ) ) &
               * 0.25_rp / gsqrt(ks,i,j,i_xyw)
       ! at KS
       ! momws at the surface is assumed to be zero
       w(ks,i,j) = momws * 0.5_rp &
                 - ( j13g(ks,i,j,i_xyz) * ( momx(ks,i,j) + momx(ks,i-1,j) ) &
                   + j23g(ks,i,j,i_xyz) * ( momy(ks,i,j) + momy(ks,i,j-1) ) ) &
                   * 0.5_rp / ( dens(ks,i,j) * gsqrt(ks,i,j,i_xyz) )
    enddo
    enddo
!OCL XFILL
    !$omp parallel do private(i,j) OMP_SCHEDULE_ collapse(2)
    do j = 1, ja
    do i = 1, ia
       w(ke,i,j) = 0.5_rp * ( momz(ke-1,i,j) ) / dens(ke,i,j)
    enddo
    enddo

!OCL XFILL
    !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
    do j = 1, ja
    do i = 2, ia
    do k = ks, ke
       u(k,i,j) = 0.5_rp * ( momx(k,i-1,j)+momx(k,i,j) ) / dens(k,i,j)
    enddo
    enddo
    enddo
!OCL XFILL
    !$omp parallel do private(j,k) OMP_SCHEDULE_ collapse(2)
    do j = 1, ja
    do k = ks, ke
       u(k,1,j) = momx(k,1,j) / dens(k,1,j)
    enddo
    enddo

!OCL XFILL
    !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
    do j = 2, ja
    do i = 1, ia
    do k = ks, ke
       v(k,i,j) = 0.5_rp * ( momy(k,i,j-1)+momy(k,i,j) ) / dens(k,i,j)
    enddo
    enddo
    enddo
!OCL XFILL
    !$omp parallel do private(i,k) OMP_SCHEDULE_ collapse(2)
    do i = 1, ia
    do k = ks, ke
       v(k,i,1) = momy(k,i,1) / dens(k,i,1)
    enddo
    enddo

    !$omp parallel do private(i,j) OMP_SCHEDULE_ collapse(2)
    do j  = 1, ja
    do i  = 1, ia
       w(   1:ks-1,i,j) = w(ks,i,j)
       u(   1:ks-1,i,j) = u(ks,i,j)
       v(   1:ks-1,i,j) = v(ks,i,j)
       w(ke+1:ka,  i,j) = w(ke,i,j)
       u(ke+1:ka,  i,j) = u(ke,i,j)
       v(ke+1:ka,  i,j) = v(ke,i,j)
    enddo
    enddo

    call comm_vars8( w(:,:,:), 1 )
    call comm_vars8( u(:,:,:), 2 )
    call comm_vars8( v(:,:,:), 3 )
    call comm_wait ( w(:,:,:), 1, .false. )
    call comm_wait ( u(:,:,:), 2, .false. )
    call comm_wait ( v(:,:,:), 3, .false. )

!OCL XFILL
    !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
    do j = 1, ja
    do i = 1, ia
    do k = ks, ke
       pott(k,i,j) = rhot(k,i,j) / dens(k,i,j)
    enddo
    enddo
    enddo

    !$omp parallel do OMP_SCHEDULE_ collapse(2) &
    !$omp private(i,j,q,rpt,qdry,rtot)
    do j = 1, ja
    do i = 1, ia
       do k = ks, ke
          do iq = 1, qa
             q(iq) = qtrc(k,i,j,iq)
          end do
          call thermodyn_qd( qdry, q(:), tracer_mass(:) )
          call thermodyn_r ( rtot, q(:), tracer_r(:), qdry )
          rpt(k) = rtot * pott(k,i,j)
       end do

       n2(ks,i,j) = grav * ( rpt(ks+1) - rpt(ks) ) / ( ( cz(ks+1,i,j) - cz(ks,i,j) ) * rpt(ks) )
       do k = ks+1,ke-1
          n2(k,i,j) = grav * ( rpt(k+1) - rpt(k-1) ) / ( ( cz(k+1,i,j) - cz(k-1,i,j) ) * rpt(k) )
       end do
       n2(ke,i,j) = grav * ( rpt(ke) - rpt(ke-1) ) / ( ( cz(ke,i,j) - cz(ke-1,i,j) ) * rpt(ke) )
    end do
    end do

    return
  end subroutine atmos_diagnostic_get


end module scale_atmos_diagnostic