MOM_diag_buffers.F90

1!> Provides buffers that can dynamically grow as needed. These are primarily intended for the
2!! diagnostics which need to store intermediate or partial states of state variables
4
5use mom_io, only : stdout, stderr
6
7! This file is part of MOM6. See LICENSE.md for the license.
8
9implicit none ; private
10
12
13type, abstract :: buffer_base
14end type buffer_base
15
16!> Holds a 2d field
17type, extends(buffer_base) :: buffer_2d
18 real, dimension(:,:), allocatable :: field !< The actual 2d field to be stored [arbitrary]
19end type buffer_2d
20
21!> Holds a 3d field
22type, extends(buffer_base) :: buffer_3d
23 real, dimension(:,:,:), allocatable :: field !< The actual 3d field to be stored [arbitrary]
24end type buffer_3d
25
26!> The base class for the diagnostic buffers in this module
27type, abstract :: diag_buffer_base ; private
28 integer :: is !< The start slot of the array i-direction
29 integer :: js !< The start slot of the array j-direction
30 integer :: ie !< The end slot of the array i-direction
31 integer :: je !< The end slot of the array j-direction
32 real :: fill_value = 0. !< Set the fill value to use when growing the buffer [arbitrary]
33
34 integer, allocatable, dimension(:) :: ids !< List of diagnostic ids whose slot corresponds to the row in the buffer
35 integer :: length = 0 !< The number of slots in the buffer
36
37 contains
38
39 procedure(a_grow), deferred :: grow !< Increase the size of the buffer
40 procedure, public :: set_fill_value !< Set the fill value to use when growing the buffer
41 procedure, public :: check_capacity_by_id !< Check the size size of the buffer and increase if necessary
42 procedure, public :: set_horizontal_extents !< Define the horizontal extents of the arrays
43 procedure, public :: mark_available !< Mark that a slot in the buffer can be reused
44 procedure, public :: grow_ids !< Increase the size of the vector storing diagnostic ids
45 procedure, public :: find_buffer_slot !< Find the slot corresponding to a specific diagnostic id
46end type diag_buffer_base
47
48!> Dynamically growing buffer for 2D arrays.
49type, extends(diag_buffer_base), public :: diag_buffer_2d ; private
50 type(buffer_2d), public, dimension(:), allocatable :: buffer !< The actual 2D buffer which will dynamically grow
51
52 contains
53
54 procedure, public :: grow => grow_2d !< Increase the size of the buffer
55 procedure, public :: store => store_2d !< Store a field in the buffer, increasing as necessary
56 procedure, public :: set_extents_from_array => set_extents_from_array_2d !< Set extents from array bounds
57end type diag_buffer_2d
58
59!> Dynamically growing buffer for 3D arrays.
60type, extends(diag_buffer_base), public :: diag_buffer_3d ; private
61 type(buffer_3d), public, dimension(:), allocatable :: buffer !< The actual 2D buffer which will dynamically grow
62 integer :: ks !< The start slot in the k-dimension
63 integer :: ke !< The last slot in the k-dimension
64
65 contains
66
67 procedure, public :: set_vertical_extent !< Set the vertical extents of the buffer
68 procedure, public :: grow => grow_3d !< Increase the size of the buffer
69 procedure, public :: store => store_3d !< Store a field in the buffer, increasing as necessary
70 procedure, public :: set_extents_from_array => set_extents_from_array_3d !< Set extents from array bounds
71end type diag_buffer_3d
72
73contains
74
75!> Signature for the grow methods on n-dimension diagnostic buffer types
76subroutine a_grow(this)
77 class(diag_buffer_base), intent(inout) :: this !< The diagnostic buffer
78end subroutine
79
80!> Set the fill value to use when growing the buffer
81subroutine set_fill_value(this, fill_value)
82 class(diag_buffer_base), intent(inout) :: this !< The diagnostic buffer
83 real, intent(in) :: fill_value !< The fill value to use when growing the buffer [arbitrary]
84
85 this%fill_value = fill_value
86end subroutine set_fill_value
87
88!> Mark a slot in the buffer as unused based on a diagnostic id. For example,
89!! the data in that slot has already been consumed and can thus be overwritten
90subroutine mark_available(this, id)
91 class(diag_buffer_base), intent(inout) :: this !< The diagnostic buffer
92 integer, intent(in) :: id !< The diagnostic id
93 integer :: slot
94
95 slot = this%find_buffer_slot(id)
96 this%ids(slot) = 0
97end subroutine mark_available
98
99!> Return the slot of the buffer corresponding to the diagnostic id
100pure function find_buffer_slot(this, id) result(slot)
101 class(diag_buffer_base), intent(in) :: this !< The diagnostic buffer
102 integer, intent(in) :: id !< The diagnostic id
103
104 integer, dimension(1) :: temp
105 integer :: slot !< The slot in the buffer corresponding to the diagnostic id
106
107 if (allocated(this%ids)) then
108 !NOTE: Alternatively could do slot = SUM(findloc(...))
109 temp = findloc(this%ids(:), id)
110 slot = temp(1)
111 else
112 slot = 0
113 endif
114
115end function find_buffer_slot
116
117!> Grow the ids array by one
118subroutine grow_ids(this)
119 class(diag_buffer_base), intent(inout) :: this !< This buffer
120
121 integer, allocatable, dimension(:) :: temp
122 integer :: n
123
124 n = this%length
125
126 allocate(temp(n+1))
127 if (n>0) temp(1:n) = this%ids(:)
128 call move_alloc(temp, this%ids)
129end subroutine grow_ids
130
131!> Check whether the id already has a slot reserved. If not, find a new empty slot and if
132!! need be, grow the buffer.
133impure function check_capacity_by_id(this, id) result(slot)
134 class(diag_buffer_base), intent(inout) :: this !< This 2d buffer
135 integer, intent(in) :: id !< The diagnostic id
136 integer :: slot
137
138 slot = this%find_buffer_slot(id)
139 if (slot==0) then
140 ! Check to see if there is an open slot
141 if (allocated(this%ids)) slot = this%find_buffer_slot(0)
142 ! If slot is still 0, then the buffer must grow
143 if (slot==0) then
144 call this%grow()
145 slot = this%length
146 endif
147 this%ids(slot) = id
148 endif
149end function check_capacity_by_id
150
151!> Set the horizontal extents of the buffer
152subroutine set_horizontal_extents(this, is, ie, js, je)
153 class(diag_buffer_base), intent(inout) :: this !< The diagnostic buffer
154 integer, intent(in) :: is !< The start slot of the array i-direction
155 integer, intent(in) :: ie !< The end slot of the array i-direction
156 integer, intent(in) :: js !< The start slot of the array j-direction
157 integer, intent(in) :: je !< The end slot of the array j-direction
158
159 this%is = is ; this%ie = ie ; this%js = js ; this%je = je
160end subroutine set_horizontal_extents
161
162!> Set the vertical extent of the buffer
163subroutine set_vertical_extent(this, ks, ke)
164 class(diag_buffer_3d), intent(inout) :: this !< The diagnostic buffer
165 integer, intent(in) :: ks !< The start slot of the array i-direction
166 integer, intent(in) :: ke !< The end slot of the array i-direction
167
168 this%ks = ks ; this%ke = ke
169end subroutine set_vertical_extent
170
171!> Set the extents of a 2D buffer from the bounds of a 2D array
172subroutine set_extents_from_array_2d(this, array, fill_value_in)
173 class(diag_buffer_2d), intent(inout) :: this !< The diagnostic buffer
174 real, dimension(:,:), intent(in) :: array !< The array whose bounds define the buffer extents
175 real, optional, intent(in) :: fill_value_in !< Optional fill value
176
177 call this%set_horizontal_extents(lbound(array,1), ubound(array,1), &
178 lbound(array,2), ubound(array,2))
179 if (present(fill_value_in)) call this%set_fill_value(fill_value_in)
180end subroutine set_extents_from_array_2d
181
182!> Set the extents of a 3D buffer from the bounds of a 3D array
183subroutine set_extents_from_array_3d(this, array, fill_value_in)
184 class(diag_buffer_3d), intent(inout) :: this !< The diagnostic buffer
185 real, dimension(:,:,:), intent(in) :: array !< The array whose bounds define the buffer extents
186 real, optional, intent(in) :: fill_value_in !< Optional fill value
187
188 call this%set_horizontal_extents(lbound(array,1), ubound(array,1), &
189 lbound(array,2), ubound(array,2))
190 call this%set_vertical_extent(lbound(array,3), ubound(array,3))
191 if (present(fill_value_in)) call this%set_fill_value(fill_value_in)
192end subroutine set_extents_from_array_3d
193
194!> Grow a 2d diagnostic buffer
195subroutine grow_2d(this)
196 class(diag_buffer_2d), intent(inout) :: this
197
198 integer :: i, n
199 integer :: is, ie, js, je
200 type(buffer_2d), dimension(:), allocatable :: new_buffer
201
202 ! Grow the ID array
203 call this%grow_ids()
204
205 is = this%is ; ie = this%ie ; js = this%js ; je = this%je
206 n = this%length
207
208 allocate(new_buffer(n+1))
209 do i=1,n
210 allocate(new_buffer(i)%field(is:ie,js:je))
211 new_buffer(i)%field(:,:) = this%buffer(i)%field(:,:)
212 enddo
213 allocate(new_buffer(n+1)%field(is:ie,js:je), source=this%fill_value)
214 call move_alloc(new_buffer, this%buffer)
215 this%length = n+1
216
217end subroutine grow_2d
218
219!> Store a 2D array into this buffer
220subroutine store_2d(this, data, id)
221 class(diag_buffer_2d), intent(inout) :: this !< This 2d buffer
222 real, dimension(:,:), intent(in) :: data !< The data to be stored in the buffer [arbitrary]
223 integer, intent(in) :: id !< The diagnostic id
224
225 integer :: slot
226
227 slot = this%check_capacity_by_id(id)
228 this%buffer(slot)%field(:,:) = data(:,:)
229end subroutine store_2d
230
231!> Grow a 2d diagnostic buffer
232subroutine grow_3d(this)
233 class(diag_buffer_3d), intent(inout) :: this
234
235 integer :: i, n
236 integer :: is, ie, js, je, ks, ke
237 type(buffer_3d), dimension(:), allocatable :: new_buffer
238
239 ! Grow the ID array
240 call this%grow_ids()
241
242 is = this%is ; ie = this%ie ; js = this%js ; je = this%je ; ks = this%ks ; ke = this%ke
243 n = this%length
244
245 allocate(new_buffer(n+1))
246 do i=1,n
247 allocate(new_buffer(i)%field(is:ie,js:je,ks:ke))
248 new_buffer(i)%field(:,:,:) = this%buffer(i)%field(:,:,:)
249 enddo
250 allocate(new_buffer(n+1)%field(is:ie,js:je,ks:ke), source=this%fill_value)
251 call move_alloc(new_buffer, this%buffer)
252 this%length = n+1
253
254end subroutine grow_3d
255
256!> Store a 3d array into this buffer
257subroutine store_3d(this, data, id)
258 class(diag_buffer_3d), intent(inout) :: this !< This 3d buffer
259 real, dimension(:,:,:), intent(in) :: data !< The data to be stored in the buffer [arbitrary]
260 integer, intent(in) :: id !< The diagnostic id
261
262 integer :: slot
263
264 ! Find the first slot in the ids array that is 0, i.e. this is a portion of the buffer that can be reused
265 slot = this%check_capacity_by_id(id)
266 this%buffer(slot)%field(:,:,:) = data(:,:,:)
267end subroutine store_3d
268
269!> Unit tests for the 2d version of the diag buffer
270function diag_buffer_unit_tests_2d(verbose) result(fail)
271 logical, intent(in) :: verbose !< If true, write results to stdout
272 logical :: fail !< True if any of the unit tests fail
273
274 fail = .false.
275 write(stdout,*) '==== MOM_diag_buffers: diag_buffers_unit_tests_2d ==='
276 fail = fail .or. new_buffer_2d()
277 fail = fail .or. grow_buffer_2d()
278 fail = fail .or. fill_value_2d()
279 fail = fail .or. store_buffer_2d()
280 fail = fail .or. reuse_buffer_2d()
281
282 contains
283
284 !> Ensure properties of a newly initialized buffer
285 function new_buffer_2d() result(local_fail)
286 type(diag_buffer_2d) :: buffer
287 logical :: local_fail !< True if any of the unit tests fail
288 local_fail = .false.
289 local_fail = local_fail .or. allocated(buffer%buffer)
290 if (verbose) write(stdout,*) "new_buffer_2d: ", local_fail
291 local_fail = local_fail .or. allocated(buffer%ids)
292 if (verbose) write(stdout,*) "new_buffer_2d: ", local_fail
293 local_fail = local_fail .or. buffer%length /= 0
294 if (verbose) write(stdout,*) "new_buffer_2d: ", local_fail
295 end function new_buffer_2d
296
297 !> Test the growing of a buffer
298 function grow_buffer_2d() result(local_fail)
299 type(diag_buffer_2d) :: buffer
300 logical :: local_fail !< True if any of the unit tests fail
301 integer, parameter :: is=1, ie=2, js=3, je=6
302 integer :: i
303
304 local_fail = .false.
305
306 call buffer%set_horizontal_extents(is=is, ie=ie, js=js, je=je)
307 ! Grow the buffer 3 times
308 do i=1,3
309 call buffer%grow()
310 local_fail = local_fail .or. (buffer%length /= i)
311 local_fail = local_fail .or. (lbound(buffer%buffer(i)%field, 1) /= is)
312 local_fail = local_fail .or. (ubound(buffer%buffer(i)%field, 1) /= ie)
313 local_fail = local_fail .or. (lbound(buffer%buffer(i)%field, 2) /= js)
314 local_fail = local_fail .or. (ubound(buffer%buffer(i)%field, 2) /= je)
315 enddo
316 if (verbose) write(stdout,*) "grow_buffer_2d: ", local_fail
317 end function grow_buffer_2d
318
319 !> Test that growing new buffer fills the array with a set fill value
320 function fill_value_2d() result(local_fail)
321 type(diag_buffer_2d) :: buffer
322 logical :: local_fail !< True if any of the unit tests fail
323 integer, parameter :: is=1, ie=2, js=3, je=6
324 real, parameter :: fill_value = -123.456
325 integer :: i
326
327
328 local_fail = .false.
329
330 call buffer%set_horizontal_extents(is=is, ie=ie, js=js, je=je)
331 call buffer%set_fill_value(fill_value)
332 ! Grow the buffer 3 times
333 call buffer%grow()
334 if (any(buffer%buffer(1)%field(:,:) /= fill_value)) local_fail = .true.
335 if (verbose) write(stdout,*) "fill_value_2d: ", local_fail
336 end function fill_value_2d
337
338 !> Test storing a buffer based on a unique id
339 function store_buffer_2d() result(local_fail)
340 type(diag_buffer_2d) :: buffer
341 logical :: local_fail !< True if any of the unit tests fail
342
343 integer, parameter :: is=1, ie=2, js=3, je=6, nlen=3
344 integer :: i, slot
345 real, allocatable, dimension(:,:,:) :: test_2d
346
347 local_fail = .false.
348
349 allocate(test_2d(nlen, is:ie, js:je))
350 call random_number(test_2d)
351 buffer%is = is
352 buffer%ie = ie
353 buffer%js = js
354 buffer%je = je
355
356 do i=1,nlen
357 call buffer%store(test_2d(i,:,:), i*3)
358 slot = buffer%find_buffer_slot(i*3)
359 local_fail = local_fail .or. any(buffer%buffer(slot)%field(:,:) /= test_2d(i,:,:))
360 enddo
361
362 if (verbose) write(stdout,*) "store_buffer_2d: ", local_fail
363 end function store_buffer_2d
364
365 !> Test the reuse of a buffer. Fill it first like store_buffer_2d. Then,
366 !! loop through again, but use the slots of the buffer in the following
367 !! order: 2, 1, 3
368 function reuse_buffer_2d() result(local_fail)
369 type(diag_buffer_2d) :: buffer
370 logical :: local_fail !< True if any of the unit tests fail
371
372 integer, parameter :: is=1, ie=2, js=3, je=6, nlen=3
373 integer :: i, new_i, id, new_id
374 real, dimension(nlen, is:ie, js:je) :: test_2d_first, test_2d_second
375 integer, dimension(nlen) :: reorder = [2,1,3]
376
377 local_fail = .false.
378 call random_number(test_2d_first)
379 call random_number(test_2d_second)
380
381 call buffer%set_horizontal_extents(is=is, ie=ie, js=js, je=je)
382
383 do i=1,nlen
384 call buffer%store(test_2d_first(i,:,:), id=i*3)
385 enddo
386
387 do i=1,nlen
388 new_i = reorder(i)
389 ! id and new_id are multiplied by primes to make sure they are unique
390 id = reorder(i)*3
391 new_id = i*7
392 call buffer%mark_available(id=reorder(i)*3)
393 call buffer%store(test_2d_second(i,:,:), id=new_id)
394 local_fail = local_fail .or. buffer%find_buffer_slot(new_id) /= new_i
395 test_2d_first(new_i,:,:) = test_2d_second(i,:,:)
396 enddo
397 local_fail = local_fail .or. any(buffer%ids /= [14, 7, 21])
398 do i=1,nlen
399 local_fail = local_fail .or. any(buffer%buffer(i)%field(:,:) /= test_2d_first(i,:,:))
400 enddo
401 if (verbose) write(stdout,*) "reuse_buffer_2d: ", local_fail
402 end function reuse_buffer_2d
403
404end function diag_buffer_unit_tests_2d
405
406!> Test the 3d version of the buffer
407function diag_buffer_unit_tests_3d(verbose) result(fail)
408 logical, intent(in) :: verbose !< If true, write results to stdout
409 logical :: fail !< True if any of the unit tests fail
410
411 fail = .false.
412 write(stdout,*) '==== MOM_diag_buffers: diag_buffers_unit_tests_3d ==='
413 fail = fail .or. new_buffer_3d()
414 fail = fail .or. grow_buffer_3d()
415 fail = fail .or. fill_value_3d()
416 fail = fail .or. store_buffer_3d()
417 fail = fail .or. reuse_buffer_3d()
418
419 contains
420
421 !> Ensure properties of a newly initialized buffer
422 function new_buffer_3d() result(local_fail)
423 type(diag_buffer_3d) :: buffer
424 logical :: local_fail !< True if any of the unit tests fail
425 local_fail = .false.
426 local_fail = local_fail .or. allocated(buffer%buffer)
427 local_fail = local_fail .or. allocated(buffer%ids)
428 local_fail = local_fail .or. buffer%length /= 0
429 if (verbose) write(stdout,*) "new_buffer_3d: ", local_fail
430 end function new_buffer_3d
431
432 !> Test the growing of a buffer
433 function grow_buffer_3d() result(local_fail)
434 type(diag_buffer_3d) :: buffer
435 logical :: local_fail !< True if any of the unit tests fail
436 integer, parameter :: is=1, ie=2, js=3, je=6, ks=1, ke=10
437 integer :: i
438
439 local_fail = .false.
440
441 call buffer%set_horizontal_extents(is=is, ie=ie, js=js, je=je)
442 call buffer%set_vertical_extent(ks=ks, ke=ke)
443 ! Grow the buffer 3 times
444 do i=1,3
445 call buffer%grow()
446 local_fail = local_fail .or. (buffer%length /= i)
447 local_fail = local_fail .or. (lbound(buffer%buffer(i)%field, 1) /= is)
448 local_fail = local_fail .or. (ubound(buffer%buffer(i)%field, 1) /= ie)
449 local_fail = local_fail .or. (lbound(buffer%buffer(i)%field, 2) /= js)
450 local_fail = local_fail .or. (ubound(buffer%buffer(i)%field, 2) /= je)
451 local_fail = local_fail .or. (lbound(buffer%buffer(i)%field, 3) /= ks)
452 local_fail = local_fail .or. (ubound(buffer%buffer(i)%field, 3) /= ke)
453 if (verbose) write(stdout,*) "grow_buffer_3d: ", local_fail
454 enddo
455 if (verbose) write(stdout,*) "grow_buffer_3d: ", local_fail
456 end function grow_buffer_3d
457
458 !> Test that growing new buffer fills the array with a set fill value
459 function fill_value_3d() result(local_fail)
460 type(diag_buffer_3d) :: buffer
461 logical :: local_fail !< True if any of the unit tests fail
462 integer, parameter :: is=1, ie=2, js=3, je=6
463 real, parameter :: fill_value = -123.456
464 integer :: i
465
466
467 local_fail = .false.
468
469 call buffer%set_horizontal_extents(is=is, ie=ie, js=js, je=je)
470 call buffer%set_fill_value(fill_value)
471 ! Grow the buffer 3 times
472 call buffer%grow()
473 if (any(buffer%buffer(1)%field(:,:,:) /= fill_value)) local_fail = .true.
474 if (verbose) write(stdout,*) "fill_value_3d: ", local_fail
475 end function fill_value_3d
476
477 !> Test storing a buffer based on a unique id
478 function store_buffer_3d() result(local_fail)
479 type(diag_buffer_3d) :: buffer
480 logical :: local_fail !< True if any of the unit tests fail
481
482 integer, parameter :: is=1, ie=2, js=3, je=6, ks=1, ke=10, nlen=3
483 integer :: i, slot
484 real, dimension(nlen,is:ie,js:je,ks:ke) :: test_3d
485
486 local_fail = .false.
487 call random_number(test_3d)
488 buffer%is = is
489 buffer%ie = ie
490 buffer%js = js
491 buffer%je = je
492 buffer%ks = ks
493 buffer%ke = ke
494
495 do i=1,nlen
496 call buffer%store(test_3d(i,:,:,:), i*3)
497 slot = buffer%find_buffer_slot(i*3)
498 local_fail = local_fail .or. any(buffer%buffer(slot)%field(:,:,:) /= test_3d(i,:,:,:))
499 enddo
500
501 if (verbose) write(stdout,*) "store_buffer_3d: ", local_fail
502 end function store_buffer_3d
503
504 !> Test the reuse of a buffer. Fill it first like store_buffer_3d. Then,
505 !! loop through again, but use the slots of the buffer in the following
506 !! order: 2, 1, 3
507 function reuse_buffer_3d() result(local_fail)
508 type(diag_buffer_3d) :: buffer
509 logical :: local_fail !< True if any of the unit tests fail
510
511 integer, parameter :: is=1, ie=2, js=3, je=6, ks=1, ke=10, nlen=3
512 integer :: i, new_i, id, new_id
513 real, dimension(nlen, is:ie, js:je, ks:ke) :: test_3d_first, test_3d_second
514 integer, dimension(nlen) :: reorder = [2,1,3]
515
516 local_fail = .false.
517 call random_number(test_3d_first)
518 call random_number(test_3d_second)
519
520 buffer%is = is
521 buffer%ie = ie
522 buffer%js = js
523 buffer%je = je
524 buffer%ks = ks
525 buffer%ke = ke
526
527 do i=1,nlen
528 call buffer%store(test_3d_first(i,:,:,:), id=i*3)
529 enddo
530
531 do i=1,nlen
532 new_i = reorder(i)
533 ! id and new_id are multiplied by primes to make sure they are unique
534 id = reorder(i)*3
535 new_id = i*7
536 call buffer%mark_available(id=reorder(i)*3)
537 call buffer%store(test_3d_second(i,:,:,:), id=new_id)
538 local_fail = local_fail .or. buffer%find_buffer_slot(new_id) /= new_i
539 test_3d_first(new_i,:,:,:) = test_3d_second(i,:,:,:)
540 enddo
541 local_fail = local_fail .or. any(buffer%ids /= [14, 7, 21])
542 do i=1,nlen
543 local_fail = local_fail .or. any(buffer%buffer(i)%field(:,:,:) /= test_3d_first(i,:,:,:))
544 enddo
545 if (verbose) write(stdout,*) "reuse_buffer_3d: ", local_fail
546 end function reuse_buffer_3d
547
548end function diag_buffer_unit_tests_3d
549
550end module mom_diag_buffers
551