19 use mom_eos,
only : gsw_sp_from_sr, gsw_pt_from_ct
32 use coupler_types_mod,
only : coupler_type_send_data
34 implicit none ;
private
36 #include <MOM_memory.h>
51 real :: mono_n2_column_fraction = 0.
54 real :: mono_n2_depth = -1.
63 real,
pointer,
dimension(:,:,:) :: &
68 real,
pointer,
dimension(:,:,:) :: &
74 real,
pointer,
dimension(:,:,:) :: h_rlay => null()
76 real,
pointer,
dimension(:,:,:) :: uh_rlay => null()
78 real,
pointer,
dimension(:,:,:) :: vh_rlay => null()
80 real,
pointer,
dimension(:,:,:) :: uhgm_rlay => null()
82 real,
pointer,
dimension(:,:,:) :: vhgm_rlay => null()
86 real,
pointer,
dimension(:,:) :: &
90 cfl_cg1_x => null(), &
94 real,
pointer,
dimension(:,:,:) :: &
98 ke_coradv => null(), &
105 ke_horvisc => null(), &
109 integer :: id_u = -1, id_v = -1, id_h = -1
110 integer :: id_e = -1, id_e_d = -1
111 integer :: id_du_dt = -1, id_dv_dt = -1
112 integer :: id_col_ht = -1, id_dh_dt = -1
113 integer :: id_ke = -1, id_dkedt = -1
114 integer :: id_pe_to_ke = -1, id_ke_coradv = -1
115 integer :: id_ke_adv = -1, id_ke_visc = -1
116 integer :: id_ke_horvisc = -1, id_ke_dia = -1
117 integer :: id_uh_rlay = -1, id_vh_rlay = -1
118 integer :: id_uhgm_rlay = -1, id_vhgm_rlay = -1
119 integer :: id_h_rlay = -1, id_rd1 = -1
120 integer :: id_rml = -1, id_rcv = -1
121 integer :: id_cg1 = -1, id_cfl_cg1 = -1
122 integer :: id_cfl_cg1_x = -1, id_cfl_cg1_y = -1
123 integer :: id_cg_ebt = -1, id_rd_ebt = -1
124 integer :: id_p_ebt = -1
125 integer :: id_temp_int = -1, id_salt_int = -1
126 integer :: id_mass_wt = -1, id_col_mass = -1
127 integer :: id_masscello = -1, id_masso = -1
128 integer :: id_volcello = -1
129 integer :: id_tpot = -1, id_sprac = -1
130 integer :: id_tob = -1, id_sob = -1
131 integer :: id_thetaoga = -1, id_soga = -1
132 integer :: id_sosga = -1, id_tosga = -1
133 integer :: id_temp_layer_ave = -1, id_salt_layer_ave = -1
134 integer :: id_pbo = -1
135 integer :: id_thkcello = -1, id_rhoinsitu = -1
136 integer :: id_rhopot0 = -1, id_rhopot2 = -1
137 integer :: id_h_pre_sync = -1
141 type(
p3d) :: var_ptr(max_fields_)
143 type(
p3d) :: deriv(max_fields_)
144 type(
p3d) :: prev_val(max_fields_)
147 integer :: nlay(max_fields_)
148 integer :: num_time_deriv = 0
150 type(group_pass_type) :: pass_ke_uv
159 integer :: id_zos = -1, id_zossq = -1
160 integer :: id_volo = -1, id_speed = -1
161 integer :: id_ssh = -1, id_ssh_ga = -1
162 integer :: id_sst = -1, id_sst_sq = -1, id_sstcon = -1
163 integer :: id_sss = -1, id_sss_sq = -1, id_sssabs = -1
164 integer :: id_ssu = -1, id_ssv = -1
167 integer :: id_fraz = -1
168 integer :: id_salt_deficit = -1
169 integer :: id_heat_pme = -1
170 integer :: id_intern_heat = -1
178 integer :: id_uhtr = -1, id_umo = -1, id_umo_2d = -1
179 integer :: id_vhtr = -1, id_vmo = -1, id_vmo_2d = -1
180 integer :: id_dynamics_h = -1, id_dynamics_h_tendency = -1
187 dt, diag_pre_sync, G, GV, US, CS, eta_bt)
191 real,
dimension(SZIB_(G),SZJ_(G),SZK_(G)), &
193 real,
dimension(SZI_(G),SZJB_(G),SZK_(G)), &
195 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)), &
197 real,
dimension(SZIB_(G),SZJ_(G),SZK_(G)), &
200 real,
dimension(SZI_(G),SZJB_(G),SZK_(G)), &
209 real,
dimension(:,:),
pointer :: p_surf
212 real,
intent(in) :: dt
217 real,
dimension(SZI_(G),SZJ_(G)), &
218 optional,
intent(in) :: eta_bt
224 integer i, j, k, is, ie, js, je, isq, ieq, jsq, jeq, nz, nkmb
226 real :: rcv(szi_(g),szj_(g),szk_(g))
227 real :: work_3d(szi_(g),szj_(g),szk_(g))
228 real :: work_2d(szi_(g),szj_(g))
231 real :: surface_field(szi_(g),szj_(g))
232 real :: pressure_1d(szi_(g))
237 real :: absurdly_small_freq2
241 real,
dimension(SZK_(G)) :: temp_layer_ave, salt_layer_ave
242 real :: thetaoga, soga, masso, tosga, sosga
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
249 absurdly_small_freq2 = 1e-34*us%T_to_s**2
252 "calculate_diagnostic_fields: Module must be initialized before used.")
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)
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)
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)
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)
271 call calculate_energy_diagnostics(u, v, h, uh, vh, adp, cdp, g, gv, us, cs)
280 if (nkmb==0 .and. nz > 1) nkmb = nz
282 if (cs%id_u > 0)
call post_data(cs%id_u, u, cs%diag)
284 if (cs%id_v > 0)
call post_data(cs%id_v, v, cs%diag)
286 if (cs%id_h > 0)
call post_data(cs%id_h, h, cs%diag)
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)
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
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
305 if (cs%id_e_D > 0)
call post_data(cs%id_e_D, cs%e_D, cs%diag)
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)
317 if (cs%id_masso > 0)
then
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
323 call post_data(cs%id_masso, masso, cs%diag)
327 if (cs%id_thkcello>0 .or. cs%id_volcello>0)
then
328 if (gv%Boussinesq)
then
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)
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)
337 if (cs%id_volcello > 0)
then
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)
345 if (
associated(p_surf))
then
347 pressure_1d(i) = p_surf(i,j)
356 pressure_1d(i) = pressure_1d(i) + 0.5*gv%H_to_Pa*h(i,j,k)
360 work_3d(:,j,k), is, ie-is+1, tv%eqn_of_state)
362 work_3d(i,j,k) = (gv%H_to_kg_m2*h(i,j,k)) / work_3d(i,j,k)
365 pressure_1d(i) = pressure_1d(i) + 0.5*gv%H_to_Pa*h(i,j,k)
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
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)
380 if (tv%T_is_conT)
then
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)
393 if (cs%id_tob > 0)
call post_data(cs%id_tob, tv%T(:,:,nz), cs%diag, mask=g%mask2dT)
397 if (tv%S_is_absS)
then
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)
410 if (cs%id_sob > 0)
call post_data(cs%id_sob, tv%S(:,:,nz), cs%diag, mask=g%mask2dT)
414 if (cs%id_thetaoga>0)
then
416 call post_data(cs%id_thetaoga, thetaoga, cs%diag)
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)
424 tosga = global_area_mean(surface_field, g)
425 call post_data(cs%id_tosga, tosga, cs%diag)
429 if (cs%id_soga>0)
then
431 call post_data(cs%id_soga, soga, cs%diag)
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)
439 sosga = global_area_mean(surface_field, g)
440 call post_data(cs%id_sosga, sosga, cs%diag)
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)
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)
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
461 if (
associated(tv%eqn_of_state))
then
462 pressure_1d(:) = tv%P_Ref
464 do k=1,nz ;
do j=js-1,je+1
466 rcv(:,j,k), is-1, ie-is+3, tv%eqn_of_state, scale=us%kg_m3_to_R)
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
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)
476 if (
associated(cs%h_Rlay))
then
481 do k=1,nkmb ;
do i=is,ie
482 cs%h_Rlay(i,j,k) = 0.0
484 do k=nkmb+1,nz ;
do i=is,ie
485 cs%h_Rlay(i,j,k) = h(i,j,k)
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
494 if (cs%id_h_Rlay > 0)
call post_data(cs%id_h_Rlay, cs%h_Rlay, cs%diag)
497 if (
associated(cs%uh_Rlay))
then
502 do k=1,nkmb ;
do i=isq,ieq
503 cs%uh_Rlay(i,j,k) = 0.0
505 do k=nkmb+1,nz ;
do i=isq,ieq
506 cs%uh_Rlay(i,j,k) = uh(i,j,k)
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
516 if (cs%id_uh_Rlay > 0)
call post_data(cs%id_uh_Rlay, cs%uh_Rlay, cs%diag)
519 if (
associated(cs%vh_Rlay))
then
524 do k=1,nkmb ;
do i=is,ie
525 cs%vh_Rlay(i,j,k) = 0.0
527 do k=nkmb+1,nz ;
do i=is,ie
528 cs%vh_Rlay(i,j,k) = vh(i,j,k)
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
537 if (cs%id_vh_Rlay > 0)
call post_data(cs%id_vh_Rlay, cs%vh_Rlay, cs%diag)
540 if (
associated(cs%uhGM_Rlay) .and.
associated(cdp%uhGM))
then
545 do k=1,nkmb ;
do i=isq,ieq
546 cs%uhGM_Rlay(i,j,k) = 0.0
548 do k=nkmb+1,nz ;
do i=isq,ieq
549 cs%uhGM_Rlay(i,j,k) = cdp%uhGM(i,j,k)
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
558 if (cs%id_uhGM_Rlay > 0)
call post_data(cs%id_uhGM_Rlay, cs%uhGM_Rlay, cs%diag)
561 if (
associated(cs%vhGM_Rlay) .and.
associated(cdp%vhGM))
then
566 do k=1,nkmb ;
do i=is,ie
567 cs%vhGM_Rlay(i,j,k) = 0.0
569 do k=nkmb+1,nz ;
do i=is,ie
570 cs%vhGM_Rlay(i,j,k) = cdp%vhGM(i,j,k)
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
579 if (cs%id_vhGM_Rlay > 0)
call post_data(cs%id_vhGM_Rlay, cs%vhGM_Rlay, cs%diag)
583 if (
associated(tv%eqn_of_state))
then
584 if (cs%id_rhopot0 > 0)
then
587 do k=1,nz ;
do j=js,je
589 rcv(:,j,k),is,ie-is+1, tv%eqn_of_state)
591 if (cs%id_rhopot0 > 0)
call post_data(cs%id_rhopot0, rcv, cs%diag)
593 if (cs%id_rhopot2 > 0)
then
594 pressure_1d(:) = 2.e7
596 do k=1,nz ;
do j=js,je
598 rcv(:,j,k),is,ie-is+1, tv%eqn_of_state)
600 if (cs%id_rhopot2 > 0)
call post_data(cs%id_rhopot2, rcv, cs%diag)
602 if (cs%id_rhoinsitu > 0)
then
607 pressure_1d(:) = pressure_1d(:) + 0.5 * h(:,j,k) * gv%H_to_Pa
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
613 if (cs%id_rhoinsitu > 0)
call post_data(cs%id_rhoinsitu, rcv, cs%diag)
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
623 do j=js,je ;
do i=is,ie
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)
636 call post_data(cs%id_Rd1, cs%Rd1, cs%diag)
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))
642 call post_data(cs%id_cfl_cg1, cs%cfl_cg1, cs%diag)
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)
648 call post_data(cs%id_cfl_cg1_x, cs%cfl_cg1_x, cs%diag)
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)
654 call post_data(cs%id_cfl_cg1_y, cs%cfl_cg1_y, cs%diag)
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)
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)
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
671 do j=js,je ;
do i=is,ie
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)
684 call post_data(cs%id_Rd_ebt, cs%Rd1, cs%diag)
695 real,
dimension(:), &
697 real,
intent(in) :: R_in
698 integer,
intent(inout) :: k
700 integer,
intent(in) :: nz
701 real,
intent(out) :: wt
702 real,
intent(out) :: wt_p
709 integer :: k_upper, k_lower, k_new, inc
712 if ((k < 1) .or. (k > nz)) k = nz/2
714 k_upper = k ; k_lower = k ; inc = 1
715 if (r_in < rlist(k))
then
717 k_lower = max(k_lower-inc, 1)
718 if ((k_lower == 1) .or. (r_in >= rlist(k_lower)))
exit
724 k_upper = min(k_upper+inc, nz)
725 if ((k_upper == nz) .or. (r_in < rlist(k_upper)))
exit
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
737 if (k_upper <= k_lower+1)
exit
738 k_new = (k_upper + k_lower) / 2
739 if (r_in < rlist(k_new))
then
751 wt = (rlist(k_upper) - r_in) / (rlist(k_upper) - rlist(k_lower))
765 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)), &
769 real,
dimension(:,:),
pointer :: p_surf
775 real,
dimension(SZI_(G), SZJ_(G)) :: &
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
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)
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)
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)
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)
823 call post_data(cs%id_col_ht, z_bot, cs%diag)
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
833 do j=g%jscB,g%jecB+1 ;
do i=g%iscB,g%iecB+1
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)
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
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
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
858 if (cs%id_col_mass > 0)
then
859 call post_data(cs%id_col_mass, mass, cs%diag)
861 if (cs%id_pbo > 0)
then
862 do j=js,je ;
do i=is,ie ; btm_pres(i,j) = 0.0 ;
enddo ;
enddo
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)
872 call post_data(cs%id_pbo, btm_pres, cs%diag)
879 subroutine calculate_energy_diagnostics(u, v, h, uh, vh, ADp, CDp, G, GV, US, CS)
882 real,
dimension(SZIB_(G),SZJ_(G),SZK_(G)), &
884 real,
dimension(SZI_(G),SZJB_(G),SZK_(G)), &
886 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)), &
888 real,
dimension(SZIB_(G),SZJ_(G),SZK_(G)), &
891 real,
dimension(SZI_(G),SZJB_(G),SZK_(G)), &
903 real :: KE_u(SZIB_(G),SZJ_(G))
904 real :: KE_v(SZI_(G),SZJB_(G))
905 real :: KE_h(SZI_(G),SZJ_(G))
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
911 do j=js-1,je ;
do i=is-1,ie
912 ke_u(i,j) = 0.0 ; ke_v(i,j) = 0.0
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
922 enddo ;
enddo ;
enddo
923 if (cs%id_KE > 0)
call post_data(cs%id_KE, cs%KE, cs%diag)
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
934 if (
associated(cs%dKE_dt))
then
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)
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)
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)
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))
952 if (cs%id_dKEdt > 0)
call post_data(cs%id_dKEdt, cs%dKE_dt, cs%diag)
955 if (
associated(cs%PE_to_KE))
then
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)
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)
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))
970 if (cs%id_PE_to_KE > 0)
call post_data(cs%id_PE_to_KE, cs%PE_to_KE, cs%diag)
973 if (
associated(cs%KE_CorAdv))
then
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)
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)
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))
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))
992 if (cs%id_KE_Coradv > 0)
call post_data(cs%id_KE_Coradv, cs%KE_Coradv, cs%diag)
995 if (
associated(cs%KE_adv))
then
999 ke_u(:,:) = 0. ; ke_v(:,:) = 0.
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)
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)
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))
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))
1020 if (cs%id_KE_adv > 0)
call post_data(cs%id_KE_adv, cs%KE_adv, cs%diag)
1023 if (
associated(cs%KE_visc))
then
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)
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)
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))
1038 if (cs%id_KE_visc > 0)
call post_data(cs%id_KE_visc, cs%KE_visc, cs%diag)
1041 if (
associated(cs%KE_horvisc))
then
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)
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)
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))
1056 if (cs%id_KE_horvisc > 0)
call post_data(cs%id_KE_horvisc, cs%KE_horvisc, cs%diag)
1059 if (
associated(cs%KE_dia))
then
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)
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)
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)))
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))
1078 if (cs%id_KE_dia > 0)
call post_data(cs%id_KE_dia, cs%KE_dia, cs%diag)
1085 integer,
intent(in),
dimension(3) :: lb
1086 real,
dimension(lb(1):,lb(2):,:),
target :: f_ptr
1088 real,
dimension(lb(1):,lb(2):,:),
target :: deriv_ptr
1100 if (.not.
associated(cs))
call mom_error(fatal, &
1101 "register_time_deriv: Module must be initialized before it is used.")
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.")
1109 m = cs%num_time_deriv+1 ; cs%num_time_deriv = m
1111 ub(:) = lb(:) + shape(f_ptr) - 1
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)))
1117 cs%var_ptr(m)%p => f_ptr
1118 cs%prev_val(m)%p(:,:,:) = f_ptr(:,:,:)
1124 real,
intent(in) :: dt
1131 integer :: i, j, k, m
1133 if (dt > 0.0)
then ; idt = 1.0/dt
1134 else ;
return ;
endif
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
1161 type(
surface),
intent(in) :: sfc_state
1162 real,
dimension(SZI_(G),SZJ_(G)), &
1166 real,
dimension(SZI_(G),SZJ_(G)) :: work_2d
1167 integer :: i, j, is, ie, js, je
1169 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
1171 if (ids%id_ssh > 0) &
1172 call post_data(ids%id_ssh, ssh, diag, mask=g%mask2dT)
1174 if (ids%id_ssu > 0) &
1175 call post_data(ids%id_ssu, sfc_state%u, diag, mask=g%mask2dCu)
1177 if (ids%id_ssv > 0) &
1178 call post_data(ids%id_ssv, sfc_state%v, diag, mask=g%mask2dCv)
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))
1185 call post_data(ids%id_speed, work_2d, diag, mask=g%mask2dT)
1200 real,
intent(in) :: dt_int
1201 type(
surface),
intent(in) :: sfc_state
1203 real,
dimension(SZI_(G),SZJ_(G)), &
1205 real,
dimension(SZI_(G),SZJ_(G)),
intent(in) :: ssh_ibc
1208 real,
dimension(SZI_(G),SZJ_(G)) :: work_2d
1209 real,
dimension(SZI_(G),SZJ_(G)) :: &
1212 real :: zos_area_mean, volo, ssh_ga
1213 integer :: i, j, is, ie, js, je
1215 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
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)
1226 if (ids%id_zos > 0 .or. ids%id_zossq > 0)
then
1228 do j=js,je ;
do i=is,ie
1229 zos(i,j) = ssh_ibc(i,j)
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
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)
1240 call post_data(ids%id_zossq, work_2d, diag, mask=g%mask2dT)
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))
1254 i_time_int = 0.0 ;
if (dt_int > 0.0) i_time_int = 1.0 / dt_int
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
1261 call post_data(ids%id_fraz, work_2d, diag, mask=g%mask2dT)
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
1269 call post_data(ids%id_salt_deficit, work_2d, diag, mask=g%mask2dT)
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)
1277 call post_data(ids%id_Heat_PmE, work_2d, diag, mask=g%mask2dT)
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)
1285 call post_data(ids%id_intern_heat, work_2d, diag, mask=g%mask2dT)
1288 if (tv%T_is_conT)
then
1290 if (ids%id_sstcon > 0)
call post_data(ids%id_sstcon, sfc_state%SST, diag, mask=g%mask2dT)
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))
1296 if (ids%id_sst > 0)
call post_data(ids%id_sst, work_2d, diag, mask=g%mask2dT)
1299 if (ids%id_sst > 0)
call post_data(ids%id_sst, sfc_state%SST, diag, mask=g%mask2dT)
1302 if (tv%S_is_absS)
then
1304 if (ids%id_sssabs > 0)
call post_data(ids%id_sssabs, sfc_state%SSS, diag, mask=g%mask2dT)
1307 do j=js,je ;
do i=is,ie
1308 work_2d(i,j) = gsw_sp_from_sr(sfc_state%SSS(i,j))
1310 if (ids%id_sss > 0)
call post_data(ids%id_sss, work_2d, diag, mask=g%mask2dT)
1313 if (ids%id_sss > 0)
call post_data(ids%id_sss, sfc_state%SSS, diag, mask=g%mask2dT)
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)
1320 call post_data(ids%id_sst_sq, work_2d, diag, mask=g%mask2dT)
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)
1326 call post_data(ids%id_sss_sq, work_2d, diag, mask=g%mask2dT)
1329 call coupler_type_send_data(sfc_state%tr_fields, get_diag_time_end(diag))
1337 diag, dt_trans, Reg)
1341 real,
dimension(SZIB_(G),SZJ_(G),SZK_(G)),
intent(in) :: uhtr
1343 real,
dimension(SZI_(G),SZJB_(G),SZK_(G)),
intent(in) :: vhtr
1345 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)), &
1350 real,
intent(in) :: dt_trans
1354 real,
dimension(SZIB_(G), SZJ_(G)) :: umo2d
1355 real,
dimension(SZI_(G), SZJB_(G)) :: vmo2d
1356 real,
dimension(SZIB_(G), SZJ_(G), SZK_(G)) :: umo
1357 real,
dimension(SZI_(G), SZJB_(G), SZK_(G)) :: vmo
1358 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)) :: h_tend
1361 real :: h_to_kg_m2_dt
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
1367 h_to_kg_m2_dt = gv%H_to_kg_m2 * us%L_to_m**2 * us%s_to_T * idt
1372 if (ids%id_umo_2d > 0)
then
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)
1379 if (ids%id_umo > 0)
then
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)
1386 if (ids%id_vmo_2d > 0)
then
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)
1393 if (ids%id_vmo > 0)
then
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)
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)
1406 if (ids%id_dynamics_h_tendency > 0)
then
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)
1422 subroutine mom_diagnostics_init(MIS, ADp, CDp, Time, G, GV, US, param_file, diag, CS, tv)
1430 type(time_type),
intent(in) :: time
1436 type(
diag_ctrl),
target,
intent(inout) :: diag
1443 real :: omega, f2_min, convert_h
1445 # include "version_variable.h"
1446 character(len=40) :: mdl =
"MOM_diagnostics"
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
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
1457 if (
associated(cs))
then
1458 call mom_error(warning,
"MOM_diagnostics_init called with an associated "// &
1459 "control structure.")
1465 use_temperature =
associated(tv%T)
1466 call get_param(param_file, mdl,
"ADIABATIC", adiabatic, default=.false., &
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.)
1480 if (gv%Boussinesq)
then
1481 thickness_units =
"m" ; flux_units =
"m3 s-1" ; convert_h = gv%H_to_m
1483 thickness_units =
"kg m-2" ; flux_units =
"kg s-1" ; convert_h = gv%H_to_kg_m2
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.)
1490 cs%id_masso = register_scalar_field(
'ocean_model',
'masso', time, &
1491 diag,
'Mass of liquid ocean',
'kg', standard_name=
'sea_water_mass')
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)
1502 cs%id_volcello = diag_get_volume_cell_measure_dm_id(diag)
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')
1509 if (tv%S_is_absS)
then
1510 cs%id_Sprac = register_diag_field(
'ocean_model',
'salt', diag%axesTL, &
1511 time,
'Salinity',
'psu')
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')
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')
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')
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')
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)
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)
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)
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)
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)
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')
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)
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)
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)
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)
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)
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)
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)
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)
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)
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)
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)
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)
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)
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)
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)
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)
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')
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
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)
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')
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')
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')
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')
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', &
1732 type(time_type),
intent(in) :: time
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')
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')
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')
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')
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)
1807 type(time_type),
intent(in) :: time
1815 character(len=48) :: thickness_units
1818 if (gv%Boussinesq)
then
1819 h_convert = gv%H_to_m
1821 h_convert = gv%H_to_kg_m2
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)
1858 type(
diag_ctrl),
target,
intent(inout) :: diag
1862 logical :: use_temperature
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.)
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.)
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.)
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.)
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.)
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.)
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.)
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.)
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')
1902 call post_data(id, g%areaT, diag, .true.)
1903 call diag_register_area_ids(diag, id_area_t=id)
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.)
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.)
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.)
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)
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.)
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.)
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.)
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.)
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.)
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.)
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.)
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.)
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.)
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.)
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.)
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.)
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.)
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.)
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.)
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', &
2005 if (id > 0)
call post_data(id, g%mask2dT, diag, .true.)
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.)
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.)
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
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)
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)
2056 if (.not.
associated(cs%dv_dt))
then
2057 call safe_alloc_ptr(cs%dv_dt,isd,ied,jsdb,jedb,nz)
2060 if (.not.
associated(cs%dh_dt))
then
2061 call safe_alloc_ptr(cs%dh_dt,isd,ied,jsd,jed,nz)
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)
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)
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)
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)
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)
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)
2122 do m=1,cs%num_time_deriv ;
deallocate(cs%prev_val(m)%p) ;
enddo