CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
ccpp_suite_simulator.F90
1
17!
18! ########################################################################################
20 use machine, only: kind_phys
21 use module_ccpp_suite_simulator, only: base_physics_process, sim_lwrad, sim_swrad, &
22 sim_pbl, sim_gwd, sim_dcnv, sim_scnv, sim_cldmp
23 implicit none
24 public ccpp_suite_simulator_run
25contains
26
27 ! ######################################################################################
28 !
29 ! SUBROUTINE ccpp_suite_simulator_run
30 !
31 ! ######################################################################################
32!! \section arg_table_ccpp_suite_simulator_run
33!! \htmlinclude ccpp_suite_simulator_run.html
34!!
35 subroutine ccpp_suite_simulator_run(do_ccpp_suite_sim, kdt, nCol, nLay, dtp, jdat, &
36 iactive_T, iactive_u, iactive_v, iactive_q, proc_start, proc_end, physics_process,&
37 in_pre_active, in_post_active, tgrs, ugrs, vgrs, qgrs, active_phys_tend, gt0, gu0,&
38 gv0, gq0, errmsg, errflg)
39
40 ! Inputs
41 logical, intent(in) :: do_ccpp_suite_sim
42 integer, intent(in) :: kdt, ncol, nlay, jdat(8), iactive_t, iactive_u, &
43 iactive_v, iactive_q
44 real(kind_phys), intent(in) :: dtp, tgrs(:,:), ugrs(:,:), vgrs(:,:), qgrs(:,:,:), &
45 active_phys_tend(:,:,:)
46 ! Outputs
47 type(base_physics_process),intent(inout) :: physics_process(:)
48 real(kind_phys), intent(inout) :: gt0(:,:), gu0(:,:), gv0(:,:), gq0(:,:)
49 character(len=*),intent(out) :: errmsg
50 integer, intent(out) :: errflg
51 integer, intent(inout) :: proc_start, proc_end
52 logical, intent(inout) :: in_pre_active, in_post_active
53
54 ! Locals
55 integer :: icol, year, month, day, hour, min, sec, iprc
56 real(kind_phys), dimension(nCol,nLay) :: gt1, gu1, gv1, dtdt, dudt, dvdt, gq1, dqdt
57
58 ! Initialize CCPP error handling variables
59 errmsg = ''
60 errflg = 0
61
62 if (.not. do_ccpp_suite_sim) return
63
64 ! Current forecast time (Data-format specific)
65 year = jdat(1)
66 month = jdat(2)
67 day = jdat(3)
68 hour = jdat(5)
69 min = jdat(6)
70 sec = jdat(7)
71
72 ! Set state at beginning of the physics timestep.
73 gt1(:,:) = tgrs(:,:)
74 gu1(:,:) = ugrs(:,:)
75 gv1(:,:) = vgrs(:,:)
76 gq1(:,:) = qgrs(:,:,1)
77 dtdt(:,:) = 0.
78 dudt(:,:) = 0.
79 dvdt(:,:) = 0.
80 dqdt(:,:) = 0.
81
82 !
83 ! Set bookeeping indices
84 !
85 if (in_pre_active) then
86 proc_start = 1
87 proc_end = max(1,physics_process(1)%iactive_scheme-1)
88 endif
89 if (in_post_active) then
90 proc_start = physics_process(1)%iactive_scheme
91 proc_end = size(physics_process)
92 endif
93
94 !
95 ! Simulate internal physics timestep evolution.
96 !
97 do iprc = proc_start,proc_end
98 do icol = 1,ncol
99
100 ! Reset locals
101 physics_process(iprc)%tend1d%T(:) = 0.
102 physics_process(iprc)%tend1d%u(:) = 0.
103 physics_process(iprc)%tend1d%v(:) = 0.
104 physics_process(iprc)%tend1d%q(:) = 0.
105
106 ! Using scheme simulator
107 ! Very simple...
108 ! Interpolate 2D data (time,level) tendency to local time.
109 ! Here the data is already on the SCM vertical coordinate.
110 !
111 ! In theory the data can be of any dimensionality and the onus falls on the
112 ! developer to extend the type "base_physics_process" to work with for their
113 ! application.
114 !
115 if (physics_process(iprc)%use_sim) then
116 if (physics_process(iprc)%name == "LWRAD") then
117 call sim_lwrad(year, month, day, hour, min, sec, physics_process(iprc))
118 endif
119 if (physics_process(iprc)%name == "SWRAD")then
120 call sim_swrad(year, month, day, hour, min, sec, physics_process(iprc))
121 endif
122 if (physics_process(iprc)%name == "GWD")then
123 call sim_gwd(year, month, day, hour, min, sec, physics_process(iprc))
124 endif
125 if (physics_process(iprc)%name == "PBL")then
126 call sim_pbl(year, month, day, hour, min, sec, physics_process(iprc))
127 endif
128 if (physics_process(iprc)%name == "SCNV")then
129 call sim_scnv(year, month, day, hour, min, sec, physics_process(iprc))
130 endif
131 if (physics_process(iprc)%name == "DCNV")then
132 call sim_dcnv(year, month, day, hour, min, sec, physics_process(iprc))
133 endif
134 if (physics_process(iprc)%name == "cldMP")then
135 call sim_cldmp(year, month, day, hour, min, sec, physics_process(iprc))
136 endif
137
138 ! Using data tendency from "active" scheme(s).
139 else
140 if (iactive_t > 0) physics_process(iprc)%tend1d%T = active_phys_tend(icol,:,iactive_t)
141 if (iactive_u > 0) physics_process(iprc)%tend1d%u = active_phys_tend(icol,:,iactive_u)
142 if (iactive_v > 0) physics_process(iprc)%tend1d%v = active_phys_tend(icol,:,iactive_v)
143 if (iactive_q > 0) physics_process(iprc)%tend1d%q = active_phys_tend(icol,:,iactive_q)
144 endif
145
146 ! Update state now? (time-split scheme)
147 if (physics_process(iprc)%time_split) then
148 gt1(icol,:) = gt1(icol,:) + (dtdt(icol,:) + physics_process(iprc)%tend1d%T)*dtp
149 gu1(icol,:) = gu1(icol,:) + (dudt(icol,:) + physics_process(iprc)%tend1d%u)*dtp
150 gv1(icol,:) = gv1(icol,:) + (dvdt(icol,:) + physics_process(iprc)%tend1d%v)*dtp
151 gq1(icol,:) = gq1(icol,:) + (dqdt(icol,:) + physics_process(iprc)%tend1d%q)*dtp
152 dtdt(icol,:) = 0.
153 dudt(icol,:) = 0.
154 dvdt(icol,:) = 0.
155 dqdt(icol,:) = 0.
156 ! Accumulate tendencies, update later? (process-split scheme)
157 else
158 dtdt(icol,:) = dtdt(icol,:) + physics_process(iprc)%tend1d%T
159 dudt(icol,:) = dudt(icol,:) + physics_process(iprc)%tend1d%u
160 dvdt(icol,:) = dvdt(icol,:) + physics_process(iprc)%tend1d%v
161 dqdt(icol,:) = dqdt(icol,:) + physics_process(iprc)%tend1d%q
162 endif
163 enddo ! END: Loop over columns
164
165 ! Print diagnostics
166 if (physics_process(iprc)%use_sim) then
167 if (physics_process(iprc)%time_split) then
168 write(*,'(a25,i2,a4,i2,a5,a10,a35)') 'CCPP suite simulator: ',iprc,' of ',proc_end,' ',physics_process(iprc)%name,'time split scheme (simulated)'
169 else
170 write(*,'(a25,i2,a4,i2,a5,a10,a35)') 'CCPP suite simulator: ',iprc,' of ',proc_end,' ',physics_process(iprc)%name,'process split scheme (simulated)'
171 endif
172 else
173 if (physics_process(iprc)%time_split) then
174 write(*,'(a25,i2,a4,i2,a5,a10,a35)') 'CCPP suite simulator: ',iprc,' of ',proc_end,' ',physics_process(iprc)%name,'time split scheme (active)'
175 else
176 write(*,'(a25,i2,a4,i2,a5,a10,a35)') 'CCPP suite simulator: ',iprc,' of ',proc_end,' ',physics_process(iprc)%name,'process split scheme (active)'
177 endif
178 write(*,'(a25,i2)') ' # prog. vars.: ',physics_process(1)%nprg_active
179 endif
180 enddo ! END: Loop over physics processes
181
182 !
183 ! Update state with accumulated tendencies (process-split only)
184 ! (Suites where active scheme is last physical process)
185 !
186 iprc = minval([iprc,proc_end])
187 if (.not. physics_process(iprc)%time_split) then
188 do icol = 1,ncol
189 gt0(icol,:) = gt1(icol,:) + dtdt(icol,:)*dtp
190 gu0(icol,:) = gu1(icol,:) + dudt(icol,:)*dtp
191 gv0(icol,:) = gv1(icol,:) + dvdt(icol,:)*dtp
192 gq0(icol,:) = gq1(icol,:) + dqdt(icol,:)*dtp
193 enddo
194 endif
195
196 !
197 ! Update bookeeping indices
198 !
199 if (in_pre_active) then
200 in_pre_active = .false.
201 in_post_active = .true.
202 endif
203
204 if (size(physics_process) == proc_end) then
205 in_pre_active = .true.
206 in_post_active = .false.
207 endif
208
209 end subroutine ccpp_suite_simulator_run
210
211end module ccpp_suite_simulator
This type contains the meta information and data for each physics process.