MOM6
mom_hor_index Module Reference

Detailed Description

Defines the horizontal index type (hor_index_type) used for providing index ranges.

The hor_index_type provides the decalarations and loop ranges for almost all data with horizontal extent.

Declarations and loop ranges should always be coded with the symmetric memory model in mind. The non-symmetric memory mode will then also work, albeit with a different (less efficient) communication pattern.

Using the hor_index_type HI:

  • declaration of h-point data is of the form h(HI%isd:HI%ied,HI%jsd:HI%jed)
  • declaration of q-point data is of the form q(HI%IsdB:HI%IedB,HI%JsdB:HI%JedB)
  • declaration of u-point data is of the form u(HI%IsdB:HI%IedB,HI%jsd:HI%jed)
  • declaration of v-point data is of the form v(HI%isd:HI%ied,HI%JsdB:HI%JedB).

For more detail explanation of horizontal indexing see Horizontal indexing and memory.

Data Types

interface  assignment(=)
 Copy the contents of one horizontal index type into another. More...
 
type  hor_index_type
 Container for horizontal index ranges for data, computational and global domains. More...
 

Functions/Subroutines

subroutine, public hor_index_init (Domain, HI, param_file, local_indexing, index_offset)
 Sets various index values in a hor_index_type. More...
 
subroutine hit_assign (HI1, HI2)
 HIT_assign copies one hor_index_type into another. It is accessed via an assignment (=) operator. More...
 

Function/Subroutine Documentation

◆ hit_assign()

subroutine mom_hor_index::hit_assign ( type(hor_index_type), intent(out)  HI1,
type(hor_index_type), intent(in)  HI2 
)
private

HIT_assign copies one hor_index_type into another. It is accessed via an assignment (=) operator.

Parameters
[out]hi1Horizontal index type to copy to
[in]hi2Horizontal index type to copy from

Definition at line 96 of file MOM_hor_index.F90.

96  type(hor_index_type), intent(out) :: HI1 !< Horizontal index type to copy to
97  type(hor_index_type), intent(in) :: HI2 !< Horizontal index type to copy from
98  ! This subroutine copies all components of the horizontal array index type
99  ! variable on the RHS (HI2) to the variable on the LHS (HI1).
100 
101  hi1%isc = hi2%isc ; hi1%iec = hi2%iec ; hi1%jsc = hi2%jsc ; hi1%jec = hi2%jec
102  hi1%isd = hi2%isd ; hi1%ied = hi2%ied ; hi1%jsd = hi2%jsd ; hi1%jed = hi2%jed
103  hi1%isg = hi2%isg ; hi1%ieg = hi2%ieg ; hi1%jsg = hi2%jsg ; hi1%jeg = hi2%jeg
104 
105  hi1%IscB = hi2%IscB ; hi1%IecB = hi2%IecB ; hi1%JscB = hi2%JscB ; hi1%JecB = hi2%JecB
106  hi1%IsdB = hi2%IsdB ; hi1%IedB = hi2%IedB ; hi1%JsdB = hi2%JsdB ; hi1%JedB = hi2%JedB
107  hi1%IsgB = hi2%IsgB ; hi1%IegB = hi2%IegB ; hi1%JsgB = hi2%JsgB ; hi1%JegB = hi2%JegB
108 
109  hi1%idg_offset = hi2%idg_offset ; hi1%jdg_offset = hi2%jdg_offset
110  hi1%symmetric = hi2%symmetric
111 

◆ hor_index_init()

subroutine, public mom_hor_index::hor_index_init ( type(mom_domain_type), intent(in)  Domain,
type(hor_index_type), intent(inout)  HI,
type(param_file_type), intent(in)  param_file,
logical, intent(in), optional  local_indexing,
integer, intent(in), optional  index_offset 
)

Sets various index values in a hor_index_type.

Parameters
[in]domainThe MOM domain from which to extract information.
[in,out]hiA horizontal index type to populate with data
[in]param_fileParameter file handle
[in]local_indexingIf true, all tracer data domains start at 1
[in]index_offsetA fixed additional offset to all indices

Definition at line 58 of file MOM_hor_index.F90.

58  type(MOM_domain_type), intent(in) :: Domain !< The MOM domain from which to extract information.
59  type(hor_index_type), intent(inout) :: HI !< A horizontal index type to populate with data
60  type(param_file_type), intent(in) :: param_file !< Parameter file handle
61  logical, optional, intent(in) :: local_indexing !< If true, all tracer data domains start at 1
62  integer, optional, intent(in) :: index_offset !< A fixed additional offset to all indices
63 
64 ! This include declares and sets the variable "version".
65 #include "version_variable.h"
66 
67  ! get_domain_extent ensures that domains start at 1 for compatibility between
68  ! static and dynamically allocated arrays.
69  call get_domain_extent(domain, hi%isc, hi%iec, hi%jsc, hi%jec, &
70  hi%isd, hi%ied, hi%jsd, hi%jed, &
71  hi%isg, hi%ieg, hi%jsg, hi%jeg, &
72  hi%idg_offset, hi%jdg_offset, hi%symmetric, &
73  local_indexing=local_indexing)
74 
75  ! Read all relevant parameters and write them to the model log.
76  call log_version(param_file, "MOM_hor_index", version, &
77  "Sets the horizontal array index types.")
78 
79  hi%IscB = hi%isc ; hi%JscB = hi%jsc
80  hi%IsdB = hi%isd ; hi%JsdB = hi%jsd
81  hi%IsgB = hi%isg ; hi%JsgB = hi%jsg
82  if (hi%symmetric) then
83  hi%IscB = hi%isc-1 ; hi%JscB = hi%jsc-1
84  hi%IsdB = hi%isd-1 ; hi%JsdB = hi%jsd-1
85  hi%IsgB = hi%isg-1 ; hi%JsgB = hi%jsg-1
86  endif
87  hi%IecB = hi%iec ; hi%JecB = hi%jec
88  hi%IedB = hi%ied ; hi%JedB = hi%jed
89  hi%IegB = hi%ieg ; hi%JegB = hi%jeg
90 

References mom_domains::get_domain_extent().

Referenced by mom_oda_driver_mod::init_oda(), mom::initialize_mom(), and mom_grid::mom_grid_init().

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