external_gwave_initialization.F90

1! This file is part of MOM6, the Modular Ocean Model version 6.
2! See the LICENSE file for licensing information.
3! SPDX-License-Identifier: Apache-2.0
4
5!> Initialization for the "external gravity wave wave" configuration
7
8use mom_error_handler, only : mom_mesg, mom_error, fatal, is_root_pe
9use mom_file_parser, only : get_param, log_version, param_file_type
10use mom_get_input, only : directories
11use mom_grid, only : ocean_grid_type
12use mom_tracer_registry, only : tracer_registry_type
16implicit none ; private
17
18#include <MOM_memory.h>
19
21
22! A note on unit descriptions in comments: MOM6 uses units that can be rescaled for dimensional
23! consistency testing. These are noted in comments with units like Z, H, L, and T, along with
24! their mks counterparts with notation like "a velocity [Z T-1 ~> m s-1]". If the units
25! vary with the Boussinesq approximation, the Boussinesq variant is given first.
26
27contains
28
29!> This subroutine initializes layer thicknesses for the external_gwave experiment.
30subroutine external_gwave_initialize_thickness(h, G, GV, US, param_file, just_read)
31 type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure.
32 type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure.
33 type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
34 real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), &
35 intent(out) :: h !< The thickness that is being initialized [Z ~> m]
36 type(param_file_type), intent(in) :: param_file !< A structure indicating the open file
37 !! to parse for model parameter values.
38 logical, intent(in) :: just_read !< If true, this call will only read
39 !! parameters without changing h.
40 ! Local variables
41 real :: eta1d(szk_(gv)+1) ! Interface height relative to the sea surface
42 ! positive upward [Z ~> m].
43 real :: ssh_anomaly_height ! Vertical height of ssh anomaly [Z ~> m]
44 real :: ssh_anomaly_width ! Lateral width of anomaly, often in [km] or [degrees_E]
45 character(len=40) :: mdl = "external_gwave_initialize_thickness" ! This subroutine's name.
46 ! This include declares and sets the variable "version".
47# include "version_variable.h"
48 integer :: i, j, k, is, ie, js, je, nz
49 real :: pi ! The ratio of the circumference of a circle to its diameter [nondim]
50 real :: xnondim ! A normalized x position [nondim]
51
52 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
53
54 if (.not.just_read) &
55 call mom_mesg(" external_gwave_initialization.F90, external_gwave_initialize_thickness: setting thickness", 5)
56
57 if (.not.just_read) call log_version(param_file, mdl, version, "")
58 call get_param(param_file, mdl, "SSH_ANOMALY_HEIGHT", ssh_anomaly_height, &
59 "The vertical displacement of the SSH anomaly. ", &
60 units="m", scale=us%m_to_Z, fail_if_missing=.not.just_read, do_not_log=just_read)
61 call get_param(param_file, mdl, "SSH_ANOMALY_WIDTH", ssh_anomaly_width, &
62 "The lateral width of the SSH anomaly. ", &
63 units=g%x_ax_unit_short, fail_if_missing=.not.just_read, do_not_log=just_read)
64
65 if (just_read) return ! All run-time parameters have been read, so return.
66
67 pi = 4.0*atan(1.0)
68 do j=g%jsc,g%jec ; do i=g%isc,g%iec
69 xnondim = (g%geoLonT(i,j)-g%west_lon-0.5*g%len_lon) / ssh_anomaly_width
70 xnondim = min(1., abs(xnondim))
71 eta1d(1) = ssh_anomaly_height * 0.5 * ( 1. + cos(pi*xnondim) ) ! Cosine bell
72 do k=2,nz
73 eta1d(k) = -g%max_depth & ! Stretch interior interfaces with SSH
74 + (eta1d(1)+g%max_depth) * ( real(nz+1-k)/real(nz) ) ! Stratification
75 enddo
76 eta1d(nz+1) = -g%max_depth ! Force bottom interface to bottom
77 do k=1,nz
78 h(i,j,k) = eta1d(k) - eta1d(k+1)
79 enddo
80 enddo ; enddo
81
83