Go to the documentation of this file.
17 implicit none ;
private
19 #include <MOM_memory.h>
34 real,
dimension(SZI_(G),SZJ_(G), SZK_(GV)), &
38 logical,
optional,
intent(in) :: just_read_params
41 real :: e0(szk_(gv)+1)
43 real :: eta1d(szk_(gv)+1)
46 real :: diskrad, rad, xcenter, xradius, lonc, latc, xoffset
49 # include "version_variable.h"
50 character(len=40) :: mdl =
"circle_obcs_initialization"
51 integer :: i, j, k, is, ie, js, je, nz
53 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
55 just_read = .false. ;
if (
present(just_read_params)) just_read = just_read_params
58 call mom_mesg(
" circle_obcs_initialization.F90, circle_obcs_initialize_thickness: setting thickness", 5)
60 if (.not.just_read)
call log_version(param_file, mdl, version,
"")
62 call get_param(param_file, mdl,
"DISK_RADIUS", diskrad, &
63 "The radius of the initially elevated disk in the "//&
64 "circle_obcs test case.", units=g%x_axis_units, &
65 fail_if_missing=.not.just_read, do_not_log=just_read)
66 call get_param(param_file, mdl,
"DISK_X_OFFSET", xoffset, &
67 "The x-offset of the initially elevated disk in the "//&
68 "circle_obcs test case.", units=g%x_axis_units, &
69 default = 0.0, do_not_log=just_read)
70 call get_param(param_file, mdl,
"DISK_IC_AMPLITUDE", ic_amp, &
71 "Initial amplitude of interface height displacements "//&
72 "in the circle_obcs test case.", &
73 units=
'm', default=5.0, scale=gv%m_to_H, do_not_log=just_read)
78 e0(k) = -g%max_depth * real(k-1) / real(nz)
82 do j=js,je ;
do i=is,ie
83 eta1d(nz+1) = -g%bathyT(i,j)
86 if (eta1d(k) < (eta1d(k+1) + gv%Angstrom_Z))
then
87 eta1d(k) = eta1d(k+1) + gv%Angstrom_Z
88 h(i,j,k) = gv%Angstrom_H
90 h(i,j,k) = gv%Z_to_H * (eta1d(k) - eta1d(k+1))
97 latc = g%south_lat + 0.5*g%len_lat
98 lonc = g%west_lon + 0.5*g%len_lon + xoffset
99 do j=js,je ;
do i=is,ie
100 rad = sqrt((g%geoLonT(i,j)-lonc)**2+(g%geoLatT(i,j)-latc)**2)/(diskrad)
103 rad = rad*(2.*asin(1.))
106 h(i,j,k) = h(i,j,k) + ic_amp * 0.5*(1.+cos(rad))
110 h(i,j,k) = h(i,j,k) - 0.5*(1.+cos(rad)) * ic_amp * real( 2*k-nz )
Configures the model for the "circle_obcs" experiment which tests Open Boundary Conditions radiating ...
Provides a transparent vertical ocean grid type and supporting routines.
subroutine, public set_up_sponge_field(sp_val, f_ptr, G, nlay, CS, sp_val_i_mean)
This subroutine stores the reference profile for the variable whose address is given by f_ptr....
An overloaded interface to log version information about modules.
Pointers to an assortment of thermodynamic fields that may be available, including potential temperat...
subroutine, public mom_mesg(message, verb, all_print)
This provides a convenient interface for writing an informative comment.
Provides subroutines for quantities specific to the equation of state.
This module contains the tracer_registry_type and the subroutines that handle registration of tracers...
A structure that can be parsed to read and document run-time parameters.
subroutine, public circle_obcs_initialize_thickness(h, G, GV, param_file, just_read_params)
This subroutine initializes layer thicknesses for the circle_obcs experiment.
An overloaded interface to read and log the values of various types of parameters.
A control structure for the equation of state.
Describes the vertical ocean grid, including unit conversion factors.
Provides transparent structures with groups of MOM6 variables and supporting routines.
The MOM6 facility to parse input files for runtime parameters.
Implements sponge regions in isopycnal mode.
Calculate the derivatives of density with temperature and salinity from T, S, and P.
Type to carry basic tracer information.
Provides the ocean grid type.
This control structure holds memory and parameters for the MOM_sponge module.
logical function, public is_root_pe()
This returns .true. if the current PE is the root PE.
subroutine, public mom_error(level, message, all_print)
This provides a convenient interface for writing an mpp_error message with run-time filter based on a...
subroutine, public initialize_sponge(Iresttime, int_height, G, param_file, CS, GV, Iresttime_i_mean, int_height_i_mean)
This subroutine determines the number of points which are within sponges in this computational domain...
Routines for error handling and I/O management.
Ocean grid type. See mom_grid for details.
Calculates density of sea water from T, S and P.