394 type(verticalgrid_type),
intent(in) :: gv
395 logical,
intent(in) :: use_fluxes
397 real,
dimension(SZIB_(G),SZJ_(G),SZK_(G)), &
400 real,
dimension(SZI_(G),SZJB_(G),SZK_(G)), &
403 real,
dimension(SZI_(G),SZJ_(G)), &
409 real,
dimension(SZIB_(G),SZJ_(G),SZK_(G)), &
410 intent(in) :: bc_accel_u
412 real,
dimension(SZI_(G),SZJB_(G),SZK_(G)), &
413 intent(in) :: bc_accel_v
415 type(forcing),
intent(in) :: fluxes
417 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)), &
420 real,
dimension(SZI_(G),SZJ_(G)), &
421 intent(in) :: eta_pf_in
424 real,
dimension(SZIB_(G),SZJ_(G),SZK_(G)), &
428 real,
dimension(SZI_(G),SZJB_(G),SZK_(G)), &
432 real,
dimension(SZIB_(G),SZJ_(G),SZK_(G)), &
433 intent(out) :: accel_layer_u
435 real,
dimension(SZI_(G),SZJB_(G),SZK_(G)), &
436 intent(out) :: accel_layer_v
438 real,
dimension(SZI_(G),SZJ_(G)), &
439 intent(inout) :: eta_out
441 real,
dimension(SZIB_(G),SZJ_(G)), &
442 intent(out) :: uhbtav
444 real,
dimension(SZI_(G),SZJB_(G)), &
445 intent(out) :: vhbtav
447 type(legacy_barotropic_cs),
pointer :: cs
449 real,
dimension(SZIB_(G),SZJ_(G),SZK_(G)), &
450 intent(in),
optional :: visc_rem_u
455 real,
dimension(SZI_(G),SZJB_(G),SZK_(G)), &
456 intent(in),
optional :: visc_rem_v
461 real,
dimension(SZI_(G),SZJ_(G)), &
462 intent(out),
optional :: etaav
464 real,
dimension(SZIB_(G),SZJ_(G)), &
465 intent(out),
optional :: uhbt_out
467 real,
dimension(SZI_(G),SZJB_(G)), &
468 intent(out),
optional :: vhbt_out
470 type(ocean_obc_type), &
471 pointer,
optional :: obc
473 type(bt_cont_type), &
474 pointer,
optional :: bt_cont
476 real,
dimension(:,:), &
477 pointer,
optional :: eta_pf_start
479 real,
dimension(:,:), &
480 pointer,
optional :: taux_bot
482 real,
dimension(:,:), &
483 pointer,
optional :: tauy_bot
485 real,
dimension(:,:,:), &
486 pointer,
optional :: uh0, u_uh0
487 real,
dimension(:,:,:), &
488 pointer,
optional :: vh0, v_vh0
559 real :: ubt_cor(szib_(g),szj_(g))
560 real :: vbt_cor(szi_(g),szjb_(g))
562 real :: wt_u(szib_(g),szj_(g),szk_(g))
563 real :: wt_v(szi_(g),szjb_(g),szk_(g))
566 real,
dimension(SZIB_(G),SZJ_(G)) :: &
569 real,
dimension(SZI_(G),SZJB_(G)) :: &
572 real,
dimension(SZI_(G),SZJ_(G)) :: &
579 real :: q(szibw_(cs),szjbw_(cs))
580 real,
dimension(SZIBW_(CS),SZJW_(CS)) :: &
609 real,
dimension(SZIW_(CS),SZJBW_(CS)) :: &
636 real,
target,
dimension(SZIW_(CS),SZJW_(CS)) :: &
640 real,
pointer,
dimension(:,:) :: &
643 real,
dimension(SZIW_(CS),SZJW_(CS)) :: &
663 type(local_bt_cont_u_type),
dimension(SZIBW_(CS),SZJW_(CS)) :: &
665 type(local_bt_cont_v_type),
dimension(SZIW_(CS),SZJBW_(CS)) :: &
697 logical :: do_hifreq_output
698 logical :: use_visc_rem, use_bt_cont
699 logical :: do_ave, find_etaav, find_pf, find_cor
700 logical :: ice_is_rigid, nonblock_setup, interp_eta_pf
701 logical :: project_velocity, add_uh0
704 real :: ice_strength = 0.0
713 real,
allocatable,
dimension(:) :: wt_vel, wt_eta, wt_accel, wt_trans, wt_accel2
714 real :: sum_wt_vel, sum_wt_eta, sum_wt_accel, sum_wt_trans
715 real :: i_sum_wt_vel, i_sum_wt_eta, i_sum_wt_accel, i_sum_wt_trans
719 logical :: apply_obcs, apply_obc_flather
720 type(bt_obc_type) :: bt_obc
722 type(memory_size_type) :: ms
723 character(len=200) :: mesg
724 integer :: pid_ubt, pid_eta, pid_e_anom, pid_etaav, pid_uhbtav, pid_ubtav
725 integer :: pid_q, pid_eta_pf, pid_dyn_coef_eta, pid_eta_src
726 integer :: pid_dcor_u, pid_datu_res, pid_tmp_u, pid_gtot_e, pid_gtot_w
727 integer :: pid_bt_rem_u, pid_datu, pid_bt_force_u, pid_cor_ref
728 integer :: pid_eta_pf_1, pid_d_eta_pf, pid_uhbt0
729 integer :: isv, iev, jsv, jev
731 integer :: isvf, ievf, jsvf, jevf, num_cycles
732 integer :: i, j, k, n
733 integer :: is, ie, js, je, nz, isq, ieq, jsq, jeq
734 integer :: isd, ied, jsd, jed, isdb, iedb, jsdb, jedb
736 if (.not.
associated(cs))
call mom_error(fatal, &
737 "legacy_btstep: Module MOM_legacy_barotropic must be initialized before it is used.")
738 if (.not.cs%split)
return 739 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
740 isq = g%IscB ; ieq = g%IecB ; jsq = g%JscB ; jeq = g%JecB
741 isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
742 isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
743 ms%isdw = cs%isdw ; ms%iedw = cs%iedw ; ms%jsdw = cs%jsdw ; ms%jedw = cs%jedw
746 use_bt_cont = .false.
747 if (
present(bt_cont)) use_bt_cont = (
associated(bt_cont))
749 interp_eta_pf = .false.
750 if (
present(eta_pf_start)) interp_eta_pf = (
associated(eta_pf_start))
752 project_velocity = cs%BT_project_velocity
756 if ((.not.use_bt_cont) .and. cs%Nonlinear_continuity .and. &
757 (cs%Nonlin_cont_update_period > 0)) stencil = 2
760 if (cs%use_wide_halos) &
761 num_cycles = min((is-cs%isdw) / stencil, (js-cs%jsdw) / stencil)
762 isvf = is - (num_cycles-1)*stencil ; ievf = ie + (num_cycles-1)*stencil
763 jsvf = js - (num_cycles-1)*stencil ; jevf = je + (num_cycles-1)*stencil
765 do_ave = query_averaging_enabled(cs%diag)
766 find_etaav =
present(etaav)
767 use_visc_rem =
present(visc_rem_u)
768 if ((use_visc_rem) .neqv.
present(visc_rem_v))
call mom_error(fatal, &
769 "btstep: Either both visc_rem_u and visc_rem_v or neither"// &
770 " one must be present in call to btstep.")
771 find_pf = (do_ave .and. ((cs%id_PFu_bt > 0) .or. (cs%id_PFv_bt > 0)))
772 find_cor = (do_ave .and. ((cs%id_Coru_bt > 0) .or. (cs%id_Corv_bt > 0)))
775 if (
present(uh0)) add_uh0 =
associated(uh0)
776 if (add_uh0 .and. .not.(
present(vh0) .and.
present(u_uh0) .and. &
777 present(v_vh0)))
call mom_error(fatal, &
778 "legacy_btstep: vh0, u_uh0, and v_vh0 must be present if uh0 is used.")
779 if (add_uh0 .and. .not.(
associated(vh0) .and.
associated(u_uh0) .and. &
780 associated(v_vh0)))
call mom_error(fatal, &
781 "legacy_btstep: vh0, u_uh0, and v_vh0 must be associated if uh0 is used.")
782 if (add_uh0 .and. use_fluxes)
call mom_error(warning, &
783 "legacy_btstep: with use_fluxes, add_uh0 does nothing!")
784 if (use_fluxes) add_uh0 = .false.
787 nonblock_setup = g%nonblocking_updates
789 if (id_clock_calc_pre > 0)
call cpu_clock_begin(id_clock_calc_pre)
791 apply_obcs = .false. ; apply_u_obcs = .false. ; apply_v_obcs = .false.
792 apply_obc_flather = .false.
793 if (
present(obc))
then ;
if (
associated(obc))
then 794 apply_u_obcs = obc%Flather_u_BCs_exist_globally .or. obc%specified_u_BCs_exist_globally
795 apply_v_obcs = obc%Flather_v_BCs_exist_globally .or. obc%specified_v_BCs_exist_globally
796 apply_obc_flather = obc%Flather_u_BCs_exist_globally .or. obc%Flather_v_BCs_exist_globally
797 apply_obcs = obc%specified_u_BCs_exist_globally .or. obc%specified_v_BCs_exist_globally .or. apply_obc_flather
799 if (apply_obc_flather .and. .not.gv%Boussinesq)
call mom_error(fatal, &
800 "legacy_btstep: Flather open boundary conditions have not yet been "// &
801 "implemented for a non-Boussinesq model.")
804 nstep = ceiling(dt/cs%dtbt - 0.0001)
805 if (is_root_pe() .and. (nstep /= cs%nstep_last))
then 806 write(mesg,
'("legacy_btstep is using a dynamic barotropic timestep of ", ES12.6, & 807 & " seconds, max ", ES12.6, ".")') (dt/nstep), cs%dtbt_max
808 call mom_mesg(mesg, 3)
810 cs%nstep_last = nstep
813 instep = 1.0 /
real(nstep)
818 do_ave = query_averaging_enabled(cs%diag)
820 do_hifreq_output = .false.
821 if ((cs%id_ubt_hifreq > 0) .or. (cs%id_vbt_hifreq > 0) .or. &
822 (cs%id_eta_hifreq > 0) .or. (cs%id_eta_pred_hifreq > 0) .or. &
823 (cs%id_uhbt_hifreq > 0) .or. (cs%id_vhbt_hifreq > 0))
then 824 do_hifreq_output = query_averaging_enabled(cs%diag, time_int_in, time_end_in)
825 if (do_hifreq_output) &
826 time_bt_start = time_end_in - set_time(int(floor(dt+0.5)))
832 if (cs%linearized_BT_PV)
then 834 do j=jsvf-2,jevf+1 ;
do i=isvf-2,ievf+1
838 do j=jsvf-1,jevf+1 ;
do i=isvf-2,ievf+1
839 dcor_u(i,j) = cs%D_u_Cor(i,j)
842 do j=jsvf-2,jevf+1 ;
do i=isvf-1,ievf+1
843 dcor_v(i,j) = cs%D_v_Cor(i,j)
846 q(:,:) = 0.0 ; dcor_u(:,:) = 0.0 ; dcor_v(:,:) = 0.0
850 do j=js,je ;
do i=is-1,ie
851 dcor_u(i,j) = 0.5 * (g%bathyT(i+1,j) + g%bathyT(i,j))
854 do j=js-1,je ;
do i=is,ie
855 dcor_v(i,j) = 0.5 * (g%bathyT(i,j+1) + g%bathyT(i,j))
858 do j=js-1,je ;
do i=is-1,ie
859 q(i,j) = 0.25 * g%CoriolisBu(i,j) * &
860 ((g%areaT(i,j) + g%areaT(i+1,j+1)) + (g%areaT(i+1,j) + g%areaT(i,j+1))) / &
861 ((g%areaT(i,j) * g%bathyT(i,j) + g%areaT(i+1,j+1) * g%bathyT(i+1,j+1)) + &
862 (g%areaT(i+1,j) * g%bathyT(i+1,j) + g%areaT(i,j+1) * g%bathyT(i,j+1)))
868 if (id_clock_calc_pre > 0)
call cpu_clock_end(id_clock_calc_pre)
869 if (id_clock_pass_pre > 0)
call cpu_clock_begin(id_clock_pass_pre)
870 if (nonblock_setup)
then 871 pid_q = pass_var_start(q, cs%BT_Domain, to_all, position=corner)
872 pid_dcor_u = pass_vector_start(dcor_u, dcor_v, cs%BT_Domain, to_all+scalar_pair)
874 call pass_var(q, cs%BT_Domain, to_all, position=corner)
875 call pass_vector(dcor_u, dcor_v, cs%BT_Domain, to_all+scalar_pair)
877 if (id_clock_pass_pre > 0)
call cpu_clock_end(id_clock_pass_pre)
878 if (id_clock_calc_pre > 0)
call cpu_clock_begin(id_clock_calc_pre)
881 if (nonblock_setup)
then 884 if ((.not.use_bt_cont) .and. cs%rescale_D_bt .and. (ievf>ie))
then 885 if (id_clock_calc_pre > 0)
call cpu_clock_end(id_clock_calc_pre)
886 if (id_clock_pass_pre > 0)
call cpu_clock_begin(id_clock_pass_pre)
887 pid_datu_res = pass_vector_start(cs%Datu_res, cs%Datv_res, cs%BT_Domain, &
889 if (id_clock_pass_pre > 0)
call cpu_clock_end(id_clock_pass_pre)
890 if (id_clock_calc_pre > 0)
call cpu_clock_begin(id_clock_calc_pre)
895 do j=cs%jsdw,cs%jedw ;
do i=cs%isdw,cs%iedw
896 gtot_e(i,j) = 0.0 ; gtot_w(i,j) = 0.0
897 gtot_n(i,j) = 0.0 ; gtot_s(i,j) = 0.0
900 if (interp_eta_pf)
then 901 eta_pf_1(i,j) = 0.0 ; d_eta_pf(i,j) = 0.0
903 p_surf_dyn(i,j) = 0.0
904 if (cs%dynamic_psurf) dyn_coef_eta(i,j) = 0.0
909 do j=cs%jsdw,cs%jedw ;
do i=cs%isdw-1,cs%iedw
910 cor_ref_u(i,j) = 0.0 ; bt_force_u(i,j) = 0.0 ; ubt(i,j) = 0.0
911 datu(i,j) = 0.0 ; bt_rem_u(i,j) = 0.0 ; uhbt0(i,j) = 0.0
913 do j=cs%jsdw-1,cs%jedw ;
do i=cs%isdw,cs%iedw
914 cor_ref_v(i,j) = 0.0 ; bt_force_v(i,j) = 0.0 ; vbt(i,j) = 0.0
915 datv(i,j) = 0.0 ; bt_rem_v(i,j) = 0.0 ; vhbt0(i,j) = 0.0
920 if (interp_eta_pf)
then 921 do j=g%jsd,g%jed ;
do i=g%isd,g%ied
922 eta(i,j) = eta_in(i,j)
923 eta_pf_1(i,j) = eta_pf_start(i,j)
924 d_eta_pf(i,j) = eta_pf_in(i,j) - eta_pf_start(i,j)
927 do j=g%jsd,g%jed ;
do i=g%isd,g%ied
928 eta(i,j) = eta_in(i,j)
929 eta_pf(i,j) = eta_pf_in(i,j)
933 if (use_visc_rem)
then 935 do k=1,nz ;
do j=js-1,je+1 ;
do i=isq-1,ieq+1
938 if (visc_rem_u(i,j,k) <= 0.0)
then ; visc_rem = 0.0
939 elseif (visc_rem_u(i,j,k) >= 1.0)
then ; visc_rem = 1.0
940 elseif (visc_rem_u(i,j,k)**2 > visc_rem_u(i,j,k) - 0.5*instep)
then 941 visc_rem = visc_rem_u(i,j,k)
942 else ; visc_rem = 1.0 - 0.5*instep/visc_rem_u(i,j,k) ;
endif 943 wt_u(i,j,k) = cs%frhatu(i,j,k) * visc_rem
944 enddo ;
enddo ;
enddo 946 do k=1,nz ;
do j=jsq-1,jeq+1 ;
do i=is-1,ie+1
948 if (visc_rem_v(i,j,k) <= 0.0)
then ; visc_rem = 0.0
949 elseif (visc_rem_v(i,j,k) >= 1.0)
then ; visc_rem = 1.0
950 elseif (visc_rem_v(i,j,k)**2 > visc_rem_v(i,j,k) - 0.5*instep)
then 951 visc_rem = visc_rem_v(i,j,k)
952 else ; visc_rem = 1.0 - 0.5*instep/visc_rem_v(i,j,k) ;
endif 953 wt_v(i,j,k) = cs%frhatv(i,j,k) * visc_rem
954 enddo ;
enddo ;
enddo 956 do k=1,nz ;
do j=js-1,je+1 ;
do i=isq-1,ieq+1
957 wt_u(i,j,k) = cs%frhatu(i,j,k)
958 enddo ;
enddo ;
enddo 959 do k=1,nz ;
do j=jsq-1,jeq+1 ;
do i=is-1,ie+1
960 wt_v(i,j,k) = cs%frhatv(i,j,k)
961 enddo ;
enddo ;
enddo 969 do j=js,je ;
do i=is-1,ie
970 gtot_e(i,j) = gtot_e(i,j) + pbce(i,j,k) * wt_u(i,j,k)
971 gtot_w(i+1,j) = gtot_w(i+1,j) + pbce(i+1,j,k) * wt_u(i,j,k)
973 do j=js-1,je ;
do i=is,ie
974 gtot_n(i,j) = gtot_n(i,j) + pbce(i,j,k) * wt_v(i,j,k)
975 gtot_s(i,j+1) = gtot_s(i,j+1) + pbce(i,j+1,k) * wt_v(i,j,k)
980 call tidal_forcing_sensitivity(g, cs%tides_CSp, det_de)
981 dgeo_de = 1.0 + det_de + cs%G_extra
983 dgeo_de = 1.0 + cs%G_extra
986 if (nonblock_setup .and. .not.cs%linearized_BT_PV)
then 987 if (id_clock_calc_pre > 0)
call cpu_clock_end(id_clock_calc_pre)
988 if (id_clock_pass_pre > 0)
call cpu_clock_begin(id_clock_pass_pre)
989 call pass_var_complete(pid_q, q, cs%BT_Domain, to_all, position=corner)
990 call pass_vector_complete(pid_dcor_u, dcor_u, dcor_v, cs%BT_Domain, to_all+scalar_pair)
991 if (id_clock_pass_pre > 0)
call cpu_clock_end(id_clock_pass_pre)
992 if (id_clock_calc_pre > 0)
call cpu_clock_begin(id_clock_calc_pre)
997 if (use_bt_cont)
then 998 call set_local_bt_cont_types(bt_cont, btcl_u, btcl_v, g, ms, cs%BT_Domain, 1+ievf-ie)
1000 if (cs%rescale_D_bt .and. (ievf>ie))
then 1003 if (id_clock_calc_pre > 0)
call cpu_clock_end(id_clock_calc_pre)
1004 if (id_clock_pass_pre > 0)
call cpu_clock_begin(id_clock_pass_pre)
1005 if (nonblock_setup)
then 1006 call pass_vector_complete(pid_datu_res, cs%Datu_res, cs%Datv_res, &
1007 cs%BT_Domain, to_all+scalar_pair)
1009 call pass_vector(cs%Datu_res, cs%Datv_res, cs%BT_Domain, to_all+scalar_pair)
1011 if (id_clock_pass_pre > 0)
call cpu_clock_end(id_clock_pass_pre)
1012 if (id_clock_calc_pre > 0)
call cpu_clock_begin(id_clock_calc_pre)
1014 if (cs%Nonlinear_continuity)
then 1015 call find_face_areas(datu, datv, g, gv, cs, ms, cs%rescale_D_bt, eta, 1)
1017 call find_face_areas(datu, datv, g, gv, cs, ms, cs%rescale_D_bt, halo=1)
1022 if (apply_obcs)
then 1023 call set_up_bt_obc(obc, eta, bt_obc, g, gv, ms, ievf-ie, use_bt_cont, &
1024 datu, datv, btcl_u, btcl_v)
1033 if (use_visc_rem)
then 1034 do j=js,je ;
do i=is-1,ie
1035 bt_force_u(i,j) = fluxes%taux(i,j) * i_rho0*cs%IDatu(i,j)*visc_rem_u(i,j,1)
1037 do j=js-1,je ;
do i=is,ie
1038 bt_force_v(i,j) = fluxes%tauy(i,j) * i_rho0*cs%IDatv(i,j)*visc_rem_v(i,j,1)
1041 do j=js,je ;
do i=is-1,ie
1042 bt_force_u(i,j) = fluxes%taux(i,j) * i_rho0 * cs%IDatu(i,j)
1044 do j=js-1,je ;
do i=is,ie
1045 bt_force_v(i,j) = fluxes%tauy(i,j) * i_rho0 * cs%IDatv(i,j)
1048 if (
present(taux_bot) .and.
present(tauy_bot))
then 1049 if (
associated(taux_bot) .and.
associated(tauy_bot))
then 1050 do j=js,je ;
do i=is-1,ie
1051 bt_force_u(i,j) = bt_force_u(i,j) - taux_bot(i,j) * i_rho0 * cs%IDatu(i,j)
1053 do j=js-1,je ;
do i=is,ie
1054 bt_force_v(i,j) = bt_force_v(i,j) - tauy_bot(i,j) * i_rho0 * cs%IDatv(i,j)
1061 do k=1,nz ;
do j=js,je ;
do i=isq,ieq
1062 bt_force_u(i,j) = bt_force_u(i,j) + wt_u(i,j,k) * bc_accel_u(i,j,k)
1063 enddo ;
enddo ;
enddo 1064 do k=1,nz ;
do j=jsq,jeq ;
do i=is,ie
1065 bt_force_v(i,j) = bt_force_v(i,j) + wt_v(i,j,k) * bc_accel_v(i,j,k)
1066 enddo ;
enddo ;
enddo 1071 do j=js,je ;
do i=is-1,ie ; uhbt(i,j) = 0.0 ; ubt(i,j) = 0.0 ;
enddo ;
enddo 1072 do j=js-1,je ;
do i=is,ie ; vhbt(i,j) = 0.0 ; vbt(i,j) = 0.0 ;
enddo ;
enddo 1073 do k=1,nz ;
do j=js,je ;
do i=is-1,ie
1074 uhbt(i,j) = uhbt(i,j) + uh0(i,j,k)
1075 ubt(i,j) = ubt(i,j) + wt_u(i,j,k) * u_uh0(i,j,k)
1076 enddo ;
enddo ;
enddo 1077 do k=1,nz ;
do j=js-1,je ;
do i=is,ie
1078 vhbt(i,j) = vhbt(i,j) + vh0(i,j,k)
1079 vbt(i,j) = vbt(i,j) + wt_v(i,j,k) * v_vh0(i,j,k)
1080 enddo ;
enddo ;
enddo 1081 if (use_bt_cont)
then 1082 do j=js,je ;
do i=is-1,ie
1083 uhbt0(i,j) = uhbt(i,j) - find_uhbt(ubt(i,j),btcl_u(i,j))
1085 do j=js-1,je ;
do i=is,ie
1086 vhbt0(i,j) = vhbt(i,j) - find_vhbt(vbt(i,j),btcl_v(i,j))
1089 do j=js,je ;
do i=is-1,ie
1090 uhbt0(i,j) = uhbt(i,j) - datu(i,j)*ubt(i,j)
1092 do j=js-1,je ;
do i=is,ie
1093 vhbt0(i,j) = vhbt(i,j) - datv(i,j)*vbt(i,j)
1099 do j=jsvf-1,jevf+1 ;
do i=isvf-2,ievf+1
1100 ubt(i,j) = 0.0 ; uhbt(i,j) = 0.0 ; u_accel_bt(i,j) = 0.0
1102 do j=jsvf-2,jevf+1 ;
do i=isvf-1,ievf+1
1103 vbt(i,j) = 0.0 ; vhbt(i,j) = 0.0 ; v_accel_bt(i,j) = 0.0
1105 if (use_fluxes)
then 1106 do k=1,nz ;
do j=js,je ;
do i=is-1,ie
1107 uhbt(i,j) = uhbt(i,j) + u_in(i,j,k)
1108 enddo ;
enddo ;
enddo 1109 do k=1,nz ;
do j=js-1,je ;
do i=is,ie
1110 vhbt(i,j) = vhbt(i,j) + v_in(i,j,k)
1111 enddo ;
enddo ;
enddo 1112 if (use_bt_cont)
then 1113 do j=js,je ;
do i=is-1,ie
1114 ubt(i,j) = uhbt_to_ubt(uhbt(i,j),btcl_u(i,j), guess=cs%ubt_IC(i,j))
1116 do j=js-1,je ;
do i=is,ie
1117 vbt(i,j) = vhbt_to_vbt(vhbt(i,j),btcl_v(i,j), guess=cs%vbt_IC(i,j))
1120 do j=js,je ;
do i=is-1,ie
1121 ubt(i,j) = 0.0 ;
if (datu(i,j) > 0.0) ubt(i,j) = uhbt(i,j) / datu(i,j)
1123 do j=js-1,je ;
do i=is,ie
1124 vbt(i,j) = 0.0 ;
if (datv(i,j) > 0.0) vbt(i,j) = vhbt(i,j) / datv(i,j)
1128 do k=1,nz ;
do j=js,je ;
do i=is-1,ie
1129 ubt(i,j) = ubt(i,j) + wt_u(i,j,k) * u_in(i,j,k)
1130 enddo ;
enddo ;
enddo 1131 do k=1,nz ;
do j=js-1,je ;
do i=is,ie
1132 vbt(i,j) = vbt(i,j) + wt_v(i,j,k) * v_in(i,j,k)
1133 enddo ;
enddo ;
enddo 1136 if (cs%gradual_BT_ICs)
then 1137 if (use_fluxes)
then 1138 if (use_bt_cont)
then 1139 do j=js,je ;
do i=is-1,ie
1140 vel_tmp = uhbt_to_ubt(cs%uhbt_IC(i,j),btcl_u(i,j), guess=cs%ubt_IC(i,j))
1141 bt_force_u(i,j) = bt_force_u(i,j) + (ubt(i,j) - vel_tmp) * idt
1144 do j=js-1,je ;
do i=is,ie
1145 vel_tmp = vhbt_to_vbt(cs%vhbt_IC(i,j),btcl_v(i,j), guess=cs%vbt_IC(i,j))
1146 bt_force_v(i,j) = bt_force_v(i,j) + (vbt(i,j) - vel_tmp) * idt
1150 do j=js,je ;
do i=is-1,ie
1151 vel_tmp = 0.0 ;
if (datu(i,j) > 0.0) vel_tmp = cs%uhbt_IC(i,j) / datu(i,j)
1152 bt_force_u(i,j) = bt_force_u(i,j) + (ubt(i,j) - vel_tmp) * idt
1155 do j=js-1,je ;
do i=is,ie
1156 vel_tmp = 0.0 ;
if (datv(i,j) > 0.0) vel_tmp = cs%vhbt_IC(i,j) / datv(i,j)
1157 bt_force_v(i,j) = bt_force_v(i,j) + (vbt(i,j) - vel_tmp) * idt
1162 do j=js,je ;
do i=is-1,ie
1163 bt_force_u(i,j) = bt_force_u(i,j) + (ubt(i,j) - cs%ubt_IC(i,j)) * idt
1164 ubt(i,j) = cs%ubt_IC(i,j)
1166 do j=js-1,je ;
do i=is,ie
1167 bt_force_v(i,j) = bt_force_v(i,j) + (vbt(i,j) - cs%vbt_IC(i,j)) * idt
1168 vbt(i,j) = cs%vbt_IC(i,j)
1173 if ((isq > is-1) .or. (jsq > js-1))
then 1176 if (id_clock_calc_pre > 0)
call cpu_clock_end(id_clock_calc_pre)
1177 if (id_clock_pass_pre > 0)
call cpu_clock_begin(id_clock_pass_pre)
1178 tmp_u(:,:) = 0.0 ; tmp_v(:,:) = 0.0
1179 do j=js,je ;
do i=isq,ieq ; tmp_u(i,j) = bt_force_u(i,j) ;
enddo ;
enddo 1180 do j=jsq,jeq ;
do i=is,ie ; tmp_v(i,j) = bt_force_v(i,j) ;
enddo ;
enddo 1181 if (nonblock_setup)
then 1182 pid_tmp_u = pass_vector_start(tmp_u, tmp_v, g%Domain)
1184 call pass_vector(tmp_u, tmp_v, g%Domain, complete=.true.)
1185 do j=jsd,jed ;
do i=isdb,iedb ; bt_force_u(i,j) = tmp_u(i,j) ;
enddo ;
enddo 1186 do j=jsdb,jedb ;
do i=isd,ied ; bt_force_v(i,j) = tmp_v(i,j) ;
enddo ;
enddo 1188 if (id_clock_pass_pre > 0)
call cpu_clock_end(id_clock_pass_pre)
1189 if (id_clock_calc_pre > 0)
call cpu_clock_begin(id_clock_calc_pre)
1192 if (nonblock_setup)
then 1193 if (id_clock_calc_pre > 0)
call cpu_clock_end(id_clock_calc_pre)
1194 if (id_clock_pass_pre > 0)
call cpu_clock_begin(id_clock_pass_pre)
1195 pid_gtot_e = pass_vector_start(gtot_e, gtot_n, cs%BT_Domain, &
1196 to_all+scalar_pair, agrid, complete=.false.)
1197 pid_gtot_w = pass_vector_start(gtot_w, gtot_s, cs%BT_Domain, &
1198 to_all+scalar_pair, agrid)
1199 if (id_clock_pass_pre > 0)
call cpu_clock_end(id_clock_pass_pre)
1200 if (id_clock_calc_pre > 0)
call cpu_clock_begin(id_clock_calc_pre)
1205 do j=jsvf-1,jevf ;
do i=isvf-1,ievf+1
1206 if (cs%Sadourny)
then 1207 amer(i-1,j) = dcor_u(i-1,j) * q(i-1,j)
1208 bmer(i,j) = dcor_u(i,j) * q(i,j)
1209 cmer(i,j+1) = dcor_u(i,j+1) * q(i,j)
1210 dmer(i-1,j+1) = dcor_u(i-1,j+1) * q(i-1,j)
1212 amer(i-1,j) = dcor_u(i-1,j) * &
1213 ((q(i,j) + q(i-1,j-1)) + q(i-1,j)) / 3.0
1214 bmer(i,j) = dcor_u(i,j) * &
1215 (q(i,j) + (q(i-1,j) + q(i,j-1))) / 3.0
1216 cmer(i,j+1) = dcor_u(i,j+1) * &
1217 (q(i,j) + (q(i-1,j) + q(i,j+1))) / 3.0
1218 dmer(i-1,j+1) = dcor_u(i-1,j+1) * &
1219 ((q(i,j) + q(i-1,j+1)) + q(i-1,j)) / 3.0
1224 do j=jsvf-1,jevf+1 ;
do i=isvf-1,ievf
1225 if (cs%Sadourny)
then 1226 azon(i,j) = dcor_v(i+1,j) * q(i,j)
1227 bzon(i,j) = dcor_v(i,j) * q(i,j)
1228 czon(i,j) = dcor_v(i,j-1) * q(i,j-1)
1229 dzon(i,j) = dcor_v(i+1,j-1) * q(i,j-1)
1231 azon(i,j) = dcor_v(i+1,j) * &
1232 (q(i,j) + (q(i+1,j) + q(i,j-1))) / 3.0
1233 bzon(i,j) = dcor_v(i,j) * &
1234 (q(i,j) + (q(i-1,j) + q(i,j-1))) / 3.0
1235 czon(i,j) = dcor_v(i,j-1) * &
1236 ((q(i,j) + q(i-1,j-1)) + q(i,j-1)) / 3.0
1237 dzon(i,j) = dcor_v(i+1,j-1) * &
1238 ((q(i,j) + q(i+1,j-1)) + q(i,j-1)) / 3.0
1244 if (use_fluxes)
then 1245 do j=js-1,je+1 ;
do i=is-1,ie ; uhbt(i,j) = 0.0 ;
enddo ;
enddo 1246 do j=js-1,je ;
do i=is-1,ie+1 ; vhbt(i,j) = 0.0 ;
enddo ;
enddo 1247 do k=1,nz ;
do j=js-1,je+1 ;
do i=is-1,ie
1248 uhbt(i,j) = uhbt(i,j) + u_cor(i,j,k)
1249 enddo ;
enddo ;
enddo 1250 do k=1,nz ;
do j=js-1,je ;
do i=is-1,ie+1
1251 vhbt(i,j) = vhbt(i,j) + v_cor(i,j,k)
1252 enddo ;
enddo ;
enddo 1253 if (use_bt_cont)
then 1254 do j=js-1,je+1 ;
do i=is-1,ie
1255 ubt_cor(i,j) = uhbt_to_ubt(uhbt(i,j), btcl_u(i,j), guess=cs%ubtav(i,j))
1257 do j=js-1,je ;
do i=is-1,ie+1
1258 vbt_cor(i,j) = vhbt_to_vbt(vhbt(i,j), btcl_v(i,j), guess=cs%vbtav(i,j))
1261 do j=js-1,je+1 ;
do i=is-1,ie
1262 ubt_cor(i,j) = 0.0 ;
if (datu(i,j)>0.0) ubt_cor(i,j) = uhbt(i,j)/datu(i,j)
1264 do j=js-1,je ;
do i=is-1,ie+1
1265 vbt_cor(i,j) = 0.0 ;
if (datv(i,j)>0.0) vbt_cor(i,j) = vhbt(i,j)/datv(i,j)
1269 do j=js-1,je+1 ;
do i=is-1,ie ; ubt_cor(i,j) = 0.0 ;
enddo ;
enddo 1270 do j=js-1,je ;
do i=is-1,ie+1 ; vbt_cor(i,j) = 0.0 ;
enddo ;
enddo 1271 do k=1,nz ;
do j=js-1,je+1 ;
do i=is-1,ie
1272 ubt_cor(i,j) = ubt_cor(i,j) + wt_u(i,j,k) * u_cor(i,j,k)
1273 enddo ;
enddo ;
enddo 1274 do k=1,nz ;
do j=js-1,je ;
do i=is-1,ie+1
1275 vbt_cor(i,j) = vbt_cor(i,j) + wt_v(i,j,k) * v_cor(i,j,k)
1276 enddo ;
enddo ;
enddo 1280 do j=js,je ;
do i=is-1,ie
1282 ((azon(i,j) * vbt_cor(i+1,j) + czon(i,j) * vbt_cor(i ,j-1)) + &
1283 (bzon(i,j) * vbt_cor(i ,j) + dzon(i,j) * vbt_cor(i+1,j-1)))
1286 do j=js-1,je ;
do i=is,ie
1287 cor_ref_v(i,j) = -1.0 * &
1288 ((amer(i-1,j) * ubt_cor(i-1,j) + cmer(i ,j+1) * ubt_cor(i ,j+1)) + &
1289 (bmer(i ,j) * ubt_cor(i ,j) + dmer(i-1,j+1) * ubt_cor(i-1,j+1)))
1293 if (id_clock_calc_pre > 0)
call cpu_clock_end(id_clock_calc_pre)
1294 if (id_clock_pass_pre > 0)
call cpu_clock_begin(id_clock_pass_pre)
1295 if (nonblock_setup)
then 1296 if ((isq > is-1) .or. (jsq > js-1))
then 1297 call pass_vector_complete(pid_tmp_u, tmp_u, tmp_v, g%Domain)
1298 do j=jsd,jed ;
do i=isdb,iedb ; bt_force_u(i,j) = tmp_u(i,j) ;
enddo ;
enddo 1299 do j=jsdb,jedb ;
do i=isd,ied ; bt_force_v(i,j) = tmp_v(i,j) ;
enddo ;
enddo 1301 call pass_vector_complete(pid_gtot_e, gtot_e, gtot_n, cs%BT_Domain, to_all+scalar_pair, agrid)
1302 call pass_vector_complete(pid_gtot_w, gtot_w, gtot_s, cs%BT_Domain, to_all+scalar_pair, agrid)
1304 call pass_vector(gtot_e, gtot_n, cs%BT_Domain, to_all+scalar_pair, agrid, complete=.false.)
1305 call pass_vector(gtot_w, gtot_s, cs%BT_Domain, to_all+scalar_pair, agrid, complete=.true.)
1309 do j=jsvf-1,jevf+1 ;
do i=isvf-1,ievf+1
1310 if (cs%ua_polarity(i,j) < 0.0)
call swap(gtot_e(i,j), gtot_w(i,j))
1311 if (cs%va_polarity(i,j) < 0.0)
call swap(gtot_n(i,j), gtot_s(i,j))
1315 if (nonblock_setup)
then 1316 if (.not.use_bt_cont) &
1317 pid_datu = pass_vector_start(datu, datv, cs%BT_Domain, to_all+scalar_pair)
1319 if (((g%isd > cs%isdw) .or. (g%jsd > cs%jsdw)) .or. (isq <= is-1) .or. (jsq <= js-1)) &
1320 pid_bt_force_u = pass_vector_start(bt_force_u, bt_force_v, cs%BT_Domain, &
1322 if (add_uh0) pid_uhbt0 = pass_vector_start(uhbt0, vhbt0, cs%BT_Domain)
1323 pid_cor_ref = pass_vector_start(cor_ref_u, cor_ref_v, cs%BT_Domain)
1325 if (id_clock_pass_pre > 0)
call cpu_clock_end(id_clock_pass_pre)
1326 if (id_clock_calc_pre > 0)
call cpu_clock_begin(id_clock_calc_pre)
1328 do j=js-1,je+1 ;
do i=is-1,ie ; bt_rem_u(i,j) = g%mask2dCu(i,j) ;
enddo ;
enddo 1329 do j=js-1,je ;
do i=is-1,ie+1 ; bt_rem_v(i,j) = g%mask2dCv(i,j) ;
enddo ;
enddo 1330 if ((use_visc_rem) .and. (cs%drag_amp > 0.0))
then 1331 do j=js-1,je+1 ;
do i=is-1,ie ; av_rem_u(i,j) = 0.0 ;
enddo ;
enddo 1332 do j=js-1,je ;
do i=is-1,ie+1 ; av_rem_v(i,j) = 0.0 ;
enddo ;
enddo 1333 do k=1,nz ;
do j=js-1,je+1 ;
do i=is-1,ie
1334 av_rem_u(i,j) = av_rem_u(i,j) + cs%frhatu(i,j,k) * visc_rem_u(i,j,k)
1335 enddo ;
enddo ;
enddo 1336 do k=1,nz ;
do j=js-1,je ;
do i=is-1,ie+1
1337 av_rem_v(i,j) = av_rem_v(i,j) + cs%frhatv(i,j,k) * visc_rem_v(i,j,k)
1338 enddo ;
enddo ;
enddo 1339 if (cs%strong_drag)
then 1340 do j=js-1,je+1 ;
do i=is-1,ie
1341 bt_rem_u(i,j) = g%mask2dCu(i,j) * cs%drag_amp * &
1342 ((nstep * av_rem_u(i,j)) / (1.0 + (nstep-1)*av_rem_u(i,j)))
1344 do j=js-1,je ;
do i=is-1,ie+1
1345 bt_rem_v(i,j) = g%mask2dCv(i,j) * cs%drag_amp * &
1346 ((nstep * av_rem_v(i,j)) / (1.0 + (nstep-1)*av_rem_v(i,j)))
1349 do j=js-1,je+1 ;
do i=is-1,ie
1351 if (g%mask2dCu(i,j) * av_rem_u(i,j) > 0.0) &
1352 bt_rem_u(i,j) = g%mask2dCu(i,j) * cs%drag_amp * (av_rem_u(i,j)**instep)
1354 do j=js-1,je ;
do i=is-1,ie+1
1356 if (g%mask2dCv(i,j) * av_rem_v(i,j) > 0.0) &
1357 bt_rem_v(i,j) = g%mask2dCv(i,j) * cs%drag_amp * (av_rem_v(i,j)**instep)
1363 if (find_etaav)
then ;
do j=jsvf-1,jevf+1 ;
do i=isvf-1,ievf+1
1364 eta_sum(i,j) = 0.0 ; eta_wtd(i,j) = 0.0
1365 enddo ;
enddo ;
else ;
do j=jsvf-1,jevf+1 ;
do i=isvf-1,ievf+1
1367 enddo ;
enddo ;
endif 1368 do j=jsvf-1,jevf+1 ;
do i=isvf-1,ievf
1369 ubt_sum(i,j) = 0.0 ; uhbt_sum(i,j) = 0.0
1370 pfu_bt_sum(i,j) = 0.0 ; coru_bt_sum(i,j) = 0.0
1371 ubt_wtd(i,j) = 0.0 ; ubt_trans(i,j) = 0.0
1373 do j=jsvf-1,jevf ;
do i=isvf-1,ievf+1
1374 vbt_sum(i,j) = 0.0 ; vhbt_sum(i,j) = 0.0
1375 pfv_bt_sum(i,j) = 0.0 ; corv_bt_sum(i,j) = 0.0
1376 vbt_wtd(i,j) = 0.0 ; vbt_trans(i,j) = 0.0
1380 do j=jsvf-1,jevf+1;
do i=isvf-1,ievf+1 ; eta_src(i,j) = 0.0 ;
enddo ;
enddo 1381 if (cs%bound_BT_corr)
then ;
do j=js,je ;
do i=is,ie
1382 if (abs(cs%eta_cor(i,j)) > dt*cs%eta_cor_bound(i,j)) &
1383 cs%eta_cor(i,j) = sign(dt*cs%eta_cor_bound(i,j),cs%eta_cor(i,j))
1384 enddo ;
enddo ;
endif 1385 do j=js,je ;
do i=is,ie
1386 eta_src(i,j) = g%mask2dT(i,j) * (instep * cs%eta_cor(i,j) + dtbt * cs%eta_source(i,j))
1389 if (cs%dynamic_psurf)
then 1390 ice_is_rigid = (
associated(fluxes%rigidity_ice_u) .and. &
1391 associated(fluxes%rigidity_ice_v))
1392 h_min_dyn = gv%m_to_H * cs%Dmin_dyn_psurf
1393 if (ice_is_rigid .and. use_bt_cont) &
1394 call bt_cont_to_face_areas(bt_cont, datu, datv, g, ms, 0, .true.)
1395 if (ice_is_rigid)
then ;
do j=js,je ;
do i=is,ie
1401 idt_max2 = 0.5 * (dgeo_de * (1.0 + 2.0*bebt)) * (g%IareaT(i,j) * &
1402 ((gtot_e(i,j) * (datu(i,j)*g%IdxCu(i,j)) + &
1403 gtot_w(i,j) * (datu(i-1,j)*g%IdxCu(i-1,j))) + &
1404 (gtot_n(i,j) * (datv(i,j)*g%IdyCv(i,j)) + &
1405 gtot_s(i,j) * (datv(i,j-1)*g%IdyCv(i,j-1)))) + &
1406 ((g%CoriolisBu(i,j)**2 + g%CoriolisBu(i-1,j-1)**2) + &
1407 (g%CoriolisBu(i-1,j)**2 + g%CoriolisBu(i,j-1)**2)))
1408 h_eff_dx2 = max(h_min_dyn * (g%IdxT(i,j)**2 + g%IdyT(i,j)**2), &
1410 ((datu(i,j)*g%IdxCu(i,j) + datu(i-1,j)*g%IdxCu(i-1,j)) + &
1411 (datv(i,j)*g%IdyCv(i,j) + datv(i,j-1)*g%IdyCv(i,j-1)) ) )
1412 dyn_coef_max = cs%const_dyn_psurf * max(0.0, 1.0 - dtbt**2 * idt_max2) / &
1413 (dtbt**2 * h_eff_dx2)
1416 ice_strength = ((fluxes%rigidity_ice_u(i,j) + fluxes%rigidity_ice_u(i-1,j)) + &
1417 (fluxes%rigidity_ice_v(i,j) + fluxes%rigidity_ice_v(i,j-1))) / &
1418 (cs%ice_strength_length**2 * dtbt)
1421 dyn_coef_eta(i,j) = min(dyn_coef_max, ice_strength * gv%H_to_m)
1422 enddo ;
enddo ;
endif 1425 if (id_clock_calc_pre > 0)
call cpu_clock_end(id_clock_calc_pre)
1426 if (id_clock_pass_pre > 0)
call cpu_clock_begin(id_clock_pass_pre)
1427 if (nonblock_setup)
then 1428 if (cs%dynamic_psurf) pid_dyn_coef_eta = &
1429 pass_var_start(dyn_coef_eta, cs%BT_Domain, complete=.false.)
1430 if (interp_eta_pf)
then 1431 pid_eta_pf_1 = pass_var_start(eta_pf_1, cs%BT_Domain, complete=.false.)
1432 pid_d_eta_pf = pass_var_start(d_eta_pf, cs%BT_Domain, complete=.false.)
1434 if ((g%isd > cs%isdw) .or. (g%jsd > cs%jsdw)) &
1435 pid_eta_pf = pass_var_start(eta_pf, cs%BT_Domain, complete=.false.)
1437 pid_eta_src = pass_var_start(eta_src, cs%BT_Domain)
1441 pid_bt_rem_u = pass_vector_start(bt_rem_u, bt_rem_v, cs%BT_Domain, &
1444 if (interp_eta_pf)
then 1445 call pass_var(eta_pf_1, cs%BT_Domain, complete=.false.)
1446 call pass_var(d_eta_pf, cs%BT_Domain, complete=.false.)
1450 if ((g%isd > cs%isdw) .or. (g%jsd > cs%jsdw)) &
1451 call pass_var(eta_pf, cs%BT_Domain, complete=.false.)
1453 if (cs%dynamic_psurf)
call pass_var(dyn_coef_eta, cs%BT_Domain, complete=.false.)
1456 call pass_vector(bt_rem_u, bt_rem_v, cs%BT_Domain, to_all+scalar_pair)
1457 if (.not.use_bt_cont) &
1458 call pass_vector(datu, datv, cs%BT_Domain, to_all+scalar_pair)
1459 if (((g%isd > cs%isdw) .or. (g%jsd > cs%jsdw)) .or. (isq <= is-1) .or. (jsq <= js-1)) &
1460 call pass_vector(bt_force_u, bt_force_v, cs%BT_Domain, complete=.false.)
1461 if (g%nonblocking_updates)
then 1462 call pass_var(eta_src, cs%BT_Domain, complete=.true.)
1463 if (add_uh0)
call pass_vector(uhbt0, vhbt0, cs%BT_Domain, complete=.false.)
1464 call pass_vector(cor_ref_u, cor_ref_v, cs%BT_Domain, complete=.true.)
1466 call pass_var(eta_src, cs%BT_Domain, complete=.false.)
1467 if (add_uh0)
call pass_vector(uhbt0, vhbt0, cs%BT_Domain, complete=.false.)
1468 call pass_vector(cor_ref_u, cor_ref_v, cs%BT_Domain, complete=.false.)
1471 if (id_clock_pass_pre > 0)
call cpu_clock_end(id_clock_pass_pre)
1472 if (id_clock_calc_pre > 0)
call cpu_clock_begin(id_clock_calc_pre)
1475 if (nonblock_setup)
then 1476 if (id_clock_calc_pre > 0)
call cpu_clock_end(id_clock_calc_pre)
1477 if (id_clock_pass_pre > 0)
call cpu_clock_begin(id_clock_pass_pre)
1479 if (.not.use_bt_cont) &
1480 call pass_vector_complete(pid_datu, datu, datv, cs%BT_Domain, to_all+scalar_pair)
1482 if (((g%isd > cs%isdw) .or. (g%jsd > cs%jsdw)) .or. (isq <= is-1) .or. (jsq <= js-1)) &
1483 call pass_vector_complete(pid_bt_force_u, bt_force_u, bt_force_v, cs%BT_Domain)
1484 if (add_uh0)
call pass_vector_complete(pid_uhbt0, uhbt0, vhbt0, cs%BT_Domain)
1485 call pass_vector_complete(pid_cor_ref, cor_ref_u, cor_ref_v, cs%BT_Domain)
1487 if (cs%dynamic_psurf) &
1488 call pass_var_complete(pid_dyn_coef_eta, dyn_coef_eta, cs%BT_Domain)
1489 if (interp_eta_pf)
then 1490 call pass_var_complete(pid_eta_pf_1, eta_pf_1, cs%BT_Domain)
1491 call pass_var_complete(pid_d_eta_pf, d_eta_pf, cs%BT_Domain)
1493 if ((g%isd > cs%isdw) .or. (g%jsd > cs%jsdw)) &
1494 call pass_var_complete(pid_eta_pf, eta_pf, cs%BT_Domain)
1496 call pass_var_complete(pid_eta_src, eta_src, cs%BT_Domain)
1499 call pass_vector_complete(pid_bt_rem_u, bt_rem_u, bt_rem_v, cs%BT_Domain,&
1502 if (id_clock_pass_pre > 0)
call cpu_clock_end(id_clock_pass_pre)
1503 if (id_clock_calc_pre > 0)
call cpu_clock_begin(id_clock_calc_pre)
1507 call uvchksum(
"BT [uv]hbt", uhbt, vhbt, cs%debug_BT_HI, haloshift=0)
1508 call uvchksum(
"BT Initial [uv]bt", ubt, vbt, cs%debug_BT_HI, haloshift=0)
1509 call hchksum(eta,
"BT Initial eta",cs%debug_BT_HI, haloshift=0, scale=gv%H_to_m)
1510 call uvchksum(
"BT BT_force_[uv]", &
1511 bt_force_u, bt_force_v, cs%debug_BT_HI, haloshift=0)
1512 if (interp_eta_pf)
then 1513 call hchksum(eta_pf_1,
"BT eta_PF_1",cs%debug_BT_HI,haloshift=0, scale=gv%H_to_m)
1514 call hchksum(d_eta_pf,
"BT d_eta_PF",cs%debug_BT_HI,haloshift=0, scale=gv%H_to_m)
1516 call hchksum(eta_pf,
"BT eta_PF",cs%debug_BT_HI,haloshift=0, scale=gv%H_to_m)
1517 call hchksum(eta_pf_in,
"BT eta_PF_in",g%HI,haloshift=0, scale=gv%H_to_m)
1519 call uvchksum(
"BT Cor_ref_[uv]", cor_ref_u, cor_ref_v, &
1520 cs%debug_BT_HI, haloshift=0)
1521 call uvchksum(
"BT [uv]hbt0", uhbt0, vhbt0, cs%debug_BT_HI, haloshift=0, scale=gv%H_to_m)
1522 if (.not. use_bt_cont)
then 1523 call uvchksum(
"BT Dat[uv]", datu, datv, &
1524 cs%debug_BT_HI,haloshift=1, scale=gv%H_to_m)
1526 call uvchksum(
"BT wt_[uv]", wt_u, wt_v, g%HI,haloshift=1)
1527 call uvchksum(
"BT frhat[uv]", cs%frhatu, cs%frhatv, g%HI, haloshift=1)
1528 call uvchksum(
"BT bc_accel_[uv]", bc_accel_u, bc_accel_v, g%HI,haloshift=0)
1529 call uvchksum(
"BT IDat[uv]", cs%IDatu, cs%IDatv, g%HI,haloshift=0)
1530 if (use_visc_rem)
then 1531 call uvchksum(
"BT visc_rem_[uv]", visc_rem_u, visc_rem_v, g%HI, haloshift=1)
1535 if (query_averaging_enabled(cs%diag))
then 1536 if (cs%id_eta_st > 0)
call post_data(cs%id_eta_st, eta(isd:ied,jsd:jed), cs%diag)
1537 if (cs%id_ubt_st > 0)
call post_data(cs%id_ubt_st, ubt(isdb:iedb,jsd:jed), cs%diag)
1538 if (cs%id_vbt_st > 0)
call post_data(cs%id_vbt_st, vbt(isd:ied,jsdb:jedb), cs%diag)
1541 if (id_clock_calc_pre > 0)
call cpu_clock_end(id_clock_calc_pre)
1542 if (id_clock_calc > 0)
call cpu_clock_begin(id_clock_calc)
1544 if (project_velocity)
then ; eta_pf_bt => eta ;
else ; eta_pf_bt => eta_pred ;
endif 1546 if (cs%dt_bt_filter >= 0.0)
then 1547 dt_filt = 0.5 * max(0.0, min(cs%dt_bt_filter, 2.0*dt))
1549 dt_filt = 0.5 * max(0.0, dt * min(-cs%dt_bt_filter, 2.0))
1551 nfilter = ceiling(dt_filt / dtbt)
1554 sum_wt_vel = 0.0 ; sum_wt_eta = 0.0 ; sum_wt_accel = 0.0 ; sum_wt_trans = 0.0
1555 allocate(wt_vel(nstep+nfilter)) ;
allocate(wt_eta(nstep+nfilter))
1556 allocate(wt_trans(nstep+nfilter+1)) ;
allocate(wt_accel(nstep+nfilter+1))
1557 allocate(wt_accel2(nstep+nfilter+1))
1558 do n=1,nstep+nfilter
1560 if ( (n==nstep) .or. (dt_filt - abs(n-nstep)*dtbt >= 0.0))
then 1561 wt_vel(n) = 1.0 ; wt_eta(n) = 1.0
1562 elseif (dtbt + dt_filt - abs(n-nstep)*dtbt > 0.0)
then 1563 wt_vel(n) = 1.0 + (dt_filt / dtbt) - abs(n-nstep) ; wt_eta(n) = wt_vel(n)
1565 wt_vel(n) = 0.0 ; wt_eta(n) = 0.0
1571 sum_wt_vel = sum_wt_vel + wt_vel(n) ; sum_wt_eta = sum_wt_eta + wt_eta(n)
1573 wt_trans(nstep+nfilter+1) = 0.0 ; wt_accel(nstep+nfilter+1) = 0.0
1574 do n=nstep+nfilter,1,-1
1575 wt_trans(n) = wt_trans(n+1) + wt_eta(n)
1576 wt_accel(n) = wt_accel(n+1) + wt_vel(n)
1577 sum_wt_accel = sum_wt_accel + wt_accel(n) ; sum_wt_trans = sum_wt_trans + wt_trans(n)
1580 i_sum_wt_vel = 1.0 / sum_wt_vel ; i_sum_wt_accel = 1.0 / sum_wt_accel
1581 i_sum_wt_eta = 1.0 / sum_wt_eta ; i_sum_wt_trans = 1.0 / sum_wt_trans
1582 do n=1,nstep+nfilter
1583 wt_vel(n) = wt_vel(n) * i_sum_wt_vel
1584 wt_accel2(n) = wt_accel(n)
1585 wt_accel(n) = wt_accel(n) * i_sum_wt_accel
1586 wt_eta(n) = wt_eta(n) * i_sum_wt_eta
1591 sum_wt_vel = 0.0 ; sum_wt_eta = 0.0 ; sum_wt_accel = 0.0 ; sum_wt_trans = 0.0
1594 isv=is ; iev=ie ; jsv=js ; jev=je
1595 do n=1,nstep+nfilter
1597 sum_wt_vel = sum_wt_vel + wt_vel(n)
1598 sum_wt_eta = sum_wt_eta + wt_eta(n)
1599 sum_wt_accel = sum_wt_accel + wt_accel2(n)
1600 sum_wt_trans = sum_wt_trans + wt_trans(n)
1602 if (cs%clip_velocity)
then 1603 do j=jsv,jev ;
do i=isv-1,iev
1604 if (abs(ubt(i,j)) > cs%maxvel)
then 1606 ubt(i,j) = sign(cs%maxvel,ubt(i,j))
1609 do j=jsv-1,jev ;
do i=isv,iev
1610 if (abs(vbt(i,j)) > cs%maxvel)
then 1612 vbt(i,j) = sign(cs%maxvel,vbt(i,j))
1617 if ((iev - stencil < ie) .or. (jev - stencil < je))
then 1618 if (id_clock_calc > 0)
call cpu_clock_end(id_clock_calc)
1619 if (id_clock_pass_step > 0)
call cpu_clock_begin(id_clock_pass_step)
1620 if (g%nonblocking_updates)
then 1621 pid_ubt = pass_vector_start(ubt, vbt, cs%BT_Domain)
1622 pid_eta = pass_var_start(eta, cs%BT_Domain)
1623 call pass_vector_complete(pid_ubt, ubt, vbt, cs%BT_Domain)
1624 call pass_var_complete(pid_eta, eta, cs%BT_Domain)
1626 call pass_var(eta, cs%BT_Domain)
1627 call pass_vector(ubt, vbt, cs%BT_Domain)
1629 isv = isvf ; iev = ievf ; jsv = jsvf ; jev = jevf
1630 if (id_clock_pass_step > 0)
call cpu_clock_end(id_clock_pass_step)
1631 if (id_clock_calc > 0)
call cpu_clock_begin(id_clock_calc)
1633 isv = isv+stencil ; iev = iev-stencil
1634 jsv = jsv+stencil ; jev = jev-stencil
1637 if ((.not.use_bt_cont) .and. cs%Nonlinear_continuity .and. &
1638 (cs%Nonlin_cont_update_period > 0))
then 1639 if ((n>1) .and. (mod(n-1,cs%Nonlin_cont_update_period) == 0)) &
1640 call find_face_areas(datu, datv, g, gv, cs, ms, cs%rescale_D_bt, eta, 1+iev-ie)
1643 if (cs%dynamic_psurf .or. .not.project_velocity)
then 1644 if (use_bt_cont)
then 1646 do j=jsv-1,jev+1 ;
do i=isv-2,iev+1
1647 uhbt(i,j) = find_uhbt(ubt(i,j),btcl_u(i,j)) + uhbt0(i,j)
1649 do j=jsv-2,jev+1 ;
do i=isv-1,iev+1
1650 vhbt(i,j) = find_vhbt(vbt(i,j),btcl_v(i,j)) + vhbt0(i,j)
1652 do j=jsv-1,jev+1 ;
do i=isv-1,iev+1
1653 eta_pred(i,j) = (eta(i,j) + eta_src(i,j)) + (dtbt * cs%IareaT(i,j)) * &
1654 ((uhbt(i-1,j) - uhbt(i,j)) + (vhbt(i,j-1) - vhbt(i,j)))
1659 do j=jsv-1,jev+1 ;
do i=isv-1,iev+1
1660 eta_pred(i,j) = (eta(i,j) + eta_src(i,j)) + (dtbt * cs%IareaT(i,j)) * &
1661 (((datu(i-1,j)*ubt(i-1,j) + uhbt0(i-1,j)) - &
1662 (datu(i,j)*ubt(i,j) + uhbt0(i,j))) + &
1663 ((datv(i,j-1)*vbt(i,j-1) + vhbt0(i,j-1)) - &
1664 (datv(i,j)*vbt(i,j) + vhbt0(i,j))))
1668 if (cs%dynamic_psurf)
then ;
do j=jsv-1,jev+1 ;
do i=isv-1,iev+1
1669 p_surf_dyn(i,j) = dyn_coef_eta(i,j) * (eta_pred(i,j) - eta(i,j))
1670 enddo ;
enddo ;
endif 1676 if (find_etaav)
then ;
do j=js,je ;
do i=is,ie
1677 eta_sum(i,j) = eta_sum(i,j) + wt_accel2(n) * eta_pf_bt(i,j)
1678 enddo ;
enddo ;
endif 1680 if (interp_eta_pf)
then 1682 do j=jsv-1,jev+1 ;
do i=isv-1,iev+1
1683 eta_pf(i,j) = eta_pf_1(i,j) + wt_end*d_eta_pf(i,j)
1687 if (apply_obc_flather)
then 1689 do j=jsv,jev ;
do i=isv-2,iev+1
1690 ubt_old(i,j) = ubt(i,j)
1693 do j=jsv-2,jev+1 ;
do i=isv,iev
1694 vbt_old(i,j) = vbt(i,j)
1710 if (mod(n+g%first_direction,2)==1)
then 1713 do j=jsv-1,jev ;
do i=isv-1,iev+1
1714 cor = -1.0*((amer(i-1,j) * ubt(i-1,j) + cmer(i,j+1) * ubt(i,j+1)) + &
1715 (bmer(i,j) * ubt(i,j) + dmer(i-1,j+1) * ubt(i-1,j+1))) - cor_ref_v(i,j)
1716 gradp = ((eta_pf_bt(i,j)-eta_pf(i,j))*gtot_n(i,j) - &
1717 (eta_pf_bt(i,j+1)-eta_pf(i,j+1))*gtot_s(i,j+1)) * &
1718 dgeo_de * cs%IdyCv(i,j)
1719 if (cs%dynamic_psurf) &
1720 gradp = gradp + (p_surf_dyn(i,j) - p_surf_dyn(i,j+1)) * cs%IdyCv(i,j)
1721 v_accel_bt(i,j) = v_accel_bt(i,j) + wt_accel(n) * (cor + gradp)
1722 if (find_pf) pfv_bt_sum(i,j) = pfv_bt_sum(i,j) + wt_accel2(n) * gradp
1723 if (find_cor) corv_bt_sum(i,j) = corv_bt_sum(i,j) + wt_accel2(n) * cor
1725 if (apply_v_obcs)
then ;
if (obc%segnum_v(i,j) /= obc_none) cycle ;
endif 1728 vbt(i,j) = bt_rem_v(i,j) * (vbt(i,j) + &
1729 dtbt * ((bt_force_v(i,j) + cor) + gradp))
1730 if (project_velocity)
then 1731 vel_trans = (1.0 + be_proj)*vbt(i,j) - be_proj*vel_prev
1733 vel_trans = (1.0-bebt)*vel_prev + bebt*vbt(i,j)
1735 if (use_bt_cont)
then 1736 vhbt(i,j) = find_vhbt(vel_trans,btcl_v(i,j)) + vhbt0(i,j)
1738 vhbt(i,j) = datv(i,j)*vel_trans + vhbt0(i,j)
1741 vbt_sum(i,j) = vbt_sum(i,j) + wt_trans(n) * vel_trans
1742 vhbt_sum(i,j) = vhbt_sum(i,j) + wt_trans(n) * vhbt(i,j)
1743 vbt_wtd(i,j) = vbt_wtd(i,j) + wt_vel(n) * vbt(i,j)
1747 do j=jsv,jev ;
do i=isv-1,iev
1748 cor = ((azon(i,j) * vbt(i+1,j) + czon(i,j) * vbt(i,j-1)) + &
1749 (bzon(i,j) * vbt(i,j) + dzon(i,j) * vbt(i+1,j-1))) - cor_ref_u(i,j)
1750 gradp = ((eta_pf_bt(i,j)-eta_pf(i,j))*gtot_e(i,j) - &
1751 (eta_pf_bt(i+1,j)-eta_pf(i+1,j))*gtot_w(i+1,j)) * &
1752 dgeo_de * cs%IdxCu(i,j)
1753 if (cs%dynamic_psurf) &
1754 gradp = gradp + (p_surf_dyn(i,j) - p_surf_dyn(i+1,j)) * cs%IdxCu(i,j)
1755 u_accel_bt(i,j) = u_accel_bt(i,j) + wt_accel(n) * (cor + gradp)
1756 if (find_pf) pfu_bt_sum(i,j) = pfu_bt_sum(i,j) + wt_accel2(n) * gradp
1757 if (find_cor) coru_bt_sum(i,j) = coru_bt_sum(i,j) + wt_accel2(n) * cor
1759 if (apply_u_obcs)
then ;
if (obc%segnum_u(i,j) /= obc_none) cycle ;
endif 1762 ubt(i,j) = bt_rem_u(i,j) * (ubt(i,j) + &
1763 dtbt * ((bt_force_u(i,j) + cor) + gradp))
1764 if (project_velocity)
then 1765 vel_trans = (1.0 + be_proj)*ubt(i,j) - be_proj*vel_prev
1767 vel_trans = (1.0-bebt)*vel_prev + bebt*ubt(i,j)
1769 if (use_bt_cont)
then 1770 uhbt(i,j) = find_uhbt(vel_trans, btcl_u(i,j)) + uhbt0(i,j)
1772 uhbt(i,j) = datu(i,j)*vel_trans + uhbt0(i,j)
1775 ubt_sum(i,j) = ubt_sum(i,j) + wt_trans(n) * vel_trans
1776 uhbt_sum(i,j) = uhbt_sum(i,j) + wt_trans(n) * uhbt(i,j)
1777 ubt_wtd(i,j) = ubt_wtd(i,j) + wt_vel(n) * ubt(i,j)
1783 do j=jsv-1,jev+1 ;
do i=isv-1,iev
1784 cor = ((azon(i,j) * vbt(i+1,j) + czon(i,j) * vbt(i,j-1)) + &
1785 (bzon(i,j) * vbt(i,j) + dzon(i,j) * vbt(i+1,j-1))) - cor_ref_u(i,j)
1786 gradp = ((eta_pf_bt(i,j)-eta_pf(i,j))*gtot_e(i,j) - &
1787 (eta_pf_bt(i+1,j)-eta_pf(i+1,j))*gtot_w(i+1,j)) * &
1788 dgeo_de * cs%IdxCu(i,j)
1789 if (cs%dynamic_psurf) &
1790 gradp = gradp + (p_surf_dyn(i,j) - p_surf_dyn(i+1,j)) * cs%IdxCu(i,j)
1791 u_accel_bt(i,j) = u_accel_bt(i,j) + wt_accel(n) * (cor + gradp)
1792 if (find_pf) pfu_bt_sum(i,j) = pfu_bt_sum(i,j) + wt_accel2(n) * gradp
1793 if (find_cor) coru_bt_sum(i,j) = coru_bt_sum(i,j) + wt_accel2(n) * cor
1795 if (apply_u_obcs)
then ;
if (obc%segnum_u(i,j) /= obc_none) cycle ;
endif 1798 ubt(i,j) = bt_rem_u(i,j) * (ubt(i,j) + &
1799 dtbt * ((bt_force_u(i,j) + cor) + gradp))
1800 if (project_velocity)
then 1801 vel_trans = (1.0 + be_proj)*ubt(i,j) - be_proj*vel_prev
1803 vel_trans = (1.0-bebt)*vel_prev + bebt*ubt(i,j)
1805 if (use_bt_cont)
then 1806 uhbt(i,j) = find_uhbt(vel_trans,btcl_u(i,j)) + uhbt0(i,j)
1808 uhbt(i,j) = datu(i,j)*vel_trans + uhbt0(i,j)
1811 ubt_sum(i,j) = ubt_sum(i,j) + wt_trans(n) * vel_trans
1812 uhbt_sum(i,j) = uhbt_sum(i,j) + wt_trans(n) * uhbt(i,j)
1813 ubt_wtd(i,j) = ubt_wtd(i,j) + wt_vel(n) * ubt(i,j)
1817 do j=jsv-1,jev ;
do i=isv,iev
1818 cor = -1.0*((amer(i-1,j) * ubt(i-1,j) + bmer(i,j) * ubt(i,j)) + &
1819 (cmer(i,j+1) * ubt(i,j+1) + dmer(i-1,j+1) * ubt(i-1,j+1))) - cor_ref_v(i,j)
1820 gradp = ((eta_pf_bt(i,j)-eta_pf(i,j))*gtot_n(i,j) - &
1821 (eta_pf_bt(i,j+1)-eta_pf(i,j+1))*gtot_s(i,j+1)) * &
1822 dgeo_de * cs%IdyCv(i,j)
1823 if (cs%dynamic_psurf) &
1824 gradp = gradp + (p_surf_dyn(i,j) - p_surf_dyn(i,j+1)) * cs%IdyCv(i,j)
1825 v_accel_bt(i,j) = v_accel_bt(i,j) + wt_accel(n) * (cor + gradp)
1826 if (find_pf) pfv_bt_sum(i,j) = pfv_bt_sum(i,j) + wt_accel2(n) * gradp
1827 if (find_cor) corv_bt_sum(i,j) = corv_bt_sum(i,j) + wt_accel2(n) * cor
1829 if (apply_v_obcs)
then ;
if (obc%segnum_v(i,j) /= obc_none) cycle ;
endif 1832 vbt(i,j) = bt_rem_v(i,j) * (vbt(i,j) + &
1833 dtbt * ((bt_force_v(i,j) + cor) + gradp))
1834 if (project_velocity)
then 1835 vel_trans = (1.0 + be_proj)*vbt(i,j) - be_proj*vel_prev
1837 vel_trans = (1.0-bebt)*vel_prev + bebt*vbt(i,j)
1839 if (use_bt_cont)
then 1840 vhbt(i,j) = find_vhbt(vel_trans,btcl_v(i,j)) + vhbt0(i,j)
1842 vhbt(i,j) = datv(i,j)*vel_trans + vhbt0(i,j)
1845 vbt_sum(i,j) = vbt_sum(i,j) + wt_trans(n) * vel_trans
1846 vhbt_sum(i,j) = vhbt_sum(i,j) + wt_trans(n) * vhbt(i,j)
1847 vbt_wtd(i,j) = vbt_wtd(i,j) + wt_vel(n) * vbt(i,j)
1852 if (apply_obcs)
then 1853 call apply_velocity_obcs(obc, ubt, vbt, uhbt, vhbt, &
1854 ubt_trans, vbt_trans, eta, ubt_old, vbt_old, bt_obc, &
1855 g, ms, iev-ie, dtbt, bebt, use_bt_cont, datu, datv, btcl_u, btcl_v, &
1857 if (apply_u_obcs)
then ;
do j=js,je ;
do i=is-1,ie
1858 if (obc%segnum_u(i,j) /= obc_none)
then 1859 ubt_sum(i,j) = ubt_sum(i,j) + wt_trans(n) * ubt_trans(i,j)
1860 uhbt_sum(i,j) = uhbt_sum(i,j) + wt_trans(n) * uhbt(i,j)
1861 ubt_wtd(i,j) = ubt_wtd(i,j) + wt_vel(n) * ubt(i,j)
1863 enddo ;
enddo ;
endif 1864 if (apply_v_obcs)
then ;
do j=js-1,je ;
do i=is,ie
1865 if (obc%segnum_v(i,j) /= obc_none)
then 1866 vbt_sum(i,j) = vbt_sum(i,j) + wt_trans(n) * vbt_trans(i,j)
1867 vhbt_sum(i,j) = vhbt_sum(i,j) + wt_trans(n) * vhbt(i,j)
1868 vbt_wtd(i,j) = vbt_wtd(i,j) + wt_vel(n) * vbt(i,j)
1870 enddo ;
enddo ;
endif 1873 if (cs%debug_bt)
then 1874 call uvchksum(
"BT [uv]hbt just after OBC", &
1875 uhbt, vhbt, cs%debug_BT_HI, haloshift=iev-ie)
1879 do j=jsv,jev ;
do i=isv,iev
1880 eta(i,j) = (eta(i,j) + eta_src(i,j)) + (dtbt * cs%IareaT(i,j)) * &
1881 ((uhbt(i-1,j) - uhbt(i,j)) + (vhbt(i,j-1) - vhbt(i,j)))
1882 eta_wtd(i,j) = eta_wtd(i,j) + eta(i,j) * wt_eta(n)
1884 if (apply_obcs)
call apply_eta_obcs(obc, eta, ubt_old, vbt_old, bt_obc, &
1885 g, ms, iev-ie, dtbt)
1887 if (do_hifreq_output)
then 1888 time_step_end = time_bt_start + set_time(int(floor(n*dtbt+0.5)))
1889 call enable_averaging(dtbt, time_step_end, cs%diag)
1890 if (cs%id_ubt_hifreq > 0)
call post_data(cs%id_ubt_hifreq, ubt(isdb:iedb,jsd:jed), cs%diag)
1891 if (cs%id_vbt_hifreq > 0)
call post_data(cs%id_vbt_hifreq, vbt(isd:ied,jsdb:jedb), cs%diag)
1892 if (cs%id_eta_hifreq > 0)
call post_data(cs%id_eta_hifreq, eta(isd:ied,jsd:jed), cs%diag)
1893 if (cs%id_uhbt_hifreq > 0)
call post_data(cs%id_uhbt_hifreq, uhbt(isdb:iedb,jsd:jed), cs%diag)
1894 if (cs%id_vhbt_hifreq > 0)
call post_data(cs%id_vhbt_hifreq, vhbt(isd:ied,jsdb:jedb), cs%diag)
1895 if (cs%id_eta_pred_hifreq > 0)
call post_data(cs%id_eta_pred_hifreq, eta_pf_bt(isd:ied,jsd:jed), cs%diag)
1898 if (cs%debug_bt)
then 1899 write(mesg,
'("BT step ",I4)') n
1900 call uvchksum(trim(mesg)//
" [uv]bt", &
1901 ubt, vbt, cs%debug_BT_HI, haloshift=iev-ie)
1902 call hchksum(eta, trim(mesg)//
" eta", cs%debug_BT_HI, haloshift=iev-ie, scale=gv%H_to_m)
1906 if (id_clock_calc > 0)
call cpu_clock_end(id_clock_calc)
1907 if (id_clock_calc_post > 0)
call cpu_clock_begin(id_clock_calc_post)
1910 if (do_hifreq_output)
call enable_averaging(time_int_in, time_end_in, cs%diag)
1912 i_sum_wt_vel = 1.0 / sum_wt_vel ; i_sum_wt_eta = 1.0 / sum_wt_eta
1913 i_sum_wt_accel = 1.0 / sum_wt_accel ; i_sum_wt_trans = 1.0 / sum_wt_trans
1915 if (find_etaav)
then ;
do j=js,je ;
do i=is,ie
1916 etaav(i,j) = eta_sum(i,j) * i_sum_wt_accel
1917 enddo ;
enddo ;
endif 1918 do j=js-1,je+1 ;
do i=is-1,ie+1 ; e_anom(i,j) = 0.0 ;
enddo ;
enddo 1919 if (interp_eta_pf)
then 1920 do j=js,je ;
do i=is,ie
1921 e_anom(i,j) = dgeo_de * (0.5 * (eta(i,j) + eta_in(i,j)) - &
1922 (eta_pf_1(i,j) + 0.5*d_eta_pf(i,j)))
1925 do j=js,je ;
do i=is,ie
1926 e_anom(i,j) = dgeo_de * (0.5 * (eta(i,j) + eta_in(i,j)) - eta_pf(i,j))
1931 do j=js,je ;
do i=is,ie
1932 eta_out(i,j) = eta_wtd(i,j) * i_sum_wt_eta
1935 if (id_clock_calc_post > 0)
call cpu_clock_end(id_clock_calc_post)
1936 if (id_clock_pass_post > 0)
call cpu_clock_begin(id_clock_pass_post)
1937 if (g%nonblocking_updates)
then 1938 pid_e_anom = pass_var_start(e_anom, g%Domain)
1940 if (find_etaav)
call pass_var(etaav, g%Domain, complete=.false.)
1941 call pass_var(e_anom, g%Domain)
1943 if (id_clock_pass_post > 0)
call cpu_clock_end(id_clock_pass_post)
1944 if (id_clock_calc_post > 0)
call cpu_clock_begin(id_clock_calc_post)
1946 do j=js,je ;
do i=is-1,ie
1947 cs%ubtav(i,j) = ubt_sum(i,j) * i_sum_wt_trans
1948 uhbtav(i,j) = uhbt_sum(i,j) * i_sum_wt_trans
1950 ubt_wtd(i,j) = ubt_wtd(i,j) * i_sum_wt_vel
1953 do j=js-1,je ;
do i=is,ie
1954 cs%vbtav(i,j) = vbt_sum(i,j) * i_sum_wt_trans
1955 vhbtav(i,j) = vhbt_sum(i,j) * i_sum_wt_trans
1957 vbt_wtd(i,j) = vbt_wtd(i,j) * i_sum_wt_vel
1960 if (
present(uhbt_out))
then 1962 if (use_bt_cont)
then ;
do j=js,je ;
do i=is-1,ie
1963 uhbt_out(i,j) = find_uhbt(ubt_wtd(i,j), btcl_u(i,j)) + uhbt0(i,j)
1964 enddo ;
enddo ;
else ;
do j=js,je ;
do i=is-1,ie
1965 uhbt_out(i,j) = ubt_wtd(i,j) * datu(i,j) + uhbt0(i,j)
1966 enddo ;
enddo ;
endif 1969 if (use_bt_cont)
then ;
do j=js-1,je ;
do i=is,ie
1970 vhbt_out(i,j) = find_vhbt(vbt_wtd(i,j), btcl_v(i,j)) + vhbt0(i,j)
1971 enddo ;
enddo ;
else ;
do j=js-1,je ;
do i=is,ie
1972 vhbt_out(i,j) = vbt_wtd(i,j) * datv(i,j) + vhbt0(i,j)
1973 enddo ;
enddo ;
endif 1976 if (id_clock_calc_post > 0)
call cpu_clock_end(id_clock_calc_post)
1977 if (id_clock_pass_post > 0)
call cpu_clock_begin(id_clock_pass_post)
1978 if (g%nonblocking_updates)
then 1979 call pass_var_complete(pid_e_anom, e_anom, g%Domain)
1981 if (find_etaav) pid_etaav = pass_var_start(etaav, g%Domain)
1982 pid_uhbtav = pass_vector_start(uhbtav, vhbtav, g%Domain, complete=.false.)
1983 pid_ubtav = pass_vector_start(cs%ubtav, cs%vbtav, g%Domain)
1985 call pass_vector(cs%ubtav, cs%vbtav, g%Domain, complete=.false.)
1986 call pass_vector(uhbtav, vhbtav, g%Domain)
1988 if (id_clock_pass_post > 0)
call cpu_clock_end(id_clock_pass_post)
1989 if (id_clock_calc_post > 0)
call cpu_clock_begin(id_clock_calc_post)
1996 do j=js,je ;
do i=is-1,ie
1997 accel_layer_u(i,j,k) = u_accel_bt(i,j) - &
1998 ((pbce(i+1,j,k) - gtot_w(i+1,j)) * e_anom(i+1,j) - &
1999 (pbce(i,j,k) - gtot_e(i,j)) * e_anom(i,j)) * cs%IdxCu(i,j)
2001 do j=js-1,je ;
do i=is,ie
2002 accel_layer_v(i,j,k) = v_accel_bt(i,j) - &
2003 ((pbce(i,j+1,k) - gtot_s(i,j+1))*e_anom(i,j+1) - &
2004 (pbce(i,j,k) - gtot_n(i,j))*e_anom(i,j)) * cs%IdyCv(i,j)
2008 if (id_clock_calc_post > 0)
call cpu_clock_end(id_clock_calc_post)
2011 if (query_averaging_enabled(cs%diag))
then 2013 do j=js,je ;
do i=is-1,ie ; cs%ubt_IC(i,j) = ubt_wtd(i,j) ;
enddo ;
enddo 2014 do j=js-1,je ;
do i=is,ie ; cs%vbt_IC(i,j) = vbt_wtd(i,j) ;
enddo ;
enddo 2015 if (
present(uhbt_out))
then 2016 do j=js,je ;
do i=is-1,ie ; cs%uhbt_IC(i,j) = uhbt_out(i,j) ;
enddo ;
enddo 2017 do j=js-1,je ;
do i=is,ie ; cs%vhbt_IC(i,j) = vhbt_out(i,j) ;
enddo ;
enddo 2018 elseif (use_bt_cont)
then 2019 do j=js,je ;
do i=is-1,ie
2020 cs%uhbt_IC(i,j) = find_uhbt(ubt_wtd(i,j), btcl_u(i,j)) + uhbt0(i,j)
2022 do j=js-1,je ;
do i=is,ie
2023 cs%vhbt_IC(i,j) = find_vhbt(vbt_wtd(i,j), btcl_v(i,j)) + vhbt0(i,j)
2026 do j=js,je ;
do i=is-1,ie
2027 cs%uhbt_IC(i,j) = ubt_wtd(i,j) * datu(i,j) + uhbt0(i,j)
2029 do j=js-1,je ;
do i=is,ie
2030 cs%vhbt_IC(i,j) = vbt_wtd(i,j) * datv(i,j) + vhbt0(i,j)
2035 if (cs%id_PFu_bt > 0)
then 2036 do j=js,je ;
do i=is-1,ie
2037 pfu_bt_sum(i,j) = pfu_bt_sum(i,j) * i_sum_wt_accel
2039 call post_data(cs%id_PFu_bt, pfu_bt_sum(isdb:iedb,jsd:jed), cs%diag)
2041 if (cs%id_PFv_bt > 0)
then 2042 do j=js-1,je ;
do i=is,ie
2043 pfv_bt_sum(i,j) = pfv_bt_sum(i,j) * i_sum_wt_accel
2045 call post_data(cs%id_PFv_bt, pfv_bt_sum(isd:ied,jsdb:jedb), cs%diag)
2047 if (cs%id_Coru_bt > 0)
then 2048 do j=js,je ;
do i=is-1,ie
2049 coru_bt_sum(i,j) = coru_bt_sum(i,j) * i_sum_wt_accel
2051 call post_data(cs%id_Coru_bt, coru_bt_sum(isdb:iedb,jsd:jed), cs%diag)
2053 if (cs%id_Corv_bt > 0)
then 2054 do j=js-1,je ;
do i=is,ie
2055 corv_bt_sum(i,j) = corv_bt_sum(i,j) * i_sum_wt_accel
2057 call post_data(cs%id_Corv_bt, corv_bt_sum(isd:ied,jsdb:jedb), cs%diag)
2059 if (cs%id_ubtforce > 0)
call post_data(cs%id_ubtforce, bt_force_u(isdb:iedb,jsd:jed), cs%diag)
2060 if (cs%id_vbtforce > 0)
call post_data(cs%id_vbtforce, bt_force_v(isd:ied,jsdb:jedb), cs%diag)
2061 if (cs%id_uaccel > 0)
call post_data(cs%id_uaccel, u_accel_bt(isdb:iedb,jsd:jed), cs%diag)
2062 if (cs%id_vaccel > 0)
call post_data(cs%id_vaccel, v_accel_bt(isd:ied,jsdb:jedb), cs%diag)
2063 if (cs%id_Datu_res > 0)
call post_data(cs%id_Datu_res, cs%Datu_res(isdb:iedb,jsd:jed), cs%diag)
2064 if (cs%id_Datv_res > 0)
call post_data(cs%id_Datv_res, cs%Datv_res(isd:ied,jsdb:jedb), cs%diag)
2066 if (cs%id_eta_cor > 0)
call post_data(cs%id_eta_cor, cs%eta_cor, cs%diag)
2067 if (cs%id_eta_bt > 0)
call post_data(cs%id_eta_bt, eta_out, cs%diag)
2068 if (cs%id_gtotn > 0)
call post_data(cs%id_gtotn, gtot_n(isd:ied,jsd:jed), cs%diag)
2069 if (cs%id_gtots > 0)
call post_data(cs%id_gtots, gtot_s(isd:ied,jsd:jed), cs%diag)
2070 if (cs%id_gtote > 0)
call post_data(cs%id_gtote, gtot_e(isd:ied,jsd:jed), cs%diag)
2071 if (cs%id_gtotw > 0)
call post_data(cs%id_gtotw, gtot_w(isd:ied,jsd:jed), cs%diag)
2072 if (cs%id_ubt > 0)
call post_data(cs%id_ubt, ubt_wtd(isdb:iedb,jsd:jed), cs%diag)
2073 if (cs%id_vbt > 0)
call post_data(cs%id_vbt, vbt_wtd(isd:ied,jsdb:jedb), cs%diag)
2074 if (cs%id_ubtav > 0)
call post_data(cs%id_ubtav, cs%ubtav, cs%diag)
2075 if (cs%id_vbtav > 0)
call post_data(cs%id_vbtav, cs%vbtav, cs%diag)
2076 if (use_visc_rem)
then 2077 if (cs%id_visc_rem_u > 0)
call post_data(cs%id_visc_rem_u, visc_rem_u, cs%diag)
2078 if (cs%id_visc_rem_v > 0)
call post_data(cs%id_visc_rem_v, visc_rem_v, cs%diag)
2081 if (cs%id_frhatu > 0)
call post_data(cs%id_frhatu, cs%frhatu, cs%diag)
2082 if (cs%id_uhbt > 0)
call post_data(cs%id_uhbt, uhbtav, cs%diag)
2083 if (cs%id_frhatv > 0)
call post_data(cs%id_frhatv, cs%frhatv, cs%diag)
2084 if (cs%id_vhbt > 0)
call post_data(cs%id_vhbt, vhbtav, cs%diag)
2086 if (cs%id_frhatu1 > 0)
call post_data(cs%id_frhatu1, cs%frhatu1, cs%diag)
2087 if (cs%id_frhatv1 > 0)
call post_data(cs%id_frhatv1, cs%frhatv1, cs%diag)
2089 if (cs%id_frhatu1 > 0) cs%frhatu1(:,:,:) = cs%frhatu(:,:,:)
2090 if (cs%id_frhatv1 > 0) cs%frhatv1(:,:,:) = cs%frhatv(:,:,:)
2093 if (g%nonblocking_updates)
then 2094 if (find_etaav)
call pass_var_complete(pid_etaav, etaav, g%Domain)
2095 call pass_vector_complete(pid_uhbtav, uhbtav, vhbtav, g%Domain)
2096 call pass_vector_complete(pid_ubtav, cs%ubtav, cs%vbtav, g%Domain)
2099 if (apply_obcs)
call destroy_bt_obc(bt_obc)
Ocean grid type. See mom_grid for details.