Line data Source code
1 : ! Computes deep and shallow convective cloud fractions
2 : ! CCPP-ized: Haipeng Lin, February 2025
3 : module convective_cloud_cover
4 :
5 : use ccpp_kinds, only: kind_phys
6 :
7 : implicit none
8 : private
9 : save
10 :
11 : ! public CCPP-compliant subroutines
12 : public :: convective_cloud_cover_init
13 : public :: convective_cloud_cover_run
14 :
15 : ! Module variables for tuning parameters
16 : real(kind_phys) :: sh1 ! Shallow convection tuning parameter 1 [fraction]
17 : real(kind_phys) :: sh2 ! Shallow convection tuning parameter 2 [m2 s kg-1]
18 : real(kind_phys) :: dp1 ! Deep convection tuning parameter 1 [fraction]
19 : real(kind_phys) :: dp2 ! Deep convection tuning parameter 2 [m2 s kg-1]
20 :
21 : contains
22 :
23 : !> \section arg_table_convective_cloud_cover_init Argument Table
24 : !! \htmlinclude convective_cloud_cover_init.html
25 1024 : subroutine convective_cloud_cover_init( &
26 : amIRoot, iulog, &
27 1024 : sh1_in, sh2_in, dp1_in, dp2_in, errmsg, errflg)
28 :
29 : logical, intent(in) :: amIRoot
30 : integer, intent(in) :: iulog ! log output unit
31 : real(kind_phys), intent(in) :: sh1_in ! Shallow convection parameter 1
32 : real(kind_phys), intent(in) :: sh2_in ! Shallow convection parameter 2
33 : real(kind_phys), intent(in) :: dp1_in ! Deep convection parameter 1
34 : real(kind_phys), intent(in) :: dp2_in ! Deep convection parameter 2
35 : character(len=512), intent(out) :: errmsg
36 : integer, intent(out) :: errflg
37 :
38 1024 : errmsg = ''
39 1024 : errflg = 0
40 :
41 : ! Set module variables from namelist
42 1024 : sh1 = sh1_in
43 1024 : sh2 = sh2_in
44 1024 : dp1 = dp1_in
45 1024 : dp2 = dp2_in
46 :
47 1024 : if(amIRoot) then
48 2 : write(iulog,*) 'tuning parameters convective_cloud_cover: dp1 ',dp1,' dp2 ',dp2,' sh1 ',sh1,' sh2 ',sh2
49 : endif
50 :
51 1024 : end subroutine convective_cloud_cover_init
52 :
53 : ! Compute convective cloud cover (deep and shallow)
54 : ! Should produce typical numbers of 20%
55 : ! Shallow and deep convective cloudiness are evaluated separately and summed
56 : ! because the processes are evaluated separately.
57 : !> \section arg_table_convective_cloud_cover_run Argument Table
58 : !! \htmlinclude convective_cloud_cover_run.html
59 70392 : subroutine convective_cloud_cover_run( &
60 : ncol, pver, &
61 : top_lev_cloudphys, &
62 70392 : use_shfrc, shfrc, &
63 140784 : cmfmc_total, cmfmc_sh, &
64 140784 : shallowcu, deepcu, concld, &
65 70392 : errmsg, errflg)
66 :
67 : ! Input arguments
68 : integer, intent(in) :: ncol
69 : integer, intent(in) :: pver
70 : integer, intent(in) :: top_lev_cloudphys ! Top vertical level for cloud physics [index]
71 :
72 : logical, intent(in) :: use_shfrc ! Use cloud area fraction provided by shallow convection? [flag]
73 : real(kind_phys), intent(in) :: shfrc(:, :) ! Input shallow cloud fraction [fraction]
74 :
75 : real(kind_phys), intent(in) :: cmfmc_total(:, :) ! atmosphere_convective_mass_flux_due_to_all_convection [kg m-2 s-1]
76 : real(kind_phys), intent(in) :: cmfmc_sh(:, :) ! atmosphere_convective_mass_flux_due_to_shallow_convection [kg m-2 s-1]
77 :
78 : ! Output arguments
79 : real(kind_phys), intent(out) :: shallowcu(:, :) ! Shallow convective cloud fraction [fraction]
80 : real(kind_phys), intent(out) :: deepcu(:, :) ! Deep convective cloud fraction [fraction]
81 : real(kind_phys), intent(out) :: concld(:, :) ! Convective cloud cover [fraction]
82 : character(len=512), intent(out) :: errmsg
83 : integer, intent(out) :: errflg
84 :
85 : ! Local variables
86 : integer :: i, k
87 :
88 70392 : errmsg = ''
89 70392 : errflg = 0
90 :
91 28436184 : concld(:ncol, :) = 0.0_kind_phys
92 1900584 : do k = top_lev_cloudphys, pver
93 28436184 : do i = 1, ncol
94 26535600 : if (.not. use_shfrc) then
95 : ! Compute the shallow convection cloud cover using the shallow convective mass flux.
96 53071200 : shallowcu(i, k) = max(0.0_kind_phys, &
97 79606800 : min(sh1*log(1.0_kind_phys + sh2*cmfmc_sh(i, k + 1)), 0.30_kind_phys))
98 : else
99 : ! If use shallow convective calculated clouds, then just assign
100 0 : shallowcu(i, k) = shfrc(i, k)
101 : end if
102 :
103 : ! REMOVECAM: This could be changed to use deep convective mass flux
104 : ! since it is independently available in CCPP, once CAM is retired
105 53071200 : deepcu(i, k) = max(0.0_kind_phys, &
106 : min(dp1 * &
107 : log(1.0_kind_phys + &
108 106142400 : dp2 * (cmfmc_total(i, k + 1) - cmfmc_sh(i, k + 1)) &
109 : ), &
110 159213600 : 0.60_kind_phys))
111 :
112 : ! Estimate of local convective cloud cover based on convective mass flux
113 : ! Modify local large-scale relative humidity to account for presence of
114 : ! convective cloud when evaluating relative humidity based layered cloud amount
115 28365792 : concld(i, k) = min(shallowcu(i, k) + deepcu(i, k), 0.80_kind_phys)
116 : end do
117 : end do
118 :
119 70392 : end subroutine convective_cloud_cover_run
120 :
121 : end module convective_cloud_cover
|