!= 乾燥対流調節 ! != Dry convective adjustment ! ! Authors:: Yoshiyuki O. Takahashi, Yasuhiro MORIKAWA, Yukiko YAMADA, Shin-ichi Takehiro ! Version:: $Id: dry_conv_adjust.f90,v 1.91 2025/09/17 22:00:00 takepiro Exp $ ! Tag Name:: $Name: $ ! Copyright:: Copyright (C) GFD Dennou Club, 2008-2025. All rights reserved. ! License:: See COPYRIGHT[link:../../../COPYRIGHT] ! module dry_conv_adjust ! != 乾燥対流調節スキーム ! != Dry convective adjustment ! ! Note that Japanese and English are described in parallel. ! !== Procedures List ! ! DryConvAdjust :: 乾燥対流調節 ! ------------ :: ------------ ! DryConvAdjust :: Dry convective adjustment ! !== NAMELIST ! ! NAMELIST#dry_conv_adjust_nml ! ! モジュール引用 ; USE statements ! ! 格子点設定 ! Grid points settings ! use gridset, only: imax, & ! 経度格子点数. ! Number of grid points in longitude & jmax, & ! 緯度格子点数. ! Number of grid points in latitude & kmax ! 鉛直層数. ! Number of vertical level ! 組成に関わる配列の設定 ! Settings of array for atmospheric composition ! use composition, only: & & ncmax ! 成分の数 ! Number of composition ! 種別型パラメタ ! Kind type parameter ! use dc_types, only: DP, & ! 倍精度実数型. Double precision. & STRING ! 文字列. Strings. ! NAMELIST ファイル入力に関するユーティリティ ! Utilities for NAMELIST file input ! use namelist_util, only: MaxNmlArySize ! NAMELIST から読み込む配列の最大サイズ. ! Maximum size of arrays loaded from NAMELIST ! メッセージ出力 ! Message output ! use dc_message, only: MessageNotify ! 宣言文 ; Declaration statements ! implicit none private ! 公開手続き ! Public procedure ! public:: DryConvAdjust public:: DryConvAdjustInit ! 公開変数 ! Public variables ! logical, save :: dry_conv_adjust_inited = .false. ! 初期設定フラグ. ! Initialization flag ! 非公開変数 ! Private variables ! real(DP), save:: CrtlRH ! 臨界相対湿度. ! Critical relative humidity integer, save:: ItrtMax ! イテレーション回数. ! Number of iteration real(DP), save:: AdjustCriterion(1:MaxNmlArySize) ! 調節を行う基準 (湿潤静的エネルギーの差の温度換算値) ! Criterion of adjustment (temperature difference ! equivalent to difference of dry static energy) logical, save :: FlagAdjustMom ! Flag for adjusting momentum logical, save :: FlagAdjustMR ! Flag for adjusting mixing ratio character(*), parameter:: module_name = 'dry_conv_adjust' ! モジュールの名称. ! Module name character(*), parameter:: version = & & '$Name: $' // & & '$Id: dry_conv_adjust.f90,v 1.9 2025/09/17 22:00:00 takepiro Exp $' ! モジュールのバージョン ! Module version contains !-------------------------------------------------------------------------------------- subroutine DryConvAdjust( & & xyz_Temp, xyz_U, xyz_V, xyzf_QMix, & ! (inout) & xyz_Press, xyr_Press & ! (in) & ) ! ! 乾燥対流調節スキームにより, 温度を調節. ! ! Adjust temperature by dry convective adjustment ! ! モジュール引用 ; USE statements ! ! 物理・数学定数設定 ! Physical and mathematical constants settings ! use constants0, only: & & GasRUniv ! $ R^{*} $ [J K-1 mol-1]. ! 普遍気体定数. Universal gas constant ! 物理定数設定 ! Physical constants settings ! use constants, only: & & Grav, & ! $ g $ [m s-2]. ! 重力加速度. ! Gravitational acceleration & GasRDry, & ! $ R $ [J kg-1 K-1]. ! 乾燥大気の気体定数. ! Gas constant of air & CpDry ! $ C_p $ [J kg-1 K-1]. ! 乾燥大気の定圧比熱. ! Specific heat of air at constant pressure ! 時刻管理 ! Time control ! use timeset, only: & & DelTime, & ! $ \Delta t $ & TimeN, & ! ステップ $ t $ の時刻. Time of step $ t $. & TimesetClockStart, TimesetClockStop ! ヒストリデータ出力 ! History data output ! use gtool_historyauto, only: HistoryAutoPut ! 宣言文 ; Declaration statements ! implicit none real(DP), intent(inout):: xyz_Temp (0:imax-1, 1:jmax, 1:kmax) ! $ T $ . 温度. Temperature real(DP), intent(inout):: xyz_U (0:imax-1, 1:jmax, 1:kmax) ! $ U $ . Eastward wind velocity real(DP), intent(inout):: xyz_V (0:imax-1, 1:jmax, 1:kmax) ! $ V $ . Northward wind velocity real(DP), intent(inout):: xyzf_QMix(0:imax-1, 1:jmax, 1:kmax, 1:ncmax) ! $ q $ . Mixing ratio real(DP), intent(in ):: xyz_Press(0:imax-1, 1:jmax, 1:kmax) ! $ p $ . 気圧 (整数レベル). ! Air pressure (full level) real(DP), intent(in ):: xyr_Press(0:imax-1, 1:jmax, 0:kmax) ! $ \hat{p} $ . 気圧 (半整数レベル). ! Air pressure (half level) ! 作業変数 ! Work variables ! real(DP):: xyz_DTempDt(0:imax-1, 1:jmax, 1:kmax) ! 温度変化率. ! Temperature tendency real(DP):: xyz_TempB (0:imax-1, 1:jmax, 1:kmax) ! 調節前の温度. ! Temperature before adjustment logical:: xy_Adjust (0:imax-1, 1:jmax) ! 今回調節されたか否か?. ! Whether it was adjusted this time or not? logical:: xy_AdjustB (0:imax-1, 1:jmax) ! 前回調節されたか否か?. ! Whether it was adjusted last time or not? real(DP):: xyz_DelPress(0:imax-1, 1:jmax, 1:kmax) ! $ \Delta p $ ! real(DP):: xyz_DelMass (0:imax-1, 1:jmax, 1:kmax) ! $ \Delta p / g $ ! real(DP):: xyr_ConvAdjustFactor(0:imax-1, 1:jmax, 0:kmax) ! $ \frac{1}{2} \frac{ R }{Cp} ! \frac{ p_{k} - p_{k+1} }{ p_{k+1/2} } $ real(DP):: TempEquivToExcEne ! Temperature equivalent to the excess dry static energy ! (Dry static energy difference devided by specific heat) real(DP):: DelTempUppLev ! k+1 番目の層における調節による温度の変化量. ! Temperature variation by adjustment at k+1 level real(DP):: DelTempLowLev ! k 番目の層における調節による温度の変化量. ! Temperature variation by adjustment at k level logical:: Adjust ! 今回全領域において一度でも調節されたか否か?. ! Whether it was adjusted even once in global ! this time or not? real(DP):: TempLowLevBefAdj ! Variables for check routine real(DP):: TempUppLevBefAdj real(DP):: ExchangeMass ! ! Mass transport real(DP):: ExchangeMassDenom ! ! Variable for mass transport calculation real(DP):: ExchangeMassLowLim ! ! Lower limit of mass transport calculation real(DP), parameter :: ExchangeMassLowLimTempDiff = 1.0d-5 ! ! Lower limit of temperature difference ! between two layers for mass transport ! calculation real(DP):: DelULowLev ! ! Eastward wind velocity change real(DP):: DelUUppLev ! ! Eastward wind velocity change real(DP):: DelVLowLev ! ! Eastward wind velocity change real(DP):: DelVUppLev ! ! Eastward wind velocity change real(DP):: f_DelQMixLowLev(1:ncmax) ! ! Mixing ratio change real(DP):: f_DelQMixUppLev(1:ncmax) ! ! Mixing ratio change integer:: i ! 経度方向に回る DO ループ用作業変数 ! Work variables for DO loop in longitude integer:: j ! 緯度方向に回る DO ループ用作業変数 ! Work variables for DO loop in latitude integer:: k ! 鉛直方向に回る DO ループ用作業変数 ! Work variables for DO loop in vertical direction integer:: itr ! イテレーション方向に回る DO ループ用作業変数 ! Work variables for DO loop in iteration direction ! 実行文 ; Executable statement ! ! 初期化確認 ! Initialization check ! if ( .not. dry_conv_adjust_inited ) then call MessageNotify( 'E', module_name, 'This module has not been initialized.' ) end if ! 計算時間計測開始 ! Start measurement of computation time ! call TimesetClockStart( module_name ) ! 調節前 "Temp" の保存 ! Store "Temp" before adjustment ! xyz_TempB = xyz_Temp ! Calculate some values used for dry convective adjustment ! do k = 1, kmax xyz_DelPress(:,:,k) = xyr_Press(:,:,k-1) - xyr_Press(:,:,k) end do xyz_DelMass = xyz_DelPress / Grav ! \frac{1}{2} \frac{ R }{Cp} \frac{ p_{k} - p_{k+1} }{ p_{k+1/2} } ! ! The value at k = 0 is not used. k = 0 xyr_ConvAdjustFactor(:,:,k) = 0.0_DP ! do k = 1, kmax-1 xyr_ConvAdjustFactor(:,:,k) = & & GasRDry / CpDry & & * ( xyz_Press(:,:,k) - xyz_Press(:,:,k+1) ) & & / xyr_Press(:,:,k) / 2.0_DP end do ! The value at k = kmax is not used. k = kmax xyr_ConvAdjustFactor(:,:,k) = 0.0d0 ! 調節 ! Adjustment ! xy_AdjustB = .true. ! 繰り返し ! Iteration ! do itr = 1, ItrtMax xy_Adjust = .false. do k = 1, kmax-1 do j = 1, jmax do i = 0, imax-1 if ( xy_AdjustB(i,j) ) then ! Temperature equivalent to the excess dry static energy ! (Dry static energy difference devided by specific heat) ! TempEquivToExcEne = & & xyz_Temp(i,j,k) - xyz_Temp(i,j,k+1) & & - xyr_ConvAdjustFactor(i,j,k) * ( xyz_Temp(i,j,k) + xyz_Temp(i,j,k+1) ) ! Check vertical gradient of dry static energy ! if ( TempEquivToExcEne > AdjustCriterion(itr) ) then DelTempUppLev = & & xyz_DelPress(i,j,k) * TempEquivToExcEne & & / ( xyr_ConvAdjustFactor(i,j,k) & & * ( xyz_DelPress(i,j,k ) & & - xyz_DelPress(i,j,k+1) ) & & + xyz_DelPress(i,j,k) + xyz_DelPress(i,j,k+1) ) DelTempLowLev = & & - xyz_DelPress(i,j,k+1) * DelTempUppLev & & / xyz_DelPress(i,j,k) !========= ! save temperature before adjustment !--------- TempLowLevBefAdj = xyz_Temp(i,j,k ) TempUppLevBefAdj = xyz_Temp(i,j,k+1) ! 温度の調節 ! Adjust temperature ! xyz_Temp(i,j,k ) = xyz_Temp(i,j,k ) + DelTempLowLev xyz_Temp(i,j,k+1) = xyz_Temp(i,j,k+1) + DelTempUppLev !========= ! check routine !--------- !!$ write( 6, * ) '=====' !!$ write( 6, * ) 'Energy difference before and after adjustment and each energy' !!$ write( 6, * ) & !!$ & ( CpDry * TempLowLevBefAdj ) & !!$ & * xyz_DelPress(i,j,k ) / Grav & !!$ & + ( CpDry * TempUppLevBefAdj ) & !!$ & * xyz_DelPress(i,j,k+1) / Grav & !!$ & - ( CpDry * xyz_Temp(i,j,k ) ) & !!$ & * xyz_DelPress(i,j,k ) / Grav & !!$ & - ( CpDry * xyz_Temp(i,j,k+1) ) & !!$ & * xyz_DelPress(i,j,k+1) / Grav, & !!$ & ( CpDry * TempLowLevBefAdj ) & !!$ & * xyz_DelPress(i,j,k ) / Grav, & !!$ & ( CpDry * TempUppLevBefAdj ) & !!$ & * xyz_DelPress(i,j,k+1) / Grav, & !!$ & ( CpDry * xyz_Temp(i,j,k ) ) & !!$ & * xyz_DelPress(i,j,k ) / Grav, & !!$ & ( CpDry * xyz_Temp(i,j,k+1) ) & !!$ & * xyz_DelPress(i,j,k+1) / Grav !!$ write( 6, * ) 'Difference of dry static energy after adjustment' !!$ write( 6, * ) & !!$ & ( CpDry * xyz_Temp(i,j,k ) ) & !!$ & - ( CpDry * xyz_Temp(i,j,k+1) ) & !!$ & - CpDry * xyr_ConvAdjustFactor(i,j,k) & !!$ & * ( xyz_Temp(i,j,k) + xyz_Temp(i,j,k+1) ), & !!$ & ( CpDry * xyz_Temp(i,j,k ) ), & !!$ & ( CpDry * xyz_Temp(i,j,k+1) ), & !!$ & - CpDry * xyr_ConvAdjustFactor(i,j,k) & !!$ & * ( xyz_Temp(i,j,k) + xyz_Temp(i,j,k+1) ) !========= ! ! Mass exchange ! Denominator ExchangeMassDenom = & & CpDry * ( TempLowLevBefAdj - TempUppLevBefAdj ) & & - GasRDry & & * ( TempLowLevBefAdj + TempUppLevBefAdj ) / 2.0_DP & & / xyr_Press(i,j,k) & & * ( xyz_Press(i,j,k) - xyz_Press(i,j,k+1) ) ExchangeMassLowLim = CpDry * ExchangeMassLowLimTempDiff ! If a static energy difference between two layers is smaller ! than a specified lower limit, momentum and mixing ratio are ! not mixed to ensure numerical stability. ! If the lower limit is zero, some calculations are unstable. ! (yot, 2013/10/02) if ( ExchangeMassDenom > ExchangeMassLowLim ) then ExchangeMass = & & - CpDry * DelTempLowLev / ExchangeMassDenom & & * xyz_DelMass(i,j,k) else ExchangeMass = 0.0_DP end if ! Limitation of amount of mass exchange not to ! reverse a gradient ExchangeMass = & & min( ExchangeMass, & & xyz_DelMass(i,j,k) * xyz_DelMass(i,j,k+1) & & / ( xyz_DelMass(i,j,k) + xyz_DelMass(i,j,k+1) ) & & ) if ( FlagAdjustMom ) then DelULowLev = & & ( xyz_U(i,j,k+1) - xyz_U(i,j,k) ) & & * ExchangeMass & & / xyz_DelMass(i,j,k ) DelUUppLev = & & - ( xyz_U(i,j,k+1) - xyz_U(i,j,k) ) & & * ExchangeMass & & / xyz_DelMass(i,j,k+1) DelVLowLev = & & ( xyz_V(i,j,k+1) - xyz_V(i,j,k) ) & & * ExchangeMass & & / xyz_DelMass(i,j,k ) DelVUppLev = & & - ( xyz_V(i,j,k+1) - xyz_V(i,j,k) ) & & * ExchangeMass & & / xyz_DelMass(i,j,k+1) xyz_U(i,j,k ) = xyz_U(i,j,k ) + DelULowLev xyz_U(i,j,k+1) = xyz_U(i,j,k+1) + DelUUppLev xyz_V(i,j,k ) = xyz_V(i,j,k ) + DelVLowLev xyz_V(i,j,k+1) = xyz_V(i,j,k+1) + DelVUppLev end if if ( FlagAdjustMR ) then f_DelQMixLowLev = & & ( xyzf_QMix(i,j,k+1,:) - xyzf_QMix(i,j,k,:) ) & & * ExchangeMass & & / xyz_DelMass(i,j,k ) f_DelQMixUppLev = & & - ( xyzf_QMix(i,j,k+1,:) - xyzf_QMix(i,j,k,:) ) & & * ExchangeMass & & / xyz_DelMass(i,j,k+1) xyzf_QMix(i,j,k ,:) = xyzf_QMix(i,j,k ,:) + f_DelQMixLowLev xyzf_QMix(i,j,k+1,:) = xyzf_QMix(i,j,k+1,:) + f_DelQMixUppLev end if ! 調節したか否か? ! Whether it was adjusted or not? ! xy_Adjust(i,j) = .true. end if end if end do end do end do Adjust = .false. do i = 0, imax-1 do j = 1, jmax xy_AdjustB(i,j) = xy_Adjust(i,j) Adjust = Adjust .or. xy_Adjust(i,j) end do end do if ( .not. Adjust ) exit end do ! 温度変化率 ! Calculate temperature tendency ! xyz_DTempDt = ( xyz_Temp - xyz_TempB ) / ( 2.0_DP * DelTime ) ! ヒストリデータ出力 ! History data output ! call HistoryAutoPut( TimeN, 'DTempDtDryConv', xyz_DTempDt ) ! 計算時間計測一時停止 ! Pause measurement of computation time ! call TimesetClockStop( module_name ) end subroutine DryConvAdjust !-------------------------------------------------------------------------------------- subroutine DryConvAdjustInit ! ! dry_conv_adjust モジュールの初期化を行います. ! NAMELIST#dry_conv_adjust_nml の読み込みはこの手続きで行われます. ! ! "dyr_conv_adjust" module is initialized. ! "NAMELIST#dry_conv_adjust_nml" is loaded in this procedure. ! ! モジュール引用 ; USE statements ! ! NAMELIST ファイル入力に関するユーティリティ ! Utilities for NAMELIST file input ! use namelist_util, only: namelist_filename, NmlutilMsg, NmlutilAryValid ! ファイル入出力補助 ! File I/O support ! use dc_iounit, only: FileOpen ! 種別型パラメタ ! Kind type parameter ! use dc_types, only: STDOUT ! 標準出力の装置番号. Unit number of standard output ! 文字列操作 ! Character handling ! use dc_string, only: StoA ! ヒストリデータ出力 ! History data output ! use gtool_historyauto, only: HistoryAutoAddVariable ! メッセージ制御 ! Message control ! use mpi_messagecntl, only : DoesOutputMPIMessage ! 宣言文 ; Declaration statements ! implicit none integer:: unit_nml ! NAMELIST ファイルオープン用装置番号. ! Unit number for NAMELIST file open integer:: iostat_nml ! NAMELIST 読み込み時の IOSTAT. ! IOSTAT of NAMELIST read ! NAMELIST 変数群 ! NAMELIST group name ! namelist /dry_conv_adjust_nml/ & & CrtlRH, ItrtMax, AdjustCriterion, & & FlagAdjustMom, FlagAdjustMR ! デフォルト値については初期化手続 "dry_conv_adjust#DryConvAdjustInit" ! のソースコードを参照のこと. ! ! Refer to source codes in the initialization procedure ! "dry_conv_adjust#DryConvAdjustInit" for the default values. ! ! 実行文 ; Executable statement ! if ( dry_conv_adjust_inited ) return ! デフォルト値の設定 ! Default values settings ! ! default values used in AGCM5 !!$ CrtlRH = 0.990_DP !!$ ItrtMax = 10 !!$ AdjustCriterion(1:ItrtMax) = & !!$ & (/ 0.01_DP, 0.02_DP, 0.02_DP, 0.05_DP, 0.05_DP, & !!$ & 0.10_DP, 0.10_DP, 0.20_DP, 0.20_DP, 0.40_DP /) ! CrtlRH = 1.0_DP ItrtMax = 10 AdjustCriterion(1:ItrtMax) = 0.0_DP FlagAdjustMom = .false. FlagAdjustMR = .false. ! NAMELIST の読み込み ! NAMELIST is input ! if ( trim(namelist_filename) /= '' ) then call FileOpen( unit_nml, & ! (out) & namelist_filename, mode = 'r' ) ! (in) rewind( unit_nml ) read( unit_nml, & ! (in) & nml = dry_conv_adjust_nml, & ! (out) & iostat = iostat_nml ) ! (out) close( unit_nml ) call NmlutilMsg( iostat_nml, module_name ) ! (in) if ( iostat_nml == 0 .AND. DoesOutputMPIMessage() ) write( STDOUT, nml = dry_conv_adjust_nml ) end if ! イテレーション回数, 不安定の許容誤差のチェック ! Check number of iteration, admissible error of unstability ! call NmlutilAryValid( module_name, & ! (in) & AdjustCriterion, 'AdjustCriterion', & ! (in) & ItrtMax, 'ItrtMax' ) ! (in) ! ヒストリデータ出力のためのへの変数登録 ! Register of variables for history data output ! call HistoryAutoAddVariable( 'DTempDtDryConv', & & (/ 'lon ', 'lat ', 'sig ', 'time' /), & & 'dry convective adjustment heating', 'K s-1' ) ! 印字 ; Print ! call MessageNotify( 'M', module_name, '----- Initialization Messages -----' ) call MessageNotify( 'M', module_name, ' CrtlRH = %f', d = (/ CrtlRH /) ) call MessageNotify( 'M', module_name, ' ItrtMax = %d', i = (/ ItrtMax /) ) call MessageNotify( 'M', module_name, ' AdjustCriterion = (/ %f /)', & & d = AdjustCriterion(1:ItrtMax), n = (/ ItrtMax /) ) call MessageNotify( 'M', module_name, ' FlagAdjustMom = %b', l = (/ FlagAdjustMom /) ) call MessageNotify( 'M', module_name, ' FlagAdjustMR = %b', l = (/ FlagAdjustMR /) ) call MessageNotify( 'M', module_name, '-- version = %c', c1 = trim(version) ) dry_conv_adjust_inited = .true. end subroutine DryConvAdjustInit !-------------------------------------------------------------------------------------- end module dry_conv_adjust