CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
GFS_physics_post.F90
1
8! ###########################################################################################
10 use machine, only : kind_phys, kind_dbl_prec, kind_sngl_prec
11 implicit none
12 public gfs_physics_post_run
13contains
14
18 subroutine gfs_physics_post_run(nCol, nLev, ntoz, ntracp100, nprocess, nprocess_summed, &
19 dtidx, is_photochem, ldiag3d, ip_physics, ip_photochem, ip_prod_loss, ip_ozmix, &
20 ip_temp, ip_overhead_ozone, do3_dt_prd, do3_dt_ozmx, do3_dt_temp, do3_dt_ohoz, &
21 dtend, errmsg, errflg)
22
23 ! Inputs
24 integer, intent(in) :: &
25 ncol, & !< Horizontal dimension
26 nlev, & !< Number of vertical layers
27 ntoz, & !< Index for ozone mixing ratio
28 ntracp100, & !< Number of tracers plus 100
29 nprocess, & !< Number of processes that cause changes in state variables
30 nprocess_summed,& !< Number of causes in dtidx per tracer summed for total physics tendency
31 ip_physics, & !< Index for process in diagnostic tendency output
32 ip_photochem, & !< Index for process in diagnostic tendency output
33 ip_prod_loss, & !< Index for process in diagnostic tendency output
34 ip_ozmix, & !< Index for process in diagnostic tendency output
35 ip_temp, & !< Index for process in diagnostic tendency output
36 ip_overhead_ozone
37 integer, intent(in), dimension(:,:) :: &
38 dtidx
39 logical, intent(in) :: &
40 ldiag3d
41 logical, intent(in), dimension(:) :: &
42 is_photochem
43
44 ! Inputs (optional)
45 real(kind=kind_phys), intent(in), dimension(:,:), pointer, optional :: &
46 do3_dt_prd, &
47 do3_dt_ozmx, &
48 do3_dt_temp, &
49 do3_dt_ohoz
50
51 ! Outputs
52 real(kind=kind_phys), intent(inout), dimension(:,:,:), optional :: &
53 dtend
54 character(len=*), intent(out) :: &
55 errmsg
56 integer, intent(out) :: &
57 errflg
58
59 ! Locals
60 integer :: idtend, ichem, iphys, itrac
61 logical :: all_true(nprocess)
62
63 ! Initialize CCPP error handling variables
64 errmsg = ''
65 errflg = 0
66
67 if(.not.ldiag3d) then
68 return
69 endif
70
71 ! #######################################################################################
72 !
73 ! Ozone physics diagnostics
74 !
75 ! #######################################################################################
76 idtend = dtidx(100+ntoz,ip_prod_loss)
77 if (idtend >= 1 .and. associated(do3_dt_prd)) then
78 dtend(:,:,idtend) = dtend(:,:,idtend) + do3_dt_prd
79 endif
80 !
81 idtend = dtidx(100+ntoz,ip_ozmix)
82 if (idtend >= 1 .and. associated(do3_dt_ozmx)) then
83 dtend(:,:,idtend) = dtend(:,:,idtend) + do3_dt_ozmx
84 endif
85 !
86 idtend = dtidx(100+ntoz,ip_temp)
87 if (idtend >= 1 .and. associated(do3_dt_temp)) then
88 dtend(:,:,idtend) = dtend(:,:,idtend) + do3_dt_temp
89 endif
90 !
91 idtend = dtidx(100+ntoz,ip_overhead_ozone)
92 if (idtend >= 1 .and. associated(do3_dt_ohoz)) then
93 dtend(:,:,idtend) = dtend(:,:,idtend) + do3_dt_ohoz
94 endif
95
96 ! #######################################################################################
97 !
98 ! Total (photochemical) tendencies.
99 !
100 ! #######################################################################################
101 itrac = ntoz+100
102 ichem = dtidx(itrac, ip_photochem)
103 if(ichem >= 1) then
104 call sum_it(ichem, itrac, is_photochem)
105 endif
106
107 ! #######################################################################################
108 !
109 ! Total (physics) tendencies
110 !
111 ! #######################################################################################
112 all_true = .true.
113 do itrac = 2,ntracp100
114 iphys = dtidx(itrac,ip_physics)
115 if(iphys >= 1) then
116 call sum_it(iphys, itrac, all_true)
117 endif
118 enddo
119
120 contains
121
123 subroutine sum_it(isum,itrac,sum_me)
124 integer, intent(in) :: isum ! third index of dtend of summary process
125 integer, intent(in) :: itrac ! tracer or state variable being summed
126 logical, intent(in) :: sum_me(nprocess) ! false = skip this process
127 logical :: first
128 integer :: idtend, iprocess
129
130 first=.true.
131 do iprocess=1,nprocess
132 if(iprocess>nprocess_summed) then
133 exit ! Don't sum up the sums.
134 else if(.not.sum_me(iprocess)) then
135 cycle ! We were asked to skip this one.
136 endif
137 idtend = dtidx(itrac,iprocess)
138 if(idtend>=1) then
139 ! This tendency was calculated for this tracer, so
140 ! accumulate it into the total tendency.
141 if(first) then
142 dtend(:,:,isum) = dtend(:,:,idtend)
143 first=.false.
144 else
145 dtend(:,:,isum) = dtend(:,:,isum) + dtend(:,:,idtend)
146 endif
147 endif
148 enddo
149 if(first) then
150 ! No tendencies were calculated, so sum is 0:
151 dtend(:,:,isum) = 0
152 endif
153 end subroutine sum_it
154 end subroutine gfs_physics_post_run
155end module gfs_physics_post