34 type(model_state_type),
target,
intent(inout) :: current_state
44 type(model_state_type),
target,
intent(inout) :: current_state
46 integer :: current_x_index, current_y_index
48 current_x_index=current_state%column_local_x
49 current_y_index=current_state%column_local_y
51 if (current_state%halo_column)
return 62 subroutine advect_q_field(current_state, current_x_index, current_y_index)
63 type(model_state_type),
target,
intent(inout) :: current_state
64 integer,
intent(in) :: current_x_index, current_y_index
68 do n=1,current_state%number_q_fields
69 do k=2,current_state%local_grid%size(z_index)-1
71 current_state%sq(n)%data(k, current_y_index, current_x_index)=&
72 current_state%sq(n)%data(k, current_y_index, current_x_index)+current_state%global_grid%configuration%horizontal%cx*&
73 0.5_default_precision*(current_state%u%data(k, current_y_index, current_x_index-1)*&
74 current_state%q(n)%data(k, current_y_index, current_x_index-1)-&
75 current_state%u%data(k, current_y_index, current_x_index)*&
76 current_state%q(n)%data(k, current_y_index, current_x_index+1))
79 current_state%sq(n)%data(k, current_y_index, current_x_index)=&
80 current_state%sq(n)%data(k, current_y_index, current_x_index)+&
81 current_state%global_grid%configuration%horizontal%cy*0.5_default_precision*&
82 (current_state%v%data(k, current_y_index-1, current_x_index)*&
83 current_state%q(n)%data(k, current_y_index-1, current_x_index)-&
84 current_state%v%data(k, current_y_index, current_x_index)*&
85 current_state%q(n)%data(k, current_y_index+1, current_x_index))
88 current_state%sq(n)%data(k, current_y_index, current_x_index)=&
89 current_state%sq(n)%data(k, current_y_index, current_x_index)+&
90 2.0_default_precision*(current_state%global_grid%configuration%vertical%tzc1(k)*&
91 current_state%w%data(k-1, current_y_index, current_x_index)*&
92 current_state%q(n)%data(k-1, current_y_index, current_x_index)-&
93 current_state%global_grid%configuration%vertical%tzc2(k)*&
94 current_state%w%data(k, current_y_index, current_x_index)*&
95 current_state%q(n)%data(k+1, current_y_index, current_x_index))
100 k=current_state%local_grid%size(z_index)
102 current_state%sq(n)%data(k, current_y_index, current_x_index)=&
103 current_state%sq(n)%data(k, current_y_index, current_x_index)+current_state%global_grid%configuration%horizontal%cx*&
104 0.5_default_precision*(current_state%u%data(k, current_y_index, current_x_index-1)*&
105 current_state%q(n)%data(k, current_y_index, current_x_index-1)-&
106 current_state%u%data(k, current_y_index, current_x_index)*&
107 current_state%q(n)%data(k, current_y_index, current_x_index+1))
110 current_state%sq(n)%data(k, current_y_index, current_x_index)=&
111 current_state%sq(n)%data(k, current_y_index, current_x_index)+current_state%global_grid%configuration%horizontal%cy*&
112 0.5_default_precision*(current_state%v%data(k, current_y_index-1, current_x_index)*&
113 current_state%q(n)%data(k, current_y_index-1, current_x_index)-&
114 current_state%v%data(k, current_y_index, current_x_index)*&
115 current_state%q(n)%data(k, current_y_index+1, current_x_index))
118 current_state%sq(n)%data(k, current_y_index, current_x_index)=&
119 current_state%sq(n)%data(k, current_y_index, current_x_index)+&
120 current_state%global_grid%configuration%vertical%tzc1(k)*2.0_default_precision*&
121 current_state%w%data(k-1, current_y_index, current_x_index)*&
122 current_state%q(n)%data(k-1, current_y_index, current_x_index)
132 subroutine advect_th_field(current_state, current_x_index, current_y_index)
133 type(model_state_type),
target,
intent(inout) :: current_state
134 integer,
intent(in) :: current_x_index, current_y_index
138 if (current_state%th%active)
then 140 do k=2,current_state%local_grid%size(z_index)-1
142 current_state%sth%data(k, current_y_index, current_x_index)= &
143 current_state%global_grid%configuration%horizontal%cx*&
144 0.5_default_precision*(current_state%u%data(k, current_y_index, current_x_index-1)*&
145 current_state%th%data(k, current_y_index, current_x_index-1)-&
146 current_state%u%data(k, current_y_index, current_x_index)*&
147 current_state%th%data(k, current_y_index, current_x_index+1))
150 current_state%sth%data(k, current_y_index, current_x_index)=current_state%sth%data(k, current_y_index, current_x_index)+&
151 current_state%global_grid%configuration%horizontal%cy*0.5_default_precision*&
152 (current_state%v%data(k, current_y_index-1, current_x_index)*&
153 current_state%th%data(k, current_y_index-1, current_x_index)-&
154 current_state%v%data(k, current_y_index, current_x_index)*&
155 current_state%th%data(k, current_y_index+1, current_x_index))
158 current_state%sth%data(k, current_y_index, current_x_index)=current_state%sth%data(k, current_y_index, current_x_index)+&
159 2.0_default_precision*(current_state%global_grid%configuration%vertical%tzc1(k)*&
160 current_state%w%data(k-1, current_y_index, current_x_index)*&
161 current_state%th%data(k-1, current_y_index, current_x_index)-&
162 current_state%global_grid%configuration%vertical%tzc2(k)*&
163 current_state%w%data(k, current_y_index, current_x_index)*&
164 current_state%th%data(k+1, current_y_index, current_x_index))
170 k=current_state%local_grid%size(z_index)
172 current_state%sth%data(k, current_y_index, current_x_index)= &
173 current_state%global_grid%configuration%horizontal%cx*&
174 0.5_default_precision*(current_state%u%data(k, current_y_index, current_x_index-1)*&
175 current_state%th%data(k, current_y_index, current_x_index-1)-&
176 current_state%u%data(k, current_y_index, current_x_index)*&
177 current_state%th%data(k, current_y_index, current_x_index+1))
180 current_state%sth%data(k, current_y_index, current_x_index)=current_state%sth%data(k, current_y_index, current_x_index)+&
181 current_state%global_grid%configuration%horizontal%cy*0.5_default_precision*&
182 (current_state%v%data(k, current_y_index-1, current_x_index)*&
183 current_state%th%data(k, current_y_index-1, current_x_index)-&
184 current_state%v%data(k, current_y_index, current_x_index)*&
185 current_state%th%data(k, current_y_index+1, current_x_index))
188 current_state%sth%data(k, current_y_index, current_x_index)=current_state%sth%data(k, current_y_index, current_x_index)+&
189 current_state%global_grid%configuration%vertical%tzc1(k)*2.0_default_precision*&
190 current_state%w%data(k-1, current_y_index, current_x_index)*current_state%th%data(k-1, current_y_index, &
202 type(model_state_type),
target,
intent(inout) :: current_state
203 integer,
intent(in) :: current_x_index, current_y_index
207 do k=2,current_state%local_grid%size(z_index)-1
209 current_state%su%data(k, current_y_index, current_x_index)=&
210 current_state%global_grid%configuration%horizontal%tcx*(current_state%u%data(k, current_y_index, current_x_index-1)*&
211 (current_state%u%data(k, current_y_index, current_x_index)+&
212 current_state%u%data(k, current_y_index, current_x_index-1))-&
213 current_state%u%data(k, current_y_index, current_x_index+1)*&
214 (current_state%u%data(k, current_y_index, current_x_index)+&
215 current_state%u%data(k, current_y_index, current_x_index+1)))
217 current_state%su%data(k, current_y_index, current_x_index)=current_state%su%data(k, current_y_index, current_x_index)+&
218 current_state%global_grid%configuration%horizontal%tcy*(current_state%u%data(k, current_y_index-1, current_x_index)*&
219 (current_state%v%data(k, current_y_index-1, current_x_index)+&
220 current_state%v%data(k, current_y_index-1, current_x_index+1))-&
221 current_state%u%data(k, current_y_index+1, current_x_index)*&
222 (current_state%v%data(k, current_y_index, current_x_index)+&
223 current_state%v%data(k, current_y_index, current_x_index+1)))
226 current_state%su%data(k, current_y_index, current_x_index)=current_state%su%data(k, current_y_index, current_x_index)+&
227 (current_state%global_grid%configuration%vertical%tzc1(k)*current_state%u%data(k-1, current_y_index, current_x_index)*&
228 (current_state%w%data(k-1, current_y_index, current_x_index)+&
229 current_state%w%data(k-1, current_y_index, current_x_index+1))-&
230 current_state%global_grid%configuration%vertical%tzc2(k)*current_state%u%data(k+1, current_y_index, current_x_index)*&
231 (current_state%w%data(k, current_y_index, current_x_index)+&
232 current_state%w%data(k, current_y_index, current_x_index+1)))
237 current_state%sv%data(k, current_y_index, current_x_index)=current_state%global_grid%configuration%horizontal%tcy*(&
238 current_state%v%data(k, current_y_index-1, current_x_index)*&
239 (current_state%v%data(k, current_y_index, current_x_index)+&
240 current_state%v%data(k, current_y_index-1, current_x_index))-&
241 current_state%v%data(k, current_y_index+1, current_x_index)*&
242 (current_state%v%data(k, current_y_index, current_x_index)+&
243 current_state%v%data(k, current_y_index+1, current_x_index)))
245 current_state%sv%data(k, current_y_index, current_x_index)=current_state%sv%data(k, current_y_index, current_x_index)+&
246 current_state%global_grid%configuration%horizontal%tcx*(current_state%v%data(k, current_y_index, current_x_index-1)*&
247 (current_state%u%data(k, current_y_index, current_x_index-1)+&
248 current_state%u%data(k, current_y_index+1, current_x_index-1))-&
249 current_state%v%data(k, current_y_index, current_x_index+1)*&
250 (current_state%u%data(k, current_y_index, current_x_index)+&
251 current_state%u%data(k, current_y_index+1, current_x_index)))
254 current_state%sv%data(k, current_y_index, current_x_index)=current_state%sv%data(k, current_y_index, current_x_index)+&
255 (current_state%global_grid%configuration%vertical%tzc1(k)*current_state%v%data(k-1, current_y_index, current_x_index)*&
256 (current_state%w%data(k-1, current_y_index, current_x_index)+&
257 current_state%w%data(k-1, current_y_index+1, current_x_index))-&
258 current_state%global_grid%configuration%vertical%tzc2(k)*current_state%v%data(k+1, current_y_index, current_x_index)*&
259 (current_state%w%data(k, current_y_index, current_x_index)+&
260 current_state%w%data(k, current_y_index+1, current_x_index)))
265 current_state%sw%data(k, current_y_index, current_x_index)=(current_state%global_grid%configuration%vertical%tzd1(k)*&
266 current_state%w%data(k-1, current_y_index, current_x_index)*&
267 (current_state%w%data(k, current_y_index, current_x_index)+&
268 current_state%w%data(k-1, current_y_index, current_x_index))-&
269 current_state%global_grid%configuration%vertical%tzd2(k)*current_state%w%data(k+1, current_y_index, current_x_index)*&
270 (current_state%w%data(k, current_y_index, current_x_index)+&
271 current_state%w%data(k+1, current_y_index, current_x_index)))
273 current_state%sw%data(k, current_y_index, current_x_index)=current_state%sw%data(k, current_y_index, current_x_index)+&
274 current_state%global_grid%configuration%horizontal%tcx*(current_state%w%data(k, current_y_index, current_x_index-1)*&
275 (current_state%u%data(k, current_y_index, current_x_index-1)+&
276 current_state%u%data(k+1, current_y_index, current_x_index-1))-&
277 current_state%w%data(k, current_y_index, current_x_index+1)*&
278 (current_state%u%data(k, current_y_index, current_x_index)+&
279 current_state%u%data(k+1, current_y_index, current_x_index)))
282 current_state%sw%data(k, current_y_index, current_x_index)=current_state%sw%data(k, current_y_index, current_x_index)+&
283 current_state%global_grid%configuration%horizontal%tcy*(current_state%w%data(k, current_y_index-1, current_x_index)*&
284 (current_state%v%data(k, current_y_index-1, current_x_index)+&
285 current_state%v%data(k+1, current_y_index-1, current_x_index))-&
286 current_state%w%data(k, current_y_index+1, current_x_index)*&
287 (current_state%v%data(k, current_y_index, current_x_index)+&
288 current_state%v%data(k+1, current_y_index, current_x_index)))
294 k=current_state%local_grid%size(z_index)
296 current_state%su%data(k, current_y_index, current_x_index)=current_state%global_grid%configuration%horizontal%tcx*&
297 (current_state%u%data(k, current_y_index, current_x_index-1)*&
298 (current_state%u%data(k, current_y_index, current_x_index)+&
299 current_state%u%data(k, current_y_index, current_x_index-1))-&
300 current_state%u%data(k, current_y_index, current_x_index+1)*&
301 (current_state%u%data(k, current_y_index, current_x_index)+&
302 current_state%u%data(k, current_y_index, current_x_index+1)))
304 current_state%su%data(k, current_y_index, current_x_index)=current_state%su%data(k, current_y_index, current_x_index)+&
305 current_state%global_grid%configuration%horizontal%tcy*(current_state%u%data(k, current_y_index-1, current_x_index)*&
306 (current_state%v%data(k, current_y_index-1, current_x_index)+&
307 current_state%v%data(k, current_y_index-1, current_x_index+1))-&
308 current_state%u%data(k, current_y_index+1, current_x_index)*&
309 (current_state%v%data(k, current_y_index, current_x_index)+&
310 current_state%v%data(k, current_y_index, current_x_index+1)))
313 current_state%su%data(k, current_y_index, current_x_index)=current_state%su%data(k, current_y_index, current_x_index)+&
314 current_state%global_grid%configuration%vertical%tzc1(k)*current_state%u%data(k-1, current_y_index, current_x_index)*&
315 (current_state%w%data(k-1, current_y_index, current_x_index)+&
316 current_state%w%data(k-1, current_y_index, current_x_index+1))
321 current_state%sv%data(k, current_y_index, current_x_index)=current_state%global_grid%configuration%horizontal%tcy*&
322 (current_state%v%data(k, current_y_index-1, current_x_index)*&
323 (current_state%v%data(k, current_y_index, current_x_index)+&
324 current_state%v%data(k, current_y_index-1, current_x_index))-&
325 current_state%v%data(k, current_y_index+1, current_x_index)*&
326 (current_state%v%data(k, current_y_index, current_x_index)+&
327 current_state%v%data(k, current_y_index+1, current_x_index)))
329 current_state%sv%data(k, current_y_index, current_x_index)=current_state%sv%data(k, current_y_index, current_x_index)+&
330 current_state%global_grid%configuration%horizontal%tcx*(current_state%v%data(k, current_y_index, current_x_index-1)*&
331 (current_state%u%data(k, current_y_index, current_x_index-1)+&
332 current_state%u%data(k, current_y_index+1, current_x_index-1))-&
333 current_state%v%data(k, current_y_index, current_x_index+1)*&
334 (current_state%u%data(k, current_y_index, current_x_index)+&
335 current_state%u%data(k, current_y_index+1, current_x_index)))
338 current_state%sv%data(k, current_y_index, current_x_index)=current_state%sv%data(k, current_y_index, current_x_index)+&
339 current_state%global_grid%configuration%vertical%tzc1(k)*current_state%v%data(k-1, current_y_index, current_x_index)*&
340 (current_state%w%data(k-1, current_y_index, current_x_index)+&
341 current_state%w%data(k-1, current_y_index+1, current_x_index))
352 character(len=*),
intent(in) :: field
354 if (len_trim(field) .ne. 0)
then 355 if (trim(field) .eq.
"pw" .or. trim(field) .eq.
"any")
then
subroutine advect_q_field(current_state, current_x_index, current_y_index)
Advects the q fields in a column.
type(component_descriptor_type) function, public pwadvection_get_descriptor()
Provides a description of this component for the core to register.
Piacsek-Williams advection scheme.
character(len=string_length) function, public options_get_string(options_database, key, index)
Retrieves a string value from the database that matches the provided key.
integer, parameter, public default_precision
MPI communication type which we use for the prognostic and calculation data.
integer, parameter, public z_index
Grid index parameters.
Contains common definitions for the data and datatypes used by MONC.
The ModelState which represents the current state of a run.
Description of a component.
subroutine advect_th_field(current_state, current_x_index, current_y_index)
Advects the theta field in a column.
Map data structure that holds string (length 20 maximum) key value pairs.
subroutine advect_flow_fields(current_state, current_x_index, current_y_index)
Advects the flow fields depending upon which fields are active in the model in a column.
Interfaces and types that MONC components must specify.
Collection data structures.
subroutine timestep_callback(current_state)
Called per column of data, this will perform Piacsek-Williams advection on the applicable fields for ...
Functionality to support the different types of grid and abstraction between global grids and local o...
Manages the options database. Contains administration functions and deduce runtime options from the c...
logical function determine_if_advection_here(field)
Parses a field string (read in from the configuration file) and determines whether this algorithm sho...
subroutine initialisation_callback(current_state)
Initialisation callback, will set up the configuration of this advection scheme.
The model state which represents the current state of a run.
integer, parameter, public y_index
integer, parameter, public x_index