CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
GFS_DCNV_generic_post.F90
1
3
5
6 contains
7
11 subroutine gfs_dcnv_generic_post_run (im, levs, lssav, ldiag3d, qdiag3d, ras, &
12 cscnv, frain, rain1, dtf, cld1d, save_u, save_v, save_t, gu0, gv0, gt0, &
13 ud_mf, dd_mf, dt_mf, con_g, npdf3d, num_p3d, ncnvcld3d, nsamftrac, &
14 rainc, cldwrk, upd_mf, dwn_mf, det_mf, dtend, dtidx, index_of_process_dcnv, &
15 index_of_temperature, index_of_x_wind, index_of_y_wind, ntqv, gq0, save_q, &
16 cnvw, cnvc, cnvw_phy_f3d, cnvc_phy_f3d, flag_for_dcnv_generic_tend, &
17 ntcw,ntiw,ntclamt,ntrw,ntsw,ntrnc,ntsnc,ntgl, &
18 ntgnc, nthl, nthnc, nthv, ntgv, ntrz, ntgz, nthz, ntsigma, ntrac,clw, &
19 satmedmf, trans_trac, errmsg, errflg)
20
21
22 use machine, only: kind_phys
23
24 implicit none
25
26 integer, intent(in) :: im, levs, nsamftrac
27 logical, intent(in) :: lssav, ldiag3d, qdiag3d, ras, cscnv
28 logical, intent(in) :: flag_for_dcnv_generic_tend
29
30 real(kind=kind_phys), intent(in) :: frain, dtf
31 real(kind=kind_phys), dimension(:), intent(in) :: rain1, cld1d
32 real(kind=kind_phys), dimension(:,:), intent(in) :: save_u, save_v, save_t
33 real(kind=kind_phys), dimension(:,:), intent(in) :: gu0, gv0, gt0
34 real(kind=kind_phys), dimension(:,:,:), intent(in) :: gq0, save_q
35 real(kind=kind_phys), dimension(:,:), intent(in) :: dd_mf, dt_mf
36 real(kind=kind_phys), dimension(:,:), intent(in), optional :: ud_mf
37 real(kind=kind_phys), intent(in) :: con_g
38 integer, intent(in) :: npdf3d, num_p3d, ncnvcld3d
39 logical, intent(in) :: satmedmf, trans_trac
40
41 real(kind=kind_phys), dimension(:), intent(inout) :: rainc, cldwrk
42 real(kind=kind_phys), dimension(:,:), intent(inout), optional :: upd_mf, dwn_mf, det_mf
43 real(kind=kind_phys), dimension(:,:), intent(inout) :: cnvw, cnvc
44
45 real(kind=kind_phys), dimension(:,:,:), intent(inout), optional :: dtend
46 integer, intent(in) :: dtidx(:,:), index_of_process_dcnv, index_of_temperature, &
47 index_of_x_wind, index_of_y_wind, ntqv
48 integer, intent(in) :: ntcw,ntiw,ntclamt,ntrw,ntsw,ntrnc,ntsnc,ntgl, &
49 ntgnc, nthl, nthnc, nthv, ntgv, ntrz, ntgz, nthz, &
50 ntsigma, ntrac
51 real(kind=kind_phys), dimension(:,:,:), intent(in) :: clw
52
53
54 real(kind=kind_phys), dimension(:,:), intent(inout), optional :: cnvw_phy_f3d, cnvc_phy_f3d
55
56 character(len=*), intent(out) :: errmsg
57 integer, intent(out) :: errflg
58
59 integer :: i, k, n, idtend, tracers
60
61 ! Initialize CCPP error handling variables
62 errmsg = ''
63 errflg = 0
64
65 if (.not. ras .and. .not. cscnv) then
66 if (npdf3d == 3 .and. num_p3d == 4) then
67 do k=1,levs
68 do i=1,im
69 cnvw_phy_f3d(i,k) = cnvw(i,k)
70 cnvc_phy_f3d(i,k) = cnvc(i,k)
71 cnvw(i,k) = 0.0
72 cnvc(i,k) = 0.0
73 enddo
74 enddo
75 elseif (npdf3d == 0 .and. ncnvcld3d == 1) then
76 do k=1,levs
77 do i=1,im
78 cnvw_phy_f3d(i,k) = cnvw(i,k)
79 cnvw(i,k) = 0.0
80 enddo
81 enddo
82 endif
83 endif ! if (.not. ras .and. .not. cscnv)
84
85 do i=1,im
86 rainc(i) = frain * rain1(i)
87 enddo
88!
89 if (lssav) then
90 do i=1,im
91 cldwrk(i) = cldwrk(i) + cld1d(i) * dtf
92 enddo
93
94 if (ldiag3d .and. flag_for_dcnv_generic_tend) then
95 idtend=dtidx(index_of_temperature,index_of_process_dcnv)
96 if(idtend>=1) then
97 dtend(:,:,idtend) = dtend(:,:,idtend) + (gt0-save_t)*frain
98 endif
99
100 idtend=dtidx(index_of_x_wind,index_of_process_dcnv)
101 if(idtend>=1) then
102 dtend(:,:,idtend) = dtend(:,:,idtend) + (gu0-save_u)*frain
103 endif
104
105 idtend=dtidx(index_of_y_wind,index_of_process_dcnv)
106 if(idtend>=1) then
107 dtend(:,:,idtend) = dtend(:,:,idtend) + (gv0-save_v)*frain
108 endif
109
110 if (cscnv .or. satmedmf .or. trans_trac .or. ras) then
111 tracers = 2
112 do n=2,ntrac
113 if ( n /= ntcw .and. n /= ntiw .and. n /= ntclamt .and. &
114 n /= ntrw .and. n /= ntsw .and. n /= ntrnc .and. &
115 n /= ntsnc .and. n /= ntgl .and. n /= ntgnc .and. &
116 n /= nthl .and. n /= nthnc .and. n /= nthv .and. &
117 n /= ntrz .and. n /= ntgz .and. n /= nthz .and. &
118 n /= ntgv .and. n /= ntsigma) then
119 tracers = tracers + 1
120 idtend = dtidx(100+n,index_of_process_dcnv)
121 if(idtend>0) then
122 dtend(:,:,idtend) = dtend(:,:,idtend) + clw(:,:,tracers)-save_q(:,:,n) * frain
123 endif
124 endif
125 enddo
126 else
127 do n=2,ntrac
128 idtend = dtidx(100+n,index_of_process_dcnv)
129 if(idtend>0) then
130 dtend(:,:,idtend) = dtend(:,:,idtend) + (gq0(:,:,n)-save_q(:,:,n))*frain
131 endif
132 enddo
133 endif
134 idtend = dtidx(100+ntqv, index_of_process_dcnv)
135 if(idtend>=1) then
136 dtend(:,:,idtend) = dtend(:,:,idtend) + (gq0(:,:,ntqv) - save_q(:,:,ntqv)) * frain
137 endif
138
139 ! convective mass fluxes
140 if(qdiag3d) then
141 do k=1,levs
142 do i=1,im
143 upd_mf(i,k) = upd_mf(i,k) + ud_mf(i,k) * (con_g*frain)
144 dwn_mf(i,k) = dwn_mf(i,k) + dd_mf(i,k) * (con_g*frain)
145 det_mf(i,k) = det_mf(i,k) + dt_mf(i,k) * (con_g*frain)
146 enddo
147 enddo
148 endif
149 endif ! if (ldiag3d)
150
151 endif ! if (lssav)
152
153 end subroutine gfs_dcnv_generic_post_run
154 end module gfs_dcnv_generic_post