Class NumDiffusion
In: util/numdiffusion.f90

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

Methods

Included Modules

gridset timeset differentiate_center2 StoreSet

Public Instance methods

Subroutine :

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

[Source]

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

    !暗黙の型宣言禁止
    implicit none
    
    ! CReSS マニュアルでは, 2 次精度中心差分の場合, Alpha < 1/8 くらいが適当と述べている. 
    ! deepconv では NuH として 500 以上の値が入っていたそうだ,
    NuH  = Alpha * ( DelX ** 2.0d0 ) / DelTimeLong
    NuV  = Alpha * ( DelZ ** 2.0d0 ) / DelTimeLong
    
    !確認
    write(*,*) "NumDiffusion_init, NuH: ", NuH
    write(*,*) "NumDiffusion_init, NuV: ", NuV

  end subroutine NumDiffusion_init
Function :
pz_NumDiffVelX(DimXMin:DimXMax, DimZMin:DimZMax) :real(8)
: 数値拡散
pz_VarX(DimXMin:DimXMax, DimZMin:DimZMax) :real(8), intent(in)
: 物理量

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

[Source]

  function pz_NumDiffVelX(pz_VarX)
    !
    ! z 方向に半格子ずれた点での数値拡散項を評価
    !
    
    !暗黙の型宣言禁止
    implicit none
    
    !変数定義
    real(8), intent(in) :: pz_VarX(DimXMin:DimXMax, DimZMin:DimZMax)
                                                    !物理量
    real(8)             :: pz_NumDiffVelX(DimXMin:DimXMax, DimZMin:DimZMax)
                                                    !数値拡散
    
    pz_NumDiffVelX = 0.0d0  !初期化
    pz_NumDiffVelX = NuH * ( pz_dx_xz( xz_dx_pz( pz_VarX ) ) ) + NuV * ( pz_dz_pr( pr_dz_pz( pz_VarX ) ) )
    
  end function pz_NumDiffVelX
Function :
xr_NumDiffVelZ(DimXMin:DimXMax, DimZMin:DimZMax) :real(8)
: 数値拡散
xr_VarZ(DimXMin:DimXMax, DimZMin:DimZMax) :real(8), intent(in)
: 物理量

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

[Source]

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

    xr_NumDiffVelZ = 0.0d0 
    xr_NumDiffVelZ = NuH * ( xr_dx_pr( pr_dx_xr( xr_VarZ ) ) ) + NuV * ( xr_dz_xz( xz_dz_xr( xr_VarZ ) ) )
    
  end function xr_NumDiffVelZ
Function :
xz_NumDiffScalar(DimXMin:DimXMax, DimZMin:DimZMax) :real(8)
: 水平方向の数値拡散
xz_Scalar(DimXMin:DimXMax, DimZMin:DimZMax) :real(8), intent(in)
: スカラー量

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

[Source]

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

    !変数定義
    real(8), intent(in) :: xz_Scalar(DimXMin:DimXMax, DimZMin:DimZMax)
                                                    !スカラー量
    real(8)             :: xz_NumDiffScalar(DimXMin:DimXMax, DimZMin:DimZMax)
                                                    !水平方向の数値拡散
    
    xz_NumDiffScalar = 0.0d0  !初期化
    xz_NumDiffScalar = NuH * (xz_dx_pz(pz_dx_xz( xz_Scalar ))) + NuV * (xz_dz_xr(xr_dz_xz( xz_Scalar ))) 
    
  end function xz_NumDiffScalar
Function :
xz_NumDiffScalar2(DimXMin:DimXMax, DimZMin:DimZMax) :real(8)
: 水平方向の数値拡散
xz_Scalar(DimXMin:DimXMax, DimZMin:DimZMax) :real(8), intent(in)
: スカラー量

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

[Source]

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

    !変数定義
    real(8), intent(in) :: xz_Scalar(DimXMin:DimXMax, DimZMin:DimZMax)
                                                    !スカラー量
    real(8)             :: xz_NumDiffScalar2(DimXMin:DimXMax, DimZMin:DimZMax)
                                                    !水平方向の数値拡散
    
    xz_NumDiffScalar2 = 0.0d0  !初期化
    xz_NumDiffScalar2 = NuH * (xz_dx_pz(pz_dx_xz( xz_Scalar ))) + NuV * (xz_dz_xr(xr_dz_xz( xz_Scalar ))) 
    
  end function xz_NumDiffScalar2

[Validate]