CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
GFS_SCNV_generic_post.F90
1
3
5
6 contains
7
11 subroutine gfs_scnv_generic_post_run (im, levs, nn, lssav, ldiag3d, qdiag3d, &
12 frain, gu0, gv0, gt0, gq0, save_u, save_v, save_t, save_q, &
13 clw, shcnvcw, rain1, npdf3d, num_p3d, ncnvcld3d, cnvc, cnvw, nsamftrac, &
14 rainc, cnvprcp, cnvprcpb, cnvw_phy_f3d, cnvc_phy_f3d, &
15 dtend, dtidx, index_of_temperature, index_of_x_wind, index_of_y_wind, &
16 index_of_process_scnv, ntqv, flag_for_scnv_generic_tend, &
17 ntcw,ntiw,ntclamt,ntrw,ntsw,ntrnc,ntsnc,ntgl,ntgnc,ntsigma, &
18 imfshalcnv, imfshalcnv_sas, imfshalcnv_samf, ntrac, &
19 cscnv, satmedmf, trans_trac, ras, errmsg, errflg)
20
21 use machine, only: kind_phys
22
23 implicit none
24
25 integer, intent(in) :: im, levs, nn, ntqv, nsamftrac
26 integer, intent(in) :: ntcw,ntiw,ntclamt,ntrw,ntsw,ntrnc,ntsnc,ntgl,ntgnc,ntsigma,ntrac
27 logical, intent(in) :: lssav, ldiag3d, qdiag3d, flag_for_scnv_generic_tend
28 real(kind=kind_phys), intent(in) :: frain
29 real(kind=kind_phys), dimension(:,:), intent(in) :: gu0, gv0, gt0
30 real(kind=kind_phys), dimension(:,:), intent(in) :: save_u, save_v, save_t
31 real(kind=kind_phys), dimension(:,:,:), intent(in) :: save_q, gq0
32
33 ! dtend only allocated if ldiag3d == .true.
34 real(kind=kind_phys), intent(inout), optional :: dtend(:,:,:)
35 integer, intent(in) :: dtidx(:,:)
36 integer, intent(in) :: index_of_temperature, index_of_x_wind, index_of_y_wind, index_of_process_scnv
37 real(kind=kind_phys), dimension(:,:,:), intent(in) :: clw
38
39 ! Post code for SAS/SAMF
40 integer, intent(in) :: npdf3d, num_p3d, ncnvcld3d
41 logical, intent(in) :: shcnvcw
42 real(kind=kind_phys), dimension(:), intent(in) :: rain1
43 real(kind=kind_phys), dimension(:, :), intent(in) :: cnvw, cnvc
44 real(kind=kind_phys), dimension(:), intent(inout) :: rainc, cnvprcp, cnvprcpb
45 ! The following arrays may not be allocated, depending on certain flags and microphysics schemes.
46 ! Since Intel 15 crashes when passing unallocated arrays to arrays defined with explicit shape,
47 ! use assumed-shape arrays. Note that Intel 18 and GNU 6.2.0-8.1.0 tolerate explicit-shape arrays
48 ! as long as these do not get used when not allocated.
49 real(kind=kind_phys), dimension(:,:), intent(inout), optional :: cnvw_phy_f3d, cnvc_phy_f3d
50 integer, intent(in) :: imfshalcnv, imfshalcnv_sas, imfshalcnv_samf
51 logical, intent(in) :: cscnv, satmedmf, trans_trac, ras
52
53 character(len=*), intent(out) :: errmsg
54 integer, intent(out) :: errflg
55
56 integer :: i, k, n, idtend, tracers
57 real(kind=kind_phys) :: tem
58
59 ! Initialize CCPP error handling variables
60 errmsg = ''
61 errflg = 0
62
63 if (imfshalcnv==imfshalcnv_sas .or. imfshalcnv==imfshalcnv_samf) then
64 do i=1,im
65 rainc(i) = rainc(i) + frain * rain1(i)
66 enddo
67! 'cnvw' and 'cnvc' are set to zero before computation starts:
68 if (shcnvcw .and. num_p3d == 4 .and. npdf3d == 3) then
69 do k=1,levs
70 do i=1,im
71 cnvw_phy_f3d(i,k) = cnvw_phy_f3d(i,k) + cnvw(i,k)
72 cnvc_phy_f3d(i,k) = cnvc_phy_f3d(i,k) + cnvc(i,k)
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_phy_f3d(i,k) + cnvw(i,k)
79 enddo
80 enddo
81 endif
82 endif
83
84 if (lssav .and. flag_for_scnv_generic_tend) then
85 if (ldiag3d) then
86 idtend = dtidx(index_of_temperature, index_of_process_scnv)
87 if(idtend>=1) then
88 dtend(:,:,idtend) = dtend(:,:,idtend) + (gt0 - save_t) * frain
89 endif
90
91 idtend = dtidx(index_of_x_wind, index_of_process_scnv)
92 if(idtend>=1) then
93 dtend(:,:,idtend) = dtend(:,:,idtend) + (gu0 - save_u) * frain
94 endif
95
96 idtend = dtidx(index_of_y_wind, index_of_process_scnv)
97 if(idtend>=1) then
98 dtend(:,:,idtend) = dtend(:,:,idtend) + (gv0 - save_v) * frain
99 endif
100
101 if (cscnv .or. satmedmf .or. trans_trac .or. ras) then
102 tracers = 2
103 do n=2,ntrac
104 if ( n /= ntcw .and. n /= ntiw .and. n /= ntclamt .and. &
105 n /= ntrw .and. n /= ntsw .and. n /= ntrnc .and. &
106 n /= ntsnc .and. n /= ntgl .and. n /= ntgnc .and. n /= ntsigma) then
107 tracers = tracers + 1
108 idtend = dtidx(100+n,index_of_process_scnv)
109 if(idtend>0) then
110 dtend(:,:,idtend) = dtend(:,:,idtend) + clw(:,:,tracers)-save_q(:,:,n) * frain
111 endif
112 endif
113 enddo
114 else
115 do n=2,ntrac
116 idtend = dtidx(100+n,index_of_process_scnv)
117 if(idtend>0) then
118 dtend(:,:,idtend) = dtend(:,:,idtend) + (gq0(:,:,n)-save_q(:,:,n))*frain
119 endif
120 enddo
121 endif
122 idtend = dtidx(100+ntqv, index_of_process_scnv)
123 if(idtend>=1) then
124 dtend(:,:,idtend) = dtend(:,:,idtend) + (gq0(:,:,ntqv) - save_q(:,:,ntqv)) * frain
125 endif
126 endif
127 endif
128
129 end subroutine gfs_scnv_generic_post_run
130
131 end module gfs_scnv_generic_post