!= Chou and Lee (1996) による短波放射モデル ! != Short wave radiation model described by Chou and Lee (1996) ! ! Authors:: Yoshiyuki O. Takahashi ! Version:: $Id: rad_CL1996.f90,v 1.4 2011/06/19 11:12:42 yot Exp $ ! Tag Name:: $Name: $ ! Copyright:: Copyright (C) GFD Dennou Club, 2008. All rights reserved. ! License:: See COPYRIGHT[link:../../../COPYRIGHT] ! module rad_CL1996 ! != Chou and Lee (1996) による短波放射モデル ! != Short wave radiation model described by Chou and Lee (1996) ! ! Note that Japanese and English are described in parallel. ! ! 短波放射モデル. ! ! This is a model of short wave radiation. ! !== References ! ! Chou, M.-D., and K.-T. Lee, ! Parameterizations for the absorption of solar radiation by water vapor and ozone, ! J. Atmos. Sci., 53, 1203-1208, 1996. ! !== Procedures List ! !!$ ! RadiationFluxDennouAGCM :: 放射フラックスの計算 !!$ ! RadiationDTempDt :: 放射フラックスによる温度変化の計算 !!$ ! RadiationFluxOutput :: 放射フラックスの出力 !!$ ! RadiationFinalize :: 終了処理 (モジュール内部の変数の割り付け解除) !!$ ! ------------ :: ------------ !!$ ! RadiationFluxDennouAGCM :: Calculate radiation flux !!$ ! RadiationDTempDt :: Calculate temperature tendency with radiation flux !!$ ! RadiationFluxOutput :: Output radiation fluxes !!$ ! RadiationFinalize :: Termination (deallocate variables in this module) ! !== NAMELIST ! !!$ ! NAMELIST#radiation_DennouAGCM_nml ! ! USE statements ! ! ! Kind type parameter ! use dc_types, only: DP, & ! Double precision. & STRING, & ! Strings. & TOKEN ! Keywords. ! メッセージ出力 ! Message output ! use dc_message, only: MessageNotify ! ! 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 ! Declaration statements ! implicit none private ! ! Public procedure ! public :: RadCL1996NumBands public :: RadCL1996UVVISParams public :: RadCL1996IRH2ONumKDFBin public :: RadCL1996IRH2OKDFParams public :: RadCL1996ScaleH2OVapMass public :: RadCL1996Init integer , parameter :: nband1 = 8 ! * 14500 to 57143 cm-1 (0.175 to 0.70 micron) integer , parameter :: nband2 = 3 ! * 2600 to 14500 cm-1 (0.70-10 micron) integer , parameter :: nkdf = 10 real(DP), save :: a_UVVIFracSolarFlux(1:nband1) real(DP), save :: a_UVVISO3AbsCoef (1:nband1) real(DP), save :: a_UVVISRayScatCoef (1:nband1) real(DP), save :: a_IRH2Okdfk (1:nkdf) ! k real(DP), save :: aa_IRH2Okdfdgi(1:nkdf,nband1+1:nband1+nband2) ! k dist. func. logical , save :: rad_cl1996_inited data rad_cl1996_inited /.false./ ! Table 1 ! Unit of k is g-1 cm2. data a_IRH2Okdfk & & / & & 0.0010 , 0.0133 , 0.0422 , 0.1334 , 0.4217 , 1.334 , 5.623 , 31.62 , 177.8 , 1000.0 & & / ! Table 1 data aa_IRH2Okdfdgi & & / & & 0.20673, 0.03497, 0.03011, 0.02260, 0.01336, 0.00696, 0.00441, 0.00115, 0.00026, 0.00000, & ! 0.7 -1.22 & 0.08236, 0.01157, 0.01133, 0.01143, 0.01240, 0.01258, 0.01381, 0.00650, 0.00244, 0.00094, & ! 1.22-2.27 & 0.01074, 0.00360, 0.00411, 0.00421, 0.00389, 0.00326, 0.00499, 0.00465, 0.00245, 0.00145 & ! 2.27-10 & / !!$ data aa_KDFParams & !!$ & / & !!$ ! (0.7-10) !!$ & 0.29983, 0.05014, 0.04555, 0.03824, 0.02965, 0.02280, 0.02321, 0.01230, 0.00515, 0.00239 & !!$ & / !!$ data aa_UVVISBandParams & !!$ & / & !!$ & , & ! 0.175-0.225 (UV-C) !!$ & , & ! 0.225-0.245, 0.260-0.280 (UV-C) !!$ & , & ! 0.245-0.260 (UV-C) !!$ & , & ! 0.280-0.295 (UV-B) !!$ & , & ! 0.295-0.310 (UV-B) !!$ & , & ! 0.310-0.320 (UV-B) !!$ & , & ! 0.320-0.400 (UV-A) !!$ & & ! 0.400-0.700 (PAR) !!$ & / ! Table 3. data a_UVVIFracSolarFlux & & / 0.00057, 0.00367, 0.00083, 0.00417, 0.00600, 0.00556, 0.05913, 0.39081 / ! Table 3. ! Unit is (cm-atm)_{STP}^{-1} data a_UVVISO3AbsCoef & & / 30.47 , 187.24 , 301.92 , 42.83 , 7.09 , 1.25 , 0.0345 , 0.0539 / ! Table 3. ! Unit is (mb)^{-1} data a_UVVISRayScatCoef & & / 0.00604, 0.00170, 0.00222, 0.00132, 0.00107, 0.00091, 0.00055, 0.00012 / character(*), parameter:: module_name = 'rad_CL1996' ! モジュールの名称. ! Module name character(*), parameter:: version = & & '$Name: $' // & & '$Id: rad_CL1996.f90,v 1.4 2011/06/19 11:12:42 yot Exp $' ! モジュールのバージョン ! Module version contains !-------------------------------------------------------------------------------------- subroutine RadCL1996NumBands( & & nbands1, nbands2 & ! (out) & ) integer, intent(out) :: nbands1 integer, intent(out) :: nbands2 ! 初期化確認 ! Initialization check ! if ( .not. rad_cl1996_inited ) then call MessageNotify( 'E', module_name, 'This module has not been initialized.' ) end if nbands1 = nband1 nbands2 = nband2 end subroutine RadCL1996NumBands !-------------------------------------------------------------------------------------- subroutine RadCL1996IRH2ONumKDFBin( & & nbin & ! (out) & ) integer, intent(out) :: nbin ! 初期化確認 ! Initialization check ! if ( .not. rad_cl1996_inited ) then call MessageNotify( 'E', module_name, 'This module has not been initialized.' ) end if nbin = nkdf end subroutine RadCL1996IRH2ONumKDFBin !-------------------------------------------------------------------------------------- subroutine RadCL1996ScaleH2OVapMass( & & xyz_Temp, xyz_DelH2OVapMass, xyz_Press, & ! (in ) & xyz_DelH2OVapMassScaled & ! (out) & ) ! USE statements ! real(DP), intent(in ):: xyz_Temp (0:imax-1, 1:jmax, 1:kmax) real(DP), intent(in ):: xyz_DelH2OVapMass (0:imax-1, 1:jmax, 1:kmax) real(DP), intent(in ):: xyz_Press (0:imax-1, 1:jmax, 1:kmax) real(DP), intent(out):: xyz_DelH2OVapMassScaled(0:imax-1, 1:jmax, 1:kmax) ! ! Work variables ! real(DP), parameter :: H2OScaleIndex = 0.8_DP real(DP), parameter :: RefPress = 300.0d2 real(DP), parameter :: RefTemp = 240.0d0 ! 初期化確認 ! Initialization check ! if ( .not. rad_cl1996_inited ) then call MessageNotify( 'E', module_name, 'This module has not been initialized.' ) end if xyz_DelH2OVapMassScaled = & & ( xyz_Press / RefPress )**H2OScaleIndex & & * exp( 0.00135_DP * ( xyz_Temp - RefTemp ) ) & & * xyz_DelH2OVapMass end subroutine RadCL1996ScaleH2OVapMass !-------------------------------------------------------------------------------------- subroutine RadCL1996IRH2OKDFParams( & & iband, ikdfbin, & ! (in ) & KDFAbsCoef, KDFWeight & ! (out) & ) ! USE statements ! integer , intent(in ):: iband integer , intent(in ):: ikdfbin real(DP), intent(out):: KDFAbsCoef real(DP), intent(out):: KDFWeight ! ! Work variables ! integer :: l integer :: m ! 初期化確認 ! Initialization check ! if ( .not. rad_cl1996_inited ) then call MessageNotify( 'E', module_name, 'This module has not been initialized.' ) end if l = iband m = ikdfbin KDFAbsCoef = a_IRH2Okdfk (m) KDFWeight = aa_IRH2Okdfdgi(m,l) end subroutine RadCL1996IRH2OKDFParams !-------------------------------------------------------------------------------------- subroutine RadCL1996UVVISParams( & & iband, & ! (in ) & UVVISFracSolarFlux, UVVISO3AbsCoef, UVVISRayScatCoef & ! (out) & ) ! USE statements ! integer , intent(in ):: iband real(DP), intent(out):: UVVISFracSolarFlux real(DP), intent(out):: UVVISO3AbsCoef real(DP), intent(out):: UVVISRayScatCoef ! ! Work variables ! integer :: l ! 初期化確認 ! Initialization check ! if ( .not. rad_cl1996_inited ) then call MessageNotify( 'E', module_name, 'This module has not been initialized.' ) end if l = iband UVVISFracSolarFlux = a_UVVIFracSolarFlux(l) UVVISO3AbsCoef = a_UVVISO3AbsCoef (l) UVVISRayScatCoef = a_UVVISRayScatCoef (l) end subroutine RadCL1996UVVISParams !-------------------------------------------------------------------------------------- subroutine RadCL1996Init !!$ ! NAMELIST ファイル入力に関するユーティリティ !!$ ! Utilities for NAMELIST file input !!$ ! !!$ use namelist_util, only: NmlutilMsg !!$ ! ファイル入出力補助 !!$ ! File I/O support !!$ ! !!$ use dc_iounit, only: FileOpen !!$ ! ヒストリデータ出力 !!$ ! History data output !!$ ! !!$ use gtool_historyauto, only: HistoryAutoAddVariable !!$ integer:: unit_nml ! NAMELIST ファイルオープン用装置番号. !!$ ! Unit number for NAMELIST file open !!$ integer:: iostat_nml ! NAMELIST 読み込み時の IOSTAT. !!$ ! IOSTAT of NAMELIST read !!$ ! NAMELIST 変数群 !!$ ! NAMELIST group name !!$ ! !!$ namelist /rad_CL1996_nml/ !& !!$ & ShortAtmosAlbedo !!$ ! !!$ ! デフォルト値については初期化手続 "rad_CL1996#RadCL1996Init" !!$ ! のソースコードを参照のこと. !!$ ! !!$ ! Refer to source codes in the initialization procedure !!$ ! "rad_LH74#RadLH74Init" for the default values. !!$ ! if ( rad_cl1996_inited ) return ! デフォルト値の設定 ! Default values settings ! !!$ ! 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 = rad_CL1996_nml, & ! (out) !!$ & iostat = iostat_nml ) ! (out) !!$ close( unit_nml ) !!$ !!$ call NmlutilMsg( iostat_nml, module_name ) ! (in) !!$ end if ! Unit is changed of k from g-1 cm2 to kg-1 m2. ! a_IRH2Okdfk = a_IRH2Okdfk * 1.0d3 * 1.0d-4 ! Convert unit of O3 absorption coefficient from (cm-atm)^{-1} to (kg m-2)^{-1} ! In order to convert unit from cm-1 atm-1 to m2 kg-1, multiply by ! 1.0d2 / 101325.0d0 * 8.31432d0 / ( 48.0d-3 ) * 273.15d0. ! a_UVVISO3AbsCoef = a_UVVISO3AbsCoef & & * 1.0d2 / 101325.0d0 * 8.31432d0 / ( 48.0d-3 ) * 273.15d0 ! Convert unit of Rayleigh scattering coefficient from (mb)^{-1} to (kg m-2)^{-1} ! memo: ! 1 mbar = 1e2 Pa ! 1e-2 mbar = 1 Pa ! ! In order to convert unit from (mb)^{-1} to (kg m-2)^{-1}, scattering ! coefficient is multiplied by 1.0d-1 * g, where g is the gravitational ! acceleration. In the old version of the code, a variable Grav in constants ! module is used for gravitational acceleration, as is shown in a line below. ! However, this variable is replaced with a constant value 9.8_DP. ! This is because the value of Grav may be changed to a value which is ! different from the Earth's value. (Of course, the scattering coefficient ! in unit of (kg m-2)^{-1} is independent on gravitational acceleration.) ! !!$ a_UVVISRayScatCoef = a_UVVISRayScatCoef * 1.0d-2 * Grav a_UVVISRayScatCoef = a_UVVISRayScatCoef * 1.0d-2 * 9.8_DP ! 印字 ; Print ! call MessageNotify( 'M', module_name, '----- Initialization Messages -----' ) !!$ call MessageNotify( 'M', module_name, 'ShortAtmosAlbedo = %f', d = (/ ShortAtmosAlbedo /) ) call MessageNotify( 'M', module_name, '-- version = %c', c1 = trim(version) ) rad_cl1996_inited = .true. end subroutine RadCL1996Init !-------------------------------------------------------------------------------------- end module rad_CL1996