MOM6
MOM_checksum_packages.F90
Go to the documentation of this file.
1 !> Provides routines that do checksums of groups of MOM variables
3 
4 ! This file is part of MOM6. See LICENSE.md for the license.
5 
6 ! This module provides several routines that do check-sums of groups
7 ! of variables in the various dynamic solver routines.
8 
9 use mom_debugging, only : hchksum, uvchksum
10 use mom_domains, only : sum_across_pes, min_across_pes, max_across_pes
12 use mom_grid, only : ocean_grid_type
16 
17 implicit none ; private
18 
21 
22 !> Write out checksums of the MOM6 state variables
24  module procedure mom_state_chksum_5arg
25  module procedure mom_state_chksum_3arg
26 end interface
27 
28 #include <MOM_memory.h>
29 
30 !> A type for storing statistica about a variable
31 type :: stats ; private
32  real :: minimum = 1.e34 !< The minimum value
33  real :: maximum = -1.e34 !< The maximum value
34  real :: average = 0. !< The average value
35 end type stats
36 
37 contains
38 
39 ! =============================================================================
40 
41 !> Write out chksums for the model's basic state variables, including transports.
42 subroutine mom_state_chksum_5arg(mesg, u, v, h, uh, vh, G, GV, US, haloshift, symmetric, vel_scale)
43  character(len=*), &
44  intent(in) :: mesg !< A message that appears on the chksum lines.
45  type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure.
46  type(verticalgrid_type), intent(in) :: GV !< The ocean's vertical grid structure.
47  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)), &
48  intent(in) :: u !< The zonal velocity [L T-1 ~> m s-1] or other units.
49  real, dimension(SZI_(G),SZJB_(G),SZK_(G)), &
50  intent(in) :: v !< The meridional velocity [L T-1 ~> m s-1] or other units.
51  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
52  intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2].
53  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)), &
54  intent(in) :: uh !< Volume flux through zonal faces = u*h*dy
55  !! [H L2 T-1 ~> m3 s-1 or kg s-1].
56  real, dimension(SZI_(G),SZJB_(G),SZK_(G)), &
57  intent(in) :: vh !< Volume flux through meridional faces = v*h*dx
58  !! [H L2 T-1 ~> m3 s-1 or kg s-1].
59  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
60  integer, optional, intent(in) :: haloshift !< The width of halos to check (default 0).
61  logical, optional, intent(in) :: symmetric !< If true, do checksums on the fully symmetric
62  !! computational domain.
63  real, optional, intent(in) :: vel_scale !< The scaling factor to convert velocities to [m s-1]
64 
65  real :: scale_vel ! The scaling factor to convert velocities to [m s-1]
66  logical :: sym
67  integer :: is, ie, js, je, nz, hs
68  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
69 
70  ! Note that for the chksum calls to be useful for reproducing across PE
71  ! counts, there must be no redundant points, so all variables use is..ie
72  ! and js...je as their extent.
73  hs = 1 ; if (present(haloshift)) hs=haloshift
74  sym = .false. ; if (present(symmetric)) sym=symmetric
75  scale_vel = us%L_T_to_m_s ; if (present(vel_scale)) scale_vel = vel_scale
76 
77  call uvchksum(mesg//" [uv]", u, v, g%HI, haloshift=hs, symmetric=sym, scale=scale_vel)
78  call hchksum(h, mesg//" h", g%HI, haloshift=hs, scale=gv%H_to_m)
79  call uvchksum(mesg//" [uv]h", uh, vh, g%HI, haloshift=hs, &
80  symmetric=sym, scale=gv%H_to_m*us%L_to_m**2*us%s_to_T)
81 end subroutine mom_state_chksum_5arg
82 
83 ! =============================================================================
84 
85 !> Write out chksums for the model's basic state variables.
86 subroutine mom_state_chksum_3arg(mesg, u, v, h, G, GV, US, haloshift, symmetric)
87  character(len=*), intent(in) :: mesg !< A message that appears on the chksum lines.
88  type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure.
89  type(verticalgrid_type), intent(in) :: GV !< The ocean's vertical grid structure.
90  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)), &
91  intent(in) :: u !< Zonal velocity [L T-1 ~> m s-1] or [m s-1].
92  real, dimension(SZI_(G),SZJB_(G),SZK_(G)), &
93  intent(in) :: v !< Meridional velocity [L T-1 ~> m s-1] or [m s-1]..
94  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
95  intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2].
96  type(unit_scale_type), optional, intent(in) :: US !< A dimensional unit scaling type, which is
97  !! used to rescale u and v if present.
98  integer, optional, intent(in) :: haloshift !< The width of halos to check (default 0).
99  logical, optional, intent(in) :: symmetric !< If true, do checksums on the fully
100  !! symmetric computational domain.
101  real :: L_T_to_m_s ! A rescaling factor for velocities [m T s-1 L-1 ~> nondim] or [nondim]
102  integer :: is, ie, js, je, nz, hs
103  logical :: sym
104 
105  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
106  l_t_to_m_s = 1.0 ; if (present(us)) l_t_to_m_s = us%L_T_to_m_s
107 
108  ! Note that for the chksum calls to be useful for reproducing across PE
109  ! counts, there must be no redundant points, so all variables use is..ie
110  ! and js...je as their extent.
111  hs=1; if (present(haloshift)) hs=haloshift
112  sym=.false.; if (present(symmetric)) sym=symmetric
113  call uvchksum(mesg//" u", u, v, g%HI, haloshift=hs, symmetric=sym, scale=l_t_to_m_s)
114  call hchksum(h, mesg//" h",g%HI, haloshift=hs, scale=gv%H_to_m)
115 end subroutine mom_state_chksum_3arg
116 
117 ! =============================================================================
118 
119 !> Write out chksums for the model's thermodynamic state variables.
120 subroutine mom_thermo_chksum(mesg, tv, G, US, haloshift)
121  character(len=*), intent(in) :: mesg !< A message that appears on the chksum lines.
122  type(thermo_var_ptrs), intent(in) :: tv !< A structure pointing to various
123  !! thermodynamic variables.
124  type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure.
125  type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
126  integer, optional, intent(in) :: haloshift !< The width of halos to check (default 0).
127 
128  integer :: is, ie, js, je, nz, hs
129  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
130  hs=1; if (present(haloshift)) hs=haloshift
131 
132  if (associated(tv%T)) call hchksum(tv%T, mesg//" T",g%HI,haloshift=hs)
133  if (associated(tv%S)) call hchksum(tv%S, mesg//" S",g%HI,haloshift=hs)
134  if (associated(tv%frazil)) call hchksum(tv%frazil, mesg//" frazil",g%HI,haloshift=hs)
135  if (associated(tv%salt_deficit)) &
136  call hchksum(tv%salt_deficit, mesg//" salt deficit",g%HI,haloshift=hs, scale=us%R_to_kg_m3*us%Z_to_m)
137 
138 end subroutine mom_thermo_chksum
139 
140 ! =============================================================================
141 
142 !> Write out chksums for the ocean surface variables.
143 subroutine mom_surface_chksum(mesg, sfc, G, haloshift, symmetric)
144  character(len=*), intent(in) :: mesg !< A message that appears on the chksum lines.
145  type(surface), intent(inout) :: sfc !< transparent ocean surface state
146  !! structure shared with the calling routine
147  !! data in this structure is intent out.
148  type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure.
149  integer, optional, intent(in) :: haloshift !< The width of halos to check (default 0).
150  logical, optional, intent(in) :: symmetric !< If true, do checksums on the fully symmetric
151  !! computational domain.
152 
153  integer :: hs
154  logical :: sym
155 
156  sym = .false. ; if (present(symmetric)) sym = symmetric
157  hs = 1 ; if (present(haloshift)) hs = haloshift
158 
159  if (allocated(sfc%SST)) call hchksum(sfc%SST, mesg//" SST",g%HI,haloshift=hs)
160  if (allocated(sfc%SSS)) call hchksum(sfc%SSS, mesg//" SSS",g%HI,haloshift=hs)
161  if (allocated(sfc%sea_lev)) call hchksum(sfc%sea_lev, mesg//" sea_lev",g%HI,haloshift=hs)
162  if (allocated(sfc%Hml)) call hchksum(sfc%Hml, mesg//" Hml",g%HI,haloshift=hs)
163  if (allocated(sfc%u) .and. allocated(sfc%v)) &
164  call uvchksum(mesg//" SSU", sfc%u, sfc%v, g%HI, haloshift=hs, symmetric=sym)
165 ! if (allocated(sfc%salt_deficit)) call hchksum(sfc%salt_deficit, mesg//" salt deficit",G%HI,haloshift=hs)
166  if (associated(sfc%frazil)) call hchksum(sfc%frazil, mesg//" frazil",g%HI,haloshift=hs)
167 
168 end subroutine mom_surface_chksum
169 
170 ! =============================================================================
171 
172 !> Write out chksums for the model's accelerations
173 subroutine mom_accel_chksum(mesg, CAu, CAv, PFu, PFv, diffu, diffv, G, GV, US, pbce, &
174  u_accel_bt, v_accel_bt, symmetric)
175  character(len=*), intent(in) :: mesg !< A message that appears on the chksum lines.
176  type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure.
177  type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure.
178  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)), &
179  intent(in) :: cau !< Zonal acceleration due to Coriolis
180  !! and momentum advection terms [L T-2 ~> m s-2].
181  real, dimension(SZI_(G),SZJB_(G),SZK_(G)), &
182  intent(in) :: cav !< Meridional acceleration due to Coriolis
183  !! and momentum advection terms [L T-2 ~> m s-2].
184  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)), &
185  intent(in) :: pfu !< Zonal acceleration due to pressure gradients
186  !! (equal to -dM/dx) [L T-2 ~> m s-2].
187  real, dimension(SZI_(G),SZJB_(G),SZK_(G)), &
188  intent(in) :: pfv !< Meridional acceleration due to pressure gradients
189  !! (equal to -dM/dy) [L T-2 ~> m s-2].
190  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)), &
191  intent(in) :: diffu !< Zonal acceleration due to convergence of the
192  !! along-isopycnal stress tensor [L T-2 ~> m s-2].
193  real, dimension(SZI_(G),SZJB_(G),SZK_(G)), &
194  intent(in) :: diffv !< Meridional acceleration due to convergence of
195  !! the along-isopycnal stress tensor [L T-2 ~> m s-2].
196  type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
197  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
198  optional, intent(in) :: pbce !< The baroclinic pressure anomaly in each layer
199  !! due to free surface height anomalies
200  !! [L2 T-2 H-1 ~> m s-2 or m4 s-2 kg-1].
201  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)), &
202  optional, intent(in) :: u_accel_bt !< The zonal acceleration from terms in the
203  !! barotropic solver [L T-2 ~> m s-2].
204  real, dimension(SZI_(G),SZJB_(G),SZK_(G)), &
205  optional, intent(in) :: v_accel_bt !< The meridional acceleration from terms in
206  !! the barotropic solver [L T-2 ~> m s-2].
207  logical, optional, intent(in) :: symmetric !< If true, do checksums on the fully symmetric
208  !! computational domain.
209 
210  integer :: is, ie, js, je, nz
211  logical :: sym
212 
213  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
214  sym=.false.; if (present(symmetric)) sym=symmetric
215 
216  ! Note that for the chksum calls to be useful for reproducing across PE
217  ! counts, there must be no redundant points, so all variables use is..ie
218  ! and js...je as their extent.
219  call uvchksum(mesg//" CA[uv]", cau, cav, g%HI, haloshift=0, symmetric=sym, scale=us%L_T2_to_m_s2)
220  call uvchksum(mesg//" PF[uv]", pfu, pfv, g%HI, haloshift=0, symmetric=sym, scale=us%L_T2_to_m_s2)
221  call uvchksum(mesg//" diffu", diffu, diffv, g%HI,haloshift=0, symmetric=sym, scale=us%L_T2_to_m_s2)
222  if (present(pbce)) &
223  call hchksum(pbce, mesg//" pbce",g%HI,haloshift=0, scale=gv%m_to_H*us%L_T_to_m_s**2)
224  if (present(u_accel_bt) .and. present(v_accel_bt)) &
225  call uvchksum(mesg//" [uv]_accel_bt", u_accel_bt, v_accel_bt, g%HI,haloshift=0, symmetric=sym, &
226  scale=us%L_T2_to_m_s2)
227 end subroutine mom_accel_chksum
228 
229 ! =============================================================================
230 
231 !> Monitor and write out statistics for the model's state variables.
232 subroutine mom_state_stats(mesg, u, v, h, Temp, Salt, G, GV, US, allowChange, permitDiminishing)
233  type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure.
234  type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure.
235  character(len=*), intent(in) :: mesg !< A message that appears on the chksum lines.
236  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)), &
237  intent(in) :: u !< The zonal velocity [L T-1 ~> m s-1].
238  real, dimension(SZI_(G),SZJB_(G),SZK_(G)), &
239  intent(in) :: v !< The meridional velocity [L T-1 ~> m s-1].
240  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
241  intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2].
242  real, pointer, dimension(:,:,:), &
243  intent(in) :: temp !< Temperature [degC].
244  real, pointer, dimension(:,:,:), &
245  intent(in) :: salt !< Salinity [ppt].
246  type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
247  logical, optional, intent(in) :: allowchange !< do not flag an error
248  !! if the statistics change.
249  logical, optional, intent(in) :: permitdiminishing !< do not flag error if the
250  !! extrema are diminishing.
251 
252  ! Local variables
253  real :: vol, dv ! The total ocean volume and its change [m3] (unscaled to permit reproducing sum).
254  real :: area ! The total ocean surface area [m2] (unscaled to permit reproducing sum).
255  real :: h_minimum ! The minimum layer thicknesses [H ~> m or kg m-2]
256  logical :: do_ts ! If true, evaluate statistics for temperature and salinity
257  type(stats) :: t, s, delt, dels
258 
259  ! NOTE: save data is not normally allowed but we use it for debugging purposes here on the
260  ! assumption we will not turn this on with threads
261  type(stats), save :: oldt, olds
262  logical, save :: firstcall = .true.
263  real, save :: oldvol ! The previous total ocean volume [m3]
264 
265  character(len=80) :: lmsg
266  integer :: is, ie, js, je, nz, i, j, k
267 
268  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
269  do_ts = associated(temp) .and. associated(salt)
270 
271  ! First collect local stats
272  area = 0. ; vol = 0.
273  do j = js, je ; do i = is, ie
274  area = area + us%L_to_m**2*g%areaT(i,j)
275  enddo ; enddo
276  t%minimum = 1.e34 ; t%maximum = -1.e34 ; t%average = 0.
277  s%minimum = 1.e34 ; s%maximum = -1.e34 ; s%average = 0.
278  h_minimum = 1.e34*gv%m_to_H
279  do k = 1, nz ; do j = js, je ; do i = is, ie
280  if (g%mask2dT(i,j)>0.) then
281  dv = us%L_to_m**2*g%areaT(i,j)*gv%H_to_m*h(i,j,k) ; vol = vol + dv
282  if (do_ts .and. h(i,j,k)>0.) then
283  t%minimum = min( t%minimum, temp(i,j,k) ) ; t%maximum = max( t%maximum, temp(i,j,k) )
284  t%average = t%average + dv*temp(i,j,k)
285  s%minimum = min( s%minimum, salt(i,j,k) ) ; s%maximum = max( s%maximum, salt(i,j,k) )
286  s%average = s%average + dv*salt(i,j,k)
287  endif
288  if (h_minimum > h(i,j,k)) h_minimum = h(i,j,k)
289  endif
290  enddo ; enddo ; enddo
291  call sum_across_pes( area ) ; call sum_across_pes( vol )
292  if (do_ts) then
293  call min_across_pes( t%minimum ) ; call max_across_pes( t%maximum ) ; call sum_across_pes( t%average )
294  call min_across_pes( s%minimum ) ; call max_across_pes( s%maximum ) ; call sum_across_pes( s%average )
295  t%average = t%average / vol ; s%average = s%average / vol
296  endif
297  if (is_root_pe()) then
298  if (.not.firstcall) then
299  dv = vol - oldvol
300  delt%minimum = t%minimum - oldt%minimum ; delt%maximum = t%maximum - oldt%maximum
301  delt%average = t%average - oldt%average
302  dels%minimum = s%minimum - olds%minimum ; dels%maximum = s%maximum - olds%maximum
303  dels%average = s%average - olds%average
304  write(lmsg(1:80),'(2(a,es12.4))') 'Mean thickness =', vol/area,' frac. delta=',dv/vol
305  call mom_mesg(lmsg//trim(mesg))
306  if (do_ts) then
307  write(lmsg(1:80),'(a,3es12.4)') 'Temp min/mean/max =',t%minimum,t%average,t%maximum
308  call mom_mesg(lmsg//trim(mesg))
309  write(lmsg(1:80),'(a,3es12.4)') 'delT min/mean/max =',delt%minimum,delt%average,delt%maximum
310  call mom_mesg(lmsg//trim(mesg))
311  write(lmsg(1:80),'(a,3es12.4)') 'Salt min/mean/max =',s%minimum,s%average,s%maximum
312  call mom_mesg(lmsg//trim(mesg))
313  write(lmsg(1:80),'(a,3es12.4)') 'delS min/mean/max =',dels%minimum,dels%average,dels%maximum
314  call mom_mesg(lmsg//trim(mesg))
315  endif
316  else
317  write(lmsg(1:80),'(a,es12.4)') 'Mean thickness =', vol/area
318  call mom_mesg(lmsg//trim(mesg))
319  if (do_ts) then
320  write(lmsg(1:80),'(a,3es12.4)') 'Temp min/mean/max =', t%minimum, t%average, t%maximum
321  call mom_mesg(lmsg//trim(mesg))
322  write(lmsg(1:80),'(a,3es12.4)') 'Salt min/mean/max =', s%minimum, s%average, s%maximum
323  call mom_mesg(lmsg//trim(mesg))
324  endif
325  endif
326  endif
327  firstcall = .false. ; oldvol = vol
328  oldt%minimum = t%minimum ; oldt%maximum = t%maximum ; oldt%average = t%average
329  olds%minimum = s%minimum ; olds%maximum = s%maximum ; olds%average = s%average
330 
331  if (do_ts .and. t%minimum<-5.0) then
332  do j = js, je ; do i = is, ie
333  if (minval(temp(i,j,:)) == t%minimum) then
334  write(0,'(a,2f12.5)') 'x,y=', g%geoLonT(i,j), g%geoLatT(i,j)
335  write(0,'(a3,3a12)') 'k','h','Temp','Salt'
336  do k = 1, nz
337  write(0,'(i3,3es12.4)') k, h(i,j,k), temp(i,j,k), salt(i,j,k)
338  enddo
339  stop 'Extremum detected'
340  endif
341  enddo ; enddo
342  endif
343 
344  if (h_minimum<0.0) then
345  do j = js, je ; do i = is, ie
346  if (minval(h(i,j,:)) == h_minimum) then
347  write(0,'(a,2f12.5)') 'x,y=',g%geoLonT(i,j),g%geoLatT(i,j)
348  write(0,'(a3,3a12)') 'k','h','Temp','Salt'
349  do k = 1, nz
350  write(0,'(i3,3es12.4)') k, h(i,j,k), temp(i,j,k), salt(i,j,k)
351  enddo
352  stop 'Negative thickness detected'
353  endif
354  enddo ; enddo
355  endif
356 
357 end subroutine mom_state_stats
358 
359 end module mom_checksum_packages
mom_checksum_packages::stats
A type for storing statistica about a variable.
Definition: MOM_checksum_packages.F90:31
mom_variables::surface
Pointers to various fields which may be used describe the surface state of MOM, and which will be ret...
Definition: MOM_variables.F90:38
mom_checksum_packages::mom_state_chksum
Write out checksums of the MOM6 state variables.
Definition: MOM_checksum_packages.F90:23
mom_verticalgrid
Provides a transparent vertical ocean grid type and supporting routines.
Definition: MOM_verticalGrid.F90:2
mom_variables::thermo_var_ptrs
Pointers to an assortment of thermodynamic fields that may be available, including potential temperat...
Definition: MOM_variables.F90:78
mom_checksum_packages::mom_state_chksum_5arg
subroutine mom_state_chksum_5arg(mesg, u, v, h, uh, vh, G, GV, US, haloshift, symmetric, vel_scale)
Write out chksums for the model's basic state variables, including transports.
Definition: MOM_checksum_packages.F90:43
mom_error_handler::mom_mesg
subroutine, public mom_mesg(message, verb, all_print)
This provides a convenient interface for writing an informative comment.
Definition: MOM_error_handler.F90:53
mom_checksum_packages::mom_accel_chksum
subroutine, public mom_accel_chksum(mesg, CAu, CAv, PFu, PFv, diffu, diffv, G, GV, US, pbce, u_accel_bt, v_accel_bt, symmetric)
Write out chksums for the model's accelerations.
Definition: MOM_checksum_packages.F90:175
mom_checksum_packages
Provides routines that do checksums of groups of MOM variables.
Definition: MOM_checksum_packages.F90:2
mom_unit_scaling::unit_scale_type
Describes various unit conversion factors.
Definition: MOM_unit_scaling.F90:14
mom_verticalgrid::verticalgrid_type
Describes the vertical ocean grid, including unit conversion factors.
Definition: MOM_verticalGrid.F90:24
mom_domains
Describes the decomposed MOM domain and has routines for communications across PEs.
Definition: MOM_domains.F90:2
mom_variables
Provides transparent structures with groups of MOM6 variables and supporting routines.
Definition: MOM_variables.F90:2
mom_checksum_packages::mom_state_stats
subroutine, public mom_state_stats(mesg, u, v, h, Temp, Salt, G, GV, US, allowChange, permitDiminishing)
Monitor and write out statistics for the model's state variables.
Definition: MOM_checksum_packages.F90:233
mom_checksum_packages::mom_state_chksum_3arg
subroutine mom_state_chksum_3arg(mesg, u, v, h, G, GV, US, haloshift, symmetric)
Write out chksums for the model's basic state variables.
Definition: MOM_checksum_packages.F90:87
mom_checksum_packages::mom_surface_chksum
subroutine, public mom_surface_chksum(mesg, sfc, G, haloshift, symmetric)
Write out chksums for the ocean surface variables.
Definition: MOM_checksum_packages.F90:144
mom_grid
Provides the ocean grid type.
Definition: MOM_grid.F90:2
mom_unit_scaling
Provides a transparent unit rescaling type to facilitate dimensional consistency testing.
Definition: MOM_unit_scaling.F90:2
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_debugging
Provides checksumming functions for debugging.
Definition: MOM_debugging.F90:7
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_checksum_packages::mom_thermo_chksum
subroutine, public mom_thermo_chksum(mesg, tv, G, US, haloshift)
Write out chksums for the model's thermodynamic state variables.
Definition: MOM_checksum_packages.F90:121