Class FillNegative
In: ../src/utils/fillnegative.f90

負の雲水量などを穴埋めするためのルーチン

 * 負となった場合には, 周囲の 8 格子点の値を削って穴埋めする.

Methods

Included Modules

dc_types dc_message gtool_historyauto timeset gridset basicset composition

Public Instance methods

Subroutine :
xyz_CDensAl(imin:imax, jmin:jmax, kmin:kmax) :real(DP), intent(inout)

[Source]

  subroutine FillNegativeDensity(xyz_CDensAl)
    
    implicit none
    
    real(DP), intent(inout)  :: xyz_CDensAl(imin:imax, jmin:jmax, kmin:kmax)
    real(DP)                 :: xyz_CDensOrig(imin:imax, jmin:jmax, kmin:kmax)
    real(DP)                 :: xyz_CDensWork(imin:imax, jmin:jmax, kmin:kmax)
    real(DP)                 :: xyz_Del(imin:imax, jmin:jmax, kmin:kmax)


    ! 初期値を保管 Store Initial Value
    !
    xyz_CDensOrig = xyz_CDensAl
    
    ! 移流によって負になった部分を埋める
    ! Negative values due to advection are corrected.
    !
    xyz_CDensWork = xyz_CDensAl
    xyz_CDensAl   = xyz_FillNegativeDensity_xyz( xyz_CDensWork ) 
    
    ! 埋め切れなかった部分をゼロにする. 
    ! Negative values mixing ratios are corrected.
    !
    xyz_CDensWork = xyz_CDensAl
    xyz_CDensAl = max( 0.0d0, xyz_CDensWork )
    
    ! 埋めた/削った量を保管
    ! Adjustion value is stored.
    !    
    xyz_Del = ( xyz_CDensAl - xyz_CDensOrig ) / DelTimeShort

    ! Output
    call HistoryAutoPut(TimeN, 'CDensFill', xyz_Del(1:nx,1:ny,1:nz))
    
  end subroutine FillNegativeDensity
Subroutine :
xyzf_QMixAl(imin:imax, jmin:jmax, kmin:kmax, ncmax) :real(DP), intent(inout)

[Source]

  subroutine FillNegativeQMix(xyzf_QMixAl)
    
    implicit none
    
    real(DP), intent(inout)  :: xyzf_QMixAl(imin:imax, jmin:jmax, kmin:kmax, ncmax)
    real(DP)                 :: xyzf_QMixOrig(imin:imax, jmin:jmax, kmin:kmax, ncmax)
    real(DP)                 :: xyzf_QMixWork(imin:imax, jmin:jmax, kmin:kmax, ncmax)
    real(DP)                 :: xyza_Del(imin:imax, jmin:jmax, kmin:kmax, ncmax)
    real(DP)                 :: DelTime
    integer                  :: l

    ! 時間刻み幅
    !
    DelTime = 2.0d0 * DelTimeLong

    ! 初期値を保管 Store Initial Value
    !
    xyzf_QMixOrig = xyzf_QMixAl
    
    ! 移流によって負になった部分を埋める
    ! Negative values due to advection are corrected.
    !
    xyzf_QMixWork = xyzf_QMixAl
    xyzf_QMixAl   = xyza_FillNegative_xyza( xyzf_QMixWork ) 
    
    ! 埋め切れなかった部分をゼロにする. 
    ! Negative values mixing ratios are corrected.
    !
    xyzf_QMixWork = xyzf_QMixAl
    xyzf_QMixAl = max( - xyzf_QMixBZ, xyzf_QMixWork )

    ! 埋めた/削った量を保管
    ! Adjustion value is stored.
    !    
    xyza_Del = (xyzf_QMixAl - xyzf_QMixOrig) / DelTime

    ! Output
    do l = 1, ncmax
      call HistoryAutoPut(TimeN, trim(SpcWetSymbol(l))//'_Fill', xyza_Del(1:nx,1:ny,1:nz,l))
    end do

  end subroutine FillNegativeQMix
Subroutine :

基準値を取得. 基本場と擾乱成分の和がゼロよりも小さくなることは許容しない

[Source]

  subroutine FillNegative_Init( )    
    !
    ! 基準値を取得. 
    ! 基本場と擾乱成分の和がゼロよりも小さくなることは許容しない
    !

    implicit none

    !入力変数
    integer  :: l

    !初期化
    allocate( xyza_Basic(imin:imax,jmin:jmax,kmin:kmax, ncmax) )
    allocate( xyz_Dens  (imin:imax,jmin:jmax,kmin:kmax) )
    allocate( xza_Basic(imin:imax,kmin:kmax, ncmax) )
    allocate( xz_Dens(imin:imax,kmin:kmax ) )
    
    !値の代入
    xyza_Basic = xyzf_QMixBZ
    xyz_Dens   = xyz_DensBZ

    xza_Basic = xyzf_QMixBZ(:,1,:,:)
    xz_Dens   = xyz_DensBZ(:,1,:)

    call HistoryAutoAddVariable( varname="CDensFill", dims=(/'x','y','z','t'/), longname='Filling Negative term of cloud density ', units='kg.m-3.s-1', xtype='float')
    do l = 1, ncmax
      call HistoryAutoAddVariable( varname=trim(SpcWetSymbol(l))//'_Fill', dims=(/'x','y','z','t'/), longname='Filling Negative term of ' //trim(SpcWetSymbol(l))//' mixing ratio', units='kg.kg-1.s-1', xtype='float')
    end do

  end subroutine FillNegative_Init
xyz_Dens
Variable :
xyz_Dens(:,:,:) :real(DP), save, allocatable
Function :
xyz_FillNegativeDensity_xyz(imin:imax, jmin:jmax, kmin:kmax) :real(8)
xyz_Var(imin:imax, jmin:jmax, kmin:kmax) :real(8), intent(inout)

主成分凝結対流用. 雲密度に対してのみ適用可能. 行なっていることは木星版と同じ.

[Source]

  function xyz_FillNegativeDensity_xyz( xyz_Var )
    ! 主成分凝結対流用. 
    ! 雲密度に対してのみ適用可能. 
    ! 行なっていることは木星版と同じ. 

    implicit none

    real(8), intent(inout) :: xyz_Var(imin:imax, jmin:jmax, kmin:kmax)
    real(8)                :: xyz_FillNegativeDensity_xyz(imin:imax, jmin:jmax, kmin:kmax)
    real(8)                :: xyz_DQFILL(imin:imax, jmin:jmax, kmin:kmax)
    real(8)                :: xyz_QSUMPN(imin:imax, jmin:jmax, kmin:kmax)
    real(8)                :: xz_Var(imin:imax, kmin:kmax)
    real(8)                :: xz_DQFILL(imin:imax, kmin:kmax)
    real(8)                :: xz_QSUMPN(imin:imax, kmin:kmax)
    real(8), parameter     :: EPS = 1.0d-60  !零での割り算を防ぐ
    integer                :: i, j, k

    if (jmin==jmax) then     

      !初期化
!      xz_QSUMPN = 0.0d0
!      xz_DQFILL = 0.0d0
      xz_Var = xyz_Var(:,1,:)

      do k = kmin+2, kmax-2
        do i = imin+2, imax-2
          xz_QSUMPN(i,k) = 1.0d0 / ( (    MAX( 0.0d0, xz_Var(i-1,k) ) + MAX( 0.0d0, xz_Var(i+1,k) ) + MAX( 0.0d0, xz_Var(i,k-1) ) + MAX( 0.0d0, xz_Var(i,k+1) ) ) * 0.75d0 + (   MAX( 0.0d0, xz_Var(i-2,k) ) + MAX( 0.0d0, xz_Var(i+2,k) ) + MAX( 0.0d0, xz_Var(i,k-2) ) + MAX( 0.0d0, xz_Var(i,k+2) ) ) * 0.25d0 + EPS ) 
        end do
      end do

      do k = kmin+2, kmax-2
        do i = imin+2, imax-2
          xz_DQFILL(i,k) = - MIN( 0.0d0, xz_Var(i,k) ) + MAX( 0.0d0, xz_Var(i,k) ) * ( ( MIN( 0.0d0, xz_Var(i-1,k) ) * xz_QSUMPN(i-1,k) + MIN( 0.0d0, xz_Var(i+1,k) ) * xz_QSUMPN(i+1,k) + MIN( 0.0d0, xz_Var(i,k-1) ) * xz_QSUMPN(i,k-1) + MIN( 0.0d0, xz_Var(i,k+1) ) * xz_QSUMPN(i,k+1) ) * 0.75d0 + ( MIN( 0.0d0, xz_Var(i-2,k) ) * xz_QSUMPN(i-2,k) + MIN( 0.0d0, xz_Var(i+2,k) ) * xz_QSUMPN(i+2,k) + MIN( 0.0d0, xz_Var(i,k-2) ) * xz_QSUMPN(i,k-2) + MIN( 0.0d0, xz_Var(i,k+2) ) * xz_QSUMPN(i,k+2) ) * 0.25d0 )
        end do
      end do
      !出力
      xyz_FillNegativeDensity_xyz(:,1,:) = xz_Var + xz_DQFILL

    else

      !初期化
!      xyz_QSUMPN = 0.0d0
!      xyz_DQFILL = 0.0d0
      
      do k = kmin+2, kmax-2
        do j = jmin+2, jmax-2
          do i = imin+2, imax-2
            xyz_QSUMPN(i,j,k) = 1.0d0 / ( (    MAX( 0.0d0, xyz_Var(i-1,j,k) ) + MAX( 0.0d0, xyz_Var(i+1,j,k) ) + MAX( 0.0d0, xyz_Var(i,j-1,k) ) + MAX( 0.0d0, xyz_Var(i,j+1,k) ) + MAX( 0.0d0, xyz_Var(i,j,k-1) ) + MAX( 0.0d0, xyz_Var(i,j,k+1) ) ) * 0.75d0 + (   MAX( 0.0d0, xyz_Var(i-2,j,k) ) + MAX( 0.0d0, xyz_Var(i+2,j,k) ) + MAX( 0.0d0, xyz_Var(i,j-2,k ) ) + MAX( 0.0d0, xyz_Var(i,j+2,k) ) + MAX( 0.0d0, xyz_Var(i,j,k-2) ) + MAX( 0.0d0, xyz_Var(i,j,k+2) ) ) * 0.25d0 + EPS ) 
          end do
        end do
      end do

      do k = kmin+2, kmax-2
        do j = jmin+2, jmax-2
          do i = imin+2, imax-2
            xyz_DQFILL(i,j,k) = - MIN( 0.0d0, xyz_Var(i,j,k) ) + MAX( 0.0d0, xyz_Var(i,j,k) ) * ( ( MIN( 0.0d0, xyz_Var(i-1,j,k) ) * xyz_QSUMPN(i-1,j,k) + MIN( 0.0d0, xyz_Var(i+1,j,k) ) * xyz_QSUMPN(i+1,j,k) + MIN( 0.0d0, xyz_Var(i,j-1,k) ) * xyz_QSUMPN(i,j-1,k) + MIN( 0.0d0, xyz_Var(i,j+1,k) ) * xyz_QSUMPN(i,j+1,k) + MIN( 0.0d0, xyz_Var(i,j,k-1) ) * xyz_QSUMPN(i,j,k-1) + MIN( 0.0d0, xyz_Var(i,j,k+1) ) * xyz_QSUMPN(i,j,k+1) ) * 0.75d0 + ( MIN( 0.0d0, xyz_Var(i-2,j,k) ) * xyz_QSUMPN(i-2,j,k) + MIN( 0.0d0, xyz_Var(i+2,j,k) ) * xyz_QSUMPN(i+2,j,k) + MIN( 0.0d0, xyz_Var(i,j-2,k) ) * xyz_QSUMPN(i,j-2,k) + MIN( 0.0d0, xyz_Var(i,j+2,k) ) * xyz_QSUMPN(i,j+2,k) + MIN( 0.0d0, xyz_Var(i,j,k-2) ) * xyz_QSUMPN(i,j,k-2) + MIN( 0.0d0, xyz_Var(i,j,k+2) ) * xyz_QSUMPN(i,j,k+2) ) * 0.25d0 )
          end do
        end do
      end do
      
      !出力
      xyz_FillNegativeDensity_xyz = xyz_Var + xyz_DQFILL

    end if
    
  end function xyz_FillNegativeDensity_xyz
xyza_Basic
Variable :
xyza_Basic(:,:,:,:) :real(DP), save, allocatable
Function :
xyza_FillNegative_xyza(imin:imax,jmin:jmax,kmin:kmax, ncmax) :real(DP)
xyza_Var(imin:imax,jmin:jmax,kmin:kmax, ncmax) :real(DP), intent(in)

[Source]

  function xyza_FillNegative_xyza( xyza_Var )
    
    implicit none

    real(DP), intent(in) :: xyza_Var(imin:imax,jmin:jmax,kmin:kmax, ncmax)
    real(DP) :: xyza_FillNegative_xyza(imin:imax,jmin:jmax,kmin:kmax, ncmax)
    real(DP) :: xyza_DQFILL(imin:imax,jmin:jmax,kmin:kmax, ncmax)
    real(DP) :: xyza_QSUMPN(imin:imax,jmin:jmax,kmin:kmax, ncmax)
    real(DP) :: xza_Var(imin:imax, kmin:kmax, ncmax)
    real(DP) :: xza_DQFILL(imin:imax, kmin:kmax, ncmax)
    real(DP) :: xza_QSUMPN(imin:imax, kmin:kmax, ncmax)
    real(DP), parameter     :: EPS = 1.0d-60  !零での割り算を防ぐ
    integer                :: i, j, k, s

    
    if (jmin==jmax) then 

      !初期化
!      xza_QSUMPN = 0.0d0
      xza_DQFILL = 0.0d0
      xza_Var(:,:,:) = xyza_Var(:,1,:,:)
      
      do s = 1, ncmax
        do k = kmin+2, kmax-2
          do i = imin+2, imax-2
            xza_QSUMPN(i,k,s) = 1.0d0 / ( (    MAX( 0.0d0, xza_Basic(i-1,k,s) + xza_Var(i-1,k,s) ) * xz_Dens(i-1,k) + MAX( 0.0d0, xza_Basic(i+1,k,s) + xza_Var(i+1,k,s) ) * xz_Dens(i+1,k) + MAX( 0.0d0, xza_Basic(i,k-1,s) + xza_Var(i,k-1,s) ) * xz_Dens(i,k-1) + MAX( 0.0d0, xza_Basic(i,k+1,s) + xza_Var(i,k+1,s) ) * xz_Dens(i,k+1) ) * 0.75d0 + (   MAX( 0.0d0, xza_Basic(i-2,k,s) + xza_Var(i-2,k,s) ) * xz_Dens(i-2,k) + MAX( 0.0d0, xza_Basic(i+2,k,s) + xza_Var(i+2,k,s) ) * xz_Dens(i+2,k) + MAX( 0.0d0, xza_Basic(i,k-2,s) + xza_Var(i,k-2,s) ) * xz_Dens(i,k-2) + MAX( 0.0d0, xza_Basic(i,k+2,s) + xza_Var(i,k+2,s) ) * xz_Dens(i,k+2) ) * 0.25d0 + EPS ) 
          end do
        end do
      end do
      
      do s = 1, ncmax
        do k = 1, nz
          do i = 1, nx
            xza_DQFILL(i,k,s) = - MIN( 0.0d0, xza_Basic(i,k,s) + xza_Var(i,k,s) ) + MAX( 0.0d0, xza_Basic(i,k,s) + xza_Var(i,k,s) ) * ( ( MIN( 0.0d0, xza_Basic(i-1,k,s) + xza_Var(i-1,k,s) ) * xza_QSUMPN(i-1,k,s) * xz_Dens(i-1,k) + MIN( 0.0d0, xza_Basic(i+1,k,s) + xza_Var(i+1,k,s) ) * xza_QSUMPN(i+1,k,s) * xz_Dens(i+1,k) + MIN( 0.0d0, xza_Basic(i,k-1,s) + xza_Var(i,k-1,s) ) * xza_QSUMPN(i,k-1,s) * xz_Dens(i,k-1) + MIN( 0.0d0, xza_Basic(i,k+1,s) + xza_Var(i,k+1,s) ) * xza_QSUMPN(i,k+1,s) * xz_Dens(i,k+1) ) * 0.75d0 + ( MIN( 0.0d0, xza_Basic(i-2,k,s) + xza_Var(i-2,k,s) ) * xza_QSUMPN(i-2,k,s) * xz_Dens(i-2,k) + MIN( 0.0d0, xza_Basic(i+2,k,s) + xza_Var(i+2,k,s) ) * xza_QSUMPN(i+2,k,s) * xz_Dens(i+2,k) + MIN( 0.0d0, xza_Basic(i,k-2,s) + xza_Var(i,k-2,s) ) * xza_QSUMPN(i,k-2,s) * xz_Dens(i,k-2) + MIN( 0.0d0, xza_Basic(i,k+2,s) + xza_Var(i,k+2,s) ) * xza_QSUMPN(i,k+2,s) * xz_Dens(i,k+2) ) * 0.25d0 )
          end do
        end do
      end do
      
      !出力
      xyza_FillNegative_xyza(:,1,:,:) = xza_Var + xza_DQFILL
      
    else
      
!      !初期化
!      xyza_QSUMPN = 0.0d0
      xyza_DQFILL = 0.0d0
      
      do s = 1, ncmax
        do k = kmin+2, kmax-2
          do j = jmin+2, jmax-2
            do i = imin+2, imax-2
              xyza_QSUMPN(i,j,k,s) = 1.0d0 / ( (    MAX( 0.0d0, xyza_Basic(i-1,j,k,s) + xyza_Var(i-1,j,k,s) ) * xyz_Dens(i-1,j,k) + MAX( 0.0d0, xyza_Basic(i+1,j,k,s) + xyza_Var(i+1,j,k,s) ) * xyz_Dens(i+1,j,k) + MAX( 0.0d0, xyza_Basic(i,j-1,k,s) + xyza_Var(i,j-1,k,s) ) * xyz_Dens(i,j-1,k) + MAX( 0.0d0, xyza_Basic(i,j+1,k,s) + xyza_Var(i,j+1,k,s) ) * xyz_Dens(i,j+1,k) + MAX( 0.0d0, xyza_Basic(i,j,k-1,s) + xyza_Var(i,j,k-1,s) ) * xyz_Dens(i,j,k-1) + MAX( 0.0d0, xyza_Basic(i,j,k+1,s) + xyza_Var(i,j,k+1,s) ) * xyz_Dens(i,j,k+1) ) * 0.75d0 + (   MAX( 0.0d0, xyza_Basic(i-2,j,k,s) + xyza_Var(i-2,j,k,s) ) * xyz_Dens(i-2,j,k) + MAX( 0.0d0, xyza_Basic(i+2,j,k,s) + xyza_Var(i+2,j,k,s) ) * xyz_Dens(i+2,j,k) + MAX( 0.0d0, xyza_Basic(i,j-2,k,s) + xyza_Var(i,j-2,k,s) ) * xyz_Dens(i,j-2,k) + MAX( 0.0d0, xyza_Basic(i,j+2,k,s) + xyza_Var(i,j+2,k,s) ) * xyz_Dens(i,j+2,k) + MAX( 0.0d0, xyza_Basic(i,j,k-2,s) + xyza_Var(i,j,k-2,s) ) * xyz_Dens(i,j,k-2) + MAX( 0.0d0, xyza_Basic(i,j,k+2,s) + xyza_Var(i,j,k+2,s) ) * xyz_Dens(i,j,k+2) ) * 0.25d0 + EPS ) 
            end do
          end do
        end do
      end do
      
      do s = 1, ncmax
        do k = 1, nz
          do j = 1, ny
            do i = 1, nx
              xyza_DQFILL(i,j,k,s) = - MIN( 0.0d0, xyza_Basic(i,j,k,s) + xyza_Var(i,j,k,s) ) + MAX( 0.0d0, xyza_Basic(i,j,k,s) + xyza_Var(i,j,k,s) ) * ( ( MIN( 0.0d0, xyza_Basic(i-1,j,k,s) + xyza_Var(i-1,j,k,s) ) * xyza_QSUMPN(i-1,j,k,s) * xyz_Dens(i-1,j,k) + MIN( 0.0d0, xyza_Basic(i+1,j,k,s) + xyza_Var(i+1,j,k,s) ) * xyza_QSUMPN(i+1,j,k,s) * xyz_Dens(i+1,j,k) + MIN( 0.0d0, xyza_Basic(i,j-1,k,s) + xyza_Var(i,j-1,k,s) ) * xyza_QSUMPN(i,j-1,k,s) * xyz_Dens(i,j-1,k) + MIN( 0.0d0, xyza_Basic(i,j+1,k,s) + xyza_Var(i,j+1,k,s) ) * xyza_QSUMPN(i,j+1,k,s) * xyz_Dens(i,j+1,k) + MIN( 0.0d0, xyza_Basic(i,j,k-1,s) + xyza_Var(i,j,k-1,s) ) * xyza_QSUMPN(i,j,k-1,s) * xyz_Dens(i,j,k-1) + MIN( 0.0d0, xyza_Basic(i,j,k+1,s) + xyza_Var(i,j,k+1,s) ) * xyza_QSUMPN(i,j,k+1,s) * xyz_Dens(i,j,k+1) ) * 0.75d0 + ( MIN( 0.0d0, xyza_Basic(i-2,j,k,s) + xyza_Var(i-2,j,k,s) ) * xyza_QSUMPN(i-2,j,k,s) * xyz_Dens(i-2,j,k) + MIN( 0.0d0, xyza_Basic(i+2,j,k,s) + xyza_Var(i+2,j,k,s) ) * xyza_QSUMPN(i+2,j,k,s) * xyz_Dens(i+2,j,k) + MIN( 0.0d0, xyza_Basic(i,j-2,k,s) + xyza_Var(i,j-2,k,s) ) * xyza_QSUMPN(i,j-2,k,s) * xyz_Dens(i,j-2,k) + MIN( 0.0d0, xyza_Basic(i,j+2,k,s) + xyza_Var(i,j+2,k,s) ) * xyza_QSUMPN(i,j+2,k,s) * xyz_Dens(i,j+2,k) + MIN( 0.0d0, xyza_Basic(i,j,k-2,s) + xyza_Var(i,j,k-2,s) ) * xyza_QSUMPN(i,j,k-2,s) * xyz_Dens(i,j,k-2) + MIN( 0.0d0, xyza_Basic(i,j,k+2,s) + xyza_Var(i,j,k+2,s) ) * xyza_QSUMPN(i,j,k+2,s) * xyz_Dens(i,j,k+2) ) * 0.25d0 )
            end do
          end do
        end do
      end do

      !出力
      xyza_FillNegative_xyza = xyza_Var + xyza_DQFILL
      
    end if

  end function xyza_FillNegative_xyza
xz_Dens
Variable :
xz_Dens(:,:) :real(DP), save, allocatable
xza_Basic
Variable :
xza_Basic(:,:,:) :real(DP), save, allocatable