| Class | HeatFlux |
| In: |
physics/heatflux.f90
|
下部境界でのフラックスの計算モジュール
| Function : | |||
| xz_HeatFluxBulk(DimXMin:DimXMax,DimZMin:DimZMax) : | real(8)
| ||
| xz_PotTemp(DimXMin:DimXMax,DimZMin:DimZMax) : | real(8), intent(in)
|
下部境界からのフラックスによる温度の変化率をバルク法に 基づいて計算する.
function xz_HeatFluxBulk( xz_PotTemp )
!
! 下部境界からのフラックスによる温度の変化率をバルク法に
! 基づいて計算する.
!
!暗黙の型宣言禁止
implicit none
!変数定義
real(8), intent(in) :: xz_PotTemp(DimXMin:DimXMax,DimZMin:DimZMax)
!温位の擾乱成分
real(8) :: xz_HeatFluxBulk(DimXMin:DimXMax,DimZMin:DimZMax)
!地表面熱フラックス
real(8) :: VelX = 100.0d0 !下層での水平速度を決めうち
!初期化
! * 全ての値をゼロに固定
xz_HeatFluxBulk = 0.0d0
!地表面熱フラックスを計算
! * 単位は K/s
! * 格子点 xz では, 物理領域の最下端の添え字は RegZMin+1
!
xz_HeatFluxBulk(:,RegZMin+1) = max( 0.0d0, - Bulk * VelX * xz_ExnerBasicZ(:,RegZMin+1) * xz_PotTemp(:,RegZMin+1) / DelZ )
end function xz_HeatFluxBulk
| Function : | |
| xz_HeatFluxDiff(DimXMin:DimXMax,DimZMin:DimZMax) : | real(8) |
| xz_PotTemp(DimXMin:DimXMax,DimZMin:DimZMax) : | real(8), intent(in) |
function xz_HeatFluxDiff( xz_PotTemp )
implicit none
real(8), intent(in) :: xz_PotTemp(DimXMin:DimXMax,DimZMin:DimZMax)
real(8) :: xz_HeatFluxDiff(DimXMin:DimXMax,DimZMin:DimZMax)
xz_HeatFluxDiff = 0.0d0
xz_HeatFluxDiff(:,RegZMin+1) = - Kappa * xz_PotTemp(:,RegZMin+1) * xz_ExnerBasicZ(:, RegZMin+1) / ( ( DelZ * 5.0d-1 ) ** 2.0d0 )
call StorePotTempFlux( xz_HeatFluxDiff )
end function xz_HeatFluxDiff
| Function : | |||
| xz_MixRtFluxBulk(DimXMin:DimXMax,DimZMin:DimZMax) : | real(8)
| ||
| xz_MixRt(DimXMin:DimXMax,DimZMin:DimZMax) : | real(8), intent(in)
|
下部境界からのフラックスによる凝結成分混合比の変化率をバルク法に 基づいて計算する.
function xz_MixRtFluxBulk( xz_MixRt )
!
! 下部境界からのフラックスによる凝結成分混合比の変化率をバルク法に
! 基づいて計算する.
!
!暗黙の型宣言禁止
implicit none
!変数定義
real(8), intent(in) :: xz_MixRt(DimXMin:DimXMax,DimZMin:DimZMax)
!温位の擾乱成分
real(8) :: xz_MixRtFluxBulk(DimXMin:DimXMax,DimZMin:DimZMax)
!地表面熱フラックス
real(8) :: VelX = 100.0d0 !下層での水平速度を決めうち
!初期化
! * 全ての値をゼロに固定
xz_MixRtFluxBulk = 0.0d0
xz_MixRtFluxBulk(:,RegZMin+1) = max( 0.0d0, - Bulk * VelX * xz_MixRt(:,RegZMin+1) / DelZ )
end function xz_MixRtFluxBulk
| Function : | |
| xz_MixRtFluxDiff(DimXMin:DimXMax,DimZMin:DimZMax) : | real(8) |
| xz_MixRt(DimXMin:DimXMax,DimZMin:DimZMax) : | real(8), intent(in) |
function xz_MixRtFluxDiff( xz_MixRt )
implicit none
real(8), intent(in) :: xz_MixRt(DimXMin:DimXMax,DimZMin:DimZMax)
real(8) :: xz_MixRtFluxDiff(DimXMin:DimXMax,DimZMin:DimZMax)
xz_MixRtFluxDiff = 0.0d0
xz_MixRtFluxDiff(:,RegZMin+1) = max( 0.0d0, - Kappa * xz_MixRt(:,RegZMin+1) / ( ( DelZ * 5.0d-1 ) ** 2.0d0 ) )
! xz_MixRtFluxDiff(:,RegZMin+1) = &
! & - Kappa * xz_MixRt(:,RegZMin+1) &
! & / ( ( DelZ * 5.0d-1 ) ** 2.0d0 )
end function xz_MixRtFluxDiff
| Function : | |
| xza_MixRtFluxDiff(DimXMin:DimXMax,DimZMin:DimZMax,SpcNum) : | real(8) |
| xza_MixRt(DimXMin:DimXMax,DimZMin:DimZMax,SpcNum) : | real(8), intent(in) |
function xza_MixRtFluxDiff( xza_MixRt )
implicit none
real(8), intent(in) :: xza_MixRt(DimXMin:DimXMax,DimZMin:DimZMax,SpcNum)
real(8) :: xza_MixRtFluxDiff(DimXMin:DimXMax,DimZMin:DimZMax,SpcNum)
xza_MixRtFluxDiff = 0.0d0
xza_MixRtFluxDiff(:,RegZMin+1,:) = - Kappa * xza_MixRt(:,RegZMin+1,:) / ( ( DelZ * 5.0d-1 ) ** 2.0d0 )
! xza_MixRtFluxDiff(:,RegZMin+1,:) = &
! & max( &
! & 0.0d0, &
! & - Kappa * xza_MixRt(:,RegZMin+1,:) &
! & / ( ( DelZ * 5.0d-1 ) ** 2.0d0 ) &
! & )
call StoreMixRtFlux( xza_MixRtFluxDiff )
end function xza_MixRtFluxDiff