Class | cloud_utils |
In: |
radiation/cloud_utils.f90
|
Note that Japanese and English are described in parallel.
�²ã����å¸���è¨å�.
In this module, the amount of cloud or cloud optical depth are set. This module is under development and is still a preliminary version.
!$ ! RadiationFluxDennouAGCM : | �¾å����������¹ã���ç®� |
!$ ! ———— : | ———— |
!$ ! RadiationFluxDennouAGCM : | Calculate radiation flux |
Subroutine : | |||||
ParentRoutineName : | character(*), intent(in)
| ||||
xyr_Press(0:imax-1, 1:jmax, 0:kmax) : | real(DP), intent(in) | ||||
xyz_TempB(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in) | ||||
xyz_QH2OVapB(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in) | ||||
xyz_QH2OLiqB(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in) | ||||
xyz_QH2OSolB(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in) | ||||
xyz_Temp(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in) | ||||
xyz_QH2OVap(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in) | ||||
xyz_QH2OLiq(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in) | ||||
xyz_QH2OSol(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in) | ||||
xy_Rain(0:imax-1, 1:jmax) : | real(DP), intent(in) | ||||
xy_Snow(0:imax-1, 1:jmax) : | real(DP), intent(in) |
subroutine CloudUtilConsChk( ParentRoutineName, xyr_Press, xyz_TempB, xyz_QH2OVapB, xyz_QH2OLiqB, xyz_QH2OSolB, xyz_Temp , xyz_QH2OVap , xyz_QH2OLiq , xyz_QH2OSol , xy_Rain, xy_Snow ) ! ���»ç��� ! Time control ! use timeset, only: DelTime ! $ \Delta t $ [s] ! ����å®��°è¨å®� ! Physical constants settings ! use constants, only: Grav, CpDry, LatentHeat, LatentHeatFusion ! $ L $ [J kg-1] . ! ��解ã�����. ! Latent heat of fusion character(*), intent(in) :: ParentRoutineName !!$ logical , intent(in) :: FlagIncludeSnowPhaseChange real(DP), intent(in) :: xyr_Press (0:imax-1, 1:jmax, 0:kmax) real(DP), intent(in) :: xyz_TempB (0:imax-1, 1:jmax, 1:kmax) real(DP), intent(in) :: xyz_QH2OVapB(0:imax-1, 1:jmax, 1:kmax) real(DP), intent(in) :: xyz_QH2OLiqB(0:imax-1, 1:jmax, 1:kmax) real(DP), intent(in) :: xyz_QH2OSolB(0:imax-1, 1:jmax, 1:kmax) real(DP), intent(in) :: xyz_Temp (0:imax-1, 1:jmax, 1:kmax) real(DP), intent(in) :: xyz_QH2OVap (0:imax-1, 1:jmax, 1:kmax) real(DP), intent(in) :: xyz_QH2OLiq (0:imax-1, 1:jmax, 1:kmax) real(DP), intent(in) :: xyz_QH2OSol (0:imax-1, 1:jmax, 1:kmax) real(DP), intent(in) :: xy_Rain (0:imax-1, 1:jmax) real(DP), intent(in) :: xy_Snow (0:imax-1, 1:jmax) ! Local variables ! real(DP) :: xyz_DelMass(0:imax-1, 1:jmax, 1:kmax) real(DP) :: xy_Val(0:imax-1, 1:jmax) real(DP) :: xy_SumB(0:imax-1, 1:jmax) real(DP) :: xy_Sum(0:imax-1, 1:jmax) real(DP) :: xy_Ratio(0:imax-1, 1:jmax) integer :: i integer :: j integer :: k do k = 1, kmax xyz_DelMass(:,:,k) = ( xyr_Press(:,:,k-1) - xyr_Press(:,:,k) ) / Grav end do xy_Sum = 0.0_DP do k = kmax, 1, -1 xy_Val = CpDry * xyz_TempB(:,:,k) + LatentHeat * xyz_QH2OVapB(:,:,k) - LatentHeatFusion * xyz_QH2OSolB(:,:,k) xy_Sum = xy_Sum + xy_Val * xyz_DelMass(:,:,k) end do xy_SumB = xy_Sum xy_Sum = 0.0_DP do k = kmax, 1, -1 xy_Val = CpDry * xyz_Temp (:,:,k) + LatentHeat * xyz_QH2OVap (:,:,k) - LatentHeatFusion * xyz_QH2OSol (:,:,k) xy_Sum = xy_Sum + xy_Val * xyz_DelMass(:,:,k) end do !!$ if ( FlagIncludeSnowPhaseChange ) then xy_Sum = xy_Sum - LatentHeatFusion * xy_Snow * 2.0_DP * DelTime !!$ end if xy_Ratio = ( xy_Sum - xy_SumB ) / ( xy_Sum + 1.0d-100 ) do j = 1, jmax do i = 0, imax-1 if ( abs( xy_Ratio(i,j) ) > 1.0d-10 ) then call MessageNotify( 'M', module_name, '%c, Modified condensate static energy is not conserved, %f.', c1 = trim(ParentRoutineName), d = (/ xy_Ratio(i,j) /) ) end if end do end do xy_Sum = 0.0_DP do k = kmax, 1, -1 xy_Val = xyz_QH2OVapB(:,:,k) + xyz_QH2OLiqB(:,:,k) + xyz_QH2OSolB(:,:,k) xy_Sum = xy_Sum + xy_Val * xyz_DelMass(:,:,k) end do xy_SumB = xy_Sum xy_Sum = 0.0_DP do k = kmax, 1, -1 xy_Val = xyz_QH2OVap (:,:,k) + xyz_QH2OLiq (:,:,k) + xyz_QH2OSol (:,:,k) xy_Sum = xy_Sum + xy_Val * xyz_DelMass(:,:,k) end do xy_Sum = xy_Sum + ( xy_Rain + xy_Snow ) * 2.0_DP * DelTime xy_Ratio = ( xy_Sum - xy_SumB ) / ( xy_Sum + 1.0d-100 ) do j = 1, jmax do i = 0, imax-1 if ( abs( xy_Ratio(i,j) ) > 1.0d-10 ) then call MessageNotify( 'M', module_name, '%c, H2O mass is not conserved, %f.', c1 = trim(ParentRoutineName), d = (/ xy_Ratio(i,j) /) ) end if end do end do end subroutine CloudUtilConsChk
Subroutine : | |
xyz_TransCloudOneLayer(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in ) |
xyz_CloudCover(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in ) |
xyrr_OverlappedCloudTrans(0:imax-1, 1:jmax, 0:kmax, 0:kmax) : | real(DP), intent(out) |
subroutine CloudUtilsCalcOverlapCloudTrans( xyz_TransCloudOneLayer, xyz_CloudCover, xyrr_OverlappedCloudTrans ) ! USE statements ! ! ���¹ã�������¼ã�¿å�ºå�� ! History data output ! use gtool_historyauto, only: HistoryAutoPut ! ���»ç��� ! Time control ! use timeset, only: TimeN, EndTime, TimesetClockStart, TimesetClockStop !!$ use sort, only : SortQuick real(DP), intent(in ) :: xyz_TransCloudOneLayer (0:imax-1, 1:jmax, 1:kmax) real(DP), intent(in ) :: xyz_CloudCover (0:imax-1, 1:jmax, 1:kmax) real(DP), intent(out) :: xyrr_OverlappedCloudTrans(0:imax-1, 1:jmax, 0:kmax, 0:kmax) real(DP) :: xyz_EffCloudCover (0:imax-1, 1:jmax, 1:kmax) real(DP) :: xyz_CloudCoverSorted (0:imax-1, 1:jmax, 1:kmax) real(DP) :: xyz_EffCloudCoverSorted (0:imax-1, 1:jmax, 1:kmax) real(DP) :: xyz_TransCloudOneLayerSorted(0:imax-1, 1:jmax, 1:kmax) real(DP) :: CloudCoverSortedCur real(DP) :: EffCloudCoverSortedCur real(DP) :: TransCloudOneLayerSortedCur integer :: KInsPos integer :: i integer :: j integer :: k integer :: kk integer :: kkk ! å®�è¡��� ; Executable statement ! ! ������確è� ! Initialization check ! if ( .not. cloud_utils_inited ) then call MessageNotify( 'E', module_name, 'This module has not been initialized.' ) end if ! Cloud optical depth ! select case ( IDCloudOverlapType ) case ( IDCloudOverlapTypeRandom ) xyz_EffCloudCover = xyz_CloudCover * ( 1.0_DP - xyz_TransCloudOneLayer ) do k = 0, kmax kk = k xyrr_OverlappedCloudTrans(:,:,k,kk) = 1.0_DP do kk = k+1, kmax xyrr_OverlappedCloudTrans(:,:,k,kk) = xyrr_OverlappedCloudTrans(:,:,k,kk-1) * ( 1.0_DP - xyz_EffCloudCover(:,:,kk) ) end do end do do k = 0, kmax do kk = 0, k-1 xyrr_OverlappedCloudTrans(:,:,k,kk) = xyrr_OverlappedCloudTrans(:,:,kk,k) end do end do case ( IDCloudOverlapTypeMaxOverlap ) ! see Chou et al. (2001) xyz_EffCloudCover = xyz_CloudCover * ( 1.0_DP - xyz_TransCloudOneLayer ) ! Original method (computationally expensive, probably) ! !!$ do k = 0, kmax !!$ kk = k !!$ xyrr_OverlappedCloudTrans(:,:,k,kk) = 1.0_DP !!$ do kk = k+1, kmax !!$ !!$ xyz_CloudCoverSorted = xyz_CloudCover !!$ xyz_EffCloudCoverSorted = xyz_EffCloudCover !!$ xyz_TransCloudOneLayerSorted = xyz_TransCloudOneLayer !!$ !!$ call SortQuick( imax, jmax, kk-k, & !!$ & xyz_CloudCoverSorted (:,:,k+1:kk), & !!$ & xyz_EffCloudCoverSorted (:,:,k+1:kk), & !!$ & xyz_TransCloudOneLayerSorted(:,:,k+1:kk) & !!$ & ) !!$ !!$ xyrr_OverlappedCloudTrans(:,:,k,kk) = 0.0_DP !!$ do kkk = k+1, kk !!$ xyrr_OverlappedCloudTrans(:,:,k,kk) = & !!$ & xyz_EffCloudCoverSorted(:,:,kkk) & !!$ & + xyrr_OverlappedCloudTrans(:,:,k,kk) & !!$ & * xyz_TransCloudOneLayerSorted(:,:,kkk) !!$ end do !!$ xyrr_OverlappedCloudTrans(:,:,k,kk) = & !!$ & 1.0_DP - xyrr_OverlappedCloudTrans(:,:,k,kk) !!$ !!$ end do !!$ end do ! Economical method (probably) ! do k = 0, kmax !!$ do kkk = 1, kmax !!$ xyz_CloudCoverSorted(:,:,kkk) = real( kmax-kkk ) / real(kmax) !!$! xyz_CloudCoverSorted(:,:,kkk) = abs( 0.55d0 - real( kmax-kkk ) / real(kmax) ) !!$ end do !!$ ! debug output !!$ if ( k == 0 ) then !!$ kk = kmax !!$ do kkk = k+1, kk !!$ write( 6, * ) kkk, xyz_CloudCoverSorted(0,jmax/2+1,kkk) !!$ end do !!$ end if xyz_CloudCoverSorted = xyz_CloudCover xyz_EffCloudCoverSorted = xyz_EffCloudCover xyz_TransCloudOneLayerSorted = xyz_TransCloudOneLayer kk = k xyrr_OverlappedCloudTrans(:,:,k,kk) = 1.0_DP do kk = k+1, kmax do j = 1, jmax do i = 0, imax-1 ! xyz_CloudCoverSorted(i,j,kk) is inserved in an appropriate position. ! KInsPos = kk loop : do kkk = k+1, kk-1 if ( xyz_CloudCoverSorted(i,j,kk) < xyz_CloudCoverSorted(i,j,kkk) ) then KInsPos = kkk exit loop end if end do loop ! values are saved CloudCoverSortedCur = xyz_CloudCoverSorted (i,j,kk) EffCloudCoverSortedCur = xyz_EffCloudCoverSorted (i,j,kk) TransCloudOneLayerSortedCur = xyz_TransCloudOneLayerSorted(i,j,kk) ! values are shifted upward to empty an array at insert position do kkk = kk, KInsPos+1, -1 xyz_CloudCoverSorted (i,j,kkk) = xyz_CloudCoverSorted (i,j,kkk-1) xyz_EffCloudCoverSorted (i,j,kkk) = xyz_EffCloudCoverSorted (i,j,kkk-1) xyz_TransCloudOneLayerSorted(i,j,kkk) = xyz_TransCloudOneLayerSorted(i,j,kkk-1) end do kkk = KInsPos xyz_CloudCoverSorted (i,j,kkk) = CloudCoverSortedCur xyz_EffCloudCoverSorted (i,j,kkk) = EffCloudCoverSortedCur xyz_TransCloudOneLayerSorted(i,j,kkk) = TransCloudOneLayerSortedCur end do end do !!$ xyz_CloudCoverSorted = xyz_CloudCover !!$ do kkk = 1, kmax !!$ xyz_CloudCoverSorted(:,:,kkk) = real( kmax-kkk ) / real(kmax) !!$ end do !!$ xyz_EffCloudCoverSorted = xyz_EffCloudCover !!$ xyz_TransCloudOneLayerSorted = xyz_TransCloudOneLayer !!$ !!$ call SortQuick( imax, jmax, kk-k, & !!$ & xyz_CloudCoverSorted (:,:,k+1:kk), & !!$ & xyz_EffCloudCoverSorted (:,:,k+1:kk), & !!$ & xyz_TransCloudOneLayerSorted(:,:,k+1:kk) & !!$ & ) !!$ ! debug output !!$ if ( ( k == 0 ) .and. ( kk == kmax-2 ) ) then !!$ do kkk = k+1, kk !!$ write( 6, * ) kkk, xyz_CloudCoverSorted(0,jmax/2+1,kkk) !!$ end do !!$ write( 6, * ) '-----' !!$ end if xyrr_OverlappedCloudTrans(:,:,k,kk) = 0.0_DP do kkk = k+1, kk xyrr_OverlappedCloudTrans(:,:,k,kk) = xyz_EffCloudCoverSorted(:,:,kkk) + xyrr_OverlappedCloudTrans(:,:,k,kk) * xyz_TransCloudOneLayerSorted(:,:,kkk) end do xyrr_OverlappedCloudTrans(:,:,k,kk) = 1.0_DP - xyrr_OverlappedCloudTrans(:,:,k,kk) end do end do do k = 0, kmax do kk = 0, k-1 xyrr_OverlappedCloudTrans(:,:,k,kk) = xyrr_OverlappedCloudTrans(:,:,kk,k) end do end do end select ! Output effective cloud cover ! !!$ call HistoryAutoPut( TimeN, 'EffCloudCover', & !!$ & 1.0_DP - xyrr_OverlappedCloudTrans(:,:,0,kmax) ) end subroutine CloudUtilsCalcOverlapCloudTrans
Subroutine : | |
ArgFlagSnow : | logical, intent(in) |
This procedure input/output NAMELIST#cloud_utils_nml .
subroutine CloudUtilsInit( ArgFlagSnow ) ! ���¡ã�¤ã���¥å�ºå��è£��� ! File I/O support ! use dc_iounit, only: FileOpen ! NAMELIST ���¡ã�¤ã���¥å�����¢ã�������¼ã���£ã������ ! Utilities for NAMELIST file input ! use namelist_util, only: namelist_filename, NmlutilMsg, NmlutilAryValid ! ���¹ã�������¼ã�¿å�ºå�� ! History data output ! use gtool_historyauto, only: HistoryAutoAddVariable ! 飽å��æ¯�湿ã����� ! Evaluate saturation specific humidity ! use saturate, only: SaturateInit ! 宣è��� ; Declaration statements ! logical, intent(in) :: ArgFlagSnow character(STRING) :: CloudOverlapType integer:: unit_nml ! NAMELIST ���¡ã�¤ã�����¼ã���³ç���ç½����. ! Unit number for NAMELIST file open integer:: iostat_nml ! NAMELIST èªã�¿è¾¼�¿æ���� IOSTAT. ! IOSTAT of NAMELIST read ! NAMELIST å¤��°ç¾¤ ! NAMELIST group name ! namelist /cloud_utils_nml/ CloudOverlapType, CCNMixRatPerUnitMass ! ! �����������¤ã���¤ã��������������ç¶� "cloud_utils#CloudUtilsInit" ! ���½ã�¼ã�¹ã�³ã�¼ã�������§ã������. ! ! Refer to source codes in the initialization procedure ! "cloud_utils#CloudUtilsInit" for the default values. ! ! å®�è¡��� ; Executable statement ! if ( cloud_utils_inited ) return FlagSnow = ArgFlagSnow ! �����������¤ã��¨å®� ! Default values settings ! CloudOverlapType = "Random" !!$ CloudOverlapType = "MaxOverlap" CCNMixRatPerUnitMass = 1.0e8_DP ! NAMELIST ����¿è¾¼�� ! NAMELIST is input ! if ( trim(namelist_filename) /= '' ) then call FileOpen( unit_nml, namelist_filename, mode = 'r' ) ! (in) rewind( unit_nml ) read( unit_nml, nml = cloud_utils_nml, iostat = iostat_nml ) ! (out) close( unit_nml ) call NmlutilMsg( iostat_nml, module_name ) ! (in) end if select case ( CloudOverlapType ) case ( 'Random' ) IDCloudOverlapType = IDCloudOverlapTypeRandom case ( 'MaxOverlap' ) IDCloudOverlapType = IDCloudOverlapTypeMaxOverlap case default call MessageNotify( 'E', module_name, 'CloudOverlapType=<%c> is not supported.', c1 = trim(CloudOverlapType) ) end select ! Initialization of modules used in this module ! ! 飽å��æ¯�湿ã����� ! Evaluate saturation specific humidity ! call SaturateInit ! ���¹ã�������¼ã�¿å�ºå�����������¸ã����°ç�»é�� ! Register of variables for history data output ! !!$ call HistoryAutoAddVariable( 'EffCloudCover', & !!$ & (/ 'lon ', 'lat ', 'time' /), & !!$ & 'effective cloud cover', '1' ) ! �°å� ; Print ! call MessageNotify( 'M', module_name, '----- Initialization Messages -----' ) call MessageNotify( 'M', module_name, 'CloudOverlapType = %c', c1 = trim(CloudOverlapType) ) call MessageNotify( 'M', module_name, 'CCNMixRatPerUnitMass = %f', d = (/ CCNMixRatPerUnitMass /) ) call MessageNotify( 'M', module_name, '-- version = %c', c1 = trim(version) ) cloud_utils_inited = .true. end subroutine CloudUtilsInit
Subroutine : | |
xyz_CloudCover(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in ) |
xyz_DelCloudOptDep(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(inout) |
subroutine CloudUtilsLocalizeCloud( xyz_CloudCover, xyz_DelCloudOptDep ) ! USE statements ! real(DP), intent(in ) :: xyz_CloudCover (0:imax-1, 1:jmax, 1:kmax) real(DP), intent(inout) :: xyz_DelCloudOptDep(0:imax-1, 1:jmax, 1:kmax) ! å®�è¡��� ; Executable statement ! ! ������確è� ! Initialization check ! if ( .not. cloud_utils_inited ) then call MessageNotify( 'E', module_name, 'This module has not been initialized.' ) end if ! Cloud optical depth is scaled by considering cloud cover less than 1. xyz_DelCloudOptDep = xyz_DelCloudOptDep / max( xyz_CloudCover, 1.0d-3 ) end subroutine CloudUtilsLocalizeCloud
Subroutine : | |
Press : | real(DP), intent(in ) |
PressLI : | real(DP), intent(in ) |
PressUI : | real(DP), intent(in ) |
PRCPArea : | real(DP), intent(in ) |
PRCPEvapArea : | real(DP), intent(in ) |
Temp : | real(DP), intent(inout) |
QH2OVap : | real(DP), intent(inout) |
SurfRainFlux : | real(DP), intent(inout) |
SurfSnowFlux : | real(DP), intent(inout) |
subroutine CloudUtilsPRCPEvap1Grid( Press, PressLI, PressUI, PRCPArea, PRCPEvapArea, Temp, QH2OVap, SurfRainFlux, SurfSnowFlux ) ! USE statements ! ! ���»ç��� ! Time control ! use timeset, only: DelTime ! $ \Delta t $ [s] ! ����å®��°è¨å®� ! Physical constants settings ! use constants, only: Grav ! $ g $ [m s-2]. ! ��������åº�. ! Gravitational acceleration real(DP), intent(in ) :: Press real(DP), intent(in ) :: PressLI real(DP), intent(in ) :: PressUI real(DP), intent(in ) :: PRCPArea real(DP), intent(in ) :: PRCPEvapArea real(DP), intent(inout) :: Temp real(DP), intent(inout) :: QH2OVap real(DP), intent(inout) :: SurfRainFlux real(DP), intent(inout) :: SurfSnowFlux ! Local variables ! real(DP) :: DelTemp real(DP) :: DelQH2OVap real(DP) :: DelMass real(DP) :: QH2OVapSat real(DP) :: PRCPFlux real(DP) :: DelPRCPFlux !!$ real(DP) :: DelQH2OVap character(STRING) :: CharPhase integer :: l ! å®�è¡��� ; Executable statement ! ! ������確è� ! Initialization check ! if ( .not. cloud_utils_inited ) then call MessageNotify( 'E', module_name, 'This module has not been initialized.' ) end if DelMass = ( PressLI - PressUI ) / Grav do l = 1, 2 select case ( l ) case ( 1 ) ! liquid CharPhase = 'liquid' PRCPFlux = SurfRainFlux case ( 2 ) ! solid CharPhase = 'solid' PRCPFlux = SurfSnowFlux case default call MessageNotify( 'E', module_name, 'Unexpected number for select case.' ) end select call CloudUtilsPRCPEvap1GridCore( ( 2.0_DP * DelTime ), CharPhase, DelMass, Press, Temp, QH2OVap, PRCPFlux, PRCPArea, PRCPEvapArea, DelPRCPFlux, DelTemp, DelQH2OVap ) select case ( l ) case ( 1 ) ! liquid SurfRainFlux = PRCPFlux + DelPRCPFlux case ( 2 ) ! solid SurfSnowFlux = PRCPFlux + DelPRCPFlux case default call MessageNotify( 'E', module_name, 'Unexpected number for select case.' ) end select QH2OVap = QH2OVap + DelQH2OVap Temp = Temp + DelTemp end do end subroutine CloudUtilsPRCPEvap1Grid
Subroutine : | |
PressLI : | real(DP), intent(in ) |
PressUI : | real(DP), intent(in ) |
Temp : | real(DP), intent(inout) |
SurfRainFlux : | real(DP), intent(inout) |
SurfSnowFlux : | real(DP), intent(inout) |
subroutine CloudUtilsPRCPStepPC1Grid( PressLI, PressUI, Temp, SurfRainFlux, SurfSnowFlux ) ! ���»ç��� ! Time control ! use timeset, only: DelTime ! $ \Delta t $ [s] ! ����å®��°è¨å®� ! Physical constants settings ! use constants, only: CpDry, Grav, LatentHeatFusion ! $ L $ [J kg-1] . ! ��解ã�����. ! Latent heat of fusion ! ����µ·æ°·ã����°ã��¨å®� ! Setting constants of snow and sea ice ! use constants_snowseaice, only: TempCondWater real(DP), intent(in ) :: PressLI real(DP), intent(in ) :: PressUI real(DP), intent(inout) :: Temp real(DP), intent(inout) :: SurfRainFlux real(DP), intent(inout) :: SurfSnowFlux ! ä½�æ¥å��� ! Work variables ! real(DP) :: DelMass real(DP) :: MassMaxFreezeRate real(DP) :: MassFreezeRate real(DP) :: MassMaxMeltRate real(DP) :: MassMeltRate ! ������確è� ! Initialization check ! if ( .not. cloud_utils_inited ) then call MessageNotify( 'E', module_name, 'This module has not been initialized.' ) end if DelMass = ( PressLI - PressUI ) / Grav ! Freezing and melting switching at temperature of TempCondWater MassMaxFreezeRate = CpDry * ( TempCondWater - Temp ) * DelMass / LatentHeatFusion / ( 2.0_DP * DelTime ) if ( MassMaxFreezeRate >= 0.0_DP ) then ! freezing if ( SurfRainFlux >= MassMaxFreezeRate ) then MassFreezeRate = MassMaxFreezeRate else MassFreezeRate = SurfRainFlux end if SurfRainFlux = SurfRainFlux - MassFreezeRate SurfSnowFlux = SurfSnowFlux + MassFreezeRate Temp = Temp + LatentHeatFusion * MassFreezeRate * 2.0_DP * DelTime / ( CpDry * DelMass ) else ! melting MassMaxMeltRate = - MassMaxFreezeRate if ( SurfSnowFlux >= MassMaxMeltRate ) then MassMeltRate = MassMaxMeltRate else MassMeltRate = SurfSnowFlux end if SurfRainFlux = SurfRainFlux + MassMeltRate SurfSnowFlux = SurfSnowFlux - MassMeltRate Temp = Temp - LatentHeatFusion * MassMeltRate * 2.0_DP * DelTime / ( CpDry * DelMass ) end if end subroutine CloudUtilsPRCPStepPC1Grid
Subroutine : | |
xyz_CloudCover(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in ) |
xyz_DelCloudOptDep(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(inout) |
subroutine CloudUtilsSmearCloudOptDep( xyz_CloudCover, xyz_DelCloudOptDep ) ! USE statements ! real(DP), intent(in ) :: xyz_CloudCover (0:imax-1, 1:jmax, 1:kmax) real(DP), intent(inout) :: xyz_DelCloudOptDep(0:imax-1, 1:jmax, 1:kmax) ! å®�è¡��� ; Executable statement ! ! ������確è� ! Initialization check ! if ( .not. cloud_utils_inited ) then call MessageNotify( 'E', module_name, 'This module has not been initialized.' ) end if ! Cloud optical depth is scaled by the way of Kiehl et al. (1994). xyz_DelCloudOptDep = xyz_DelCloudOptDep * xyz_CloudCover**1.5_DP end subroutine CloudUtilsSmearCloudOptDep
Variable : | |||
CCNMixRatPerUnitMass : | real(DP), save
|
Subroutine : | |
TimeStep : | real(DP) , intent(in ) |
CharPhase : | character(*), intent(in ) |
DelMass : | real(DP) , intent(in ) |
Press : | real(DP) , intent(in ) |
Temp : | real(DP) , intent(in ) |
QH2OVap : | real(DP) , intent(in ) |
PRCPFlux : | real(DP) , intent(in ) |
PRCPArea : | real(DP) , intent(in ) |
PRCPEvapArea : | real(DP) , intent(in ) |
DelPRCPFlux : | real(DP) , intent(out) |
DelTemp : | real(DP) , intent(out) |
DelQH2OVap : | real(DP) , intent(out) |
subroutine CloudUtilsPRCPEvap1GridCore( TimeStep, CharPhase, DelMass, Press, Temp, QH2OVap, PRCPFlux, PRCPArea, PRCPEvapArea, DelPRCPFlux, DelTemp, DelQH2OVap ) ! �����»æ�°å¦å®��°è¨å®� ! Physical and mathematical constants settings ! use constants0, only: PI ! $ \pi $ . ! �����. Circular constant ! ����å®��°è¨å®� ! Physical constants settings ! use constants, only: Grav, CpDry, GasRDry, LatentHeat, LatentHeatFusion, EpsV ! $ \epsilon_v $ . ! æ°´è�¸æ���å��閴�. ! Molecular weight of water vapor ! 飽å��æ¯�湿ã����� ! Evaluate saturation specific humidity ! use saturate, only: CalcQVapSatOnLiq, CalcDQVapSatDTempOnLiq, CalcQVapSatOnSol, CalcDQVapSatDTempOnSol real(DP) , intent(in ) :: TimeStep character(*), intent(in ) :: CharPhase real(DP) , intent(in ) :: DelMass real(DP) , intent(in ) :: Press real(DP) , intent(in ) :: Temp real(DP) , intent(in ) :: QH2OVap real(DP) , intent(in ) :: PRCPFlux real(DP) , intent(in ) :: PRCPArea real(DP) , intent(in ) :: PRCPEvapArea real(DP) , intent(out) :: DelPRCPFlux real(DP) , intent(out) :: DelTemp real(DP) , intent(out) :: DelQH2OVap ! Parameters for evaporation of rain real(DP), parameter :: DensWater = 1.0d3 ! rho_w real(DP) :: CCNND ! number density of CCN, N0 or Nt (m-3) real(DP), parameter :: H2OVapDiffCoef = 1.0d-5 ! Kd real(DP), parameter :: PRCPFallVel0 = 10.0_DP ! m s-1 real(DP) :: PRCPFallVelRatio real(DP) :: PRCPFallVelFactor real(DP) :: PRCPFallVel ! m s-1 real(DP) :: Dens ! rho real(DP) :: DensPRCP ! (rho q_r) real(DP) :: DelZ real(DP) :: FactorF real(DP) :: FactorG real(DP) :: FactorH real(DP) :: FactorI real(DP) :: LatentHeatSubl real(DP) :: LatentHeatLocal real(DP) :: VirTemp real(DP) :: QH2OVapSat real(DP) :: DQH2OVapSatDTemp real(DP) :: QH2OVapSatA real(DP) :: TempN real(DP) :: QH2OVapN real(DP) :: PRCPN real(DP) :: TempA real(DP) :: QH2OVapA real(DP) :: PRCPA real(DP) :: DelPRCP real(DP), parameter :: DelTempThreshold = 1.0e-2_DP integer, parameter :: ItrMax = 100 real(DP) :: a_DelTemp(ItrMax) integer :: itr LatentHeatSubl = LatentHeat + LatentHeatFusion select case ( CharPhase ) case ( 'liquid' ) ! for liquid water PRCPFallVelRatio = 1.0_DP LatentHeatLocal = LatentHeat case ( 'solid' ) ! for solid water (ice) PRCPFallVelRatio = 0.5_DP LatentHeatLocal = LatentHeatSubl case default call MessageNotify( 'E', module_name, 'Unexpected character for select case.' ) end select VirTemp = Temp * ( 1.0_DP + ((( 1.0_DP / EpsV ) - 1.0_DP ) * QH2OVap) ) Dens = Press / ( GasRDry * VirTemp ) DelZ = DelMass / Dens PRCPFallVel = PRCPFallVel0 * PRCPFallVelRatio DensPRCP = ( PRCPFlux / ( PRCPArea + 1.0e-10_DP ) ) / PRCPFallVel ! cloud condensation CCNND = CCNMixRatPerUnitMass * Dens FactorF = Dens * H2OVapDiffCoef * PRCPEvapArea * ( 48.0_DP * ( PI * CCNND )**2 / DensWater * DensPRCP )**(1.0_DP/3.0_DP) FactorG = DelZ * TimeStep * FactorF FactorH = FactorG / DelMass FactorI = LatentHeatLocal / CpDry * FactorH TempN = Temp QH2OVapN = QH2OVap PRCPN = PRCPFlux * TimeStep loop_evap : do itr = 1, ItrMax select case ( CharPhase ) case ( 'liquid' ) ! for liquid water QH2OVapSat = CalcQVapSatOnLiq( TempN, Press ) DQH2OVapSatDTemp = CalcDQVapSatDTempOnLiq( TempN, QH2OVapSat ) case ( 'solid' ) ! for solid water (ice) QH2OVapSat = CalcQVapSatOnSol( TempN, Press ) DQH2OVapSatDTemp = CalcDQVapSatDTempOnSol( TempN, QH2OVapSat ) case default call MessageNotify( 'E', module_name, 'Unexpected character for select case.' ) end select DelTemp = - FactorI * ( QH2OVapSat - QH2OVapN ) / ( 1.0_DP + FactorH + FactorI * DQH2OVapSatDTemp ) DelQH2OVap = - CpDry * DelTemp / LatentHeatLocal DelPRCP = - DelQH2OVap * DelMass if ( ( PRCPN + DelPRCP ) >= 0.0_DP ) then ! part of precipitation is evaporated ! nothing to do else ! all precipitation is evaporated DelPRCP = - PRCPN DelQH2OVap = - DelPRCP / DelMass DelTemp = - LatentHeatLocal * DelQH2OVap / CpDry end if if ( abs( DelTemp ) < DelTempThreshold ) exit loop_evap PRCPA = PRCPN + DelPRCP TempA = TempN + DelTemp QH2OVapA = QH2OVapN + DelQH2OVap PRCPN = PRCPA TempN = TempA QH2OVapN = QH2OVapA a_DelTemp(itr) = DelTemp end do loop_evap if ( itr > 100 ) then write( 6, * ) a_DelTemp call MessageNotify( 'E', module_name, 'Evaporation loop is not convergent, %d, %f.', i = (/ itr /), d = (/ DelTemp /) ) end if DelPRCPFlux = DelPRCP / TimeStep end subroutine CloudUtilsPRCPEvap1GridCore
Subroutine : | |
CharPhase : | character(*), intent(in ) |
DelMass : | real(DP) , intent(in ) |
Press : | real(DP) , intent(in ) |
QH2OVap : | real(DP) , intent(in ) |
QH2OVapSat : | real(DP) , intent(in ) |
VirTemp : | real(DP) , intent(in ) |
PRCP : | real(DP) , intent(in ) |
PRCPArea : | real(DP) , intent(in ) |
PRCPEvapArea : | real(DP) , intent(in ) |
DelPRCPFlux : | real(DP) , intent(out) |
subroutine CloudUtilsPRCPEvap1GridCoreExp( CharPhase, DelMass, Press, QH2OVap, QH2OVapSat, VirTemp, PRCP, PRCPArea, PRCPEvapArea, DelPRCPFlux ) ! �����»æ�°å¦å®��°è¨å®� ! Physical and mathematical constants settings ! use constants0, only: PI ! $ \pi $ . ! �����. Circular constant ! ����å®��°è¨å®� ! Physical constants settings ! use constants, only: Grav, GasRDry ! $ R $ [J kg-1 K-1]. ! ä¹¾ç�¥å¤§æ°����ä½�å®���. ! Gas constant of air character(*), intent(in ) :: CharPhase real(DP) , intent(in ) :: DelMass real(DP) , intent(in ) :: Press real(DP) , intent(in ) :: QH2OVap real(DP) , intent(in ) :: QH2OVapSat real(DP) , intent(in ) :: VirTemp real(DP) , intent(in ) :: PRCP real(DP) , intent(in ) :: PRCPArea real(DP) , intent(in ) :: PRCPEvapArea real(DP) , intent(out) :: DelPRCPFlux ! Parameters for evaporation of rain real(DP), parameter :: DensWater = 1.0d3 ! rho_w ! Values below are from Kessler (1969) real(DP), parameter :: PRCPFallVelFactor0 = 130.0d0 ! K real(DP), parameter :: MedianDiameterFactor = 3.67d0 ! C' real(DP), parameter :: PRCPDistFactor = 1.0d7 ! N0 real(DP), parameter :: PRCPEvapRatUnitDiamFactor = 2.24d3 ! C real(DP), parameter :: H2OVapDiffCoef = 1.0d-5 ! Kd real(DP), parameter :: PRCPFallVelSimple0 = 10.0d0 ! m s-1 real(DP) :: PRCPFallVelRatio real(DP) :: PRCPFallVelFactor real(DP) :: Dens0 ! rho_0 real(DP) :: V00 ! V_{00} real(DP) :: PRCPEvapFactor real(DP) :: Dens ! rho real(DP) :: DensPRCP ! (rho q_r) !!$ real(DP) :: RainArea !!$ ! a_p !!$ real(DP) :: RainEvapArea !!$ ! A = max( a_p - a, 0 ) !!$ real(DP) :: xy_CloudCover (0:imax-1, 1:jmax) !!$ ! a real(DP) :: PRCPEvapRate real(DP) :: DelZ select case ( CharPhase ) case ( 'liquid' ) ! for liquid water PRCPFallVelRatio = 1.0_DP case ( 'solid' ) ! for solid water (ice) PRCPFallVelRatio = 0.5_DP case ( 'mixture' ) ! for mixture, this is only for test PRCPFallVelRatio = 1.0_DP case default call MessageNotify( 'E', module_name, 'Unexpected character for select case.' ) end select ! ! Values for evaporation of rain Dens = Press / ( GasRDry * VirTemp ) DelZ = DelMass / Dens if ( .false. ) then ! ECMWF version PRCPFallVelFactor = PRCPFallVelFactor0 * PRCPFallVelRatio ! Parameters for evaporation of rain Dens0 = 1013.0d2 / ( GasRDry * 300.0_DP ) V00 = PRCPFallVelFactor * sqrt( MedianDiameterFactor ) / ( PI * DensWater * PRCPDistFactor )**(1.0d0/8.0d0) PRCPEvapFactor = PRCPEvapRatUnitDiamFactor * 1.429624558860304d0 * H2OVapDiffCoef * PRCPDistFactor**(7.0d0/20.0d0) / ( PI * DensWater )**(13.0d0/20.0d0) !!$ RainArea = RainArea !!$ xy_CloudCover = CloudCover !!$ xy_RainEvapArea = max( xy_RainArea - xy_CloudCover, 0.0_DP ) !!$ RainEvapArea = RainEvapArea DensPRCP = ( PRCP / ( PRCPArea + 1.0d-10 ) / ( V00 * sqrt( Dens0 / Dens ) ) )**(8.0d0/9.0d0) PRCPEvapRate = Dens * PRCPEvapArea * PRCPEvapFactor * max( QH2OVapSat - QH2OVap, 0.0_DP ) * DensPRCP**(13.0d0/20.0d0) else ! simple version V00 = PRCPFallVelSimple0 * PRCPFallVelRatio PRCPEvapFactor = ( 48.0_DP * ( PI * PRCPDistFactor )**2 / DensWater )**(1.0_DP/3.0_DP) * H2OVapDiffCoef DensPRCP = ( PRCP / ( PRCPArea + 1.0d-10 ) ) / V00 PRCPEvapRate = Dens * PRCPEvapArea * PRCPEvapFactor * max( QH2OVapSat - QH2OVap, 0.0_DP ) * DensPRCP**(1.0_DP/3.0_DP) end if ! PRCPEvapRate (kg m-3 s-1) ! DelZ (m) ! DelPRCPFlux (kg m-2 s-1) DelPRCPFlux = PRCPEvapRate * DelZ DelPRCPFlux = min( DelPRCPFlux, PRCP ) end subroutine CloudUtilsPRCPEvap1GridCoreExp
Variable : | |||
cloud_utils_inited = .false. : | logical, save
|
Constant : | |||
module_name = ‘cloud_utils‘ : | character(*), parameter
|
Constant : | |||
version = ’$Name: $’ // ’$Id: cloud_utils.f90,v 1.7 2015/02/11 11:55:19 yot Exp $’ : | character(*), parameter
|