Skip to content

Commit ba5b5ba

Browse files
alperaltuntasalperaltuntas
authored
Merge pull request #22 from gustavo-marques/check_fluxes
Updates IOB fluxes and adds option to write restart files
2 parents 369d5d6 + f5a8483 commit ba5b5ba

2 files changed

Lines changed: 116 additions & 58 deletions

File tree

‎config_src/mct_driver/coupler_indices.F90‎

Lines changed: 58 additions & 45 deletions
Original file line numberDiff line numberDiff line change
@@ -1,19 +1,19 @@
11
module coupler_indices
22

33
! MCT types
4-
use mct_mod, only : mct_aVect
4+
use mct_mod, only : mct_aVect
55
! MCT fucntions
6-
use mct_mod, only : mct_avect_indexra, mct_aVect_init, mct_aVect_clean
7-
use seq_flds_mod, only : seq_flds_x2o_fields, seq_flds_o2x_fields
8-
use seq_flds_mod, only : seq_flds_i2o_per_cat, ice_ncat
6+
use mct_mod, only : mct_avect_indexra, mct_aVect_init, mct_aVect_clean
7+
use seq_flds_mod, only : seq_flds_x2o_fields, seq_flds_o2x_fields
8+
use seq_flds_mod, only : seq_flds_i2o_per_cat, ice_ncat
99

1010
! MOM types
11-
use MOM_grid, only : ocean_grid_type
11+
use MOM_grid, only : ocean_grid_type
1212
use MOM_surface_forcing, only: ice_ocean_boundary_type
1313
! MOM functions
14-
use MOM_domains, only : pass_var, AGRID
15-
use ocean_model_mod, only : ocean_public_type
16-
14+
use MOM_domains, only : pass_var, AGRID
15+
use ocean_model_mod, only : ocean_public_type
16+
use MOM_error_handler, only : MOM_error, FATAL
1717
implicit none
1818

1919
private
@@ -52,8 +52,8 @@ module coupler_indices
5252
integer :: x2o_Foxx_taux ! zonal wind stress (taux) (W/m2 )
5353
integer :: x2o_Foxx_tauy ! meridonal wind stress (tauy) (W/m2 )
5454
integer :: x2o_Foxx_swnet ! net short-wave heat flux (W/m2 )
55-
integer :: x2o_Foxx_sen ! sensible heat flux (W/m2 )
56-
integer :: x2o_Foxx_lat
55+
integer :: x2o_Foxx_sen ! sensible heat flux (W/m2)
56+
integer :: x2o_Foxx_lat ! latent heat flux (W/m2)
5757
integer :: x2o_Foxx_lwup ! longwave radiation (up) (W/m2 )
5858
integer :: x2o_Faxa_lwdn ! longwave radiation (down) (W/m2 )
5959
integer :: x2o_Fioi_melth ! heat flux from snow & ice melt (W/m2 )
@@ -116,6 +116,7 @@ subroutine coupler_indices_init(ind)
116116
call mct_aVect_init(x2o, rList=seq_flds_x2o_fields, lsize=1)
117117
call mct_aVect_init(o2x, rList=seq_flds_o2x_fields, lsize=1)
118118

119+
! ocean to coupler
119120
ind%o2x_So_t = mct_avect_indexra(o2x,'So_t')
120121
ind%o2x_So_u = mct_avect_indexra(o2x,'So_u')
121122
ind%o2x_So_v = mct_avect_indexra(o2x,'So_v')
@@ -127,14 +128,15 @@ subroutine coupler_indices_init(ind)
127128
ind%o2x_Fioo_q = mct_avect_indexra(o2x,'Fioo_q')
128129
ind%o2x_Faoo_fco2_ocn = mct_avect_indexra(o2x,'Faoo_fco2_ocn',perrWith='quiet')
129130
ind%o2x_Faoo_fdms_ocn = mct_avect_indexra(o2x,'Faoo_fdms_ocn',perrWith='quiet')
131+
132+
! coupler to ocean
130133
ind%x2o_Si_ifrac = mct_avect_indexra(x2o,'Si_ifrac')
131134
ind%x2o_Sa_pslv = mct_avect_indexra(x2o,'Sa_pslv')
132135
ind%x2o_So_duu10n = mct_avect_indexra(x2o,'So_duu10n')
133136
! QL, 150526, from wav
134137
ind%x2o_Sw_lamult = mct_avect_indexra(x2o,'Sw_lamult')
135138
ind%x2o_Sw_ustokes = mct_avect_indexra(x2o,'Sw_ustokes')
136139
ind%x2o_Sw_vstokes = mct_avect_indexra(x2o,'Sw_vstokes')
137-
138140
ind%x2o_Foxx_tauy = mct_avect_indexra(x2o,'Foxx_tauy')
139141
ind%x2o_Foxx_taux = mct_avect_indexra(x2o,'Foxx_taux')
140142
ind%x2o_Foxx_swnet = mct_avect_indexra(x2o,'Foxx_swnet')
@@ -226,7 +228,8 @@ subroutine ocn_export(ind, ocn_public, grid, o2x)
226228
do i=grid%isc,grid%iec
227229
n = n+1
228230
ig = i + grid%idg_offset
229-
o2x(ind%o2x_So_t, n) = ocn_public%t_surf(ig,jg) * grid%mask2dT(i,j)
231+
! surface temperature in Kelvin
232+
o2x(ind%o2x_So_t, n) = ocn_public%t_surf(ig,jg) * grid%mask2dT(i,j)
230233
o2x(ind%o2x_So_s, n) = ocn_public%s_surf(ig,jg) * grid%mask2dT(i,j)
231234
o2x(ind%o2x_So_u, n) = ocn_public%u_surf(ig,jg) * grid%mask2dT(i,j)
232235
o2x(ind%o2x_So_v, n) = ocn_public%v_surf(ig,jg) * grid%mask2dT(i,j)
@@ -292,15 +295,18 @@ subroutine ocn_export(ind, ocn_public, grid, o2x)
292295
end subroutine ocn_export
293296

294297

295-
subroutine fill_ice_ocean_bnd(ice_ocean_boundary, grid, x2o_o, ind)
296-
type(ice_ocean_boundary_type), intent(inout) :: ice_ocean_boundary !< A type for the ice ocean boundary
297-
type(ocean_grid_type), intent(in) :: grid
298-
!type(mct_aVect), intent(in) :: x2o_o
299-
real(kind=8), intent(in) :: x2o_o(:,:)
300-
type(cpl_indices), intent(inout) :: ind
298+
subroutine fill_ice_ocean_bnd(ice_ocean_boundary, grid, x2o_o, ind, sw_decomp, c1, c2, c3, c4)
299+
type(ice_ocean_boundary_type), intent(inout) :: ice_ocean_boundary !< A type for the ice ocean boundary
300+
type(ocean_grid_type), intent(in) :: grid
301+
!type(mct_aVect), intent(in) :: x2o_o
302+
real(kind=8), intent(in) :: x2o_o(:,:)
303+
type(cpl_indices), intent(inout) :: ind
304+
logical, intent(in) :: sw_decomp !< controls if shortwave is decomposed
305+
!! into four components
306+
real(kind=8), intent(in), optional :: c1, c2, c3, c4 !< Coeffs. used in the shortwave decomposition
301307

302308
! local variables
303-
integer :: i, j, k, ig, jg
309+
integer :: i, j, k, ig, jg !< grid indices
304310

305311
! variable that are not in ice_ocean_boundary:
306312
! latent (x2o_Foxx_lat)
@@ -324,45 +330,52 @@ subroutine fill_ice_ocean_bnd(ice_ocean_boundary, grid, x2o_o, ind)
324330

325331
! need wind_stress_multiplier?
326332

327-
! Copy from x2o to ice_ocean_boundary. ice_ocean_boundary uses global indexing with no halos.
328-
write(*,*) 'max. k is:', (grid%jec-grid%jsc+1) * (grid%iec-grid%isc+1)
329-
! zonal wind stress (taux)
330-
write(*,*) 'taux', SIZE(x2o_o(ind%x2o_Foxx_taux,:))
331-
write(*,*) 'ice_ocean_boundary%u_flux', SIZE(ice_ocean_boundary%u_flux(:,:))
332333
k = 0
333334
do j = grid%jsc, grid%jec
334335
jg = j + grid%jdg_offset
335336
do i = grid%isc, grid%iec
336337
k = k + 1 ! Increment position within gindex
337338
ig = i + grid%idg_offset
339+
! sea-level pressure (Pa)
340+
ice_ocean_boundary%p(ig,jg) = x2o_o(ind%x2o_Sa_pslv,k)
338341
! zonal wind stress (taux)
339-
ice_ocean_boundary%u_flux(ig,jg) = 0.0 ! x20_o(ind%x2o_Foxx_taux,k)
342+
ice_ocean_boundary%u_flux(ig,jg) = x2o_o(ind%x2o_Foxx_taux,k)
340343
! meridional wind stress (tauy)
341-
ice_ocean_boundary%v_flux(ig,jg) = 0.0 ! x20_o(ind%x2o_Foxx_tauy,k)
342-
! sensible heat flux
343-
ice_ocean_boundary%t_flux(ig,jg) = 0.0 ! x20_o(ind%x2o_Foxx_sen,k)
344+
ice_ocean_boundary%v_flux(ig,jg) = x2o_o(ind%x2o_Foxx_tauy,k)
345+
! sensible heat flux (W/m2)
346+
ice_ocean_boundary%t_flux(ig,jg) = -x2o_o(ind%x2o_Foxx_sen,k)
344347
! salt flux
345-
ice_ocean_boundary%salt_flux(ig,jg) = 0.0 ! x20_o(ind%x2o_Fioi_salt,k)
346-
! heat flux from snow & ice melt
347-
ice_ocean_boundary%calving_hflx(ig,jg) = 0.0 ! x20_o(ind%x2o_Fioi_melth,k)
348+
ice_ocean_boundary%salt_flux(ig,jg) = -x2o_o(ind%x2o_Fioi_salt,k)
349+
! heat content from frozen runoff
350+
ice_ocean_boundary%calving_hflx(ig,jg) = 0.0
348351
! snow melt flux
349-
ice_ocean_boundary%fprec(ig,jg) = 0.0 ! x20_o(ind%x2o_Fioi_meltw,k)
352+
!ice_ocean_boundary%fprec(ig,jg) = x2o_o(ind%x2o_Fioi_meltw,k)
350353
! river runoff flux
351-
ice_ocean_boundary%runoff(ig,jg) = 0.0 ! x20_o(ind%x2o_Foxx_rofl,k)
354+
ice_ocean_boundary%runoff(ig,jg) = x2o_o(ind%x2o_Foxx_rofl,k)
352355
! ice runoff flux
353-
ice_ocean_boundary%calving(ig,jg) = 0.0 ! x20_o(ind%x2o_Foxx_rofi,k)
356+
ice_ocean_boundary%calving(ig,jg) = x2o_o(ind%x2o_Foxx_rofi,k)
354357
! liquid precipitation (rain)
355-
ice_ocean_boundary%lprec(ig,jg) = 0.0 ! x20_o(ind%x2o_Faxa_rain,k)
356-
! froze precipitation (snow)
357-
ice_ocean_boundary%fprec(ig,jg) = 0.0 ! x20_o(ind%x2o_Faxa_snow,k)
358-
!!!!!!! LONGWAVE NEEDS TO BE FIXED !!!!!!!
359-
! longwave radiation (up)
360-
ice_ocean_boundary%lw_flux(ig,jg) = 0.0 ! x20_o(k,ind%x2o_Foxx_lwup)
361-
! longwave radiation (down)
362-
ice_ocean_boundary%lw_flux(ig,jg) = 0.0 ! x20_o(k,ind%x2o_Faxa_lwdn)
363-
!!!!!!! SHORTWAVE NEEDS TO BE COMBINED !!!!!!!
364-
! net short-wave heat flux
365-
ice_ocean_boundary%u_flux(ig,jg) = 0.0 ! x20_o(k,ind%x2o_Foxx_swnet)
358+
ice_ocean_boundary%lprec(ig,jg) = x2o_o(ind%x2o_Faxa_rain,k)
359+
! frozen precipitation (snow)
360+
ice_ocean_boundary%fprec(ig,jg) = x2o_o(ind%x2o_Faxa_snow,k)
361+
! evaporation flux, MOM6 calls q_flux specific humidity (kg/m2/s)
362+
ice_ocean_boundary%q_flux(ig,jg) = -x2o_o(ind%x2o_Foxx_evap,k)
363+
! longwave radiation, sum up and down (W/m2)
364+
ice_ocean_boundary%lw_flux(ig,jg) = x2o_o(ind%x2o_Faxa_lwdn,k) + x2o_o(ind%x2o_Foxx_lwup,k)
365+
if (sw_decomp) then
366+
! Use runtime coefficients to decompose net short-wave heat flux into 4 components
367+
! 1) visible, direct shortwave (W/m2)
368+
ice_ocean_boundary%sw_flux_vis_dir(ig,jg) = x2o_o(ind%x2o_Foxx_swnet,k)*c1
369+
! 2) visible, diffuse shortwave (W/m2)
370+
ice_ocean_boundary%sw_flux_vis_dif(ig,jg) = x2o_o(ind%x2o_Foxx_swnet,k)*c2
371+
! 3) near-IR, direct shortwave (W/m2)
372+
ice_ocean_boundary%sw_flux_nir_dir(ig,jg) = x2o_o(ind%x2o_Foxx_swnet,k)*c3
373+
! 4) near-IR, diffuse shortwave (W/m2)
374+
ice_ocean_boundary%sw_flux_nir_dif(ig,jg) = x2o_o(ind%x2o_Foxx_swnet,k)*c4
375+
else
376+
call MOM_error(FATAL,"fill_ice_ocean_bnd: this option has not been implemented yet."// &
377+
"Shortwave must be decomposed using coeffs. c1, c2, c3, c4.");
378+
endif
366379
enddo
367380
enddo
368381

‎config_src/mct_driver/ocn_comp_mct.F90‎

Lines changed: 58 additions & 13 deletions
Original file line numberDiff line numberDiff line change
@@ -26,17 +26,18 @@ module ocn_comp_mct
2626

2727
! MOM6 modules
2828
use ocean_model_mod, only: ocean_state_type, ocean_public_type, ocean_model_init_sfc
29-
use ocean_model_mod, only: ocean_model_init, get_state_pointers
29+
use ocean_model_mod, only: ocean_model_init, get_state_pointers, ocean_model_restart
3030
use ocean_model_mod, only: ice_ocean_boundary_type, update_ocean_model
3131
use MOM_domains, only: MOM_infra_init, num_pes, root_pe, pe_here
3232
use MOM_grid, only: ocean_grid_type, get_global_grid_size
3333
use MOM_variables, only: surface
3434
use MOM_error_handler, only: MOM_error, FATAL, is_root_pe
3535
use MOM_time_manager, only: time_type, set_date, set_time, set_calendar_type, NOLEAP
36+
use MOM_file_parser, only: get_param, log_version, param_file_type
37+
use MOM_get_input, only: Get_MOM_Input, directories
3638
use coupler_indices, only: coupler_indices_init, cpl_indices
3739
use coupler_indices, only: ocn_export, fill_ice_ocean_bnd
3840

39-
4041
! By default make data private
4142
implicit none; private
4243
! Public member functions
@@ -54,8 +55,10 @@ module ocn_comp_mct
5455
type(surface), pointer :: ocn_surface => NULL() !< The ocean surface state
5556
type(ice_ocean_boundary_type) :: ice_ocean_boundary !< The ice ocean boundary type
5657
type(seq_infodata_type), pointer :: infodata !< The input info type
57-
type(cpl_indices), public :: ind !< Variable IDs
58-
58+
type(cpl_indices), public :: ind !< Variable IDs
59+
! runtime params
60+
logical :: sw_decomp !< Controls whether shortwave is decomposed into four components
61+
real :: c1, c2, c3, c4 !< Coeffs. used in the shortwave decomposition
5962
end type MCT_MOM_Data
6063

6164
type(MCT_MOM_Data) :: glb !< global structure
@@ -88,6 +91,11 @@ subroutine ocn_init_mct( EClock, cdata_o, x2o_o, o2x_o, NLFilename )
8891
logical :: ldiag_cpl = .false.
8992
integer :: isc, iec, jsc, jec, ni, nj !< Indices for the start and end of the domain
9093
!! in the x and y dir., respectively.
94+
! runi-time params
95+
type(param_file_type) :: param_file !< A structure to parse for run-time parameters
96+
type(directories) :: dirs_tmp !< A structure containing several relevant directory paths
97+
character(len=40) :: mdl = "ocn_comp_mct" !< This module's name.
98+
9199
! mct variables (these are local for now)
92100
integer :: MOM_MCT_ID
93101
type(mct_gsMap), pointer :: MOM_MCT_gsMap => NULL() !< 2d, points to cdata
@@ -194,6 +202,34 @@ subroutine ocn_init_mct( EClock, cdata_o, x2o_o, o2x_o, NLFilename )
194202
glb%ocn_public%pelist(:) = (/(i,i=pe0,pe0+npes)/)
195203
! \todo Set other bits of glb$ocn_public
196204

205+
! This include declares and sets the variable "version".
206+
! read useful runtime params
207+
call get_MOM_Input(param_file, dirs_tmp, check_params=.false.)
208+
!call log_version(param_file, mdl, version, "")
209+
call get_param(param_file, mdl, "SW_DECOMP", glb%sw_decomp, &
210+
"If True, read coeffs c1, c2, c3 and c4 and decompose" // &
211+
"the net shortwave radiation (SW) into four components:\n" // &
212+
"visible, direct shortwave = c1 * SW \n" // &
213+
"visible, diffuse shortwave = c2 * SW \n" // &
214+
"near-IR, direct shortwave = c3 * SW \n" // &
215+
"near-IR, diffuse shortwave = c4 * SW", default=.true.)
216+
if (glb%sw_decomp) then
217+
call get_param(param_file, mdl, "SW_c1", glb%c1, &
218+
"Coeff. used to convert net shortwave rad. into \n"//&
219+
"visible, direct shortwave.", units="nondim", default=0.285)
220+
call get_param(param_file, mdl, "SW_c2", glb%c2, &
221+
"Coeff. used to convert net shortwave rad. into \n"//&
222+
"visible, diffuse shortwave.", units="nondim", default=0.285)
223+
call get_param(param_file, mdl, "SW_c3", glb%c3, &
224+
"Coeff. used to convert net shortwave rad. into \n"//&
225+
"near-IR, direct shortwave.", units="nondim", default=0.215)
226+
call get_param(param_file, mdl, "SW_c4", glb%c4, &
227+
"Coeff. used to convert net shortwave rad. into \n"//&
228+
"near-IR, diffuse shortwave.", units="nondim", default=0.215)
229+
else
230+
glb%c1 = 0.0; glb%c2 = 0.0; glb%c3 = 0.0; glb%c4 = 0.0
231+
endif
232+
197233
! Initialize the MOM6 model
198234
call ocean_model_init(glb%ocn_public, glb%ocn_state, time_init, time_in)
199235

@@ -250,8 +286,8 @@ subroutine ocn_init_mct( EClock, cdata_o, x2o_o, o2x_o, NLFilename )
250286
! \todo Need interface to get dt from MOM6
251287
ncouple_per_day = seconds_in_day / ocn_cpl_dt
252288
mom_cpl_dt = seconds_in_day / ncouple_per_day
253-
if (mom_cpl_dt /= ocn_cpl_dt) then ! \todo this check is trivial for now.
254-
write(*,*) 'ERROR pop_cpl_dt and ocn_cpl_dt must be identical'
289+
if (mom_cpl_dt /= ocn_cpl_dt) then
290+
write(*,*) 'ERROR mom_cpl_dt and ocn_cpl_dt must be identical'
255291
call exit(0)
256292
end if
257293

@@ -325,6 +361,7 @@ subroutine ocn_run_mct( EClock, cdata_o, x2o_o, o2x_o)
325361
type(time_type) :: time_start !< Start of coupled time interval to pass to MOM6
326362
type(time_type) :: coupling_timestep !< Coupled time interval to pass to MOM6
327363
character(len=128) :: err_msg
364+
character(len=32) :: timestamp !< Name of intermediate restart file
328365

329366
! Compute the time at the start of this coupling interval
330367
call ESMF_ClockGet(EClock, PrevTime=time_start_ESMF, rc=rc)
@@ -351,26 +388,34 @@ subroutine ocn_run_mct( EClock, cdata_o, x2o_o, o2x_o)
351388
endif
352389

353390
! Translate the coupling time interval
354-
call ESMF_ClockGet(EClock, TimeStep=ocn_cpl_interval, rc=rc)
391+
call ESMF_ClockGet(EClock, TimeStep=ocn_cpl_interval, rc=rc)
355392
call ESMF_TimeIntervalGet(ocn_cpl_interval, yy=year, mm=month, d=day, s=seconds, sn=seconds_n, sd=seconds_d, rc=rc)
356393
coupling_timestep = set_time(seconds, days=day, err_msg=err_msg)
357394

358395
! set (actually, get from mct) the cdata pointers:
359396
! \todo this was done in _init_, is it needed again. Does this infodata need to be in glb%?
360397
call seq_cdata_setptrs(cdata_o, infodata=glb%infodata)
361398

362-
! Check alarms for flag to write restart at end of day
363-
write_restart_at_eod = seq_timemgr_RestartAlarmIsOn(EClock)
364-
! \todo Let MOM6 know to write restart...
365-
if (debug .and. is_root_pe()) write(6,*) 'ocn_run_mct, write_restart_at_eod=', write_restart_at_eod
366-
367399
! fill ice ocean boundary
368-
call fill_ice_ocean_bnd(glb%ice_ocean_boundary, glb%grid, x2o_o%rattr, glb%ind)
400+
call fill_ice_ocean_bnd(glb%ice_ocean_boundary, glb%grid, x2o_o%rattr, glb%ind, glb%sw_decomp, &
401+
glb%c1, glb%c2, glb%c3, glb%c4)
369402
if (debug .and. is_root_pe()) write(6,*) 'fill_ice_ocean_bnd'
370403

371404
call update_ocean_model(glb%ice_ocean_boundary, glb%ocn_state, glb%ocn_public, &
372405
time_start, coupling_timestep)
373406

407+
!--- write out intermediate restart file when needed.
408+
! Check alarms for flag to write restart at end of day
409+
write_restart_at_eod = seq_timemgr_RestartAlarmIsOn(EClock)
410+
if (debug .and. is_root_pe()) write(6,*) 'ocn_run_mct, write_restart_at_eod=', write_restart_at_eod
411+
412+
if (write_restart_at_eod) then
413+
!timestamp = date_to_string(EClock)
414+
! \todo add time stamp to ocean_model_restart
415+
!call ocean_model_restart(glb%ocn_state, timestamp)
416+
call ocean_model_restart(glb%ocn_state)
417+
endif
418+
374419
end subroutine ocn_run_mct
375420

376421
!> Finalizes MOM6

0 commit comments

Comments
 (0)
⚡