Line data Source code
1 : !=======================================================================
2 : !
3 : ! Main driver for time stepping of CICE.
4 : !
5 : ! authors Elizabeth C. Hunke, LANL
6 : ! Philip W. Jones, LANL
7 : ! William H. Lipscomb, LANL
8 : !
9 : ! 2006 ECH: moved exit timeLoop to prevent execution of unnecessary timestep
10 : ! 2006 ECH: Streamlined for efficiency
11 : ! 2006 ECH: Converted to free source form (F90)
12 : ! 2007 BPB: Modified Delta-Eddington shortwave interface
13 : ! 2008 ECH: moved ESMF code to its own driver
14 :
15 : module CICE_RunMod
16 :
17 : use ice_kinds_mod
18 : use ice_fileunits, only: nu_diag
19 : use ice_arrays_column, only: oceanmixed_ice
20 : use ice_constants, only: c0, c1
21 : use ice_constants, only: field_loc_center, field_type_scalar
22 : use ice_exit, only: abort_ice
23 : use icepack_intfc, only: icepack_warnings_flush, icepack_warnings_aborted
24 : use icepack_intfc, only: icepack_max_iso, icepack_max_aero
25 : use icepack_intfc, only: icepack_query_parameters
26 : use icepack_intfc, only: icepack_query_tracer_flags, icepack_query_tracer_sizes
27 :
28 : implicit none
29 : private
30 : public :: CICE_Run, ice_step
31 :
32 : !=======================================================================
33 :
34 : contains
35 :
36 : !=======================================================================
37 : !
38 : ! This is the main driver routine for advancing CICE forward in time.
39 : !
40 : ! author Elizabeth C. Hunke, LANL
41 : ! Philip W. Jones, LANL
42 : ! William H. Lipscomb, LANL
43 :
44 2800 : subroutine CICE_Run
45 :
46 : use ice_calendar, only: istep, istep1, time, dt, stop_now, calendar
47 : use ice_forcing, only: get_forcing_atmo, get_forcing_ocn, &
48 : get_wave_spec
49 : use ice_forcing_bgc, only: get_forcing_bgc, get_atm_bgc, &
50 : fiso_default, faero_default
51 : use ice_flux, only: init_flux_atm, init_flux_ocn
52 : use ice_timers, only: ice_timer_start, ice_timer_stop, &
53 : timer_couple, timer_step
54 : logical (kind=log_kind) :: &
55 : tr_iso, tr_aero, tr_zaero, skl_bgc, z_tracers, wave_spec, tr_fsd
56 : character(len=*), parameter :: subname = '(CICE_Run)'
57 :
58 : !--------------------------------------------------------------------
59 : ! initialize error code and step timer
60 : !--------------------------------------------------------------------
61 :
62 2800 : call ice_timer_start(timer_step) ! start timing entire run
63 :
64 : call icepack_query_parameters(skl_bgc_out=skl_bgc, &
65 : z_tracers_out=z_tracers, &
66 2800 : wave_spec_out=wave_spec)
67 : call icepack_query_tracer_flags(tr_iso_out=tr_iso, &
68 : tr_aero_out=tr_aero, &
69 : tr_zaero_out=tr_zaero, &
70 2800 : tr_fsd_out=tr_fsd)
71 2800 : call icepack_warnings_flush(nu_diag)
72 2800 : if (icepack_warnings_aborted()) call abort_ice(error_message=subname, &
73 0 : file=__FILE__, line=__LINE__)
74 :
75 : #ifndef CICE_IN_NEMO
76 : !--------------------------------------------------------------------
77 : ! timestep loop
78 : !--------------------------------------------------------------------
79 :
80 602648 : timeLoop: do
81 : #endif
82 :
83 605448 : call ice_step
84 :
85 605448 : istep = istep + 1 ! update time step counters
86 605448 : istep1 = istep1 + 1
87 605448 : time = time + dt ! determine the time and date
88 :
89 605448 : call calendar(time) ! at the end of the timestep
90 :
91 : #ifndef CICE_IN_NEMO
92 605448 : if (stop_now >= 1) exit timeLoop
93 : #endif
94 :
95 602648 : call ice_timer_start(timer_couple) ! atm/ocn coupling
96 :
97 : ! for now, wave_spectrum is constant in time
98 : ! if (tr_fsd .and. wave_spec) call get_wave_spec ! wave spectrum in ice
99 602648 : call get_forcing_atmo ! atmospheric forcing from data
100 602648 : call get_forcing_ocn(dt) ! ocean forcing from data
101 :
102 : ! isotopes
103 602648 : if (tr_iso) call fiso_default ! default values
104 : ! aerosols
105 : ! if (tr_aero) call faero_data ! data file
106 : ! if (tr_zaero) call fzaero_data ! data file (gx1)
107 602648 : if (tr_aero .or. tr_zaero) call faero_default ! default values
108 :
109 602648 : if (skl_bgc .or. z_tracers) call get_forcing_bgc ! biogeochemistry
110 602648 : if (z_tracers) call get_atm_bgc ! biogeochemistry
111 :
112 602648 : call init_flux_atm ! Initialize atmosphere fluxes sent to coupler
113 602648 : call init_flux_ocn ! initialize ocean fluxes sent to coupler
114 :
115 605448 : call ice_timer_stop(timer_couple) ! atm/ocn coupling
116 :
117 : #ifndef CICE_IN_NEMO
118 : enddo timeLoop
119 : #endif
120 :
121 : !--------------------------------------------------------------------
122 : ! end of timestep loop
123 : !--------------------------------------------------------------------
124 :
125 2800 : call ice_timer_stop(timer_step) ! end timestepping loop timer
126 :
127 2800 : end subroutine CICE_Run
128 :
129 : !=======================================================================
130 : !
131 : ! Calls drivers for physics components, some initialization, and output
132 : !
133 : ! author Elizabeth C. Hunke, LANL
134 : ! William H. Lipscomb, LANL
135 :
136 605448 : subroutine ice_step
137 :
138 : use ice_boundary, only: ice_HaloUpdate
139 : use ice_calendar, only: dt, dt_dyn, ndtd, diagfreq, write_restart, istep
140 : use ice_diagnostics, only: init_mass_diags, runtime_diags
141 : use ice_diagnostics_bgc, only: hbrine_diags, zsal_diags, bgc_diags
142 : use ice_domain, only: halo_info, nblocks
143 : use ice_dyn_eap, only: write_restart_eap
144 : use ice_dyn_shared, only: kdyn, kridge
145 : use ice_flux, only: scale_factor, init_history_therm, &
146 : daidtt, daidtd, dvidtt, dvidtd, dagedtt, dagedtd
147 : use ice_history, only: accum_hist
148 : use ice_history_bgc, only: init_history_bgc
149 : use ice_restart, only: final_restart
150 : use ice_restart_column, only: write_restart_age, write_restart_FY, &
151 : write_restart_lvl, write_restart_pond_cesm, write_restart_pond_lvl, &
152 : write_restart_pond_topo, write_restart_aero, write_restart_fsd, &
153 : write_restart_iso, write_restart_bgc, write_restart_hbrine
154 : use ice_restart_driver, only: dumpfile
155 : use ice_restoring, only: restore_ice, ice_HaloRestore
156 : use ice_step_mod, only: prep_radiation, step_therm1, step_therm2, &
157 : update_state, step_dyn_horiz, step_dyn_ridge, step_radiation, &
158 : biogeochemistry, save_init, step_dyn_wave
159 : use ice_timers, only: ice_timer_start, ice_timer_stop, &
160 : timer_diags, timer_column, timer_thermo, timer_bound, &
161 : timer_hist, timer_readwrite
162 :
163 : integer (kind=int_kind) :: &
164 : iblk , & ! block index
165 : k , & ! dynamics supercycling index
166 : ktherm ! thermodynamics is off when ktherm = -1
167 :
168 : real (kind=dbl_kind) :: &
169 53064 : offset ! d(age)/dt time offset
170 :
171 : logical (kind=log_kind) :: &
172 : tr_iage, tr_FY, tr_lvl, tr_fsd, &
173 : tr_pond_cesm, tr_pond_lvl, tr_pond_topo, tr_brine, tr_iso, tr_aero, &
174 : calc_Tsfc, skl_bgc, solve_zsal, z_tracers, wave_spec
175 :
176 : character(len=*), parameter :: subname = '(ice_step)'
177 :
178 : call icepack_query_parameters(calc_Tsfc_out=calc_Tsfc, skl_bgc_out=skl_bgc, &
179 : solve_zsal_out=solve_zsal, z_tracers_out=z_tracers, ktherm_out=ktherm, &
180 605448 : wave_spec_out=wave_spec)
181 : call icepack_query_tracer_flags(tr_iage_out=tr_iage, tr_FY_out=tr_FY, &
182 : tr_lvl_out=tr_lvl, tr_pond_cesm_out=tr_pond_cesm, tr_pond_lvl_out=tr_pond_lvl, &
183 : tr_pond_topo_out=tr_pond_topo, tr_brine_out=tr_brine, tr_aero_out=tr_aero, &
184 605448 : tr_iso_out=tr_iso, tr_fsd_out=tr_fsd)
185 605448 : call icepack_warnings_flush(nu_diag)
186 605448 : if (icepack_warnings_aborted()) call abort_ice(error_message=subname, &
187 0 : file=__FILE__, line=__LINE__)
188 :
189 : !-----------------------------------------------------------------
190 : ! restoring on grid boundaries
191 : !-----------------------------------------------------------------
192 :
193 605448 : if (restore_ice) call ice_HaloRestore
194 :
195 : !-----------------------------------------------------------------
196 : ! initialize diagnostics and save initial state values
197 : !-----------------------------------------------------------------
198 :
199 605448 : call ice_timer_start(timer_diags) ! diagnostics/history
200 605448 : call init_mass_diags ! diagnostics per timestep
201 605448 : call init_history_therm
202 605448 : call init_history_bgc
203 605448 : call ice_timer_stop(timer_diags) ! diagnostics/history
204 :
205 605448 : call ice_timer_start(timer_column) ! column physics
206 605448 : call ice_timer_start(timer_thermo) ! thermodynamics
207 :
208 605448 : call save_init
209 :
210 338520 : !$OMP PARALLEL DO PRIVATE(iblk)
211 1440840 : do iblk = 1, nblocks
212 :
213 1779360 : if (ktherm >= 0) then
214 :
215 : !-----------------------------------------------------------------
216 : ! scale radiation fields
217 : !-----------------------------------------------------------------
218 :
219 1173264 : if (calc_Tsfc) call prep_radiation (iblk)
220 :
221 : !-----------------------------------------------------------------
222 : ! thermodynamics and biogeochemistry
223 : !-----------------------------------------------------------------
224 :
225 1173264 : call step_therm1 (dt, iblk) ! vertical thermodynamics
226 1173264 : call biogeochemistry (dt, iblk) ! biogeochemistry
227 1173264 : call step_therm2 (dt, iblk) ! ice thickness distribution thermo
228 :
229 : endif ! ktherm > 0
230 :
231 : enddo ! iblk
232 : !$OMP END PARALLEL DO
233 :
234 : ! clean up, update tendency diagnostics
235 605448 : offset = dt
236 605448 : call update_state (dt, daidtt, dvidtt, dagedtt, offset)
237 :
238 605448 : call ice_timer_stop(timer_thermo) ! thermodynamics
239 605448 : call ice_timer_stop(timer_column) ! column physics
240 :
241 : !-----------------------------------------------------------------
242 : ! dynamics, transport, ridging
243 : !-----------------------------------------------------------------
244 :
245 : ! wave fracture of the floe size distribution
246 : ! note this is called outside of the dynamics subcycling loop
247 605448 : if (tr_fsd .and. wave_spec) call step_dyn_wave(dt)
248 :
249 1212432 : do k = 1, ndtd
250 :
251 : ! momentum, stress, transport
252 606984 : call step_dyn_horiz (dt_dyn)
253 :
254 : ! ridging
255 340056 : !$OMP PARALLEL DO PRIVATE(iblk)
256 1440840 : do iblk = 1, nblocks
257 1780896 : if (kridge > 0) call step_dyn_ridge (dt_dyn, ndtd, iblk)
258 : enddo
259 : !$OMP END PARALLEL DO
260 :
261 : ! clean up, update tendency diagnostics
262 606984 : offset = c0
263 1212432 : call update_state (dt_dyn, daidtd, dvidtd, dagedtd, offset)
264 :
265 : enddo
266 :
267 : !-----------------------------------------------------------------
268 : ! albedo, shortwave radiation
269 : !-----------------------------------------------------------------
270 :
271 605448 : call ice_timer_start(timer_column) ! column physics
272 605448 : call ice_timer_start(timer_thermo) ! thermodynamics
273 :
274 : !MHRI: CHECK THIS OMP
275 338520 : !$OMP PARALLEL DO PRIVATE(iblk)
276 1440840 : do iblk = 1, nblocks
277 :
278 1173912 : if (ktherm >= 0) call step_radiation (dt, iblk)
279 :
280 : !-----------------------------------------------------------------
281 : ! get ready for coupling and the next time step
282 : !-----------------------------------------------------------------
283 :
284 1779360 : call coupling_prep (iblk)
285 :
286 : enddo ! iblk
287 : !$OMP END PARALLEL DO
288 :
289 605448 : call ice_timer_start(timer_bound)
290 : call ice_HaloUpdate (scale_factor, halo_info, &
291 605448 : field_loc_center, field_type_scalar)
292 605448 : call ice_timer_stop(timer_bound)
293 :
294 605448 : call ice_timer_stop(timer_thermo) ! thermodynamics
295 605448 : call ice_timer_stop(timer_column) ! column physics
296 :
297 : !-----------------------------------------------------------------
298 : ! write data
299 : !-----------------------------------------------------------------
300 :
301 605448 : call ice_timer_start(timer_diags) ! diagnostics
302 605448 : if (mod(istep,diagfreq) == 0) then
303 32472 : call runtime_diags(dt) ! log file
304 32472 : if (solve_zsal) call zsal_diags
305 32472 : if (skl_bgc .or. z_tracers) call bgc_diags
306 32472 : if (tr_brine) call hbrine_diags
307 : endif
308 605448 : call ice_timer_stop(timer_diags) ! diagnostics
309 :
310 605448 : call ice_timer_start(timer_hist) ! history
311 605448 : call accum_hist (dt) ! history file
312 605448 : call ice_timer_stop(timer_hist) ! history
313 :
314 605448 : call ice_timer_start(timer_readwrite) ! reading/writing
315 605448 : if (write_restart == 1) then
316 4146 : call dumpfile ! core variables for restarting
317 4146 : if (tr_iage) call write_restart_age
318 4146 : if (tr_FY) call write_restart_FY
319 4146 : if (tr_lvl) call write_restart_lvl
320 4146 : if (tr_pond_cesm) call write_restart_pond_cesm
321 4146 : if (tr_pond_lvl) call write_restart_pond_lvl
322 4146 : if (tr_pond_topo) call write_restart_pond_topo
323 4146 : if (tr_fsd) call write_restart_fsd
324 4146 : if (tr_iso) call write_restart_iso
325 4146 : if (tr_aero) call write_restart_aero
326 4146 : if (solve_zsal .or. skl_bgc .or. z_tracers) &
327 921 : call write_restart_bgc
328 4146 : if (tr_brine) call write_restart_hbrine
329 4146 : if (kdyn == 2) call write_restart_eap
330 4146 : call final_restart
331 : endif
332 :
333 605448 : call ice_timer_stop(timer_readwrite) ! reading/writing
334 :
335 605448 : end subroutine ice_step
336 :
337 : !=======================================================================
338 : !
339 : ! Prepare for coupling
340 : !
341 : ! authors: Elizabeth C. Hunke, LANL
342 :
343 4326744 : subroutine coupling_prep (iblk)
344 :
345 : use ice_arrays_column, only: alvdfn, alidfn, alvdrn, alidrn, &
346 : albicen, albsnon, albpndn, apeffn, fzsal_g, fzsal, snowfracn
347 : use ice_blocks, only: nx_block, ny_block, get_block, block
348 : use ice_domain, only: blocks_ice
349 : use ice_calendar, only: dt, nstreams
350 : use ice_domain_size, only: ncat
351 : use ice_flux, only: alvdf, alidf, alvdr, alidr, albice, albsno, &
352 : albpnd, albcnt, apeff_ai, fpond, fresh, l_mpond_fresh, &
353 : alvdf_ai, alidf_ai, alvdr_ai, alidr_ai, fhocn_ai, &
354 : fresh_ai, fsalt_ai, fsalt, &
355 : fswthru_ai, fhocn, scale_factor, snowfrac, &
356 : fswthru, fswthru_vdr, fswthru_vdf, fswthru_idr, fswthru_idf, &
357 : swvdr, swidr, swvdf, swidf, Tf, Tair, Qa, strairxT, strairyT, &
358 : fsens, flat, fswabs, flwout, evap, Tref, Qref, &
359 : scale_fluxes, frzmlt_init, frzmlt
360 : use ice_flux_bgc, only: faero_ocn, fiso_ocn, Qref_iso, fiso_evap, &
361 : fzsal_ai, fzsal_g_ai, flux_bio, flux_bio_ai
362 : use ice_grid, only: tmask
363 : use ice_state, only: aicen, aice
364 : #ifdef CICE_IN_NEMO
365 : use ice_state, only: aice_init
366 : use ice_flux, only: flatn_f, fsurfn_f
367 : #endif
368 : use ice_step_mod, only: ocean_mixed_layer
369 : use ice_timers, only: timer_couple, ice_timer_start, ice_timer_stop
370 :
371 : integer (kind=int_kind), intent(in) :: &
372 : iblk ! block index
373 :
374 : ! local variables
375 :
376 : integer (kind=int_kind) :: &
377 : ilo,ihi,jlo,jhi, & ! beginning and end of physical domain
378 : n , & ! thickness category index
379 : i,j , & ! horizontal indices
380 : k , & ! tracer index
381 : nbtrcr !
382 :
383 : type (block) :: &
384 : this_block ! block information for current block
385 :
386 : logical (kind=log_kind) :: &
387 : calc_Tsfc !
388 :
389 : real (kind=dbl_kind) :: &
390 292440 : cszn , & ! counter for history averaging
391 292440 : puny , & !
392 292440 : rhofresh , & !
393 292440 : netsw ! flag for shortwave radiation presence
394 :
395 : character(len=*), parameter :: subname = '(coupling_prep)'
396 :
397 4326744 : call icepack_query_parameters(puny_out=puny, rhofresh_out=rhofresh)
398 4326744 : call icepack_query_tracer_sizes(nbtrcr_out=nbtrcr)
399 4326744 : call icepack_query_parameters(calc_Tsfc_out=calc_Tsfc)
400 4326744 : call icepack_warnings_flush(nu_diag)
401 4326744 : if (icepack_warnings_aborted()) call abort_ice(error_message=subname, &
402 0 : file=__FILE__, line=__LINE__)
403 :
404 : !-----------------------------------------------------------------
405 : ! Save current value of frzmlt for diagnostics.
406 : ! Update mixed layer with heat and radiation from ice.
407 : !-----------------------------------------------------------------
408 :
409 122734032 : do j = 1, ny_block
410 1425314616 : do i = 1, nx_block
411 1420987872 : frzmlt_init (i,j,iblk) = frzmlt(i,j,iblk)
412 : enddo
413 : enddo
414 :
415 4326744 : call ice_timer_start(timer_couple) ! atm/ocn coupling
416 :
417 4326744 : if (oceanmixed_ice) &
418 4326744 : call ocean_mixed_layer (dt,iblk) ! ocean surface fluxes and sst
419 :
420 : !-----------------------------------------------------------------
421 : ! Aggregate albedos
422 : !-----------------------------------------------------------------
423 :
424 122734032 : do j = 1, ny_block
425 1425314616 : do i = 1, nx_block
426 1302580584 : alvdf(i,j,iblk) = c0
427 1302580584 : alidf(i,j,iblk) = c0
428 1302580584 : alvdr(i,j,iblk) = c0
429 1302580584 : alidr(i,j,iblk) = c0
430 :
431 1302580584 : albice(i,j,iblk) = c0
432 1302580584 : albsno(i,j,iblk) = c0
433 1302580584 : albpnd(i,j,iblk) = c0
434 1302580584 : apeff_ai(i,j,iblk) = c0
435 1302580584 : snowfrac(i,j,iblk) = c0
436 :
437 : ! for history averaging
438 1302580584 : cszn = c0
439 1302580584 : netsw = swvdr(i,j,iblk)+swidr(i,j,iblk)+swvdf(i,j,iblk)+swidf(i,j,iblk)
440 1302580584 : if (netsw > puny) cszn = c1
441 3453448824 : do n = 1, nstreams
442 3335041536 : albcnt(i,j,iblk,n) = albcnt(i,j,iblk,n) + cszn
443 : enddo
444 : enddo
445 : enddo
446 :
447 4326744 : this_block = get_block(blocks_ice(iblk),iblk)
448 4326744 : ilo = this_block%ilo
449 4326744 : ihi = this_block%ihi
450 4326744 : jlo = this_block%jlo
451 4326744 : jhi = this_block%jhi
452 :
453 25204104 : do n = 1, ncat
454 551730504 : do j = jlo, jhi
455 5377160400 : do i = ilo, ihi
456 5356283040 : if (aicen(i,j,n,iblk) > puny) then
457 :
458 0 : alvdf(i,j,iblk) = alvdf(i,j,iblk) &
459 769514868 : + alvdfn(i,j,n,iblk)*aicen(i,j,n,iblk)
460 0 : alidf(i,j,iblk) = alidf(i,j,iblk) &
461 769514868 : + alidfn(i,j,n,iblk)*aicen(i,j,n,iblk)
462 0 : alvdr(i,j,iblk) = alvdr(i,j,iblk) &
463 769514868 : + alvdrn(i,j,n,iblk)*aicen(i,j,n,iblk)
464 0 : alidr(i,j,iblk) = alidr(i,j,iblk) &
465 769514868 : + alidrn(i,j,n,iblk)*aicen(i,j,n,iblk)
466 :
467 0 : netsw = swvdr(i,j,iblk) + swidr(i,j,iblk) &
468 769514868 : + swvdf(i,j,iblk) + swidf(i,j,iblk)
469 769514868 : if (netsw > puny) then ! sun above horizon
470 0 : albice(i,j,iblk) = albice(i,j,iblk) &
471 354945742 : + albicen(i,j,n,iblk)*aicen(i,j,n,iblk)
472 0 : albsno(i,j,iblk) = albsno(i,j,iblk) &
473 354945742 : + albsnon(i,j,n,iblk)*aicen(i,j,n,iblk)
474 0 : albpnd(i,j,iblk) = albpnd(i,j,iblk) &
475 354945742 : + albpndn(i,j,n,iblk)*aicen(i,j,n,iblk)
476 : endif
477 :
478 0 : apeff_ai(i,j,iblk) = apeff_ai(i,j,iblk) & ! for history
479 769514868 : + apeffn(i,j,n,iblk)*aicen(i,j,n,iblk)
480 0 : snowfrac(i,j,iblk) = snowfrac(i,j,iblk) & ! for history
481 769514868 : + snowfracn(i,j,n,iblk)*aicen(i,j,n,iblk)
482 :
483 : endif ! aicen > puny
484 : enddo
485 : enddo
486 : enddo
487 :
488 122734032 : do j = 1, ny_block
489 1425314616 : do i = 1, nx_block
490 :
491 : !-----------------------------------------------------------------
492 : ! reduce fresh by fpond for coupling
493 : !-----------------------------------------------------------------
494 :
495 1302580584 : if (l_mpond_fresh) then
496 0 : fpond(i,j,iblk) = fpond(i,j,iblk) * rhofresh/dt
497 0 : fresh(i,j,iblk) = fresh(i,j,iblk) - fpond(i,j,iblk)
498 : endif
499 :
500 : !----------------------------------------------------------------
501 : ! Store grid box mean albedos and fluxes before scaling by aice
502 : !----------------------------------------------------------------
503 :
504 1302580584 : alvdf_ai (i,j,iblk) = alvdf (i,j,iblk)
505 1302580584 : alidf_ai (i,j,iblk) = alidf (i,j,iblk)
506 1302580584 : alvdr_ai (i,j,iblk) = alvdr (i,j,iblk)
507 1302580584 : alidr_ai (i,j,iblk) = alidr (i,j,iblk)
508 1302580584 : fresh_ai (i,j,iblk) = fresh (i,j,iblk)
509 1302580584 : fsalt_ai (i,j,iblk) = fsalt (i,j,iblk)
510 1302580584 : fhocn_ai (i,j,iblk) = fhocn (i,j,iblk)
511 1302580584 : fswthru_ai(i,j,iblk) = fswthru(i,j,iblk)
512 1302580584 : fzsal_ai (i,j,iblk) = fzsal (i,j,iblk)
513 1302580584 : fzsal_g_ai(i,j,iblk) = fzsal_g(i,j,iblk)
514 :
515 1302580584 : if (nbtrcr > 0) then
516 3190426416 : do k = 1, nbtrcr
517 3190426416 : flux_bio_ai (i,j,k,iblk) = flux_bio (i,j,k,iblk)
518 : enddo
519 : endif
520 :
521 : !-----------------------------------------------------------------
522 : ! Save net shortwave for scaling factor in scale_factor
523 : !-----------------------------------------------------------------
524 0 : scale_factor(i,j,iblk) = &
525 0 : swvdr(i,j,iblk)*(c1 - alvdr_ai(i,j,iblk)) &
526 0 : + swvdf(i,j,iblk)*(c1 - alvdf_ai(i,j,iblk)) &
527 0 : + swidr(i,j,iblk)*(c1 - alidr_ai(i,j,iblk)) &
528 1420987872 : + swidf(i,j,iblk)*(c1 - alidf_ai(i,j,iblk))
529 :
530 : enddo
531 : enddo
532 :
533 : !-----------------------------------------------------------------
534 : ! Divide fluxes by ice area
535 : ! - the CESM coupler assumes fluxes are per unit ice area
536 : ! - also needed for global budget in diagnostics
537 : !-----------------------------------------------------------------
538 :
539 : call scale_fluxes (nx_block, ny_block, &
540 0 : tmask (:,:,iblk), nbtrcr, &
541 : icepack_max_aero, &
542 0 : aice (:,:,iblk), Tf (:,:,iblk), &
543 0 : Tair (:,:,iblk), Qa (:,:,iblk), &
544 0 : strairxT (:,:,iblk), strairyT(:,:,iblk), &
545 0 : fsens (:,:,iblk), flat (:,:,iblk), &
546 0 : fswabs (:,:,iblk), flwout (:,:,iblk), &
547 0 : evap (:,:,iblk), &
548 0 : Tref (:,:,iblk), Qref (:,:,iblk), &
549 0 : fresh (:,:,iblk), fsalt (:,:,iblk), &
550 0 : fhocn (:,:,iblk), &
551 0 : fswthru (:,:,iblk), &
552 0 : fswthru_vdr (:,:,iblk), &
553 0 : fswthru_vdf (:,:,iblk), &
554 0 : fswthru_idr (:,:,iblk), &
555 0 : fswthru_idf (:,:,iblk), &
556 0 : faero_ocn(:,:,:,iblk), &
557 0 : alvdr (:,:,iblk), alidr (:,:,iblk), &
558 0 : alvdf (:,:,iblk), alidf (:,:,iblk), &
559 0 : fzsal (:,:,iblk), fzsal_g (:,:,iblk), &
560 0 : flux_bio (:,:,1:nbtrcr,iblk), &
561 0 : Qref_iso =Qref_iso (:,:,:,iblk), &
562 0 : fiso_evap=fiso_evap(:,:,:,iblk), &
563 4326744 : fiso_ocn =fiso_ocn (:,:,:,iblk))
564 :
565 : #ifdef CICE_IN_NEMO
566 : !echmod - comment this out for efficiency, if .not. calc_Tsfc
567 : if (.not. calc_Tsfc) then
568 :
569 : !---------------------------------------------------------------
570 : ! If surface fluxes were provided, conserve these fluxes at ice
571 : ! free points by passing to ocean.
572 : !---------------------------------------------------------------
573 :
574 : call sfcflux_to_ocn &
575 : (nx_block, ny_block, &
576 : tmask (:,:,iblk), aice_init(:,:,iblk), &
577 : fsurfn_f (:,:,:,iblk), flatn_f(:,:,:,iblk), &
578 : fresh (:,:,iblk), fhocn (:,:,iblk))
579 : endif
580 : !echmod
581 : #endif
582 4326744 : call ice_timer_stop(timer_couple) ! atm/ocn coupling
583 :
584 4326744 : end subroutine coupling_prep
585 :
586 : #ifdef CICE_IN_NEMO
587 :
588 : !=======================================================================
589 : !
590 : ! If surface heat fluxes are provided to CICE instead of CICE calculating
591 : ! them internally (i.e. .not. calc_Tsfc), then these heat fluxes can
592 : ! be provided at points which do not have ice. (This is could be due to
593 : ! the heat fluxes being calculated on a lower resolution grid or the
594 : ! heat fluxes not recalculated at every CICE timestep.) At ice free points,
595 : ! conserve energy and water by passing these fluxes to the ocean.
596 : !
597 : ! author: A. McLaren, Met Office
598 :
599 : subroutine sfcflux_to_ocn(nx_block, ny_block, &
600 : tmask, aice, &
601 : fsurfn_f, flatn_f, &
602 : fresh, fhocn)
603 :
604 : use ice_domain_size, only: ncat
605 :
606 : integer (kind=int_kind), intent(in) :: &
607 : nx_block, ny_block ! block dimensions
608 :
609 : logical (kind=log_kind), dimension (nx_block,ny_block), intent(in) :: &
610 : tmask ! land/boundary mask, thickness (T-cell)
611 :
612 : real (kind=dbl_kind), dimension(nx_block,ny_block), intent(in):: &
613 : aice ! initial ice concentration
614 :
615 : real (kind=dbl_kind), dimension(nx_block,ny_block,ncat), intent(in) :: &
616 : fsurfn_f, & ! net surface heat flux (provided as forcing)
617 : flatn_f ! latent heat flux (provided as forcing)
618 :
619 : real (kind=dbl_kind), dimension(nx_block,ny_block), intent(inout):: &
620 : fresh , & ! fresh water flux to ocean (kg/m2/s)
621 : fhocn ! actual ocn/ice heat flx (W/m**2)
622 :
623 :
624 : ! local variables
625 : integer (kind=int_kind) :: &
626 : i, j, n ! horizontal indices
627 :
628 : real (kind=dbl_kind) :: &
629 : puny, & !
630 : rLsub ! 1/Lsub
631 :
632 : character(len=*), parameter :: subname = '(sfcflux_to_ocn)'
633 :
634 : call icepack_query_parameters(puny_out=puny)
635 : call icepack_warnings_flush(nu_diag)
636 : if (icepack_warnings_aborted()) call abort_ice(error_message=subname, &
637 : file=__FILE__, line=__LINE__)
638 : rLsub = c1 / Lsub
639 :
640 : do n = 1, ncat
641 : do j = 1, ny_block
642 : do i = 1, nx_block
643 : if (tmask(i,j) .and. aice(i,j) <= puny) then
644 : fhocn(i,j) = fhocn(i,j) &
645 : + fsurfn_f(i,j,n) + flatn_f(i,j,n)
646 : fresh(i,j) = fresh(i,j) &
647 : + flatn_f(i,j,n) * rLsub
648 : endif
649 : enddo ! i
650 : enddo ! j
651 : enddo ! n
652 :
653 :
654 : end subroutine sfcflux_to_ocn
655 :
656 : #endif
657 :
658 : !=======================================================================
659 :
660 : end module CICE_RunMod
661 :
662 : !=======================================================================
|