18 use random_numbers_mod,
only: initializerandomnumberstream, getrandomnumbers, randomnumberstream
20 implicit none ;
private
22 #include <MOM_memory.h>
37 real,
dimension(G%isd:G%ied,G%jsd:G%jed), &
40 real,
intent(in) :: max_depth
48 # include "version_variable.h"
49 character(len=40) :: mdl =
"Neverland_initialize_topography"
50 integer :: i, j, is, ie, js, je, isd, ied, jsd, jed
51 real :: nl_roughness_amp, nl_top_amp
52 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
53 isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
55 call mom_mesg(
" Neverland_initialization.F90, Neverland_initialize_topography: setting topography", 5)
58 call get_param(param_file, mdl,
"NL_ROUGHNESS_AMP", nl_roughness_amp, &
59 "Amplitude of wavy signal in bathymetry.", default=0.05)
60 call get_param(param_file, mdl,
"NL_CONTINENT_AMP", nl_top_amp, &
61 "Scale factor for topography - 0.0 for no continents.", default=1.0)
66 do j=js,je ;
do i=is,ie
67 x = (g%geoLonT(i,j)-g%west_lon) / g%len_lon
68 y =( g%geoLatT(i,j)-g%south_lat) / g%len_lat
70 d(i,j) = 1.0 - 1.1 *
spike(y-1,0.12) - 1.1 *
spike(y,0.12) - &
72 (1.2 *
spike(x,0.2) + 1.2 *
spike(x-1.0,0.2)) *
spike(min(0.0,y-0.3),0.2) &
73 + 1.2 *
spike(x-0.5,0.2) *
spike(min(0.0,y-0.55),0.2) &
74 + 1.2 * (
spike(x,0.12) +
spike(x-1,0.12)) *
spike(max(0.0,y-0.06),0.12) &
79 - nl_roughness_amp * cos(14*pi*x) * sin(14*pi*y) &
80 - nl_roughness_amp * cos(20*pi*x) * cos(20*pi*y)
81 if (d(i,j) < 0.0) d(i,j) = 0.0
82 d(i,j) = d(i,j) * max_depth
90 real ,
intent(in) :: x
91 real ,
intent(in) :: l
95 cosbell = 0.5 * (1 + cos(pi*min(abs(x/l),1.0)))
99 real function spike(x, L)
101 real ,
intent(in) :: x
102 real ,
intent(in) :: l
106 spike = (1 - sin(pi*min(abs(x/l),0.5)))
117 real,
intent(out),
dimension(SZI_(G),SZJ_(G),SZK_(GV)) :: h
122 type(
eos_type),
pointer :: eqn_of_state
124 real,
intent(in) :: p_ref
127 real :: e0(szk_(g)+1)
129 real,
dimension(SZK_(G)) :: h_profile
135 type(randomnumberstream) :: rns
136 character(len=40) :: mdl =
"Neverland_initialize_thickness"
137 integer :: i, j, k, k1, is, ie, js, je, nz, itt
139 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
141 call mom_mesg(
" Neverland_initialization.F90, Neverland_initialize_thickness: setting thickness", 5)
142 call get_param(param_file, mdl,
"INIT_THICKNESS_PROFILE", h_profile, &
143 "Profile of initial layer thicknesses.", units=
"m", scale=us%m_to_Z, &
144 fail_if_missing=.true.)
145 call get_param(param_file, mdl,
"NL_THICKNESS_PERT_AMP", pert_amp, &
146 "Amplitude of finite scale perturbations as fraction of depth.", &
147 units=
"nondim", default=0.)
148 call get_param(param_file, mdl,
"NL_THICKNESS_NOISE_AMP", h_noise, &
149 "Amplitude of noise to scale layer by.", units=
"nondim", default=0.)
154 e0(k+1) = e0(k) - h_profile(k)
157 do j=js,je ;
do i=is,ie
158 e_interface = -g%bathyT(i,j)
160 h(i,j,k) = gv%Z_to_H * (e0(k) - e_interface)
161 x=(g%geoLonT(i,j)-g%west_lon)/g%len_lon
162 y=(g%geoLatT(i,j)-g%south_lat)/g%len_lat
163 r1=sqrt((x-0.7)**2+(y-0.2)**2)
164 r2=sqrt((x-0.3)**2+(y-0.25)**2)
165 h(i,j,k) = h(i,j,k) + pert_amp * (e0(k) - e0(nz+1)) * gv%Z_to_H * &
167 if (h_noise /= 0.)
then
168 rns = initializerandomnumberstream( int( 4096*(x + (y+1.)) ) )
169 call getrandomnumbers(rns, noise)
170 noise = h_noise * 2. * ( noise - 0.5 )
171 h(i,j,k) = ( 1. + noise ) * h(i,j,k)
173 h(i,j,k) = max( gv%Angstrom_H, h(i,j,k) )
174 e_interface = e_interface + gv%H_to_Z * h(i,j,k)
176 h(i,j,1) = gv%Z_to_H * (e0(1) - e_interface)
177 h(i,j,1) = max( gv%Angstrom_H, h(i,j,1) )