CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
cu_c3_driver_post.F90
1
3
5
6 implicit none
7
8 private
9
10 public :: cu_c3_driver_post_run
11
12 contains
13
18 subroutine cu_c3_driver_post_run (im, km, t, q, prevst, prevsq, cactiv, cactiv_m, conv_act, conv_act_m, dt, garea, raincv, maxupmf, refl_10cm, errmsg, errflg)
19
20 use machine, only: kind_phys
21
22 implicit none
23
24 ! Interface variables
25 integer, intent(in) :: im, km
26 real(kind_phys), intent(in) :: t(:,:)
27 real(kind_phys), intent(in) :: q(:,:)
28 real(kind_phys), dimension(:),intent(in) :: garea
29 real(kind_phys), intent(out), optional :: prevst(:,:)
30 real(kind_phys), intent(out), optional :: prevsq(:,:)
31 integer, intent(in), optional :: cactiv(:)
32 integer, intent(in), optional :: cactiv_m(:)
33 real(kind_phys), intent(out), optional :: conv_act(:)
34 real(kind_phys), intent(out), optional :: conv_act_m(:)
35 ! for Radar reflectivity
36 real(kind_phys), intent(in) :: dt
37 real(kind_phys), intent(in) :: raincv(:)
38 real(kind_phys), intent(in), optional :: maxupmf(:)
39 real(kind_phys), intent(inout) :: refl_10cm(:,:)
40 character(len=*), intent(out) :: errmsg
41!$acc declare copyin(t,q,cactiv,cactiv_m) copyout(prevst,prevsq,conv_act,conv_act_m)
42 integer, intent(out) :: errflg
43
44 ! Local variables
45 real(kind_phys), parameter :: dbzmin=-10.0
46 real(kind_phys) :: cuprate
47 real(kind_phys) :: ze, ze_conv, dbz_sum
48 integer :: i, k
49
50 ! Initialize CCPP error handling variables
51 errmsg = ''
52 errflg = 0
53
54!$acc kernels
55 prevst(:,:) = t(:,:)
56 prevsq(:,:) = q(:,:)
57
58 do i = 1, im
59 if (cactiv(i).gt.0) then
60 conv_act(i) = conv_act(i)+1.0
61 else
62 conv_act(i)=0.0
63 endif
64 if (cactiv_m(i).gt.0) then
65 conv_act_m(i) = conv_act_m(i)+1.0
66 else
67 conv_act_m(i)=0.0
68 endif
69 ! reflectivity parameterization for parameterized convection (reference:Unipost MDLFLD.f)
70 ze = 0.0
71 ze_conv = 0.0
72 dbz_sum = 0.0
73 cuprate = 1.e3*raincv(i) * 3600.0 / dt ! cu precip rate (mm/h)
74 if(cuprate .lt. 0.05) cuprate=0.
75 ze_conv = 300.0 * cuprate**1.5
76 if (maxupmf(i).gt.0.1 .and. cuprate.gt.0.) then
77 do k = 1, km
78 ze = 10._kind_phys ** (0.1 * refl_10cm(i,k))
79 dbz_sum = max(dbzmin, 10.0 * log10(ze + ze_conv))
80 refl_10cm(i,k) = dbz_sum
81 enddo
82 endif
83 enddo
84!$acc end kernels
85
86 end subroutine cu_c3_driver_post_run
87
88end module cu_c3_driver_post
subroutine, public cu_c3_driver_post_run(im, km, t, q, prevst, prevsq, cactiv, cactiv_m, conv_act, conv_act_m, dt, garea, raincv, maxupmf, refl_10cm, errmsg, errflg)