12 use mom_cpu_clock,
only : cpu_clock_id, cpu_clock_begin, cpu_clock_end
14 use mom_cpu_clock,
only : clock_module_driver, clock_module, clock_routine
20 use mom_domains,
only : to_south, to_west, to_all, cgrid_ne, scalar_pair
21 use mom_domains,
only : to_north, to_east, omit_corners
34 use mom_time_manager,
only :
operator(-),
operator(>),
operator(*),
operator(/)
37 use mom_barotropic,
only : barotropic_init, btstep, btcalc, bt_mass_source
65 implicit none ;
private
67 #include <MOM_memory.h>
71 real allocable_,
dimension(NIMEMB_PTR_,NJMEM_,NKMEM_) :: &
72 cau, & !< CAu = f*v - u.grad(u) [m s-2]
73 pfu, & !< PFu = -dM/dx [m s-2]
76 real allocable_,
dimension(NIMEM_,NJMEMB_PTR_,NKMEM_) :: &
77 cav, & !< CAv = -f*u - u.grad(v) [m s-2]
78 pfv, & !< PFv = -dM/dy [m s-2]
81 real allocable_,
dimension(NIMEMB_PTR_,NJMEM_,NKMEM_) :: visc_rem_u
87 real allocable_,
dimension(NIMEMB_PTR_,NJMEM_,NKMEM_) :: u_accel_bt
91 real allocable_,
dimension(NIMEM_,NJMEMB_PTR_,NKMEM_) :: visc_rem_v
97 real allocable_,
dimension(NIMEM_,NJMEMB_PTR_,NKMEM_) :: v_accel_bt
103 real allocable_,
dimension(NIMEM_,NJMEM_) :: eta
106 real allocable_,
dimension(NIMEMB_PTR_,NJMEM_,NKMEM_) :: u_av
109 real allocable_,
dimension(NIMEM_,NJMEMB_PTR_,NKMEM_) :: v_av
112 real allocable_,
dimension(NIMEM_,NJMEM_,NKMEM_) :: h_av
114 real allocable_,
dimension(NIMEM_,NJMEM_) :: eta_pf
116 real allocable_,
dimension(NIMEMB_PTR_,NJMEM_) :: uhbt
119 real allocable_,
dimension(NIMEM_,NJMEMB_PTR_) :: vhbt
122 real allocable_,
dimension(NIMEM_,NJMEM_,NKMEM_) :: pbce
126 real,
pointer,
dimension(:,:) :: taux_bot => null()
127 real,
pointer,
dimension(:,:) :: tauy_bot => null()
134 logical :: bt_use_layer_fluxes
137 logical :: split_bottom_stress
152 logical :: module_is_initialized = .false.
155 integer :: id_uh = -1, id_vh = -1
156 integer :: id_umo = -1, id_vmo = -1
157 integer :: id_umo_2d = -1, id_vmo_2d = -1
158 integer :: id_pfu = -1, id_pfv = -1
159 integer :: id_cau = -1, id_cav = -1
162 integer :: id_uav = -1, id_vav = -1
163 integer :: id_u_bt_accel = -1, id_v_bt_accel = -1
197 type(
ale_cs),
pointer :: ale_csp => null()
206 type(group_pass_type) :: pass_eta
207 type(group_pass_type) :: pass_visc_rem
208 type(group_pass_type) :: pass_uvp
209 type(group_pass_type) :: pass_hp_uv
210 type(group_pass_type) :: pass_uv
211 type(group_pass_type) :: pass_h
212 type(group_pass_type) :: pass_av_uvh
217 public step_mom_dyn_split_rk2
218 public register_restarts_dyn_split_rk2
219 public initialize_dyn_split_rk2
220 public end_dyn_split_rk2
223 integer :: id_clock_cor, id_clock_pres, id_clock_vertvisc
224 integer :: id_clock_horvisc, id_clock_mom_update
225 integer :: id_clock_continuity, id_clock_thick_diff
226 integer :: id_clock_btstep, id_clock_btcalc, id_clock_btforce
227 integer :: id_clock_pass, id_clock_pass_init
233 subroutine step_mom_dyn_split_rk2(u, v, h, tv, visc, &
234 Time_local, dt, forces, p_surf_begin, p_surf_end, &
235 uh, vh, uhtr, vhtr, eta_av, &
236 G, GV, US, CS, calc_dtbt, VarMix, MEKE, thickness_diffuse_CSp, Waves)
240 real,
dimension(SZIB_(G),SZJ_(G),SZK_(G)), &
241 target,
intent(inout) :: u
242 real,
dimension(SZI_(G),SZJB_(G),SZK_(G)), &
243 target,
intent(inout) :: v
244 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)), &
248 type(time_type),
intent(in) :: time_local
249 real,
intent(in) :: dt
251 real,
dimension(:,:),
pointer :: p_surf_begin
253 real,
dimension(:,:),
pointer :: p_surf_end
255 real,
dimension(SZIB_(G),SZJ_(G),SZK_(G)), &
256 target,
intent(inout) :: uh
258 real,
dimension(SZI_(G),SZJB_(G),SZK_(G)), &
259 target,
intent(inout) :: vh
261 real,
dimension(SZIB_(G),SZJ_(G),SZK_(G)), &
262 intent(inout) :: uhtr
264 real,
dimension(SZI_(G),SZJB_(G),SZK_(G)), &
265 intent(inout) :: vhtr
267 real,
dimension(SZI_(G),SZJ_(G)),
intent(out) :: eta_av
270 logical,
intent(in) :: calc_dtbt
281 real,
dimension(SZIB_(G),SZJ_(G),SZK_(G)) :: up
282 real,
dimension(SZI_(G),SZJB_(G),SZK_(G)) :: vp
283 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)) :: hp
285 real,
dimension(SZIB_(G),SZJ_(G),SZK_(G)) :: u_bc_accel
286 real,
dimension(SZI_(G),SZJB_(G),SZK_(G)) :: v_bc_accel
290 real,
dimension(SZIB_(G),SZJ_(G),SZK_(G)),
target :: uh_in
291 real,
dimension(SZI_(G),SZJB_(G),SZK_(G)),
target :: vh_in
295 real,
dimension(SZIB_(G),SZJ_(G)) :: uhbt_out
296 real,
dimension(SZI_(G),SZJB_(G)) :: vhbt_out
300 real,
dimension(SZI_(G),SZJ_(G)) :: eta_pred
304 real,
dimension(SZIB_(G),SZJ_(G),SZK_(G)),
target :: u_adj
305 real,
dimension(SZI_(G),SZJB_(G),SZK_(G)),
target :: v_adj
310 real,
dimension(SZIB_(G),SZJ_(G),SZK_(G)) :: u_old_rad_obc
311 real,
dimension(SZI_(G),SZJB_(G),SZK_(G)) :: v_old_rad_obc
316 real,
pointer,
dimension(:,:) :: &
317 p_surf => null(), eta_pf_start => null(), &
318 taux_bot => null(), tauy_bot => null(), &
321 real,
pointer,
dimension(:,:,:) :: &
322 uh_ptr => null(), u_ptr => null(), vh_ptr => null(), v_ptr => null(), &
323 u_init => null(), v_init => null(), &
328 logical :: dyn_p_surf
329 logical :: bt_cont_bt_thick
333 logical :: showcalltree, sym
335 integer :: i, j, k, is, ie, js, je, isq, ieq, jsq, jeq, nz
336 integer :: cont_stencil
337 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
338 isq = g%IscB ; ieq = g%IecB ; jsq = g%JscB ; jeq = g%JecB
339 u_av => cs%u_av ; v_av => cs%v_av ; h_av => cs%h_av ; eta => cs%eta
342 sym=.false.;
if (g%Domain%symmetric) sym=.true.
344 showcalltree = calltree_showquery()
345 if (showcalltree)
call calltree_enter(
"step_MOM_dyn_split_RK2(), MOM_dynamics_split_RK2.F90")
349 do j=g%jsd,g%jed ;
do i=g%isdB,g%iedB ; up(i,j,k) = 0.0 ;
enddo ;
enddo
350 do j=g%jsdB,g%jedB ;
do i=g%isd,g%ied ; vp(i,j,k) = 0.0 ;
enddo ;
enddo
351 do j=g%jsd,g%jed ;
do i=g%isd,g%ied ; hp(i,j,k) = h(i,j,k) ;
enddo ;
enddo
355 call updatecfltruncationvalue(time_local, cs%vertvisc_CSp)
358 call mom_state_chksum(
"Start predictor ", u, v, h, uh, vh, g, gv, symmetric=sym)
363 dyn_p_surf =
associated(p_surf_begin) .and.
associated(p_surf_end)
366 call safe_alloc_ptr(eta_pf_start,g%isd,g%ied,g%jsd,g%jed)
367 eta_pf_start(:,:) = 0.0
369 p_surf => forces%p_surf
372 if (
associated(cs%OBC))
then
373 if (cs%debug_OBC)
call open_boundary_test_extern_h(g, cs%OBC, h)
375 do k=1,nz ;
do j=g%jsd,g%jed ;
do i=g%IsdB,g%IedB
376 u_old_rad_obc(i,j,k) = u_av(i,j,k)
377 enddo ;
enddo ;
enddo
378 do k=1,nz ;
do j=g%JsdB,g%JedB ;
do i=g%isd,g%ied
379 v_old_rad_obc(i,j,k) = v_av(i,j,k)
380 enddo ;
enddo ;
enddo
383 bt_cont_bt_thick = .false.
384 if (
associated(cs%BT_cont)) bt_cont_bt_thick = &
385 (
allocated(cs%BT_cont%h_u) .and.
allocated(cs%BT_cont%h_v))
387 if (cs%split_bottom_stress)
then
388 taux_bot => cs%taux_bot ; tauy_bot => cs%tauy_bot
393 cont_stencil = continuity_stencil(cs%continuity_CSp)
396 call cpu_clock_begin(id_clock_pass)
398 call create_group_pass(cs%pass_visc_rem, cs%visc_rem_u, cs%visc_rem_v, g%Domain, &
399 to_all+scalar_pair, cgrid_ne, halo=max(1,cont_stencil))
409 call cpu_clock_end(id_clock_pass)
415 if (cs%begw == 0.0)
call enable_averaging(dt, time_local, cs%diag)
416 call cpu_clock_begin(id_clock_pres)
417 call pressureforce(h, tv, cs%PFu, cs%PFv, g, gv, us, cs%PressureForce_CSp, &
418 cs%ALE_CSp, p_surf, cs%pbce, cs%eta_PF)
420 pa_to_eta = 1.0 / gv%H_to_Pa
422 do j=jsq,jeq+1 ;
do i=isq,ieq+1
423 eta_pf_start(i,j) = cs%eta_PF(i,j) - pa_to_eta * &
424 (p_surf_begin(i,j) - p_surf_end(i,j))
427 call cpu_clock_end(id_clock_pres)
428 call disable_averaging(cs%diag)
429 if (showcalltree)
call calltree_waypoint(
"done with PressureForce (step_MOM_dyn_split_RK2)")
431 if (
associated(cs%OBC)) then;
if (cs%OBC%update_OBC)
then
432 call update_obc_data(cs%OBC, g, gv, us, tv, h, cs%update_OBC_CSp, time_local)
434 if (
associated(cs%OBC) .and. cs%debug_OBC) &
435 call open_boundary_zero_normal_flow(cs%OBC, g, cs%PFu, cs%PFv)
437 if (g%nonblocking_updates) &
438 call start_group_pass(cs%pass_eta, g%Domain, clock=id_clock_pass)
441 call cpu_clock_begin(id_clock_cor)
442 call coradcalc(u_av, v_av, h_av, uh, vh, cs%CAu, cs%CAv, cs%OBC, cs%ADp, &
443 g, gv, us, cs%CoriolisAdv_CSp)
444 call cpu_clock_end(id_clock_cor)
445 if (showcalltree)
call calltree_waypoint(
"done with CorAdCalc (step_MOM_dyn_split_RK2)")
448 call cpu_clock_begin(id_clock_btforce)
451 do j=js,je ;
do i=isq,ieq
452 u_bc_accel(i,j,k) = (cs%Cau(i,j,k) + cs%PFu(i,j,k)) + us%s_to_T*cs%diffu(i,j,k)
454 do j=jsq,jeq ;
do i=is,ie
455 v_bc_accel(i,j,k) = (cs%Cav(i,j,k) + cs%PFv(i,j,k)) + us%s_to_T*cs%diffv(i,j,k)
458 if (
associated(cs%OBC))
then
459 call open_boundary_zero_normal_flow(cs%OBC, g, u_bc_accel, v_bc_accel)
461 call cpu_clock_end(id_clock_btforce)
464 call mom_accel_chksum(
"pre-btstep accel", cs%CAu, cs%CAv, cs%PFu, cs%PFv, &
465 cs%diffu, cs%diffv, g, gv, us, cs%pbce, u_bc_accel, v_bc_accel, &
470 call check_redundant(
"pre-btstep u_bc_accel ", u_bc_accel, v_bc_accel, g)
473 call cpu_clock_begin(id_clock_vertvisc)
476 do j=js,je ;
do i=isq,ieq
477 up(i,j,k) = g%mask2dCu(i,j) * (u(i,j,k) + dt * u_bc_accel(i,j,k))
479 do j=jsq,jeq ;
do i=is,ie
480 vp(i,j,k) = g%mask2dCv(i,j) * (v(i,j,k) + dt * v_bc_accel(i,j,k))
484 call enable_averaging(dt, time_local, cs%diag)
485 call set_viscous_ml(u, v, h, tv, forces, visc, dt, g, gv, us, &
487 call disable_averaging(cs%diag)
490 call uvchksum(
"before vertvisc: up", up, vp, g%HI, haloshift=0, symmetric=sym)
492 call vertvisc_coef(up, vp, h, forces, visc, dt, g, gv, us, cs%vertvisc_CSp, cs%OBC)
493 call vertvisc_remnant(visc, cs%visc_rem_u, cs%visc_rem_v, dt, g, gv, us, cs%vertvisc_CSp)
494 call cpu_clock_end(id_clock_vertvisc)
495 if (showcalltree)
call calltree_waypoint(
"done with vertvisc_coef (step_MOM_dyn_split_RK2)")
498 call cpu_clock_begin(id_clock_pass)
499 if (g%nonblocking_updates)
then
500 call complete_group_pass(cs%pass_eta, g%Domain)
501 call start_group_pass(cs%pass_visc_rem, g%Domain)
503 call do_group_pass(cs%pass_eta, g%Domain)
504 call do_group_pass(cs%pass_visc_rem, g%Domain)
506 call cpu_clock_end(id_clock_pass)
508 call cpu_clock_begin(id_clock_btcalc)
510 if (.not.bt_cont_bt_thick) &
511 call btcalc(h, g, gv, cs%barotropic_CSp, obc=cs%OBC)
512 call bt_mass_source(h, eta, .true., g, gv, cs%barotropic_CSp)
513 call cpu_clock_end(id_clock_btcalc)
515 if (g%nonblocking_updates) &
516 call complete_group_pass(cs%pass_visc_rem, g%Domain, clock=id_clock_pass)
519 if (
associated(cs%BT_cont) .or. cs%BT_use_layer_fluxes)
then
520 call cpu_clock_begin(id_clock_continuity)
521 call continuity(u, v, h, hp, uh_in, vh_in, dt, g, gv, us, &
522 cs%continuity_CSp, obc=cs%OBC, visc_rem_u=cs%visc_rem_u, &
523 visc_rem_v=cs%visc_rem_v, bt_cont=cs%BT_cont)
524 call cpu_clock_end(id_clock_continuity)
525 if (bt_cont_bt_thick)
then
526 call btcalc(h, g, gv, cs%barotropic_CSp, cs%BT_cont%h_u, cs%BT_cont%h_v, &
529 if (showcalltree)
call calltree_waypoint(
"done with continuity[BT_cont] (step_MOM_dyn_split_RK2)")
532 if (cs%BT_use_layer_fluxes)
then
533 uh_ptr => uh_in; vh_ptr => vh_in; u_ptr => u; v_ptr => v
536 u_init => u ; v_init => v
537 call cpu_clock_begin(id_clock_btstep)
538 if (calc_dtbt)
call set_dtbt(g, gv, us, cs%barotropic_CSp, eta, cs%pbce)
539 if (showcalltree)
call calltree_enter(
"btstep(), MOM_barotropic.F90")
541 call btstep(u, v, eta, dt, u_bc_accel, v_bc_accel, forces, cs%pbce, cs%eta_PF, &
542 u_av, v_av, cs%u_accel_bt, cs%v_accel_bt, eta_pred, cs%uhbt, cs%vhbt, &
543 g, gv, us, cs%barotropic_CSp, cs%visc_rem_u, cs%visc_rem_v, &
544 obc=cs%OBC, bt_cont=cs%BT_cont, eta_pf_start=eta_pf_start, &
545 taux_bot=taux_bot, tauy_bot=tauy_bot, &
546 uh0=uh_ptr, vh0=vh_ptr, u_uh0=u_ptr, v_vh0=v_ptr)
547 if (showcalltree)
call calltree_leave(
"btstep()")
548 call cpu_clock_end(id_clock_btstep)
552 call cpu_clock_begin(id_clock_mom_update)
556 do j=jsq,jeq ;
do i=is,ie
557 vp(i,j,k) = g%mask2dCv(i,j) * (v_init(i,j,k) + dt_pred * &
558 (v_bc_accel(i,j,k) + cs%v_accel_bt(i,j,k)))
560 do j=js,je ;
do i=isq,ieq
561 up(i,j,k) = g%mask2dCu(i,j) * (u_init(i,j,k) + dt_pred * &
562 (u_bc_accel(i,j,k) + cs%u_accel_bt(i,j,k)))
565 call cpu_clock_end(id_clock_mom_update)
568 call uvchksum(
"Predictor 1 [uv]", up, vp, g%HI, haloshift=0, symmetric=sym)
569 call hchksum(h,
"Predictor 1 h", g%HI, haloshift=1, scale=gv%H_to_m)
570 call uvchksum(
"Predictor 1 [uv]h", uh, vh, g%HI,haloshift=2, &
571 symmetric=sym, scale=gv%H_to_m)
573 call mom_accel_chksum(
"Predictor accel", cs%CAu, cs%CAv, cs%PFu, cs%PFv, &
574 cs%diffu, cs%diffv, g, gv, us, cs%pbce, cs%u_accel_bt, cs%v_accel_bt, symmetric=sym)
575 call mom_state_chksum(
"Predictor 1 init", u_init, v_init, h, uh, vh, g, gv, haloshift=2, &
583 call cpu_clock_begin(id_clock_vertvisc)
585 call uvchksum(
"0 before vertvisc: [uv]p", up, vp, g%HI,haloshift=0, symmetric=sym)
587 call vertvisc_coef(up, vp, h, forces, visc, dt_pred, g, gv, us, cs%vertvisc_CSp, &
589 call vertvisc(up, vp, h, forces, visc, dt_pred, cs%OBC, cs%ADp, cs%CDp, g, &
590 gv, us, cs%vertvisc_CSp, cs%taux_bot, cs%tauy_bot, waves=waves)
591 if (showcalltree)
call calltree_waypoint(
"done with vertvisc (step_MOM_dyn_split_RK2)")
592 if (g%nonblocking_updates)
then
593 call cpu_clock_end(id_clock_vertvisc)
594 call start_group_pass(cs%pass_uvp, g%Domain, clock=id_clock_pass)
595 call cpu_clock_begin(id_clock_vertvisc)
597 call vertvisc_remnant(visc, cs%visc_rem_u, cs%visc_rem_v, dt_pred, g, gv, us, cs%vertvisc_CSp)
598 call cpu_clock_end(id_clock_vertvisc)
600 call do_group_pass(cs%pass_visc_rem, g%Domain, clock=id_clock_pass)
601 if (g%nonblocking_updates)
then
602 call complete_group_pass(cs%pass_uvp, g%Domain, clock=id_clock_pass)
604 call do_group_pass(cs%pass_uvp, g%Domain, clock=id_clock_pass)
609 call cpu_clock_begin(id_clock_continuity)
610 call continuity(up, vp, h, hp, uh, vh, dt, g, gv, us, cs%continuity_CSp, &
611 cs%uhbt, cs%vhbt, cs%OBC, cs%visc_rem_u, cs%visc_rem_v, &
612 u_av, v_av, bt_cont=cs%BT_cont)
613 call cpu_clock_end(id_clock_continuity)
614 if (showcalltree)
call calltree_waypoint(
"done with continuity (step_MOM_dyn_split_RK2)")
616 call do_group_pass(cs%pass_hp_uv, g%Domain, clock=id_clock_pass)
618 if (
associated(cs%OBC))
then
621 call uvchksum(
"Pre OBC avg [uv]", u_av, v_av, g%HI, haloshift=1, symmetric=sym)
623 call radiation_open_bdry_conds(cs%OBC, u_av, u_old_rad_obc, v_av, v_old_rad_obc, g, dt_pred)
626 call uvchksum(
"Post OBC avg [uv]", u_av, v_av, g%HI, haloshift=1, symmetric=sym)
632 if (g%nonblocking_updates)
then
633 call start_group_pass(cs%pass_av_uvh, g%Domain, clock=id_clock_pass)
638 do k=1,nz ;
do j=js-2,je+2 ;
do i=is-2,ie+2
639 h_av(i,j,k) = 0.5*(h(i,j,k) + hp(i,j,k))
640 enddo ;
enddo ;
enddo
643 call enable_averaging(dt, time_local, cs%diag)
649 call cpu_clock_begin(id_clock_btcalc)
650 call bt_mass_source(hp, eta_pred, .false., g, gv, cs%barotropic_CSp)
651 call cpu_clock_end(id_clock_btcalc)
653 if (cs%begw /= 0.0)
then
658 do k=1,nz ;
do j=js-1,je+1 ;
do i=is-1,ie+1
659 hp(i,j,k) = (1.0-cs%begw)*h(i,j,k) + cs%begw*hp(i,j,k)
660 enddo ;
enddo ;
enddo
664 call cpu_clock_begin(id_clock_pres)
665 call pressureforce(hp, tv, cs%PFu, cs%PFv, g, gv, us, cs%PressureForce_CSp, &
666 cs%ALE_CSp, p_surf, cs%pbce, cs%eta_PF)
667 call cpu_clock_end(id_clock_pres)
668 if (showcalltree)
call calltree_waypoint(
"done with PressureForce[hp=(1-b).h+b.h] (step_MOM_dyn_split_RK2)")
671 if (g%nonblocking_updates) &
672 call complete_group_pass(cs%pass_av_uvh, g%Domain, clock=id_clock_pass)
674 if (bt_cont_bt_thick)
then
675 call btcalc(h, g, gv, cs%barotropic_CSp, cs%BT_cont%h_u, cs%BT_cont%h_v, &
677 if (showcalltree)
call calltree_waypoint(
"done with btcalc[BT_cont_BT_thick] (step_MOM_dyn_split_RK2)")
681 call mom_state_chksum(
"Predictor ", up, vp, hp, uh, vh, g, gv, symmetric=sym)
682 call uvchksum(
"Predictor avg [uv]", u_av, v_av, g%HI, haloshift=1, symmetric=sym)
683 call hchksum(h_av,
"Predictor avg h", g%HI, haloshift=0, scale=gv%H_to_m)
690 call cpu_clock_begin(id_clock_horvisc)
691 call horizontal_viscosity(u_av, v_av, h_av, cs%diffu, cs%diffv, &
692 meke, varmix, g, gv, us, cs%hor_visc_CSp, &
693 obc=cs%OBC, bt=cs%barotropic_CSp)
694 call cpu_clock_end(id_clock_horvisc)
695 if (showcalltree)
call calltree_waypoint(
"done with horizontal_viscosity (step_MOM_dyn_split_RK2)")
698 call cpu_clock_begin(id_clock_cor)
699 call coradcalc(u_av, v_av, h_av, uh, vh, cs%CAu, cs%CAv, cs%OBC, cs%ADp, &
700 g, gv, us, cs%CoriolisAdv_CSp)
701 call cpu_clock_end(id_clock_cor)
702 if (showcalltree)
call calltree_waypoint(
"done with CorAdCalc (step_MOM_dyn_split_RK2)")
707 call cpu_clock_begin(id_clock_btforce)
710 do j=js,je ;
do i=isq,ieq
711 u_bc_accel(i,j,k) = (cs%Cau(i,j,k) + cs%PFu(i,j,k)) + us%s_to_T*cs%diffu(i,j,k)
713 do j=jsq,jeq ;
do i=is,ie
714 v_bc_accel(i,j,k) = (cs%Cav(i,j,k) + cs%PFv(i,j,k)) + us%s_to_T*cs%diffv(i,j,k)
717 if (
associated(cs%OBC))
then
718 call open_boundary_zero_normal_flow(cs%OBC, g, u_bc_accel, v_bc_accel)
720 call cpu_clock_end(id_clock_btforce)
723 call mom_accel_chksum(
"corr pre-btstep accel", cs%CAu, cs%CAv, cs%PFu, cs%PFv, &
724 cs%diffu, cs%diffv, g, gv, us, cs%pbce, u_bc_accel, v_bc_accel, &
728 call check_redundant(
"corr pre-btstep CS%diff ", cs%diffu, cs%diffv, g)
729 call check_redundant(
"corr pre-btstep u_bc_accel ", u_bc_accel, v_bc_accel, g)
734 call cpu_clock_begin(id_clock_btstep)
735 if (cs%BT_use_layer_fluxes)
then
736 uh_ptr => uh ; vh_ptr => vh ; u_ptr => u_av ; v_ptr => v_av
739 if (showcalltree)
call calltree_enter(
"btstep(), MOM_barotropic.F90")
741 call btstep(u, v, eta, dt, u_bc_accel, v_bc_accel, forces, cs%pbce, &
742 cs%eta_PF, u_av, v_av, cs%u_accel_bt, cs%v_accel_bt, &
743 eta_pred, cs%uhbt, cs%vhbt, g, gv, us, cs%barotropic_CSp, &
744 cs%visc_rem_u, cs%visc_rem_v, etaav=eta_av, obc=cs%OBC, &
745 bt_cont = cs%BT_cont, eta_pf_start=eta_pf_start, &
746 taux_bot=taux_bot, tauy_bot=tauy_bot, &
747 uh0=uh_ptr, vh0=vh_ptr, u_uh0=u_ptr, v_vh0=v_ptr)
748 do j=js,je ;
do i=is,ie ; eta(i,j) = eta_pred(i,j) ;
enddo ;
enddo
749 call cpu_clock_end(id_clock_btstep)
750 if (showcalltree)
call calltree_leave(
"btstep()")
757 call cpu_clock_begin(id_clock_mom_update)
760 do j=js,je ;
do i=isq,ieq
761 u(i,j,k) = g%mask2dCu(i,j) * (u_init(i,j,k) + dt * &
762 (u_bc_accel(i,j,k) + cs%u_accel_bt(i,j,k)))
764 do j=jsq,jeq ;
do i=is,ie
765 v(i,j,k) = g%mask2dCv(i,j) * (v_init(i,j,k) + dt * &
766 (v_bc_accel(i,j,k) + cs%v_accel_bt(i,j,k)))
769 call cpu_clock_end(id_clock_mom_update)
772 call uvchksum(
"Corrector 1 [uv]", u, v, g%HI,haloshift=0, symmetric=sym)
773 call hchksum(h,
"Corrector 1 h", g%HI, haloshift=2, scale=gv%H_to_m)
774 call uvchksum(
"Corrector 1 [uv]h", uh, vh, g%HI, haloshift=2, &
775 symmetric=sym, scale=gv%H_to_m)
777 call mom_accel_chksum(
"Corrector accel", cs%CAu, cs%CAv, cs%PFu, cs%PFv, &
778 cs%diffu, cs%diffv, g, gv, us, cs%pbce, cs%u_accel_bt, cs%v_accel_bt, &
784 call cpu_clock_begin(id_clock_vertvisc)
785 call vertvisc_coef(u, v, h, forces, visc, dt, g, gv, us, cs%vertvisc_CSp, cs%OBC)
786 call vertvisc(u, v, h, forces, visc, dt, cs%OBC, cs%ADp, cs%CDp, g, gv, us, &
787 cs%vertvisc_CSp, cs%taux_bot, cs%tauy_bot,waves=waves)
788 if (g%nonblocking_updates)
then
789 call cpu_clock_end(id_clock_vertvisc)
790 call start_group_pass(cs%pass_uv, g%Domain, clock=id_clock_pass)
791 call cpu_clock_begin(id_clock_vertvisc)
793 call vertvisc_remnant(visc, cs%visc_rem_u, cs%visc_rem_v, dt, g, gv, us, cs%vertvisc_CSp)
794 call cpu_clock_end(id_clock_vertvisc)
795 if (showcalltree)
call calltree_waypoint(
"done with vertvisc (step_MOM_dyn_split_RK2)")
799 do k=1,nz ;
do j=js-2,je+2 ;
do i=is-2,ie+2
800 h_av(i,j,k) = h(i,j,k)
801 enddo ;
enddo ;
enddo
803 call do_group_pass(cs%pass_visc_rem, g%Domain, clock=id_clock_pass)
804 if (g%nonblocking_updates)
then
805 call complete_group_pass(cs%pass_uv, g%Domain, clock=id_clock_pass)
807 call do_group_pass(cs%pass_uv, g%Domain, clock=id_clock_pass)
813 call cpu_clock_begin(id_clock_continuity)
814 call continuity(u, v, h, h, uh, vh, dt, g, gv, us, &
815 cs%continuity_CSp, cs%uhbt, cs%vhbt, cs%OBC, &
816 cs%visc_rem_u, cs%visc_rem_v, u_av, v_av)
817 call cpu_clock_end(id_clock_continuity)
818 call do_group_pass(cs%pass_h, g%Domain, clock=id_clock_pass)
821 call diag_update_remap_grids(cs%diag)
822 if (showcalltree)
call calltree_waypoint(
"done with continuity (step_MOM_dyn_split_RK2)")
824 if (g%nonblocking_updates)
then
825 call start_group_pass(cs%pass_av_uvh, g%Domain, clock=id_clock_pass)
827 call do_group_pass(cs%pass_av_uvh, g%domain, clock=id_clock_pass)
830 if (
associated(cs%OBC))
then
831 call radiation_open_bdry_conds(cs%OBC, u, u_old_rad_obc, v, v_old_rad_obc, g, dt)
836 do k=1,nz ;
do j=js-2,je+2 ;
do i=is-2,ie+2
837 h_av(i,j,k) = 0.5*(h_av(i,j,k) + h(i,j,k))
838 enddo ;
enddo ;
enddo
840 if (g%nonblocking_updates) &
841 call complete_group_pass(cs%pass_av_uvh, g%Domain, clock=id_clock_pass)
845 do j=js-2,je+2 ;
do i=isq-2,ieq+2
846 uhtr(i,j,k) = uhtr(i,j,k) + uh(i,j,k)*dt
848 do j=jsq-2,jeq+2 ;
do i=is-2,ie+2
849 vhtr(i,j,k) = vhtr(i,j,k) + vh(i,j,k)*dt
858 if (cs%id_PFu > 0)
call post_data(cs%id_PFu, cs%PFu, cs%diag)
859 if (cs%id_PFv > 0)
call post_data(cs%id_PFv, cs%PFv, cs%diag)
860 if (cs%id_CAu > 0)
call post_data(cs%id_CAu, cs%CAu, cs%diag)
861 if (cs%id_CAv > 0)
call post_data(cs%id_CAv, cs%CAv, cs%diag)
864 if (cs%id_uh > 0)
call post_data(cs%id_uh , uh, cs%diag)
865 if (cs%id_vh > 0)
call post_data(cs%id_vh , vh, cs%diag)
866 if (cs%id_uav > 0)
call post_data(cs%id_uav, u_av, cs%diag)
867 if (cs%id_vav > 0)
call post_data(cs%id_vav, v_av, cs%diag)
868 if (cs%id_u_BT_accel > 0)
call post_data(cs%id_u_BT_accel, cs%u_accel_bt, cs%diag)
869 if (cs%id_v_BT_accel > 0)
call post_data(cs%id_v_BT_accel, cs%v_accel_bt, cs%diag)
873 call uvchksum(
"Corrector avg [uv]", u_av, v_av, g%HI,haloshift=1, symmetric=sym)
874 call hchksum(h_av,
"Corrector avg h", g%HI, haloshift=1, scale=gv%H_to_m)
878 if (showcalltree)
call calltree_leave(
"step_MOM_dyn_split_RK2()")
880 end subroutine step_mom_dyn_split_rk2
885 subroutine register_restarts_dyn_split_rk2(HI, GV, param_file, CS, restart_CS, uh, vh)
891 real,
dimension(SZIB_(HI),SZJ_(HI),SZK_(GV)), &
892 target,
intent(inout) :: uh
893 real,
dimension(SZI_(HI),SZJB_(HI),SZK_(GV)), &
894 target,
intent(inout) :: vh
897 character(len=40) :: mdl =
"MOM_dynamics_split_RK2"
898 character(len=48) :: thickness_units, flux_units
900 integer :: isd, ied, jsd, jed, nz, isdb, iedb, jsdb, jedb
901 isd = hi%isd ; ied = hi%ied ; jsd = hi%jsd ; jed = hi%jed ; nz = gv%ke
902 isdb = hi%IsdB ; iedb = hi%IedB ; jsdb = hi%JsdB ; jedb = hi%JedB
905 if (
associated(cs))
then
906 call mom_error(warning,
"register_restarts_dyn_split_RK2 called with an associated "// &
907 "control structure.")
912 alloc_(cs%diffu(isdb:iedb,jsd:jed,nz)) ; cs%diffu(:,:,:) = 0.0
913 alloc_(cs%diffv(isd:ied,jsdb:jedb,nz)) ; cs%diffv(:,:,:) = 0.0
914 alloc_(cs%CAu(isdb:iedb,jsd:jed,nz)) ; cs%CAu(:,:,:) = 0.0
915 alloc_(cs%CAv(isd:ied,jsdb:jedb,nz)) ; cs%CAv(:,:,:) = 0.0
916 alloc_(cs%PFu(isdb:iedb,jsd:jed,nz)) ; cs%PFu(:,:,:) = 0.0
917 alloc_(cs%PFv(isd:ied,jsdb:jedb,nz)) ; cs%PFv(:,:,:) = 0.0
919 alloc_(cs%eta(isd:ied,jsd:jed)) ; cs%eta(:,:) = 0.0
920 alloc_(cs%u_av(isdb:iedb,jsd:jed,nz)) ; cs%u_av(:,:,:) = 0.0
921 alloc_(cs%v_av(isd:ied,jsdb:jedb,nz)) ; cs%v_av(:,:,:) = 0.0
922 alloc_(cs%h_av(isd:ied,jsd:jed,nz)) ; cs%h_av(:,:,:) = gv%Angstrom_H
924 thickness_units = get_thickness_units(gv)
925 flux_units = get_flux_units(gv)
927 if (gv%Boussinesq)
then
928 vd = var_desc(
"sfc",thickness_units,
"Free surface Height",
'h',
'1')
930 vd = var_desc(
"p_bot",thickness_units,
"Bottom Pressure",
'h',
'1')
934 vd = var_desc(
"u2",
"m s-1",
"Auxiliary Zonal velocity",
'u',
'L')
937 vd = var_desc(
"v2",
"m s-1",
"Auxiliary Meridional velocity",
'v',
'L')
940 vd = var_desc(
"h2",thickness_units,
"Auxiliary Layer Thickness",
'h',
'L')
943 vd = var_desc(
"uh",flux_units,
"Zonal thickness flux",
'u',
'L')
946 vd = var_desc(
"vh",flux_units,
"Meridional thickness flux",
'v',
'L')
949 vd = var_desc(
"diffu",
"m s-2",
"Zonal horizontal viscous acceleration",
'u',
'L')
952 vd = var_desc(
"diffv",
"m s-2",
"Meridional horizontal viscous acceleration",
'v',
'L')
955 call register_barotropic_restarts(hi, gv, param_file, cs%barotropic_CSp, &
958 end subroutine register_restarts_dyn_split_rk2
962 subroutine initialize_dyn_split_rk2(u, v, h, uh, vh, eta, Time, G, GV, US, param_file, &
963 diag, CS, restart_CS, dt, Accel_diag, Cont_diag, MIS, &
964 VarMix, MEKE, thickness_diffuse_CSp, &
965 OBC, update_OBC_CSp, ALE_CSp, setVisc_CSp, &
966 visc, dirs, ntrunc, calc_dtbt)
970 real,
dimension(SZIB_(G),SZJ_(G),SZK_(G)), &
972 real,
dimension(SZI_(G),SZJB_(G),SZK_(G)), &
974 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)) ,
intent(inout) :: h
975 real,
dimension(SZIB_(G),SZJ_(G),SZK_(G)), &
976 target,
intent(inout) :: uh
977 real,
dimension(SZI_(G),SZJB_(G),SZK_(G)), &
978 target,
intent(inout) :: vh
979 real,
dimension(SZI_(G),SZJ_(G)),
intent(inout) :: eta
980 type(time_type),
target,
intent(in) :: time
982 type(
diag_ctrl),
target,
intent(inout) :: diag
985 real,
intent(in) :: dt
999 type(
ale_cs),
pointer :: ale_csp
1003 integer,
target,
intent(inout) :: ntrunc
1006 logical,
intent(out) :: calc_dtbt
1009 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)) :: h_tmp
1010 character(len=40) :: mdl =
"MOM_dynamics_split_RK2"
1011 character(len=48) :: thickness_units, flux_units, eta_rest_name
1017 type(group_pass_type) :: pass_av_h_uvh
1018 logical :: use_tides, debug_truncations
1020 integer :: i, j, k, is, ie, js, je, isd, ied, jsd, jed, nz
1021 integer :: isdb, iedb, jsdb, jedb
1022 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
1023 isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
1024 isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
1026 if (.not.
associated(cs))
call mom_error(fatal, &
1027 "initialize_dyn_split_RK2 called with an unassociated control structure.")
1028 if (cs%module_is_initialized)
then
1029 call mom_error(warning,
"initialize_dyn_split_RK2 called with a control "// &
1030 "structure that has already been initialized.")
1033 cs%module_is_initialized = .true.
1037 call get_param(param_file, mdl,
"TIDES", use_tides, &
1038 "If true, apply tidal momentum forcing.", default=.false.)
1039 call get_param(param_file, mdl,
"BE", cs%be, &
1040 "If SPLIT is true, BE determines the relative weighting "//&
1041 "of a 2nd-order Runga-Kutta baroclinic time stepping "//&
1042 "scheme (0.5) and a backward Euler scheme (1) that is "//&
1043 "used for the Coriolis and inertial terms. BE may be "//&
1044 "from 0.5 to 1, but instability may occur near 0.5. "//&
1045 "BE is also applicable if SPLIT is false and USE_RK2 "//&
1046 "is true.", units=
"nondim", default=0.6)
1047 call get_param(param_file, mdl,
"BEGW", cs%begw, &
1048 "If SPLIT is true, BEGW is a number from 0 to 1 that "//&
1049 "controls the extent to which the treatment of gravity "//&
1050 "waves is forward-backward (0) or simulated backward "//&
1051 "Euler (1). 0 is almost always used. "//&
1052 "If SPLIT is false and USE_RK2 is true, BEGW can be "//&
1053 "between 0 and 0.5 to damp gravity waves.", &
1054 units=
"nondim", default=0.0)
1056 call get_param(param_file, mdl,
"SPLIT_BOTTOM_STRESS", cs%split_bottom_stress, &
1057 "If true, provide the bottom stress calculated by the "//&
1058 "vertical viscosity to the barotropic solver.", default=.false.)
1059 call get_param(param_file, mdl,
"BT_USE_LAYER_FLUXES", cs%BT_use_layer_fluxes, &
1060 "If true, use the summed layered fluxes plus an "//&
1061 "adjustment due to the change in the barotropic velocity "//&
1062 "in the barotropic continuity equation.", default=.true.)
1063 call get_param(param_file, mdl,
"DEBUG", cs%debug, &
1064 "If true, write out verbose debugging data.", &
1065 default=.false., debuggingparam=.true.)
1066 call get_param(param_file, mdl,
"DEBUG_OBC", cs%debug_OBC, default=.false.)
1067 call get_param(param_file, mdl,
"DEBUG_TRUNCATIONS", debug_truncations, &
1070 allocate(cs%taux_bot(isdb:iedb,jsd:jed)) ; cs%taux_bot(:,:) = 0.0
1071 allocate(cs%tauy_bot(isd:ied,jsdb:jedb)) ; cs%tauy_bot(:,:) = 0.0
1073 alloc_(cs%uhbt(isdb:iedb,jsd:jed)) ; cs%uhbt(:,:) = 0.0
1074 alloc_(cs%vhbt(isd:ied,jsdb:jedb)) ; cs%vhbt(:,:) = 0.0
1075 alloc_(cs%visc_rem_u(isdb:iedb,jsd:jed,nz)) ; cs%visc_rem_u(:,:,:) = 0.0
1076 alloc_(cs%visc_rem_v(isd:ied,jsdb:jedb,nz)) ; cs%visc_rem_v(:,:,:) = 0.0
1077 alloc_(cs%eta_PF(isd:ied,jsd:jed)) ; cs%eta_PF(:,:) = 0.0
1078 alloc_(cs%pbce(isd:ied,jsd:jed,nz)) ; cs%pbce(:,:,:) = 0.0
1080 alloc_(cs%u_accel_bt(isdb:iedb,jsd:jed,nz)) ; cs%u_accel_bt(:,:,:) = 0.0
1081 alloc_(cs%v_accel_bt(isd:ied,jsdb:jedb,nz)) ; cs%v_accel_bt(:,:,:) = 0.0
1083 mis%diffu => cs%diffu
1084 mis%diffv => cs%diffv
1090 mis%u_accel_bt => cs%u_accel_bt
1091 mis%v_accel_bt => cs%v_accel_bt
1095 cs%ADp => accel_diag
1097 accel_diag%diffu => cs%diffu
1098 accel_diag%diffv => cs%diffv
1099 accel_diag%PFu => cs%PFu
1100 accel_diag%PFv => cs%PFv
1101 accel_diag%CAu => cs%CAu
1102 accel_diag%CAv => cs%CAv
1108 call continuity_init(time, g, gv, param_file, diag, cs%continuity_CSp)
1109 call coriolisadv_init(time, g, param_file, diag, cs%ADp, cs%CoriolisAdv_CSp)
1110 if (use_tides)
call tidal_forcing_init(time, g, param_file, cs%tides_CSp)
1111 call pressureforce_init(time, g, gv, us, param_file, diag, cs%PressureForce_CSp, &
1113 call hor_visc_init(time, g, us, param_file, diag, cs%hor_visc_CSp, meke)
1114 call vertvisc_init(mis, time, g, gv, us, param_file, diag, cs%ADp, dirs, &
1115 ntrunc, cs%vertvisc_CSp)
1116 if (.not.
associated(setvisc_csp))
call mom_error(fatal, &
1117 "initialize_dyn_split_RK2 called with setVisc_CSp unassociated.")
1118 cs%set_visc_CSp => setvisc_csp
1119 call updatecfltruncationvalue(time, cs%vertvisc_CSp, &
1120 activate=is_new_run(restart_cs) )
1122 if (
associated(ale_csp)) cs%ALE_CSp => ale_csp
1123 if (
associated(obc)) cs%OBC => obc
1124 if (
associated(update_obc_csp)) cs%update_OBC_CSp => update_obc_csp
1126 eta_rest_name =
"sfc" ;
if (.not.gv%Boussinesq) eta_rest_name =
"p_bot"
1133 if (gv%Boussinesq)
then
1134 do j=js,je ;
do i=is,ie ; cs%eta(i,j) = -gv%Z_to_H * g%bathyT(i,j) ;
enddo ;
enddo
1136 do k=1,nz ;
do j=js,je ;
do i=is,ie
1137 cs%eta(i,j) = cs%eta(i,j) + h(i,j,k)
1138 enddo ;
enddo ;
enddo
1139 elseif ((gv%m_to_H_restart /= 0.0) .and. (gv%m_to_H_restart /= gv%m_to_H))
then
1140 h_rescale = gv%m_to_H / gv%m_to_H_restart
1141 do j=js,je ;
do i=is,ie ; cs%eta(i,j) = h_rescale * cs%eta(i,j) ;
enddo ;
enddo
1144 do j=js,je ;
do i=is,ie ; eta(i,j) = cs%eta(i,j) ;
enddo ;
enddo
1146 call barotropic_init(u, v, h, cs%eta, time, g, gv, us, param_file, diag, &
1147 cs%barotropic_CSp, restart_cs, calc_dtbt, cs%BT_cont, &
1152 call horizontal_viscosity(u, v, h, cs%diffu, cs%diffv, meke, varmix, &
1153 g, gv, us, cs%hor_visc_CSp, &
1154 obc=cs%OBC, bt=cs%barotropic_CSp)
1157 cs%u_av(:,:,:) = u(:,:,:)
1158 cs%v_av(:,:,:) = v(:,:,:)
1164 h_tmp(:,:,:) = h(:,:,:)
1165 call continuity(u, v, h, h_tmp, uh, vh, dt, g, gv, us, cs%continuity_CSp, obc=cs%OBC)
1166 call pass_var(h_tmp, g%Domain, clock=id_clock_pass_init)
1167 cs%h_av(:,:,:) = 0.5*(h(:,:,:) + h_tmp(:,:,:))
1170 cs%h_av(:,:,:) = h(:,:,:)
1171 elseif ((gv%m_to_H_restart /= 0.0) .and. (gv%m_to_H_restart /= gv%m_to_H))
then
1172 h_rescale = gv%m_to_H / gv%m_to_H_restart
1173 do k=1,nz ;
do j=js,je ;
do i=is,ie ; cs%h_av(i,j,k) = h_rescale * cs%h_av(i,j,k) ;
enddo ;
enddo ;
enddo
1175 if ((gv%m_to_H_restart /= 0.0) .and. (gv%m_to_H_restart /= gv%m_to_H))
then
1176 uh_rescale = gv%m_to_H / gv%m_to_H_restart
1177 do k=1,nz ;
do j=js,je ;
do i=g%IscB,g%IecB ; uh(i,j,k) = uh_rescale * uh(i,j,k) ;
enddo ;
enddo ;
enddo
1178 do k=1,nz ;
do j=g%JscB,g%JecB ;
do i=is,ie ; vh(i,j,k) = uh_rescale * vh(i,j,k) ;
enddo ;
enddo ;
enddo
1182 call cpu_clock_begin(id_clock_pass_init)
1186 call do_group_pass(pass_av_h_uvh, g%Domain)
1187 call cpu_clock_end(id_clock_pass_init)
1189 flux_units = get_flux_units(gv)
1190 h_convert = gv%H_to_m ;
if (.not.gv%Boussinesq) h_convert = gv%H_to_kg_m2
1191 cs%id_uh = register_diag_field(
'ocean_model',
'uh', diag%axesCuL, time, &
1192 'Zonal Thickness Flux', flux_units, y_cell_method=
'sum', v_extensive=.true., &
1193 conversion=h_convert)
1194 cs%id_vh = register_diag_field(
'ocean_model',
'vh', diag%axesCvL, time, &
1195 'Meridional Thickness Flux', flux_units, x_cell_method=
'sum', v_extensive=.true., &
1196 conversion=h_convert)
1198 cs%id_CAu = register_diag_field(
'ocean_model',
'CAu', diag%axesCuL, time, &
1199 'Zonal Coriolis and Advective Acceleration',
'm s-2')
1200 cs%id_CAv = register_diag_field(
'ocean_model',
'CAv', diag%axesCvL, time, &
1201 'Meridional Coriolis and Advective Acceleration',
'm s-2')
1202 cs%id_PFu = register_diag_field(
'ocean_model',
'PFu', diag%axesCuL, time, &
1203 'Zonal Pressure Force Acceleration',
'm s-2')
1204 cs%id_PFv = register_diag_field(
'ocean_model',
'PFv', diag%axesCvL, time, &
1205 'Meridional Pressure Force Acceleration',
'm s-2')
1207 cs%id_uav = register_diag_field(
'ocean_model',
'uav', diag%axesCuL, time, &
1208 'Barotropic-step Averaged Zonal Velocity',
'm s-1')
1209 cs%id_vav = register_diag_field(
'ocean_model',
'vav', diag%axesCvL, time, &
1210 'Barotropic-step Averaged Meridional Velocity',
'm s-1')
1212 cs%id_u_BT_accel = register_diag_field(
'ocean_model',
'u_BT_accel', diag%axesCuL, time, &
1213 'Barotropic Anomaly Zonal Acceleration',
'm s-1')
1214 cs%id_v_BT_accel = register_diag_field(
'ocean_model',
'v_BT_accel', diag%axesCvL, time, &
1215 'Barotropic Anomaly Meridional Acceleration',
'm s-1')
1217 id_clock_cor = cpu_clock_id(
'(Ocean Coriolis & mom advection)', grain=clock_module)
1218 id_clock_continuity = cpu_clock_id(
'(Ocean continuity equation)', grain=clock_module)
1219 id_clock_pres = cpu_clock_id(
'(Ocean pressure force)', grain=clock_module)
1220 id_clock_vertvisc = cpu_clock_id(
'(Ocean vertical viscosity)', grain=clock_module)
1221 id_clock_horvisc = cpu_clock_id(
'(Ocean horizontal viscosity)', grain=clock_module)
1222 id_clock_mom_update = cpu_clock_id(
'(Ocean momentum increments)', grain=clock_module)
1223 id_clock_pass = cpu_clock_id(
'(Ocean message passing)', grain=clock_module)
1224 id_clock_pass_init = cpu_clock_id(
'(Ocean init message passing)', grain=clock_routine)
1225 id_clock_btcalc = cpu_clock_id(
'(Ocean barotropic mode calc)', grain=clock_module)
1226 id_clock_btstep = cpu_clock_id(
'(Ocean barotropic mode stepping)', grain=clock_module)
1227 id_clock_btforce = cpu_clock_id(
'(Ocean barotropic forcing calc)', grain=clock_module)
1229 end subroutine initialize_dyn_split_rk2
1233 subroutine end_dyn_split_rk2(CS)
1236 dealloc_(cs%diffu) ; dealloc_(cs%diffv)
1237 dealloc_(cs%CAu) ; dealloc_(cs%CAv)
1238 dealloc_(cs%PFu) ; dealloc_(cs%PFv)
1240 if (
associated(cs%taux_bot))
deallocate(cs%taux_bot)
1241 if (
associated(cs%tauy_bot))
deallocate(cs%tauy_bot)
1242 dealloc_(cs%uhbt) ; dealloc_(cs%vhbt)
1243 dealloc_(cs%u_accel_bt) ; dealloc_(cs%v_accel_bt)
1244 dealloc_(cs%visc_rem_u) ; dealloc_(cs%visc_rem_v)
1246 dealloc_(cs%eta) ; dealloc_(cs%eta_PF) ; dealloc_(cs%pbce)
1247 dealloc_(cs%h_av) ; dealloc_(cs%u_av) ; dealloc_(cs%v_av)
1249 call dealloc_bt_cont_type(cs%BT_cont)
1252 end subroutine end_dyn_split_rk2