MOM6
MOM_spatial_means.F90
Go to the documentation of this file.
1 !> Functions and routines to take area, volume, mass-weighted, layerwise, zonal or meridional means
3 
4 ! This file is part of MOM6. See LICENSE.md for the license.
5 
6 use mom_coms, only : efp_type, operator(+), operator(-), assignment(=)
8 use mom_coms, only : reproducing_sum
10 use mom_error_handler, only : mom_error, note, warning, fatal, is_root_pe
12 use mom_grid, only : ocean_grid_type
14 
15 implicit none ; private
16 
17 #include <MOM_memory.h>
18 
21 public :: global_area_integral
24 
25 contains
26 
27 !> Return the global area mean of a variable. This uses reproducing sums.
28 function global_area_mean(var, G, scale)
29  type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
30  real, dimension(SZI_(G), SZJ_(G)), intent(in) :: var !< The variable to average
31  real, optional, intent(in) :: scale !< A rescaling factor for the variable
32 
33  real, dimension(SZI_(G), SZJ_(G)) :: tmpforsumming
34  real :: global_area_mean
35 
36  real :: scalefac ! An overall scaling factor for the areas and variable.
37  integer :: i, j, is, ie, js, je
38  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
39 
40  scalefac = g%US%L_to_m**2 ; if (present(scale)) scalefac = g%US%L_to_m**2*scale
41 
42  tmpforsumming(:,:) = 0.
43  do j=js,je ; do i=is,ie
44  tmpforsumming(i,j) = var(i,j) * (scalefac * g%areaT(i,j) * g%mask2dT(i,j))
45  enddo ; enddo
46  global_area_mean = reproducing_sum(tmpforsumming) * g%IareaT_global
47 
48 end function global_area_mean
49 
50 !> Return the global area integral of a variable. This uses reproducing sums.
51 function global_area_integral(var, G, scale)
52  type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
53  real, dimension(SZI_(G), SZJ_(G)), intent(in) :: var !< The variable to integrate
54  real, optional, intent(in) :: scale !< A rescaling factor for the variable
55  real, dimension(SZI_(G), SZJ_(G)) :: tmpforsumming
56  real :: global_area_integral
57 
58  real :: scalefac ! An overall scaling factor for the areas and variable.
59  integer :: i, j, is, ie, js, je
60  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
61 
62  scalefac = g%US%L_to_m**2 ; if (present(scale)) scalefac = g%US%L_to_m**2*scale
63 
64  tmpforsumming(:,:) = 0.
65  do j=js,je ; do i=is, ie
66  tmpforsumming(i,j) = var(i,j) * (scalefac * g%areaT(i,j) * g%mask2dT(i,j))
67  enddo ; enddo
68  global_area_integral = reproducing_sum(tmpforsumming)
69 
70 end function global_area_integral
71 
72 !> Return the layerwise global thickness-weighted mean of a variable. This uses reproducing sums.
73 function global_layer_mean(var, h, G, GV, scale)
74  type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
75  type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure
76  real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), intent(in) :: var !< The variable to average
77  real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2]
78  real, optional, intent(in) :: scale !< A rescaling factor for the variable
79  real, dimension(SZK_(GV)) :: global_layer_mean
80 
81  real, dimension(SZI_(G), SZJ_(G), SZK_(GV)) :: tmpforsumming, weight
82  real, dimension(SZK_(GV)) :: scalarij, weightij
83  real, dimension(SZK_(GV)) :: global_temp_scalar, global_weight_scalar
84  real :: scalefac ! A scaling factor for the variable.
85  integer :: i, j, k, is, ie, js, je, nz
86  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
87 
88  scalefac = 1.0 ; if (present(scale)) scalefac = scale
89  tmpforsumming(:,:,:) = 0. ; weight(:,:,:) = 0.
90 
91  do k=1,nz ; do j=js,je ; do i=is,ie
92  weight(i,j,k) = (gv%H_to_m * h(i,j,k)) * (g%US%L_to_m**2*g%areaT(i,j) * g%mask2dT(i,j))
93  tmpforsumming(i,j,k) = scalefac * var(i,j,k) * weight(i,j,k)
94  enddo ; enddo ; enddo
95 
96  global_temp_scalar = reproducing_sum(tmpforsumming,sums=scalarij)
97  global_weight_scalar = reproducing_sum(weight,sums=weightij)
98 
99  do k=1, nz
100  global_layer_mean(k) = scalarij(k) / weightij(k)
101  enddo
102 
103 end function global_layer_mean
104 
105 !> Find the global thickness-weighted mean of a variable. This uses reproducing sums.
106 function global_volume_mean(var, h, G, GV, scale)
107  type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
108  type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure
109  real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), &
110  intent(in) :: var !< The variable being averaged
111  real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), &
112  intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2]
113  real, optional, intent(in) :: scale !< A rescaling factor for the variable
114  real :: global_volume_mean !< The thickness-weighted average of var
115 
116  real :: scalefac ! A scaling factor for the variable.
117  real :: weight_here
118  real, dimension(SZI_(G), SZJ_(G)) :: tmpforsumming, sum_weight
119  integer :: i, j, k, is, ie, js, je, nz
120  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
121 
122  scalefac = 1.0 ; if (present(scale)) scalefac = scale
123  tmpforsumming(:,:) = 0. ; sum_weight(:,:) = 0.
124 
125  do k=1,nz ; do j=js,je ; do i=is,ie
126  weight_here = (gv%H_to_m * h(i,j,k)) * (g%US%L_to_m**2*g%areaT(i,j) * g%mask2dT(i,j))
127  tmpforsumming(i,j) = tmpforsumming(i,j) + scalefac * var(i,j,k) * weight_here
128  sum_weight(i,j) = sum_weight(i,j) + weight_here
129  enddo ; enddo ; enddo
130  global_volume_mean = (reproducing_sum(tmpforsumming)) / &
131  (reproducing_sum(sum_weight))
132 
133 end function global_volume_mean
134 
135 
136 !> Find the global mass-weighted integral of a variable. This uses reproducing sums.
137 function global_mass_integral(h, G, GV, var, on_PE_only, scale)
138  type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
139  type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure
140  real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), &
141  intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2]
142  real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), &
143  optional, intent(in) :: var !< The variable being integrated
144  logical, optional, intent(in) :: on_pe_only !< If present and true, the sum is only
145  !! done on the local PE, and it is _not_ order invariant.
146  real, optional, intent(in) :: scale !< A rescaling factor for the variable
147  real :: global_mass_integral !< The mass-weighted integral of var (or 1) in
148  !! kg times the units of var
149 
150  real, dimension(SZI_(G), SZJ_(G)) :: tmpforsumming
151  real :: scalefac ! An overall scaling factor for the areas and variable.
152  logical :: global_sum
153  integer :: i, j, k, is, ie, js, je, nz
154  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
155 
156  scalefac = g%US%L_to_m**2 ; if (present(scale)) scalefac = g%US%L_to_m**2*scale
157  tmpforsumming(:,:) = 0.0
158 
159  if (present(var)) then
160  do k=1,nz ; do j=js,je ; do i=is,ie
161  tmpforsumming(i,j) = tmpforsumming(i,j) + var(i,j,k) * &
162  ((gv%H_to_kg_m2 * h(i,j,k)) * (scalefac*g%areaT(i,j) * g%mask2dT(i,j)))
163  enddo ; enddo ; enddo
164  else
165  do k=1,nz ; do j=js,je ; do i=is,ie
166  tmpforsumming(i,j) = tmpforsumming(i,j) + &
167  ((gv%H_to_kg_m2 * h(i,j,k)) * (scalefac*g%areaT(i,j) * g%mask2dT(i,j)))
168  enddo ; enddo ; enddo
169  endif
170  global_sum = .true. ; if (present(on_pe_only)) global_sum = .not.on_pe_only
171  if (global_sum) then
172  global_mass_integral = reproducing_sum(tmpforsumming)
173  else
175  do j=js,je ; do i=is,ie
176  global_mass_integral = global_mass_integral + tmpforsumming(i,j)
177  enddo ; enddo
178  endif
179 
180 end function global_mass_integral
181 
182 
183 !> Determine the global mean of a field along rows of constant i, returning it
184 !! in a 1-d array using the local indexing. This uses reproducing sums.
185 subroutine global_i_mean(array, i_mean, G, mask, scale, tmp_scale)
186  type(ocean_grid_type), intent(inout) :: g !< The ocean's grid structure
187  real, dimension(SZI_(G),SZJ_(G)), intent(in) :: array !< The variable being averaged
188  real, dimension(SZJ_(G)), intent(out) :: i_mean !< Global mean of array along its i-axis
189  real, dimension(SZI_(G),SZJ_(G)), &
190  optional, intent(in) :: mask !< An array used for weighting the i-mean
191  real, optional, intent(in) :: scale !< A rescaling factor for the output variable
192  real, optional, intent(in) :: tmp_scale !< A rescaling factor for the internal
193  !! calculations that is removed from the output
194 
195  ! Local variables
196  type(efp_type), allocatable, dimension(:) :: asum, mask_sum
197  real :: scalefac ! A scaling factor for the variable.
198  real :: unscale ! A factor for undoing any internal rescaling before output.
199  real :: mask_sum_r
200  integer :: is, ie, js, je, idg_off, jdg_off
201  integer :: i, j
202 
203  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
204  idg_off = g%idg_offset ; jdg_off = g%jdg_offset
205 
206  scalefac = 1.0 ; if (present(scale)) scalefac = scale
207  unscale = 1.0
208  if (present(tmp_scale)) then ; if (tmp_scale /= 0.0) then
209  scalefac = scalefac * tmp_scale ; unscale = 1.0 / tmp_scale
210  endif ; endif
212 
213  allocate(asum(g%jsg:g%jeg))
214  if (present(mask)) then
215  allocate(mask_sum(g%jsg:g%jeg))
216 
217  do j=g%jsg,g%jeg
218  asum(j) = real_to_efp(0.0) ; mask_sum(j) = real_to_efp(0.0)
219  enddo
220 
221  do i=is,ie ; do j=js,je
222  asum(j+jdg_off) = asum(j+jdg_off) + real_to_efp(scalefac*array(i,j)*mask(i,j))
223  mask_sum(j+jdg_off) = mask_sum(j+jdg_off) + real_to_efp(mask(i,j))
224  enddo ; enddo
225 
226  if (query_efp_overflow_error()) call mom_error(fatal, &
227  "global_i_mean overflow error occurred before sums across PEs.")
228 
229  call efp_list_sum_across_pes(asum(g%jsg:g%jeg), g%jeg-g%jsg+1)
230  call efp_list_sum_across_pes(mask_sum(g%jsg:g%jeg), g%jeg-g%jsg+1)
231 
232  if (query_efp_overflow_error()) call mom_error(fatal, &
233  "global_i_mean overflow error occurred during sums across PEs.")
234 
235  do j=js,je
236  mask_sum_r = efp_to_real(mask_sum(j+jdg_off))
237  if (mask_sum_r == 0.0 ) then ; i_mean(j) = 0.0 ; else
238  i_mean(j) = efp_to_real(asum(j+jdg_off)) / mask_sum_r
239  endif
240  enddo
241 
242  deallocate(mask_sum)
243  else
244  do j=g%jsg,g%jeg ; asum(j) = real_to_efp(0.0) ; enddo
245 
246  do i=is,ie ; do j=js,je
247  asum(j+jdg_off) = asum(j+jdg_off) + real_to_efp(scalefac*array(i,j))
248  enddo ; enddo
249 
250  if (query_efp_overflow_error()) call mom_error(fatal, &
251  "global_i_mean overflow error occurred before sum across PEs.")
252 
253  call efp_list_sum_across_pes(asum(g%jsg:g%jeg), g%jeg-g%jsg+1)
254 
255  if (query_efp_overflow_error()) call mom_error(fatal, &
256  "global_i_mean overflow error occurred during sum across PEs.")
257 
258  do j=js,je
259  i_mean(j) = efp_to_real(asum(j+jdg_off)) / real(g%ieg-g%isg+1)
260  enddo
261  endif
262 
263  if (unscale /= 1.0) then ; do j=js,je ; i_mean(j) = unscale*i_mean(j) ; enddo ; endif
264 
265  deallocate(asum)
266 
267 end subroutine global_i_mean
268 
269 !> Determine the global mean of a field along rows of constant j, returning it
270 !! in a 1-d array using the local indexing. This uses reproducing sums.
271 subroutine global_j_mean(array, j_mean, G, mask, scale, tmp_scale)
272  type(ocean_grid_type), intent(inout) :: g !< The ocean's grid structure
273  real, dimension(SZI_(G),SZJ_(G)), intent(in) :: array !< The variable being averaged
274  real, dimension(SZI_(G)), intent(out) :: j_mean !< Global mean of array along its j-axis
275  real, dimension(SZI_(G),SZJ_(G)), &
276  optional, intent(in) :: mask !< An array used for weighting the j-mean
277  real, optional, intent(in) :: scale !< A rescaling factor for the output variable
278  real, optional, intent(in) :: tmp_scale !< A rescaling factor for the internal
279  !! calculations that is removed from the output
280 
281  ! Local variables
282  type(efp_type), allocatable, dimension(:) :: asum, mask_sum
283  real :: mask_sum_r
284  real :: scalefac ! A scaling factor for the variable.
285  real :: unscale ! A factor for undoing any internal rescaling before output.
286  integer :: is, ie, js, je, idg_off, jdg_off
287  integer :: i, j
288 
289  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
290  idg_off = g%idg_offset ; jdg_off = g%jdg_offset
291 
292  scalefac = 1.0 ; if (present(scale)) scalefac = scale
293  unscale = 1.0
294  if (present(tmp_scale)) then ; if (tmp_scale /= 0.0) then
295  scalefac = scalefac * tmp_scale ; unscale = 1.0 / tmp_scale
296  endif ; endif
298 
299  allocate(asum(g%isg:g%ieg))
300  if (present(mask)) then
301  allocate (mask_sum(g%isg:g%ieg))
302 
303  do i=g%isg,g%ieg
304  asum(i) = real_to_efp(0.0) ; mask_sum(i) = real_to_efp(0.0)
305  enddo
306 
307  do i=is,ie ; do j=js,je
308  asum(i+idg_off) = asum(i+idg_off) + real_to_efp(scalefac*array(i,j)*mask(i,j))
309  mask_sum(i+idg_off) = mask_sum(i+idg_off) + real_to_efp(mask(i,j))
310  enddo ; enddo
311 
312  if (query_efp_overflow_error()) call mom_error(fatal, &
313  "global_j_mean overflow error occurred before sums across PEs.")
314 
315  call efp_list_sum_across_pes(asum(g%isg:g%ieg), g%ieg-g%isg+1)
316  call efp_list_sum_across_pes(mask_sum(g%isg:g%ieg), g%ieg-g%isg+1)
317 
318  if (query_efp_overflow_error()) call mom_error(fatal, &
319  "global_j_mean overflow error occurred during sums across PEs.")
320 
321  do i=is,ie
322  mask_sum_r = efp_to_real(mask_sum(i+idg_off))
323  if (mask_sum_r == 0.0 ) then ; j_mean(i) = 0.0 ; else
324  j_mean(i) = efp_to_real(asum(i+idg_off)) / mask_sum_r
325  endif
326  enddo
327 
328  deallocate(mask_sum)
329  else
330  do i=g%isg,g%ieg ; asum(i) = real_to_efp(0.0) ; enddo
331 
332  do i=is,ie ; do j=js,je
333  asum(i+idg_off) = asum(i+idg_off) + real_to_efp(scalefac*array(i,j))
334  enddo ; enddo
335 
336  if (query_efp_overflow_error()) call mom_error(fatal, &
337  "global_j_mean overflow error occurred before sum across PEs.")
338 
339  call efp_list_sum_across_pes(asum(g%isg:g%ieg), g%ieg-g%isg+1)
340 
341  if (query_efp_overflow_error()) call mom_error(fatal, &
342  "global_j_mean overflow error occurred during sum across PEs.")
343 
344  do i=is,ie
345  j_mean(i) = efp_to_real(asum(i+idg_off)) / real(g%jeg-g%jsg+1)
346  enddo
347  endif
348 
349  if (unscale /= 1.0) then ; do i=is,ie ; j_mean(i) = unscale*j_mean(i) ; enddo ; endif
350 
351  deallocate(asum)
352 
353 end subroutine global_j_mean
354 
355 !> Adjust 2d array such that area mean is zero without moving the zero contour
356 subroutine adjust_area_mean_to_zero(array, G, scaling, unit_scale)
357  type(ocean_grid_type), intent(in) :: g !< Grid structure
358  real, dimension(SZI_(G),SZJ_(G)), intent(inout) :: array !< 2D array to be adjusted
359  real, optional, intent(out) :: scaling !< The scaling factor used
360  real, optional, intent(in) :: unit_scale !< A rescaling factor for the variable
361  ! Local variables
362  real, dimension(SZI_(G), SZJ_(G)) :: posvals, negvals, areaxposvals, areaxnegvals
363  integer :: i,j
364  real :: scalefac ! A scaling factor for the variable.
365  real :: i_scalefac ! The Adcroft reciprocal of scalefac
366  real :: areaintposvals, areaintnegvals, posscale, negscale
367 
368  scalefac = 1.0 ; if (present(unit_scale)) scalefac = unit_scale
369  i_scalefac = 0.0 ; if (scalefac /= 0.0) i_scalefac = 1.0 / scalefac
370 
371  areaxposvals(:,:) = 0.
372  areaxnegvals(:,:) = 0.
373 
374  do j=g%jsc,g%jec ; do i=g%isc,g%iec
375  posvals(i,j) = max(0., scalefac*array(i,j))
376  areaxposvals(i,j) = g%US%L_to_m**2*g%areaT(i,j) * posvals(i,j)
377  negvals(i,j) = min(0., scalefac*array(i,j))
378  areaxnegvals(i,j) = g%US%L_to_m**2*g%areaT(i,j) * negvals(i,j)
379  enddo ; enddo
380 
381  areaintposvals = reproducing_sum( areaxposvals )
382  areaintnegvals = reproducing_sum( areaxnegvals )
383 
384  posscale = 0.0 ; negscale = 0.0
385  if ((areaintposvals>0.).and.(areaintnegvals<0.)) then ! Only adjust if possible
386  if (areaintposvals>-areaintnegvals) then ! Scale down positive values
387  posscale = - areaintnegvals / areaintposvals
388  do j=g%jsc,g%jec ; do i=g%isc,g%iec
389  array(i,j) = ((posscale * posvals(i,j)) + negvals(i,j)) * i_scalefac
390  enddo ; enddo
391  elseif (areaintposvals<-areaintnegvals) then ! Scale down negative values
392  negscale = - areaintposvals / areaintnegvals
393  do j=g%jsc,g%jec ; do i=g%isc,g%iec
394  array(i,j) = (posvals(i,j) + (negscale * negvals(i,j))) * i_scalefac
395  enddo ; enddo
396  endif
397  endif
398  if (present(scaling)) scaling = posscale - negscale
399 
400 end subroutine adjust_area_mean_to_zero
401 
402 end module mom_spatial_means
mom_spatial_means::global_volume_mean
real function, public global_volume_mean(var, h, G, GV, scale)
Find the global thickness-weighted mean of a variable. This uses reproducing sums.
Definition: MOM_spatial_means.F90:107
mom_spatial_means::global_i_mean
subroutine, public global_i_mean(array, i_mean, G, mask, scale, tmp_scale)
Determine the global mean of a field along rows of constant i, returning it in a 1-d array using the ...
Definition: MOM_spatial_means.F90:186
mom_spatial_means
Functions and routines to take area, volume, mass-weighted, layerwise, zonal or meridional means.
Definition: MOM_spatial_means.F90:2
mom_spatial_means::adjust_area_mean_to_zero
subroutine, public adjust_area_mean_to_zero(array, G, scaling, unit_scale)
Adjust 2d array such that area mean is zero without moving the zero contour.
Definition: MOM_spatial_means.F90:357
mom_verticalgrid
Provides a transparent vertical ocean grid type and supporting routines.
Definition: MOM_verticalGrid.F90:2
mom_coms::efp_type
The Extended Fixed Point (EFP) type provides a public interface for doing sums and taking differences...
Definition: MOM_coms.F90:64
mom_file_parser::log_version
An overloaded interface to log version information about modules.
Definition: MOM_file_parser.F90:109
mom_coms::query_efp_overflow_error
logical function, public query_efp_overflow_error()
Returns the status of the module's error flag.
Definition: MOM_coms.F90:603
mom_file_parser::param_file_type
A structure that can be parsed to read and document run-time parameters.
Definition: MOM_file_parser.F90:54
mom_coms::efp_to_real
real function, public efp_to_real(EFP1)
Return the real number that an extended-fixed-point number corresponds with.
Definition: MOM_coms.F90:650
mom_file_parser::get_param
An overloaded interface to read and log the values of various types of parameters.
Definition: MOM_file_parser.F90:102
mom_coms::real_to_efp
type(efp_type) function, public real_to_efp(val, overflow)
Return the extended-fixed-point number that a real number corresponds with.
Definition: MOM_coms.F90:674
mom_coms
Interfaces to non-domain-oriented communication subroutines, including the MOM6 reproducing sums faci...
Definition: MOM_coms.F90:3
mom_verticalgrid::verticalgrid_type
Describes the vertical ocean grid, including unit conversion factors.
Definition: MOM_verticalGrid.F90:24
mom_spatial_means::global_mass_integral
real function, public global_mass_integral(h, G, GV, var, on_PE_only, scale)
Find the global mass-weighted integral of a variable. This uses reproducing sums.
Definition: MOM_spatial_means.F90:138
mom_spatial_means::global_area_mean
real function, public global_area_mean(var, G, scale)
Return the global area mean of a variable. This uses reproducing sums.
Definition: MOM_spatial_means.F90:29
mom_file_parser
The MOM6 facility to parse input files for runtime parameters.
Definition: MOM_file_parser.F90:2
mom_grid
Provides the ocean grid type.
Definition: MOM_grid.F90:2
mom_spatial_means::global_area_integral
real function, public global_area_integral(var, G, scale)
Return the global area integral of a variable. This uses reproducing sums.
Definition: MOM_spatial_means.F90:52
mom_error_handler::is_root_pe
logical function, public is_root_pe()
This returns .true. if the current PE is the root PE.
Definition: MOM_error_handler.F90:44
mom_error_handler::mom_error
subroutine, public mom_error(level, message, all_print)
This provides a convenient interface for writing an mpp_error message with run-time filter based on a...
Definition: MOM_error_handler.F90:72
mom_coms::reproducing_sum
Find an accurate and order-invariant sum of distributed 2d or 3d fields.
Definition: MOM_coms.F90:54
mom_coms::reset_efp_overflow_error
subroutine, public reset_efp_overflow_error()
Reset the module's error flag to false.
Definition: MOM_coms.F90:609
mom_spatial_means::global_j_mean
subroutine, public global_j_mean(array, j_mean, G, mask, scale, tmp_scale)
Determine the global mean of a field along rows of constant j, returning it in a 1-d array using the ...
Definition: MOM_spatial_means.F90:272
mom_error_handler
Routines for error handling and I/O management.
Definition: MOM_error_handler.F90:2
mom_grid::ocean_grid_type
Ocean grid type. See mom_grid for details.
Definition: MOM_grid.F90:26
mom_spatial_means::global_layer_mean
real function, dimension(gv %ke), public global_layer_mean(var, h, G, GV, scale)
Return the layerwise global thickness-weighted mean of a variable. This uses reproducing sums.
Definition: MOM_spatial_means.F90:74
mom_coms::efp_list_sum_across_pes
subroutine, public efp_list_sum_across_pes(EFPs, nval, errors)
This subroutine does a sum across PEs of a list of EFP variables, returning the sums in place,...
Definition: MOM_coms.F90:698