scale_spnudge.F90

Go to the documentation of this file.

#include "scalelib.h"
module scale_spnudge
  !-----------------------------------------------------------------------------
  !
  !++ used modules
  !
  use scale_precision
  use scale_io
  use scale_prc
 
  !-----------------------------------------------------------------------------
  implicit none
  private
  !-----------------------------------------------------------------------------
  !
  !++ Public procedure
  !
  public :: spnudge_setup
 
  !-----------------------------------------------------------------------------
  !
  !++ Public variables
  !
  logical, public :: spnudge_uv         = .false.
  logical, public :: spnudge_uv_divfree = .false.
  integer, public :: spnudge_uv_lm      = 3
  integer, public :: spnudge_uv_mm      = 3
 
  logical, public :: spnudge_pt      = .false.
  integer, public :: spnudge_pt_lm   = 3
  integer, public :: spnudge_pt_mm   = 3
 
  logical, public :: spnudge_qv      = .false.
  integer, public :: spnudge_qv_lm   = 3
  integer, public :: spnudge_qv_mm   = 3
 
  real(rp), allocatable, public :: spnudge_u_alpha(:,:,:)
  real(rp), allocatable, public :: spnudge_v_alpha(:,:,:)
  real(rp), allocatable, public :: spnudge_pt_alpha(:,:,:)
  real(rp), allocatable, public :: spnudge_qv_alpha(:,:,:)
 
  !-----------------------------------------------------------------------------
  !
  !++ Private variables
  !
 
contains
 
  subroutine spnudge_setup(KA, KS, KE, IA, IS, IE, JA, JS, JE)
    use scale_dft, only: &
       dft_setup
    use scale_atmos_grid_cartesc_real, only: &
       real_cz => atmos_grid_cartesc_real_cz, &
       real_czuy => atmos_grid_cartesc_real_czuy, &
       real_czxv => atmos_grid_cartesc_real_czxv
    implicit none
    integer, intent(in) :: ka, ks, ke
    integer, intent(in) :: ia, is, ie
    integer, intent(in) :: ja, js, je
 
    real(rp) :: spnudge_uv_tau    = 0.0_rp
    real(rp) :: spnudge_pt_tau = 0.0_rp
    real(rp) :: spnudge_qv_tau = 0.0_rp
 
    real(rp) :: spnudge_level1 = 0.0_rp
    real(rp) :: spnudge_level2 = 0.0_rp
    real(rp) :: spnudge_level3 = 1e10_rp
    real(rp) :: spnudge_level4 = 1e10_rp
 
    namelist /param_spnudge/ &
      spnudge_uv, &
      spnudge_uv_divfree, &
      spnudge_uv_lm, &
      spnudge_uv_mm, &
      spnudge_uv_tau, &
      spnudge_pt, &
      spnudge_pt_lm, &
      spnudge_pt_mm, &
      spnudge_pt_tau, &
      spnudge_qv, &
      spnudge_qv_lm, &
      spnudge_qv_mm, &
      spnudge_qv_tau, &
      spnudge_level1, &
      spnudge_level2, &
      spnudge_level3, &
      spnudge_level4
 
    real(rp) :: uv_alpha, pt_alpha, qv_alpha
 
 
    integer :: k, i, j
    integer :: ierr
    !---------------------------------------------------------------------------
 
    log_newline
    log_info("SPNUDGE_setup",*) 'Setup'
 
    !--- read namelist
    rewind(io_fid_conf)
    read(io_fid_conf,nml=param_spnudge,iostat=ierr)
    if( ierr < 0 ) then !--- missing
       log_info("SPNUDGE_setup",*) 'Not found namelist. Default used.'
    elseif( ierr > 0 ) then !--- fatal error
       log_error("SPNUDGE_setup",*) 'Not appropriate names in namelist PARAM_SPNUDGE. Check!'
       call prc_abort
    endif
    log_nml(param_spnudge)
 
    if ( spnudge_level1 > spnudge_level2 ) then
       log_error("SPNUDGE_setup",*) 'SPNUDGE_level1 must be lowere or equal to SPNUDGE_level2'
       call prc_abort
    end if
    if ( spnudge_level2 > spnudge_level3 ) then
       log_error("SPNUDGE_setup",*) 'SPNUDGE_level2 must be lowere or equal to SPNUDGE_level3'
       call prc_abort
    end if
    if ( spnudge_level3 > spnudge_level4 ) then
       log_error("SPNUDGE_setup",*) 'SPNUDGE_level3 must be lowere or equal to SPNUDGE_level4'
       call prc_abort
    end if
 
    if ( spnudge_uv .or. spnudge_pt .or. spnudge_qv ) then
       log_warn("SPNUDGE_setup",*) 'Spectrul nudging is still experimental'
    end if
 
    allocate( spnudge_u_alpha(ka,ia,ja) )
    allocate( spnudge_v_alpha(ka,ia,ja) )
    allocate( spnudge_pt_alpha(ka,ia,ja) )
    allocate( spnudge_qv_alpha(ka,ia,ja) )
 
    if ( spnudge_uv .and. spnudge_uv_tau > 0.0_rp ) then
       uv_alpha = 1.0_rp / spnudge_uv_tau
 
       !$omp parallel do
       do j = js, je
       do i = is, ie
       do k = ks, ke
          if ( real_czuy(k,i,j) < spnudge_level1 ) then
             spnudge_u_alpha(k,i,j) = 0.0_rp
          elseif ( real_czuy(k,i,j) < spnudge_level2 ) then
             spnudge_u_alpha(k,i,j) = ( real_czuy(k,i,j) - spnudge_level1 ) / ( spnudge_level2 - spnudge_level1 ) * uv_alpha
          elseif ( real_czuy(k,i,j) < spnudge_level3 ) then
             spnudge_u_alpha(k,i,j) = uv_alpha
          elseif ( real_czuy(k,i,j) < spnudge_level4 ) then
             spnudge_u_alpha(k,i,j) = ( real_czuy(k,i,j) - spnudge_level3 ) / ( spnudge_level4 - spnudge_level3 ) * uv_alpha
          else
             spnudge_u_alpha(k,i,j) = 0.0_rp
          endif
       enddo
       enddo
       enddo
 
       !$omp parallel do
       do j = js, je
       do i = is, ie
       do k = ks, ke
          if ( real_czxv(k,i,j) < spnudge_level1 ) then
             spnudge_v_alpha(k,i,j) = 0.0_rp
          elseif ( real_czxv(k,i,j) < spnudge_level2 ) then
             spnudge_v_alpha(k,i,j) = ( real_czxv(k,i,j) - spnudge_level1 ) / ( spnudge_level2 - spnudge_level1 ) * uv_alpha
          elseif ( real_czxv(k,i,j) < spnudge_level3 ) then
             spnudge_v_alpha(k,i,j) = uv_alpha
          elseif ( real_czxv(k,i,j) < spnudge_level4 ) then
             spnudge_v_alpha(k,i,j) = ( real_czxv(k,i,j) - spnudge_level3 ) / ( spnudge_level4 - spnudge_level3 ) * uv_alpha
          else
             spnudge_v_alpha(k,i,j) = 0.0_rp
          endif
       enddo
       enddo
       enddo
 
    else
 
       !$omp parallel do
       do j = js, je
       do i = is, ie
       do k = ks, ke
          spnudge_u_alpha(k,i,j) = 0.0_rp
          spnudge_v_alpha(k,i,j) = 0.0_rp
       end do
       end do
       end do
 
    endif
 
    if ( spnudge_pt .and. spnudge_pt_tau > 0.0_rp ) then
        pt_alpha = 1.0_rp / spnudge_pt_tau
 
        !$omp parallel do
        do j = js, je
        do i = is, ie
        do k = ks, ke
           if ( real_cz(k,i,j) < spnudge_level1 ) then
              spnudge_pt_alpha(k,i,j) = 0.0_rp
           elseif ( real_cz(k,i,j) < spnudge_level2 ) then
              spnudge_pt_alpha(k,i,j) = ( real_cz(k,i,j) - spnudge_level1 ) / ( spnudge_level2 - spnudge_level1 ) * pt_alpha
           elseif ( real_cz(k,i,j) < spnudge_level3 ) then
              spnudge_pt_alpha(k,i,j) = pt_alpha
           elseif ( real_cz(k,i,j) < spnudge_level4 ) then
              spnudge_pt_alpha(k,i,j) = ( real_cz(k,i,j) - spnudge_level3 ) / ( spnudge_level4 - spnudge_level3 ) * pt_alpha
           else
              spnudge_pt_alpha(k,i,j) = 0.0_rp
           endif
        enddo
        enddo
        enddo
 
     else
 
        !$omp parallel do
        do j = js, je
        do i = is, ie
        do k = ks, ke
           spnudge_pt_alpha(k,i,j) = 0.0_rp
        end do
        end do
        end do
 
    endif
 
    if ( spnudge_qv .and. spnudge_qv_tau > 0.0_rp ) then
        qv_alpha = 1.0_rp / spnudge_qv_tau
 
        !$omp parallel do
        do j = js, je
        do i = is, ie
        do k = ks, ke
           if ( real_cz(k,i,j) < spnudge_level1 ) then
              spnudge_qv_alpha(k,i,j) = 0.0_rp
           elseif ( real_cz(k,i,j) < spnudge_level2 ) then
              spnudge_qv_alpha(k,i,j) = ( real_cz(k,i,j) - spnudge_level1 ) / ( spnudge_level2 - spnudge_level1 ) * qv_alpha
           elseif ( real_cz(k,i,j) < spnudge_level3 ) then
              spnudge_qv_alpha(k,i,j) = qv_alpha
           elseif ( real_cz(k,i,j) < spnudge_level4 ) then
              spnudge_qv_alpha(k,i,j) = ( real_cz(k,i,j) - spnudge_level3 ) / ( spnudge_level4 - spnudge_level3 ) * qv_alpha
           else
              spnudge_qv_alpha(k,i,j) = 0.0_rp
           endif
        enddo
        enddo
        enddo
 
     else
 
        !$omp parallel do
        do j = js, je
        do i = is, ie
        do k = ks, ke
           spnudge_qv_alpha(k,i,j) = 0.0_rp
        end do
        end do
        end do
 
     end if
 
     call dft_setup( ka, ks, ke, ia, is, ie, ja, js, je, &
                     max( spnudge_uv_lm, spnudge_pt_lm, spnudge_qv_lm ), &
                     max( spnudge_uv_mm, spnudge_pt_mm, spnudge_qv_mm ) )
 
     return
   end subroutine spnudge_setup
 
end module scale_spnudge