MOM_driver.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
5program mom6
6
7!********+*********+*********+*********+*********+*********+*********+**
8!* *
9!* The Modular Ocean Model, version 6 *
10!* MOM6 *
11!* *
12!* By Alistair Adcroft, Stephen Griffies and Robert Hallberg *
13!* *
14!* This file is the ocean-only driver for Version 6 of the Modular *
15!* Ocean Model (MOM). A separate ocean interface for use with *
16!* coupled models is provided in ocean_model_MOM.F90. These two *
17!* drivers are kept in separate directories for convenience of code *
18!* selection during compiling. This file orchestrates the calls to *
19!* the MOM initialization routines, to the subroutine that steps *
20!* the model, and coordinates the output and saving restarts. A *
21!* description of all of the files that constitute MOM is found in *
22!* the comments at the beginning of MOM.F90. The arguments of each *
23!* subroutine are described where the subroutine is defined. *
24!* *
25!* Macros written all in capital letters are defined in MOM_memory.h. *
26!* *
27!********+*********+*********+*********+*********+*********+*********+**
28
29 use mom_cpu_clock, only : cpu_clock_id, cpu_clock_begin, cpu_clock_end
30 use mom_cpu_clock, only : clock_component
31 use mom_data_override, only : data_override_init
33 use mom_diag_manager_infra, only : diag_manager_set_time_end_infra
37 use mom, only : step_offline
38 use mom, only : save_mom_restart
39 use mom_coms, only : set_pelist
40 use mom_domains, only : mom_infra_init, mom_infra_end, set_mom_thread_affinity
41 use mom_ensemble_manager, only : ensemble_manager_init, get_ensemble_size
42 use mom_ensemble_manager, only : ensemble_pelist_setup
43 use mom_error_handler, only : mom_error, mom_mesg, warning, fatal, is_root_pe
45 use mom_file_parser, only : read_param, get_param, log_param, log_version, param_file_type
50 use mom_grid, only : ocean_grid_type
56 use mom_interpolate, only : time_interp_external_init
57 use mom_io, only : file_exists, open_ascii_file, close_file
58 use mom_io, only : check_nml_error, io_infra_init, io_infra_end
59 use mom_io, only : append_file, readonly_file
63 use mom_time_manager, only : time_type, set_date, get_date, real_to_time, time_to_real
64 use mom_time_manager, only : operator(+), operator(-), operator(*), operator(/)
65 use mom_time_manager, only : operator(>), operator(<), operator(>=)
66 use mom_time_manager, only : increment_date, set_calendar_type, month_name
67 use mom_time_manager, only : julian, gregorian, noleap, thirty_day_months, no_calendar
70 use mom_variables, only : surface
76
77 implicit none
78
79#include <MOM_memory.h>
80
81 ! A structure with the driving mechanical surface forces
82 type(mech_forcing) :: forces
83 ! A structure containing pointers to the thermodynamic forcing fields
84 ! at the ocean surface.
85 type(forcing) :: fluxes
86 ! A structure containing pointers to the ocean surface state fields.
87 type(surface) :: sfc_state
88
89 ! A pointer to a structure containing metrics and related information.
90 type(ocean_grid_type), pointer :: grid => null()
91 type(verticalgrid_type), pointer :: gv => null()
92 ! A pointer to a structure containing dimensional unit scaling factors.
93 type(unit_scale_type), pointer :: us => null()
94
95 ! If .true., use the ice shelf model for part of the domain.
96 logical :: use_ice_shelf = .false.
97
98 ! If .true., use surface wave coupling
99 logical :: use_waves = .false.
100
101 ! This is .true. if incremental restart files may be saved.
102 logical :: permit_incr_restart = .true.
103
104 ! nmax is the number of iterations after which to stop so that the
105 ! simulation does not exceed its CPU time limit. nmax is determined by
106 ! evaluating the CPU time used between successive calls to write_cputime.
107 ! Initially it is set to be very large.
108 integer :: nmax=2000000000
109
110 ! A structure containing several relevant directory paths.
111 type(directories) :: dirs
112
113 ! A suite of time types for use by MOM
114 type(time_type), target :: time ! A copy of the ocean model's time.
115 ! Other modules can set pointers to this and
116 ! change it to manage diagnostics.
117 type(time_type) :: master_time ! The ocean model's master clock. No other
118 ! modules are ever given access to this.
119 type(time_type) :: time1 ! The value of the ocean model's time at the
120 ! start of a call to step_MOM.
121 type(time_type) :: start_time ! The start time of the simulation.
122 type(time_type) :: segment_start_time ! The start time of this run segment.
123 type(time_type) :: time_end ! End time for the segment or experiment.
124 type(time_type) :: restart_time ! The next time to write restart files.
125 type(time_type) :: time_step_ocean ! A time_type version of dt_forcing.
126 logical :: segment_start_time_set ! True if segment_start_time has been set to a valid value.
127
128 real :: elapsed_time = 0.0 ! Elapsed time in this run [T ~> s].
129 logical :: elapsed_time_master ! If true, elapsed time is used to set the model's master
130 ! clock (Time). This is needed if Time_step_ocean is not
131 ! an exact representation of dt_forcing.
132 real :: dt_forcing ! The coupling time step [T ~> s].
133 real :: dt ! The nominal baroclinic dynamics time step [T ~> s].
134 integer :: ntstep ! The number of baroclinic dynamics time steps within dt_forcing.
135 real :: dt_therm ! The thermodynamic timestep [T ~> s]
136 real :: dt_dyn ! The actual dynamic timestep used [T ~> s]. The value of dt_dyn
137 ! is chosen so that dt_forcing is an integer multiple of dt_dyn.
138 real :: dtdia ! The diabatic timestep [T ~> s]
139 real :: t_elapsed_seg ! The elapsed time in this run segment [T ~> s]
140 integer :: n, ns, n_max, nts, n_last_thermo
141 logical :: diabatic_first, single_step_call, initialize_smb
142 type(time_type) :: time2, time_chg ! Temporary time variables
143
144 integer :: restart_control ! An integer that is bit-tested to determine whether
145 ! incremental restart files are saved and whether they
146 ! have a time stamped name. +1 (bit 0) for generic
147 ! files and +2 (bit 1) for time-stamped files. A
148 ! restart file is saved at the end of a run segment
149 ! unless Restart_control is negative.
150
151 real :: time_unit ! The time unit for the following input fields [s].
152 type(time_type) :: restint ! The time between saves of the restart file.
153 type(time_type) :: daymax ! The final day of the simulation.
154
155 integer :: cpu_steps ! The number of steps between writing CPU time.
156 integer :: date(6) ! Possibly the start date of this run segment.
157 type(param_file_type) :: param_file ! The structure indicating the file(s)
158 ! containing all run-time parameters.
159
160 integer :: calendar_type=-1 ! A coded integer indicating the calendar type.
161
162 integer :: unit, io_status, ierr
163 integer :: initclock, mainclock, termclock
164
165 logical :: debug ! If true, write verbose checksums for debugging purposes.
166 logical :: offline_tracer_mode ! If false, use the model in prognostic mode where
167 ! the barotropic and baroclinic dynamics, thermodynamics,
168 ! etc. are stepped forward integrated in time.
169 ! If true, then all of the above are bypassed with all
170 ! fields necessary to integrate only the tracer advection
171 ! and diffusion equation are read in from files stored from
172 ! a previous integration of the prognostic model
173
174 type(mom_control_struct) :: mom_csp !< The control structure with all the MOM6 internal types,
175 !! parameters and variables
176 type(tracer_flow_control_cs), pointer :: &
177 tracer_flow_csp => null() !< A pointer to the tracer flow control structure
178 type(surface_forcing_cs), pointer :: surface_forcing_csp => null()
179 type(write_cputime_cs), pointer :: write_cpu_csp => null()
180 type(ice_shelf_cs), pointer :: ice_shelf_csp => null()
181 logical :: override_shelf_fluxes !< If true, and shelf dynamics are active,
182 !! the data_override feature is enabled (only for MOSAIC grid types)
183 type(wave_parameters_cs), pointer :: waves_csp => null()
184 type(diag_ctrl), pointer :: &
185 diag => null() !< A pointer to the diagnostic regulatory structure
186 !-----------------------------------------------------------------------
187
188 character(len=4), parameter :: vers_num = 'v2.0'
189 ! This include declares and sets the variable "version".
190# include "version_variable.h"
191 character(len=40) :: mod_name = "MOM_main (MOM_driver)" ! This module's name.
192
193 ! These are the variables that might be read via the namelist capability.
194 integer :: date_init(6)=0 ! The start date of the whole simulation.
195 character(len=16) :: calendar = 'julian' ! The name of the calendar type.
196 integer :: years=0, months=0, days=0 ! These may determine the segment run
197 integer :: hours=0, minutes=0, seconds=0 ! length, if read from a namelist.
198 integer :: ocean_nthreads = 1
199 logical :: use_hyper_thread = .false.
200 namelist /ocean_solo_nml/ date_init, calendar, months, days, hours, minutes, seconds, &
201 ocean_nthreads, use_hyper_thread
202
203 !=====================================================================
204
205 call write_cputime_start_clock(write_cpu_csp)
206
207 call mom_infra_init() ; call io_infra_init()
208
209 !allocate(forces,fluxes,sfc_state)
210
211 ! Initialize the ensemble manager based on settings in input.nml(ensemble.nml).
212 call initialize_ocean_only_ensembles()
213
214 ! These clocks are on the global pelist.
215 initclock = cpu_clock_id( 'Initialization' )
216 mainclock = cpu_clock_id( 'Main loop' )
217 termclock = cpu_clock_id( 'Termination' )
218 call cpu_clock_begin(initclock)
219
220 call mom_mesg('======== Model being driven by MOM_driver ========', 2)
221 call calltree_waypoint("Program MOM_main, MOM_driver.F90")
222
223 if (file_exists('input.nml')) then
224 ! Provide for namelist specification of the run length and calendar data.
225 call open_ascii_file(unit, 'input.nml', action=readonly_file)
226 read(unit, ocean_solo_nml, iostat=io_status)
227 call close_file(unit)
228 ierr = check_nml_error(io_status,'ocean_solo_nml')
229 if (is_root_pe() .and. (years+months+days+hours+minutes+seconds > 0)) write(*,ocean_solo_nml)
230 endif
231
232 ! This call sets the number and affinity of threads with openMP.
233 !$ call set_MOM_thread_affinity(ocean_nthreads, use_hyper_thread)
234
235 ! This call is required to initiate dirs%restart_input_dir for ocean_solo.res
236 ! The contents of dirs will be reread in initialize_MOM.
237 call get_mom_input(dirs=dirs)
238
239 segment_start_time_set = .false.
240 ! Read ocean_solo restart, which can override settings from the namelist.
241 if (file_exists(trim(dirs%restart_input_dir)//'ocean_solo.res')) then
242 date(:) = -1
243 call open_ascii_file(unit, trim(dirs%restart_input_dir)//'ocean_solo.res', action=readonly_file)
244 read(unit,*) calendar_type
245 read(unit,*) date_init
246 read(unit,*) date
247 call close_file(unit)
248
249 call set_calendar_type(calendar_type)
250 if (sum(date) >= 0) then
251 ! In this case, the segment starts at a time fixed by ocean_solo.res
252 segment_start_time = set_date(date(1), date(2), date(3), date(4), date(5), date(6))
253 segment_start_time_set = .true.
254 endif
255 else
256 calendar = uppercase(calendar)
257 if (calendar(1:6) == 'JULIAN') then ; calendar_type = julian
258 elseif (calendar(1:9) == 'GREGORIAN') then ; calendar_type = gregorian
259 elseif (calendar(1:6) == 'NOLEAP') then ; calendar_type = noleap
260 elseif (calendar(1:10)=='THIRTY_DAY') then ; calendar_type = thirty_day_months
261 elseif (calendar(1:11)=='NO_CALENDAR') then ; calendar_type = no_calendar
262 elseif (calendar(1:1) /= ' ') then
263 call mom_error(fatal,'MOM_driver: Invalid namelist value '//trim(calendar)//' for calendar')
264 else
265 call mom_error(fatal,'MOM_driver: No namelist value for calendar')
266 endif
267 call set_calendar_type(calendar_type)
268 endif
269
270
271 if (sum(date_init) > 0) then
272 start_time = set_date(date_init(1), date_init(2), date_init(3), &
273 date_init(4), date_init(5), date_init(6))
274 else
275 start_time = real_to_time(0.0)
276 endif
277
278 call time_interp_external_init()
279
280 ! Call initialize MOM with an optional Ice Shelf CS which, if present triggers
281 ! initialization of ice shelf parameters and arrays.
282 if (segment_start_time_set) then
283 ! In this case, the segment starts at a time fixed by ocean_solo.res
284 time = segment_start_time
285 call initialize_mom(time, start_time, param_file, dirs, mom_csp, &
286 segment_start_time, offline_tracer_mode=offline_tracer_mode, &
287 diag_ptr=diag, tracer_flow_csp=tracer_flow_csp, ice_shelf_csp=ice_shelf_csp, &
288 waves_csp=waves_csp)
289 else
290 ! In this case, the segment starts at a time read from the MOM restart file
291 ! or is left at Start_time by MOM_initialize.
292 time = start_time
293 call initialize_mom(time, start_time, param_file, dirs, mom_csp, &
294 offline_tracer_mode=offline_tracer_mode, diag_ptr=diag, &
295 tracer_flow_csp=tracer_flow_csp, ice_shelf_csp=ice_shelf_csp, waves_csp=waves_csp)
296 endif
297
298 call get_mom_state_elements(mom_csp, g=grid, gv=gv, us=us, c_p_scaled=fluxes%C_p)
299 master_time = time
300 use_ice_shelf = associated(ice_shelf_csp)
301
302 if (use_ice_shelf) then
303 ! These arrays are not initialized in most solo cases, but are needed
304 ! when using an ice shelf
305 call initialize_ice_shelf_fluxes(ice_shelf_csp, grid, us, fluxes)
306 call initialize_ice_shelf_forces(ice_shelf_csp, grid, us, forces)
307 call ice_shelf_query(ice_shelf_csp, grid, data_override_shelf_fluxes=override_shelf_fluxes)
308 if (override_shelf_fluxes) call data_override_init(ocean_domain_in=grid%domain%mpp_domain)
309 call get_param(param_file, mod_name, "INITIALIZE_ICE_SHEET_SMB", &
310 initialize_smb, "Read in a constant SMB for the ice sheet", default=.false.)
311 if (initialize_smb) call initialize_ice_smb(fluxes%shelf_sfc_mass_flux, grid, us, param_file)
312 endif
313
314
315 call calltree_waypoint("done initialize_MOM")
316
317 call extract_surface_state(mom_csp, sfc_state)
318
319 if (use_ice_shelf .and. allocated(sfc_state%frazil)) &
320 call adjust_ice_sheet_frazil(sfc_state, fluxes, ice_shelf_csp)
321
322 call surface_forcing_init(time, grid, us, param_file, diag, &
323 surface_forcing_csp, tracer_flow_csp)
324 call calltree_waypoint("done surface_forcing_init")
325
326
327 call get_param(param_file, mod_name, "USE_WAVES", use_waves, &
328 "If true, enables surface wave modules.",default=.false.)
329 ! MOM_wave_interface_init is called regardless of the value of USE_WAVES because
330 ! it also initializes statistical waves.
331 call mom_wave_interface_init(time, grid, gv, us, param_file, waves_csp, diag)
332
333 segment_start_time = time
334 elapsed_time = 0.0
335
336 ! Read all relevant parameters and write them to the model log.
337 call log_version(param_file, mod_name, version, "")
338 call get_param(param_file, mod_name, "DT", dt, &
339 units="s", scale=us%s_to_T, fail_if_missing=.true.)
340 call get_param(param_file, mod_name, "DT_FORCING", dt_forcing, &
341 "The time step for changing forcing, coupling with other "//&
342 "components, or potentially writing certain diagnostics. "//&
343 "The default value is given by DT.", &
344 units="s", default=us%T_to_s*dt, scale=us%s_to_T)
345 if (offline_tracer_mode) then
346 call get_param(param_file, mod_name, "DT_OFFLINE", dt_forcing, &
347 "Length of time between reading in of input fields", &
348 units="s", scale=us%s_to_T, fail_if_missing=.true.)
349 dt = dt_forcing
350 endif
351 ntstep = max(1,ceiling(dt_forcing/dt - 0.001))
352
353 time_step_ocean = real_to_time(dt_forcing, unscale=us%T_to_s)
354 elapsed_time_master = (abs(dt_forcing - time_to_real(time_step_ocean, scale=us%s_to_T)) > 1.0e-12*dt_forcing)
355 if (elapsed_time_master) &
356 call mom_mesg("Using real elapsed time for the master clock.", 2)
357
358 ! Determine the segment end time, either from the namelist file or parsed input file.
359 ! Note that Time_unit always is in [s].
360 call get_param(param_file, mod_name, "TIMEUNIT", time_unit, &
361 "The time unit for DAYMAX, ENERGYSAVEDAYS, and RESTINT.", &
362 units="s", default=86400.0)
363 if (years+months+days+hours+minutes+seconds > 0) then
364 time_end = increment_date(time, years, months, days, hours, minutes, seconds)
365 call mom_mesg('Segment run length determined from ocean_solo_nml.', 2)
366 call get_param(param_file, mod_name, "DAYMAX", daymax, timeunit=time_unit, &
367 default=time_end, do_not_log=.true.)
368 call log_param(param_file, mod_name, "DAYMAX", daymax, &
369 "The final time of the whole simulation, in units of "//&
370 "TIMEUNIT seconds. This also sets the potential end "//&
371 "time of the present run segment if the end time is "//&
372 "not set via ocean_solo_nml in input.nml.", &
373 timeunit=time_unit)
374 else
375 call get_param(param_file, mod_name, "DAYMAX", daymax, &
376 "The final time of the whole simulation, in units of "//&
377 "TIMEUNIT seconds. This also sets the potential end "//&
378 "time of the present run segment if the end time is "//&
379 "not set via ocean_solo_nml in input.nml.", &
380 timeunit=time_unit, fail_if_missing=.true.)
381 time_end = daymax
382 endif
383
384 call diag_manager_set_time_end_infra(time_end)
385
386 call get_param(param_file, mod_name, "SINGLE_STEPPING_CALL", single_step_call, &
387 "If true, advance the state of MOM with a single step "//&
388 "including both dynamics and thermodynamics. If false "//&
389 "the two phases are advanced with separate calls.", default=.true.)
390 call get_param(param_file, mod_name, "DT_THERM", dt_therm, &
391 "The thermodynamic and tracer advection time step. "//&
392 "Ideally DT_THERM should be an integer multiple of DT "//&
393 "and less than the forcing or coupling time-step, unless "//&
394 "THERMO_SPANS_COUPLING is true, in which case DT_THERM "//&
395 "can be an integer multiple of the coupling timestep. By "//&
396 "default DT_THERM is set to DT.", &
397 units="s", default=us%T_to_s*dt, scale=us%s_to_T)
398 call get_param(param_file, mod_name, "DIABATIC_FIRST", diabatic_first, &
399 "If true, apply diabatic and thermodynamic processes, "//&
400 "including buoyancy forcing and mass gain or loss, "//&
401 "before stepping the dynamics forward.", default=.false.)
402
403
404 if (time >= time_end) call mom_error(fatal, &
405 "MOM_driver: The run has been started at or after the end time of the run.")
406
407 call get_param(param_file, mod_name, "RESTART_CONTROL", restart_control, &
408 "An integer whose bits encode which restart files are "//&
409 "written. Add 2 (bit 1) for a time-stamped file, and odd "//&
410 "(bit 0) for a non-time-stamped file. A non-time-stamped "//&
411 "restart file is saved at the end of the run segment "//&
412 "for any non-negative value.", default=1)
413 call get_param(param_file, mod_name, "RESTINT", restint, &
414 "The interval between saves of the restart file in units "//&
415 "of TIMEUNIT. Use 0 (the default) to not save "//&
416 "incremental restart files at all.", default=real_to_time(0.0), &
417 timeunit=time_unit)
418 call get_param(param_file, mod_name, "WRITE_CPU_STEPS", cpu_steps, &
419 "The number of coupled timesteps between writing the cpu "//&
420 "time. If this is not positive, do not check cpu time, and "//&
421 "the segment run-length can not be set via an elapsed CPU time.", &
422 default=1000)
423 call get_param(param_file, "MOM", "DEBUG", debug, &
424 "If true, write out verbose debugging data.", &
425 default=.false., debuggingparam=.true.)
426
427 call log_param(param_file, mod_name, "ELAPSED TIME AS MASTER", elapsed_time_master)
428
429 if (cpu_steps > 0) &
430 call mom_write_cputime_init(param_file, dirs%output_directory, start_time, &
431 write_cpu_csp)
432
433 ! Close the param_file. No further parsing of input is possible after this.
434 call close_param_file(param_file)
436
437 ! Write out a time stamp file.
438 if (is_root_pe() .and. (calendar_type /= no_calendar)) call write_time_stamp_file(time)
439
440 if (cpu_steps > 0) call write_cputime(time, 0, write_cpu_csp)
441
442 if (((.not.btest(restart_control,1)) .and. (.not.btest(restart_control,0))) &
443 .or. (restart_control < 0)) permit_incr_restart = .false.
444
445 if (restint > real_to_time(0.0)) then
446 ! restart_time is the next integral multiple of restint.
447 restart_time = start_time + restint * &
448 (1 + ((time + time_step_ocean) - start_time) / restint)
449 else
450 ! Set the time so late that there is no intermediate restart.
451 restart_time = time_end + time_step_ocean
452 permit_incr_restart = .false.
453 endif
454
455 call cpu_clock_end(initclock) !end initialization
456
457 call cpu_clock_begin(mainclock) !begin main loop
458
459 ns = 1
460 do while ((ns < nmax) .and. (time < time_end))
461 call calltree_enter("Main loop, MOM_driver.F90",ns)
462
463 ! Set the forcing for the next steps.
464 if (.not. offline_tracer_mode) then
465 call set_forcing(sfc_state, forces, fluxes, time, time_step_ocean, grid, us, &
466 surface_forcing_csp)
467 endif
468 if (debug) then
469 call mom_mech_forcing_chksum("After set forcing", forces, grid, us, haloshift=0)
470 call mom_forcing_chksum("After set forcing", fluxes, grid, us, haloshift=0)
471 endif
472
473 if (use_ice_shelf) then
474 call shelf_calc_flux(sfc_state, fluxes, time, dt_forcing, ice_shelf_csp)
475 call add_shelf_forces(grid, us, ice_shelf_csp, forces, external_call=.true.)
476 endif
477 fluxes%fluxes_used = .false.
478 fluxes%dt_buoy_accum = dt_forcing
479
480 if (use_waves) then
481 call update_surface_waves(grid, gv, us, time, time_step_ocean, waves_csp)
482 endif
483
484 if (ns==1) then
485 call finish_mom_initialization(time, dirs, mom_csp)
486 endif
487
488 ! This call steps the model over a time dt_forcing.
489 time1 = master_time ; time = master_time
490 if (offline_tracer_mode) then
491 call step_offline(forces, fluxes, sfc_state, time1, dt_forcing, mom_csp)
492 elseif (single_step_call) then
493 call step_mom(forces, fluxes, sfc_state, time1, dt_forcing, mom_csp, waves=waves_csp)
494 else
495 n_max = 1 ; if (dt_forcing > dt) n_max = ceiling(dt_forcing/dt - 0.001)
496 dt_dyn = dt_forcing / real(n_max)
497
498 nts = max(1,min(n_max,floor(dt_therm/dt_dyn + 0.001)))
499 n_last_thermo = 0
500
501 time2 = time1 ; t_elapsed_seg = 0.0
502 do n=1,n_max
503 if (diabatic_first) then
504 if (modulo(n-1,nts)==0) then
505 dtdia = dt_dyn*min(ntstep,n_max-(n-1))
506 call step_mom(forces, fluxes, sfc_state, time2, dtdia, mom_csp, &
507 do_dynamics=.false., do_thermodynamics=.true., &
508 start_cycle=(n==1), end_cycle=.false., cycle_length=dt_forcing)
509 endif
510
511 call step_mom(forces, fluxes, sfc_state, time2, dt_dyn, mom_csp, &
512 do_dynamics=.true., do_thermodynamics=.false., &
513 start_cycle=.false., end_cycle=(n==n_max), cycle_length=dt_forcing)
514 else
515 call step_mom(forces, fluxes, sfc_state, time2, dt_dyn, mom_csp, &
516 do_dynamics=.true., do_thermodynamics=.false., &
517 start_cycle=(n==1), end_cycle=.false., cycle_length=dt_forcing)
518
519 if ((modulo(n,nts)==0) .or. (n==n_max)) then
520 dtdia = dt_dyn*(n - n_last_thermo)
521 ! Back up Time2 to the start of the thermodynamic segment.
522 if (n > n_last_thermo+1) &
523 time2 = time2 - real_to_time((dtdia - dt_dyn), unscale=us%T_to_s)
524 call step_mom(forces, fluxes, sfc_state, time2, dtdia, mom_csp, &
525 do_dynamics=.false., do_thermodynamics=.true., &
526 start_cycle=.false., end_cycle=(n==n_max), cycle_length=dt_forcing)
527 n_last_thermo = n
528 endif
529 endif
530
531 t_elapsed_seg = t_elapsed_seg + dt_dyn
532 time2 = time1 + real_to_time(t_elapsed_seg, unscale=us%T_to_s)
533 enddo
534 endif
535
536! Time = Time + Time_step_ocean
537! This is here to enable fractional-second time steps.
538 elapsed_time = elapsed_time + dt_forcing
539 if (elapsed_time > 2.0e9*us%s_to_T) then
540 ! This is here to ensure that the conversion from a real to an integer can be accurately
541 ! represented in long runs (longer than ~63 years). It will also ensure that elapsed time
542 ! does not lose resolution of order the timetype's resolution, provided that the timestep and
543 ! tick are larger than 10-5 seconds. If a clock with a finer resolution is used, a smaller
544 ! value would be required.
545 time_chg = real_to_time(elapsed_time, unscale=us%T_to_s)
546 segment_start_time = segment_start_time + time_chg
547 elapsed_time = elapsed_time - time_to_real(time_chg, scale=us%s_to_T)
548 endif
549 if (elapsed_time_master) then
550 master_time = segment_start_time + real_to_time(elapsed_time, unscale=us%T_to_s)
551 else
552 master_time = master_time + time_step_ocean
553 endif
554 time = master_time
555
556 if (cpu_steps > 0) then ; if (mod(ns, cpu_steps) == 0) then
557 call write_cputime(time, ns+ntstep-1, write_cpu_csp, nmax)
558 endif ; endif
559
560 call mech_forcing_diags(forces, dt_forcing, grid, time, diag, surface_forcing_csp%handles)
561
562 if (.not. offline_tracer_mode) then
563 if (fluxes%fluxes_used) then
564 call forcing_diagnostics(fluxes, sfc_state, grid, us, time, &
565 diag, surface_forcing_csp%handles)
566 else
567 call mom_error(fatal, "The solo MOM_driver is not yet set up to handle "//&
568 "thermodynamic time steps that are longer than the coupling timestep.")
569 endif
570 endif
571
572! See if it is time to write out a restart file - timestamped or not.
573 if ((permit_incr_restart) .and. (fluxes%fluxes_used) .and. &
574 (time + (time_step_ocean/2) > restart_time)) then
575 if (btest(restart_control,1)) then
576 call save_mom_restart(mom_csp, dirs%restart_output_dir, time, grid, &
577 time_stamped=.true., gv=gv)
578 call forcing_save_restart(surface_forcing_csp, grid, time, &
579 dirs%restart_output_dir, .true.)
580 if (use_ice_shelf) call ice_shelf_save_restart(ice_shelf_csp, time, &
581 dirs%restart_output_dir, .true.)
582 endif
583 if (btest(restart_control,0)) then
584 call save_mom_restart(mom_csp, dirs%restart_output_dir, time, grid, gv=gv)
585 call forcing_save_restart(surface_forcing_csp, grid, time, &
586 dirs%restart_output_dir)
587 if (use_ice_shelf) call ice_shelf_save_restart(ice_shelf_csp, time, &
588 dirs%restart_output_dir)
589 endif
590 restart_time = restart_time + restint
591 endif
592
593 ns = ns + ntstep
594 call calltree_leave("Main loop")
595 enddo
596
597 call cpu_clock_end(mainclock)
598 call cpu_clock_begin(termclock)
599 if (restart_control>=0) then
600 if (.not.mom_state_is_synchronized(mom_csp)) &
601 call mom_error(warning, "End of MOM_main reached with inconsistent "//&
602 "dynamics and advective times. Additional restart fields "//&
603 "that have not been coded yet would be required for reproducibility.")
604 if (.not.fluxes%fluxes_used .and. .not.offline_tracer_mode) call mom_error(fatal, &
605 "End of MOM_main reached with unused buoyancy fluxes. "//&
606 "For conservation, the ocean restart files can only be "//&
607 "created after the buoyancy forcing is applied.")
608
609 call save_mom_restart(mom_csp, dirs%restart_output_dir, time, grid, gv=gv)
610 if (use_ice_shelf) call ice_shelf_save_restart(ice_shelf_csp, time, &
611 dirs%restart_output_dir)
612
613 ! Write the ocean solo restart file.
614 call write_ocean_solo_res(time, start_time, calendar_type, &
615 trim(dirs%restart_output_dir)//'ocean_solo.res')
616 endif
617
618 if (is_root_pe()) then
619 call open_ascii_file(unit, "exitcode")
620 if (time < daymax) then
621 write(unit,*) 9
622 else
623 write(unit,*) 0
624 endif
625 call close_file(unit)
626 endif
627
628 call calltree_waypoint("End MOM_main")
629 if (use_ice_shelf) call ice_shelf_end(ice_shelf_csp)
630 call diag_mediator_end(time, diag, end_diag_manager=.true.)
631 if (cpu_steps > 0) call write_cputime(time, ns-1, write_cpu_csp, call_end=.true.)
632 call cpu_clock_end(termclock)
633
634 call mom_end(mom_csp)
635
636 ! This closes out the infrastructure, including clocks, I/O and message passing communicators.
637 call io_infra_end() ; call mom_infra_end()
638
639contains
640
641!> Write out the ocean solo restart file to the indicated path.
642subroutine write_ocean_solo_res(Time, Start_time, calendar, file_path)
643 type(time_type), intent(in) :: Time !< The current model time.
644 type(time_type), intent(in) :: Start_Time !< The start time of the simulation.
645 integer, intent(in) :: calendar !< A coded integer indicating the calendar type.
646 character(len=*), intent(in) :: file_path !< The full path and name of the restart file
647
648 ! Local variables
649 integer :: unit
650 integer :: yr, mon, day, hr, mins, sec ! Temp variables for writing the date.
651
652 if (.not.is_root_pe()) return
653
654 call open_ascii_file(unit, trim(file_path))
655 write(unit, '(i6,8x,a)') calendar, &
656 '(Calendar: no_calendar=0, thirty_day_months=1, julian=2, gregorian=3, noleap=4)'
657
658 call get_date(start_time, yr, mon, day, hr, mins, sec)
659 write(unit, '(6i6,8x,a)') yr, mon, day, hr, mins, sec, &
660 'Model start time: year, month, day, hour, minute, second'
661 call get_date(time, yr, mon, day, hr, mins, sec)
662 write(unit, '(6i6,8x,a)') yr, mon, day, hr, mins, sec, &
663 'Current model time: year, month, day, hour, minute, second'
664 call close_file(unit)
665end subroutine write_ocean_solo_res
666
667
668!> Write out an ascii time stamp file with the model time, following FMS conventions.
669subroutine write_time_stamp_file(Time)
670 type(time_type), intent(in) :: Time !< The current model time.
671 ! Local variables
672 integer :: unit
673 integer :: yr, mon, day, hr, mins, sec ! Temp variables for writing the date.
674 character(len=9) :: month ! The name of the month
675
676 if (.not.is_root_pe()) return
677
678 call open_ascii_file(unit, 'time_stamp.out', action=append_file)
679 call get_date(time, yr, mon, day, hr, mins, sec)
680 month = month_name(mon)
681 write(unit,'(6i4,2x,a3)') yr, mon, day, hr, mins, sec, month(1:3)
682 call get_date(time_end, yr, mon, day, hr, mins, sec)
683 month = month_name(mon)
684 write(unit,'(6i4,2x,a3)') yr, mon, day, hr, mins, sec, month(1:3)
685 call close_file(unit)
686end subroutine write_time_stamp_file
687
688!> Initialize the ensemble manager. If there are no settings for ensemble_size
689!! in input.nml(ensemble.nml), these should not do anything. In coupled
690!! configurations, this all occurs in the external driver.
691subroutine initialize_ocean_only_ensembles()
692 integer, dimension(:), allocatable :: ocean_PElist
693 integer, dimension(0) :: atm_PElist, land_PElist, ice_PElist
694 integer :: ensemble_size, nPEs_per, ensemble_info(6)
695
696 call ensemble_manager_init() ; ensemble_info(:) = get_ensemble_size()
697 ensemble_size = ensemble_info(1) ; npes_per = ensemble_info(2)
698 if (ensemble_size > 1) then ! There are multiple ensemble members.
699 allocate(ocean_pelist(npes_per))
700 call ensemble_pelist_setup(.true., 0, npes_per, 0, 0, atm_pelist, ocean_pelist, &
701 land_pelist, ice_pelist)
702 call set_pelist(ocean_pelist)
703 deallocate(ocean_pelist)
704 endif
705end subroutine initialize_ocean_only_ensembles
706
707end program mom6