CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
cs_conv_aw_adj.F90
1
3
8
9 implicit none
10
11 private
12
13 public :: cs_conv_aw_adj_run
14
15 contains
16
22!\section gen_cs_conv_aw_adj_run CPT cs_conv_aw_adj_run General Algorithm
23 subroutine cs_conv_aw_adj_run(im, levs, do_cscnv, do_aw, do_shoc, &
24 ntrac, ntcw, ntclamt, nncl, con_g, sigmafrac, &
25 gt0, gq0, save_t, save_q, prsi, cldfrac, subcldfrac, &
26 prcp, imp_physics, imp_physics_mg, errmsg, errflg)
27
28 use machine, only: kind_phys
29
30 implicit none
31
32! --- interface variables
33 integer, intent(in) :: im, levs
34 logical, intent(in) :: do_cscnv, do_aw, do_shoc
35 integer, intent(in) :: ntrac, ntcw, ntclamt, nncl
36 real(kind_phys), intent(in) :: con_g
37 real(kind_phys), dimension(:,:), intent(inout) :: sigmafrac
38 real(kind_phys), dimension(:,:), intent(inout) :: gt0
39 real(kind_phys), dimension(:,:,:), intent(inout) :: gq0
40 real(kind_phys), dimension(:,:), intent(in) :: save_t
41 real(kind_phys), dimension(:,:,:), intent(in) :: save_q
42 real(kind_phys), dimension(:,:), intent(in) :: prsi
43 real(kind_phys), dimension(:,:), intent(inout), optional :: cldfrac
44 real(kind_phys), dimension(:,:), intent(inout), optional :: subcldfrac
45 real(kind_phys), dimension(:), intent(inout) :: prcp
46 integer, intent(in ) :: imp_physics, imp_physics_mg
47 character(len=*), intent( out) :: errmsg
48 integer, intent( out) :: errflg
49
50! --- local variables
51 real(kind=kind_phys), dimension(im) :: temrain1
52 real(kind=kind_phys) :: tem1, tem2
53 real(kind=kind_phys) :: onebg
54 integer :: i, k, n
55
56! --- initialize CCPP error handling variables
57 errmsg = ''
58 errflg = 0
59
60 !if (do_cscnv .and. do_aw) then
61
62 onebg = 1.0_kind_phys/con_g
63
64! Arakawa-Wu adjustment of large-scale microphysics tendencies:
65! reduce by factor of (1-sigma)
66! these are microphysics increments. We want to keep (1-sigma) of the increment,
67! we will remove sigma*increment from final values
68! fsigma = 0. ! don't apply any AW correction, in addition comment next line
69! fsigma = sigmafrac
70
71! adjust sfc rainrate for conservation
72! vertically integrate reduction of water increments, reduce precip by that amount
73
74 temrain1(:) = 0.0
75 do k = 1,levs
76 do i = 1,im
77 tem1 = sigmafrac(i,k)
78 gt0(i,k) = gt0(i,k) - tem1 * (gt0(i,k)-save_t(i,k))
79 tem2 = tem1 * (gq0(i,k,1)-save_q(i,k,1))
80 gq0(i,k,1) = gq0(i,k,1) - tem2
81 temrain1(i) = temrain1(i) - (prsi(i,k)-prsi(i,k+1)) * tem2 * onebg
82 enddo
83 enddo
84! add convective clouds if shoc is true and not MG microphysics
85 if (do_shoc .and. imp_physics /= imp_physics_mg) then
86 do k = 1,levs
87 do i = 1,im
88 subcldfrac(i,k) = min(1.0, subcldfrac(i,k) + sigmafrac(i,k))
89 enddo
90 enddo
91 endif
92 !
93 do n=ntcw,ntcw+nncl-1
94 do k = 1,levs
95 do i = 1,im
96 tem1 = sigmafrac(i,k) * (gq0(i,k,n)-save_q(i,k,n))
97 gq0(i,k,n) = gq0(i,k,n) - tem1
98 temrain1(i) = temrain1(i) - (prsi(i,k)-prsi(i,k+1)) * tem1 * onebg
99 enddo
100 enddo
101 enddo
102 !
103 do i = 1,im
104 prcp(i) = max(prcp(i) - temrain1(i)*0.001, 0.0_kind_phys)
105 enddo
106
107 !endif
108
109 return
110
111 end subroutine cs_conv_aw_adj_run
112
113
114end module cs_conv_aw_adj
subroutine, public cs_conv_aw_adj_run(im, levs, do_cscnv, do_aw, do_shoc, ntrac, ntcw, ntclamt, nncl, con_g, sigmafrac, gt0, gq0, save_t, save_q, prsi, cldfrac, subcldfrac, prcp, imp_physics, imp_physics_mg, errmsg, errflg)
This subroutine adjusts surface rainrate for conservation.