!= Module NumDiffusion_3D ! ! Authors:: SUGIYAMA Ko-ichiro, ODAKA Masatsugu ! Version:: $Id: numdiffusion_3d.f90,v 1.4 2008-06-19 17:06:10 odakker Exp $ ! Tag Name:: $Name: arare4-20080627 $ ! Copyright:: Copyright (C) GFD Dennou Club, 2006. All rights reserved. ! License:: See COPYRIGHT[link:../../COPYRIGHT] ! !== Overview ! !数値拡散項の計算モジュール ! !== Error Handling ! !== Bugs ! !== Note ! !== Future Plans ! ! module NumDiffusion_3d ! !数値拡散項の計算モジュール !現在は 2 次精度中心差分を利用 ! !モジュール読み込み use dc_types, only: DP use dc_message, only: MessageNotify use gridset_3d, only: DimXMin, &! x 方向の配列の下限 & DimXMax, &! x 方向の配列の上限 & DimYMin, &! y 方向の配列の下限 & DimYMax, &! y 方向の配列の上限 & DimZMin, &! z 方向の配列の下限 & DimZMax, &! z 方向の配列の上限 & SpcNum, &! & x_dx, &! X 方向の格子点間隔 & y_dy, &! X 方向の格子点間隔 & z_dz ! Z 方向の格子点間隔 use timeset, only: DelTimeLong !長い時間ステップ use xyz_deriv_module, only: xyz_dx_pyz, xyz_dy_xqz, xyz_dz_xyr, & & pyz_dx_xyz, xqz_dy_xyz, xyr_dz_xyz, & & xyr_dx_pyr, pyz_dy_pqz, pyz_dz_pyr, & & pyr_dx_xyr, pqz_dy_pyz, pyr_dz_pyz, & & xqz_dx_pqz, xyr_dy_xqr, xqz_dz_xqr, & & pqz_dx_xqz, xqr_dy_xyr, xqr_dz_xqz use StorePotTemp_3d, only: StorePotTempDiff use StoreMixRt_3d, only: StoreMixRtDiff !暗黙の型宣言禁止 implicit none !属性の指定 private !関数を public にする public NumDiffusion_Init public xyz_NumDiffScalar public xyz_NumDiffKm public xyza_NumDiffScalar public pyz_NumDiffVelX public xqz_NumDiffVelY public xyr_NumDiffVelZ !変数の定義 real(DP) :: NuH = 0.0d0 !数値粘性の係数 (水平方向) real(DP) :: NuV = 0.0d0 !数値粘性の係数 (鉛直方向) real(DP) :: Alpha = 1.0d-4 save NuH, NuV, Alpha contains !!!-------------------------------------------------------------------!!! 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 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 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 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 !!!-------------------------------------------------------------------!!! 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(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(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 end module NumDiffusion_3d