MOM6
coord_hycom Module Reference

Detailed Description

Regrid columns for the HyCOM coordinate.

Data Types

type  hycom_cs
 Control structure containing required parameters for the HyCOM coordinate. More...
 

Functions/Subroutines

subroutine, public init_coord_hycom (CS, nk, coordinateResolution, target_density, interp_CS, rho_scale)
 Initialise a hycom_CS with pointers to parameters. More...
 
subroutine, public end_coord_hycom (CS)
 This subroutine deallocates memory in the control structure for the coord_hycom module. More...
 
subroutine, public set_hycom_params (CS, max_interface_depths, max_layer_thickness, interp_CS)
 This subroutine can be used to set the parameters for the coord_hycom module. More...
 
subroutine, public build_hycom1_column (CS, eqn_of_state, nz, depth, h, T, S, p_col, z_col, z_col_new, zScale, h_neglect, h_neglect_edge)
 Build a HyCOM coordinate column. More...
 

Function/Subroutine Documentation

◆ build_hycom1_column()

subroutine, public coord_hycom::build_hycom1_column ( type(hycom_cs), intent(in)  CS,
type(eos_type), pointer  eqn_of_state,
integer, intent(in)  nz,
real, intent(in)  depth,
real, dimension(nz), intent(in)  h,
real, dimension(nz), intent(in)  T,
real, dimension(nz), intent(in)  S,
real, dimension(nz), intent(in)  p_col,
real, dimension(nz+1), intent(in)  z_col,
real, dimension(cs%nk+1), intent(inout)  z_col_new,
real, intent(in), optional  zScale,
real, intent(in), optional  h_neglect,
real, intent(in), optional  h_neglect_edge 
)

Build a HyCOM coordinate column.

Parameters
[in]csCoordinate control structure
eqn_of_stateEquation of state structure
[in]nzNumber of levels
[in]depthDepth of ocean bottom (positive [H ~> m or kg m-2])
[in]tTemperature of column [degC]
[in]sSalinity of column [ppt]
[in]hLayer thicknesses, in [m] or [H ~> m or kg m-2]
[in]p_colLayer pressure [Pa]
[in]z_colInterface positions relative to the surface [H ~> m or kg m-2]
[in,out]z_col_newAbsolute positions of interfaces
[in]zscaleScaling factor from the input thicknesses in [m] to desired units for zInterface, perhaps m_to_H.
[in]h_neglectA negligibly small width for the purpose of cell reconstructions in the same units as h.
[in]h_neglect_edgeA negligibly small width for the purpose of edge value calculations in the same units as h0.

Definition at line 105 of file coord_hycom.F90.

105  type(hycom_CS), intent(in) :: CS !< Coordinate control structure
106  type(EOS_type), pointer :: eqn_of_state !< Equation of state structure
107  integer, intent(in) :: nz !< Number of levels
108  real, intent(in) :: depth !< Depth of ocean bottom (positive [H ~> m or kg m-2])
109  real, dimension(nz), intent(in) :: T !< Temperature of column [degC]
110  real, dimension(nz), intent(in) :: S !< Salinity of column [ppt]
111  real, dimension(nz), intent(in) :: h !< Layer thicknesses, in [m] or [H ~> m or kg m-2]
112  real, dimension(nz), intent(in) :: p_col !< Layer pressure [Pa]
113  real, dimension(nz+1), intent(in) :: z_col !< Interface positions relative to the surface [H ~> m or kg m-2]
114  real, dimension(CS%nk+1), intent(inout) :: z_col_new !< Absolute positions of interfaces
115  real, optional, intent(in) :: zScale !< Scaling factor from the input thicknesses in [m]
116  !! to desired units for zInterface, perhaps m_to_H.
117  real, optional, intent(in) :: h_neglect !< A negligibly small width for the
118  !! purpose of cell reconstructions
119  !! in the same units as h.
120  real, optional, intent(in) :: h_neglect_edge !< A negligibly small width
121  !! for the purpose of edge value calculations
122  !! in the same units as h0.
123 
124  ! Local variables
125  integer :: k
126  real, dimension(nz) :: rho_col ! Layer densities in a column [R ~> kg m-3]
127  real, dimension(CS%nk) :: h_col_new ! New layer thicknesses
128  real :: z_scale
129  real :: stretching ! z* stretching, converts z* to z.
130  real :: nominal_z ! Nominal depth of interface when using z* [Z ~> m]
131  real :: hNew
132  logical :: maximum_depths_set ! If true, the maximum depths of interface have been set.
133  logical :: maximum_h_set ! If true, the maximum layer thicknesses have been set.
134 
135  maximum_depths_set = allocated(cs%max_interface_depths)
136  maximum_h_set = allocated(cs%max_layer_thickness)
137 
138  z_scale = 1.0 ; if (present(zscale)) z_scale = zscale
139 
140  ! Work bottom recording potential density
141  call calculate_density(t, s, p_col, rho_col, 1, nz, eqn_of_state, scale=cs%kg_m3_to_R)
142  ! This ensures the potential density profile is monotonic
143  ! although not necessarily single valued.
144  do k = nz-1, 1, -1
145  rho_col(k) = min( rho_col(k), rho_col(k+1) )
146  enddo
147 
148  ! Interpolates for the target interface position with the rho_col profile
149  ! Based on global density profile, interpolate to generate a new grid
150  call build_and_interpolate_grid(cs%interp_CS, rho_col, nz, h(:), z_col, &
151  cs%target_density, cs%nk, h_col_new, z_col_new, h_neglect, h_neglect_edge)
152 
153  ! Sweep down the interfaces and make sure that the interface is at least
154  ! as deep as a nominal target z* grid
155  nominal_z = 0.
156  stretching = z_col(nz+1) / depth ! Stretches z* to z
157  do k = 2, cs%nk+1
158  nominal_z = nominal_z + (z_scale * cs%coordinateResolution(k-1)) * stretching
159  z_col_new(k) = max( z_col_new(k), nominal_z )
160  z_col_new(k) = min( z_col_new(k), z_col(nz+1) )
161  enddo
162 
163  if (maximum_depths_set .and. maximum_h_set) then ; do k=2,cs%nk
164  ! The loop bounds are 2 & nz so the top and bottom interfaces do not move.
165  ! Recall that z_col_new is positive downward.
166  z_col_new(k) = min(z_col_new(k), cs%max_interface_depths(k), &
167  z_col_new(k-1) + cs%max_layer_thickness(k-1))
168  enddo ; elseif (maximum_depths_set) then ; do k=2,cs%nk
169  z_col_new(k) = min(z_col_new(k), cs%max_interface_depths(k))
170  enddo ; elseif (maximum_h_set) then ; do k=2,cs%nk
171  z_col_new(k) = min(z_col_new(k), z_col_new(k-1) + cs%max_layer_thickness(k-1))
172  enddo ; endif

References regrid_interp::build_and_interpolate_grid().

Referenced by mom_regridding::build_grid_hycom1().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ end_coord_hycom()

subroutine, public coord_hycom::end_coord_hycom ( type(hycom_cs), pointer  CS)

This subroutine deallocates memory in the control structure for the coord_hycom module.

Parameters
csCoordinate control structure

Definition at line 65 of file coord_hycom.F90.

65  type(hycom_CS), pointer :: CS !< Coordinate control structure
66 
67  ! nothing to do
68  if (.not. associated(cs)) return
69  deallocate(cs%coordinateResolution)
70  deallocate(cs%target_density)
71  if (allocated(cs%max_interface_depths)) deallocate(cs%max_interface_depths)
72  if (allocated(cs%max_layer_thickness)) deallocate(cs%max_layer_thickness)
73  deallocate(cs)

Referenced by mom_regridding::end_regridding().

Here is the caller graph for this function:

◆ init_coord_hycom()

subroutine, public coord_hycom::init_coord_hycom ( type(hycom_cs), pointer  CS,
integer, intent(in)  nk,
real, dimension(nk), intent(in)  coordinateResolution,
real, dimension(nk+1), intent(in)  target_density,
type(interp_cs_type), intent(in)  interp_CS,
real, intent(in), optional  rho_scale 
)

Initialise a hycom_CS with pointers to parameters.

Parameters
csUnassociated pointer to hold the control structure
[in]nkNumber of layers in generated grid
[in]coordinateresolutionNominal near-surface resolution [m]
[in]target_densityInterface target densities [R ~> kg m-3]
[in]interp_csControls for interpolation
[in]rho_scaleA dimensional scaling factor for target_density

Definition at line 43 of file coord_hycom.F90.

43  type(hycom_CS), pointer :: CS !< Unassociated pointer to hold the control structure
44  integer, intent(in) :: nk !< Number of layers in generated grid
45  real, dimension(nk), intent(in) :: coordinateResolution !< Nominal near-surface resolution [m]
46  real, dimension(nk+1),intent(in) :: target_density !< Interface target densities [R ~> kg m-3]
47  type(interp_CS_type), intent(in) :: interp_CS !< Controls for interpolation
48  real, optional, intent(in) :: rho_scale !< A dimensional scaling factor for target_density
49 
50  if (associated(cs)) call mom_error(fatal, "init_coord_hycom: CS already associated!")
51  allocate(cs)
52  allocate(cs%coordinateResolution(nk))
53  allocate(cs%target_density(nk+1))
54 
55  cs%nk = nk
56  cs%coordinateResolution(:) = coordinateresolution(:)
57  cs%target_density(:) = target_density(:)
58  cs%interp_CS = interp_cs
59  cs%kg_m3_to_R = 1.0 ; if (present(rho_scale)) cs%kg_m3_to_R = rho_scale
60 

References mom_error_handler::mom_error().

Referenced by mom_regridding::initcoord().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ set_hycom_params()

subroutine, public coord_hycom::set_hycom_params ( type(hycom_cs), pointer  CS,
real, dimension(:), intent(in), optional  max_interface_depths,
real, dimension(:), intent(in), optional  max_layer_thickness,
type(interp_cs_type), intent(in), optional  interp_CS 
)

This subroutine can be used to set the parameters for the coord_hycom module.

Parameters
csCoordinate control structure
[in]max_interface_depthsMaximum depths of interfaces in m
[in]max_layer_thicknessMaximum thicknesses of layers in m
[in]interp_csControls for interpolation

Definition at line 78 of file coord_hycom.F90.

78  type(hycom_CS), pointer :: CS !< Coordinate control structure
79  real, dimension(:), optional, intent(in) :: max_interface_depths !< Maximum depths of interfaces in m
80  real, dimension(:), optional, intent(in) :: max_layer_thickness !< Maximum thicknesses of layers in m
81  type(interp_CS_type), optional, intent(in) :: interp_CS !< Controls for interpolation
82 
83  if (.not. associated(cs)) call mom_error(fatal, "set_hycom_params: CS not associated")
84 
85  if (present(max_interface_depths)) then
86  if (size(max_interface_depths) /= cs%nk+1) &
87  call mom_error(fatal, "set_hycom_params: max_interface_depths inconsistent size")
88  allocate(cs%max_interface_depths(cs%nk+1))
89  cs%max_interface_depths(:) = max_interface_depths(:)
90  endif
91 
92  if (present(max_layer_thickness)) then
93  if (size(max_layer_thickness) /= cs%nk) &
94  call mom_error(fatal, "set_hycom_params: max_layer_thickness inconsistent size")
95  allocate(cs%max_layer_thickness(cs%nk))
96  cs%max_layer_thickness(:) = max_layer_thickness(:)
97  endif
98 
99  if (present(interp_cs)) cs%interp_CS = interp_cs

References mom_error_handler::mom_error().

Referenced by mom_regridding::set_regrid_max_depths(), mom_regridding::set_regrid_max_thickness(), and mom_regridding::set_regrid_params().

Here is the call graph for this function:
Here is the caller graph for this function: