regrid_consts.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!> Contains constants for interpreting input parameters that control regridding.
6module regrid_consts
7
8use mom_error_handler, only : mom_error, fatal
10
11implicit none ; public
12
13! List of regridding types. These should be consecutive and starting at 1.
14! This allows them to be used as array indices.
15integer, parameter :: regridding_layer = 1 !< Layer mode identifier
16integer, parameter :: regridding_zstar = 2 !< z* coordinates identifier
17integer, parameter :: regridding_rho = 3 !< Density coordinates identifier
18integer, parameter :: regridding_sigma = 4 !< Sigma coordinates identifier
19integer, parameter :: regridding_arbitrary = 5 !< Arbitrary coordinates identifier
20integer, parameter :: regridding_hycom1 = 6 !< Simple HyCOM coordinates without BBL
21integer, parameter :: regridding_sigma_shelf_zstar = 8 !< Identifiered for z* coordinates at the bottom,
22 !! sigma-near the top
23integer, parameter :: regridding_adaptive = 9 !< Adaptive coordinate mode identifier
24integer, parameter :: regridding_hybgen = 10 !< Hybgen coordinates identifier
25
26character(len=*), parameter :: regridding_layer_string = "LAYER" !< Layer string
27character(len=*), parameter :: regridding_zstar_string_old = "Z*" !< z* string (legacy name)
28character(len=*), parameter :: regridding_zstar_string = "ZSTAR" !< z* string
29character(len=*), parameter :: regridding_rho_string = "RHO" !< Rho string
30character(len=*), parameter :: regridding_sigma_string = "SIGMA" !< Sigma string
31character(len=*), parameter :: regridding_arbitrary_string = "ARB" !< Arbitrary coordinates
32character(len=*), parameter :: regridding_hycom1_string = "HYCOM1" !< Hycom string
33character(len=*), parameter :: regridding_hybgen_string = "HYBGEN" !< Hybgen string
34character(len=*), parameter :: regridding_sigma_shelf_zstar_string = "SIGMA_SHELF_ZSTAR" !< Hybrid z*/sigma
35character(len=*), parameter :: regridding_adaptive_string = "ADAPTIVE" !< Adaptive coordinate string
36character(len=*), parameter :: default_coordinate_mode = regridding_layer_string !< Default coordinate mode
37
38!> Returns a string with the coordinate units associated with the coordinate mode.
39interface coordinateunits
40 module procedure coordinateunitsi
41 module procedure coordinateunitss
42end interface
43
44!> Returns true if the coordinate is dependent on the state density, returns false otherwise.
45interface state_dependent
46 module procedure state_dependent_char
47 module procedure state_dependent_int
48end interface
49
50contains
51
52!> Parse a string parameter specifying the coordinate mode and
53!! return the appropriate enumerated integer
54function coordinatemode(string)
55 integer :: coordinatemode !< Enumerated integer indicating coordinate mode
56 character(len=*), intent(in) :: string !< String to indicate coordinate mode
57 select case ( uppercase(trim(string)) )
58 case (trim(regridding_layer_string)); coordinatemode = regridding_layer
59 case (trim(regridding_zstar_string)); coordinatemode = regridding_zstar
60 case (trim(regridding_zstar_string_old)); coordinatemode = regridding_zstar
61 case (trim(regridding_rho_string)); coordinatemode = regridding_rho
62 case (trim(regridding_sigma_string)); coordinatemode = regridding_sigma
63 case (trim(regridding_hycom1_string)); coordinatemode = regridding_hycom1
64 case (trim(regridding_hybgen_string)); coordinatemode = regridding_hybgen
65 case (trim(regridding_arbitrary_string)); coordinatemode = regridding_arbitrary
66 case (trim(regridding_sigma_shelf_zstar_string)); coordinatemode = regridding_sigma_shelf_zstar
67 case (trim(regridding_adaptive_string)); coordinatemode = regridding_adaptive
68 case default ; call mom_error(fatal, "coordinateMode: "//&
69 "Unrecognized choice of coordinate ("//trim(string)//").")
70 end select
71end function coordinatemode
72
73!> Returns a string with the coordinate units associated with the
74!! enumerated integer,
75function coordinateunitsi(coordMode)
76 character(len=16) :: coordinateunitsi !< Units of coordinate
77 integer, intent(in) :: coordmode !< Coordinate mode
78 select case ( coordmode )
79 case (regridding_layer); coordinateunitsi = "kg m^-3"
80 case (regridding_zstar); coordinateunitsi = "m"
81 case (regridding_sigma_shelf_zstar); coordinateunitsi = "m"
82 case (regridding_rho); coordinateunitsi = "kg m^-3"
83 case (regridding_sigma); coordinateunitsi = "Non-dimensional"
84 case (regridding_hycom1); coordinateunitsi = "m"
85 case (regridding_hybgen); coordinateunitsi = "m"
86 case (regridding_adaptive); coordinateunitsi = "m"
87 case default ; call mom_error(fatal, "coordinateUnts: "//&
88 "Unrecognized coordinate mode.")
89 end select
90end function coordinateunitsi
91
92!> Returns a string with the coordinate units associated with the
93!! string defining the coordinate mode.
94function coordinateunitss(string)
95 character(len=16) :: coordinateunitss !< Units of coordinate
96 character(len=*), intent(in) :: string !< Coordinate mode
97 integer :: coordmode
98 coordmode = coordinatemode(string)
100end function coordinateunitss
101
102!> Returns true if the coordinate is dependent on the state density, returns false otherwise.
103logical function state_dependent_char(string)
104 character(len=*), intent(in) :: string !< String to indicate coordinate mode
105
107
108end function state_dependent_char
109
110!> Returns true if the coordinate is dependent on the state density, returns false otherwise.
111logical function state_dependent_int(mode)
112 integer, intent(in) :: mode !< Coordinate mode
113 select case ( mode )
114 case (regridding_layer); state_dependent_int = .true.
115 case (regridding_zstar); state_dependent_int = .false.
116 case (regridding_sigma_shelf_zstar); state_dependent_int = .false.
117 case (regridding_rho); state_dependent_int = .true.
118 case (regridding_sigma); state_dependent_int = .false.
119 case (regridding_hycom1); state_dependent_int = .true.
120 case (regridding_hybgen); state_dependent_int = .true.
121 case (regridding_adaptive); state_dependent_int = .true.
122 case default ; call mom_error(fatal, "state_dependent: "//&
123 "Unrecognized choice of coordinate.")
124 end select
125end function state_dependent_int
126
127end module regrid_consts