Class NumDiffusion_3d
In: util/numdiffusion_3d.f90

数値拡散項の計算モジュール 現在は 2 次精度中心差分を利用

Methods

Included Modules

dc_types dc_message gridset_3d timeset xyz_deriv_module StorePotTemp_3d StoreMixRt_3d

Public Instance methods

Subroutine :

NumDiffusion モジュールの初期化ルーチン

[Source]

  subroutine NumDiffusion_init()
    !
    ! NumDiffusion モジュールの初期化ルーチン
    !

    !暗黙の型宣言禁止
    implicit none

    real(DP) :: DelXMin, DelYMin, DelZMin

    DelXMin = minval(x_dx)
    DelYMin = minval(y_dy)
    DelZMin = minval(z_dz)
    
    ! CReSS マニュアルでは, 2 次精度中心差分の場合, Alpha < 1/8 くらいが適当と述べている. 
    ! deepconv では NuH として 500 以上の値が入っていたそうだ,
    NuH  = Alpha * ( SQRT(DelXMin*DelYMin) ** 2.0d0 ) / DelTimeLong
    NuV  = Alpha * ( DelZMin ** 2.0d0 ) / DelTimeLong
    
    !確認

    call MessageNotify( "M", "NumDiffusion_init", "NuH = %f", d=(/NuH/) )
    call MessageNotify( "M", "NumDiffusion_init", "NuV = %f", d=(/NuV/) )

  end subroutine NumDiffusion_init
Function :
pyz_NumDiffVelX(DimXMin:DimXMax,DimYMin:DImYMax,DimZMin:DimZMax) :real(DP)
: 数値拡散
pyz_VarX(DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax) :real(DP), intent(in)
: 物理量

z 方向に半格子ずれた点での数値拡散項を評価

[Source]

  function pyz_NumDiffVelX(pyz_VarX)
    !
    ! z 方向に半格子ずれた点での数値拡散項を評価
    !
    
    !暗黙の型宣言禁止
    implicit none
    
    !変数定義
    real(DP), intent(in) :: pyz_VarX (DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax)
                                                    !物理量
    real(DP)             :: pyz_NumDiffVelX (DimXMin:DimXMax,DimYMin:DImYMax,DimZMin:DimZMax)
                                                    !数値拡散
    
    pyz_NumDiffVelX = 0.0d0  !初期化
    pyz_NumDiffVelX = NuH * ( pyz_dx_xyz( xyz_dx_pyz( pyz_VarX ) ) ) + NuH * ( pyz_dy_pqz( pqz_dy_pyz( pyz_VarX ) ) ) + NuV * ( pyz_dz_pyr( pyr_dz_pyz( pyz_VarX ) ) )
    
  end function pyz_NumDiffVelX
Function :
xqz_NumDiffVelY(DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax) :real(DP)
: 数値拡散
xqz_VarY(DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax) :real(DP), intent(in)
: 物理量

z 方向に半格子ずれた点での数値拡散項を評価

[Source]

  function xqz_NumDiffVelY(xqz_VarY)
    !
    ! z 方向に半格子ずれた点での数値拡散項を評価
    !
    
    !暗黙の型宣言禁止
    implicit none
    
    !変数定義
    real(DP), intent(in) :: xqz_VarY (DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax)
                                                    !物理量
    real(DP)             :: xqz_NumDiffVelY (DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax)
                                                    !数値拡散
    
!   xqz_NumDiffVelY = 0.0d0  !初期化
    xqz_NumDiffVelY = NuH * ( xqz_dx_pqz( pqz_dx_xqz( xqz_VarY ) ) ) + NuH * ( xqz_dy_xyz( xyz_dy_xqz( xqz_VarY ) ) ) + NuV * ( xqz_dz_xqr( xqr_dz_xqz( xqz_VarY ) ) )
    
  end function xqz_NumDiffVelY
Function :
xyr_NumDiffVelZ(DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax) :real(DP)
: 数値拡散
xyr_VarZ(DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax) :real(DP), intent(in)
: 物理量

x 方向に半格子ずれた点での数値拡散項を評価

[Source]

  function xyr_NumDiffVelZ(xyr_VarZ)
    !
    ! x 方向に半格子ずれた点での数値拡散項を評価
    !
    
    !暗黙の型宣言禁止
    implicit none
    
    !変数定義
    real(DP), intent(in) :: xyr_VarZ(DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax)
                                                    !物理量
    real(DP)             :: xyr_NumDiffVelZ(DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax)
                                                    !数値拡散

!   xyr_NumDiffVelZ = 0.0d0 
    xyr_NumDiffVelZ = NuH * ( xyr_dx_pyr( pyr_dx_xyr( xyr_VarZ ) ) ) + NuH * ( xyr_dy_xqr( xqr_dy_xyr( xyr_VarZ ) ) ) + NuV * ( xyr_dz_xyz( xyz_dz_xyr( xyr_VarZ ) ) )
    
  end function xyr_NumDiffVelZ
Function :
xyz_NumDiffKm(DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax) :real(8)
: 水平方向の数値拡散
xyz_Scalar(DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax) :real(8), intent(in)
: スカラー量

x, z 方向に半格子ずれた点での数値拡散項を評価

[Source]

  function xyz_NumDiffKm(xyz_Scalar)
    !
    ! x, z 方向に半格子ずれた点での数値拡散項を評価
    !
    
    !暗黙の型宣言禁止
    implicit none

    !変数定義
    real(8), intent(in) :: xyz_Scalar(DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax)
                                                    !スカラー量
    real(8)             :: xyz_NumDiffKm(DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax)
                                                    !水平方向の数値拡散
    
    xyz_NumDiffKm = NuH * (xyz_dx_pyz(pyz_dx_xyz( xyz_Scalar ))) + NuH * (xyz_dy_xqz(xqz_dy_xyz( xyz_Scalar ))) + NuV * (xyz_dz_xyr(xyr_dz_xyz( xyz_Scalar ))) 
    
  end function xyz_NumDiffKm
Function :
xyz_NumDiffScalar(DimXMin:DimXMax,DimYMin:DImYMax,DimZMin:DimZMax) :real(DP)
: 水平方向の数値拡散
xyz_Scalar(DimXMin:DimXMax,DimYMin:DImYMax,DimZMin:DimZMax) :real(DP), intent(in)
: スカラー量

x, y, z 方向に半格子ずれた点での数値拡散項を評価

[Source]

  function xyz_NumDiffScalar(xyz_Scalar)
    !
    ! x, y, z 方向に半格子ずれた点での数値拡散項を評価
    !
    
    !暗黙の型宣言禁止
    implicit none

    !変数定義
    real(DP), intent(in) :: xyz_Scalar (DimXMin:DimXMax,DimYMin:DImYMax,DimZMin:DimZMax)
                                                    !スカラー量
    real(DP)             :: xyz_NumDiffScalar (DimXMin:DimXMax,DimYMin:DImYMax,DimZMin:DimZMax)
                                                    !水平方向の数値拡散
    
!    xz_NumDiffScalar = 0.0d0  !初期化
    xyz_NumDiffScalar = NuH * (xyz_dx_pyz(pyz_dx_xyz( xyz_Scalar ))) + NuH * (xyz_dy_xqz(xqz_dy_xyz( xyz_Scalar ))) + NuV * (xyz_dz_xyr(xyr_dz_xyz( xyz_Scalar ))) 
    
    call StorePotTempDiff( xyz_NumDiffScalar )

  end function xyz_NumDiffScalar
Function :
xyza_NumDiffScalar(DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax,SpcNum) :real(DP)
: 水平方向の数値拡散
xyza_Scalar(DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax,SpcNum) :real(DP), intent(in)
: スカラー量

x, z 方向に半格子ずれた点での数値拡散項を評価

[Source]

  function xyza_NumDiffScalar(xyza_Scalar)
    !
    ! x, z 方向に半格子ずれた点での数値拡散項を評価
    !
    
    !暗黙の型宣言禁止
    implicit none!

    !変数定義
    real(DP), intent(in) :: xyza_Scalar (DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax,SpcNum)
                                                    !スカラー量
    real(DP)             :: xyza_NumDiffScalar (DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax,SpcNum)
                                                    !水平方向の数値拡散
    integer             :: s
    

    do s = 1, SpcNum
      xyza_NumDiffScalar(:,:,:,s) = NuH * (xyz_dx_pyz(pyz_dx_xyz( xyza_Scalar(:,:,:,s) ))) + NuH * (xyz_dy_xqz(xqz_dy_xyz( xyza_Scalar(:,:,:,s) ))) + NuV * (xyz_dz_xyr(xyr_dz_xyz( xyza_Scalar(:,:,:,s) ))) 
    end do

    call StoreMixRtDiff( xyza_NumDiffScalar )

  end function xyza_NumDiffScalar