coord_sigma.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!> Regrid columns for the sigma coordinate
6module coord_sigma
7
8use mom_error_handler, only : mom_error, fatal
9
10implicit none ; private
11
12!> Control structure containing required parameters for the sigma coordinate
13type, public :: sigma_cs ; private
14
15 !> Number of levels
16 integer :: nk
17
18 !> Minimum thickness allowed for layers [H ~> m or kg m-2]
19 real :: min_thickness
20
21 !> Target coordinate resolution [nondim]
22 real, allocatable, dimension(:) :: coordinateresolution
23end type sigma_cs
24
26
27contains
28
29!> Initialise a sigma_CS with pointers to parameters
30subroutine init_coord_sigma(CS, nk, coordinateResolution)
31 type(sigma_cs), pointer :: cs !< Unassociated pointer to hold the control structure
32 integer, intent(in) :: nk !< Number of layers in the grid
33 real, dimension(:), intent(in) :: coordinateresolution !< Nominal coordinate resolution [nondim]
34
35 if (associated(cs)) call mom_error(fatal, "init_coord_sigma: CS already associated!")
36 allocate(cs)
37 allocate(cs%coordinateResolution(nk))
38
39 cs%nk = nk
40 cs%coordinateResolution = coordinateresolution
41end subroutine init_coord_sigma
42
43!> This subroutine deallocates memory in the control structure for the coord_sigma module
44subroutine end_coord_sigma(CS)
45 type(sigma_cs), pointer :: cs !< Coordinate control structure
46
47 ! nothing to do
48 if (.not. associated(cs)) return
49 deallocate(cs%coordinateResolution)
50 deallocate(cs)
51end subroutine end_coord_sigma
52
53!> This subroutine can be used to set the parameters for the coord_sigma module
54subroutine set_sigma_params(CS, min_thickness)
55 type(sigma_cs), pointer :: cs !< Coordinate control structure
56 real, optional, intent(in) :: min_thickness !< Minimum allowed thickness [H ~> m or kg m-2]
57
58 if (.not. associated(cs)) call mom_error(fatal, "set_sigma_params: CS not associated")
59
60 if (present(min_thickness)) cs%min_thickness = min_thickness
61end subroutine set_sigma_params
62
63
64!> Build a sigma coordinate column
65subroutine build_sigma_column(CS, depth, totalThickness, zInterface)
66 type(sigma_cs), intent(in) :: cs !< Coordinate control structure
67 real, intent(in) :: depth !< Depth of ocean bottom (positive [H ~> m or kg m-2])
68 real, intent(in) :: totalthickness !< Column thickness (positive [H ~> m or kg m-2])
69 real, dimension(CS%nk+1), intent(inout) :: zinterface !< Absolute positions of interfaces [H ~> m or kg m-2]
70
71 ! Local variables
72 integer :: k
73
74 zinterface(cs%nk+1) = -depth
75 do k = cs%nk,1,-1
76 zinterface(k) = zinterface(k+1) + (totalthickness * cs%coordinateResolution(k))
77 ! Adjust interface position to accommodate inflating layers
78 ! without disturbing the interface above
79 if (zinterface(k) < (zinterface(k+1) + cs%min_thickness)) then
80 zinterface(k) = zinterface(k+1) + cs%min_thickness
81 endif
82 enddo
83end subroutine build_sigma_column
84
85end module coord_sigma