MOM6
mom_diagnostics Module Reference

Detailed Description

Calculates any requested diagnostic quantities that are not calculated in the various subroutines. Diagnostic quantities are requested by allocating them memory.

Data Types

type  diagnostics_cs
 The control structure for the MOM_diagnostics module. More...
 
type  surface_diag_ids
 A structure with diagnostic IDs of the surface and integrated variables. More...
 
type  transport_diag_ids
 A structure with diagnostic IDs of mass transport related diagnostics. More...
 

Functions/Subroutines

subroutine, public calculate_diagnostic_fields (u, v, h, uh, vh, tv, ADp, CDp, p_surf, dt, diag_pre_sync, G, GV, US, CS, eta_bt)
 Diagnostics not more naturally calculated elsewhere are computed here. More...
 
subroutine find_weights (Rlist, R_in, k, nz, wt, wt_p)
 This subroutine finds the location of R_in in an increasing ordered list, Rlist, returning as k the element such that Rlist(k) <= R_in < Rlist(k+1), and where wt and wt_p are the linear weights that should be assigned to elements k and k+1. More...
 
subroutine calculate_vertical_integrals (h, tv, p_surf, G, GV, US, CS)
 This subroutine calculates vertical integrals of several tracers, along with the mass-weight of these tracers, the total column mass, and the carefully calculated column height. More...
 
subroutine calculate_energy_diagnostics (u, v, h, uh, vh, ADp, CDp, G, GV, US, CS)
 This subroutine calculates terms in the mechanical energy budget. More...
 
subroutine, public register_time_deriv (lb, f_ptr, deriv_ptr, CS)
 This subroutine registers fields to calculate a diagnostic time derivative. More...
 
subroutine calculate_derivs (dt, G, CS)
 This subroutine calculates all registered time derivatives. More...
 
subroutine, public post_surface_dyn_diags (IDs, G, diag, sfc_state, ssh)
 This routine posts diagnostics of various dynamic ocean surface quantities, including velocities, speed and sea surface height, at the time the ocean state is reported back to the caller. More...
 
subroutine, public post_surface_thermo_diags (IDs, G, GV, US, diag, dt_int, sfc_state, tv, ssh, ssh_ibc)
 This routine posts diagnostics of various ocean surface and integrated quantities at the time the ocean state is reported back to the caller. More...
 
subroutine, public post_transport_diagnostics (G, GV, US, uhtr, vhtr, h, IDs, diag_pre_dyn, diag, dt_trans, Reg)
 This routine posts diagnostics of the transports, including the subgridscale contributions. More...
 
subroutine, public mom_diagnostics_init (MIS, ADp, CDp, Time, G, GV, US, param_file, diag, CS, tv)
 This subroutine registers various diagnostics and allocates space for fields that other diagnostis depend upon. More...
 
subroutine, public register_surface_diags (Time, G, US, IDs, diag, tv)
 Register diagnostics of the surface state and integrated quantities. More...
 
subroutine, public register_transport_diags (Time, G, GV, US, IDs, diag)
 Register certain diagnostics related to transports. More...
 
subroutine, public write_static_fields (G, GV, US, tv, diag)
 Offers the static fields in the ocean grid type for output via the diag_manager. More...
 
subroutine set_dependent_diagnostics (MIS, ADp, CDp, G, CS)
 This subroutine sets up diagnostics upon which other diagnostics depend. More...
 
subroutine, public mom_diagnostics_end (CS, ADp)
 Deallocate memory associated with the diagnostics module. More...
 

Function/Subroutine Documentation

◆ calculate_derivs()

subroutine mom_diagnostics::calculate_derivs ( real, intent(in)  dt,
type(ocean_grid_type), intent(inout)  G,
type(diagnostics_cs), intent(inout)  CS 
)
private

This subroutine calculates all registered time derivatives.

Parameters
[in]dtThe time interval over which differences occur [T ~> s].
[in,out]gThe ocean's grid structure.
[in,out]csControl structure returned by previous call to diagnostics_init.

Definition at line 1124 of file MOM_diagnostics.F90.

1124  real, intent(in) :: dt !< The time interval over which differences occur [T ~> s].
1125  type(ocean_grid_type), intent(inout) :: G !< The ocean's grid structure.
1126  type(diagnostics_CS), intent(inout) :: CS !< Control structure returned by previous call to
1127  !! diagnostics_init.
1128 
1129 ! This subroutine calculates all registered time derivatives.
1130  real :: Idt ! The inverse timestep [T-1 ~> s-1]
1131  integer :: i, j, k, m
1132 
1133  if (dt > 0.0) then ; idt = 1.0/dt
1134  else ; return ; endif
1135 
1136  ! Because the field is unknown, its grid index bounds are also unknown.
1137  ! Additionally, two of the fields (dudt, dvdt) require calculation of spatial
1138  ! derivatives when computing d(KE)/dt. This raises issues in non-symmetric
1139  ! mode, where the symmetric boundaries (west, south) may not be updated.
1140 
1141  ! For this reason, we explicitly loop from isc-1:iec and jsc-1:jec, in order
1142  ! to force boundary value updates, even though it may not be strictly valid
1143  ! for all fields. Note this assumes a halo, and that it has been updated.
1144 
1145  do m=1,cs%num_time_deriv
1146  do k=1,cs%nlay(m) ; do j=g%jsc-1,g%jec ; do i=g%isc-1,g%iec
1147  cs%deriv(m)%p(i,j,k) = (cs%var_ptr(m)%p(i,j,k) - cs%prev_val(m)%p(i,j,k)) * idt
1148  cs%prev_val(m)%p(i,j,k) = cs%var_ptr(m)%p(i,j,k)
1149  enddo ; enddo ; enddo
1150  enddo
1151 

Referenced by calculate_diagnostic_fields().

Here is the caller graph for this function:

◆ calculate_diagnostic_fields()

subroutine, public mom_diagnostics::calculate_diagnostic_fields ( real, dimension( g %isdb: g %iedb, g %jsd: g %jed, g %ke), intent(in)  u,
real, dimension( g %isd: g %ied, g %jsdb: g %jedb, g %ke), intent(in)  v,
real, dimension( g %isd: g %ied, g %jsd: g %jed, g %ke), intent(in)  h,
real, dimension( g %isdb: g %iedb, g %jsd: g %jed, g %ke), intent(in)  uh,
real, dimension( g %isd: g %ied, g %jsdb: g %jedb, g %ke), intent(in)  vh,
type(thermo_var_ptrs), intent(in)  tv,
type(accel_diag_ptrs), intent(in)  ADp,
type(cont_diag_ptrs), intent(in)  CDp,
real, dimension(:,:), pointer  p_surf,
real, intent(in)  dt,
type(diag_grid_storage), intent(in)  diag_pre_sync,
type(ocean_grid_type), intent(inout)  G,
type(verticalgrid_type), intent(in)  GV,
type(unit_scale_type), intent(in)  US,
type(diagnostics_cs), intent(inout)  CS,
real, dimension( g %isd: g %ied, g %jsd: g %jed), intent(in), optional  eta_bt 
)

Diagnostics not more naturally calculated elsewhere are computed here.

Parameters
[in,out]gThe ocean's grid structure.
[in]gvThe ocean's vertical grid structure.
[in]usA dimensional unit scaling type
[in]uThe zonal velocity [L T-1 ~> m s-1].
[in]vThe meridional velocity [L T-1 ~> m s-1].
[in]hLayer thicknesses [H ~> m or kg m-2].
[in]uhTransport through zonal faces = u*h*dy,
[in]vhTransport through meridional faces = v*h*dx,
[in]tvA structure pointing to various thermodynamic variables.
[in]adpstructure with pointers to accelerations in momentum equation.
[in]cdpstructure with pointers to terms in continuity equation.
p_surfA pointer to the surface pressure [Pa]. If p_surf is not associated, it is the same as setting the surface pressure to 0.
[in]dtThe time difference since the last call to this subroutine [T ~> s].
[in]diag_pre_syncTarget grids from previous timestep
[in,out]csControl structure returned by a previous call to diagnostics_init.
[in]eta_btAn optional barotropic

Definition at line 188 of file MOM_diagnostics.F90.

188  type(ocean_grid_type), intent(inout) :: G !< The ocean's grid structure.
189  type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure.
190  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
191  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)), &
192  intent(in) :: u !< The zonal velocity [L T-1 ~> m s-1].
193  real, dimension(SZI_(G),SZJB_(G),SZK_(G)), &
194  intent(in) :: v !< The meridional velocity [L T-1 ~> m s-1].
195  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
196  intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2].
197  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)), &
198  intent(in) :: uh !< Transport through zonal faces = u*h*dy,
199  !! [H L2 T-1 ~> m3 s-1 or kg s-1].
200  real, dimension(SZI_(G),SZJB_(G),SZK_(G)), &
201  intent(in) :: vh !< Transport through meridional faces = v*h*dx,
202  !! [H L2 T-1 ~> m3 s-1 or kg s-1].
203  type(thermo_var_ptrs), intent(in) :: tv !< A structure pointing to various
204  !! thermodynamic variables.
205  type(accel_diag_ptrs), intent(in) :: ADp !< structure with pointers to
206  !! accelerations in momentum equation.
207  type(cont_diag_ptrs), intent(in) :: CDp !< structure with pointers to
208  !! terms in continuity equation.
209  real, dimension(:,:), pointer :: p_surf !< A pointer to the surface pressure [Pa].
210  !! If p_surf is not associated, it is the same
211  !! as setting the surface pressure to 0.
212  real, intent(in) :: dt !< The time difference since the last
213  !! call to this subroutine [T ~> s].
214  type(diag_grid_storage), intent(in) :: diag_pre_sync !< Target grids from previous timestep
215  type(diagnostics_CS), intent(inout) :: CS !< Control structure returned by a
216  !! previous call to diagnostics_init.
217  real, dimension(SZI_(G),SZJ_(G)), &
218  optional, intent(in) :: eta_bt !< An optional barotropic
219  !! variable that gives the "correct" free surface height (Boussinesq) or total water column
220  !! mass per unit area (non-Boussinesq). This is used to dilate the layer thicknesses when
221  !! calculating interface heights [H ~> m or kg m-2].
222 
223  ! Local variables
224  integer i, j, k, is, ie, js, je, Isq, Ieq, Jsq, Jeq, nz, nkmb
225 
226  real :: Rcv(SZI_(G),SZJ_(G),SZK_(G)) ! Coordinate variable potential density [R ~> kg m-3].
227  real :: work_3d(SZI_(G),SZJ_(G),SZK_(G)) ! A 3-d temporary work array.
228  real :: work_2d(SZI_(G),SZJ_(G)) ! A 2-d temporary work array.
229 
230  ! tmp array for surface properties
231  real :: surface_field(SZI_(G),SZJ_(G))
232  real :: pressure_1d(SZI_(G)) ! Temporary array for pressure when calling EOS
233  real :: wt, wt_p
234 
235  real :: f2_h ! Squared Coriolis parameter at to h-points [T-2 ~> s-2]
236  real :: mag_beta ! Magnitude of the gradient of f [T-1 L-1 ~> s-1 m-1]
237  real :: absurdly_small_freq2 ! Srequency squared used to avoid division by 0 [T-2 ~> s-2]
238 
239  integer :: k_list
240 
241  real, dimension(SZK_(G)) :: temp_layer_ave, salt_layer_ave
242  real :: thetaoga, soga, masso, tosga, sosga
243 
244  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
245  isq = g%IscB ; ieq = g%IecB ; jsq = g%JscB ; jeq = g%JecB
246  nz = g%ke ; nkmb = gv%nk_rho_varies
247 
248  ! This value is roughly (pi / (the age of the universe) )^2.
249  absurdly_small_freq2 = 1e-34*us%T_to_s**2
250 
251  if (loc(cs)==0) call mom_error(fatal, &
252  "calculate_diagnostic_fields: Module must be initialized before used.")
253 
254  call calculate_derivs(dt, g, cs)
255 
256  if (dt > 0.0) then
257  call diag_save_grids(cs%diag)
258  call diag_copy_storage_to_diag(cs%diag, diag_pre_sync)
259 
260  if (cs%id_h_pre_sync > 0) &
261  call post_data(cs%id_h_pre_sync, diag_pre_sync%h_state, cs%diag, alt_h = diag_pre_sync%h_state)
262 
263  if (cs%id_du_dt>0) call post_data(cs%id_du_dt, cs%du_dt, cs%diag, alt_h = diag_pre_sync%h_state)
264 
265  if (cs%id_dv_dt>0) call post_data(cs%id_dv_dt, cs%dv_dt, cs%diag, alt_h = diag_pre_sync%h_state)
266 
267  if (cs%id_dh_dt>0) call post_data(cs%id_dh_dt, cs%dh_dt, cs%diag, alt_h = diag_pre_sync%h_state)
268 
269  call diag_restore_grids(cs%diag)
270 
271  call calculate_energy_diagnostics(u, v, h, uh, vh, adp, cdp, g, gv, us, cs)
272  endif
273 
274  ! smg: is the following robust to ALE? It seems a bit opaque.
275  ! If the model is NOT in isopycnal mode then nkmb=0. But we need all the
276  ! following diagnostics to treat all layers as variable density, so we set
277  ! nkmb = nz, on the expectation that loops nkmb+1,nz will not iterate.
278  ! This behavior is ANSI F77 but some compiler options can force at least
279  ! one iteration that would break the following one-line workaround!
280  if (nkmb==0 .and. nz > 1) nkmb = nz
281 
282  if (cs%id_u > 0) call post_data(cs%id_u, u, cs%diag)
283 
284  if (cs%id_v > 0) call post_data(cs%id_v, v, cs%diag)
285 
286  if (cs%id_h > 0) call post_data(cs%id_h, h, cs%diag)
287 
288  if (associated(cs%e)) then
289  call find_eta(h, tv, g, gv, us, cs%e, eta_bt)
290  if (cs%id_e > 0) call post_data(cs%id_e, cs%e, cs%diag)
291  endif
292 
293  if (associated(cs%e_D)) then
294  if (associated(cs%e)) then
295  do k=1,nz+1 ; do j=js,je ; do i=is,ie
296  cs%e_D(i,j,k) = cs%e(i,j,k) + g%bathyT(i,j)
297  enddo ; enddo ; enddo
298  else
299  call find_eta(h, tv, g, gv, us, cs%e_D, eta_bt)
300  do k=1,nz+1 ; do j=js,je ; do i=is,ie
301  cs%e_D(i,j,k) = cs%e_D(i,j,k) + g%bathyT(i,j)
302  enddo ; enddo ; enddo
303  endif
304 
305  if (cs%id_e_D > 0) call post_data(cs%id_e_D, cs%e_D, cs%diag)
306  endif
307 
308  ! mass per area of grid cell (for Bouss, use Rho0)
309  if (cs%id_masscello > 0) then
310  do k=1,nz; do j=js,je ; do i=is,ie
311  work_3d(i,j,k) = gv%H_to_kg_m2*h(i,j,k)
312  enddo ; enddo ; enddo
313  call post_data(cs%id_masscello, work_3d, cs%diag)
314  endif
315 
316  ! mass of liquid ocean (for Bouss, use Rho0). The reproducing sum requires the use of MKS units.
317  if (cs%id_masso > 0) then
318  work_2d(:,:) = 0.0
319  do k=1,nz ; do j=js,je ; do i=is,ie
320  work_2d(i,j) = work_2d(i,j) + (gv%H_to_kg_m2*h(i,j,k)) * us%L_to_m**2*g%areaT(i,j)
321  enddo ; enddo ; enddo
322  masso = reproducing_sum(work_2d)
323  call post_data(cs%id_masso, masso, cs%diag)
324  endif
325 
326  ! diagnose thickness/volumes of grid cells [m]
327  if (cs%id_thkcello>0 .or. cs%id_volcello>0) then
328  if (gv%Boussinesq) then ! thkcello = h for Boussinesq
329  if (cs%id_thkcello > 0) then ; if (gv%H_to_m == 1.0) then
330  call post_data(cs%id_thkcello, h, cs%diag)
331  else
332  do k=1,nz; do j=js,je ; do i=is,ie
333  work_3d(i,j,k) = gv%H_to_m*h(i,j,k)
334  enddo ; enddo ; enddo
335  call post_data(cs%id_thkcello, work_3d, cs%diag)
336  endif ; endif
337  if (cs%id_volcello > 0) then ! volcello = h*area for Boussinesq
338  do k=1,nz; do j=js,je ; do i=is,ie
339  work_3d(i,j,k) = ( gv%H_to_m*h(i,j,k) ) * us%L_to_m**2*g%areaT(i,j)
340  enddo ; enddo ; enddo
341  call post_data(cs%id_volcello, work_3d, cs%diag)
342  endif
343  else ! thkcello = dp/(rho*g) for non-Boussinesq
344  do j=js,je
345  if (associated(p_surf)) then ! Pressure loading at top of surface layer [Pa]
346  do i=is,ie
347  pressure_1d(i) = p_surf(i,j)
348  enddo
349  else
350  do i=is,ie
351  pressure_1d(i) = 0.0
352  enddo
353  endif
354  do k=1,nz ! Integrate vertically downward for pressure
355  do i=is,ie ! Pressure for EOS at the layer center [Pa]
356  pressure_1d(i) = pressure_1d(i) + 0.5*gv%H_to_Pa*h(i,j,k)
357  enddo
358  ! Store in-situ density [kg m-3] in work_3d
359  call calculate_density(tv%T(:,j,k),tv%S(:,j,k), pressure_1d, &
360  work_3d(:,j,k), is, ie-is+1, tv%eqn_of_state)
361  do i=is,ie ! Cell thickness = dz = dp/(g*rho) (meter); store in work_3d
362  work_3d(i,j,k) = (gv%H_to_kg_m2*h(i,j,k)) / work_3d(i,j,k)
363  enddo
364  do i=is,ie ! Pressure for EOS at the bottom interface [Pa]
365  pressure_1d(i) = pressure_1d(i) + 0.5*gv%H_to_Pa*h(i,j,k)
366  enddo
367  enddo ! k
368  enddo ! j
369  if (cs%id_thkcello > 0) call post_data(cs%id_thkcello, work_3d, cs%diag)
370  if (cs%id_volcello > 0) then
371  do k=1,nz; do j=js,je ; do i=is,ie ! volcello = dp/(rho*g)*area for non-Boussinesq
372  work_3d(i,j,k) = us%L_to_m**2*g%areaT(i,j) * work_3d(i,j,k)
373  enddo ; enddo ; enddo
374  call post_data(cs%id_volcello, work_3d, cs%diag)
375  endif
376  endif
377  endif
378 
379  ! Calculate additional, potentially derived temperature diagnostics
380  if (tv%T_is_conT) then
381  ! Internal T&S variables are conservative temperature & absolute salinity,
382  ! so they need to converted to potential temperature and practical salinity
383  ! for some diagnostics using TEOS-10 function calls.
384  if ((cs%id_Tpot > 0) .or. (cs%id_tob > 0)) then
385  do k=1,nz ; do j=js,je ; do i=is,ie
386  work_3d(i,j,k) = gsw_pt_from_ct(tv%S(i,j,k),tv%T(i,j,k))
387  enddo ; enddo ; enddo
388  if (cs%id_Tpot > 0) call post_data(cs%id_Tpot, work_3d, cs%diag)
389  if (cs%id_tob > 0) call post_data(cs%id_tob, work_3d(:,:,nz), cs%diag, mask=g%mask2dT)
390  endif
391  else
392  ! Internal T&S variables are potential temperature & practical salinity
393  if (cs%id_tob > 0) call post_data(cs%id_tob, tv%T(:,:,nz), cs%diag, mask=g%mask2dT)
394  endif
395 
396  ! Calculate additional, potentially derived salinity diagnostics
397  if (tv%S_is_absS) then
398  ! Internal T&S variables are conservative temperature & absolute salinity,
399  ! so they need to converted to potential temperature and practical salinity
400  ! for some diagnostics using TEOS-10 function calls.
401  if ((cs%id_Sprac > 0) .or. (cs%id_sob > 0)) then
402  do k=1,nz ; do j=js,je ; do i=is,ie
403  work_3d(i,j,k) = gsw_sp_from_sr(tv%S(i,j,k))
404  enddo ; enddo ; enddo
405  if (cs%id_Sprac > 0) call post_data(cs%id_Sprac, work_3d, cs%diag)
406  if (cs%id_sob > 0) call post_data(cs%id_sob, work_3d(:,:,nz), cs%diag, mask=g%mask2dT)
407  endif
408  else
409  ! Internal T&S variables are potential temperature & practical salinity
410  if (cs%id_sob > 0) call post_data(cs%id_sob, tv%S(:,:,nz), cs%diag, mask=g%mask2dT)
411  endif
412 
413  ! volume mean potential temperature
414  if (cs%id_thetaoga>0) then
415  thetaoga = global_volume_mean(tv%T, h, g, gv)
416  call post_data(cs%id_thetaoga, thetaoga, cs%diag)
417  endif
418 
419  ! area mean SST
420  if (cs%id_tosga > 0) then
421  do j=js,je ; do i=is,ie
422  surface_field(i,j) = tv%T(i,j,1)
423  enddo ; enddo
424  tosga = global_area_mean(surface_field, g)
425  call post_data(cs%id_tosga, tosga, cs%diag)
426  endif
427 
428  ! volume mean salinity
429  if (cs%id_soga>0) then
430  soga = global_volume_mean(tv%S, h, g, gv)
431  call post_data(cs%id_soga, soga, cs%diag)
432  endif
433 
434  ! area mean SSS
435  if (cs%id_sosga > 0) then
436  do j=js,je ; do i=is,ie
437  surface_field(i,j) = tv%S(i,j,1)
438  enddo ; enddo
439  sosga = global_area_mean(surface_field, g)
440  call post_data(cs%id_sosga, sosga, cs%diag)
441  endif
442 
443  ! layer mean potential temperature
444  if (cs%id_temp_layer_ave>0) then
445  temp_layer_ave = global_layer_mean(tv%T, h, g, gv)
446  call post_data(cs%id_temp_layer_ave, temp_layer_ave, cs%diag)
447  endif
448 
449  ! layer mean salinity
450  if (cs%id_salt_layer_ave>0) then
451  salt_layer_ave = global_layer_mean(tv%S, h, g, gv)
452  call post_data(cs%id_salt_layer_ave, salt_layer_ave, cs%diag)
453  endif
454 
455  call calculate_vertical_integrals(h, tv, p_surf, g, gv, us, cs)
456 
457  if ((cs%id_Rml > 0) .or. (cs%id_Rcv > 0) .or. associated(cs%h_Rlay) .or. &
458  associated(cs%uh_Rlay) .or. associated(cs%vh_Rlay) .or. &
459  associated(cs%uhGM_Rlay) .or. associated(cs%vhGM_Rlay)) then
460 
461  if (associated(tv%eqn_of_state)) then
462  pressure_1d(:) = tv%P_Ref
463  !$OMP parallel do default(shared)
464  do k=1,nz ; do j=js-1,je+1
465  call calculate_density(tv%T(:,j,k), tv%S(:,j,k), pressure_1d, &
466  rcv(:,j,k), is-1, ie-is+3, tv%eqn_of_state, scale=us%kg_m3_to_R)
467  enddo ; enddo
468  else ! Rcv should not be used much in this case, so fill in sensible values.
469  do k=1,nz ; do j=js-1,je+1 ; do i=is-1,ie+1
470  rcv(i,j,k) = gv%Rlay(k)
471  enddo ; enddo ; enddo
472  endif
473  if (cs%id_Rml > 0) call post_data(cs%id_Rml, rcv, cs%diag)
474  if (cs%id_Rcv > 0) call post_data(cs%id_Rcv, rcv, cs%diag)
475 
476  if (associated(cs%h_Rlay)) then
477  k_list = nz/2
478 !$OMP parallel do default(none) shared(is,ie,js,je,nz,nkmb,CS,Rcv,h,GV) &
479 !$OMP private(wt,wt_p) firstprivate(k_list)
480  do j=js,je
481  do k=1,nkmb ; do i=is,ie
482  cs%h_Rlay(i,j,k) = 0.0
483  enddo ; enddo
484  do k=nkmb+1,nz ; do i=is,ie
485  cs%h_Rlay(i,j,k) = h(i,j,k)
486  enddo ; enddo
487  do k=1,nkmb ; do i=is,ie
488  call find_weights(gv%Rlay, rcv(i,j,k), k_list, nz, wt, wt_p)
489  cs%h_Rlay(i,j,k_list) = cs%h_Rlay(i,j,k_list) + h(i,j,k)*wt
490  cs%h_Rlay(i,j,k_list+1) = cs%h_Rlay(i,j,k_list+1) + h(i,j,k)*wt_p
491  enddo ; enddo
492  enddo
493 
494  if (cs%id_h_Rlay > 0) call post_data(cs%id_h_Rlay, cs%h_Rlay, cs%diag)
495  endif
496 
497  if (associated(cs%uh_Rlay)) then
498  k_list = nz/2
499 !$OMP parallel do default(none) shared(Isq,Ieq,js,je,nz,nkmb,Rcv,CS,GV,uh) &
500 !$OMP private(wt,wt_p) firstprivate(k_list)
501  do j=js,je
502  do k=1,nkmb ; do i=isq,ieq
503  cs%uh_Rlay(i,j,k) = 0.0
504  enddo ; enddo
505  do k=nkmb+1,nz ; do i=isq,ieq
506  cs%uh_Rlay(i,j,k) = uh(i,j,k)
507  enddo ; enddo
508  k_list = nz/2
509  do k=1,nkmb ; do i=isq,ieq
510  call find_weights(gv%Rlay, 0.5*(rcv(i,j,k)+rcv(i+1,j,k)), k_list, nz, wt, wt_p)
511  cs%uh_Rlay(i,j,k_list) = cs%uh_Rlay(i,j,k_list) + uh(i,j,k)*wt
512  cs%uh_Rlay(i,j,k_list+1) = cs%uh_Rlay(i,j,k_list+1) + uh(i,j,k)*wt_p
513  enddo ; enddo
514  enddo
515 
516  if (cs%id_uh_Rlay > 0) call post_data(cs%id_uh_Rlay, cs%uh_Rlay, cs%diag)
517  endif
518 
519  if (associated(cs%vh_Rlay)) then
520  k_list = nz/2
521 !$OMP parallel do default(none) shared(Jsq,Jeq,is,ie,nz,nkmb,Rcv,CS,GV,vh) &
522 !$OMP private(wt,wt_p) firstprivate(k_list)
523  do j=jsq,jeq
524  do k=1,nkmb ; do i=is,ie
525  cs%vh_Rlay(i,j,k) = 0.0
526  enddo ; enddo
527  do k=nkmb+1,nz ; do i=is,ie
528  cs%vh_Rlay(i,j,k) = vh(i,j,k)
529  enddo ; enddo
530  do k=1,nkmb ; do i=is,ie
531  call find_weights(gv%Rlay, 0.5*(rcv(i,j,k)+rcv(i,j+1,k)), k_list, nz, wt, wt_p)
532  cs%vh_Rlay(i,j,k_list) = cs%vh_Rlay(i,j,k_list) + vh(i,j,k)*wt
533  cs%vh_Rlay(i,j,k_list+1) = cs%vh_Rlay(i,j,k_list+1) + vh(i,j,k)*wt_p
534  enddo ; enddo
535  enddo
536 
537  if (cs%id_vh_Rlay > 0) call post_data(cs%id_vh_Rlay, cs%vh_Rlay, cs%diag)
538  endif
539 
540  if (associated(cs%uhGM_Rlay) .and. associated(cdp%uhGM)) then
541  k_list = nz/2
542 !$OMP parallel do default(none) shared(Isq,Ieq,js,je,nz,nkmb,Rcv,CDP,CS,GV) &
543 !$OMP private(wt,wt_p) firstprivate(k_list)
544  do j=js,je
545  do k=1,nkmb ; do i=isq,ieq
546  cs%uhGM_Rlay(i,j,k) = 0.0
547  enddo ; enddo
548  do k=nkmb+1,nz ; do i=isq,ieq
549  cs%uhGM_Rlay(i,j,k) = cdp%uhGM(i,j,k)
550  enddo ; enddo
551  do k=1,nkmb ; do i=isq,ieq
552  call find_weights(gv%Rlay, 0.5*(rcv(i,j,k)+rcv(i+1,j,k)), k_list, nz, wt, wt_p)
553  cs%uhGM_Rlay(i,j,k_list) = cs%uhGM_Rlay(i,j,k_list) + cdp%uhGM(i,j,k)*wt
554  cs%uhGM_Rlay(i,j,k_list+1) = cs%uhGM_Rlay(i,j,k_list+1) + cdp%uhGM(i,j,k)*wt_p
555  enddo ; enddo
556  enddo
557 
558  if (cs%id_uhGM_Rlay > 0) call post_data(cs%id_uhGM_Rlay, cs%uhGM_Rlay, cs%diag)
559  endif
560 
561  if (associated(cs%vhGM_Rlay) .and. associated(cdp%vhGM)) then
562  k_list = nz/2
563 !$OMP parallel do default(none) shared(is,ie,Jsq,Jeq,nz,nkmb,CS,CDp,Rcv,GV) &
564 !$OMP private(wt,wt_p) firstprivate(k_list)
565  do j=jsq,jeq
566  do k=1,nkmb ; do i=is,ie
567  cs%vhGM_Rlay(i,j,k) = 0.0
568  enddo ; enddo
569  do k=nkmb+1,nz ; do i=is,ie
570  cs%vhGM_Rlay(i,j,k) = cdp%vhGM(i,j,k)
571  enddo ; enddo
572  do k=1,nkmb ; do i=is,ie
573  call find_weights(gv%Rlay, 0.5*(rcv(i,j,k)+rcv(i,j+1,k)), k_list, nz, wt, wt_p)
574  cs%vhGM_Rlay(i,j,k_list) = cs%vhGM_Rlay(i,j,k_list) + cdp%vhGM(i,j,k)*wt
575  cs%vhGM_Rlay(i,j,k_list+1) = cs%vhGM_Rlay(i,j,k_list+1) + cdp%vhGM(i,j,k)*wt_p
576  enddo ; enddo
577  enddo
578 
579  if (cs%id_vhGM_Rlay > 0) call post_data(cs%id_vhGM_Rlay, cs%vhGM_Rlay, cs%diag)
580  endif
581  endif
582 
583  if (associated(tv%eqn_of_state)) then
584  if (cs%id_rhopot0 > 0) then
585  pressure_1d(:) = 0.
586 !$OMP parallel do default(none) shared(tv,Rcv,is,ie,js,je,nz,pressure_1d)
587  do k=1,nz ; do j=js,je
588  call calculate_density(tv%T(:,j,k),tv%S(:,j,k),pressure_1d, &
589  rcv(:,j,k),is,ie-is+1, tv%eqn_of_state)
590  enddo ; enddo
591  if (cs%id_rhopot0 > 0) call post_data(cs%id_rhopot0, rcv, cs%diag)
592  endif
593  if (cs%id_rhopot2 > 0) then
594  pressure_1d(:) = 2.e7 ! 2000 dbars
595 !$OMP parallel do default(none) shared(tv,Rcv,is,ie,js,je,nz,pressure_1d)
596  do k=1,nz ; do j=js,je
597  call calculate_density(tv%T(:,j,k),tv%S(:,j,k),pressure_1d, &
598  rcv(:,j,k),is,ie-is+1, tv%eqn_of_state)
599  enddo ; enddo
600  if (cs%id_rhopot2 > 0) call post_data(cs%id_rhopot2, rcv, cs%diag)
601  endif
602  if (cs%id_rhoinsitu > 0) then
603 !$OMP parallel do default(none) shared(tv,Rcv,is,ie,js,je,nz,pressure_1d,h,GV)
604  do j=js,je
605  pressure_1d(:) = 0. ! Start at p=0 Pa at surface
606  do k=1,nz
607  pressure_1d(:) = pressure_1d(:) + 0.5 * h(:,j,k) * gv%H_to_Pa ! Pressure in middle of layer k
608  call calculate_density(tv%T(:,j,k),tv%S(:,j,k),pressure_1d, &
609  rcv(:,j,k),is,ie-is+1, tv%eqn_of_state)
610  pressure_1d(:) = pressure_1d(:) + 0.5 * h(:,j,k) * gv%H_to_Pa ! Pressure at bottom of layer k
611  enddo
612  enddo
613  if (cs%id_rhoinsitu > 0) call post_data(cs%id_rhoinsitu, rcv, cs%diag)
614  endif
615  endif
616 
617  if ((cs%id_cg1>0) .or. (cs%id_Rd1>0) .or. (cs%id_cfl_cg1>0) .or. &
618  (cs%id_cfl_cg1_x>0) .or. (cs%id_cfl_cg1_y>0)) then
619  call wave_speed(h, tv, g, gv, us, cs%cg1, cs%wave_speed_CSp)
620  if (cs%id_cg1>0) call post_data(cs%id_cg1, cs%cg1, cs%diag)
621  if (cs%id_Rd1>0) then
622  !$OMP parallel do default(shared) private(f2_h,mag_beta)
623  do j=js,je ; do i=is,ie
624  ! Blend the equatorial deformation radius with the standard one.
625  f2_h = absurdly_small_freq2 + 0.25 * &
626  ((g%CoriolisBu(i,j)**2 + g%CoriolisBu(i-1,j-1)**2) + &
627  (g%CoriolisBu(i-1,j)**2 + g%CoriolisBu(i,j-1)**2))
628  mag_beta = sqrt(0.5 * ( &
629  (((g%CoriolisBu(i,j)-g%CoriolisBu(i-1,j)) * g%IdxCv(i,j))**2 + &
630  ((g%CoriolisBu(i,j-1)-g%CoriolisBu(i-1,j-1)) * g%IdxCv(i,j-1))**2) + &
631  (((g%CoriolisBu(i,j)-g%CoriolisBu(i,j-1)) * g%IdyCu(i,j))**2 + &
632  ((g%CoriolisBu(i-1,j)-g%CoriolisBu(i-1,j-1)) * g%IdyCu(i-1,j))**2) ))
633  cs%Rd1(i,j) = cs%cg1(i,j) / sqrt(f2_h + cs%cg1(i,j) * mag_beta)
634 
635  enddo ; enddo
636  call post_data(cs%id_Rd1, cs%Rd1, cs%diag)
637  endif
638  if (cs%id_cfl_cg1>0) then
639  do j=js,je ; do i=is,ie
640  cs%cfl_cg1(i,j) = (dt*cs%cg1(i,j)) * (g%IdxT(i,j) + g%IdyT(i,j))
641  enddo ; enddo
642  call post_data(cs%id_cfl_cg1, cs%cfl_cg1, cs%diag)
643  endif
644  if (cs%id_cfl_cg1_x>0) then
645  do j=js,je ; do i=is,ie
646  cs%cfl_cg1_x(i,j) = (dt*cs%cg1(i,j)) * g%IdxT(i,j)
647  enddo ; enddo
648  call post_data(cs%id_cfl_cg1_x, cs%cfl_cg1_x, cs%diag)
649  endif
650  if (cs%id_cfl_cg1_y>0) then
651  do j=js,je ; do i=is,ie
652  cs%cfl_cg1_y(i,j) = (dt*cs%cg1(i,j)) * g%IdyT(i,j)
653  enddo ; enddo
654  call post_data(cs%id_cfl_cg1_y, cs%cfl_cg1_y, cs%diag)
655  endif
656  endif
657  if ((cs%id_cg_ebt>0) .or. (cs%id_Rd_ebt>0) .or. (cs%id_p_ebt>0)) then
658  if (cs%id_p_ebt>0) then
659  call wave_speed(h, tv, g, gv, us, cs%cg1, cs%wave_speed_CSp, use_ebt_mode=.true., &
660  mono_n2_column_fraction=cs%mono_N2_column_fraction, &
661  mono_n2_depth=cs%mono_N2_depth, modal_structure=cs%p_ebt)
662  call post_data(cs%id_p_ebt, cs%p_ebt, cs%diag)
663  else
664  call wave_speed(h, tv, g, gv, us, cs%cg1, cs%wave_speed_CSp, use_ebt_mode=.true., &
665  mono_n2_column_fraction=cs%mono_N2_column_fraction, &
666  mono_n2_depth=cs%mono_N2_depth)
667  endif
668  if (cs%id_cg_ebt>0) call post_data(cs%id_cg_ebt, cs%cg1, cs%diag)
669  if (cs%id_Rd_ebt>0) then
670  !$OMP parallel do default(shared) private(f2_h,mag_beta)
671  do j=js,je ; do i=is,ie
672  ! Blend the equatorial deformation radius with the standard one.
673  f2_h = absurdly_small_freq2 + 0.25 * &
674  ((g%CoriolisBu(i,j)**2 + g%CoriolisBu(i-1,j-1)**2) + &
675  (g%CoriolisBu(i-1,j)**2 + g%CoriolisBu(i,j-1)**2))
676  mag_beta = sqrt(0.5 * ( &
677  (((g%CoriolisBu(i,j)-g%CoriolisBu(i-1,j)) * g%IdxCv(i,j))**2 + &
678  ((g%CoriolisBu(i,j-1)-g%CoriolisBu(i-1,j-1)) * g%IdxCv(i,j-1))**2) + &
679  (((g%CoriolisBu(i,j)-g%CoriolisBu(i,j-1)) * g%IdyCu(i,j))**2 + &
680  ((g%CoriolisBu(i-1,j)-g%CoriolisBu(i-1,j-1)) * g%IdyCu(i-1,j))**2) ))
681  cs%Rd1(i,j) = cs%cg1(i,j) / sqrt(f2_h + cs%cg1(i,j) * mag_beta)
682 
683  enddo ; enddo
684  call post_data(cs%id_Rd_ebt, cs%Rd1, cs%diag)
685  endif
686  endif
687 

References calculate_derivs(), calculate_energy_diagnostics(), calculate_vertical_integrals(), mom_diag_mediator::diag_copy_storage_to_diag(), mom_diag_mediator::diag_restore_grids(), mom_diag_mediator::diag_save_grids(), find_weights(), mom_spatial_means::global_volume_mean(), mom_error_handler::mom_error(), and mom_wave_speed::wave_speed().

Here is the call graph for this function:

◆ calculate_energy_diagnostics()

subroutine mom_diagnostics::calculate_energy_diagnostics ( real, dimension( g %isdb: g %iedb, g %jsd: g %jed, g %ke), intent(in)  u,
real, dimension( g %isd: g %ied, g %jsdb: g %jedb, g %ke), intent(in)  v,
real, dimension( g %isd: g %ied, g %jsd: g %jed, g %ke), intent(in)  h,
real, dimension( g %isdb: g %iedb, g %jsd: g %jed, g %ke), intent(in)  uh,
real, dimension( g %isd: g %ied, g %jsdb: g %jedb, g %ke), intent(in)  vh,
type(accel_diag_ptrs), intent(in)  ADp,
type(cont_diag_ptrs), intent(in)  CDp,
type(ocean_grid_type), intent(inout)  G,
type(verticalgrid_type), intent(in)  GV,
type(unit_scale_type), intent(in)  US,
type(diagnostics_cs), intent(inout)  CS 
)
private

This subroutine calculates terms in the mechanical energy budget.

Parameters
[in,out]gThe ocean's grid structure.
[in]gvThe ocean's vertical grid structure.
[in]uThe zonal velocity [L T-1 ~> m s-1].
[in]vThe meridional velocity [L T-1 ~> m s-1].
[in]hLayer thicknesses [H ~> m or kg m-2].
[in]uhTransport through zonal faces=u*h*dy,
[in]vhTransport through merid faces=v*h*dx,
[in]adpStructure pointing to accelerations in momentum equation.
[in]cdpStructure pointing to terms in continuity equations.
[in]usA dimensional unit scaling type
[in,out]csControl structure returned by a previous call to diagnostics_init.

Definition at line 880 of file MOM_diagnostics.F90.

880  type(ocean_grid_type), intent(inout) :: G !< The ocean's grid structure.
881  type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure.
882  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)), &
883  intent(in) :: u !< The zonal velocity [L T-1 ~> m s-1].
884  real, dimension(SZI_(G),SZJB_(G),SZK_(G)), &
885  intent(in) :: v !< The meridional velocity [L T-1 ~> m s-1].
886  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
887  intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2].
888  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)), &
889  intent(in) :: uh !< Transport through zonal faces=u*h*dy,
890  !! [H L2 T-1 ~> m3 s-1 or kg s-1].
891  real, dimension(SZI_(G),SZJB_(G),SZK_(G)), &
892  intent(in) :: vh !< Transport through merid faces=v*h*dx,
893  !! [H L2 T-1 ~> m3 s-1 or kg s-1].
894  type(accel_diag_ptrs), intent(in) :: ADp !< Structure pointing to accelerations in momentum equation.
895  type(cont_diag_ptrs), intent(in) :: CDp !< Structure pointing to terms in continuity equations.
896  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
897  type(diagnostics_CS), intent(inout) :: CS !< Control structure returned by a previous call to
898  !! diagnostics_init.
899 
900 ! This subroutine calculates terms in the mechanical energy budget.
901 
902  ! Local variables
903  real :: KE_u(SZIB_(G),SZJ_(G))
904  real :: KE_v(SZI_(G),SZJB_(G))
905  real :: KE_h(SZI_(G),SZJ_(G))
906 
907  integer :: i, j, k, is, ie, js, je, Isq, Ieq, Jsq, Jeq, nz
908  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
909  isq = g%IscB ; ieq = g%IecB ; jsq = g%JscB ; jeq = g%JecB
910 
911  do j=js-1,je ; do i=is-1,ie
912  ke_u(i,j) = 0.0 ; ke_v(i,j) = 0.0
913  enddo ; enddo
914 
915  if (associated(cs%KE)) then
916  do k=1,nz ; do j=js,je ; do i=is,ie
917  cs%KE(i,j,k) = ((u(i,j,k) * u(i,j,k) + u(i-1,j,k) * u(i-1,j,k)) &
918  + (v(i,j,k) * v(i,j,k) + v(i,j-1,k) * v(i,j-1,k))) * 0.25
919  ! DELETE THE FOLLOWING... Make this 0 to test the momentum balance,
920  ! or a huge number to test the continuity balance.
921  ! CS%KE(i,j,k) *= 1e20
922  enddo ; enddo ; enddo
923  if (cs%id_KE > 0) call post_data(cs%id_KE, cs%KE, cs%diag)
924  endif
925 
926  if (.not.g%symmetric) then
927  if (associated(cs%dKE_dt) .OR. associated(cs%PE_to_KE) .OR. associated(cs%KE_CorAdv) .OR. &
928  associated(cs%KE_adv) .OR. associated(cs%KE_visc) .OR. associated(cs%KE_horvisc).OR. &
929  associated(cs%KE_dia) ) then
930  call create_group_pass(cs%pass_KE_uv, ke_u, ke_v, g%Domain, to_north+to_east)
931  endif
932  endif
933 
934  if (associated(cs%dKE_dt)) then
935  do k=1,nz
936  do j=js,je ; do i=isq,ieq
937  ke_u(i,j) = uh(i,j,k) * g%dxCu(i,j) * cs%du_dt(i,j,k)
938  enddo ; enddo
939  do j=jsq,jeq ; do i=is,ie
940  ke_v(i,j) = vh(i,j,k) * g%dyCv(i,j) * cs%dv_dt(i,j,k)
941  enddo ; enddo
942  do j=js,je ; do i=is,ie
943  ke_h(i,j) = cs%KE(i,j,k) * cs%dh_dt(i,j,k)
944  enddo ; enddo
945  if (.not.g%symmetric) &
946  call do_group_pass(cs%pass_KE_uv, g%domain)
947  do j=js,je ; do i=is,ie
948  cs%dKE_dt(i,j,k) = ke_h(i,j) + 0.5 * g%IareaT(i,j) &
949  * (ke_u(i,j) + ke_u(i-1,j) + ke_v(i,j) + ke_v(i,j-1))
950  enddo ; enddo
951  enddo
952  if (cs%id_dKEdt > 0) call post_data(cs%id_dKEdt, cs%dKE_dt, cs%diag)
953  endif
954 
955  if (associated(cs%PE_to_KE)) then
956  do k=1,nz
957  do j=js,je ; do i=isq,ieq
958  ke_u(i,j) = uh(i,j,k) * g%dxCu(i,j) * adp%PFu(i,j,k)
959  enddo ; enddo
960  do j=jsq,jeq ; do i=is,ie
961  ke_v(i,j) = vh(i,j,k) * g%dyCv(i,j) * adp%PFv(i,j,k)
962  enddo ; enddo
963  if (.not.g%symmetric) &
964  call do_group_pass(cs%pass_KE_uv, g%domain)
965  do j=js,je ; do i=is,ie
966  cs%PE_to_KE(i,j,k) = 0.5 * g%IareaT(i,j) &
967  * (ke_u(i,j) + ke_u(i-1,j) + ke_v(i,j) + ke_v(i,j-1))
968  enddo ; enddo
969  enddo
970  if (cs%id_PE_to_KE > 0) call post_data(cs%id_PE_to_KE, cs%PE_to_KE, cs%diag)
971  endif
972 
973  if (associated(cs%KE_CorAdv)) then
974  do k=1,nz
975  do j=js,je ; do i=isq,ieq
976  ke_u(i,j) = uh(i,j,k) * g%dxCu(i,j) * adp%CAu(i,j,k)
977  enddo ; enddo
978  do j=jsq,jeq ; do i=is,ie
979  ke_v(i,j) = vh(i,j,k) * g%dyCv(i,j) * adp%CAv(i,j,k)
980  enddo ; enddo
981  do j=js,je ; do i=is,ie
982  ke_h(i,j) = -cs%KE(i,j,k) * g%IareaT(i,j) &
983  * (uh(i,j,k) - uh(i-1,j,k) + vh(i,j,k) - vh(i,j-1,k))
984  enddo ; enddo
985  if (.not.g%symmetric) &
986  call do_group_pass(cs%pass_KE_uv, g%domain)
987  do j=js,je ; do i=is,ie
988  cs%KE_CorAdv(i,j,k) = ke_h(i,j) + 0.5 * g%IareaT(i,j) &
989  * (ke_u(i,j) + ke_u(i-1,j) + ke_v(i,j) + ke_v(i,j-1))
990  enddo ; enddo
991  enddo
992  if (cs%id_KE_Coradv > 0) call post_data(cs%id_KE_Coradv, cs%KE_Coradv, cs%diag)
993  endif
994 
995  if (associated(cs%KE_adv)) then
996  ! NOTE: All terms in KE_adv are multipled by -1, which can easily produce
997  ! negative zeros and may signal a reproducibility issue over land.
998  ! We resolve this by re-initializing and only evaluating over water points.
999  ke_u(:,:) = 0. ; ke_v(:,:) = 0.
1000  do k=1,nz
1001  do j=js,je ; do i=isq,ieq
1002  if (g%mask2dCu(i,j) /= 0.) &
1003  ke_u(i,j) = uh(i,j,k) * g%dxCu(i,j) * adp%gradKEu(i,j,k)
1004  enddo ; enddo
1005  do j=jsq,jeq ; do i=is,ie
1006  if (g%mask2dCv(i,j) /= 0.) &
1007  ke_v(i,j) = vh(i,j,k) * g%dyCv(i,j) * adp%gradKEv(i,j,k)
1008  enddo ; enddo
1009  do j=js,je ; do i=is,ie
1010  ke_h(i,j) = -cs%KE(i,j,k) * g%IareaT(i,j) &
1011  * (uh(i,j,k) - uh(i-1,j,k) + vh(i,j,k) - vh(i,j-1,k))
1012  enddo ; enddo
1013  if (.not.g%symmetric) &
1014  call do_group_pass(cs%pass_KE_uv, g%domain)
1015  do j=js,je ; do i=is,ie
1016  cs%KE_adv(i,j,k) = ke_h(i,j) + 0.5 * g%IareaT(i,j) &
1017  * (ke_u(i,j) + ke_u(i-1,j) + ke_v(i,j) + ke_v(i,j-1))
1018  enddo ; enddo
1019  enddo
1020  if (cs%id_KE_adv > 0) call post_data(cs%id_KE_adv, cs%KE_adv, cs%diag)
1021  endif
1022 
1023  if (associated(cs%KE_visc)) then
1024  do k=1,nz
1025  do j=js,je ; do i=isq,ieq
1026  ke_u(i,j) = uh(i,j,k) * g%dxCu(i,j) * adp%du_dt_visc(i,j,k)
1027  enddo ; enddo
1028  do j=jsq,jeq ; do i=is,ie
1029  ke_v(i,j) = vh(i,j,k) * g%dyCv(i,j) * adp%dv_dt_visc(i,j,k)
1030  enddo ; enddo
1031  if (.not.g%symmetric) &
1032  call do_group_pass(cs%pass_KE_uv, g%domain)
1033  do j=js,je ; do i=is,ie
1034  cs%KE_visc(i,j,k) = 0.5 * g%IareaT(i,j) &
1035  * (ke_u(i,j) + ke_u(i-1,j) + ke_v(i,j) + ke_v(i,j-1))
1036  enddo ; enddo
1037  enddo
1038  if (cs%id_KE_visc > 0) call post_data(cs%id_KE_visc, cs%KE_visc, cs%diag)
1039  endif
1040 
1041  if (associated(cs%KE_horvisc)) then
1042  do k=1,nz
1043  do j=js,je ; do i=isq,ieq
1044  ke_u(i,j) = uh(i,j,k) * g%dxCu(i,j) * adp%diffu(i,j,k)
1045  enddo ; enddo
1046  do j=jsq,jeq ; do i=is,ie
1047  ke_v(i,j) = vh(i,j,k) * g%dyCv(i,j) * adp%diffv(i,j,k)
1048  enddo ; enddo
1049  if (.not.g%symmetric) &
1050  call do_group_pass(cs%pass_KE_uv, g%domain)
1051  do j=js,je ; do i=is,ie
1052  cs%KE_horvisc(i,j,k) = 0.5 * g%IareaT(i,j) &
1053  * (ke_u(i,j) + ke_u(i-1,j) + ke_v(i,j) + ke_v(i,j-1))
1054  enddo ; enddo
1055  enddo
1056  if (cs%id_KE_horvisc > 0) call post_data(cs%id_KE_horvisc, cs%KE_horvisc, cs%diag)
1057  endif
1058 
1059  if (associated(cs%KE_dia)) then
1060  do k=1,nz
1061  do j=js,je ; do i=isq,ieq
1062  ke_u(i,j) = uh(i,j,k) * g%dxCu(i,j) * adp%du_dt_dia(i,j,k)
1063  enddo ; enddo
1064  do j=jsq,jeq ; do i=is,ie
1065  ke_v(i,j) = vh(i,j,k) * g%dyCv(i,j) * adp%dv_dt_dia(i,j,k)
1066  enddo ; enddo
1067  do j=js,je ; do i=is,ie
1068  ke_h(i,j) = cs%KE(i,j,k) &
1069  * (us%T_to_s * (cdp%diapyc_vel(i,j,k) - cdp%diapyc_vel(i,j,k+1)))
1070  enddo ; enddo
1071  if (.not.g%symmetric) &
1072  call do_group_pass(cs%pass_KE_uv, g%domain)
1073  do j=js,je ; do i=is,ie
1074  cs%KE_dia(i,j,k) = ke_h(i,j) + 0.5 * g%IareaT(i,j) &
1075  * (ke_u(i,j) + ke_u(i-1,j) + ke_v(i,j) + ke_v(i,j-1))
1076  enddo ; enddo
1077  enddo
1078  if (cs%id_KE_dia > 0) call post_data(cs%id_KE_dia, cs%KE_dia, cs%diag)
1079  endif
1080 

Referenced by calculate_diagnostic_fields().

Here is the caller graph for this function:

◆ calculate_vertical_integrals()

subroutine mom_diagnostics::calculate_vertical_integrals ( real, dimension( g %isd: g %ied, g %jsd: g %jed, g %ke), intent(in)  h,
type(thermo_var_ptrs), intent(in)  tv,
real, dimension(:,:), pointer  p_surf,
type(ocean_grid_type), intent(inout)  G,
type(verticalgrid_type), intent(in)  GV,
type(unit_scale_type), intent(in)  US,
type(diagnostics_cs), intent(inout)  CS 
)
private

This subroutine calculates vertical integrals of several tracers, along with the mass-weight of these tracers, the total column mass, and the carefully calculated column height.

Parameters
[in,out]gThe ocean's grid structure.
[in]gvThe ocean's vertical grid structure.
[in]usA dimensional unit scaling type
[in]hLayer thicknesses [H ~> m or kg m-2].
[in]tvA structure pointing to various thermodynamic variables.
p_surfA pointer to the surface pressure [Pa]. If p_surf is not associated, it is the same as setting the surface pressure to 0.
[in,out]csControl structure returned by a previous call to diagnostics_init.

Definition at line 762 of file MOM_diagnostics.F90.

762  type(ocean_grid_type), intent(inout) :: G !< The ocean's grid structure.
763  type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure.
764  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
765  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
766  intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2].
767  type(thermo_var_ptrs), intent(in) :: tv !< A structure pointing to various
768  !! thermodynamic variables.
769  real, dimension(:,:), pointer :: p_surf !< A pointer to the surface pressure [Pa].
770  !! If p_surf is not associated, it is the same
771  !! as setting the surface pressure to 0.
772  type(diagnostics_CS), intent(inout) :: CS !< Control structure returned by a
773  !! previous call to diagnostics_init.
774 
775  real, dimension(SZI_(G), SZJ_(G)) :: &
776  z_top, & ! Height of the top of a layer or the ocean [Z ~> m].
777  z_bot, & ! Height of the bottom of a layer (for id_mass) or the
778  ! (positive) depth of the ocean (for id_col_ht) [Z ~> m].
779  mass, & ! integrated mass of the water column [kg m-2]. For
780  ! non-Boussinesq models this is rho*dz. For Boussinesq
781  ! models, this is either the integral of in-situ density
782  ! (rho*dz for col_mass) or reference density (Rho_0*dz for mass_wt).
783  btm_pres,&! The pressure at the ocean bottom, or CMIP variable 'pbo'.
784  ! This is the column mass multiplied by gravity plus the pressure
785  ! at the ocean surface [Pa].
786  dpress, & ! Change in hydrostatic pressure across a layer [Pa].
787  tr_int ! vertical integral of a tracer times density,
788  ! (Rho_0 in a Boussinesq model) [TR kg m-2].
789  real :: IG_Earth ! Inverse of gravitational acceleration [s2 m-1].
790 
791  integer :: i, j, k, is, ie, js, je, nz
792  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
793 
794  if (cs%id_mass_wt > 0) then
795  do j=js,je ; do i=is,ie ; mass(i,j) = 0.0 ; enddo ; enddo
796  do k=1,nz ; do j=js,je ; do i=is,ie
797  mass(i,j) = mass(i,j) + gv%H_to_kg_m2*h(i,j,k)
798  enddo ; enddo ; enddo
799  call post_data(cs%id_mass_wt, mass, cs%diag)
800  endif
801 
802  if (cs%id_temp_int > 0) then
803  do j=js,je ; do i=is,ie ; tr_int(i,j) = 0.0 ; enddo ; enddo
804  do k=1,nz ; do j=js,je ; do i=is,ie
805  tr_int(i,j) = tr_int(i,j) + (gv%H_to_kg_m2*h(i,j,k))*tv%T(i,j,k)
806  enddo ; enddo ; enddo
807  call post_data(cs%id_temp_int, tr_int, cs%diag)
808  endif
809 
810  if (cs%id_salt_int > 0) then
811  do j=js,je ; do i=is,ie ; tr_int(i,j) = 0.0 ; enddo ; enddo
812  do k=1,nz ; do j=js,je ; do i=is,ie
813  tr_int(i,j) = tr_int(i,j) + (gv%H_to_kg_m2*h(i,j,k))*tv%S(i,j,k)
814  enddo ; enddo ; enddo
815  call post_data(cs%id_salt_int, tr_int, cs%diag)
816  endif
817 
818  if (cs%id_col_ht > 0) then
819  call find_eta(h, tv, g, gv, us, z_top)
820  do j=js,je ; do i=is,ie
821  z_bot(i,j) = z_top(i,j) + g%bathyT(i,j)
822  enddo ; enddo
823  call post_data(cs%id_col_ht, z_bot, cs%diag)
824  endif
825 
826  ! NOTE: int_density_z expects z_top and z_btm values from [ij]sq to [ij]eq+1
827  if (cs%id_col_mass > 0 .or. cs%id_pbo > 0) then
828  do j=js,je ; do i=is,ie ; mass(i,j) = 0.0 ; enddo ; enddo
829  if (gv%Boussinesq) then
830  if (associated(tv%eqn_of_state)) then
831  ig_earth = 1.0 / gv%mks_g_Earth
832 ! do j=js,je ; do i=is,ie ; z_bot(i,j) = -P_SURF(i,j)/GV%H_to_Pa ; enddo ; enddo
833  do j=g%jscB,g%jecB+1 ; do i=g%iscB,g%iecB+1
834  z_bot(i,j) = 0.0
835  enddo ; enddo
836  do k=1,nz
837  do j=g%jscB,g%jecB+1 ; do i=g%iscB,g%iecB+1
838  z_top(i,j) = z_bot(i,j)
839  z_bot(i,j) = z_top(i,j) - gv%H_to_Z*h(i,j,k)
840  enddo ; enddo
841  call int_density_dz(tv%T(:,:,k), tv%S(:,:,k), &
842  z_top, z_bot, 0.0, us%R_to_kg_m3*gv%Rho0, gv%mks_g_Earth*us%Z_to_m, &
843  g%HI, g%HI, tv%eqn_of_state, dpress)
844  do j=js,je ; do i=is,ie
845  mass(i,j) = mass(i,j) + dpress(i,j) * ig_earth
846  enddo ; enddo
847  enddo
848  else
849  do k=1,nz ; do j=js,je ; do i=is,ie
850  mass(i,j) = mass(i,j) + (gv%H_to_m*us%R_to_kg_m3*gv%Rlay(k))*h(i,j,k)
851  enddo ; enddo ; enddo
852  endif
853  else
854  do k=1,nz ; do j=js,je ; do i=is,ie
855  mass(i,j) = mass(i,j) + gv%H_to_kg_m2*h(i,j,k)
856  enddo ; enddo ; enddo
857  endif
858  if (cs%id_col_mass > 0) then
859  call post_data(cs%id_col_mass, mass, cs%diag)
860  endif
861  if (cs%id_pbo > 0) then
862  do j=js,je ; do i=is,ie ; btm_pres(i,j) = 0.0 ; enddo ; enddo
863  ! 'pbo' is defined as the sea water pressure at the sea floor
864  ! pbo = (mass * g) + p_surf
865  ! where p_surf is the sea water pressure at sea water surface.
866  do j=js,je ; do i=is,ie
867  btm_pres(i,j) = mass(i,j) * gv%mks_g_Earth
868  if (associated(p_surf)) then
869  btm_pres(i,j) = btm_pres(i,j) + p_surf(i,j)
870  endif
871  enddo ; enddo
872  call post_data(cs%id_pbo, btm_pres, cs%diag)
873  endif
874  endif
875 

Referenced by calculate_diagnostic_fields().

Here is the caller graph for this function:

◆ find_weights()

subroutine mom_diagnostics::find_weights ( real, dimension(:), intent(in)  Rlist,
real, intent(in)  R_in,
integer, intent(inout)  k,
integer, intent(in)  nz,
real, intent(out)  wt,
real, intent(out)  wt_p 
)
private

This subroutine finds the location of R_in in an increasing ordered list, Rlist, returning as k the element such that Rlist(k) <= R_in < Rlist(k+1), and where wt and wt_p are the linear weights that should be assigned to elements k and k+1.

Parameters
[in]rlistThe list of target densities [R ~> kg m-3]
[in]r_inThe density being inserted into Rlist [R ~> kg m-3]
[in,out]kThe value of k such that Rlist(k) <= R_in < Rlist(k+1) The input value is a first guess
[in]nzThe number of layers in Rlist
[out]wtThe weight of layer k for interpolation, nondim
[out]wt_pThe weight of layer k+1 for interpolation, nondim

Definition at line 695 of file MOM_diagnostics.F90.

695  real, dimension(:), &
696  intent(in) :: Rlist !< The list of target densities [R ~> kg m-3]
697  real, intent(in) :: R_in !< The density being inserted into Rlist [R ~> kg m-3]
698  integer, intent(inout) :: k !< The value of k such that Rlist(k) <= R_in < Rlist(k+1)
699  !! The input value is a first guess
700  integer, intent(in) :: nz !< The number of layers in Rlist
701  real, intent(out) :: wt !< The weight of layer k for interpolation, nondim
702  real, intent(out) :: wt_p !< The weight of layer k+1 for interpolation, nondim
703 
704  ! This subroutine finds location of R_in in an increasing ordered
705  ! list, Rlist, returning as k the element such that
706  ! Rlist(k) <= R_in < Rlist(k+1), and where wt and wt_p are the linear
707  ! weights that should be assigned to elements k and k+1.
708 
709  integer :: k_upper, k_lower, k_new, inc
710 
711  ! First, bracket the desired point.
712  if ((k < 1) .or. (k > nz)) k = nz/2
713 
714  k_upper = k ; k_lower = k ; inc = 1
715  if (r_in < rlist(k)) then
716  do
717  k_lower = max(k_lower-inc, 1)
718  if ((k_lower == 1) .or. (r_in >= rlist(k_lower))) exit
719  k_upper = k_lower
720  inc = inc*2
721  enddo
722  else
723  do
724  k_upper = min(k_upper+inc, nz)
725  if ((k_upper == nz) .or. (r_in < rlist(k_upper))) exit
726  k_lower = k_upper
727  inc = inc*2
728  enddo
729  endif
730 
731  if ((k_lower == 1) .and. (r_in <= rlist(k_lower))) then
732  k = 1 ; wt = 1.0 ; wt_p = 0.0
733  elseif ((k_upper == nz) .and. (r_in >= rlist(k_upper))) then
734  k = nz-1 ; wt = 0.0 ; wt_p = 1.0
735  else
736  do
737  if (k_upper <= k_lower+1) exit
738  k_new = (k_upper + k_lower) / 2
739  if (r_in < rlist(k_new)) then
740  k_upper = k_new
741  else
742  k_lower = k_new
743  endif
744  enddo
745 
746 ! Uncomment this as a code check
747 ! if ((R_in < Rlist(k_lower)) .or. (R_in >= Rlist(k_upper)) .or. (k_upper-k_lower /= 1)) &
748 ! write (*,*) "Error: ",R_in," is not between R(",k_lower,") = ", &
749 ! Rlist(k_lower)," and R(",k_upper,") = ",Rlist(k_upper),"."
750  k = k_lower
751  wt = (rlist(k_upper) - r_in) / (rlist(k_upper) - rlist(k_lower))
752  wt_p = 1.0 - wt
753 
754  endif
755 

Referenced by calculate_diagnostic_fields().

Here is the caller graph for this function:

◆ mom_diagnostics_end()

subroutine, public mom_diagnostics::mom_diagnostics_end ( type(diagnostics_cs), pointer  CS,
type(accel_diag_ptrs), intent(inout)  ADp 
)

Deallocate memory associated with the diagnostics module.

Parameters
csControl structure returned by a previous call to diagnostics_init.
[in,out]adpstructure with pointers to accelerations in momentum equation.

Definition at line 2088 of file MOM_diagnostics.F90.

2088  type(diagnostics_CS), pointer :: CS !< Control structure returned by a
2089  !! previous call to diagnostics_init.
2090  type(accel_diag_ptrs), intent(inout) :: ADp !< structure with pointers to
2091  !! accelerations in momentum equation.
2092  integer :: m
2093 
2094  if (associated(cs%e)) deallocate(cs%e)
2095  if (associated(cs%e_D)) deallocate(cs%e_D)
2096  if (associated(cs%KE)) deallocate(cs%KE)
2097  if (associated(cs%dKE_dt)) deallocate(cs%dKE_dt)
2098  if (associated(cs%PE_to_KE)) deallocate(cs%PE_to_KE)
2099  if (associated(cs%KE_Coradv)) deallocate(cs%KE_Coradv)
2100  if (associated(cs%KE_adv)) deallocate(cs%KE_adv)
2101  if (associated(cs%KE_visc)) deallocate(cs%KE_visc)
2102  if (associated(cs%KE_horvisc)) deallocate(cs%KE_horvisc)
2103  if (associated(cs%KE_dia)) deallocate(cs%KE_dia)
2104  if (associated(cs%dv_dt)) deallocate(cs%dv_dt)
2105  if (associated(cs%dh_dt)) deallocate(cs%dh_dt)
2106  if (associated(cs%du_dt)) deallocate(cs%du_dt)
2107  if (associated(cs%h_Rlay)) deallocate(cs%h_Rlay)
2108  if (associated(cs%uh_Rlay)) deallocate(cs%uh_Rlay)
2109  if (associated(cs%vh_Rlay)) deallocate(cs%vh_Rlay)
2110  if (associated(cs%uhGM_Rlay)) deallocate(cs%uhGM_Rlay)
2111  if (associated(cs%vhGM_Rlay)) deallocate(cs%vhGM_Rlay)
2112 
2113  if (associated(adp%gradKEu)) deallocate(adp%gradKEu)
2114  if (associated(adp%gradKEu)) deallocate(adp%gradKEu)
2115  if (associated(adp%du_dt_visc)) deallocate(adp%du_dt_visc)
2116  if (associated(adp%dv_dt_visc)) deallocate(adp%dv_dt_visc)
2117  if (associated(adp%du_dt_dia)) deallocate(adp%du_dt_dia)
2118  if (associated(adp%dv_dt_dia)) deallocate(adp%dv_dt_dia)
2119  if (associated(adp%du_other)) deallocate(adp%du_other)
2120  if (associated(adp%dv_other)) deallocate(adp%dv_other)
2121 
2122  do m=1,cs%num_time_deriv ; deallocate(cs%prev_val(m)%p) ; enddo
2123 
2124  deallocate(cs)
2125 

◆ mom_diagnostics_init()

subroutine, public mom_diagnostics::mom_diagnostics_init ( type(ocean_internal_state), intent(in)  MIS,
type(accel_diag_ptrs), intent(inout)  ADp,
type(cont_diag_ptrs), intent(inout)  CDp,
type(time_type), intent(in)  Time,
type(ocean_grid_type), intent(in)  G,
type(verticalgrid_type), intent(in)  GV,
type(unit_scale_type), intent(in)  US,
type(param_file_type), intent(in)  param_file,
type(diag_ctrl), intent(inout), target  diag,
type(diagnostics_cs), pointer  CS,
type(thermo_var_ptrs), intent(in)  tv 
)

This subroutine registers various diagnostics and allocates space for fields that other diagnostis depend upon.

Parameters
[in]misFor "MOM Internal State" a set of pointers to the fields and accelerations that make up the ocean's internal physical state.
[in,out]adpStructure with pointers to momentum equation terms.
[in,out]cdpStructure with pointers to continuity equation terms.
[in]timeCurrent model time.
[in]gThe ocean's grid structure.
[in]gvThe ocean's vertical grid structure.
[in]usA dimensional unit scaling type
[in]param_fileA structure to parse for run-time parameters.
[in,out]diagStructure to regulate diagnostic output.
csPointer set to point to control structure for this module.
[in]tvA structure pointing to various thermodynamic variables.

Definition at line 1423 of file MOM_diagnostics.F90.

1423  type(ocean_internal_state), intent(in) :: MIS !< For "MOM Internal State" a set of pointers to
1424  !! the fields and accelerations that make up the
1425  !! ocean's internal physical state.
1426  type(accel_diag_ptrs), intent(inout) :: ADp !< Structure with pointers to momentum equation
1427  !! terms.
1428  type(cont_diag_ptrs), intent(inout) :: CDp !< Structure with pointers to continuity
1429  !! equation terms.
1430  type(time_type), intent(in) :: Time !< Current model time.
1431  type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure.
1432  type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure.
1433  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
1434  type(param_file_type), intent(in) :: param_file !< A structure to parse for run-time
1435  !! parameters.
1436  type(diag_ctrl), target, intent(inout) :: diag !< Structure to regulate diagnostic output.
1437  type(diagnostics_CS), pointer :: CS !< Pointer set to point to control structure
1438  !! for this module.
1439  type(thermo_var_ptrs), intent(in) :: tv !< A structure pointing to various
1440  !! thermodynamic variables.
1441 
1442  ! Local variables
1443  real :: omega, f2_min, convert_H
1444  ! This include declares and sets the variable "version".
1445 # include "version_variable.h"
1446  character(len=40) :: mdl = "MOM_diagnostics" ! This module's name.
1447  character(len=48) :: thickness_units, flux_units
1448  logical :: use_temperature, adiabatic
1449  integer :: isd, ied, jsd, jed, IsdB, IedB, JsdB, JedB, nz, nkml, nkbl
1450  integer :: is, ie, js, je, Isq, Ieq, Jsq, Jeq, i, j
1451 
1452  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
1453  isq = g%IscB ; ieq = g%IecB ; jsq = g%JscB ; jeq = g%JecB
1454  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed ; nz = g%ke
1455  isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
1456 
1457  if (associated(cs)) then
1458  call mom_error(warning, "MOM_diagnostics_init called with an associated "// &
1459  "control structure.")
1460  return
1461  endif
1462  allocate(cs)
1463 
1464  cs%diag => diag
1465  use_temperature = associated(tv%T)
1466  call get_param(param_file, mdl, "ADIABATIC", adiabatic, default=.false., &
1467  do_not_log=.true.)
1468 
1469  ! Read all relevant parameters and write them to the model log.
1470  call log_version(param_file, mdl, version)
1471  call get_param(param_file, mdl, "DIAG_EBT_MONO_N2_COLUMN_FRACTION", cs%mono_N2_column_fraction, &
1472  "The lower fraction of water column over which N2 is limited as monotonic "// &
1473  "for the purposes of calculating the equivalent barotropic wave speed.", &
1474  units='nondim', default=0.)
1475  call get_param(param_file, mdl, "DIAG_EBT_MONO_N2_DEPTH", cs%mono_N2_depth, &
1476  "The depth below which N2 is limited as monotonic for the "// &
1477  "purposes of calculating the equivalent barotropic wave speed.", &
1478  units='m', scale=us%m_to_Z, default=-1.)
1479 
1480  if (gv%Boussinesq) then
1481  thickness_units = "m" ; flux_units = "m3 s-1" ; convert_h = gv%H_to_m
1482  else
1483  thickness_units = "kg m-2" ; flux_units = "kg s-1" ; convert_h = gv%H_to_kg_m2
1484  endif
1485 
1486  cs%id_masscello = register_diag_field('ocean_model', 'masscello', diag%axesTL,&
1487  time, 'Mass per unit area of liquid ocean grid cell', 'kg m-2', &
1488  standard_name='sea_water_mass_per_unit_area', v_extensive=.true.)
1489 
1490  cs%id_masso = register_scalar_field('ocean_model', 'masso', time, &
1491  diag, 'Mass of liquid ocean', 'kg', standard_name='sea_water_mass')
1492 
1493  cs%id_thkcello = register_diag_field('ocean_model', 'thkcello', diag%axesTL, time, &
1494  long_name = 'Cell Thickness', standard_name='cell_thickness', units='m', v_extensive=.true.)
1495  cs%id_h_pre_sync = register_diag_field('ocean_model', 'h_pre_sync', diag%axesTL, time, &
1496  long_name = 'Cell thickness from the previous timestep', units='m', &
1497  v_extensive=.true., conversion=gv%H_to_m)
1498 
1499  ! Note that CS%id_volcello would normally be registered here but because it is a "cell measure" and
1500  ! must be registered first. We earlier stored the handle of volcello but need it here for posting
1501  ! by this module.
1502  cs%id_volcello = diag_get_volume_cell_measure_dm_id(diag)
1503 
1504  if (use_temperature) then
1505  if (tv%T_is_conT) then
1506  cs%id_Tpot = register_diag_field('ocean_model', 'temp', diag%axesTL, &
1507  time, 'Potential Temperature', 'degC')
1508  endif
1509  if (tv%S_is_absS) then
1510  cs%id_Sprac = register_diag_field('ocean_model', 'salt', diag%axesTL, &
1511  time, 'Salinity', 'psu')
1512  endif
1513 
1514  cs%id_tob = register_diag_field('ocean_model','tob', diag%axesT1, time, &
1515  long_name='Sea Water Potential Temperature at Sea Floor', &
1516  standard_name='sea_water_potential_temperature_at_sea_floor', units='degC')
1517  cs%id_sob = register_diag_field('ocean_model','sob',diag%axesT1, time, &
1518  long_name='Sea Water Salinity at Sea Floor', &
1519  standard_name='sea_water_salinity_at_sea_floor', units='psu')
1520 
1521  cs%id_temp_layer_ave = register_diag_field('ocean_model', 'temp_layer_ave', &
1522  diag%axesZL, time, 'Layer Average Ocean Temperature', 'degC')
1523  cs%id_salt_layer_ave = register_diag_field('ocean_model', 'salt_layer_ave', &
1524  diag%axesZL, time, 'Layer Average Ocean Salinity', 'psu')
1525 
1526  cs%id_thetaoga = register_scalar_field('ocean_model', 'thetaoga', &
1527  time, diag, 'Global Mean Ocean Potential Temperature', 'degC',&
1528  standard_name='sea_water_potential_temperature')
1529  cs%id_soga = register_scalar_field('ocean_model', 'soga', &
1530  time, diag, 'Global Mean Ocean Salinity', 'psu', &
1531  standard_name='sea_water_salinity')
1532 
1533  cs%id_tosga = register_scalar_field('ocean_model', 'sst_global', time, diag,&
1534  long_name='Global Area Average Sea Surface Temperature', &
1535  units='degC', standard_name='sea_surface_temperature', &
1536  cmor_field_name='tosga', cmor_standard_name='sea_surface_temperature', &
1537  cmor_long_name='Sea Surface Temperature')
1538  cs%id_sosga = register_scalar_field('ocean_model', 'sss_global', time, diag,&
1539  long_name='Global Area Average Sea Surface Salinity', &
1540  units='psu', standard_name='sea_surface_salinity', &
1541  cmor_field_name='sosga', cmor_standard_name='sea_surface_salinity', &
1542  cmor_long_name='Sea Surface Salinity')
1543  endif
1544 
1545  cs%id_u = register_diag_field('ocean_model', 'u', diag%axesCuL, time, &
1546  'Zonal velocity', 'm s-1', conversion=us%L_T_to_m_s, cmor_field_name='uo', &
1547  cmor_standard_name='sea_water_x_velocity', cmor_long_name='Sea Water X Velocity')
1548  cs%id_v = register_diag_field('ocean_model', 'v', diag%axesCvL, time, &
1549  'Meridional velocity', 'm s-1', conversion=us%L_T_to_m_s, cmor_field_name='vo', &
1550  cmor_standard_name='sea_water_y_velocity', cmor_long_name='Sea Water Y Velocity')
1551  cs%id_h = register_diag_field('ocean_model', 'h', diag%axesTL, time, &
1552  'Layer Thickness', thickness_units, v_extensive=.true., conversion=convert_h)
1553 
1554  cs%id_e = register_diag_field('ocean_model', 'e', diag%axesTi, time, &
1555  'Interface Height Relative to Mean Sea Level', 'm', conversion=us%Z_to_m)
1556  if (cs%id_e>0) call safe_alloc_ptr(cs%e,isd,ied,jsd,jed,nz+1)
1557 
1558  cs%id_e_D = register_diag_field('ocean_model', 'e_D', diag%axesTi, time, &
1559  'Interface Height above the Seafloor', 'm', conversion=us%Z_to_m)
1560  if (cs%id_e_D>0) call safe_alloc_ptr(cs%e_D,isd,ied,jsd,jed,nz+1)
1561 
1562  cs%id_Rml = register_diag_field('ocean_model', 'Rml', diag%axesTL, time, &
1563  'Mixed Layer Coordinate Potential Density', 'kg m-3', conversion=us%R_to_kg_m3)
1564 
1565  cs%id_Rcv = register_diag_field('ocean_model', 'Rho_cv', diag%axesTL, time, &
1566  'Coordinate Potential Density', 'kg m-3', conversion=us%R_to_kg_m3)
1567 
1568  cs%id_rhopot0 = register_diag_field('ocean_model', 'rhopot0', diag%axesTL, time, &
1569  'Potential density referenced to surface', 'kg m-3')
1570  cs%id_rhopot2 = register_diag_field('ocean_model', 'rhopot2', diag%axesTL, time, &
1571  'Potential density referenced to 2000 dbar', 'kg m-3')
1572  cs%id_rhoinsitu = register_diag_field('ocean_model', 'rhoinsitu', diag%axesTL, time, &
1573  'In situ density', 'kg m-3')
1574 
1575  cs%id_du_dt = register_diag_field('ocean_model', 'dudt', diag%axesCuL, time, &
1576  'Zonal Acceleration', 'm s-2', conversion=us%L_T2_to_m_s2)
1577  if ((cs%id_du_dt>0) .and. .not.associated(cs%du_dt)) then
1578  call safe_alloc_ptr(cs%du_dt,isdb,iedb,jsd,jed,nz)
1579  call register_time_deriv(lbound(mis%u), mis%u, cs%du_dt, cs)
1580  endif
1581 
1582  cs%id_dv_dt = register_diag_field('ocean_model', 'dvdt', diag%axesCvL, time, &
1583  'Meridional Acceleration', 'm s-2', conversion=us%L_T2_to_m_s2)
1584  if ((cs%id_dv_dt>0) .and. .not.associated(cs%dv_dt)) then
1585  call safe_alloc_ptr(cs%dv_dt,isd,ied,jsdb,jedb,nz)
1586  call register_time_deriv(lbound(mis%v), mis%v, cs%dv_dt, cs)
1587  endif
1588 
1589  cs%id_dh_dt = register_diag_field('ocean_model', 'dhdt', diag%axesTL, time, &
1590  'Thickness tendency', trim(thickness_units)//" s-1", conversion=convert_h*us%s_to_T, v_extensive=.true.)
1591  if ((cs%id_dh_dt>0) .and. .not.associated(cs%dh_dt)) then
1592  call safe_alloc_ptr(cs%dh_dt,isd,ied,jsd,jed,nz)
1593  call register_time_deriv(lbound(mis%h), mis%h, cs%dh_dt, cs)
1594  endif
1595 
1596  ! layer thickness variables
1597  !if (GV%nk_rho_varies > 0) then
1598  cs%id_h_Rlay = register_diag_field('ocean_model', 'h_rho', diag%axesTL, time, &
1599  'Layer thicknesses in pure potential density coordinates', thickness_units, &
1600  conversion=convert_h)
1601  if (cs%id_h_Rlay>0) call safe_alloc_ptr(cs%h_Rlay,isd,ied,jsd,jed,nz)
1602 
1603  cs%id_uh_Rlay = register_diag_field('ocean_model', 'uh_rho', diag%axesCuL, time, &
1604  'Zonal volume transport in pure potential density coordinates', flux_units, &
1605  conversion=us%L_to_m**2*us%s_to_T*convert_h)
1606  if (cs%id_uh_Rlay>0) call safe_alloc_ptr(cs%uh_Rlay,isdb,iedb,jsd,jed,nz)
1607 
1608  cs%id_vh_Rlay = register_diag_field('ocean_model', 'vh_rho', diag%axesCvL, time, &
1609  'Meridional volume transport in pure potential density coordinates', flux_units, &
1610  conversion=us%L_to_m**2*us%s_to_T*convert_h)
1611  if (cs%id_vh_Rlay>0) call safe_alloc_ptr(cs%vh_Rlay,isd,ied,jsdb,jedb,nz)
1612 
1613  cs%id_uhGM_Rlay = register_diag_field('ocean_model', 'uhGM_rho', diag%axesCuL, time, &
1614  'Zonal volume transport due to interface height diffusion in pure potential '//&
1615  'density coordinates', flux_units, conversion=us%L_to_m**2*us%s_to_T*convert_h)
1616  if (cs%id_uhGM_Rlay>0) call safe_alloc_ptr(cs%uhGM_Rlay,isdb,iedb,jsd,jed,nz)
1617 
1618  cs%id_vhGM_Rlay = register_diag_field('ocean_model', 'vhGM_rho', diag%axesCvL, time, &
1619  'Meridional volume transport due to interface height diffusion in pure potential '//&
1620  'density coordinates', flux_units, conversion=us%L_to_m**2*us%s_to_T*convert_h)
1621  if (cs%id_vhGM_Rlay>0) call safe_alloc_ptr(cs%vhGM_Rlay,isd,ied,jsdb,jedb,nz)
1622  !endif
1623 
1624 
1625  ! terms in the kinetic energy budget
1626  cs%id_KE = register_diag_field('ocean_model', 'KE', diag%axesTL, time, &
1627  'Layer kinetic energy per unit mass', 'm2 s-2', &
1628  conversion=us%L_T_to_m_s**2)
1629  if (cs%id_KE>0) call safe_alloc_ptr(cs%KE,isd,ied,jsd,jed,nz)
1630 
1631  cs%id_dKEdt = register_diag_field('ocean_model', 'dKE_dt', diag%axesTL, time, &
1632  'Kinetic Energy Tendency of Layer', 'm3 s-3', &
1633  conversion=gv%H_to_m*(us%L_T_to_m_s**2)*us%s_to_T)
1634  if (cs%id_dKEdt>0) call safe_alloc_ptr(cs%dKE_dt,isd,ied,jsd,jed,nz)
1635 
1636  cs%id_PE_to_KE = register_diag_field('ocean_model', 'PE_to_KE', diag%axesTL, time, &
1637  'Potential to Kinetic Energy Conversion of Layer', 'm3 s-3', &
1638  conversion=gv%H_to_m*(us%L_T_to_m_s**2)*us%s_to_T)
1639  if (cs%id_PE_to_KE>0) call safe_alloc_ptr(cs%PE_to_KE,isd,ied,jsd,jed,nz)
1640 
1641  cs%id_KE_Coradv = register_diag_field('ocean_model', 'KE_Coradv', diag%axesTL, time, &
1642  'Kinetic Energy Source from Coriolis and Advection', 'm3 s-3', &
1643  conversion=gv%H_to_m*(us%L_T_to_m_s**2)*us%s_to_T)
1644  if (cs%id_KE_Coradv>0) call safe_alloc_ptr(cs%KE_Coradv,isd,ied,jsd,jed,nz)
1645 
1646  cs%id_KE_adv = register_diag_field('ocean_model', 'KE_adv', diag%axesTL, time, &
1647  'Kinetic Energy Source from Advection', 'm3 s-3', &
1648  conversion=gv%H_to_m*(us%L_T_to_m_s**2)*us%s_to_T)
1649  if (cs%id_KE_adv>0) call safe_alloc_ptr(cs%KE_adv,isd,ied,jsd,jed,nz)
1650 
1651  cs%id_KE_visc = register_diag_field('ocean_model', 'KE_visc', diag%axesTL, time, &
1652  'Kinetic Energy Source from Vertical Viscosity and Stresses', 'm3 s-3', &
1653  conversion=gv%H_to_m*(us%L_T_to_m_s**2)*us%s_to_T)
1654  if (cs%id_KE_visc>0) call safe_alloc_ptr(cs%KE_visc,isd,ied,jsd,jed,nz)
1655 
1656  cs%id_KE_horvisc = register_diag_field('ocean_model', 'KE_horvisc', diag%axesTL, time, &
1657  'Kinetic Energy Source from Horizontal Viscosity', 'm3 s-3', &
1658  conversion=gv%H_to_m*(us%L_T_to_m_s**2)*us%s_to_T)
1659  if (cs%id_KE_horvisc>0) call safe_alloc_ptr(cs%KE_horvisc,isd,ied,jsd,jed,nz)
1660 
1661  if (.not. adiabatic) then
1662  cs%id_KE_dia = register_diag_field('ocean_model', 'KE_dia', diag%axesTL, time, &
1663  'Kinetic Energy Source from Diapycnal Diffusion', 'm3 s-3', &
1664  conversion=gv%H_to_m*(us%L_T_to_m_s**2)*us%s_to_T)
1665  if (cs%id_KE_dia>0) call safe_alloc_ptr(cs%KE_dia,isd,ied,jsd,jed,nz)
1666  endif
1667 
1668  ! gravity wave CFLs
1669  cs%id_cg1 = register_diag_field('ocean_model', 'cg1', diag%axesT1, time, &
1670  'First baroclinic gravity wave speed', 'm s-1', conversion=us%L_T_to_m_s)
1671  cs%id_Rd1 = register_diag_field('ocean_model', 'Rd1', diag%axesT1, time, &
1672  'First baroclinic deformation radius', 'm', conversion=us%L_to_m)
1673  cs%id_cfl_cg1 = register_diag_field('ocean_model', 'CFL_cg1', diag%axesT1, time, &
1674  'CFL of first baroclinic gravity wave = dt*cg1*(1/dx+1/dy)', 'nondim')
1675  cs%id_cfl_cg1_x = register_diag_field('ocean_model', 'CFL_cg1_x', diag%axesT1, time, &
1676  'i-component of CFL of first baroclinic gravity wave = dt*cg1*/dx', 'nondim')
1677  cs%id_cfl_cg1_y = register_diag_field('ocean_model', 'CFL_cg1_y', diag%axesT1, time, &
1678  'j-component of CFL of first baroclinic gravity wave = dt*cg1*/dy', 'nondim')
1679  cs%id_cg_ebt = register_diag_field('ocean_model', 'cg_ebt', diag%axesT1, time, &
1680  'Equivalent barotropic gravity wave speed', 'm s-1', conversion=us%L_T_to_m_s)
1681  cs%id_Rd_ebt = register_diag_field('ocean_model', 'Rd_ebt', diag%axesT1, time, &
1682  'Equivalent barotropic deformation radius', 'm', conversion=us%L_to_m)
1683  cs%id_p_ebt = register_diag_field('ocean_model', 'p_ebt', diag%axesTL, time, &
1684  'Equivalent barotropic modal strcuture', 'nondim')
1685 
1686  if ((cs%id_cg1>0) .or. (cs%id_Rd1>0) .or. (cs%id_cfl_cg1>0) .or. &
1687  (cs%id_cfl_cg1_x>0) .or. (cs%id_cfl_cg1_y>0) .or. &
1688  (cs%id_cg_ebt>0) .or. (cs%id_Rd_ebt>0) .or. (cs%id_p_ebt>0)) then
1689  call wave_speed_init(cs%wave_speed_CSp)
1690  call safe_alloc_ptr(cs%cg1,isd,ied,jsd,jed)
1691  if (cs%id_Rd1>0) call safe_alloc_ptr(cs%Rd1,isd,ied,jsd,jed)
1692  if (cs%id_Rd_ebt>0) call safe_alloc_ptr(cs%Rd1,isd,ied,jsd,jed)
1693  if (cs%id_cfl_cg1>0) call safe_alloc_ptr(cs%cfl_cg1,isd,ied,jsd,jed)
1694  if (cs%id_cfl_cg1_x>0) call safe_alloc_ptr(cs%cfl_cg1_x,isd,ied,jsd,jed)
1695  if (cs%id_cfl_cg1_y>0) call safe_alloc_ptr(cs%cfl_cg1_y,isd,ied,jsd,jed)
1696  if (cs%id_p_ebt>0) call safe_alloc_ptr(cs%p_ebt,isd,ied,jsd,jed,nz)
1697  endif
1698 
1699  cs%id_mass_wt = register_diag_field('ocean_model', 'mass_wt', diag%axesT1, time, &
1700  'The column mass for calculating mass-weighted average properties', 'kg m-2')
1701 
1702  if (use_temperature) then
1703  cs%id_temp_int = register_diag_field('ocean_model', 'temp_int', diag%axesT1, time, &
1704  'Density weighted column integrated potential temperature', 'degC kg m-2', &
1705  cmor_field_name='opottempmint', &
1706  cmor_long_name='integral_wrt_depth_of_product_of_sea_water_density_and_potential_temperature',&
1707  cmor_standard_name='Depth integrated density times potential temperature')
1708 
1709  cs%id_salt_int = register_diag_field('ocean_model', 'salt_int', diag%axesT1, time, &
1710  'Density weighted column integrated salinity', 'psu kg m-2', &
1711  cmor_field_name='somint', &
1712  cmor_long_name='integral_wrt_depth_of_product_of_sea_water_density_and_salinity',&
1713  cmor_standard_name='Depth integrated density times salinity')
1714  endif
1715 
1716  cs%id_col_mass = register_diag_field('ocean_model', 'col_mass', diag%axesT1, time, &
1717  'The column integrated in situ density', 'kg m-2')
1718 
1719  cs%id_col_ht = register_diag_field('ocean_model', 'col_height', diag%axesT1, time, &
1720  'The height of the water column', 'm', conversion=us%Z_to_m)
1721  cs%id_pbo = register_diag_field('ocean_model', 'pbo', diag%axesT1, time, &
1722  long_name='Sea Water Pressure at Sea Floor', standard_name='sea_water_pressure_at_sea_floor', &
1723  units='Pa')
1724 
1725  call set_dependent_diagnostics(mis, adp, cdp, g, cs)
1726 

References mom_error_handler::mom_error(), register_time_deriv(), set_dependent_diagnostics(), and mom_wave_speed::wave_speed_init().

Here is the call graph for this function:

◆ post_surface_dyn_diags()

subroutine, public mom_diagnostics::post_surface_dyn_diags ( type(surface_diag_ids), intent(in)  IDs,
type(ocean_grid_type), intent(in)  G,
type(diag_ctrl), intent(in)  diag,
type(surface), intent(in)  sfc_state,
real, dimension(szi_(g),szj_(g)), intent(in)  ssh 
)

This routine posts diagnostics of various dynamic ocean surface quantities, including velocities, speed and sea surface height, at the time the ocean state is reported back to the caller.

Parameters
[in]idsA structure with the diagnostic IDs.
[in]gocean grid structure
[in]diagregulates diagnostic output
[in]sfc_statestructure describing the ocean surface state
[in]sshTime mean surface height without corrections for ice displacement [m]

Definition at line 1158 of file MOM_diagnostics.F90.

1158  type(surface_diag_IDs), intent(in) :: IDs !< A structure with the diagnostic IDs.
1159  type(ocean_grid_type), intent(in) :: G !< ocean grid structure
1160  type(diag_ctrl), intent(in) :: diag !< regulates diagnostic output
1161  type(surface), intent(in) :: sfc_state !< structure describing the ocean surface state
1162  real, dimension(SZI_(G),SZJ_(G)), &
1163  intent(in) :: ssh !< Time mean surface height without corrections for ice displacement [m]
1164 
1165  ! Local variables
1166  real, dimension(SZI_(G),SZJ_(G)) :: work_2d ! A 2-d work array
1167  integer :: i, j, is, ie, js, je
1168 
1169  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
1170 
1171  if (ids%id_ssh > 0) &
1172  call post_data(ids%id_ssh, ssh, diag, mask=g%mask2dT)
1173 
1174  if (ids%id_ssu > 0) &
1175  call post_data(ids%id_ssu, sfc_state%u, diag, mask=g%mask2dCu)
1176 
1177  if (ids%id_ssv > 0) &
1178  call post_data(ids%id_ssv, sfc_state%v, diag, mask=g%mask2dCv)
1179 
1180  if (ids%id_speed > 0) then
1181  do j=js,je ; do i=is,ie
1182  work_2d(i,j) = sqrt(0.5*(sfc_state%u(i-1,j)**2 + sfc_state%u(i,j)**2) + &
1183  0.5*(sfc_state%v(i,j-1)**2 + sfc_state%v(i,j)**2))
1184  enddo ; enddo
1185  call post_data(ids%id_speed, work_2d, diag, mask=g%mask2dT)
1186  endif
1187 

◆ post_surface_thermo_diags()

subroutine, public mom_diagnostics::post_surface_thermo_diags ( type(surface_diag_ids), intent(in)  IDs,
type(ocean_grid_type), intent(in)  G,
type(verticalgrid_type), intent(in)  GV,
type(unit_scale_type), intent(in)  US,
type(diag_ctrl), intent(in)  diag,
real, intent(in)  dt_int,
type(surface), intent(in)  sfc_state,
type(thermo_var_ptrs), intent(in)  tv,
real, dimension(szi_(g),szj_(g)), intent(in)  ssh,
real, dimension(szi_(g),szj_(g)), intent(in)  ssh_ibc 
)

This routine posts diagnostics of various ocean surface and integrated quantities at the time the ocean state is reported back to the caller.

Parameters
[in]idsA structure with the diagnostic IDs.
[in]gocean grid structure
[in]gvocean vertical grid structure
[in]usA dimensional unit scaling type
[in]diagregulates diagnostic output
[in]dt_inttotal time step associated with these diagnostics [T ~> s].
[in]sfc_statestructure describing the ocean surface state
[in]tvA structure pointing to various thermodynamic variables
[in]sshTime mean surface height without corrections for ice displacement [m]
[in]ssh_ibcTime mean surface height with corrections for ice displacement and the inverse barometer [m]

Definition at line 1195 of file MOM_diagnostics.F90.

1195  type(surface_diag_IDs), intent(in) :: IDs !< A structure with the diagnostic IDs.
1196  type(ocean_grid_type), intent(in) :: G !< ocean grid structure
1197  type(verticalGrid_type), intent(in) :: GV !< ocean vertical grid structure
1198  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
1199  type(diag_ctrl), intent(in) :: diag !< regulates diagnostic output
1200  real, intent(in) :: dt_int !< total time step associated with these diagnostics [T ~> s].
1201  type(surface), intent(in) :: sfc_state !< structure describing the ocean surface state
1202  type(thermo_var_ptrs), intent(in) :: tv !< A structure pointing to various thermodynamic variables
1203  real, dimension(SZI_(G),SZJ_(G)), &
1204  intent(in) :: ssh !< Time mean surface height without corrections for ice displacement [m]
1205  real, dimension(SZI_(G),SZJ_(G)), intent(in) :: ssh_ibc !< Time mean surface height with corrections
1206  !! for ice displacement and the inverse barometer [m]
1207 
1208  real, dimension(SZI_(G),SZJ_(G)) :: work_2d ! A 2-d work array
1209  real, dimension(SZI_(G),SZJ_(G)) :: &
1210  zos ! dynamic sea lev (zero area mean) from inverse-barometer adjusted ssh [m]
1211  real :: I_time_int ! The inverse of the time interval [T-1 ~> s-1].
1212  real :: zos_area_mean, volo, ssh_ga
1213  integer :: i, j, is, ie, js, je
1214 
1215  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
1216 
1217  ! area mean SSH
1218  if (ids%id_ssh_ga > 0) then
1219  ssh_ga = global_area_mean(ssh, g)
1220  call post_data(ids%id_ssh_ga, ssh_ga, diag)
1221  endif
1222 
1223  ! post the dynamic sea level, zos, and zossq.
1224  ! zos is ave_ssh with sea ice inverse barometer removed,
1225  ! and with zero global area mean.
1226  if (ids%id_zos > 0 .or. ids%id_zossq > 0) then
1227  zos(:,:) = 0.0
1228  do j=js,je ; do i=is,ie
1229  zos(i,j) = ssh_ibc(i,j)
1230  enddo ; enddo
1231  zos_area_mean = global_area_mean(zos, g)
1232  do j=js,je ; do i=is,ie
1233  zos(i,j) = zos(i,j) - g%mask2dT(i,j)*zos_area_mean
1234  enddo ; enddo
1235  if (ids%id_zos > 0) call post_data(ids%id_zos, zos, diag, mask=g%mask2dT)
1236  if (ids%id_zossq > 0) then
1237  do j=js,je ; do i=is,ie
1238  work_2d(i,j) = zos(i,j)*zos(i,j)
1239  enddo ; enddo
1240  call post_data(ids%id_zossq, work_2d, diag, mask=g%mask2dT)
1241  endif
1242  endif
1243 
1244  ! post total volume of the liquid ocean
1245  if (ids%id_volo > 0) then
1246  do j=js,je ; do i=is,ie
1247  work_2d(i,j) = g%mask2dT(i,j)*(ssh(i,j) + us%Z_to_m*g%bathyT(i,j))
1248  enddo ; enddo
1249  volo = global_area_integral(work_2d, g)
1250  call post_data(ids%id_volo, volo, diag)
1251  endif
1252 
1253  ! Use Adcroft's rule of reciprocals; it does the right thing here.
1254  i_time_int = 0.0 ; if (dt_int > 0.0) i_time_int = 1.0 / dt_int
1255 
1256  ! post time-averaged rate of frazil formation
1257  if (associated(tv%frazil) .and. (ids%id_fraz > 0)) then
1258  do j=js,je ; do i=is,ie
1259  work_2d(i,j) = tv%frazil(i,j) * i_time_int
1260  enddo ; enddo
1261  call post_data(ids%id_fraz, work_2d, diag, mask=g%mask2dT)
1262  endif
1263 
1264  ! post time-averaged salt deficit
1265  if (associated(tv%salt_deficit) .and. (ids%id_salt_deficit > 0)) then
1266  do j=js,je ; do i=is,ie
1267  work_2d(i,j) = tv%salt_deficit(i,j) * i_time_int
1268  enddo ; enddo
1269  call post_data(ids%id_salt_deficit, work_2d, diag, mask=g%mask2dT)
1270  endif
1271 
1272  ! post temperature of P-E+R
1273  if (associated(tv%TempxPmE) .and. (ids%id_Heat_PmE > 0)) then
1274  do j=js,je ; do i=is,ie
1275  work_2d(i,j) = tv%TempxPmE(i,j) * (tv%C_p * i_time_int)
1276  enddo ; enddo
1277  call post_data(ids%id_Heat_PmE, work_2d, diag, mask=g%mask2dT)
1278  endif
1279 
1280  ! post geothermal heating or internal heat source/sinks
1281  if (associated(tv%internal_heat) .and. (ids%id_intern_heat > 0)) then
1282  do j=js,je ; do i=is,ie
1283  work_2d(i,j) = tv%internal_heat(i,j) * (tv%C_p * i_time_int)
1284  enddo ; enddo
1285  call post_data(ids%id_intern_heat, work_2d, diag, mask=g%mask2dT)
1286  endif
1287 
1288  if (tv%T_is_conT) then
1289  ! Internal T&S variables are conservative temperature & absolute salinity
1290  if (ids%id_sstcon > 0) call post_data(ids%id_sstcon, sfc_state%SST, diag, mask=g%mask2dT)
1291  ! Use TEOS-10 function calls convert T&S diagnostics from conservative temp
1292  ! to potential temperature.
1293  do j=js,je ; do i=is,ie
1294  work_2d(i,j) = gsw_pt_from_ct(sfc_state%SSS(i,j),sfc_state%SST(i,j))
1295  enddo ; enddo
1296  if (ids%id_sst > 0) call post_data(ids%id_sst, work_2d, diag, mask=g%mask2dT)
1297  else
1298  ! Internal T&S variables are potential temperature & practical salinity
1299  if (ids%id_sst > 0) call post_data(ids%id_sst, sfc_state%SST, diag, mask=g%mask2dT)
1300  endif
1301 
1302  if (tv%S_is_absS) then
1303  ! Internal T&S variables are conservative temperature & absolute salinity
1304  if (ids%id_sssabs > 0) call post_data(ids%id_sssabs, sfc_state%SSS, diag, mask=g%mask2dT)
1305  ! Use TEOS-10 function calls convert T&S diagnostics from absolute salinity
1306  ! to practical salinity.
1307  do j=js,je ; do i=is,ie
1308  work_2d(i,j) = gsw_sp_from_sr(sfc_state%SSS(i,j))
1309  enddo ; enddo
1310  if (ids%id_sss > 0) call post_data(ids%id_sss, work_2d, diag, mask=g%mask2dT)
1311  else
1312  ! Internal T&S variables are potential temperature & practical salinity
1313  if (ids%id_sss > 0) call post_data(ids%id_sss, sfc_state%SSS, diag, mask=g%mask2dT)
1314  endif
1315 
1316  if (ids%id_sst_sq > 0) then
1317  do j=js,je ; do i=is,ie
1318  work_2d(i,j) = sfc_state%SST(i,j)*sfc_state%SST(i,j)
1319  enddo ; enddo
1320  call post_data(ids%id_sst_sq, work_2d, diag, mask=g%mask2dT)
1321  endif
1322  if (ids%id_sss_sq > 0) then
1323  do j=js,je ; do i=is,ie
1324  work_2d(i,j) = sfc_state%SSS(i,j)*sfc_state%SSS(i,j)
1325  enddo ; enddo
1326  call post_data(ids%id_sss_sq, work_2d, diag, mask=g%mask2dT)
1327  endif
1328 
1329  call coupler_type_send_data(sfc_state%tr_fields, get_diag_time_end(diag))
1330 

References mom_spatial_means::global_area_integral().

Here is the call graph for this function:

◆ post_transport_diagnostics()

subroutine, public mom_diagnostics::post_transport_diagnostics ( type(ocean_grid_type), intent(inout)  G,
type(verticalgrid_type), intent(in)  GV,
type(unit_scale_type), intent(in)  US,
real, dimension( g %isdb: g %iedb, g %jsd: g %jed, g %ke), intent(in)  uhtr,
real, dimension( g %isd: g %ied, g %jsdb: g %jedb, g %ke), intent(in)  vhtr,
real, dimension( g %isd: g %ied, g %jsd: g %jed, g %ke), intent(in)  h,
type(transport_diag_ids), intent(in)  IDs,
type(diag_grid_storage), intent(inout)  diag_pre_dyn,
type(diag_ctrl), intent(inout)  diag,
real, intent(in)  dt_trans,
type(tracer_registry_type), pointer  Reg 
)

This routine posts diagnostics of the transports, including the subgridscale contributions.

Parameters
[in,out]gocean grid structure
[in]gvocean vertical grid structure
[in]usA dimensional unit scaling type
[in]uhtrAccumulated zonal thickness fluxes used to advect tracers [H L2 ~> m3 or kg]
[in]vhtrAccumulated meridional thickness fluxes used to advect tracers [H L2 ~> m3 or kg]
[in]hThe updated layer thicknesses [H ~> m or kg m-2]
[in]idsA structure with the diagnostic IDs.
[in,out]diag_pre_dynStored grids from before dynamics
[in,out]diagregulates diagnostic output
[in]dt_transtotal time step associated with the transports [T ~> s].
regPointer to the tracer registry

Definition at line 1338 of file MOM_diagnostics.F90.

1338  type(ocean_grid_type), intent(inout) :: G !< ocean grid structure
1339  type(verticalGrid_type), intent(in) :: GV !< ocean vertical grid structure
1340  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
1341  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)), intent(in) :: uhtr !< Accumulated zonal thickness fluxes
1342  !! used to advect tracers [H L2 ~> m3 or kg]
1343  real, dimension(SZI_(G),SZJB_(G),SZK_(G)), intent(in) :: vhtr !< Accumulated meridional thickness fluxes
1344  !! used to advect tracers [H L2 ~> m3 or kg]
1345  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
1346  intent(in) :: h !< The updated layer thicknesses [H ~> m or kg m-2]
1347  type(transport_diag_IDs), intent(in) :: IDs !< A structure with the diagnostic IDs.
1348  type(diag_grid_storage), intent(inout) :: diag_pre_dyn !< Stored grids from before dynamics
1349  type(diag_ctrl), intent(inout) :: diag !< regulates diagnostic output
1350  real, intent(in) :: dt_trans !< total time step associated with the transports [T ~> s].
1351  type(tracer_registry_type), pointer :: Reg !< Pointer to the tracer registry
1352 
1353  ! Local variables
1354  real, dimension(SZIB_(G), SZJ_(G)) :: umo2d ! Diagnostics of integrated mass transport [kg s-1]
1355  real, dimension(SZI_(G), SZJB_(G)) :: vmo2d ! Diagnostics of integrated mass transport [kg s-1]
1356  real, dimension(SZIB_(G), SZJ_(G), SZK_(G)) :: umo ! Diagnostics of layer mass transport [kg s-1]
1357  real, dimension(SZI_(G), SZJB_(G), SZK_(G)) :: vmo ! Diagnostics of layer mass transport [kg s-1]
1358  real, dimension(SZI_(G),SZJ_(G),SZK_(G)) :: h_tend ! Change in layer thickness due to dynamics
1359  ! [H s-1 ~> m s-1 or kg m-2 s-1].
1360  real :: Idt ! The inverse of the time interval [T-1 ~> s-1]
1361  real :: H_to_kg_m2_dt ! A conversion factor from accumulated transports to fluxes
1362  ! [kg L-2 H-1 T-1 ~> kg m-3 s-1 or s-1].
1363  integer :: i, j, k, is, ie, js, je, nz
1364  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
1365 
1366  idt = 1. / dt_trans
1367  h_to_kg_m2_dt = gv%H_to_kg_m2 * us%L_to_m**2 * us%s_to_T * idt
1368 
1369  call diag_save_grids(diag)
1370  call diag_copy_storage_to_diag(diag, diag_pre_dyn)
1371 
1372  if (ids%id_umo_2d > 0) then
1373  umo2d(:,:) = 0.0
1374  do k=1,nz ; do j=js,je ; do i=is-1,ie
1375  umo2d(i,j) = umo2d(i,j) + uhtr(i,j,k) * h_to_kg_m2_dt
1376  enddo ; enddo ; enddo
1377  call post_data(ids%id_umo_2d, umo2d, diag)
1378  endif
1379  if (ids%id_umo > 0) then
1380  ! Convert to kg/s.
1381  do k=1,nz ; do j=js,je ; do i=is-1,ie
1382  umo(i,j,k) = uhtr(i,j,k) * h_to_kg_m2_dt
1383  enddo ; enddo ; enddo
1384  call post_data(ids%id_umo, umo, diag, alt_h = diag_pre_dyn%h_state)
1385  endif
1386  if (ids%id_vmo_2d > 0) then
1387  vmo2d(:,:) = 0.0
1388  do k=1,nz ; do j=js-1,je ; do i=is,ie
1389  vmo2d(i,j) = vmo2d(i,j) + vhtr(i,j,k) * h_to_kg_m2_dt
1390  enddo ; enddo ; enddo
1391  call post_data(ids%id_vmo_2d, vmo2d, diag)
1392  endif
1393  if (ids%id_vmo > 0) then
1394  ! Convert to kg/s.
1395  do k=1,nz ; do j=js-1,je ; do i=is,ie
1396  vmo(i,j,k) = vhtr(i,j,k) * h_to_kg_m2_dt
1397  enddo ; enddo ; enddo
1398  call post_data(ids%id_vmo, vmo, diag, alt_h = diag_pre_dyn%h_state)
1399  endif
1400 
1401  if (ids%id_uhtr > 0) call post_data(ids%id_uhtr, uhtr, diag, alt_h = diag_pre_dyn%h_state)
1402  if (ids%id_vhtr > 0) call post_data(ids%id_vhtr, vhtr, diag, alt_h = diag_pre_dyn%h_state)
1403  if (ids%id_dynamics_h > 0) call post_data(ids%id_dynamics_h, diag_pre_dyn%h_state, diag, &
1404  alt_h = diag_pre_dyn%h_state)
1405  ! Post the change in thicknesses
1406  if (ids%id_dynamics_h_tendency > 0) then
1407  h_tend(:,:,:) = 0.
1408  do k=1,nz ; do j=js,je ; do i=is,ie
1409  h_tend(i,j,k) = (h(i,j,k) - diag_pre_dyn%h_state(i,j,k))*idt
1410  enddo ; enddo ; enddo
1411  call post_data(ids%id_dynamics_h_tendency, h_tend, diag, alt_h = diag_pre_dyn%h_state)
1412  endif
1413 
1414  call post_tracer_transport_diagnostics(g, gv, reg, diag_pre_dyn%h_state, diag)
1415 
1416  call diag_restore_grids(diag)
1417 

References mom_diag_mediator::diag_copy_storage_to_diag(), mom_diag_mediator::diag_restore_grids(), mom_diag_mediator::diag_save_grids(), and mom_tracer_registry::post_tracer_transport_diagnostics().

Here is the call graph for this function:

◆ register_surface_diags()

subroutine, public mom_diagnostics::register_surface_diags ( type(time_type), intent(in)  Time,
type(ocean_grid_type), intent(in)  G,
type(unit_scale_type), intent(in)  US,
type(surface_diag_ids), intent(inout)  IDs,
type(diag_ctrl), intent(inout)  diag,
type(thermo_var_ptrs), intent(in)  tv 
)

Register diagnostics of the surface state and integrated quantities.

Parameters
[in]timecurrent model time
[in]gocean grid structure
[in]usA dimensional unit scaling type
[in,out]idsA structure with the diagnostic IDs.
[in,out]diagregulates diagnostic output
[in]tvA structure pointing to various thermodynamic variables

Definition at line 1732 of file MOM_diagnostics.F90.

1732  type(time_type), intent(in) :: Time !< current model time
1733  type(ocean_grid_type), intent(in) :: G !< ocean grid structure
1734  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
1735  type(surface_diag_IDs), intent(inout) :: IDs !< A structure with the diagnostic IDs.
1736  type(diag_ctrl), intent(inout) :: diag !< regulates diagnostic output
1737  type(thermo_var_ptrs), intent(in) :: tv !< A structure pointing to various thermodynamic variables
1738 
1739  ! Vertically integrated, budget, and surface state diagnostics
1740  ids%id_volo = register_scalar_field('ocean_model', 'volo', time, diag,&
1741  long_name='Total volume of liquid ocean', units='m3', &
1742  standard_name='sea_water_volume')
1743  ids%id_zos = register_diag_field('ocean_model', 'zos', diag%axesT1, time,&
1744  standard_name = 'sea_surface_height_above_geoid', &
1745  long_name= 'Sea surface height above geoid', units='m')
1746  ids%id_zossq = register_diag_field('ocean_model', 'zossq', diag%axesT1, time,&
1747  standard_name='square_of_sea_surface_height_above_geoid', &
1748  long_name='Square of sea surface height above geoid', units='m2')
1749  ids%id_ssh = register_diag_field('ocean_model', 'SSH', diag%axesT1, time, &
1750  'Sea Surface Height', 'm')
1751  ids%id_ssh_ga = register_scalar_field('ocean_model', 'ssh_ga', time, diag,&
1752  long_name='Area averaged sea surface height', units='m', &
1753  standard_name='area_averaged_sea_surface_height')
1754  ids%id_ssu = register_diag_field('ocean_model', 'SSU', diag%axesCu1, time, &
1755  'Sea Surface Zonal Velocity', 'm s-1')
1756  ids%id_ssv = register_diag_field('ocean_model', 'SSV', diag%axesCv1, time, &
1757  'Sea Surface Meridional Velocity', 'm s-1')
1758  ids%id_speed = register_diag_field('ocean_model', 'speed', diag%axesT1, time, &
1759  'Sea Surface Speed', 'm s-1')
1760 
1761  if (associated(tv%T)) then
1762  ids%id_sst = register_diag_field('ocean_model', 'SST', diag%axesT1, time, &
1763  'Sea Surface Temperature', 'degC', cmor_field_name='tos', &
1764  cmor_long_name='Sea Surface Temperature', &
1765  cmor_standard_name='sea_surface_temperature')
1766  ids%id_sst_sq = register_diag_field('ocean_model', 'SST_sq', diag%axesT1, time, &
1767  'Sea Surface Temperature Squared', 'degC2', cmor_field_name='tossq', &
1768  cmor_long_name='Square of Sea Surface Temperature ', &
1769  cmor_standard_name='square_of_sea_surface_temperature')
1770  ids%id_sss = register_diag_field('ocean_model', 'SSS', diag%axesT1, time, &
1771  'Sea Surface Salinity', 'psu', cmor_field_name='sos', &
1772  cmor_long_name='Sea Surface Salinity', &
1773  cmor_standard_name='sea_surface_salinity')
1774  ids%id_sss_sq = register_diag_field('ocean_model', 'SSS_sq', diag%axesT1, time, &
1775  'Sea Surface Salinity Squared', 'psu', cmor_field_name='sossq', &
1776  cmor_long_name='Square of Sea Surface Salinity ', &
1777  cmor_standard_name='square_of_sea_surface_salinity')
1778  if (tv%T_is_conT) then
1779  ids%id_sstcon = register_diag_field('ocean_model', 'conSST', diag%axesT1, time, &
1780  'Sea Surface Conservative Temperature', 'Celsius')
1781  endif
1782  if (tv%S_is_absS) then
1783  ids%id_sssabs = register_diag_field('ocean_model', 'absSSS', diag%axesT1, time, &
1784  'Sea Surface Absolute Salinity', 'g kg-1')
1785  endif
1786  if (associated(tv%frazil)) then
1787  ids%id_fraz = register_diag_field('ocean_model', 'frazil', diag%axesT1, time, &
1788  'Heat from frazil formation', 'W m-2', conversion=us%s_to_T, cmor_field_name='hfsifrazil', &
1789  cmor_standard_name='heat_flux_into_sea_water_due_to_frazil_ice_formation', &
1790  cmor_long_name='Heat Flux into Sea Water due to Frazil Ice Formation')
1791  endif
1792  endif
1793 
1794  ids%id_salt_deficit = register_diag_field('ocean_model', 'salt_deficit', diag%axesT1, time, &
1795  'Salt sink in ocean due to ice flux', &
1796  'psu m-2 s-1', conversion=g%US%R_to_kg_m3*g%US%Z_to_m*us%s_to_T)
1797  ids%id_Heat_PmE = register_diag_field('ocean_model', 'Heat_PmE', diag%axesT1, time, &
1798  'Heat flux into ocean from mass flux into ocean', &
1799  'W m-2', conversion=g%US%R_to_kg_m3*g%US%Z_to_m*us%s_to_T)
1800  ids%id_intern_heat = register_diag_field('ocean_model', 'internal_heat', diag%axesT1, time,&
1801  'Heat flux into ocean from geothermal or other internal sources', 'W m-2', conversion=us%s_to_T)
1802 

◆ register_time_deriv()

subroutine, public mom_diagnostics::register_time_deriv ( integer, dimension(3), intent(in)  lb,
real, dimension(lb(1):,lb(2):,:), target  f_ptr,
real, dimension(lb(1):,lb(2):,:), target  deriv_ptr,
type(diagnostics_cs), pointer  CS 
)

This subroutine registers fields to calculate a diagnostic time derivative.

Parameters
[in]lbLower index bound of f_ptr
f_ptrTime derivative operand
deriv_ptrTime derivative of f_ptr
csControl structure returned by previous call to diagnostics_init.

Definition at line 1085 of file MOM_diagnostics.F90.

1085  integer, intent(in), dimension(3) :: lb !< Lower index bound of f_ptr
1086  real, dimension(lb(1):,lb(2):,:), target :: f_ptr
1087  !< Time derivative operand
1088  real, dimension(lb(1):,lb(2):,:), target :: deriv_ptr
1089  !< Time derivative of f_ptr
1090  type(diagnostics_CS), pointer :: CS !< Control structure returned by previous call to
1091  !! diagnostics_init.
1092 
1093  ! This subroutine registers fields to calculate a diagnostic time derivative.
1094  ! NOTE: Lower bound is required for grid indexing in calculate_derivs().
1095  ! We assume that the vertical axis is 1-indexed.
1096 
1097  integer :: m !< New index of deriv_ptr in CS%deriv
1098  integer :: ub(3) !< Upper index bound of f_ptr, based on shape.
1099 
1100  if (.not.associated(cs)) call mom_error(fatal, &
1101  "register_time_deriv: Module must be initialized before it is used.")
1102 
1103  if (cs%num_time_deriv >= max_fields_) then
1104  call mom_error(warning,"MOM_diagnostics: Attempted to register more than " // &
1105  "MAX_FIELDS_ diagnostic time derivatives via register_time_deriv.")
1106  return
1107  endif
1108 
1109  m = cs%num_time_deriv+1 ; cs%num_time_deriv = m
1110 
1111  ub(:) = lb(:) + shape(f_ptr) - 1
1112 
1113  cs%nlay(m) = size(f_ptr, 3)
1114  cs%deriv(m)%p => deriv_ptr
1115  allocate(cs%prev_val(m)%p(lb(1):ub(1), lb(2):ub(2), cs%nlay(m)))
1116 
1117  cs%var_ptr(m)%p => f_ptr
1118  cs%prev_val(m)%p(:,:,:) = f_ptr(:,:,:)
1119 

References mom_error_handler::mom_error().

Referenced by mom_diagnostics_init(), and set_dependent_diagnostics().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ register_transport_diags()

subroutine, public mom_diagnostics::register_transport_diags ( type(time_type), intent(in)  Time,
type(ocean_grid_type), intent(in)  G,
type(verticalgrid_type), intent(in)  GV,
type(unit_scale_type), intent(in)  US,
type(transport_diag_ids), intent(inout)  IDs,
type(diag_ctrl), intent(inout)  diag 
)

Register certain diagnostics related to transports.

Parameters
[in]timecurrent model time
[in]gocean grid structure
[in]gvocean vertical grid structure
[in]usA dimensional unit scaling type
[in,out]idsA structure with the diagnostic IDs.
[in,out]diagregulates diagnostic output

Definition at line 1807 of file MOM_diagnostics.F90.

1807  type(time_type), intent(in) :: Time !< current model time
1808  type(ocean_grid_type), intent(in) :: G !< ocean grid structure
1809  type(verticalGrid_type), intent(in) :: GV !< ocean vertical grid structure
1810  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
1811  type(transport_diag_IDs), intent(inout) :: IDs !< A structure with the diagnostic IDs.
1812  type(diag_ctrl), intent(inout) :: diag !< regulates diagnostic output
1813 
1814  real :: H_convert
1815  character(len=48) :: thickness_units
1816 
1817  thickness_units = get_thickness_units(gv)
1818  if (gv%Boussinesq) then
1819  h_convert = gv%H_to_m
1820  else
1821  h_convert = gv%H_to_kg_m2
1822  endif
1823 
1824  ! Diagnostics related to tracer and mass transport
1825  ids%id_uhtr = register_diag_field('ocean_model', 'uhtr', diag%axesCuL, time, &
1826  'Accumulated zonal thickness fluxes to advect tracers', 'kg', &
1827  y_cell_method='sum', v_extensive=.true., conversion=h_convert*us%L_to_m**2)
1828  ids%id_vhtr = register_diag_field('ocean_model', 'vhtr', diag%axesCvL, time, &
1829  'Accumulated meridional thickness fluxes to advect tracers', 'kg', &
1830  x_cell_method='sum', v_extensive=.true., conversion=h_convert*us%L_to_m**2)
1831  ids%id_umo = register_diag_field('ocean_model', 'umo', &
1832  diag%axesCuL, time, 'Ocean Mass X Transport', 'kg s-1', &
1833  standard_name='ocean_mass_x_transport', y_cell_method='sum', v_extensive=.true.)
1834  ids%id_vmo = register_diag_field('ocean_model', 'vmo', &
1835  diag%axesCvL, time, 'Ocean Mass Y Transport', 'kg s-1', &
1836  standard_name='ocean_mass_y_transport', x_cell_method='sum', v_extensive=.true.)
1837  ids%id_umo_2d = register_diag_field('ocean_model', 'umo_2d', &
1838  diag%axesCu1, time, 'Ocean Mass X Transport Vertical Sum', 'kg s-1', &
1839  standard_name='ocean_mass_x_transport_vertical_sum', y_cell_method='sum')
1840  ids%id_vmo_2d = register_diag_field('ocean_model', 'vmo_2d', &
1841  diag%axesCv1, time, 'Ocean Mass Y Transport Vertical Sum', 'kg s-1', &
1842  standard_name='ocean_mass_y_transport_vertical_sum', x_cell_method='sum')
1843  ids%id_dynamics_h = register_diag_field('ocean_model','dynamics_h', &
1844  diag%axesTl, time, 'Change in layer thicknesses due to horizontal dynamics', &
1845  'm s-1', v_extensive=.true., conversion=gv%H_to_m)
1846  ids%id_dynamics_h_tendency = register_diag_field('ocean_model','dynamics_h_tendency', &
1847  diag%axesTl, time, 'Change in layer thicknesses due to horizontal dynamics', &
1848  'm s-1', v_extensive=.true., conversion=gv%H_to_m*us%s_to_T)
1849 

References mom_verticalgrid::get_thickness_units().

Here is the call graph for this function:

◆ set_dependent_diagnostics()

subroutine mom_diagnostics::set_dependent_diagnostics ( type(ocean_internal_state), intent(in)  MIS,
type(accel_diag_ptrs), intent(inout)  ADp,
type(cont_diag_ptrs), intent(inout)  CDp,
type(ocean_grid_type), intent(in)  G,
type(diagnostics_cs), pointer  CS 
)
private

This subroutine sets up diagnostics upon which other diagnostics depend.

Parameters
[in]misFor "MOM Internal State" a set of pointers to the fields and accelerations making up ocean internal physical state.
[in,out]adpStructure pointing to accelerations in momentum equation.
[in,out]cdpStructure pointing to terms in continuity equation.
[in]gThe ocean's grid structure.
csPointer to the control structure for this module.

Definition at line 2029 of file MOM_diagnostics.F90.

2029  type(ocean_internal_state), intent(in) :: MIS !< For "MOM Internal State" a set of pointers to
2030  !! the fields and accelerations making up ocean
2031  !! internal physical state.
2032  type(accel_diag_ptrs), intent(inout) :: ADp !< Structure pointing to accelerations in
2033  !! momentum equation.
2034  type(cont_diag_ptrs), intent(inout) :: CDp !< Structure pointing to terms in continuity
2035  !! equation.
2036  type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure.
2037  type(diagnostics_CS), pointer :: CS !< Pointer to the control structure for this
2038  !! module.
2039 
2040 ! This subroutine sets up diagnostics upon which other diagnostics depend.
2041  integer :: isd, ied, jsd, jed, IsdB, IedB, JsdB, JedB, nz
2042  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed ; nz = g%ke
2043  isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
2044 
2045  if (associated(cs%dKE_dt) .or. associated(cs%PE_to_KE) .or. &
2046  associated(cs%KE_CorAdv) .or. associated(cs%KE_adv) .or. &
2047  associated(cs%KE_visc) .or. associated(cs%KE_horvisc) .or. &
2048  associated(cs%KE_dia)) &
2049  call safe_alloc_ptr(cs%KE,isd,ied,jsd,jed,nz)
2050 
2051  if (associated(cs%dKE_dt)) then
2052  if (.not.associated(cs%du_dt)) then
2053  call safe_alloc_ptr(cs%du_dt,isdb,iedb,jsd,jed,nz)
2054  call register_time_deriv(lbound(mis%u), mis%u, cs%du_dt, cs)
2055  endif
2056  if (.not.associated(cs%dv_dt)) then
2057  call safe_alloc_ptr(cs%dv_dt,isd,ied,jsdb,jedb,nz)
2058  call register_time_deriv(lbound(mis%v), mis%v, cs%dv_dt, cs)
2059  endif
2060  if (.not.associated(cs%dh_dt)) then
2061  call safe_alloc_ptr(cs%dh_dt,isd,ied,jsd,jed,nz)
2062  call register_time_deriv(lbound(mis%h), mis%h, cs%dh_dt, cs)
2063  endif
2064  endif
2065 
2066  if (associated(cs%KE_adv)) then
2067  call safe_alloc_ptr(adp%gradKEu,isdb,iedb,jsd,jed,nz)
2068  call safe_alloc_ptr(adp%gradKEv,isd,ied,jsdb,jedb,nz)
2069  endif
2070 
2071  if (associated(cs%KE_visc)) then
2072  call safe_alloc_ptr(adp%du_dt_visc,isdb,iedb,jsd,jed,nz)
2073  call safe_alloc_ptr(adp%dv_dt_visc,isd,ied,jsdb,jedb,nz)
2074  endif
2075 
2076  if (associated(cs%KE_dia)) then
2077  call safe_alloc_ptr(adp%du_dt_dia,isdb,iedb,jsd,jed,nz)
2078  call safe_alloc_ptr(adp%dv_dt_dia,isd,ied,jsdb,jedb,nz)
2079  endif
2080 
2081  if (associated(cs%uhGM_Rlay)) call safe_alloc_ptr(cdp%uhGM,isdb,iedb,jsd,jed,nz)
2082  if (associated(cs%vhGM_Rlay)) call safe_alloc_ptr(cdp%vhGM,isd,ied,jsdb,jedb,nz)
2083 

References register_time_deriv().

Referenced by mom_diagnostics_init().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ write_static_fields()

subroutine, public mom_diagnostics::write_static_fields ( type(ocean_grid_type), intent(in)  G,
type(verticalgrid_type), intent(in)  GV,
type(unit_scale_type), intent(in)  US,
type(thermo_var_ptrs), intent(in)  tv,
type(diag_ctrl), intent(inout), target  diag 
)

Offers the static fields in the ocean grid type for output via the diag_manager.

Parameters
[in]gocean grid structure
[in]gvocean vertical grid structure
[in]usA dimensional unit scaling type
[in]tvA structure pointing to various thermodynamic variables
[in,out]diagregulates diagnostic output

Definition at line 1854 of file MOM_diagnostics.F90.

1854  type(ocean_grid_type), intent(in) :: G !< ocean grid structure
1855  type(verticalGrid_type), intent(in) :: GV !< ocean vertical grid structure
1856  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
1857  type(thermo_var_ptrs), intent(in) :: tv !< A structure pointing to various thermodynamic variables
1858  type(diag_ctrl), target, intent(inout) :: diag !< regulates diagnostic output
1859 
1860  ! Local variables
1861  integer :: id
1862  logical :: use_temperature
1863 
1864  id = register_static_field('ocean_model', 'geolat', diag%axesT1, &
1865  'Latitude of tracer (T) points', 'degrees_north')
1866  if (id > 0) call post_data(id, g%geoLatT, diag, .true.)
1867 
1868  id = register_static_field('ocean_model', 'geolon', diag%axesT1, &
1869  'Longitude of tracer (T) points', 'degrees_east')
1870  if (id > 0) call post_data(id, g%geoLonT, diag, .true.)
1871 
1872  id = register_static_field('ocean_model', 'geolat_c', diag%axesB1, &
1873  'Latitude of corner (Bu) points', 'degrees_north', interp_method='none')
1874  if (id > 0) call post_data(id, g%geoLatBu, diag, .true.)
1875 
1876  id = register_static_field('ocean_model', 'geolon_c', diag%axesB1, &
1877  'Longitude of corner (Bu) points', 'degrees_east', interp_method='none')
1878  if (id > 0) call post_data(id, g%geoLonBu, diag, .true.)
1879 
1880  id = register_static_field('ocean_model', 'geolat_v', diag%axesCv1, &
1881  'Latitude of meridional velocity (Cv) points', 'degrees_north', interp_method='none')
1882  if (id > 0) call post_data(id, g%geoLatCv, diag, .true.)
1883 
1884  id = register_static_field('ocean_model', 'geolon_v', diag%axesCv1, &
1885  'Longitude of meridional velocity (Cv) points', 'degrees_east', interp_method='none')
1886  if (id > 0) call post_data(id, g%geoLonCv, diag, .true.)
1887 
1888  id = register_static_field('ocean_model', 'geolat_u', diag%axesCu1, &
1889  'Latitude of zonal velocity (Cu) points', 'degrees_north', interp_method='none')
1890  if (id > 0) call post_data(id, g%geoLatCu, diag, .true.)
1891 
1892  id = register_static_field('ocean_model', 'geolon_u', diag%axesCu1, &
1893  'Longitude of zonal velocity (Cu) points', 'degrees_east', interp_method='none')
1894  if (id > 0) call post_data(id, g%geoLonCu, diag, .true.)
1895 
1896  id = register_static_field('ocean_model', 'area_t', diag%axesT1, &
1897  'Surface area of tracer (T) cells', 'm2', conversion=us%L_to_m**2, &
1898  cmor_field_name='areacello', cmor_standard_name='cell_area', &
1899  cmor_long_name='Ocean Grid-Cell Area', &
1900  x_cell_method='sum', y_cell_method='sum', area_cell_method='sum')
1901  if (id > 0) then
1902  call post_data(id, g%areaT, diag, .true.)
1903  call diag_register_area_ids(diag, id_area_t=id)
1904  endif
1905 
1906  id = register_static_field('ocean_model', 'area_u', diag%axesCu1, &
1907  'Surface area of x-direction flow (U) cells', 'm2', conversion=us%L_to_m**2, &
1908  cmor_field_name='areacello_cu', cmor_standard_name='cell_area', &
1909  cmor_long_name='Ocean Grid-Cell Area', &
1910  x_cell_method='sum', y_cell_method='sum', area_cell_method='sum')
1911  if (id > 0) call post_data(id, g%areaCu, diag, .true.)
1912 
1913  id = register_static_field('ocean_model', 'area_v', diag%axesCv1, &
1914  'Surface area of y-direction flow (V) cells', 'm2', conversion=us%L_to_m**2, &
1915  cmor_field_name='areacello_cv', cmor_standard_name='cell_area', &
1916  cmor_long_name='Ocean Grid-Cell Area', &
1917  x_cell_method='sum', y_cell_method='sum', area_cell_method='sum')
1918  if (id > 0) call post_data(id, g%areaCv, diag, .true.)
1919 
1920  id = register_static_field('ocean_model', 'area_q', diag%axesB1, &
1921  'Surface area of B-grid flow (Q) cells', 'm2', conversion=us%L_to_m**2, &
1922  cmor_field_name='areacello_bu', cmor_standard_name='cell_area', &
1923  cmor_long_name='Ocean Grid-Cell Area', &
1924  x_cell_method='sum', y_cell_method='sum', area_cell_method='sum')
1925  if (id > 0) call post_data(id, g%areaBu, diag, .true.)
1926 
1927  id = register_static_field('ocean_model', 'depth_ocean', diag%axesT1, &
1928  'Depth of the ocean at tracer points', 'm', &
1929  standard_name='sea_floor_depth_below_geoid', &
1930  cmor_field_name='deptho', cmor_long_name='Sea Floor Depth', &
1931  cmor_standard_name='sea_floor_depth_below_geoid', area=diag%axesT1%id_area, &
1932  x_cell_method='mean', y_cell_method='mean', area_cell_method='mean', &
1933  conversion=us%Z_to_m)
1934  if (id > 0) call post_data(id, g%bathyT, diag, .true., mask=g%mask2dT)
1935 
1936  id = register_static_field('ocean_model', 'wet', diag%axesT1, &
1937  '0 if land, 1 if ocean at tracer points', 'none', area=diag%axesT1%id_area)
1938  if (id > 0) call post_data(id, g%mask2dT, diag, .true.)
1939 
1940  id = register_static_field('ocean_model', 'wet_c', diag%axesB1, &
1941  '0 if land, 1 if ocean at corner (Bu) points', 'none', interp_method='none')
1942  if (id > 0) call post_data(id, g%mask2dBu, diag, .true.)
1943 
1944  id = register_static_field('ocean_model', 'wet_u', diag%axesCu1, &
1945  '0 if land, 1 if ocean at zonal velocity (Cu) points', 'none', interp_method='none')
1946  if (id > 0) call post_data(id, g%mask2dCu, diag, .true.)
1947 
1948  id = register_static_field('ocean_model', 'wet_v', diag%axesCv1, &
1949  '0 if land, 1 if ocean at meridional velocity (Cv) points', 'none', interp_method='none')
1950  if (id > 0) call post_data(id, g%mask2dCv, diag, .true.)
1951 
1952  id = register_static_field('ocean_model', 'Coriolis', diag%axesB1, &
1953  'Coriolis parameter at corner (Bu) points', 's-1', interp_method='none', conversion=us%s_to_T)
1954  if (id > 0) call post_data(id, g%CoriolisBu, diag, .true.)
1955 
1956  id = register_static_field('ocean_model', 'dxt', diag%axesT1, &
1957  'Delta(x) at thickness/tracer points (meter)', 'm', interp_method='none', conversion=us%L_to_m)
1958  if (id > 0) call post_data(id, g%dxT, diag, .true.)
1959 
1960  id = register_static_field('ocean_model', 'dyt', diag%axesT1, &
1961  'Delta(y) at thickness/tracer points (meter)', 'm', interp_method='none', conversion=us%L_to_m)
1962  if (id > 0) call post_data(id, g%dyT, diag, .true.)
1963 
1964  id = register_static_field('ocean_model', 'dxCu', diag%axesCu1, &
1965  'Delta(x) at u points (meter)', 'm', interp_method='none', conversion=us%L_to_m)
1966  if (id > 0) call post_data(id, g%dxCu, diag, .true.)
1967 
1968  id = register_static_field('ocean_model', 'dyCu', diag%axesCu1, &
1969  'Delta(y) at u points (meter)', 'm', interp_method='none', conversion=us%L_to_m)
1970  if (id > 0) call post_data(id, g%dyCu, diag, .true.)
1971 
1972  id = register_static_field('ocean_model', 'dxCv', diag%axesCv1, &
1973  'Delta(x) at v points (meter)', 'm', interp_method='none', conversion=us%L_to_m)
1974  if (id > 0) call post_data(id, g%dxCv, diag, .true.)
1975 
1976  id = register_static_field('ocean_model', 'dyCv', diag%axesCv1, &
1977  'Delta(y) at v points (meter)', 'm', interp_method='none', conversion=us%L_to_m)
1978  if (id > 0) call post_data(id, g%dyCv, diag, .true.)
1979 
1980  id = register_static_field('ocean_model', 'dyCuo', diag%axesCu1, &
1981  'Open meridional grid spacing at u points (meter)', 'm', interp_method='none', conversion=us%L_to_m)
1982  if (id > 0) call post_data(id, g%dy_Cu, diag, .true.)
1983 
1984  id = register_static_field('ocean_model', 'dxCvo', diag%axesCv1, &
1985  'Open zonal grid spacing at v points (meter)', 'm', interp_method='none', conversion=us%L_to_m)
1986  if (id > 0) call post_data(id, g%dx_Cv, diag, .true.)
1987 
1988  id = register_static_field('ocean_model', 'sin_rot', diag%axesT1, &
1989  'sine of the clockwise angle of the ocean grid north to true north', 'none')
1990  if (id > 0) call post_data(id, g%sin_rot, diag, .true.)
1991 
1992  id = register_static_field('ocean_model', 'cos_rot', diag%axesT1, &
1993  'cosine of the clockwise angle of the ocean grid north to true north', 'none')
1994  if (id > 0) call post_data(id, g%cos_rot, diag, .true.)
1995 
1996 
1997  ! This static diagnostic is from CF 1.8, and is the fraction of a cell
1998  ! covered by ocean, given as a percentage (poorly named).
1999  id = register_static_field('ocean_model', 'area_t_percent', diag%axesT1, &
2000  'Percentage of cell area covered by ocean', '%', &
2001  cmor_field_name='sftof', cmor_standard_name='SeaAreaFraction', &
2002  cmor_long_name='Sea Area Fraction', &
2003  x_cell_method='mean', y_cell_method='mean', area_cell_method='mean', &
2004  conversion=100.0)
2005  if (id > 0) call post_data(id, g%mask2dT, diag, .true.)
2006 
2007 
2008  id = register_static_field('ocean_model','Rho_0', diag%axesNull, &
2009  'mean ocean density used with the Boussinesq approximation', &
2010  'kg m-3', cmor_field_name='rhozero', conversion=us%R_to_kg_m3, &
2011  cmor_standard_name='reference_sea_water_density_for_boussinesq_approximation', &
2012  cmor_long_name='reference sea water density for boussinesq approximation')
2013  if (id > 0) call post_data(id, gv%Rho0, diag, .true.)
2014 
2015  use_temperature = associated(tv%T)
2016  if (use_temperature) then
2017  id = register_static_field('ocean_model','C_p', diag%axesNull, &
2018  'heat capacity of sea water', 'J kg-1 K-1', cmor_field_name='cpocean', &
2019  cmor_standard_name='specific_heat_capacity_of_sea_water', &
2020  cmor_long_name='specific_heat_capacity_of_sea_water')
2021  if (id > 0) call post_data(id, tv%C_p, diag, .true.)
2022  endif
2023