16 implicit none ;
private
18 #include <MOM_memory.h>
21 public bfb_initialize_sponges_southonly
30 logical :: first_call = .true.
38 subroutine bfb_set_coord(Rlay, g_prime, GV, param_file, eqn_of_state)
39 real,
dimension(NKMEM_),
intent(out) :: rlay
40 real,
dimension(NKMEM_),
intent(out) :: g_prime
44 type(
eos_type),
pointer :: eqn_of_state
47 real :: drho_dt, sst_s, t_bot, rho_top, rho_bot
49 character(len=40) :: mdl =
"BFB_set_coord"
51 call get_param(param_file, mdl,
"DRHO_DT", drho_dt, &
52 "Rate of change of density with temperature.", &
53 units=
"kg m-3 K-1", default=-0.2)
54 call get_param(param_file, mdl,
"SST_S", sst_s, &
55 "SST at the suothern edge of the domain.", units=
"C", default=20.0)
56 call get_param(param_file, mdl,
"T_BOT", t_bot, &
57 "Bottom Temp", units=
"C", default=5.0)
58 rho_top = gv%rho0 + drho_dt*sst_s
59 rho_bot = gv%rho0 + drho_dt*t_bot
63 rlay(k) = (rho_bot - rho_top)/(nz-1)*real(k-1) + rho_top
65 g_prime(k) = (rlay(k) - rlay(k-1)) * gv%g_Earth/gv%rho0
67 g_prime(k) = gv%g_Earth
73 if (first_call)
call write_bfb_log(param_file)
75 end subroutine bfb_set_coord
79 subroutine bfb_initialize_sponges_southonly(G, GV, US, use_temperature, tv, param_file, CSp, h)
83 logical,
intent(in) :: use_temperature
88 real,
dimension(NIMEM_, NJMEM_, NKMEM_), &
92 real :: eta(szi_(g),szj_(g),szk_(g)+1)
93 real :: idamp(szi_(g),szj_(g))
96 real :: damp, e_dense, damp_new, slat, wlon, lenlat, lenlon, nlat
97 character(len=40) :: mdl =
"BFB_initialize_sponges_southonly"
98 integer :: i, j, k, is, ie, js, je, isd, ied, jsd, jed, nz
100 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
101 isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
103 eta(:,:,:) = 0.0 ; idamp(:,:) = 0.0
111 call get_param(param_file, mdl,
"MINIMUM_DEPTH", min_depth, &
112 "The minimum depth of the ocean.", units=
"m", default=0.0, scale=us%m_to_Z)
114 call get_param(param_file, mdl,
"SOUTHLAT", slat, &
115 "The southern latitude of the domain.", units=
"degrees")
116 call get_param(param_file, mdl,
"LENLAT", lenlat, &
117 "The latitudinal length of the domain.", units=
"degrees")
118 call get_param(param_file, mdl,
"WESTLON", wlon, &
119 "The western longitude of the domain.", units=
"degrees", default=0.0)
120 call get_param(param_file, mdl,
"LENLON", lenlon, &
121 "The longitudinal length of the domain.", units=
"degrees")
123 do k=1,nz ; h0(k) = -g%max_depth * real(k-1) / real(nz) ;
enddo
128 do i=is,ie;
do j=js,je
129 if (g%geoLatT(i,j) < slat+2.0)
then ; damp = 1.0
130 elseif (g%geoLatT(i,j) < slat+4.0)
then
131 damp_new = 1.0*(slat+4.0-g%geoLatT(i,j))/2.0
139 do k = 1,nz; eta(i,j,k) = h0(k);
enddo
153 eta(i,j,nz+1) = -g%max_depth
155 if (g%bathyT(i,j) > min_depth)
then
156 idamp(i,j) = damp/86400.0
157 else ; idamp(i,j) = 0.0 ;
endif
162 call initialize_sponge(idamp, eta, g, param_file, csp, gv)
168 if (first_call)
call write_bfb_log(param_file)
170 end subroutine bfb_initialize_sponges_southonly
173 subroutine write_bfb_log(param_file)
179 #include "version_variable.h"
180 character(len=40) :: mdl =
"BFB_initialization"
185 end subroutine write_bfb_log