MONC
pwadvection.F90
Go to the documentation of this file.
1 
6  use collections_mod, only : map_type
7  use state_mod, only : model_state_type
8  use grids_mod, only : z_index, y_index, x_index
9 implicit none
10 
11 #ifndef TEST_MODE
12  private
13 #endif
14 
16 
17  logical :: l_toplevel=.true.
18 
20 contains
21 
25  pwadvection_get_descriptor%name="pw_advection"
26  pwadvection_get_descriptor%version=0.1
29  end function pwadvection_get_descriptor
30 
33  subroutine initialisation_callback(current_state)
34  type(model_state_type), target, intent(inout) :: current_state
35 
36  advect_flow=determine_if_advection_here(options_get_string(current_state%options_database, "advection_flow_fields"))
37  advect_th=determine_if_advection_here(options_get_string(current_state%options_database, "advection_theta_field"))
38  advect_q=determine_if_advection_here(options_get_string(current_state%options_database, "advection_q_fields"))
39  end subroutine initialisation_callback
40 
43  subroutine timestep_callback(current_state)
44  type(model_state_type), target, intent(inout) :: current_state
45 
46  integer :: current_x_index, current_y_index
47 
48  current_x_index=current_state%column_local_x
49  current_y_index=current_state%column_local_y
50 
51  if (current_state%halo_column) return
52 
53  if (advect_flow) call advect_flow_fields(current_state, current_x_index, current_y_index)
54  if (advect_th) call advect_th_field(current_state, current_x_index, current_y_index)
55  if (advect_q) call advect_q_field(current_state, current_x_index, current_y_index)
56  end subroutine timestep_callback
57 
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
65 
66  integer :: k, n
67 
68  do n=1,current_state%number_q_fields
69  do k=2,current_state%local_grid%size(z_index)-1
70 #ifdef U_ACTIVE
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))
77 #endif
78 #ifdef V_ACTIVE
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))
86 #endif
87 #ifdef W_ACTIVE
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))
96 #endif
97  end do
98 
99  if (l_toplevel)then
100  k=current_state%local_grid%size(z_index)
101 #ifdef U_ACTIVE
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))
108 #endif
109 #ifdef V_ACTIVE
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))
116 #endif
117 #ifdef W_ACTIVE
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)
123 #endif
124  end if
125  end do
126  end subroutine advect_q_field
127 
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
135 
136  integer :: k
137 
138  if (current_state%th%active) then
139 
140  do k=2,current_state%local_grid%size(z_index)-1
141 #ifdef U_ACTIVE
142  current_state%sth%data(k, current_y_index, current_x_index)= & !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))
148 #endif
149 #ifdef V_ACTIVE
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))
156 #endif
157 #ifdef W_ACTIVE
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))
165 #endif
166  end do
167 
168  if (l_toplevel)then
169 
170  k=current_state%local_grid%size(z_index)
171 #ifdef U_ACTIVE
172  current_state%sth%data(k, current_y_index, current_x_index)= & !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))
178 #endif
179 #ifdef V_ACTIVE
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))
186 #endif
187 #ifdef W_ACTIVE
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, &
191  current_x_index)
192 #endif
193  end if
194  end if
195  end subroutine advect_th_field
196 
201  subroutine advect_flow_fields(current_state, current_x_index, current_y_index)
202  type(model_state_type), target, intent(inout) :: current_state
203  integer, intent(in) :: current_x_index, current_y_index
204 
205  integer :: k
206 
207  do k=2,current_state%local_grid%size(z_index)-1
208 #ifdef U_ACTIVE
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)))
216 #ifdef V_ACTIVE
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)))
224 #endif
225 #ifdef W_ACTIVE
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)))
233 #endif
234 #endif
235 
236 #ifdef V_ACTIVE
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)))
244 #ifdef U_ACTIVE
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)))
252 #endif
253 #ifdef W_ACTIVE
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)))
261 #endif
262 #endif
263 
264 #ifdef W_ACTIVE
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)))
272 #ifdef U_ACTIVE
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)))
280 #endif
281 #ifdef V_ACTIVE
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)))
289 #endif
290 #endif
291  end do
292 
293  if (l_toplevel)then
294  k=current_state%local_grid%size(z_index)
295 #ifdef U_ACTIVE
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)))
303 #ifdef V_ACTIVE
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)))
311 #endif
312 #ifdef W_ACTIVE
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))
317 #endif
318 #endif
319 
320 #ifdef V_ACTIVE
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)))
328 #ifdef U_ACTIVE
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)))
336 #endif
337 #ifdef W_ACTIVE
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))
342 #endif
343 #endif
344  end if
345  end subroutine advect_flow_fields
346 
351  logical function determine_if_advection_here(field)
352  character(len=*), intent(in) :: field
353 
354  if (len_trim(field) .ne. 0) then
355  if (trim(field) .eq. "pw" .or. trim(field) .eq. "any") then
357  else
359  end if
360  else
362  end if
363  end function determine_if_advection_here
364 end module pwadvection_mod
subroutine advect_q_field(current_state, current_x_index, current_y_index)
Advects the q fields in a column.
Definition: pwadvection.F90:63
type(component_descriptor_type) function, public pwadvection_get_descriptor()
Provides a description of this component for the core to register.
Definition: pwadvection.F90:25
Piacsek-Williams advection scheme.
Definition: pwadvection.F90:2
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.
Definition: datadefn.F90:17
integer, parameter, public z_index
Grid index parameters.
Definition: grids.F90:14
Contains common definitions for the data and datatypes used by MONC.
Definition: datadefn.F90:2
The ModelState which represents the current state of a run.
Definition: state.F90:39
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.
Definition: collections.F90:86
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.
Definition: collections.F90:7
logical l_toplevel
Definition: pwadvection.F90:17
subroutine timestep_callback(current_state)
Called per column of data, this will perform Piacsek-Williams advection on the applicable fields for ...
Definition: pwadvection.F90:44
Functionality to support the different types of grid and abstraction between global grids and local o...
Definition: grids.F90:5
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...
logical advect_flow
Definition: pwadvection.F90:15
subroutine initialisation_callback(current_state)
Initialisation callback, will set up the configuration of this advection scheme.
Definition: pwadvection.F90:34
The model state which represents the current state of a run.
Definition: state.F90:2
integer, parameter, public y_index
Definition: grids.F90:14
integer, parameter, public x_index
Definition: grids.F90:14