MOM6
ocn_cpl_indices Module Reference

Data Types

type  cpl_indices_type
 Structure with indices needed for MCT attribute vectors. More...
 

Functions/Subroutines

subroutine, public cpl_indices_init (ind)
 Determines attribute vector indices. More...
 

Function/Subroutine Documentation

◆ cpl_indices_init()

subroutine, public ocn_cpl_indices::cpl_indices_init ( type(cpl_indices_type), intent(inout)  ind)

Determines attribute vector indices.

Parameters
[in,out]indStructure with coupler indices and vectors

Definition at line 85 of file ocn_cpl_indices.F90.

85  type(cpl_indices_type), intent(inout) :: ind !< Structure with coupler indices and vectors
86 
87  ! Local Variables
88  type(mct_aVect) :: o2x !< Array with ocean to coupler data
89  type(mct_aVect) :: x2o !< Array with coupler to ocean data
90  integer :: ncat !< Thickness category index
91  character(len=2) :: cncat !< Character version of ncat
92  integer :: ncol !< Column index
93  integer :: mcog_ncols !< Number of ice thickness categories?
94  integer :: lmcog_flds_sent !< Used to convert per thickness category fields?
95 
96  ! create temporary attribute vectors
97  call mct_avect_init(x2o, rlist=seq_flds_x2o_fields, lsize=1)
98  call mct_avect_init(o2x, rlist=seq_flds_o2x_fields, lsize=1)
99 
100  ! ocean to coupler
101  ind%o2x_So_t = mct_avect_indexra(o2x,'So_t')
102  ind%o2x_So_u = mct_avect_indexra(o2x,'So_u')
103  ind%o2x_So_v = mct_avect_indexra(o2x,'So_v')
104  ind%o2x_So_s = mct_avect_indexra(o2x,'So_s')
105  ind%o2x_So_dhdx = mct_avect_indexra(o2x,'So_dhdx')
106  ind%o2x_So_dhdy = mct_avect_indexra(o2x,'So_dhdy')
107  ind%o2x_So_bldepth = mct_avect_indexra(o2x,'So_bldepth')
108  ind%o2x_Fioo_q = mct_avect_indexra(o2x,'Fioo_q')
109  ind%o2x_Faoo_fco2_ocn = mct_avect_indexra(o2x,'Faoo_fco2_ocn',perrwith='quiet')
110  ind%o2x_Faoo_fdms_ocn = mct_avect_indexra(o2x,'Faoo_fdms_ocn',perrwith='quiet')
111 
112  ! coupler to ocean
113  ind%x2o_Si_ifrac = mct_avect_indexra(x2o,'Si_ifrac')
114  ind%x2o_Sa_pslv = mct_avect_indexra(x2o,'Sa_pslv')
115  ind%x2o_So_duu10n = mct_avect_indexra(x2o,'So_duu10n')
116  ind%x2o_Sw_lamult = mct_avect_indexra(x2o,'Sw_lamult')
117  ind%x2o_Sw_ustokes = mct_avect_indexra(x2o,'Sw_ustokes')
118  ind%x2o_Sw_vstokes = mct_avect_indexra(x2o,'Sw_vstokes')
119  ind%x2o_Foxx_tauy = mct_avect_indexra(x2o,'Foxx_tauy')
120  ind%x2o_Foxx_taux = mct_avect_indexra(x2o,'Foxx_taux')
121  ind%x2o_Foxx_swnet = mct_avect_indexra(x2o,'Foxx_swnet')
122  ind%x2o_Foxx_lat = mct_avect_indexra(x2o,'Foxx_lat')
123  ind%x2o_Foxx_sen = mct_avect_indexra(x2o,'Foxx_sen')
124  ind%x2o_Foxx_lwup = mct_avect_indexra(x2o,'Foxx_lwup')
125  ind%x2o_Faxa_lwdn = mct_avect_indexra(x2o,'Faxa_lwdn')
126  ind%x2o_Faxa_swvdr = mct_avect_indexra(x2o,'Faxa_swvdr',perrwith='quiet')
127  ind%x2o_Faxa_swvdf = mct_avect_indexra(x2o,'Faxa_swvdf',perrwith='quiet')
128  ind%x2o_Faxa_swndr = mct_avect_indexra(x2o,'Faxa_swndr',perrwith='quiet')
129  ind%x2o_Faxa_swndf = mct_avect_indexra(x2o,'Faxa_swndf',perrwith='quiet')
130  ind%x2o_Fioi_melth = mct_avect_indexra(x2o,'Fioi_melth')
131  ind%x2o_Fioi_meltw = mct_avect_indexra(x2o,'Fioi_meltw')
132  ind%x2o_Fioi_salt = mct_avect_indexra(x2o,'Fioi_salt')
133  ind%x2o_Fioi_bcpho = mct_avect_indexra(x2o,'Fioi_bcpho')
134  ind%x2o_Fioi_bcphi = mct_avect_indexra(x2o,'Fioi_bcphi')
135  ind%x2o_Fioi_flxdst = mct_avect_indexra(x2o,'Fioi_flxdst')
136  ind%x2o_Faxa_prec = mct_avect_indexra(x2o,'Faxa_prec')
137  ind%x2o_Faxa_snow = mct_avect_indexra(x2o,'Faxa_snow')
138  ind%x2o_Faxa_rain = mct_avect_indexra(x2o,'Faxa_rain')
139  ind%x2o_Foxx_evap = mct_avect_indexra(x2o,'Foxx_evap')
140  ind%x2o_Foxx_rofl = mct_avect_indexra(x2o,'Foxx_rofl')
141  ind%x2o_Foxx_rofi = mct_avect_indexra(x2o,'Foxx_rofi')
142  ind%x2o_Faxa_bcphidry = mct_avect_indexra(x2o,'Faxa_bcphidry')
143  ind%x2o_Faxa_bcphodry = mct_avect_indexra(x2o,'Faxa_bcphodry')
144  ind%x2o_Faxa_bcphiwet = mct_avect_indexra(x2o,'Faxa_bcphiwet')
145  ind%x2o_Faxa_ocphidry = mct_avect_indexra(x2o,'Faxa_ocphidry')
146  ind%x2o_Faxa_ocphodry = mct_avect_indexra(x2o,'Faxa_ocphodry')
147  ind%x2o_Faxa_ocphiwet = mct_avect_indexra(x2o,'Faxa_ocphiwet')
148  ind%x2o_Faxa_dstdry1 = mct_avect_indexra(x2o,'Faxa_dstdry1')
149  ind%x2o_Faxa_dstdry2 = mct_avect_indexra(x2o,'Faxa_dstdry2')
150  ind%x2o_Faxa_dstdry3 = mct_avect_indexra(x2o,'Faxa_dstdry3')
151  ind%x2o_Faxa_dstdry4 = mct_avect_indexra(x2o,'Faxa_dstdry4')
152  ind%x2o_Faxa_dstwet1 = mct_avect_indexra(x2o,'Faxa_dstwet1')
153  ind%x2o_Faxa_dstwet2 = mct_avect_indexra(x2o,'Faxa_dstwet2')
154  ind%x2o_Faxa_dstwet3 = mct_avect_indexra(x2o,'Faxa_dstwet3')
155  ind%x2o_Faxa_dstwet4 = mct_avect_indexra(x2o,'Faxa_dstwet4')
156  ind%x2o_Sa_co2prog = mct_avect_indexra(x2o,'Sa_co2prog',perrwith='quiet')
157  ind%x2o_Sa_co2diag = mct_avect_indexra(x2o,'Sa_co2diag',perrwith='quiet')
158 
159  ! optional per thickness category fields
160  ! convert cpl indices to mcog column indices
161  ! this implementation only handles columns due to ice thickness categories
162  lmcog_flds_sent = seq_flds_i2o_per_cat
163 
164  if (seq_flds_i2o_per_cat) then
165  mcog_ncols = ice_ncat+1
166  allocate(ind%x2o_frac_col(mcog_ncols))
167  allocate(ind%x2o_fracr_col(mcog_ncols))
168  allocate(ind%x2o_qsw_fracr_col(mcog_ncols))
169  ncol = 1
170  ind%x2o_frac_col(ncol) = mct_avect_indexra(x2o,'Sf_afrac')
171  ind%x2o_fracr_col(ncol) = mct_avect_indexra(x2o,'Sf_afracr')
172  ind%x2o_qsw_fracr_col(ncol) = mct_avect_indexra(x2o,'Foxx_swnet_afracr')
173 
174  do ncat = 1, ice_ncat
175  write(cncat,'(i2.2)') ncat
176  ncol = ncat+1
177  ind%x2o_frac_col(ncol) = mct_avect_indexra(x2o,'Si_ifrac_'//cncat)
178  ind%x2o_fracr_col(ncol) = ind%x2o_frac_col(ncol)
179  ind%x2o_qsw_fracr_col(ncol) = mct_avect_indexra(x2o,'PFioi_swpen_ifrac_'//cncat)
180  enddo
181  else
182  mcog_ncols = 1
183  endif
184 
185  call mct_avect_clean(x2o)
186  call mct_avect_clean(o2x)
187 

Referenced by ocn_comp_mct::ocn_init_mct().

Here is the caller graph for this function: