MOM6
mom_cap_methods.F90
Go to the documentation of this file.
1 !> Contains import/export methods for both NEMS and CMEPS.
3 
4 use esmf, only: esmf_clock, esmf_clockget, esmf_time, esmf_timeget
5 use esmf, only: esmf_timeinterval, esmf_timeintervalget
6 use esmf, only: esmf_state, esmf_stateget
7 use esmf, only: esmf_field, esmf_fieldget, esmf_fieldcreate
8 use esmf, only: esmf_gridcomp, esmf_mesh, esmf_grid, esmf_gridcreate
9 use esmf, only: esmf_distgrid, esmf_distgridcreate
10 use esmf, only: esmf_kind_r8, esmf_success, esmf_logfounderror
11 use esmf, only: esmf_logerr_passthru, esmf_logmsg_info, esmf_logwrite
12 use esmf, only: esmf_logseterror, esmf_rc_mem_allocate
13 use esmf, only: esmf_stateitem_flag, esmf_stateitem_notfound
14 use esmf, only: esmf_geomtype_flag, esmf_geomtype_grid, esmf_geomtype_mesh
15 use esmf, only: esmf_rc_val_outofrange, esmf_index_delocal, esmf_meshloc_element
16 use esmf, only: esmf_typekind_r8
17 use esmf, only: operator(/=), operator(==)
20 use mom_grid, only: ocean_grid_type
21 use mom_domains, only: pass_var
22 use mpp_domains_mod, only: mpp_get_compute_domain
23 
24 ! By default make data private
25 implicit none; private
26 
27 ! Public member functions
28 public :: mom_set_geomtype
29 public :: mom_import
30 public :: mom_export
31 
32 private :: state_getimport
33 private :: state_setexport
34 
35 !> Get field pointer
36 interface state_getfldptr
37  module procedure state_getfldptr_1d
38  module procedure state_getfldptr_2d
39 end interface
40 
41 integer :: import_cnt = 0!< used to skip using the import state
42  !! at the first count for cesm
43 type(esmf_geomtype_flag) :: geomtype !< SMF type describing type of
44  !! geometry (mesh or grid)
45 
46 contains
47 
48 !> Sets module variable geometry type
49 subroutine mom_set_geomtype(geomtype_in)
50  type(esmf_geomtype_flag), intent(in) :: geomtype_in !< ESMF type describing type of
51  !! geometry (mesh or grid)
52 
53  geomtype = geomtype_in
54 
55 end subroutine mom_set_geomtype
56 
57 !> This function has a few purposes:
58 !! (1) it imports surface fluxes using data from the mediator; and
59 !! (2) it can apply restoring in SST and SSS.
60 subroutine mom_import(ocean_public, ocean_grid, importState, ice_ocean_boundary, rc)
61  type(ocean_public_type) , intent(in) :: ocean_public !< Ocean surface state
62  type(ocean_grid_type) , intent(in) :: ocean_grid !< Ocean model grid
63  type(esmf_state) , intent(inout) :: importstate !< incoming data from mediator
64  type(ice_ocean_boundary_type) , intent(inout) :: ice_ocean_boundary !< Ocean boundary forcing
65  integer , intent(inout) :: rc !< Return code
66 
67  ! Local Variables
68  integer :: i, j, ig, jg, n
69  integer :: isc, iec, jsc, jec
70  character(len=128) :: fldname
71  real(esmf_kind_r8), allocatable :: taux(:,:)
72  real(esmf_kind_r8), allocatable :: tauy(:,:)
73  character(len=*) , parameter :: subname = '(mom_import)'
74 
75  rc = esmf_success
76 
77  ! -------
78  ! import_cnt is used to skip using the import state at the first count for cesm
79  ! -------
80 
81  ! The following are global indices without halos
82  call mpp_get_compute_domain(ocean_public%domain, isc, iec, jsc, jec)
83 
84  !----
85  ! surface height pressure
86  !----
87  call state_getimport(importstate, 'inst_pres_height_surface', &
88  isc, iec, jsc, jec, ice_ocean_boundary%p, rc=rc)
89  if (esmf_logfounderror(rctocheck=rc, msg=esmf_logerr_passthru, &
90  line=__line__, &
91  file=__file__)) &
92  return ! bail out
93 
94  !----
95  ! near-IR, direct shortwave (W/m2)
96  !----
97  call state_getimport(importstate, 'mean_net_sw_ir_dir_flx', &
98  isc, iec, jsc, jec, ice_ocean_boundary%sw_flux_nir_dir, rc=rc)
99  if (esmf_logfounderror(rctocheck=rc, msg=esmf_logerr_passthru, &
100  line=__line__, &
101  file=__file__)) &
102  return ! bail out
103 
104  !----
105  ! near-IR, diffuse shortwave (W/m2)
106  !----
107  call state_getimport(importstate, 'mean_net_sw_ir_dif_flx', &
108  isc, iec, jsc, jec, ice_ocean_boundary%sw_flux_nir_dif, rc=rc)
109  if (esmf_logfounderror(rctocheck=rc, msg=esmf_logerr_passthru, &
110  line=__line__, &
111  file=__file__)) &
112  return ! bail out
113 
114  !----
115  ! visible, direct shortwave (W/m2)
116  !----
117  call state_getimport(importstate, 'mean_net_sw_vis_dir_flx', &
118  isc, iec, jsc, jec, ice_ocean_boundary%sw_flux_vis_dir, rc=rc)
119  if (esmf_logfounderror(rctocheck=rc, msg=esmf_logerr_passthru, &
120  line=__line__, &
121  file=__file__)) &
122  return ! bail out
123 
124  !----
125  ! visible, diffuse shortwave (W/m2)
126  !----
127  call state_getimport(importstate, 'mean_net_sw_vis_dif_flx', &
128  isc, iec, jsc, jec, ice_ocean_boundary%sw_flux_vis_dif, rc=rc)
129  if (esmf_logfounderror(rctocheck=rc, msg=esmf_logerr_passthru, &
130  line=__line__, &
131  file=__file__)) &
132  return ! bail out
133 
134  ! -------
135  ! Net longwave radiation (W/m2)
136  ! -------
137  call state_getimport(importstate, 'mean_net_lw_flx', &
138  isc, iec, jsc, jec, ice_ocean_boundary%lw_flux, rc=rc)
139  if (esmf_logfounderror(rctocheck=rc, msg=esmf_logerr_passthru, &
140  line=__line__, &
141  file=__file__)) &
142  return ! bail out
143 
144  !----
145  ! zonal and meridional surface stress
146  !----
147  allocate (taux(isc:iec,jsc:jec))
148  allocate (tauy(isc:iec,jsc:jec))
149 
150  call state_getimport(importstate, 'mean_zonal_moment_flx', isc, iec, jsc, jec, taux, rc=rc)
151  if (esmf_logfounderror(rctocheck=rc, msg=esmf_logerr_passthru, &
152  line=__line__, &
153  file=__file__)) &
154  return ! bail out
155  call state_getimport(importstate, 'mean_merid_moment_flx', isc, iec, jsc, jec, tauy, rc=rc)
156  if (esmf_logfounderror(rctocheck=rc, msg=esmf_logerr_passthru, &
157  line=__line__, &
158  file=__file__)) &
159  return ! bail out
160 
161 
162  ! rotate taux and tauy from true zonal/meridional to local coordinates
163  do j = jsc, jec
164  jg = j + ocean_grid%jsc - jsc
165  do i = isc, iec
166  ig = i + ocean_grid%isc - isc
167  ice_ocean_boundary%u_flux(i,j) = ocean_grid%cos_rot(ig,jg)*taux(i,j) &
168  - ocean_grid%sin_rot(ig,jg)*tauy(i,j)
169  ice_ocean_boundary%v_flux(i,j) = ocean_grid%cos_rot(ig,jg)*tauy(i,j) &
170  + ocean_grid%sin_rot(ig,jg)*taux(i,j)
171  enddo
172  enddo
173 
174  deallocate(taux, tauy)
175 
176  !----
177  ! sensible heat flux (W/m2)
178  !----
179  call state_getimport(importstate, 'mean_sensi_heat_flx', &
180  isc, iec, jsc, jec, ice_ocean_boundary%t_flux, rc=rc)
181  if (esmf_logfounderror(rctocheck=rc, msg=esmf_logerr_passthru, &
182  line=__line__, &
183  file=__file__)) &
184  return ! bail out
185 
186  !----
187  ! evaporation flux (W/m2)
188  !----
189  call state_getimport(importstate, 'mean_evap_rate', &
190  isc, iec, jsc, jec, ice_ocean_boundary%q_flux, rc=rc)
191  if (esmf_logfounderror(rctocheck=rc, msg=esmf_logerr_passthru, &
192  line=__line__, &
193  file=__file__)) &
194  return ! bail out
195 
196  !----
197  ! liquid precipitation (rain)
198  !----
199  call state_getimport(importstate, 'mean_prec_rate', &
200  isc, iec, jsc, jec, ice_ocean_boundary%lprec, rc=rc)
201  if (esmf_logfounderror(rctocheck=rc, msg=esmf_logerr_passthru, &
202  line=__line__, &
203  file=__file__)) &
204  return ! bail out
205 
206  !----
207  ! frozen precipitation (snow)
208  !----
209  call state_getimport(importstate, 'mean_fprec_rate', &
210  isc, iec, jsc, jec, ice_ocean_boundary%fprec, rc=rc)
211  if (esmf_logfounderror(rctocheck=rc, msg=esmf_logerr_passthru, &
212  line=__line__, &
213  file=__file__)) &
214  return ! bail out
215 
216  !----
217  ! mass and heat content of liquid and frozen runoff
218  !----
219  ! Note - preset values to 0, if field does not exist in importState, then will simply return
220  ! and preset value will be used
221 
222  ! liquid runoff
223  ice_ocean_boundary%lrunoff (:,:) = 0._esmf_kind_r8
224  call state_getimport(importstate, 'Foxx_rofl', &
225  isc, iec, jsc, jec, ice_ocean_boundary%lrunoff,rc=rc)
226  if (esmf_logfounderror(rctocheck=rc, msg=esmf_logerr_passthru, &
227  line=__line__, &
228  file=__file__)) &
229  return ! bail out
230 
231  ! ice runoff
232  ice_ocean_boundary%frunoff (:,:) = 0._esmf_kind_r8
233  call state_getimport(importstate, 'Foxx_rofi', &
234  isc, iec, jsc, jec, ice_ocean_boundary%frunoff,rc=rc)
235  if (esmf_logfounderror(rctocheck=rc, msg=esmf_logerr_passthru, &
236  line=__line__, &
237  file=__file__)) &
238  return ! bail out
239 
240  ! heat content of lrunoff
241  ice_ocean_boundary%lrunoff_hflx(:,:) = 0._esmf_kind_r8
242  call state_getimport(importstate, 'mean_runoff_heat_flx', &
243  isc, iec, jsc, jec, ice_ocean_boundary%lrunoff_hflx, rc=rc)
244  if (esmf_logfounderror(rctocheck=rc, msg=esmf_logerr_passthru, &
245  line=__line__, &
246  file=__file__)) &
247  return ! bail out
248 
249  ! heat content of frunoff
250  ice_ocean_boundary%frunoff_hflx(:,:) = 0._esmf_kind_r8
251  call state_getimport(importstate, 'mean_calving_heat_flx', &
252  isc, iec, jsc, jec, ice_ocean_boundary%frunoff_hflx, rc=rc)
253  if (esmf_logfounderror(rctocheck=rc, msg=esmf_logerr_passthru, &
254  line=__line__, &
255  file=__file__)) &
256  return ! bail out
257 
258  !----
259  ! salt flux from ice
260  !----
261  ice_ocean_boundary%salt_flux(:,:) = 0._esmf_kind_r8
262  call state_getimport(importstate, 'mean_salt_rate', &
263  isc, iec, jsc, jec, ice_ocean_boundary%salt_flux,rc=rc)
264  if (esmf_logfounderror(rctocheck=rc, msg=esmf_logerr_passthru, &
265  line=__line__, &
266  file=__file__)) &
267  return ! bail out
268 
269  ! !----
270  ! ! snow&ice melt heat flux (W/m^2)
271  ! !----
272  ice_ocean_boundary%seaice_melt_heat(:,:) = 0._esmf_kind_r8
273  call state_getimport(importstate, 'net_heat_flx_to_ocn', &
274  isc, iec, jsc, jec, ice_ocean_boundary%seaice_melt_heat,rc=rc)
275  if (esmf_logfounderror(rctocheck=rc, msg=esmf_logerr_passthru, &
276  line=__line__, &
277  file=__file__)) &
278  return ! bail out
279 
280  ! !----
281  ! ! snow&ice melt water flux (W/m^2)
282  ! !----
283  ice_ocean_boundary%seaice_melt(:,:) = 0._esmf_kind_r8
284  call state_getimport(importstate, 'mean_fresh_water_to_ocean_rate', &
285  isc, iec, jsc, jec, ice_ocean_boundary%seaice_melt,rc=rc)
286  if (esmf_logfounderror(rctocheck=rc, msg=esmf_logerr_passthru, &
287  line=__line__, &
288  file=__file__)) &
289  return ! bail out
290 
291  !----
292  ! mass of overlying ice
293  !----
294  ! Note - preset values to 0, if field does not exist in importState, then will simply return
295  ! and preset value will be used
296 
297  ice_ocean_boundary%mi(:,:) = 0._esmf_kind_r8
298  call state_getimport(importstate, 'mass_of_overlying_ice', &
299  isc, iec, jsc, jec, ice_ocean_boundary%mi, rc=rc)
300  if (esmf_logfounderror(rctocheck=rc, msg=esmf_logerr_passthru, &
301  line=__line__, &
302  file=__file__)) &
303  return ! bail out
304 
305 end subroutine mom_import
306 
307 !> Maps outgoing ocean data to ESMF State
308 subroutine mom_export(ocean_public, ocean_grid, ocean_state, exportState, clock, rc)
309  type(ocean_public_type) , intent(in) :: ocean_public !< Ocean surface state
310  type(ocean_grid_type) , intent(in) :: ocean_grid !< Ocean model grid
311  type(ocean_state_type) , pointer :: ocean_state !< Ocean state
312  type(esmf_state) , intent(inout) :: exportstate !< outgoing data
313  type(esmf_clock) , intent(in) :: clock !< ESMF clock
314  integer , intent(inout) :: rc !< Return code
315 
316  ! Local variables
317  integer :: i, j, ig, jg ! indices
318  integer :: isc, iec, jsc, jec ! indices
319  integer :: iloc, jloc ! indices
320  integer :: iglob, jglob ! indices
321  integer :: n
322  integer :: icount
323  real :: slp_l, slp_r, slp_c
324  real :: slope, u_min, u_max
325  integer :: day, secs
326  type(esmf_timeinterval) :: timestep
327  integer :: dt_int
328  real :: inv_dt_int !< The inverse of coupling time interval in s-1.
329  type(esmf_stateitem_flag) :: itemflag
330  real(esmf_kind_r8), allocatable :: omask(:,:)
331  real(esmf_kind_r8), allocatable :: melt_potential(:,:)
332  real(esmf_kind_r8), allocatable :: ocz(:,:), ocm(:,:)
333  real(esmf_kind_r8), allocatable :: ocz_rot(:,:), ocm_rot(:,:)
334  real(esmf_kind_r8), allocatable :: ssh(:,:)
335  real(esmf_kind_r8), allocatable :: dhdx(:,:), dhdy(:,:)
336  real(esmf_kind_r8), allocatable :: dhdx_rot(:,:), dhdy_rot(:,:)
337  character(len=*) , parameter :: subname = '(mom_export)'
338 
339  rc = esmf_success
340 
341  call esmf_clockget( clock, timestep=timestep, rc=rc)
342  if (esmf_logfounderror(rctocheck=rc, msg=esmf_logerr_passthru, &
343  line=__line__, &
344  file=__file__)) &
345  return ! bail out
346 
347  call esmf_timeintervalget( timestep, s=dt_int, rc=rc )
348  if (esmf_logfounderror(rctocheck=rc, msg=esmf_logerr_passthru, &
349  line=__line__, &
350  file=__file__)) &
351  return ! bail out
352 
353  ! Use Adcroft's rule of reciprocals; it does the right thing here.
354  if (real(dt_int) > 0.0) then
355  inv_dt_int = 1.0 / real(dt_int)
356  else
357  inv_dt_int = 0.0
358  endif
359 
360  !----------------
361  ! Copy from ocean_public to exportstate.
362  !----------------
363 
364  call mpp_get_compute_domain(ocean_public%domain, isc, iec, jsc, jec)
365 
366  ! -------
367  ! ocean mask
368  ! -------
369 
370  allocate(omask(isc:iec, jsc:jec))
371  do j = jsc, jec
372  jg = j + ocean_grid%jsc - jsc
373  do i = isc, iec
374  ig = i + ocean_grid%isc - isc
375  omask(i,j) = nint(ocean_grid%mask2dT(ig,jg))
376  enddo
377  enddo
378 
379  call state_setexport(exportstate, 'ocean_mask', &
380  isc, iec, jsc, jec, omask, ocean_grid, rc=rc)
381  if (esmf_logfounderror(rctocheck=rc, msg=esmf_logerr_passthru, &
382  line=__line__, &
383  file=__file__)) &
384  return ! bail out
385 
386  deallocate(omask)
387 
388  ! -------
389  ! Sea surface temperature
390  ! -------
391  call state_setexport(exportstate, 'sea_surface_temperature', &
392  isc, iec, jsc, jec, ocean_public%t_surf, ocean_grid, rc=rc)
393  if (esmf_logfounderror(rctocheck=rc, msg=esmf_logerr_passthru, &
394  line=__line__, &
395  file=__file__)) &
396  return ! bail out
397 
398  ! -------
399  ! Sea surface salinity
400  ! -------
401  call state_setexport(exportstate, 's_surf', &
402  isc, iec, jsc, jec, ocean_public%s_surf, ocean_grid, rc=rc)
403  if (esmf_logfounderror(rctocheck=rc, msg=esmf_logerr_passthru, &
404  line=__line__, &
405  file=__file__)) &
406  return ! bail out
407 
408  ! -------
409  ! zonal and meridional currents
410  ! -------
411 
412  ! rotate ocn current from tripolar grid back to lat/lon grid x,y => latlon (CCW)
413  ! "ocean_grid" has halos and uses local indexing.
414 
415  allocate(ocz(isc:iec, jsc:jec))
416  allocate(ocm(isc:iec, jsc:jec))
417  allocate(ocz_rot(isc:iec, jsc:jec))
418  allocate(ocm_rot(isc:iec, jsc:jec))
419 
420  do j = jsc, jec
421  jg = j + ocean_grid%jsc - jsc
422  do i = isc, iec
423  ig = i + ocean_grid%isc - isc
424  ocz(i,j) = ocean_public%u_surf(i,j)
425  ocm(i,j) = ocean_public%v_surf(i,j)
426  ocz_rot(i,j) = ocean_grid%cos_rot(ig,jg)*ocz(i,j) + ocean_grid%sin_rot(ig,jg)*ocm(i,j)
427  ocm_rot(i,j) = ocean_grid%cos_rot(ig,jg)*ocm(i,j) - ocean_grid%sin_rot(ig,jg)*ocz(i,j)
428  enddo
429  enddo
430 
431  call state_setexport(exportstate, 'ocn_current_zonal', &
432  isc, iec, jsc, jec, ocz_rot, ocean_grid, rc=rc)
433  if (esmf_logfounderror(rctocheck=rc, msg=esmf_logerr_passthru, &
434  line=__line__, &
435  file=__file__)) &
436  return ! bail out
437 
438  call state_setexport(exportstate, 'ocn_current_merid', &
439  isc, iec, jsc, jec, ocm_rot, ocean_grid, rc=rc)
440  if (esmf_logfounderror(rctocheck=rc, msg=esmf_logerr_passthru, &
441  line=__line__, &
442  file=__file__)) &
443  return ! bail out
444 
445  deallocate(ocz, ocm, ocz_rot, ocm_rot)
446 
447  ! -------
448  ! Boundary layer depth
449  ! -------
450  call esmf_stateget(exportstate, 'So_bldepth', itemflag, rc=rc)
451  if (itemflag /= esmf_stateitem_notfound) then
452  call state_setexport(exportstate, 'So_bldepth', &
453  isc, iec, jsc, jec, ocean_public%obld, ocean_grid, rc=rc)
454  if (esmf_logfounderror(rctocheck=rc, msg=esmf_logerr_passthru, &
455  line=__line__, &
456  file=__file__)) &
457  return ! bail out
458  endif
459 
460  ! -------
461  ! Freezing melting potential
462  ! -------
463  ! melt_potential, defined positive for T>Tfreeze, so need to change sign
464  ! Convert from J/m^2 to W/m^2 and make sure Melt_potential is always <= 0
465 
466  allocate(melt_potential(isc:iec, jsc:jec))
467 
468  do j = jsc,jec
469  do i = isc,iec
470  if (ocean_public%frazil(i,j) > 0.0) then
471  melt_potential(i,j) = ocean_public%frazil(i,j) * inv_dt_int
472  else
473  melt_potential(i,j) = -ocean_public%melt_potential(i,j) * inv_dt_int
474  if (melt_potential(i,j) > 0.0) melt_potential(i,j) = 0.0
475  endif
476  enddo
477  enddo
478 
479  call state_setexport(exportstate, 'freezing_melting_potential', &
480  isc, iec, jsc, jec, melt_potential, ocean_grid, rc=rc)
481  if (esmf_logfounderror(rctocheck=rc, msg=esmf_logerr_passthru, &
482  line=__line__, &
483  file=__file__)) &
484  return ! bail out
485 
486  deallocate(melt_potential)
487 
488  ! -------
489  ! Sea level
490  ! -------
491  call esmf_stateget(exportstate, 'sea_level', itemflag, rc=rc)
492  if (itemflag /= esmf_stateitem_notfound) then
493  call state_setexport(exportstate, 'sea_level', &
494  isc, iec, jsc, jec, ocean_public%sea_lev, ocean_grid, rc=rc)
495  if (esmf_logfounderror(rctocheck=rc, msg=esmf_logerr_passthru, &
496  line=__line__, &
497  file=__file__)) &
498  return ! bail out
499  endif
500 
501  !----------------
502  ! Sea-surface zonal and meridional slopes
503  !----------------
504 
505  allocate(ssh(ocean_grid%isd:ocean_grid%ied,ocean_grid%jsd:ocean_grid%jed)) ! local indices with halos
506  allocate(dhdx(isc:iec, jsc:jec)) !global indices without halos
507  allocate(dhdy(isc:iec, jsc:jec)) !global indices without halos
508  allocate(dhdx_rot(isc:iec, jsc:jec)) !global indices without halos
509  allocate(dhdy_rot(isc:iec, jsc:jec)) !global indices without halos
510 
511  ssh = 0.0_esmf_kind_r8
512  dhdx = 0.0_esmf_kind_r8
513  dhdy = 0.0_esmf_kind_r8
514 
515  ! Make a copy of ssh in order to do a halo update (ssh has local indexing with halos)
516  do j = ocean_grid%jsc, ocean_grid%jec
517  jloc = j + ocean_grid%jdg_offset
518  do i = ocean_grid%isc,ocean_grid%iec
519  iloc = i + ocean_grid%idg_offset
520  ssh(i,j) = ocean_public%sea_lev(iloc,jloc)
521  enddo
522  enddo
523 
524  ! Update halo of ssh so we can calculate gradients (local indexing)
525  call pass_var(ssh, ocean_grid%domain)
526 
527  ! d/dx ssh
528  ! This is a simple second-order difference
529  ! dhdx(i,j) = 0.5 * (ssh(i+1,j) - ssh(i-1,j)) * ocean_grid%US%m_to_L*ocean_grid%IdxT(i,j) * ocean_grid%mask2dT(ig,jg)
530 
531  do jglob = jsc, jec
532  j = jglob + ocean_grid%jsc - jsc
533  do iglob = isc,iec
534  i = iglob + ocean_grid%isc - isc
535  ! This is a PLM slope which might be less prone to the A-grid null mode
536  slp_l = (ssh(i,j) - ssh(i-1,j)) * ocean_grid%mask2dCu(i-1,j)
537  if (ocean_grid%mask2dCu(i-1,j)==0.) slp_l = 0.
538  slp_r = (ssh(i+1,j) - ssh(i,j)) * ocean_grid%mask2dCu(i,j)
539  if (ocean_grid%mask2dCu(i+1,j)==0.) slp_r = 0.
540  slp_c = 0.5 * (slp_l + slp_r)
541  if ( (slp_l * slp_r) > 0.0 ) then
542  ! This limits the slope so that the edge values are bounded by the
543  ! two cell averages spanning the edge.
544  u_min = min( ssh(i-1,j), ssh(i,j), ssh(i+1,j) )
545  u_max = max( ssh(i-1,j), ssh(i,j), ssh(i+1,j) )
546  slope = sign( min( abs(slp_c), 2.*min( ssh(i,j) - u_min, u_max - ssh(i,j) ) ), slp_c )
547  else
548  ! Extrema in the mean values require a PCM reconstruction avoid generating
549  ! larger extreme values.
550  slope = 0.0
551  endif
552  dhdx(iglob,jglob) = slope * ocean_grid%US%m_to_L*ocean_grid%IdxT(i,j) * ocean_grid%mask2dT(i,j)
553  if (ocean_grid%mask2dT(i,j)==0.) dhdx(iglob,jglob) = 0.0
554  enddo
555  enddo
556 
557  ! d/dy ssh
558  ! This is a simple second-order difference
559  ! dhdy(i,j) = 0.5 * (ssh(i,j+1) - ssh(i,j-1)) * ocean_grid%US%m_to_L*ocean_grid%IdyT(i,j) * ocean_grid%mask2dT(ig,jg)
560 
561  do jglob = jsc, jec
562  j = jglob + ocean_grid%jsc - jsc
563  do iglob = isc,iec
564  i = iglob + ocean_grid%isc - isc
565  ! This is a PLM slope which might be less prone to the A-ocean_grid null mode
566  slp_l = ssh(i,j) - ssh(i,j-1) * ocean_grid%mask2dCv(i,j-1)
567  if (ocean_grid%mask2dCv(i,j-1)==0.) slp_l = 0.
568  slp_r = ssh(i,j+1) - ssh(i,j) * ocean_grid%mask2dCv(i,j)
569  if (ocean_grid%mask2dCv(i,j+1)==0.) slp_r = 0.
570  slp_c = 0.5 * (slp_l + slp_r)
571  if ((slp_l * slp_r) > 0.0) then
572  ! This limits the slope so that the edge values are bounded by the
573  ! two cell averages spanning the edge.
574  u_min = min( ssh(i,j-1), ssh(i,j), ssh(i,j+1) )
575  u_max = max( ssh(i,j-1), ssh(i,j), ssh(i,j+1) )
576  slope = sign( min( abs(slp_c), 2.*min( ssh(i,j) - u_min, u_max - ssh(i,j) ) ), slp_c )
577  else
578  ! Extrema in the mean values require a PCM reconstruction avoid generating
579  ! larger extreme values.
580  slope = 0.0
581  endif
582  dhdy(iglob,jglob) = slope * ocean_grid%US%m_to_L*ocean_grid%IdyT(i,j) * ocean_grid%mask2dT(i,j)
583  if (ocean_grid%mask2dT(i,j)==0.) dhdy(iglob,jglob) = 0.0
584  enddo
585  enddo
586 
587  ! rotate slopes from tripolar grid back to lat/lon grid, x,y => latlon (CCW)
588  ! "ocean_grid" uses has halos and uses local indexing.
589 
590  do j = jsc, jec
591  jg = j + ocean_grid%jsc - jsc
592  do i = isc, iec
593  ig = i + ocean_grid%isc - isc
594  dhdx_rot(i,j) = ocean_grid%cos_rot(ig,jg)*dhdx(i,j) + ocean_grid%sin_rot(ig,jg)*dhdy(i,j)
595  dhdy_rot(i,j) = ocean_grid%cos_rot(ig,jg)*dhdy(i,j) - ocean_grid%sin_rot(ig,jg)*dhdx(i,j)
596  enddo
597  enddo
598 
599  call state_setexport(exportstate, 'sea_surface_slope_zonal', &
600  isc, iec, jsc, jec, dhdx_rot, ocean_grid, rc=rc)
601  if (esmf_logfounderror(rctocheck=rc, msg=esmf_logerr_passthru, &
602  line=__line__, &
603  file=__file__)) &
604  return ! bail out
605 
606  call state_setexport(exportstate, 'sea_surface_slope_merid', &
607  isc, iec, jsc, jec, dhdy_rot, ocean_grid, rc=rc)
608  if (esmf_logfounderror(rctocheck=rc, msg=esmf_logerr_passthru, &
609  line=__line__, &
610  file=__file__)) &
611  return ! bail out
612 
613  deallocate(ssh, dhdx, dhdy, dhdx_rot, dhdy_rot)
614 
615 end subroutine mom_export
616 
617 !> Get field pointer 1D
618 subroutine state_getfldptr_1d(State, fldname, fldptr, rc)
619  type(esmf_state) , intent(in) :: State !< ESMF state
620  character(len=*) , intent(in) :: fldname !< Field name
621  real(ESMF_KIND_R8), pointer , intent(in) :: fldptr(:)!< Pointer to the 1D field
622  integer, optional , intent(out) :: rc !< Return code
623 
624  ! local variables
625  type(esmf_field) :: lfield
626  integer :: lrc
627  character(len=*),parameter :: subname='(MOM_cap:State_GetFldPtr)'
628 
629  call esmf_stateget(state, itemname=trim(fldname), field=lfield, rc=lrc)
630  if (esmf_logfounderror(rctocheck=rc, msg=esmf_logerr_passthru, &
631  line=__line__, &
632  file=__file__)) &
633  return ! bail out
634  call esmf_fieldget(lfield, farrayptr=fldptr, rc=lrc)
635  if (esmf_logfounderror(rctocheck=rc, msg=esmf_logerr_passthru, &
636  line=__line__, &
637  file=__file__)) &
638  return ! bail out
639 
640  if (present(rc)) rc = lrc
641 
642 end subroutine state_getfldptr_1d
643 
644 !> Get field pointer 2D
645 subroutine state_getfldptr_2d(State, fldname, fldptr, rc)
646  type(esmf_state) , intent(in) :: State !< ESMF state
647  character(len=*) , intent(in) :: fldname !< Field name
648  real(ESMF_KIND_R8), pointer , intent(in) :: fldptr(:,:)!< Pointer to the 2D field
649  integer, optional , intent(out) :: rc !< Return code
650 
651  ! local variables
652  type(esmf_field) :: lfield
653  integer :: lrc
654  character(len=*),parameter :: subname='(MOM_cap:State_GetFldPtr)'
655 
656  call esmf_stateget(state, itemname=trim(fldname), field=lfield, rc=lrc)
657  if (esmf_logfounderror(rctocheck=rc, msg=esmf_logerr_passthru, &
658  line=__line__, &
659  file=__file__)) &
660  return ! bail out
661  call esmf_fieldget(lfield, farrayptr=fldptr, rc=lrc)
662  if (esmf_logfounderror(rctocheck=rc, msg=esmf_logerr_passthru, &
663  line=__line__, &
664  file=__file__)) &
665  return ! bail out
666 
667  if (present(rc)) rc = lrc
668 
669 end subroutine state_getfldptr_2d
670 
671 !> Map import state field to output array
672 subroutine state_getimport(state, fldname, isc, iec, jsc, jec, output, do_sum, rc)
673  type(esmf_state) , intent(in) :: state !< ESMF state
674  character(len=*) , intent(in) :: fldname !< Field name
675  integer , intent(in) :: isc !< The start i-index of cell centers within
676  !! the computational domain
677  integer , intent(in) :: iec !< The end i-index of cell centers within the
678  !! computational domain
679  integer , intent(in) :: jsc !< The start j-index of cell centers within
680  !! the computational domain
681  integer , intent(in) :: jec !< The end j-index of cell centers within
682  !! the computational domain
683  real (esmf_kind_r8) , intent(inout) :: output(isc:iec,jsc:jec)!< Output 2D array
684  logical, optional , intent(in) :: do_sum !< If true, sums the data
685  integer , intent(out) :: rc !< Return code
686 
687  ! local variables
688  type(esmf_stateitem_flag) :: itemflag
689  integer :: n, i, j, i1, j1
690  integer :: lbnd1,lbnd2
691  real(esmf_kind_r8), pointer :: dataptr1d(:)
692  real(esmf_kind_r8), pointer :: dataptr2d(:,:)
693  character(len=*) , parameter :: subname='(MOM_cap_methods:state_getimport)'
694  ! ----------------------------------------------
695 
696  rc = esmf_success
697 
698  call esmf_stateget(state, trim(fldname), itemflag, rc=rc)
699  if (itemflag /= esmf_stateitem_notfound) then
700 
701  if (geomtype == esmf_geomtype_mesh) then
702 
703  ! get field pointer
704  call state_getfldptr(state, trim(fldname), dataptr1d, rc)
705  if (esmf_logfounderror(rctocheck=rc, msg=esmf_logerr_passthru, &
706  line=__line__, &
707  file=__file__)) &
708  return ! bail out
709 
710  ! determine output array
711  n = 0
712  do j = jsc,jec
713  do i = isc,iec
714  n = n + 1
715  if (present(do_sum)) then
716  output(i,j) = output(i,j) + dataptr1d(n)
717  else
718  output(i,j) = dataptr1d(n)
719  endif
720  enddo
721  enddo
722 
723  else if (geomtype == esmf_geomtype_grid) then
724 
725  call state_getfldptr(state, trim(fldname), dataptr2d, rc)
726  if (esmf_logfounderror(rctocheck=rc, msg=esmf_logerr_passthru, &
727  line=__line__, &
728  file=__file__)) &
729  return ! bail out
730 
731  lbnd1 = lbound(dataptr2d,1)
732  lbnd2 = lbound(dataptr2d,2)
733 
734  do j = jsc, jec
735  j1 = j + lbnd2 - jsc
736  do i = isc, iec
737  i1 = i + lbnd1 - isc
738  if (present(do_sum)) then
739  output(i,j) = output(i,j) + dataptr2d(i1,j1)
740  else
741  output(i,j) = dataptr2d(i1,j1)
742  endif
743  enddo
744  enddo
745 
746  endif
747 
748  endif
749 
750 end subroutine state_getimport
751 
752 !> Map input array to export state
753 subroutine state_setexport(state, fldname, isc, iec, jsc, jec, input, ocean_grid, rc)
754  type(esmf_state) , intent(inout) :: state !< ESMF state
755  character(len=*) , intent(in) :: fldname !< Field name
756  integer , intent(in) :: isc !< The start i-index of cell centers within
757  !! the computational domain
758  integer , intent(in) :: iec !< The end i-index of cell centers within the
759  !! computational domain
760  integer , intent(in) :: jsc !< The start j-index of cell centers within
761  !! the computational domain
762  integer , intent(in) :: jec !< The end j-index of cell centers within
763  !! the computational domain
764  real (esmf_kind_r8) , intent(in) :: input(isc:iec,jsc:jec)!< Input 2D array
765  type(ocean_grid_type) , intent(in) :: ocean_grid !< Ocean horizontal grid
766  integer , intent(out) :: rc !< Return code
767 
768  ! local variables
769  type(esmf_stateitem_flag) :: itemflag
770  integer :: n, i, j, i1, j1, ig,jg
771  integer :: lbnd1,lbnd2
772  real(esmf_kind_r8), pointer :: dataptr1d(:)
773  real(esmf_kind_r8), pointer :: dataptr2d(:,:)
774  character(len=*) , parameter :: subname='(MOM_cap_methods:state_setexport)'
775  ! ----------------------------------------------
776 
777  rc = esmf_success
778 
779  ! Indexing notes:
780  ! input array from "ocean_public" uses local indexing without halos
781  ! mask from "ocean_grid" uses local indexing with halos
782 
783  call esmf_stateget(state, trim(fldname), itemflag, rc=rc)
784  if (itemflag /= esmf_stateitem_notfound) then
785 
786  if (geomtype == esmf_geomtype_mesh) then
787 
788  call state_getfldptr(state, trim(fldname), dataptr1d, rc)
789  if (esmf_logfounderror(rctocheck=rc, msg=esmf_logerr_passthru, &
790  line=__line__, &
791  file=__file__)) &
792  return ! bail out
793 
794  n = 0
795  do j = jsc, jec
796  jg = j + ocean_grid%jsc - jsc
797  do i = isc, iec
798  ig = i + ocean_grid%isc - isc
799  n = n+1
800  dataptr1d(n) = input(i,j) * ocean_grid%mask2dT(ig,jg)
801  enddo
802  enddo
803 
804  else if (geomtype == esmf_geomtype_grid) then
805 
806  call state_getfldptr(state, trim(fldname), dataptr2d, rc)
807  if (esmf_logfounderror(rctocheck=rc, msg=esmf_logerr_passthru, &
808  line=__line__, &
809  file=__file__)) &
810  return ! bail out
811 
812  lbnd1 = lbound(dataptr2d,1)
813  lbnd2 = lbound(dataptr2d,2)
814 
815  do j = jsc, jec
816  j1 = j + lbnd2 - jsc
817  jg = j + ocean_grid%jsc - jsc
818  do i = isc, iec
819  i1 = i + lbnd1 - isc
820  ig = i + ocean_grid%isc - isc
821  dataptr2d(i1,j1) = input(i,j) * ocean_grid%mask2dT(ig,jg)
822  enddo
823  enddo
824 
825  endif
826 
827  endif
828 
829 end subroutine state_setexport
830 
831 end module mom_cap_methods
mom_cap_methods::state_getfldptr
Get field pointer.
Definition: mom_cap_methods.F90:36
mom_cap_methods::import_cnt
integer import_cnt
used to skip using the import state at the first count for cesm
Definition: mom_cap_methods.F90:41
mom_cap_methods::state_setexport
subroutine, private state_setexport(state, fldname, isc, iec, jsc, jec, input, ocean_grid, rc)
Map input array to export state.
Definition: mom_cap_methods.F90:754
mom_surface_forcing_nuopc::ice_ocean_boundary_type
Structure corresponding to forcing, but with the elements, units, and conventions that exactly confor...
Definition: mom_surface_forcing_nuopc.F90:155
mom_cap_methods::geomtype
type(esmf_geomtype_flag) geomtype
SMF type describing type of geometry (mesh or grid)
Definition: mom_cap_methods.F90:43
mom_cap_methods::state_getfldptr_1d
subroutine state_getfldptr_1d(State, fldname, fldptr, rc)
Get field pointer 1D.
Definition: mom_cap_methods.F90:619
mom_surface_forcing_nuopc
Converts the input ESMF data (import data) to a MOM-specific data type (surface_forcing_CS).
Definition: mom_surface_forcing_nuopc.F90:2
mom_domains::pass_var
Do a halo update on an array.
Definition: MOM_domains.F90:49
mom_cap_methods::mom_set_geomtype
subroutine, public mom_set_geomtype(geomtype_in)
Sets module variable geometry type.
Definition: mom_cap_methods.F90:50
mom_domains
Describes the decomposed MOM domain and has routines for communications across PEs.
Definition: MOM_domains.F90:2
mom_ocean_model_nuopc::ocean_state_type
The ocean_state_type contains all information about the state of the ocean, with a format that is pri...
Definition: mom_ocean_model_nuopc.F90:133
mom_cap_methods::state_getimport
subroutine, private state_getimport(state, fldname, isc, iec, jsc, jec, output, do_sum, rc)
Map import state field to output array.
Definition: mom_cap_methods.F90:673
mom_cap_methods::mom_export
subroutine, public mom_export(ocean_public, ocean_grid, ocean_state, exportState, clock, rc)
Maps outgoing ocean data to ESMF State.
Definition: mom_cap_methods.F90:309
mom_grid
Provides the ocean grid type.
Definition: MOM_grid.F90:2
mom_cap_methods::mom_import
subroutine, public mom_import(ocean_public, ocean_grid, importState, ice_ocean_boundary, rc)
This function has a few purposes: (1) it imports surface fluxes using data from the mediator; and (2)...
Definition: mom_cap_methods.F90:61
mom_ocean_model_nuopc::ocean_public_type
This type is used for communication with other components via the FMS coupler. The element names and ...
Definition: mom_ocean_model_nuopc.F90:86
mom_cap_methods
Contains import/export methods for both NEMS and CMEPS.
Definition: mom_cap_methods.F90:2
mom_cap_methods::state_getfldptr_2d
subroutine state_getfldptr_2d(State, fldname, fldptr, rc)
Get field pointer 2D.
Definition: mom_cap_methods.F90:646
mom_grid::ocean_grid_type
Ocean grid type. See mom_grid for details.
Definition: MOM_grid.F90:26
mom_ocean_model_nuopc
Top-level module for the MOM6 ocean model in coupled mode.
Definition: mom_ocean_model_nuopc.F90:2