MOM6
Functions/Subroutines
mom_diag_vkernels Module Reference

Detailed Description

Provides kernels for single-column interpolation, re-integration (re-mapping of integrated quantities) and intensive-variable remapping in the vertical.

Functions/Subroutines

subroutine, public interpolate_column (nsrc, h_src, u_src, ndest, h_dest, missing_value, u_dest)
 Linearly interpolate interface data, u_src, from grid h_src to a grid h_dest. More...
 
subroutine, public reintegrate_column (nsrc, h_src, uh_src, ndest, h_dest, missing_value, uh_dest)
 Conservatively calculate integrated data, uh_dest, on grid h_dest, from layer-integrated data, uh_src, on grid h_src. More...
 
logical function, public diag_vkernels_unit_tests (verbose)
 Returns true if any unit tests for module MOM_diag_vkernels fail. More...
 
logical function test_interp (verbose, missing_value, msg, nsrc, h_src, u_src, ndest, h_dest, u_true)
 Returns true if a test of interpolate_column() produces the wrong answer. More...
 
logical function test_reintegrate (verbose, missing_value, msg, nsrc, h_src, uh_src, ndest, h_dest, uh_true)
 Returns true if a test of reintegrate_column() produces the wrong answer. More...
 

Function/Subroutine Documentation

◆ diag_vkernels_unit_tests()

logical function, public mom_diag_vkernels::diag_vkernels_unit_tests ( logical, intent(in)  verbose)

Returns true if any unit tests for module MOM_diag_vkernels fail.

Parameters
[in]verboseIf true, write results to stdout

Definition at line 169 of file MOM_diag_vkernels.F90.

References test_interp(), and test_reintegrate().

Referenced by mom_unit_tests::unit_tests().

169  logical, intent(in) :: verbose !< If true, write results to stdout
170  ! Local variables
171  real, parameter :: mv=-9.999999999e9 ! Value to use for vanished layers
172  logical :: fail, v
173 
174  v = verbose
175 
176  write(0,*) '==== MOM_diag_kernels: diag_vkernels_unit_tests =========='
177  if (v) write(0,*) '- - - - - - - - - - interpolation tests - - - - - - - - -'
178 
179  fail = test_interp(v,mv,'Identity: 3 layer', &
180  3, (/1.,2.,3./), (/1.,2.,3.,4./), &
181  3, (/1.,2.,3./), (/1.,2.,3.,4./) )
182  diag_vkernels_unit_tests = fail
183 
184  fail = test_interp(v,mv,'A: 3 layer to 2', &
185  3, (/1.,1.,1./), (/1.,2.,3.,4./), &
186  2, (/1.5,1.5/), (/1.,2.5,4./) )
187  diag_vkernels_unit_tests = diag_vkernels_unit_tests .or. fail
188 
189  fail = test_interp(v,mv,'B: 2 layer to 3', &
190  2, (/1.5,1.5/), (/1.,4.,7./), &
191  3, (/1.,1.,1./), (/1.,3.,5.,7./) )
192  diag_vkernels_unit_tests = diag_vkernels_unit_tests .or. fail
193 
194  fail = test_interp(v,mv,'C: 3 layer (vanished middle) to 2', &
195  3, (/1.,0.,2./), (/1.,2.,2.,3./), &
196  2, (/1.,2./), (/1.,2.,3./) )
197  diag_vkernels_unit_tests = diag_vkernels_unit_tests .or. fail
198 
199  fail = test_interp(v,mv,'D: 3 layer (deep) to 3', &
200  3, (/1.,2.,3./), (/1.,2.,4.,7./), &
201  2, (/2.,2./), (/1.,3.,5./) )
202  diag_vkernels_unit_tests = diag_vkernels_unit_tests .or. fail
203 
204  fail = test_interp(v,mv,'E: 3 layer to 3 (deep)', &
205  3, (/1.,2.,4./), (/1.,2.,4.,8./), &
206  3, (/2.,3.,4./), (/1.,3.,6.,8./) )
207  diag_vkernels_unit_tests = diag_vkernels_unit_tests .or. fail
208 
209  fail = test_interp(v,mv,'F: 3 layer to 4 with vanished top/botton', &
210  3, (/1.,2.,4./), (/1.,2.,4.,8./), &
211  4, (/0.,2.,5.,0./), (/mv,1.,3.,8.,mv/) )
212  diag_vkernels_unit_tests = diag_vkernels_unit_tests .or. fail
213 
214  fail = test_interp(v,mv,'Fs: 3 layer to 4 with vanished top/botton (shallow)', &
215  3, (/1.,2.,4./), (/1.,2.,4.,8./), &
216  4, (/0.,2.,4.,0./), (/mv,1.,3.,7.,mv/) )
217  diag_vkernels_unit_tests = diag_vkernels_unit_tests .or. fail
218 
219  fail = test_interp(v,mv,'Fd: 3 layer to 4 with vanished top/botton (deep)', &
220  3, (/1.,2.,4./), (/1.,2.,4.,8./), &
221  4, (/0.,2.,6.,0./), (/mv,1.,3.,8.,mv/) )
222  diag_vkernels_unit_tests = diag_vkernels_unit_tests .or. fail
223 
224  if (v) write(0,*) '- - - - - - - - - - reintegration tests - - - - - - - - -'
225 
226  fail = test_reintegrate(v,mv,'Identity: 3 layer', &
227  3, (/1.,2.,3./), (/-5.,2.,1./), &
228  3, (/1.,2.,3./), (/-5.,2.,1./) )
229  diag_vkernels_unit_tests = diag_vkernels_unit_tests .or. fail
230 
231  fail = test_reintegrate(v,mv,'A: 3 layer to 2', &
232  3, (/2.,2.,2./), (/-5.,2.,1./), &
233  2, (/3.,3./), (/-4.,2./) )
234  diag_vkernels_unit_tests = diag_vkernels_unit_tests .or. fail
235 
236  fail = test_reintegrate(v,mv,'A: 3 layer to 2 (deep)', &
237  3, (/2.,2.,2./), (/-5.,2.,1./), &
238  2, (/3.,4./), (/-4.,2./) )
239  diag_vkernels_unit_tests = diag_vkernels_unit_tests .or. fail
240 
241  fail = test_reintegrate(v,mv,'A: 3 layer to 2 (shallow)', &
242  3, (/2.,2.,2./), (/-5.,2.,1./), &
243  2, (/3.,2./), (/-4.,1.5/) )
244  diag_vkernels_unit_tests = diag_vkernels_unit_tests .or. fail
245 
246  fail = test_reintegrate(v,mv,'B: 3 layer to 4 with vanished top/bottom', &
247  3, (/2.,2.,2./), (/-5.,2.,1./), &
248  4, (/0.,3.,3.,0./), (/0.,-4.,2.,0./) )
249  diag_vkernels_unit_tests = diag_vkernels_unit_tests .or. fail
250 
251  fail = test_reintegrate(v,mv,'C: 3 layer to 4 with vanished top//middle/bottom', &
252  3, (/2.,2.,2./), (/-5.,2.,1./), &
253  5, (/0.,3.,0.,3.,0./), (/0.,-4.,0.,2.,0./) )
254  diag_vkernels_unit_tests = diag_vkernels_unit_tests .or. fail
255 
256  fail = test_reintegrate(v,mv,'D: 3 layer to 3 (vanished)', &
257  3, (/2.,2.,2./), (/-5.,2.,1./), &
258  3, (/0.,0.,0./), (/0.,0.,0./) )
259  diag_vkernels_unit_tests = diag_vkernels_unit_tests .or. fail
260 
261  fail = test_reintegrate(v,mv,'D: 3 layer (vanished) to 3', &
262  3, (/0.,0.,0./), (/-5.,2.,1./), &
263  3, (/2.,2.,2./), (/mv, mv, mv/) )
264  diag_vkernels_unit_tests = diag_vkernels_unit_tests .or. fail
265 
266  fail = test_reintegrate(v,mv,'D: 3 layer (vanished) to 3 (vanished)', &
267  3, (/0.,0.,0./), (/-5.,2.,1./), &
268  3, (/0.,0.,0./), (/mv, mv, mv/) )
269  diag_vkernels_unit_tests = diag_vkernels_unit_tests .or. fail
270 
271  fail = test_reintegrate(v,mv,'D: 3 layer (vanished) to 3 (vanished)', &
272  3, (/0.,0.,0./), (/0.,0.,0./), &
273  3, (/0.,0.,0./), (/mv, mv, mv/) )
274  diag_vkernels_unit_tests = diag_vkernels_unit_tests .or. fail
275 
276  if (.not. fail) write(*,*) 'Pass'
277 
Here is the call graph for this function:
Here is the caller graph for this function:

◆ interpolate_column()

subroutine, public mom_diag_vkernels::interpolate_column ( integer, intent(in)  nsrc,
real, dimension(nsrc), intent(in)  h_src,
real, dimension(nsrc+1), intent(in)  u_src,
integer, intent(in)  ndest,
real, dimension(ndest), intent(in)  h_dest,
real, intent(in)  missing_value,
real, dimension(ndest+1), intent(inout)  u_dest 
)

Linearly interpolate interface data, u_src, from grid h_src to a grid h_dest.

Parameters
[in]nsrcNumber of source cells
[in]h_srcThickness of source cells
[in]u_srcValues at source cell interfaces
[in]ndestNumber of destination cells
[in]h_destThickness of destination cells
[in]missing_valueValue to assign in vanished cells
[in,out]u_destInterpolated value at destination cell interfaces

Definition at line 17 of file MOM_diag_vkernels.F90.

Referenced by mom_ale::ale_offline_inputs(), test_interp(), and mom_diag_remap::vertically_interpolate_diag_field().

17  integer, intent(in) :: nsrc !< Number of source cells
18  real, dimension(nsrc), intent(in) :: h_src !< Thickness of source cells
19  real, dimension(nsrc+1), intent(in) :: u_src !< Values at source cell interfaces
20  integer, intent(in) :: ndest !< Number of destination cells
21  real, dimension(ndest), intent(in) :: h_dest !< Thickness of destination cells
22  real, intent(in) :: missing_value !< Value to assign in vanished cells
23  real, dimension(ndest+1), intent(inout) :: u_dest !< Interpolated value at destination cell interfaces
24  ! Local variables
25  real :: x_dest ! Relative position of target interface
26  real :: dh ! Source cell thickness
27  real :: u1, u2 ! Values to interpolate between
28  real :: weight_a, weight_b ! Weights for interpolation
29  integer :: k_src, k_dest ! Index of cell in src and dest columns
30  logical :: still_vanished ! Used for figuring out what to mask as missing
31 
32  ! Initial values for the loop
33  still_vanished = .true.
34 
35  ! The following forces the "do while" loop to do one cycle that will set u1, u2, dh.
36  k_src = 0
37  dh = 0.
38  x_dest = 0.
39 
40  do k_dest=1, ndest+1
41  do while (dh<=x_dest .and. k_src<nsrc)
42  ! Move positions pointers forward until the interval 0 .. dh spans x_dest.
43  x_dest = x_dest - dh
44  k_src = k_src + 1
45  dh = h_src(k_src) ! Thickness of layer k_src
46 
47  ! Values that will be used for the interpolation.
48  u1 = u_src(k_src) ! Value on left of source cell
49  u2 = u_src(k_src+1) ! Value on right of source cell
50  enddo
51 
52  if (dh>0.) then
53  weight_a = max(0., ( dh - x_dest ) / dh) ! Weight of u1
54  weight_b = min(1., x_dest / dh) ! Weight of u2
55  u_dest(k_dest) = weight_a * u1 + weight_b * u2 ! Linear interpolation between u1 and u2
56  else
57  u_dest(k_dest) = 0.5 * ( u1 + u2 ) ! For a vanished layer we need to do something reasonable...
58  endif
59 
60  ! Mask vanished layers at the surface which would be under an ice-shelf.
61  ! TODO: Need to figure out what to do for an isopycnal coordinate diagnostic that could
62  ! also have vanished layers at the surface.
63  if (k_dest<=ndest) then
64  x_dest = x_dest + h_dest(k_dest) ! Position of interface k_dest+1
65  if (still_vanished .and. h_dest(k_dest)==0.) then
66  ! When the layer k_dest is vanished and all layers above are also vanished, the k_dest
67  ! interface value should be missing.
68  u_dest(k_dest) = missing_value
69  else
70  still_vanished = .false.
71  endif
72  endif
73 
74  enddo
75 
76  ! Mask vanished layers on topography
77  still_vanished = .true.
78  do k_dest=ndest, 1, -1
79  if (still_vanished .and. h_dest(k_dest)==0.) then
80  ! When the layer k_dest is vanished and all layers below are also vanished, the k_dest+1
81  ! interface value should be missing.
82  u_dest(k_dest+1) = missing_value
83  else
84  exit
85  endif
86  enddo
87 
Here is the caller graph for this function:

◆ reintegrate_column()

subroutine, public mom_diag_vkernels::reintegrate_column ( integer, intent(in)  nsrc,
real, dimension(nsrc), intent(in)  h_src,
real, dimension(nsrc), intent(in)  uh_src,
integer, intent(in)  ndest,
real, dimension(ndest), intent(in)  h_dest,
real, intent(in)  missing_value,
real, dimension(ndest), intent(inout)  uh_dest 
)

Conservatively calculate integrated data, uh_dest, on grid h_dest, from layer-integrated data, uh_src, on grid h_src.

Parameters
[in]nsrcNumber of source cells
[in]h_srcThickness of source cells
[in]uh_srcValues at source cell interfaces
[in]ndestNumber of destination cells
[in]h_destThickness of destination cells
[in]missing_valueValue to assign in vanished cells
[in,out]uh_destInterpolated value at destination cell interfaces

Definition at line 92 of file MOM_diag_vkernels.F90.

Referenced by mom_ale::ale_offline_inputs(), test_reintegrate(), and mom_diag_remap::vertically_reintegrate_diag_field().

92  integer, intent(in) :: nsrc !< Number of source cells
93  real, dimension(nsrc), intent(in) :: h_src !< Thickness of source cells
94  real, dimension(nsrc), intent(in) :: uh_src !< Values at source cell interfaces
95  integer, intent(in) :: ndest !< Number of destination cells
96  real, dimension(ndest), intent(in) :: h_dest !< Thickness of destination cells
97  real, intent(in) :: missing_value !< Value to assign in vanished cells
98  real, dimension(ndest), intent(inout) :: uh_dest !< Interpolated value at destination cell interfaces
99  ! Local variables
100  real :: x_dest ! Relative position of target interface
101  real :: h_src_rem, h_dest_rem, dh ! Incremental thicknesses
102  real :: uh_src_rem, uh_dest_rem, duh ! Incremental amounts of stuff
103  integer :: k_src, k_dest ! Index of cell in src and dest columns
104  integer :: iter
105  logical :: src_ran_out, src_exists
106 
107  uh_dest(:) = missing_value
108 
109  k_src = 0
110  k_dest = 0
111  h_dest_rem = 0.
112  h_src_rem = 0.
113  src_ran_out = .false.
114  src_exists = .false.
115 
116  do while(.true.)
117  if (h_src_rem==0. .and. k_src<nsrc) then
118  ! Supply is empty so move to the next source cell
119  k_src = k_src + 1
120  h_src_rem = h_src(k_src)
121  uh_src_rem = uh_src(k_src)
122  if (h_src_rem==0.) cycle
123  src_exists = .true. ! This stops us masking out the entire column
124  endif
125  if (h_dest_rem==0. .and. k_dest<ndest) then
126  ! Sink has no capacity so move to the next destination cell
127  k_dest = k_dest + 1
128  h_dest_rem = h_dest(k_dest)
129  uh_dest(k_dest) = 0.
130  if (h_dest_rem==0.) cycle
131  endif
132  if (k_src==nsrc .and. h_src_rem==0.) then
133  if (src_ran_out) exit ! This is the second time implying there is no more src
134  src_ran_out = .true.
135  cycle
136  endif
137  duh = 0.
138  if (h_src_rem<h_dest_rem) then
139  ! The source cell is fully within the destination cell
140  dh = h_src_rem
141  if (dh>0.) duh = uh_src_rem
142  h_src_rem = 0.
143  uh_src_rem = 0.
144  h_dest_rem = max(0., h_dest_rem - dh)
145  elseif (h_src_rem>h_dest_rem) then
146  ! Only part of the source cell can be used up
147  dh = h_dest_rem
148  duh = (dh / h_src_rem) * uh_src_rem
149  h_src_rem = max(0., h_src_rem - dh)
150  uh_src_rem = uh_src_rem - duh
151  h_dest_rem = 0.
152  else ! h_src_rem==h_dest_rem
153  ! The source cell exactly fits the destination cell
154  duh = uh_src_rem
155  h_src_rem = 0.
156  uh_src_rem = 0.
157  h_dest_rem = 0.
158  endif
159  uh_dest(k_dest) = uh_dest(k_dest) + duh
160  if (k_dest==ndest .and. (k_src==nsrc .or. h_dest_rem==0.)) exit
161  enddo
162 
163  if (.not. src_exists) uh_dest(1:ndest) = missing_value
164 
Here is the caller graph for this function:

◆ test_interp()

logical function mom_diag_vkernels::test_interp ( logical, intent(in)  verbose,
real, intent(in)  missing_value,
character(len=*), intent(in)  msg,
integer, intent(in)  nsrc,
real, dimension(nsrc), intent(in)  h_src,
real, dimension(nsrc+1), intent(in)  u_src,
integer, intent(in)  ndest,
real, dimension(ndest), intent(in)  h_dest,
real, dimension(ndest+1), intent(in)  u_true 
)
private

Returns true if a test of interpolate_column() produces the wrong answer.

Parameters
[in]verboseIf true, write results to stdout
[in]missing_valueValue to indicate missing data
[in]msgMessage to label test
[in]nsrcNumber of source cells
[in]h_srcThickness of source cells
[in]u_srcValues at source cell interfaces
[in]ndestNumber of destination cells
[in]h_destThickness of destination cells
[in]u_trueCorrect value at destination cell interfaces

Definition at line 282 of file MOM_diag_vkernels.F90.

References interpolate_column().

Referenced by diag_vkernels_unit_tests().

282  logical, intent(in) :: verbose !< If true, write results to stdout
283  real, intent(in) :: missing_value !< Value to indicate missing data
284  character(len=*), intent(in) :: msg !< Message to label test
285  integer, intent(in) :: nsrc !< Number of source cells
286  real, dimension(nsrc), intent(in) :: h_src !< Thickness of source cells
287  real, dimension(nsrc+1), intent(in) :: u_src !< Values at source cell interfaces
288  integer, intent(in) :: ndest !< Number of destination cells
289  real, dimension(ndest), intent(in) :: h_dest !< Thickness of destination cells
290  real, dimension(ndest+1), intent(in) :: u_true !< Correct value at destination cell interfaces
291  ! Local variables
292  real, dimension(ndest+1) :: u_dest ! Interpolated value at destination cell interfaces
293  integer :: k
294  real :: error
295  logical :: print_results
296 
297  ! Interpolate from src to dest
298  call interpolate_column(nsrc, h_src, u_src, ndest, h_dest, missing_value, u_dest)
299 
300  test_interp = .false.
301  do k=1,ndest+1
302  if (u_dest(k)/=u_true(k)) test_interp = .true.
303  enddo
304  if (verbose .or. test_interp) then
305  write(0,'(2a)') ' Test: ',msg
306  write(0,'(a3,3(a24))') 'k','u_result','u_true','error'
307  do k=1,ndest+1
308  error = u_dest(k)-u_true(k)
309  if (error==0.) then
310  write(0,'(i3,3(1pe24.16))') k,u_dest(k),u_true(k),u_dest(k)-u_true(k)
311  else
312  write(0,'(i3,3(1pe24.16),x,a)') k,u_dest(k),u_true(k),u_dest(k)-u_true(k),'<--- WRONG!'
313  endif
314  enddo
315  endif
Here is the call graph for this function:
Here is the caller graph for this function:

◆ test_reintegrate()

logical function mom_diag_vkernels::test_reintegrate ( logical, intent(in)  verbose,
real, intent(in)  missing_value,
character(len=*), intent(in)  msg,
integer, intent(in)  nsrc,
real, dimension(nsrc), intent(in)  h_src,
real, dimension(nsrc), intent(in)  uh_src,
integer, intent(in)  ndest,
real, dimension(ndest), intent(in)  h_dest,
real, dimension(ndest), intent(in)  uh_true 
)
private

Returns true if a test of reintegrate_column() produces the wrong answer.

Parameters
[in]verboseIf true, write results to stdout
[in]missing_valueValue to indicate missing data
[in]msgMessage to label test
[in]nsrcNumber of source cells
[in]h_srcThickness of source cells
[in]uh_srcValues of source cell stuff
[in]ndestNumber of destination cells
[in]h_destThickness of destination cells
[in]uh_trueCorrect value of destination cell stuff

Definition at line 320 of file MOM_diag_vkernels.F90.

References reintegrate_column().

Referenced by diag_vkernels_unit_tests().

320  logical, intent(in) :: verbose !< If true, write results to stdout
321  real, intent(in) :: missing_value !< Value to indicate missing data
322  character(len=*), intent(in) :: msg !< Message to label test
323  integer, intent(in) :: nsrc !< Number of source cells
324  real, dimension(nsrc), intent(in) :: h_src !< Thickness of source cells
325  real, dimension(nsrc), intent(in) :: uh_src !< Values of source cell stuff
326  integer, intent(in) :: ndest !< Number of destination cells
327  real, dimension(ndest), intent(in) :: h_dest !< Thickness of destination cells
328  real, dimension(ndest), intent(in) :: uh_true !< Correct value of destination cell stuff
329  ! Local variables
330  real, dimension(ndest) :: uh_dest ! Reintegrated value on destination cells
331  integer :: k
332  real :: error
333  logical :: print_results
334 
335  ! Interpolate from src to dest
336  call reintegrate_column(nsrc, h_src, uh_src, ndest, h_dest, missing_value, uh_dest)
337 
338  test_reintegrate = .false.
339  do k=1,ndest
340  if (uh_dest(k)/=uh_true(k)) test_reintegrate = .true.
341  enddo
342  if (verbose .or. test_reintegrate) then
343  write(0,'(2a)') ' Test: ',msg
344  write(0,'(a3,3(a24))') 'k','uh_result','uh_true','error'
345  do k=1,ndest
346  error = uh_dest(k)-uh_true(k)
347  if (error==0.) then
348  write(0,'(i3,3(1pe24.16))') k,uh_dest(k),uh_true(k),uh_dest(k)-uh_true(k)
349  else
350  write(0,'(i3,3(1pe24.16),x,a)') k,uh_dest(k),uh_true(k),uh_dest(k)-uh_true(k),'<--- WRONG!'
351  endif
352  enddo
353  endif
Here is the call graph for this function:
Here is the caller graph for this function: