7 use mom_cpu_clock,
only : cpu_clock_id, cpu_clock_begin, cpu_clock_end, clock_routine
21 use coupler_types_mod,
only : coupler_2d_bc_type, coupler_type_spawn
22 use coupler_types_mod,
only : coupler_type_increment_data, coupler_type_initialized
23 use coupler_types_mod,
only : coupler_type_copy_data, coupler_type_destructor
25 implicit none ;
private
27 #include <MOM_memory.h>
53 real,
pointer,
dimension(:,:) :: &
55 ustar_gustless => null()
59 real,
pointer,
dimension(:,:) :: &
63 real,
pointer,
dimension(:,:) :: &
65 sw_vis_dir => null(), &
66 sw_vis_dif => null(), &
67 sw_nir_dir => null(), &
68 sw_nir_dif => null(), &
72 real,
pointer,
dimension(:,:) :: &
75 seaice_melt_heat => null(), &
79 real,
pointer,
dimension(:,:) :: &
80 latent_evap_diag => null(), &
81 latent_fprec_diag => null(), &
82 latent_frunoff_diag => null()
85 real,
pointer,
dimension(:,:) :: &
92 seaice_melt => null(), &
93 netmassin => null(), &
94 netmassout => null(), &
98 real,
pointer,
dimension(:,:) :: &
99 heat_content_cond => null(), &
100 heat_content_lprec => null(), &
101 heat_content_icemelt => null(), &
103 heat_content_fprec => null(), &
104 heat_content_vprec => null(), &
105 heat_content_lrunoff => null(), &
106 heat_content_frunoff => null(), &
107 heat_content_massout => null(), &
108 heat_content_massin => null()
111 real,
pointer,
dimension(:,:) :: &
112 salt_flux => null(), &
113 salt_flux_in => null(), &
114 salt_flux_added => null()
118 real,
pointer,
dimension(:,:) :: p_surf_full => null()
121 real,
pointer,
dimension(:,:) :: p_surf => null()
124 real,
pointer,
dimension(:,:) :: p_surf_ssh => null()
128 logical :: accumulate_p_surf = .false.
134 real,
pointer,
dimension(:,:) :: &
135 tke_tidal => null(), &
136 ustar_tidal => null()
139 real,
pointer,
dimension(:,:) :: &
140 ustar_berg => null(), &
141 area_berg => null(), &
145 real,
pointer,
dimension(:,:) :: ustar_shelf => null()
147 real,
pointer,
dimension(:,:) :: frac_shelf_h => null()
151 real,
pointer,
dimension(:,:) :: iceshelf_melt => null()
155 real :: vprecglobaladj = 0.
156 real :: saltfluxglobaladj = 0.
157 real :: netfwglobaladj = 0.
158 real :: vprecglobalscl = 0.
159 real :: saltfluxglobalscl = 0.
160 real :: netfwglobalscl = 0.
162 logical :: fluxes_used = .true.
164 real :: dt_buoy_accum = -1.0
172 type(coupler_2d_bc_type) :: tr_fluxes
178 integer :: num_msg = 0
179 integer :: max_msg = 2
189 real,
pointer,
dimension(:,:) :: &
193 net_mass_src => null()
196 real,
pointer,
dimension(:,:) :: p_surf_full => null()
199 real,
pointer,
dimension(:,:) :: p_surf => null()
202 real,
pointer,
dimension(:,:) :: p_surf_ssh => null()
208 real,
pointer,
dimension(:,:) :: &
209 area_berg => null(), &
213 real,
pointer,
dimension(:,:) :: frac_shelf_u => null()
216 real,
pointer,
dimension(:,:) :: frac_shelf_v => null()
219 real,
pointer,
dimension(:,:) :: &
220 rigidity_ice_u => null(), &
221 rigidity_ice_v => null()
222 real :: dt_force_accum = -1.0
224 logical :: net_mass_src_set = .false.
225 logical :: accumulate_p_surf = .false.
229 logical :: accumulate_rigidity = .false.
233 logical :: initialized = .false.
241 integer :: id_prcme = -1, id_evap = -1
242 integer :: id_precip = -1, id_vprec = -1
243 integer :: id_lprec = -1, id_fprec = -1
244 integer :: id_lrunoff = -1, id_frunoff = -1
245 integer :: id_net_massout = -1, id_net_massin = -1
246 integer :: id_massout_flux = -1, id_massin_flux = -1
247 integer :: id_seaice_melt = -1
250 integer :: id_total_prcme = -1, id_total_evap = -1
251 integer :: id_total_precip = -1, id_total_vprec = -1
252 integer :: id_total_lprec = -1, id_total_fprec = -1
253 integer :: id_total_lrunoff = -1, id_total_frunoff = -1
254 integer :: id_total_net_massout = -1, id_total_net_massin = -1
255 integer :: id_total_seaice_melt = -1
258 integer :: id_prcme_ga = -1, id_evap_ga = -1
259 integer :: id_lprec_ga = -1, id_fprec_ga= -1
260 integer :: id_precip_ga = -1, id_vprec_ga= -1
263 integer :: id_net_heat_coupler = -1, id_net_heat_surface = -1
264 integer :: id_sens = -1, id_lwlatsens = -1
265 integer :: id_sw = -1, id_lw = -1
266 integer :: id_sw_vis = -1, id_sw_nir = -1
267 integer :: id_lat_evap = -1, id_lat_frunoff = -1
268 integer :: id_lat = -1, id_lat_fprec = -1
269 integer :: id_heat_content_lrunoff= -1, id_heat_content_frunoff = -1
270 integer :: id_heat_content_lprec = -1, id_heat_content_fprec = -1
271 integer :: id_heat_content_cond = -1, id_heat_content_surfwater= -1
272 integer :: id_heat_content_vprec = -1, id_heat_content_massout = -1
273 integer :: id_heat_added = -1, id_heat_content_massin = -1
274 integer :: id_hfrainds = -1, id_hfrunoffds = -1
275 integer :: id_seaice_melt_heat = -1, id_heat_content_icemelt = -1
278 integer :: id_total_net_heat_coupler = -1, id_total_net_heat_surface = -1
279 integer :: id_total_sens = -1, id_total_lwlatsens = -1
280 integer :: id_total_sw = -1, id_total_lw = -1
281 integer :: id_total_lat_evap = -1, id_total_lat_frunoff = -1
282 integer :: id_total_lat = -1, id_total_lat_fprec = -1
283 integer :: id_total_heat_content_lrunoff= -1, id_total_heat_content_frunoff = -1
284 integer :: id_total_heat_content_lprec = -1, id_total_heat_content_fprec = -1
285 integer :: id_total_heat_content_cond = -1, id_total_heat_content_surfwater= -1
286 integer :: id_total_heat_content_vprec = -1, id_total_heat_content_massout = -1
287 integer :: id_total_heat_added = -1, id_total_heat_content_massin = -1
288 integer :: id_total_seaice_melt_heat = -1, id_total_heat_content_icemelt = -1
291 integer :: id_net_heat_coupler_ga = -1, id_net_heat_surface_ga = -1
292 integer :: id_sens_ga = -1, id_lwlatsens_ga = -1
293 integer :: id_sw_ga = -1, id_lw_ga = -1
294 integer :: id_lat_ga = -1
297 integer :: id_saltflux = -1
298 integer :: id_saltfluxin = -1
299 integer :: id_saltfluxadded = -1
301 integer :: id_total_saltflux = -1
302 integer :: id_total_saltfluxin = -1
303 integer :: id_total_saltfluxadded = -1
305 integer :: id_vprecglobaladj = -1
306 integer :: id_vprecglobalscl = -1
307 integer :: id_saltfluxglobaladj = -1
308 integer :: id_saltfluxglobalscl = -1
309 integer :: id_netfwglobaladj = -1
310 integer :: id_netfwglobalscl = -1
313 integer :: id_taux = -1
314 integer :: id_tauy = -1
315 integer :: id_ustar = -1
317 integer :: id_psurf = -1
318 integer :: id_tke_tidal = -1
319 integer :: id_buoy = -1
322 integer :: id_ustar_berg = -1
323 integer :: id_area_berg = -1
324 integer :: id_mass_berg = -1
327 integer :: id_ustar_ice_cover = -1
328 integer :: id_frac_ice_cover = -1
331 integer :: id_clock_forcing = -1
341 subroutine extractfluxes1d(G, GV, US, fluxes, optics, nsw, j, dt_in_T, &
342 FluxRescaleDepth, useRiverHeatContent, useCalvingHeatContent, &
343 h, T, netMassInOut, netMassOut, net_heat, net_salt, pen_SW_bnd, tv, &
344 aggregate_FW, nonpenSW, netmassInOut_rate, net_Heat_Rate, &
345 net_salt_rate, pen_sw_bnd_Rate, skip_diags)
350 type(
forcing),
intent(inout) :: fluxes
353 integer,
intent(in) :: nsw
354 integer,
intent(in) :: j
355 real,
intent(in) :: dt_in_t
356 real,
intent(in) :: fluxrescaledepth
358 logical,
intent(in) :: useriverheatcontent
359 logical,
intent(in) :: usecalvingheatcontent
360 real,
dimension(SZI_(G),SZK_(G)), &
362 real,
dimension(SZI_(G),SZK_(G)), &
364 real,
dimension(SZI_(G)),
intent(out) :: netmassinout
367 real,
dimension(SZI_(G)),
intent(out) :: netmassout
371 real,
dimension(SZI_(G)),
intent(out) :: net_heat
378 real,
dimension(SZI_(G)),
intent(out) :: net_salt
381 real,
dimension(max(1,nsw),G%isd:G%ied),
intent(out) :: pen_sw_bnd
390 logical,
intent(in) :: aggregate_fw
391 real,
dimension(SZI_(G)), &
392 optional,
intent(out) :: nonpensw
395 real,
dimension(SZI_(G)), &
396 optional,
intent(out) :: net_heat_rate
398 real,
dimension(SZI_(G)), &
399 optional,
intent(out) :: net_salt_rate
401 real,
dimension(SZI_(G)), &
402 optional,
intent(out) :: netmassinout_rate
404 real,
dimension(max(1,nsw),G%isd:G%ied), &
405 optional,
intent(out) :: pen_sw_bnd_rate
407 logical,
optional,
intent(in) :: skip_diags
410 real :: htot(szi_(g))
411 real :: pen_sw_tot(szi_(g))
412 real :: pen_sw_tot_rate(szi_(g))
418 real :: rz_t_to_w_m2_degc
423 logical :: calculate_diags
424 character(len=200) :: mesg
425 integer :: is, ie, nz, i, k, n
427 logical :: do_nhr, do_nsr, do_nmior, do_pswbr
437 if (
present(net_heat_rate)) do_nhr = .true.
438 if (
present(net_salt_rate)) do_nsr = .true.
439 if (
present(netmassinout_rate)) do_nmior = .true.
440 if (
present(pen_sw_bnd_rate)) do_pswbr = .true.
443 ih_limit = 0.0 ;
if (fluxrescaledepth > 0.0) ih_limit = 1.0 / fluxrescaledepth
444 rz_t_to_w_m2_degc = fluxes%C_p*us%R_to_kg_m3*us%Z_to_m*us%s_to_T
445 i_cp = 1.0 / fluxes%C_p
446 w_m2_to_h_t = 1.0 / (us%s_to_T * gv%H_to_kg_m2 * fluxes%C_p)
448 rzcp_to_h = 1.0 / (gv%H_to_RZ * fluxes%C_p)
450 is = g%isc ; ie = g%iec ; nz = g%ke
452 calculate_diags = .true.
453 if (
present(skip_diags)) calculate_diags = .not. skip_diags
458 "mismatch in the number of bands of shortwave radiation in MOM_forcing_type extract_fluxes.")
461 if (.not.
associated(fluxes%sw))
call mom_error(fatal, &
462 "MOM_forcing_type extractFluxes1d: fluxes%sw is not associated.")
464 if (.not.
associated(fluxes%lw))
call mom_error(fatal, &
465 "MOM_forcing_type extractFluxes1d: fluxes%lw is not associated.")
467 if (.not.
associated(fluxes%latent))
call mom_error(fatal, &
468 "MOM_forcing_type extractFluxes1d: fluxes%latent is not associated.")
470 if (.not.
associated(fluxes%sens))
call mom_error(fatal, &
471 "MOM_forcing_type extractFluxes1d: fluxes%sens is not associated.")
473 if (.not.
associated(fluxes%evap))
call mom_error(fatal, &
474 "MOM_forcing_type extractFluxes1d: No evaporation defined.")
476 if (.not.
associated(fluxes%vprec))
call mom_error(fatal, &
477 "MOM_forcing_type extractFluxes1d: fluxes%vprec not defined.")
479 if ((.not.
associated(fluxes%lprec)) .or. &
480 (.not.
associated(fluxes%fprec)))
call mom_error(fatal, &
481 "MOM_forcing_type extractFluxes1d: No precipitation defined.")
483 do i=is,ie ; htot(i) = h(i,1) ;
enddo
484 do k=2,nz ;
do i=is,ie ; htot(i) = htot(i) + h(i,k) ;
enddo ;
enddo
493 scale = 1.0 ;
if ((ih_limit > 0.0) .and. (htot(i)*ih_limit < 1.0)) scale = htot(i)*ih_limit
500 pen_sw_bnd(n,i) = w_m2_to_h_t*scale*dt_in_t * max(0.0, pen_sw_bnd(n,i))
501 pen_sw_tot(i) = pen_sw_tot(i) + pen_sw_bnd(n,i)
504 pen_sw_bnd(1,i) = 0.0
508 pen_sw_tot_rate(i) = 0.0
511 pen_sw_bnd_rate(n,i) = w_m2_to_h_t*scale * max(0.0, pen_sw_bnd_rate(n,i))
512 pen_sw_tot_rate(i) = pen_sw_tot_rate(i) + pen_sw_bnd_rate(n,i)
515 pen_sw_bnd_rate(1,i) = 0.0
520 netmassinout(i) = dt_in_t * (scale * &
521 (((((( fluxes%lprec(i,j) &
522 + fluxes%fprec(i,j) ) &
523 + fluxes%evap(i,j) ) &
524 + fluxes%lrunoff(i,j) ) &
525 + fluxes%vprec(i,j) ) &
526 + fluxes%seaice_melt(i,j)) &
527 + fluxes%frunoff(i,j) ))
530 netmassinout_rate(i) = (scale * &
531 (((((( fluxes%lprec(i,j) &
532 + fluxes%fprec(i,j) ) &
533 + fluxes%evap(i,j) ) &
534 + fluxes%lrunoff(i,j) ) &
535 + fluxes%vprec(i,j) ) &
536 + fluxes%seaice_melt(i,j)) &
537 + fluxes%frunoff(i,j) ))
545 if (.not.gv%Boussinesq .and.
associated(fluxes%salt_flux))
then
546 netmassinout(i) = netmassinout(i) + dt_in_t * (scale * fluxes%salt_flux(i,j))
547 if (do_nmior) netmassinout_rate(i) = netmassinout_rate(i) + &
548 (scale * fluxes%salt_flux(i,j))
558 if (fluxes%evap(i,j) < 0.0) netmassout(i) = netmassout(i) + fluxes%evap(i,j)
563 if (fluxes%lprec(i,j) < 0.0) netmassout(i) = netmassout(i) + fluxes%lprec(i,j)
566 if (fluxes%seaice_melt(i,j) < 0.0) netmassout(i) = netmassout(i) + fluxes%seaice_melt(i,j)
570 if (fluxes%vprec(i,j) < 0.0) netmassout(i) = netmassout(i) + fluxes%vprec(i,j)
572 netmassout(i) = dt_in_t * scale * netmassout(i)
575 netmassinout(i) = gv%RZ_to_H * netmassinout(i)
576 if (do_nmior) netmassinout_rate(i) = gv%RZ_to_H * netmassinout_rate(i)
577 netmassout(i) = gv%RZ_to_H * netmassout(i)
583 if (
associated(fluxes%seaice_melt_heat))
then
584 net_heat(i) = scale * dt_in_t * w_m2_to_h_t * &
585 ( fluxes%sw(i,j) + ((fluxes%lw(i,j) + fluxes%latent(i,j)) + fluxes%sens(i,j) + &
586 fluxes%seaice_melt_heat(i,j)) )
588 if (do_nhr) net_heat_rate(i) = scale * w_m2_to_h_t * &
589 ( fluxes%sw(i,j) + ((fluxes%lw(i,j) + fluxes%latent(i,j)) + fluxes%sens(i,j) + &
590 fluxes%seaice_melt_heat(i,j)))
592 net_heat(i) = scale * dt_in_t * w_m2_to_h_t * &
593 ( fluxes%sw(i,j) + ((fluxes%lw(i,j) + fluxes%latent(i,j)) + fluxes%sens(i,j)) )
595 if (do_nhr) net_heat_rate(i) = scale * w_m2_to_h_t * &
596 ( fluxes%sw(i,j) + ((fluxes%lw(i,j) + fluxes%latent(i,j)) + fluxes%sens(i,j)) )
600 if (
associated(fluxes%heat_added))
then
601 net_heat(i) = net_heat(i) + (scale * (dt_in_t * w_m2_to_h_t)) * fluxes%heat_added(i,j)
602 if (do_nhr) net_heat_rate(i) = net_heat_rate(i) + (scale * (w_m2_to_h_t)) * fluxes%heat_added(i,j)
607 if (useriverheatcontent)
then
609 net_heat(i) = (net_heat(i) + (scale*(dt_in_t * rzcp_to_h)) * fluxes%heat_content_lrunoff(i,j)) - &
610 (gv%RZ_to_H * (scale * dt_in_t)) * fluxes%lrunoff(i,j) * t(i,1)
616 if (calculate_diags .and.
associated(tv%TempxPmE))
then
617 tv%TempxPmE(i,j) = tv%TempxPmE(i,j) + (scale * dt_in_t) * &
618 (i_cp*fluxes%heat_content_lrunoff(i,j) - fluxes%lrunoff(i,j)*t(i,1))
624 if (usecalvingheatcontent)
then
626 net_heat(i) = net_heat(i) + (scale*(dt_in_t * rzcp_to_h)) * fluxes%heat_content_frunoff(i,j) - &
627 (gv%RZ_to_H * (scale * dt_in_t)) * fluxes%frunoff(i,j) * t(i,1)
633 if (calculate_diags .and.
associated(tv%TempxPmE))
then
634 tv%TempxPmE(i,j) = tv%TempxPmE(i,j) + (scale * dt_in_t) * &
635 (i_cp*fluxes%heat_content_frunoff(i,j) - fluxes%frunoff(i,j)*t(i,1))
655 if (fluxes%num_msg < fluxes%max_msg)
then
656 if (pen_sw_tot(i) > 1.000001 * w_m2_to_h_t*scale*dt_in_t*fluxes%sw(i,j))
then
657 fluxes%num_msg = fluxes%num_msg + 1
658 write(mesg,
'("Penetrating shortwave of ",1pe17.10, &
659 &" exceeds total shortwave of ",1pe17.10,&
660 &" at ",1pg11.4,"E, "1pg11.4,"N.")') &
661 pen_sw_tot(i),w_m2_to_h_t*scale*dt_in_t * fluxes%sw(i,j),&
662 g%geoLonT(i,j),g%geoLatT(i,j)
669 net_heat(i) = net_heat(i) - pen_sw_tot(i)
671 if (do_nhr) net_heat_rate(i) = net_heat_rate(i) - pen_sw_tot_rate(i)
674 if (
present(nonpensw))
then
675 nonpensw(i) = scale * dt_in_t * w_m2_to_h_t * fluxes%sw(i,j) - pen_sw_tot(i)
680 if (do_nsr) net_salt_rate(i) = 0.0
684 if (
associated(fluxes%salt_flux))
then
685 net_salt(i) = (scale * dt_in_t * (1000.0 * fluxes%salt_flux(i,j))) * gv%RZ_to_H
687 if (do_nsr) net_salt_rate(i) = (scale * 1. * (1000.0 * fluxes%salt_flux(i,j))) * gv%RZ_to_H
691 if (calculate_diags)
then
694 if (
associated(fluxes%salt_flux))
then
696 if (calculate_diags) fluxes%netSalt(i,j) = us%kg_m3_to_R*us%m_to_Z*us%T_to_s*net_salt(i)
701 if (
associated(fluxes%heat_content_massin))
then
702 if (aggregate_fw)
then
703 if (netmassinout(i) > 0.0)
then
704 fluxes%heat_content_massin(i,j) = -fluxes%C_p * netmassout(i) * t(i,1) * gv%H_to_RZ / dt_in_t
706 fluxes%heat_content_massin(i,j) = fluxes%C_p * ( netmassinout(i) - netmassout(i) ) * &
707 t(i,1) * gv%H_to_RZ / dt_in_t
710 fluxes%heat_content_massin(i,j) = 0.
716 if (
associated(fluxes%heat_content_massout))
then
717 if (aggregate_fw)
then
718 if (netmassinout(i) > 0.0)
then
719 fluxes%heat_content_massout(i,j) = fluxes%C_p * netmassout(i) * t(i,1) * gv%H_to_RZ / dt_in_t
721 fluxes%heat_content_massout(i,j) = -fluxes%C_p * ( netmassinout(i) - netmassout(i) ) * &
722 t(i,1) * gv%H_to_RZ / dt_in_t
725 fluxes%heat_content_massout(i,j) = 0.0
736 if (
associated(fluxes%heat_content_lprec))
then
737 if (fluxes%lprec(i,j) > 0.0)
then
738 fluxes%heat_content_lprec(i,j) = fluxes%C_p*fluxes%lprec(i,j)*t(i,1)
740 fluxes%heat_content_lprec(i,j) = 0.0
747 if (
associated(fluxes%heat_content_fprec))
then
748 if (fluxes%fprec(i,j) > 0.0)
then
749 fluxes%heat_content_fprec(i,j) = fluxes%C_p*fluxes%fprec(i,j)*t(i,1)
751 fluxes%heat_content_fprec(i,j) = 0.0
756 if (
associated(fluxes%heat_content_icemelt))
then
757 if (fluxes%seaice_melt(i,j) > 0.0)
then
758 fluxes%heat_content_icemelt(i,j) = fluxes%C_p*fluxes%seaice_melt(i,j)*t(i,1)
760 fluxes%heat_content_icemelt(i,j) = 0.0
767 if (
associated(fluxes%heat_content_vprec))
then
768 if (fluxes%vprec(i,j) > 0.0)
then
769 fluxes%heat_content_vprec(i,j) = fluxes%C_p*fluxes%vprec(i,j)*t(i,1)
771 fluxes%heat_content_vprec(i,j) = 0.0
781 if (
associated(fluxes%heat_content_cond))
then
782 if (fluxes%evap(i,j) > 0.0)
then
783 fluxes%heat_content_cond(i,j) = fluxes%C_p*fluxes%evap(i,j)*t(i,1)
785 fluxes%heat_content_cond(i,j) = 0.0
790 if (.not. useriverheatcontent)
then
791 if (
associated(fluxes%lrunoff) .and.
associated(fluxes%heat_content_lrunoff))
then
792 fluxes%heat_content_lrunoff(i,j) = fluxes%C_p*fluxes%lrunoff(i,j)*t(i,1)
797 if (.not. usecalvingheatcontent)
then
798 if (
associated(fluxes%frunoff) .and.
associated(fluxes%heat_content_frunoff))
then
799 fluxes%heat_content_frunoff(i,j) = fluxes%C_p*fluxes%frunoff(i,j)*t(i,1)
813 subroutine extractfluxes2d(G, GV, US, fluxes, optics, nsw, dt_in_T, FluxRescaleDepth, &
814 useRiverHeatContent, useCalvingHeatContent, h, T, &
815 netMassInOut, netMassOut, net_heat, Net_salt, Pen_SW_bnd, tv, &
821 type(
forcing),
intent(inout) :: fluxes
823 integer,
intent(in) :: nsw
824 real,
intent(in) :: dt_in_t
825 real,
intent(in) :: fluxrescaledepth
827 logical,
intent(in) :: useriverheatcontent
828 logical,
intent(in) :: usecalvingheatcontent
829 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)), &
831 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)), &
833 real,
dimension(SZI_(G),SZJ_(G)),
intent(out) :: netmassinout
836 real,
dimension(SZI_(G),SZJ_(G)),
intent(out) :: netmassout
839 real,
dimension(SZI_(G),SZJ_(G)),
intent(out) :: net_heat
846 real,
dimension(SZI_(G),SZJ_(G)),
intent(out) :: net_salt
848 real,
dimension(max(1,nsw),G%isd:G%ied,G%jsd:G%jed),
intent(out) :: pen_sw_bnd
856 logical,
intent(in) :: aggregate_fw
865 fluxrescaledepth, useriverheatcontent, usecalvingheatcontent,&
866 h(:,j,:), t(:,j,:), netmassinout(:,j), netmassout(:,j), &
867 net_heat(:,j), net_salt(:,j), pen_sw_bnd(:,:,j), tv, aggregate_fw)
877 subroutine calculatebuoyancyflux1d(G, GV, US, fluxes, optics, nsw, h, Temp, Salt, tv, j, &
878 buoyancyFlux, netHeatMinusSW, netSalt, skip_diags)
882 type(
forcing),
intent(inout) :: fluxes
884 integer,
intent(in) :: nsw
886 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)),
intent(in) :: h
887 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)),
intent(in) :: temp
888 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)),
intent(in) :: salt
890 integer,
intent(in) :: j
891 real,
dimension(SZI_(G),SZK_(G)+1),
intent(inout) :: buoyancyflux
892 real,
dimension(SZI_(G)),
intent(inout) :: netheatminussw
894 real,
dimension(SZI_(G)),
intent(inout) :: netsalt
896 logical,
optional,
intent(in) :: skip_diags
899 integer :: start, npts, k
900 real,
parameter :: dt = 1.
901 real,
dimension(SZI_(G)) :: neth
902 real,
dimension(SZI_(G)) :: netevap
904 real,
dimension(SZI_(G)) :: netheat
905 real,
dimension(max(nsw,1), SZI_(G)) :: penswbnd
907 real,
dimension(SZI_(G)) :: pressure
908 real,
dimension(SZI_(G)) :: drhodt
909 real,
dimension(SZI_(G)) :: drhods
910 real,
dimension(SZI_(G),SZK_(G)+1) :: netpen
913 logical :: useriverheatcontent
914 logical :: usecalvingheatcontent
915 real :: depthbeforescalingfluxes
918 real :: h_limit_fluxes
921 useriverheatcontent = .false.
922 usecalvingheatcontent = .false.
924 depthbeforescalingfluxes = max( gv%Angstrom_H, 1.e-30*gv%m_to_H )
926 gorho = (gv%g_Earth * gv%H_to_Z*us%T_to_s) / gv%Rho0
927 start = 1 + g%isc - g%isd
928 npts = 1 + g%iec - g%isc
930 h_limit_fluxes = depthbeforescalingfluxes
939 call extractfluxes1d(g, gv, us, fluxes, optics, nsw, j, dt*us%s_to_T, &
940 depthbeforescalingfluxes, useriverheatcontent, usecalvingheatcontent, &
941 h(:,j,:), temp(:,j,:), neth, netevap, netheatminussw, &
942 netsalt, penswbnd, tv, .false., skip_diags=skip_diags)
946 call sumswoverbands(g, gv, us, h(:,j,:),
optics_nbands(optics), optics, j, dt*us%s_to_T, &
947 h_limit_fluxes, .true., penswbnd, netpen)
951 drhodt, drhods, start, npts, tv%eqn_of_state, scale=us%kg_m3_to_R)
954 netsalt(g%isc:g%iec) = netsalt(g%isc:g%iec) - salt(g%isc:g%iec,j,1) * neth(g%isc:g%iec)
959 netheat(g%isc:g%iec) = netheatminussw(g%isc:g%iec) + netpen(g%isc:g%iec,1)
962 buoyancyflux(g%isc:g%iec,1) = - gorho * ( drhods(g%isc:g%iec) * netsalt(g%isc:g%iec) + &
963 drhodt(g%isc:g%iec) * netheat(g%isc:g%iec) )
966 buoyancyflux(g%isc:g%iec,k) = - gorho * ( drhodt(g%isc:g%iec) * netpen(g%isc:g%iec,k) )
975 buoyancyFlux, netHeatMinusSW, netSalt, skip_diags)
979 type(
forcing),
intent(inout) :: fluxes
981 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)),
intent(in) :: h
982 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)),
intent(in) :: temp
983 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)),
intent(in) :: salt
985 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)+1),
intent(inout) :: buoyancyflux
986 real,
dimension(SZI_(G),SZJ_(G)),
optional,
intent(inout) :: netheatminussw
988 real,
dimension(SZI_(G),SZJ_(G)),
optional,
intent(inout) :: netsalt
990 logical,
optional,
intent(in) :: skip_diags
993 real,
dimension( SZI_(G) ) :: nett
994 real,
dimension( SZI_(G) ) :: nets
997 nett(g%isc:g%iec) = 0. ; nets(g%isc:g%iec) = 0.
1002 tv, j, buoyancyflux(:,j,:), nett, nets, skip_diags=skip_diags)
1003 if (
present(netheatminussw)) netheatminussw(g%isc:g%iec,j) = nett(g%isc:g%iec)
1004 if (
present(netsalt)) netsalt(g%isc:g%iec,j) = nets(g%isc:g%iec)
1012 character(len=*),
intent(in) :: mesg
1013 type(
forcing),
intent(in) :: fluxes
1016 integer,
optional,
intent(in) :: haloshift
1018 real :: rz_t_conversion
1019 integer :: is, ie, js, je, nz, hshift
1020 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
1022 hshift = 1 ;
if (
present(haloshift)) hshift = haloshift
1023 rz_t_conversion = us%R_to_kg_m3*us%Z_to_m*us%s_to_T
1028 if (
associated(fluxes%ustar)) &
1029 call hchksum(fluxes%ustar, mesg//
" fluxes%ustar",g%HI, haloshift=hshift, scale=us%Z_to_m*us%s_to_T)
1030 if (
associated(fluxes%buoy)) &
1031 call hchksum(fluxes%buoy, mesg//
" fluxes%buoy ",g%HI, haloshift=hshift, scale=us%L_to_m**2*us%s_to_T**3)
1032 if (
associated(fluxes%sw)) &
1033 call hchksum(fluxes%sw, mesg//
" fluxes%sw",g%HI,haloshift=hshift)
1034 if (
associated(fluxes%sw_vis_dir)) &
1035 call hchksum(fluxes%sw_vis_dir, mesg//
" fluxes%sw_vis_dir",g%HI,haloshift=hshift)
1036 if (
associated(fluxes%sw_vis_dif)) &
1037 call hchksum(fluxes%sw_vis_dif, mesg//
" fluxes%sw_vis_dif",g%HI,haloshift=hshift)
1038 if (
associated(fluxes%sw_nir_dir)) &
1039 call hchksum(fluxes%sw_nir_dir, mesg//
" fluxes%sw_nir_dir",g%HI,haloshift=hshift)
1040 if (
associated(fluxes%sw_nir_dif)) &
1041 call hchksum(fluxes%sw_nir_dif, mesg//
" fluxes%sw_nir_dif",g%HI,haloshift=hshift)
1042 if (
associated(fluxes%lw)) &
1043 call hchksum(fluxes%lw, mesg//
" fluxes%lw",g%HI,haloshift=hshift)
1044 if (
associated(fluxes%latent)) &
1045 call hchksum(fluxes%latent, mesg//
" fluxes%latent",g%HI,haloshift=hshift)
1046 if (
associated(fluxes%latent_evap_diag)) &
1047 call hchksum(fluxes%latent_evap_diag, mesg//
" fluxes%latent_evap_diag",g%HI,haloshift=hshift)
1048 if (
associated(fluxes%latent_fprec_diag)) &
1049 call hchksum(fluxes%latent_fprec_diag, mesg//
" fluxes%latent_fprec_diag",g%HI,haloshift=hshift)
1050 if (
associated(fluxes%latent_frunoff_diag)) &
1051 call hchksum(fluxes%latent_frunoff_diag, mesg//
" fluxes%latent_frunoff_diag",g%HI,haloshift=hshift)
1052 if (
associated(fluxes%sens)) &
1053 call hchksum(fluxes%sens, mesg//
" fluxes%sens",g%HI,haloshift=hshift)
1054 if (
associated(fluxes%evap)) &
1055 call hchksum(fluxes%evap, mesg//
" fluxes%evap",g%HI,haloshift=hshift, scale=rz_t_conversion)
1056 if (
associated(fluxes%lprec)) &
1057 call hchksum(fluxes%lprec, mesg//
" fluxes%lprec",g%HI,haloshift=hshift, scale=rz_t_conversion)
1058 if (
associated(fluxes%fprec)) &
1059 call hchksum(fluxes%fprec, mesg//
" fluxes%fprec",g%HI,haloshift=hshift, scale=rz_t_conversion)
1060 if (
associated(fluxes%vprec)) &
1061 call hchksum(fluxes%vprec, mesg//
" fluxes%vprec",g%HI,haloshift=hshift, scale=rz_t_conversion)
1062 if (
associated(fluxes%seaice_melt)) &
1063 call hchksum(fluxes%seaice_melt, mesg//
" fluxes%seaice_melt",g%HI,haloshift=hshift, scale=rz_t_conversion)
1064 if (
associated(fluxes%seaice_melt_heat)) &
1065 call hchksum(fluxes%seaice_melt_heat, mesg//
" fluxes%seaice_melt_heat",g%HI,haloshift=hshift)
1066 if (
associated(fluxes%p_surf)) &
1067 call hchksum(fluxes%p_surf, mesg//
" fluxes%p_surf",g%HI,haloshift=hshift)
1068 if (
associated(fluxes%salt_flux)) &
1069 call hchksum(fluxes%salt_flux, mesg//
" fluxes%salt_flux",g%HI,haloshift=hshift, scale=rz_t_conversion)
1070 if (
associated(fluxes%TKE_tidal)) &
1071 call hchksum(fluxes%TKE_tidal, mesg//
" fluxes%TKE_tidal",g%HI,haloshift=hshift, &
1072 scale=us%R_to_kg_m3**3*us%Z_to_m**3*us%s_to_T)
1073 if (
associated(fluxes%ustar_tidal)) &
1074 call hchksum(fluxes%ustar_tidal, mesg//
" fluxes%ustar_tidal",g%HI,haloshift=hshift, scale=us%Z_to_m*us%s_to_T)
1075 if (
associated(fluxes%lrunoff)) &
1076 call hchksum(fluxes%lrunoff, mesg//
" fluxes%lrunoff",g%HI,haloshift=hshift, scale=rz_t_conversion)
1077 if (
associated(fluxes%frunoff)) &
1078 call hchksum(fluxes%frunoff, mesg//
" fluxes%frunoff",g%HI,haloshift=hshift, scale=rz_t_conversion)
1079 if (
associated(fluxes%heat_content_lrunoff)) &
1080 call hchksum(fluxes%heat_content_lrunoff, mesg//
" fluxes%heat_content_lrunoff", g%HI, &
1081 haloshift=hshift, scale=rz_t_conversion)
1082 if (
associated(fluxes%heat_content_frunoff)) &
1083 call hchksum(fluxes%heat_content_frunoff, mesg//
" fluxes%heat_content_frunoff", g%HI, &
1084 haloshift=hshift, scale=rz_t_conversion)
1085 if (
associated(fluxes%heat_content_lprec)) &
1086 call hchksum(fluxes%heat_content_lprec, mesg//
" fluxes%heat_content_lprec", g%HI, &
1087 haloshift=hshift, scale=rz_t_conversion)
1088 if (
associated(fluxes%heat_content_fprec)) &
1089 call hchksum(fluxes%heat_content_fprec, mesg//
" fluxes%heat_content_fprec", g%HI, &
1090 haloshift=hshift, scale=rz_t_conversion)
1091 if (
associated(fluxes%heat_content_icemelt)) &
1092 call hchksum(fluxes%heat_content_icemelt, mesg//
" fluxes%heat_content_icemelt", g%HI, &
1093 haloshift=hshift, scale=rz_t_conversion)
1094 if (
associated(fluxes%heat_content_cond)) &
1095 call hchksum(fluxes%heat_content_cond, mesg//
" fluxes%heat_content_cond", g%HI, &
1096 haloshift=hshift, scale=rz_t_conversion)
1097 if (
associated(fluxes%heat_content_massout)) &
1098 call hchksum(fluxes%heat_content_massout, mesg//
" fluxes%heat_content_massout", g%HI, &
1099 haloshift=hshift, scale=rz_t_conversion)
1104 character(len=*),
intent(in) :: mesg
1108 integer,
optional,
intent(in) :: haloshift
1110 integer :: is, ie, js, je, nz, hshift
1111 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
1113 hshift=1;
if (
present(haloshift)) hshift=haloshift
1118 if (
associated(forces%taux) .and.
associated(forces%tauy)) &
1119 call uvchksum(mesg//
" forces%tau[xy]", forces%taux, forces%tauy, g%HI, &
1120 haloshift=hshift, symmetric=.true., scale=us%R_to_kg_m3*us%L_T_to_m_s**2*us%Z_to_L)
1121 if (
associated(forces%p_surf)) &
1122 call hchksum(forces%p_surf, mesg//
" forces%p_surf",g%HI,haloshift=hshift)
1123 if (
associated(forces%ustar)) &
1124 call hchksum(forces%ustar, mesg//
" forces%ustar",g%HI,haloshift=hshift, scale=us%Z_to_m*us%s_to_T)
1125 if (
associated(forces%rigidity_ice_u) .and.
associated(forces%rigidity_ice_v)) &
1126 call uvchksum(mesg//
" forces%rigidity_ice_[uv]", forces%rigidity_ice_u, &
1127 forces%rigidity_ice_v, g%HI, haloshift=hshift, symmetric=.true.)
1135 character(len=*),
intent(in) :: mesg
1136 integer,
intent(in) :: i
1137 integer,
intent(in) :: j
1139 write(0,
'(2a)')
'MOM_forcing_type, forcing_SinglePointPrint: Called from ',mesg
1140 write(0,
'(a,2es15.3)')
'MOM_forcing_type, forcing_SinglePointPrint: lon,lat = ',g%geoLonT(i,j),g%geoLatT(i,j)
1141 call locmsg(forces%taux,
'taux')
1142 call locmsg(forces%tauy,
'tauy')
1146 subroutine locmsg(array,aname)
1147 real,
dimension(:,:),
pointer :: array
1148 character(len=*) :: aname
1150 if (
associated(array))
then
1151 write(0,
'(3a,es15.3)')
'MOM_forcing_type, mech_forcing_SinglePointPrint: ',trim(aname),
' = ',array(i,j)
1153 write(0,
'(4a)')
'MOM_forcing_type, mech_forcing_SinglePointPrint: ',trim(aname),
' is not associated.'
1163 character(len=*),
intent(in) :: mesg
1164 integer,
intent(in) :: i
1165 integer,
intent(in) :: j
1167 write(0,
'(2a)')
'MOM_forcing_type, forcing_SinglePointPrint: Called from ',mesg
1168 write(0,
'(a,2es15.3)')
'MOM_forcing_type, forcing_SinglePointPrint: lon,lat = ',g%geoLonT(i,j),g%geoLatT(i,j)
1169 call locmsg(fluxes%ustar,
'ustar')
1170 call locmsg(fluxes%buoy,
'buoy')
1171 call locmsg(fluxes%sw,
'sw')
1172 call locmsg(fluxes%sw_vis_dir,
'sw_vis_dir')
1173 call locmsg(fluxes%sw_vis_dif,
'sw_vis_dif')
1174 call locmsg(fluxes%sw_nir_dir,
'sw_nir_dir')
1175 call locmsg(fluxes%sw_nir_dif,
'sw_nir_dif')
1176 call locmsg(fluxes%lw,
'lw')
1177 call locmsg(fluxes%latent,
'latent')
1178 call locmsg(fluxes%latent_evap_diag,
'latent_evap_diag')
1179 call locmsg(fluxes%latent_fprec_diag,
'latent_fprec_diag')
1180 call locmsg(fluxes%latent_frunoff_diag,
'latent_frunoff_diag')
1181 call locmsg(fluxes%sens,
'sens')
1182 call locmsg(fluxes%evap,
'evap')
1183 call locmsg(fluxes%lprec,
'lprec')
1184 call locmsg(fluxes%fprec,
'fprec')
1185 call locmsg(fluxes%vprec,
'vprec')
1186 call locmsg(fluxes%seaice_melt,
'seaice_melt')
1187 call locmsg(fluxes%seaice_melt_heat,
'seaice_melt_heat')
1188 call locmsg(fluxes%p_surf,
'p_surf')
1189 call locmsg(fluxes%salt_flux,
'salt_flux')
1190 call locmsg(fluxes%TKE_tidal,
'TKE_tidal')
1191 call locmsg(fluxes%ustar_tidal,
'ustar_tidal')
1192 call locmsg(fluxes%lrunoff,
'lrunoff')
1193 call locmsg(fluxes%frunoff,
'frunoff')
1194 call locmsg(fluxes%heat_content_lrunoff,
'heat_content_lrunoff')
1195 call locmsg(fluxes%heat_content_frunoff,
'heat_content_frunoff')
1196 call locmsg(fluxes%heat_content_lprec,
'heat_content_lprec')
1197 call locmsg(fluxes%heat_content_fprec,
'heat_content_fprec')
1198 call locmsg(fluxes%heat_content_icemelt,
'heat_content_icemelt')
1199 call locmsg(fluxes%heat_content_vprec,
'heat_content_vprec')
1200 call locmsg(fluxes%heat_content_cond,
'heat_content_cond')
1201 call locmsg(fluxes%heat_content_cond,
'heat_content_massout')
1205 subroutine locmsg(array,aname)
1206 real,
dimension(:,:),
pointer :: array
1207 character(len=*) :: aname
1209 if (
associated(array))
then
1210 write(0,
'(3a,es15.3)')
'MOM_forcing_type, forcing_SinglePointPrint: ',trim(aname),
' = ',array(i,j)
1212 write(0,
'(4a)')
'MOM_forcing_type, forcing_SinglePointPrint: ',trim(aname),
' is not associated.'
1221 type(time_type),
intent(in) :: time
1224 logical,
intent(in) :: use_temperature
1226 logical,
optional,
intent(in) :: use_berg_fluxes
1229 handles%id_clock_forcing=cpu_clock_id(
'(Ocean forcing diagnostics)', grain=clock_routine)
1232 handles%id_taux = register_diag_field(
'ocean_model',
'taux', diag%axesCu1, time, &
1233 'Zonal surface stress from ocean interactions with atmos and ice', &
1234 'Pa', conversion=us%R_to_kg_m3*us%L_T_to_m_s**2*us%Z_to_L, &
1235 standard_name=
'surface_downward_x_stress', cmor_field_name=
'tauuo', &
1236 cmor_units=
'N m-2', cmor_long_name=
'Surface Downward X Stress', &
1237 cmor_standard_name=
'surface_downward_x_stress')
1239 handles%id_tauy = register_diag_field(
'ocean_model',
'tauy', diag%axesCv1, time, &
1240 'Meridional surface stress ocean interactions with atmos and ice', &
1241 'Pa', conversion=us%R_to_kg_m3*us%L_T_to_m_s**2*us%Z_to_L, &
1242 standard_name=
'surface_downward_y_stress', cmor_field_name=
'tauvo', &
1243 cmor_units=
'N m-2', cmor_long_name=
'Surface Downward Y Stress', &
1244 cmor_standard_name=
'surface_downward_y_stress')
1246 handles%id_ustar = register_diag_field(
'ocean_model',
'ustar', diag%axesT1, time, &
1247 'Surface friction velocity = [(gustiness + tau_magnitude)/rho0]^(1/2)', &
1248 'm s-1', conversion=us%Z_to_m*us%s_to_T)
1250 if (
present(use_berg_fluxes))
then
1251 if (use_berg_fluxes)
then
1252 handles%id_ustar_berg = register_diag_field(
'ocean_model',
'ustar_berg', diag%axesT1, time, &
1253 'Friction velocity below iceberg ',
'm s-1', conversion=us%Z_to_m*us%s_to_T)
1255 handles%id_area_berg = register_diag_field(
'ocean_model',
'area_berg', diag%axesT1, time, &
1256 'Area of grid cell covered by iceberg ',
'm2 m-2')
1258 handles%id_mass_berg = register_diag_field(
'ocean_model',
'mass_berg', diag%axesT1, time, &
1259 'Mass of icebergs ',
'kg m-2')
1261 handles%id_ustar_ice_cover = register_diag_field(
'ocean_model',
'ustar_ice_cover', diag%axesT1, time, &
1262 'Friction velocity below iceberg and ice shelf together',
'm s-1', conversion=us%Z_to_m*us%s_to_T)
1264 handles%id_frac_ice_cover = register_diag_field(
'ocean_model',
'frac_ice_cover', diag%axesT1, time, &
1265 'Area of grid cell below iceberg and ice shelf together ',
'm2 m-2')
1269 handles%id_psurf = register_diag_field(
'ocean_model',
'p_surf', diag%axesT1, time, &
1270 'Pressure at ice-ocean or atmosphere-ocean interface',
'Pa', cmor_field_name=
'pso', &
1271 cmor_long_name=
'Sea Water Pressure at Sea Water Surface', &
1272 cmor_standard_name=
'sea_water_pressure_at_sea_water_surface')
1274 handles%id_TKE_tidal = register_diag_field(
'ocean_model',
'TKE_tidal', diag%axesT1, time, &
1275 'Tidal source of BBL mixing',
'W m-2', conversion=us%R_to_kg_m3*us%Z_to_m**3*us%s_to_T**3)
1277 if (.not. use_temperature)
then
1278 handles%id_buoy = register_diag_field(
'ocean_model',
'buoy', diag%axesT1, time, &
1279 'Buoyancy forcing',
'm2 s-3', conversion=us%L_to_m**2*us%s_to_T**3)
1287 handles%id_prcme = register_diag_field(
'ocean_model',
'PRCmE', diag%axesT1, time, &
1288 'Net surface water flux (precip+melt+lrunoff+ice calving-evap)',
'kg m-2 s-1', &
1289 standard_name=
'water_flux_into_sea_water', cmor_field_name=
'wfo', &
1290 cmor_standard_name=
'water_flux_into_sea_water',cmor_long_name=
'Water Flux Into Sea Water')
1293 handles%id_evap = register_diag_field(
'ocean_model',
'evap', diag%axesT1, time, &
1294 'Evaporation/condensation at ocean surface (evaporation is negative)', &
1295 'kg m-2 s-1', conversion=us%R_to_kg_m3*us%Z_to_m*us%s_to_T, &
1296 standard_name=
'water_evaporation_flux', cmor_field_name=
'evs', &
1297 cmor_standard_name=
'water_evaporation_flux', &
1298 cmor_long_name=
'Water Evaporation Flux Where Ice Free Ocean over Sea')
1301 handles%id_seaice_melt = register_diag_field(
'ocean_model',
'seaice_melt', &
1302 diag%axesT1, time,
'water flux to ocean from snow/sea ice melting(> 0) or formation(< 0)', &
1303 'kg m-2 s-1', conversion=us%R_to_kg_m3*us%Z_to_m*us%s_to_T, &
1304 standard_name=
'water_flux_into_sea_water_due_to_sea_ice_thermodynamics', &
1305 cmor_field_name=
'fsitherm', &
1306 cmor_standard_name=
'water_flux_into_sea_water_due_to_sea_ice_thermodynamics',&
1307 cmor_long_name=
'water flux to ocean from sea ice melt(> 0) or form(< 0)')
1309 handles%id_precip = register_diag_field(
'ocean_model',
'precip', diag%axesT1, time, &
1310 'Liquid + frozen precipitation into ocean',
'kg m-2 s-1')
1313 handles%id_fprec = register_diag_field(
'ocean_model',
'fprec', diag%axesT1, time, &
1314 'Frozen precipitation into ocean', &
1315 units=
'kg m-2 s-1', conversion=us%R_to_kg_m3*us%Z_to_m*us%s_to_T, &
1316 standard_name=
'snowfall_flux', cmor_field_name=
'prsn', &
1317 cmor_standard_name=
'snowfall_flux', cmor_long_name=
'Snowfall Flux where Ice Free Ocean over Sea')
1319 handles%id_lprec = register_diag_field(
'ocean_model',
'lprec', diag%axesT1, time, &
1320 'Liquid precipitation into ocean', &
1321 units=
'kg m-2 s-1', conversion=us%R_to_kg_m3*us%Z_to_m*us%s_to_T, &
1322 standard_name=
'rainfall_flux', &
1323 cmor_field_name=
'prlq', cmor_standard_name=
'rainfall_flux', &
1324 cmor_long_name=
'Rainfall Flux where Ice Free Ocean over Sea')
1326 handles%id_vprec = register_diag_field(
'ocean_model',
'vprec', diag%axesT1, time, &
1327 'Virtual liquid precip into ocean due to SSS restoring', &
1328 units=
'kg m-2 s-1', conversion=us%R_to_kg_m3*us%Z_to_m*us%s_to_T)
1330 handles%id_frunoff = register_diag_field(
'ocean_model',
'frunoff', diag%axesT1, time, &
1331 'Frozen runoff (calving) and iceberg melt into ocean', &
1332 units=
'kg m-2 s-1', conversion=us%R_to_kg_m3*us%Z_to_m*us%s_to_T, &
1333 standard_name=
'water_flux_into_sea_water_from_icebergs', &
1334 cmor_field_name=
'ficeberg', &
1335 cmor_standard_name=
'water_flux_into_sea_water_from_icebergs', &
1336 cmor_long_name=
'Water Flux into Seawater from Icebergs')
1338 handles%id_lrunoff = register_diag_field(
'ocean_model',
'lrunoff', diag%axesT1, time, &
1339 'Liquid runoff (rivers) into ocean', &
1340 units=
'kg m-2 s-1', conversion=us%R_to_kg_m3*us%Z_to_m*us%s_to_T, &
1341 standard_name=
'water_flux_into_sea_water_from_rivers', cmor_field_name=
'friver', &
1342 cmor_standard_name=
'water_flux_into_sea_water_from_rivers', &
1343 cmor_long_name=
'Water Flux into Sea Water From Rivers')
1345 handles%id_net_massout = register_diag_field(
'ocean_model',
'net_massout', diag%axesT1, time, &
1346 'Net mass leaving the ocean due to evaporation, seaice formation',
'kg m-2 s-1')
1349 handles%id_net_massin = register_diag_field(
'ocean_model',
'net_massin', diag%axesT1, time, &
1350 'Net mass entering ocean due to precip, runoff, ice melt',
'kg m-2 s-1')
1353 handles%id_massout_flux = register_diag_field(
'ocean_model',
'massout_flux', diag%axesT1, time, &
1354 'Net mass flux of freshwater out of the ocean (used in the boundary flux calculation)', &
1355 'kg m-2', conversion=diag%GV%H_to_kg_m2)
1358 handles%id_massin_flux = register_diag_field(
'ocean_model',
'massin_flux', diag%axesT1, time, &
1359 'Net mass flux of freshwater into the ocean (used in boundary flux calculation)',
'kg m-2')
1365 handles%id_total_prcme = register_scalar_field(
'ocean_model',
'total_PRCmE', time, diag, &
1366 long_name=
'Area integrated net surface water flux (precip+melt+liq runoff+ice calving-evap)',&
1367 units=
'kg s-1', standard_name=
'water_flux_into_sea_water_area_integrated', &
1368 cmor_field_name=
'total_wfo', &
1369 cmor_standard_name=
'water_flux_into_sea_water_area_integrated', &
1370 cmor_long_name=
'Water Transport Into Sea Water Area Integrated')
1372 handles%id_total_evap = register_scalar_field(
'ocean_model',
'total_evap', time, diag,&
1373 long_name=
'Area integrated evap/condense at ocean surface', &
1374 units=
'kg s-1', standard_name=
'water_evaporation_flux_area_integrated', &
1375 cmor_field_name=
'total_evs', &
1376 cmor_standard_name=
'water_evaporation_flux_area_integrated', &
1377 cmor_long_name=
'Evaporation Where Ice Free Ocean over Sea Area Integrated')
1380 handles%id_total_seaice_melt = register_scalar_field(
'ocean_model',
'total_icemelt', time, diag, &
1381 long_name=
'Area integrated sea ice melt (>0) or form (<0)', units=
'kg s-1', &
1382 standard_name=
'water_flux_into_sea_water_due_to_sea_ice_thermodynamics_area_integrated', &
1383 cmor_field_name=
'total_fsitherm', &
1384 cmor_standard_name=
'water_flux_into_sea_water_due_to_sea_ice_thermodynamics_area_integrated', &
1385 cmor_long_name=
'Water Melt/Form from Sea Ice Area Integrated')
1387 handles%id_total_precip = register_scalar_field(
'ocean_model',
'total_precip', time, diag, &
1388 long_name=
'Area integrated liquid+frozen precip into ocean', units=
'kg s-1')
1390 handles%id_total_fprec = register_scalar_field(
'ocean_model',
'total_fprec', time, diag,&
1391 long_name=
'Area integrated frozen precip into ocean', units=
'kg s-1', &
1392 standard_name=
'snowfall_flux_area_integrated', &
1393 cmor_field_name=
'total_prsn', &
1394 cmor_standard_name=
'snowfall_flux_area_integrated', &
1395 cmor_long_name=
'Snowfall Flux where Ice Free Ocean over Sea Area Integrated')
1397 handles%id_total_lprec = register_scalar_field(
'ocean_model',
'total_lprec', time, diag,&
1398 long_name=
'Area integrated liquid precip into ocean', units=
'kg s-1', &
1399 standard_name=
'rainfall_flux_area_integrated', &
1400 cmor_field_name=
'total_pr', &
1401 cmor_standard_name=
'rainfall_flux_area_integrated', &
1402 cmor_long_name=
'Rainfall Flux where Ice Free Ocean over Sea Area Integrated')
1404 handles%id_total_vprec = register_scalar_field(
'ocean_model',
'total_vprec', time, diag, &
1405 long_name=
'Area integrated virtual liquid precip due to SSS restoring', units=
'kg s-1')
1407 handles%id_total_frunoff = register_scalar_field(
'ocean_model',
'total_frunoff', time, diag, &
1408 long_name=
'Area integrated frozen runoff (calving) & iceberg melt into ocean', units=
'kg s-1',&
1409 cmor_field_name=
'total_ficeberg', &
1410 cmor_standard_name=
'water_flux_into_sea_water_from_icebergs_area_integrated', &
1411 cmor_long_name=
'Water Flux into Seawater from Icebergs Area Integrated')
1413 handles%id_total_lrunoff = register_scalar_field(
'ocean_model',
'total_lrunoff', time, diag,&
1414 long_name=
'Area integrated liquid runoff into ocean', units=
'kg s-1', &
1415 cmor_field_name=
'total_friver', &
1416 cmor_standard_name=
'water_flux_into_sea_water_from_rivers_area_integrated', &
1417 cmor_long_name=
'Water Flux into Sea Water From Rivers Area Integrated')
1419 handles%id_total_net_massout = register_scalar_field(
'ocean_model',
'total_net_massout', time, diag, &
1420 long_name=
'Area integrated mass leaving ocean due to evap and seaice form', units=
'kg s-1')
1422 handles%id_total_net_massin = register_scalar_field(
'ocean_model',
'total_net_massin', time, diag, &
1423 long_name=
'Area integrated mass entering ocean due to predip, runoff, ice melt', units=
'kg s-1')
1428 handles%id_prcme_ga = register_scalar_field(
'ocean_model',
'PRCmE_ga', time, diag, &
1429 long_name=
'Area averaged net surface water flux (precip+melt+liq runoff+ice calving-evap)',&
1430 units=
'kg m-2 s-1', standard_name=
'water_flux_into_sea_water_area_averaged', &
1431 cmor_field_name=
'ave_wfo', &
1432 cmor_standard_name=
'rainfall_flux_area_averaged', &
1433 cmor_long_name=
'Water Transport Into Sea Water Area Averaged')
1435 handles%id_evap_ga = register_scalar_field(
'ocean_model',
'evap_ga', time, diag,&
1436 long_name=
'Area averaged evap/condense at ocean surface', &
1437 units=
'kg m-2 s-1', standard_name=
'water_evaporation_flux_area_averaged', &
1438 cmor_field_name=
'ave_evs', &
1439 cmor_standard_name=
'water_evaporation_flux_area_averaged', &
1440 cmor_long_name=
'Evaporation Where Ice Free Ocean over Sea Area Averaged')
1442 handles%id_lprec_ga = register_scalar_field(
'ocean_model',
'lprec_ga', time, diag,&
1443 long_name=
'Area integrated liquid precip into ocean', units=
'kg m-2 s-1', &
1444 standard_name=
'rainfall_flux_area_averaged', &
1445 cmor_field_name=
'ave_pr', &
1446 cmor_standard_name=
'rainfall_flux_area_averaged', &
1447 cmor_long_name=
'Rainfall Flux where Ice Free Ocean over Sea Area Averaged')
1449 handles%id_fprec_ga = register_scalar_field(
'ocean_model',
'fprec_ga', time, diag,&
1450 long_name=
'Area integrated frozen precip into ocean', units=
'kg m-2 s-1', &
1451 standard_name=
'snowfall_flux_area_averaged', &
1452 cmor_field_name=
'ave_prsn', &
1453 cmor_standard_name=
'snowfall_flux_area_averaged', &
1454 cmor_long_name=
'Snowfall Flux where Ice Free Ocean over Sea Area Averaged')
1456 handles%id_precip_ga = register_scalar_field(
'ocean_model',
'precip_ga', time, diag, &
1457 long_name=
'Area averaged liquid+frozen precip into ocean', units=
'kg m-2 s-1')
1459 handles%id_vprec_ga = register_scalar_field(
'ocean_model',
'vrec_ga', time, diag, &
1460 long_name=
'Area averaged virtual liquid precip due to SSS restoring', units=
'kg m-2 s-1')
1465 handles%id_heat_content_frunoff = register_diag_field(
'ocean_model',
'heat_content_frunoff', &
1466 diag%axesT1, time,
'Heat content (relative to 0C) of solid runoff into ocean', &
1467 'W m-2', conversion=us%R_to_kg_m3*us%Z_to_m*us%s_to_T, &
1468 standard_name=
'temperature_flux_due_to_solid_runoff_expressed_as_heat_flux_into_sea_water')
1470 handles%id_heat_content_lrunoff = register_diag_field(
'ocean_model',
'heat_content_lrunoff', &
1471 diag%axesT1, time,
'Heat content (relative to 0C) of liquid runoff into ocean', &
1472 'W m-2', conversion=us%R_to_kg_m3*us%Z_to_m*us%s_to_T, &
1473 standard_name=
'temperature_flux_due_to_runoff_expressed_as_heat_flux_into_sea_water')
1475 handles%id_hfrunoffds = register_diag_field(
'ocean_model',
'hfrunoffds', &
1476 diag%axesT1, time,
'Heat content (relative to 0C) of liquid+solid runoff into ocean', &
1477 'W m-2', conversion=us%R_to_kg_m3*us%Z_to_m*us%s_to_T, &
1478 standard_name=
'temperature_flux_due_to_runoff_expressed_as_heat_flux_into_sea_water')
1480 handles%id_heat_content_lprec = register_diag_field(
'ocean_model',
'heat_content_lprec', &
1481 diag%axesT1,time,
'Heat content (relative to 0degC) of liquid precip entering ocean', &
1482 'W m-2', conversion=us%R_to_kg_m3*us%Z_to_m*us%s_to_T)
1484 handles%id_heat_content_fprec = register_diag_field(
'ocean_model',
'heat_content_fprec',&
1485 diag%axesT1,time,
'Heat content (relative to 0degC) of frozen prec entering ocean',&
1486 'W m-2', conversion=us%R_to_kg_m3*us%Z_to_m*us%s_to_T)
1488 handles%id_heat_content_icemelt = register_diag_field(
'ocean_model',
'heat_content_icemelt',&
1489 diag%axesT1,time,
'Heat content (relative to 0degC) of water flux due to sea ice melting/freezing',&
1490 'W m-2', conversion=us%R_to_kg_m3*us%Z_to_m*us%s_to_T)
1492 handles%id_heat_content_vprec = register_diag_field(
'ocean_model',
'heat_content_vprec', &
1493 diag%axesT1,time,
'Heat content (relative to 0degC) of virtual precip entering ocean',&
1494 'W m-2', conversion=us%R_to_kg_m3*us%Z_to_m*us%s_to_T)
1496 handles%id_heat_content_cond = register_diag_field(
'ocean_model',
'heat_content_cond', &
1497 diag%axesT1,time,
'Heat content (relative to 0degC) of water condensing into ocean',&
1498 'W m-2', conversion=us%R_to_kg_m3*us%Z_to_m*us%s_to_T)
1500 handles%id_hfrainds = register_diag_field(
'ocean_model',
'hfrainds', &
1501 diag%axesT1,time,
'Heat content (relative to 0degC) of liquid+frozen precip entering ocean', &
1502 'W m-2', conversion=us%R_to_kg_m3*us%Z_to_m*us%s_to_T, &
1503 standard_name=
'temperature_flux_due_to_rainfall_expressed_as_heat_flux_into_sea_water',&
1504 cmor_long_name=
'Heat Content (relative to 0degC) of Liquid + Frozen Precipitation')
1506 handles%id_heat_content_surfwater = register_diag_field(
'ocean_model',
'heat_content_surfwater',&
1507 diag%axesT1, time, &
1508 'Heat content (relative to 0degC) of net water crossing ocean surface (frozen+liquid)', &
1509 'W m-2', conversion=us%R_to_kg_m3*us%Z_to_m*us%s_to_T)
1511 handles%id_heat_content_massout = register_diag_field(
'ocean_model',
'heat_content_massout', &
1512 diag%axesT1, time,
'Heat content (relative to 0degC) of net mass leaving ocean ocean via evap and ice form',&
1513 'W m-2', conversion=us%R_to_kg_m3*us%Z_to_m*us%s_to_T, &
1514 cmor_field_name=
'hfevapds', &
1515 cmor_standard_name=
'temperature_flux_due_to_evaporation_expressed_as_heat_flux_out_of_sea_water', &
1516 cmor_long_name=
'Heat Content (relative to 0degC) of Water Leaving Ocean via Evaporation and Ice Formation')
1518 handles%id_heat_content_massin = register_diag_field(
'ocean_model',
'heat_content_massin', &
1519 diag%axesT1, time,
'Heat content (relative to 0degC) of net mass entering ocean ocean',&
1520 'W m-2', conversion=us%R_to_kg_m3*us%Z_to_m*us%s_to_T)
1522 handles%id_net_heat_coupler = register_diag_field(
'ocean_model',
'net_heat_coupler', &
1523 diag%axesT1,time,
'Surface ocean heat flux from SW+LW+latent+sensible+seaice_melt_heat (via the coupler)',&
1526 handles%id_net_heat_surface = register_diag_field(
'ocean_model',
'net_heat_surface',diag%axesT1, time, &
1527 'Surface ocean heat flux from SW+LW+lat+sens+mass transfer+frazil+restore+seaice_melt_heat or '// &
1528 'flux adjustments',&
1530 standard_name=
'surface_downward_heat_flux_in_sea_water', cmor_field_name=
'hfds', &
1531 cmor_standard_name=
'surface_downward_heat_flux_in_sea_water', &
1532 cmor_long_name=
'Surface ocean heat flux from SW+LW+latent+sensible+masstransfer+frazil+seaice_melt_heat')
1534 handles%id_sw = register_diag_field(
'ocean_model',
'SW', diag%axesT1, time, &
1535 'Shortwave radiation flux into ocean',
'W m-2', &
1536 standard_name=
'net_downward_shortwave_flux_at_sea_water_surface', &
1537 cmor_field_name=
'rsntds', &
1538 cmor_standard_name=
'net_downward_shortwave_flux_at_sea_water_surface', &
1539 cmor_long_name=
'Net Downward Shortwave Radiation at Sea Water Surface')
1540 handles%id_sw_vis = register_diag_field(
'ocean_model',
'sw_vis', diag%axesT1, time, &
1541 'Shortwave radiation direct and diffuse flux into the ocean in the visible band', &
1543 handles%id_sw_nir = register_diag_field(
'ocean_model',
'sw_nir', diag%axesT1, time, &
1544 'Shortwave radiation direct and diffuse flux into the ocean in the near-infrared band', &
1547 handles%id_LwLatSens = register_diag_field(
'ocean_model',
'LwLatSens', diag%axesT1, time, &
1548 'Combined longwave, latent, and sensible heating at ocean surface',
'W m-2')
1550 handles%id_lw = register_diag_field(
'ocean_model',
'LW', diag%axesT1, time, &
1551 'Longwave radiation flux into ocean',
'W m-2', &
1552 standard_name=
'surface_net_downward_longwave_flux', &
1553 cmor_field_name=
'rlntds', &
1554 cmor_standard_name=
'surface_net_downward_longwave_flux', &
1555 cmor_long_name=
'Surface Net Downward Longwave Radiation')
1557 handles%id_lat = register_diag_field(
'ocean_model',
'latent', diag%axesT1, time, &
1558 'Latent heat flux into ocean due to fusion and evaporation (negative means ocean heat loss)', &
1559 'W m-2', cmor_field_name=
'hflso', &
1560 cmor_standard_name=
'surface_downward_latent_heat_flux', &
1561 cmor_long_name=
'Surface Downward Latent Heat Flux due to Evap + Melt Snow/Ice')
1563 handles%id_lat_evap = register_diag_field(
'ocean_model',
'latent_evap', diag%axesT1, time, &
1564 'Latent heat flux into ocean due to evaporation/condensation',
'W m-2')
1566 handles%id_lat_fprec = register_diag_field(
'ocean_model',
'latent_fprec_diag', diag%axesT1, time,&
1567 'Latent heat flux into ocean due to melting of frozen precipitation',
'W m-2', &
1568 cmor_field_name=
'hfsnthermds', &
1569 cmor_standard_name=
'heat_flux_into_sea_water_due_to_snow_thermodynamics', &
1570 cmor_long_name=
'Latent Heat to Melt Frozen Precipitation')
1572 handles%id_lat_frunoff = register_diag_field(
'ocean_model',
'latent_frunoff', diag%axesT1, time, &
1573 'Latent heat flux into ocean due to melting of icebergs',
'W m-2', &
1574 cmor_field_name=
'hfibthermds', &
1575 cmor_standard_name=
'heat_flux_into_sea_water_due_to_iceberg_thermodynamics', &
1576 cmor_long_name=
'Latent Heat to Melt Frozen Runoff/Iceberg')
1578 handles%id_sens = register_diag_field(
'ocean_model',
'sensible', diag%axesT1, time,&
1579 'Sensible heat flux into ocean',
'W m-2', &
1580 standard_name=
'surface_downward_sensible_heat_flux', &
1581 cmor_field_name=
'hfsso', &
1582 cmor_standard_name=
'surface_downward_sensible_heat_flux', &
1583 cmor_long_name=
'Surface Downward Sensible Heat Flux')
1585 handles%id_seaice_melt_heat = register_diag_field(
'ocean_model',
'seaice_melt_heat', diag%axesT1, time,&
1586 'Heat flux into ocean due to snow and sea ice melt/freeze',
'W m-2', &
1587 standard_name=
'snow_ice_melt_heat_flux', &
1589 cmor_standard_name=
'snow_ice_melt_heat_flux', &
1590 cmor_long_name=
'Heat flux into ocean from snow and sea ice melt')
1592 handles%id_heat_added = register_diag_field(
'ocean_model',
'heat_added', diag%axesT1, time, &
1593 'Flux Adjustment or restoring surface heat flux into ocean',
'W m-2')
1599 handles%id_total_heat_content_frunoff = register_scalar_field(
'ocean_model', &
1600 'total_heat_content_frunoff', time, diag, &
1601 long_name=
'Area integrated heat content (relative to 0C) of solid runoff', &
1602 units=
'W', cmor_field_name=
'total_hfsolidrunoffds', &
1603 cmor_standard_name= &
1604 'temperature_flux_due_to_solid_runoff_expressed_as_heat_flux_into_sea_water_area_integrated',&
1606 'Temperature Flux due to Solid Runoff Expressed as Heat Flux into Sea Water Area Integrated')
1608 handles%id_total_heat_content_lrunoff = register_scalar_field(
'ocean_model', &
1609 'total_heat_content_lrunoff', time, diag, &
1610 long_name=
'Area integrated heat content (relative to 0C) of liquid runoff', &
1611 units=
'W', cmor_field_name=
'total_hfrunoffds', &
1612 cmor_standard_name= &
1613 'temperature_flux_due_to_runoff_expressed_as_heat_flux_into_sea_water_area_integrated',&
1615 'Temperature Flux due to Runoff Expressed as Heat Flux into Sea Water Area Integrated')
1617 handles%id_total_heat_content_lprec = register_scalar_field(
'ocean_model', &
1618 'total_heat_content_lprec', time, diag, &
1619 long_name=
'Area integrated heat content (relative to 0C) of liquid precip', &
1620 units=
'W', cmor_field_name=
'total_hfrainds', &
1621 cmor_standard_name= &
1622 'temperature_flux_due_to_rainfall_expressed_as_heat_flux_into_sea_water_area_integrated',&
1624 'Temperature Flux due to Rainfall Expressed as Heat Flux into Sea Water Area Integrated')
1626 handles%id_total_heat_content_fprec = register_scalar_field(
'ocean_model', &
1627 'total_heat_content_fprec', time, diag, &
1628 long_name=
'Area integrated heat content (relative to 0C) of frozen precip',&
1631 handles%id_total_heat_content_icemelt = register_scalar_field(
'ocean_model', &
1632 'total_heat_content_icemelt', time, diag,long_name= &
1633 'Area integrated heat content (relative to 0C) of water flux due sea ice melting/freezing', &
1636 handles%id_total_heat_content_vprec = register_scalar_field(
'ocean_model', &
1637 'total_heat_content_vprec', time, diag, &
1638 long_name=
'Area integrated heat content (relative to 0C) of virtual precip',&
1641 handles%id_total_heat_content_cond = register_scalar_field(
'ocean_model', &
1642 'total_heat_content_cond', time, diag, &
1643 long_name=
'Area integrated heat content (relative to 0C) of condensate',&
1646 handles%id_total_heat_content_surfwater = register_scalar_field(
'ocean_model', &
1647 'total_heat_content_surfwater', time, diag, &
1648 long_name=
'Area integrated heat content (relative to 0C) of water crossing surface',&
1651 handles%id_total_heat_content_massout = register_scalar_field(
'ocean_model', &
1652 'total_heat_content_massout', time, diag, &
1653 long_name=
'Area integrated heat content (relative to 0C) of water leaving ocean', &
1655 cmor_field_name=
'total_hfevapds', &
1656 cmor_standard_name= &
1657 'temperature_flux_due_to_evaporation_expressed_as_heat_flux_out_of_sea_water_area_integrated',&
1658 cmor_long_name=
'Heat Flux Out of Sea Water due to Evaporating Water Area Integrated')
1660 handles%id_total_heat_content_massin = register_scalar_field(
'ocean_model', &
1661 'total_heat_content_massin', time, diag, &
1662 long_name=
'Area integrated heat content (relative to 0C) of water entering ocean',&
1665 handles%id_total_net_heat_coupler = register_scalar_field(
'ocean_model', &
1666 'total_net_heat_coupler', time, diag, &
1667 long_name=
'Area integrated surface heat flux from SW+LW+latent+sensible+seaice_melt_heat (via the coupler)',&
1670 handles%id_total_net_heat_surface = register_scalar_field(
'ocean_model', &
1671 'total_net_heat_surface', time, diag, &
1672 long_name=
'Area integrated surface heat flux from SW+LW+lat+sens+mass+frazil+restore or flux adjustments', &
1674 cmor_field_name=
'total_hfds', &
1675 cmor_standard_name=
'surface_downward_heat_flux_in_sea_water_area_integrated', &
1677 'Surface Ocean Heat Flux from SW+LW+latent+sensible+mass transfer+frazil Area Integrated')
1679 handles%id_total_sw = register_scalar_field(
'ocean_model', &
1680 'total_sw', time, diag, &
1681 long_name=
'Area integrated net downward shortwave at sea water surface', &
1683 cmor_field_name=
'total_rsntds', &
1684 cmor_standard_name=
'net_downward_shortwave_flux_at_sea_water_surface_area_integrated',&
1686 'Net Downward Shortwave Radiation at Sea Water Surface Area Integrated')
1688 handles%id_total_LwLatSens = register_scalar_field(
'ocean_model',&
1689 'total_LwLatSens', time, diag, &
1690 long_name=
'Area integrated longwave+latent+sensible heating',&
1693 handles%id_total_lw = register_scalar_field(
'ocean_model', &
1694 'total_lw', time, diag, &
1695 long_name=
'Area integrated net downward longwave at sea water surface', &
1697 cmor_field_name=
'total_rlntds', &
1698 cmor_standard_name=
'surface_net_downward_longwave_flux_area_integrated',&
1700 'Surface Net Downward Longwave Radiation Area Integrated')
1702 handles%id_total_lat = register_scalar_field(
'ocean_model', &
1703 'total_lat', time, diag, &
1704 long_name=
'Area integrated surface downward latent heat flux', &
1706 cmor_field_name=
'total_hflso', &
1707 cmor_standard_name=
'surface_downward_latent_heat_flux_area_integrated',&
1709 'Surface Downward Latent Heat Flux Area Integrated')
1711 handles%id_total_lat_evap = register_scalar_field(
'ocean_model', &
1712 'total_lat_evap', time, diag, &
1713 long_name=
'Area integrated latent heat flux due to evap/condense',&
1716 handles%id_total_lat_fprec = register_scalar_field(
'ocean_model', &
1717 'total_lat_fprec', time, diag, &
1718 long_name=
'Area integrated latent heat flux due to melting frozen precip', &
1720 cmor_field_name=
'total_hfsnthermds', &
1721 cmor_standard_name=
'heat_flux_into_sea_water_due_to_snow_thermodynamics_area_integrated',&
1723 'Latent Heat to Melt Frozen Precipitation Area Integrated')
1725 handles%id_total_lat_frunoff = register_scalar_field(
'ocean_model', &
1726 'total_lat_frunoff', time, diag, &
1727 long_name=
'Area integrated latent heat flux due to melting icebergs', &
1729 cmor_field_name=
'total_hfibthermds', &
1730 cmor_standard_name=
'heat_flux_into_sea_water_due_to_iceberg_thermodynamics_area_integrated',&
1732 'Heat Flux into Sea Water due to Iceberg Thermodynamics Area Integrated')
1734 handles%id_total_sens = register_scalar_field(
'ocean_model', &
1735 'total_sens', time, diag, &
1736 long_name=
'Area integrated downward sensible heat flux', &
1738 cmor_field_name=
'total_hfsso', &
1739 cmor_standard_name=
'surface_downward_sensible_heat_flux_area_integrated',&
1741 'Surface Downward Sensible Heat Flux Area Integrated')
1743 handles%id_total_heat_added = register_scalar_field(
'ocean_model',&
1744 'total_heat_adjustment', time, diag, &
1745 long_name=
'Area integrated surface heat flux from restoring and/or flux adjustment', &
1748 handles%id_total_seaice_melt_heat = register_scalar_field(
'ocean_model',&
1749 'total_seaice_melt_heat', time, diag, &
1750 long_name=
'Area integrated surface heat flux from snow and sea ice melt', &
1756 handles%id_net_heat_coupler_ga = register_scalar_field(
'ocean_model', &
1757 'net_heat_coupler_ga', time, diag, &
1758 long_name=
'Area averaged surface heat flux from SW+LW+latent+sensible+seaice_melt_heat (via the coupler)',&
1761 handles%id_net_heat_surface_ga = register_scalar_field(
'ocean_model', &
1762 'net_heat_surface_ga', time, diag, long_name= &
1763 'Area averaged surface heat flux from SW+LW+lat+sens+mass+frazil+restore+seaice_melt_heat or flux adjustments', &
1765 cmor_field_name=
'ave_hfds', &
1766 cmor_standard_name=
'surface_downward_heat_flux_in_sea_water_area_averaged', &
1768 'Surface Ocean Heat Flux from SW+LW+latent+sensible+mass transfer+frazil Area Averaged')
1770 handles%id_sw_ga = register_scalar_field(
'ocean_model', &
1771 'sw_ga', time, diag, &
1772 long_name=
'Area averaged net downward shortwave at sea water surface', &
1774 cmor_field_name=
'ave_rsntds', &
1775 cmor_standard_name=
'net_downward_shortwave_flux_at_sea_water_surface_area_averaged',&
1777 'Net Downward Shortwave Radiation at Sea Water Surface Area Averaged')
1779 handles%id_LwLatSens_ga = register_scalar_field(
'ocean_model',&
1780 'LwLatSens_ga', time, diag, &
1781 long_name=
'Area averaged longwave+latent+sensible heating',&
1784 handles%id_lw_ga = register_scalar_field(
'ocean_model', &
1785 'lw_ga', time, diag, &
1786 long_name=
'Area averaged net downward longwave at sea water surface', &
1788 cmor_field_name=
'ave_rlntds', &
1789 cmor_standard_name=
'surface_net_downward_longwave_flux_area_averaged',&
1791 'Surface Net Downward Longwave Radiation Area Averaged')
1793 handles%id_lat_ga = register_scalar_field(
'ocean_model', &
1794 'lat_ga', time, diag, &
1795 long_name=
'Area averaged surface downward latent heat flux', &
1797 cmor_field_name=
'ave_hflso', &
1798 cmor_standard_name=
'surface_downward_latent_heat_flux_area_averaged',&
1800 'Surface Downward Latent Heat Flux Area Averaged')
1802 handles%id_sens_ga = register_scalar_field(
'ocean_model', &
1803 'sens_ga', time, diag, &
1804 long_name=
'Area averaged downward sensible heat flux', &
1806 cmor_field_name=
'ave_hfsso', &
1807 cmor_standard_name=
'surface_downward_sensible_heat_flux_area_averaged',&
1809 'Surface Downward Sensible Heat Flux Area Averaged')
1815 handles%id_saltflux = register_diag_field(
'ocean_model',
'salt_flux', diag%axesT1, time,&
1816 'Net salt flux into ocean at surface (restoring + sea-ice)', &
1817 units=
'kg m-2 s-1', conversion=us%R_to_kg_m3*us%Z_to_m*us%s_to_T, &
1818 cmor_field_name=
'sfdsi', cmor_standard_name=
'downward_sea_ice_basal_salt_flux', &
1819 cmor_long_name=
'Downward Sea Ice Basal Salt Flux')
1821 handles%id_saltFluxIn = register_diag_field(
'ocean_model',
'salt_flux_in', diag%axesT1, time, &
1822 'Salt flux into ocean at surface from coupler', &
1823 units=
'kg m-2 s-1', conversion=us%R_to_kg_m3*us%Z_to_m*us%s_to_T)
1825 handles%id_saltFluxAdded = register_diag_field(
'ocean_model',
'salt_flux_added', &
1826 diag%axesT1,time,
'Salt flux into ocean at surface due to restoring or flux adjustment', &
1827 units=
'kg m-2 s-1', conversion=us%R_to_kg_m3*us%Z_to_m*us%s_to_T)
1829 handles%id_saltFluxGlobalAdj = register_scalar_field(
'ocean_model', &
1830 'salt_flux_global_restoring_adjustment', time, diag, &
1831 'Adjustment needed to balance net global salt flux into ocean at surface', &
1834 handles%id_vPrecGlobalAdj = register_scalar_field(
'ocean_model', &
1835 'vprec_global_adjustment', time, diag, &
1836 'Adjustment needed to adjust net vprec into ocean to zero', &
1839 handles%id_netFWGlobalAdj = register_scalar_field(
'ocean_model', &
1840 'net_fresh_water_global_adjustment', time, diag, &
1841 'Adjustment needed to adjust net fresh water into ocean to zero',&
1844 handles%id_saltFluxGlobalScl = register_scalar_field(
'ocean_model', &
1845 'salt_flux_global_restoring_scaling', time, diag, &
1846 'Scaling applied to balance net global salt flux into ocean at surface', &
1849 handles%id_vPrecGlobalScl = register_scalar_field(
'ocean_model',&
1850 'vprec_global_scaling', time, diag, &
1851 'Scaling applied to adjust net vprec into ocean to zero', &
1854 handles%id_netFWGlobalScl = register_scalar_field(
'ocean_model', &
1855 'net_fresh_water_global_scaling', time, diag, &
1856 'Scaling applied to adjust net fresh water into ocean to zero', &
1862 handles%id_total_saltflux = register_scalar_field(
'ocean_model', &
1863 'total_salt_flux', time, diag, &
1864 long_name=
'Area integrated surface salt flux', units=
'kg', &
1865 cmor_field_name=
'total_sfdsi', &
1866 cmor_units=
'kg s-1', &
1867 cmor_standard_name=
'downward_sea_ice_basal_salt_flux_area_integrated',&
1868 cmor_long_name=
'Downward Sea Ice Basal Salt Flux Area Integrated')
1870 handles%id_total_saltFluxIn = register_scalar_field(
'ocean_model',
'total_salt_Flux_In', &
1871 time, diag, long_name=
'Area integrated surface salt flux at surface from coupler', units=
'kg')
1873 handles%id_total_saltFluxAdded = register_scalar_field(
'ocean_model',
'total_salt_Flux_Added', &
1874 time, diag, long_name=
'Area integrated surface salt flux due to restoring or flux adjustment', units=
'kg')
1885 type(
forcing),
intent(inout) :: fluxes
1888 real,
intent(out) :: wt2
1903 type(
forcing),
intent(inout) :: fluxes
1906 real,
intent(out) :: wt2
1915 integer :: i, j, is, ie, js, je, isq, ieq, jsq, jeq, i0, j0
1916 integer :: isd, ied, jsd, jed, isdb, iedb, jsdb, jedb, isr, ier, jsr, jer
1917 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
1918 isq = g%IscB ; ieq = g%IecB ; jsq = g%JscB ; jeq = g%JecB
1919 isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
1920 isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
1923 if (fluxes%dt_buoy_accum < 0)
call mom_error(fatal,
"fluxes_accumulate: "//&
1924 "fluxes must be initialzed before it can be augmented.")
1927 wt1 = fluxes%dt_buoy_accum / (fluxes%dt_buoy_accum + flux_tmp%dt_buoy_accum)
1929 fluxes%dt_buoy_accum = fluxes%dt_buoy_accum + flux_tmp%dt_buoy_accum
1933 if (
present(forces))
then
1934 do j=js,je ;
do i=is,ie
1935 fluxes%p_surf(i,j) = forces%p_surf(i,j)
1936 fluxes%p_surf_full(i,j) = forces%p_surf_full(i,j)
1938 fluxes%ustar(i,j) = wt1*fluxes%ustar(i,j) + wt2*forces%ustar(i,j)
1941 do j=js,je ;
do i=is,ie
1942 fluxes%p_surf(i,j) = flux_tmp%p_surf(i,j)
1943 fluxes%p_surf_full(i,j) = flux_tmp%p_surf_full(i,j)
1945 fluxes%ustar(i,j) = wt1*fluxes%ustar(i,j) + wt2*flux_tmp%ustar(i,j)
1950 do j=js,je ;
do i=is,ie
1953 fluxes%ustar_gustless(i,j) = flux_tmp%ustar_gustless(i,j)
1955 fluxes%evap(i,j) = wt1*fluxes%evap(i,j) + wt2*flux_tmp%evap(i,j)
1956 fluxes%lprec(i,j) = wt1*fluxes%lprec(i,j) + wt2*flux_tmp%lprec(i,j)
1957 fluxes%fprec(i,j) = wt1*fluxes%fprec(i,j) + wt2*flux_tmp%fprec(i,j)
1958 fluxes%vprec(i,j) = wt1*fluxes%vprec(i,j) + wt2*flux_tmp%vprec(i,j)
1959 fluxes%lrunoff(i,j) = wt1*fluxes%lrunoff(i,j) + wt2*flux_tmp%lrunoff(i,j)
1960 fluxes%frunoff(i,j) = wt1*fluxes%frunoff(i,j) + wt2*flux_tmp%frunoff(i,j)
1961 fluxes%seaice_melt(i,j) = wt1*fluxes%seaice_melt(i,j) + wt2*flux_tmp%seaice_melt(i,j)
1962 fluxes%sw(i,j) = wt1*fluxes%sw(i,j) + wt2*flux_tmp%sw(i,j)
1963 fluxes%sw_vis_dir(i,j) = wt1*fluxes%sw_vis_dir(i,j) + wt2*flux_tmp%sw_vis_dir(i,j)
1964 fluxes%sw_vis_dif(i,j) = wt1*fluxes%sw_vis_dif(i,j) + wt2*flux_tmp%sw_vis_dif(i,j)
1965 fluxes%sw_nir_dir(i,j) = wt1*fluxes%sw_nir_dir(i,j) + wt2*flux_tmp%sw_nir_dir(i,j)
1966 fluxes%sw_nir_dif(i,j) = wt1*fluxes%sw_nir_dif(i,j) + wt2*flux_tmp%sw_nir_dif(i,j)
1967 fluxes%lw(i,j) = wt1*fluxes%lw(i,j) + wt2*flux_tmp%lw(i,j)
1968 fluxes%latent(i,j) = wt1*fluxes%latent(i,j) + wt2*flux_tmp%latent(i,j)
1969 fluxes%sens(i,j) = wt1*fluxes%sens(i,j) + wt2*flux_tmp%sens(i,j)
1971 fluxes%salt_flux(i,j) = wt1*fluxes%salt_flux(i,j) + wt2*flux_tmp%salt_flux(i,j)
1973 if (
associated(fluxes%heat_added) .and.
associated(flux_tmp%heat_added))
then
1974 do j=js,je ;
do i=is,ie
1975 fluxes%heat_added(i,j) = wt1*fluxes%heat_added(i,j) + wt2*flux_tmp%heat_added(i,j)
1979 if (
associated(fluxes%heat_content_cond) .and.
associated(flux_tmp%heat_content_cond))
then
1980 do j=js,je ;
do i=is,ie
1981 fluxes%heat_content_cond(i,j) = wt1*fluxes%heat_content_cond(i,j) + wt2*flux_tmp%heat_content_cond(i,j)
1984 if (
associated(fluxes%heat_content_lprec) .and.
associated(flux_tmp%heat_content_lprec))
then
1985 do j=js,je ;
do i=is,ie
1986 fluxes%heat_content_lprec(i,j) = wt1*fluxes%heat_content_lprec(i,j) + wt2*flux_tmp%heat_content_lprec(i,j)
1989 if (
associated(fluxes%heat_content_fprec) .and.
associated(flux_tmp%heat_content_fprec))
then
1990 do j=js,je ;
do i=is,ie
1991 fluxes%heat_content_fprec(i,j) = wt1*fluxes%heat_content_fprec(i,j) + wt2*flux_tmp%heat_content_fprec(i,j)
1994 if (
associated(fluxes%heat_content_icemelt) .and.
associated(flux_tmp%heat_content_icemelt))
then
1995 do j=js,je ;
do i=is,ie
1996 fluxes%heat_content_icemelt(i,j) = wt1*fluxes%heat_content_icemelt(i,j) + wt2*flux_tmp%heat_content_icemelt(i,j)
1999 if (
associated(fluxes%heat_content_vprec) .and.
associated(flux_tmp%heat_content_vprec))
then
2000 do j=js,je ;
do i=is,ie
2001 fluxes%heat_content_vprec(i,j) = wt1*fluxes%heat_content_vprec(i,j) + wt2*flux_tmp%heat_content_vprec(i,j)
2004 if (
associated(fluxes%heat_content_lrunoff) .and.
associated(flux_tmp%heat_content_lrunoff))
then
2005 do j=js,je ;
do i=is,ie
2006 fluxes%heat_content_lrunoff(i,j) = wt1*fluxes%heat_content_lrunoff(i,j) + wt2*flux_tmp%heat_content_lrunoff(i,j)
2009 if (
associated(fluxes%heat_content_frunoff) .and.
associated(flux_tmp%heat_content_frunoff))
then
2010 do j=js,je ;
do i=is,ie
2011 fluxes%heat_content_frunoff(i,j) = wt1*fluxes%heat_content_frunoff(i,j) + wt2*flux_tmp%heat_content_frunoff(i,j)
2014 if (
associated(fluxes%heat_content_icemelt) .and.
associated(flux_tmp%heat_content_icemelt))
then
2015 do j=js,je ;
do i=is,ie
2016 fluxes%heat_content_icemelt(i,j) = wt1*fluxes%heat_content_icemelt(i,j) + wt2*flux_tmp%heat_content_icemelt(i,j)
2020 if (
associated(fluxes%ustar_shelf) .and.
associated(flux_tmp%ustar_shelf))
then
2021 do i=isd,ied ;
do j=jsd,jed
2022 fluxes%ustar_shelf(i,j) = flux_tmp%ustar_shelf(i,j)
2025 if (
associated(fluxes%iceshelf_melt) .and.
associated(flux_tmp%iceshelf_melt))
then
2026 do i=isd,ied ;
do j=jsd,jed
2027 fluxes%iceshelf_melt(i,j) = flux_tmp%iceshelf_melt(i,j)
2030 if (
associated(fluxes%frac_shelf_h) .and.
associated(flux_tmp%frac_shelf_h))
then
2031 do i=isd,ied ;
do j=jsd,jed
2032 fluxes%frac_shelf_h(i,j) = flux_tmp%frac_shelf_h(i,j)
2036 if (coupler_type_initialized(fluxes%tr_fluxes) .and. &
2037 coupler_type_initialized(flux_tmp%tr_fluxes)) &
2038 call coupler_type_increment_data(flux_tmp%tr_fluxes, fluxes%tr_fluxes, &
2039 scale_factor=wt2, scale_prev=wt1)
2047 type(
forcing),
intent(inout) :: fluxes
2049 logical,
optional,
intent(in) :: skip_pres
2051 real :: taux2, tauy2
2053 integer :: i, j, is, ie, js, je
2054 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
2056 do_pres = .true. ;
if (
present(skip_pres)) do_pres = .not.skip_pres
2058 if (
associated(forces%ustar) .and.
associated(fluxes%ustar))
then
2059 do j=js,je ;
do i=is,ie
2060 fluxes%ustar(i,j) = forces%ustar(i,j)
2065 if (
associated(forces%p_surf) .and.
associated(fluxes%p_surf))
then
2066 do j=js,je ;
do i=is,ie
2067 fluxes%p_surf(i,j) = forces%p_surf(i,j)
2071 if (
associated(forces%p_surf_full) .and.
associated(fluxes%p_surf_full))
then
2072 do j=js,je ;
do i=is,ie
2073 fluxes%p_surf_full(i,j) = forces%p_surf_full(i,j)
2077 if (
associated(forces%p_surf_SSH, forces%p_surf_full))
then
2078 fluxes%p_surf_SSH => fluxes%p_surf_full
2079 elseif (
associated(forces%p_surf_SSH, forces%p_surf))
then
2080 fluxes%p_surf_SSH => fluxes%p_surf
2090 type(
forcing),
intent(inout) :: fluxes
2093 real,
intent(in) :: rho0
2096 real :: taux2, tauy2
2098 integer :: i, j, is, ie, js, je
2099 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
2101 irho0 = us%L_to_Z / rho0
2103 if (
associated(forces%taux) .and.
associated(forces%tauy) .and. &
2104 associated(fluxes%ustar_gustless))
then
2105 do j=js,je ;
do i=is,ie
2107 if ((g%mask2dCu(i-1,j) + g%mask2dCu(i,j)) > 0) &
2108 taux2 = (g%mask2dCu(i-1,j) * forces%taux(i-1,j)**2 + &
2109 g%mask2dCu(i,j) * forces%taux(i,j)**2) / &
2110 (g%mask2dCu(i-1,j) + g%mask2dCu(i,j))
2112 if ((g%mask2dCv(i,j-1) + g%mask2dCv(i,j)) > 0) &
2113 tauy2 = (g%mask2dCv(i,j-1) * forces%tauy(i,j-1)**2 + &
2114 g%mask2dCv(i,j) * forces%tauy(i,j)**2) / &
2115 (g%mask2dCv(i,j-1) + g%mask2dCv(i,j))
2117 fluxes%ustar_gustless(i,j) = sqrt(us%L_to_Z * sqrt(taux2 + tauy2) / rho0)
2134 if (
associated(forces%net_mass_src)) &
2145 real,
dimension(SZI_(G),SZJ_(G)),
intent(out) :: net_mass_src
2148 real :: rz_t_conversion
2149 integer :: i, j, is, ie, js, je
2150 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
2152 rz_t_conversion = us%R_to_kg_m3*us%Z_to_m*us%s_to_T
2154 net_mass_src(:,:) = 0.0
2155 if (
associated(fluxes%lprec))
then ;
do j=js,je ;
do i=is,ie
2156 net_mass_src(i,j) = net_mass_src(i,j) + rz_t_conversion*fluxes%lprec(i,j)
2157 enddo ;
enddo ;
endif
2158 if (
associated(fluxes%fprec))
then ;
do j=js,je ;
do i=is,ie
2159 net_mass_src(i,j) = net_mass_src(i,j) + rz_t_conversion*fluxes%fprec(i,j)
2160 enddo ;
enddo ;
endif
2161 if (
associated(fluxes%vprec))
then ;
do j=js,je ;
do i=is,ie
2162 net_mass_src(i,j) = net_mass_src(i,j) + rz_t_conversion*fluxes%vprec(i,j)
2163 enddo ;
enddo ;
endif
2164 if (
associated(fluxes%lrunoff))
then ;
do j=js,je ;
do i=is,ie
2165 net_mass_src(i,j) = net_mass_src(i,j) + rz_t_conversion*fluxes%lrunoff(i,j)
2166 enddo ;
enddo ;
endif
2167 if (
associated(fluxes%frunoff))
then ;
do j=js,je ;
do i=is,ie
2168 net_mass_src(i,j) = net_mass_src(i,j) + rz_t_conversion*fluxes%frunoff(i,j)
2169 enddo ;
enddo ;
endif
2170 if (
associated(fluxes%evap))
then ;
do j=js,je ;
do i=is,ie
2171 net_mass_src(i,j) = net_mass_src(i,j) + rz_t_conversion*fluxes%evap(i,j)
2172 enddo ;
enddo ;
endif
2173 if (
associated(fluxes%seaice_melt))
then ;
do j=js,je ;
do i=is,ie
2174 net_mass_src(i,j) = net_mass_src(i,j) + rz_t_conversion*fluxes%seaice_melt(i,j)
2175 enddo ;
enddo ;
endif
2186 real :: taux2, tauy2
2187 integer :: i, j, is, ie, js, je
2188 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
2190 if (
associated(forces%ustar) .and.
associated(fluxes%ustar))
then
2191 do j=js,je ;
do i=is,ie
2192 forces%ustar(i,j) = fluxes%ustar(i,j)
2202 real,
intent(in) :: dt
2204 type(time_type),
intent(in) :: time_end
2208 integer :: i,j,is,ie,js,je
2210 call cpu_clock_begin(handles%id_clock_forcing)
2212 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
2216 if ((handles%id_taux > 0) .and.
associated(forces%taux)) &
2217 call post_data(handles%id_taux, forces%taux, diag)
2219 if ((handles%id_tauy > 0) .and.
associated(forces%tauy)) &
2220 call post_data(handles%id_tauy, forces%tauy, diag)
2222 if ((handles%id_mass_berg > 0) .and.
associated(forces%mass_berg)) &
2223 call post_data(handles%id_mass_berg, forces%mass_berg, diag)
2225 if ((handles%id_area_berg > 0) .and.
associated(forces%area_berg)) &
2226 call post_data(handles%id_area_berg, forces%area_berg, diag)
2231 call cpu_clock_end(handles%id_clock_forcing)
2239 type(
surface),
intent(in) :: sfc_state
2243 type(time_type),
intent(in) :: time_end
2248 real,
dimension(SZI_(G),SZJ_(G)) :: res
2249 real :: total_transport
2252 real :: rz_t_conversion
2255 integer :: i,j,is,ie,js,je
2257 call cpu_clock_begin(handles%id_clock_forcing)
2260 rz_t_conversion = us%R_to_kg_m3*us%Z_to_m*us%s_to_T
2261 i_dt = 1.0 / (us%T_to_s*fluxes%dt_buoy_accum)
2263 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
2270 if (handles%id_prcme > 0 .or. handles%id_total_prcme > 0 .or. handles%id_prcme_ga > 0)
then
2271 do j=js,je ;
do i=is,ie
2273 if (
associated(fluxes%lprec)) res(i,j) = res(i,j) + rz_t_conversion*fluxes%lprec(i,j)
2274 if (
associated(fluxes%fprec)) res(i,j) = res(i,j) + rz_t_conversion*fluxes%fprec(i,j)
2276 if (
associated(fluxes%evap)) res(i,j) = res(i,j) + rz_t_conversion*fluxes%evap(i,j)
2277 if (
associated(fluxes%lrunoff)) res(i,j) = res(i,j) + rz_t_conversion*fluxes%lrunoff(i,j)
2278 if (
associated(fluxes%frunoff)) res(i,j) = res(i,j) + rz_t_conversion*fluxes%frunoff(i,j)
2279 if (
associated(fluxes%vprec)) res(i,j) = res(i,j) + rz_t_conversion*fluxes%vprec(i,j)
2280 if (
associated(fluxes%seaice_melt)) res(i,j) = res(i,j) + rz_t_conversion*fluxes%seaice_melt(i,j)
2282 if (handles%id_prcme > 0)
call post_data(handles%id_prcme, res, diag)
2283 if (handles%id_total_prcme > 0)
then
2285 call post_data(handles%id_total_prcme, total_transport, diag)
2287 if (handles%id_prcme_ga > 0)
then
2289 call post_data(handles%id_prcme_ga, ave_flux, diag)
2293 if (handles%id_net_massout > 0 .or. handles%id_total_net_massout > 0)
then
2294 do j=js,je ;
do i=is,ie
2296 if (
associated(fluxes%lprec))
then
2297 if (fluxes%lprec(i,j) < 0.0) res(i,j) = res(i,j) + rz_t_conversion*fluxes%lprec(i,j)
2299 if (
associated(fluxes%vprec))
then
2300 if (fluxes%vprec(i,j) < 0.0) res(i,j) = res(i,j) + rz_t_conversion*fluxes%vprec(i,j)
2302 if (
associated(fluxes%evap))
then
2303 if (fluxes%evap(i,j) < 0.0) res(i,j) = res(i,j) + rz_t_conversion*fluxes%evap(i,j)
2305 if (
associated(fluxes%seaice_melt))
then
2306 if (fluxes%seaice_melt(i,j) < 0.0) &
2307 res(i,j) = res(i,j) + rz_t_conversion*fluxes%seaice_melt(i,j)
2310 if (handles%id_net_massout > 0)
call post_data(handles%id_net_massout, res, diag)
2311 if (handles%id_total_net_massout > 0)
then
2313 call post_data(handles%id_total_net_massout, total_transport, diag)
2317 if (handles%id_massout_flux > 0 .and.
associated(fluxes%netMassOut)) &
2318 call post_data(handles%id_massout_flux,fluxes%netMassOut,diag)
2320 if (handles%id_net_massin > 0 .or. handles%id_total_net_massin > 0)
then
2321 do j=js,je ;
do i=is,ie
2324 if (
associated(fluxes%fprec)) &
2325 res(i,j) = res(i,j) + rz_t_conversion*fluxes%fprec(i,j)
2326 if (
associated(fluxes%lrunoff)) &
2327 res(i,j) = res(i,j) + rz_t_conversion*fluxes%lrunoff(i,j)
2328 if (
associated(fluxes%frunoff)) &
2329 res(i,j) = res(i,j) + rz_t_conversion*fluxes%frunoff(i,j)
2331 if (
associated(fluxes%lprec))
then
2332 if (fluxes%lprec(i,j) > 0.0) res(i,j) = res(i,j) + rz_t_conversion*fluxes%lprec(i,j)
2334 if (
associated(fluxes%vprec))
then
2335 if (fluxes%vprec(i,j) > 0.0) res(i,j) = res(i,j) + rz_t_conversion*fluxes%vprec(i,j)
2338 if (
associated(fluxes%evap))
then
2339 if (fluxes%evap(i,j) > 0.0) res(i,j) = res(i,j) + rz_t_conversion*fluxes%evap(i,j)
2341 if (
associated(fluxes%seaice_melt))
then
2342 if (fluxes%seaice_melt(i,j) > 0.0) &
2343 res(i,j) = res(i,j) + rz_t_conversion*fluxes%seaice_melt(i,j)
2346 if (handles%id_net_massin > 0)
call post_data(handles%id_net_massin, res, diag)
2347 if (handles%id_total_net_massin > 0)
then
2349 call post_data(handles%id_total_net_massin, total_transport, diag)
2353 if (handles%id_massin_flux > 0 .and.
associated(fluxes%netMassIn)) &
2354 call post_data(handles%id_massin_flux,fluxes%netMassIn,diag)
2356 if ((handles%id_evap > 0) .and.
associated(fluxes%evap)) &
2357 call post_data(handles%id_evap, fluxes%evap, diag)
2358 if ((handles%id_total_evap > 0) .and.
associated(fluxes%evap))
then
2360 call post_data(handles%id_total_evap, total_transport, diag)
2362 if ((handles%id_evap_ga > 0) .and.
associated(fluxes%evap))
then
2364 call post_data(handles%id_evap_ga, ave_flux, diag)
2367 if (
associated(fluxes%lprec) .and.
associated(fluxes%fprec))
then
2368 do j=js,je ;
do i=is,ie
2369 res(i,j) = rz_t_conversion* (fluxes%lprec(i,j) + fluxes%fprec(i,j))
2371 if (handles%id_precip > 0)
call post_data(handles%id_precip, res, diag)
2372 if (handles%id_total_precip > 0)
then
2374 call post_data(handles%id_total_precip, total_transport, diag)
2376 if (handles%id_precip_ga > 0)
then
2378 call post_data(handles%id_precip_ga, ave_flux, diag)
2382 if (
associated(fluxes%lprec))
then
2383 if (handles%id_lprec > 0)
call post_data(handles%id_lprec, fluxes%lprec, diag)
2384 if (handles%id_total_lprec > 0)
then
2386 call post_data(handles%id_total_lprec, total_transport, diag)
2388 if (handles%id_lprec_ga > 0)
then
2390 call post_data(handles%id_lprec_ga, ave_flux, diag)
2394 if (
associated(fluxes%fprec))
then
2395 if (handles%id_fprec > 0)
call post_data(handles%id_fprec, fluxes%fprec, diag)
2396 if (handles%id_total_fprec > 0)
then
2398 call post_data(handles%id_total_fprec, total_transport, diag)
2400 if (handles%id_fprec_ga > 0)
then
2402 call post_data(handles%id_fprec_ga, ave_flux, diag)
2406 if (
associated(fluxes%vprec))
then
2407 if (handles%id_vprec > 0)
call post_data(handles%id_vprec, fluxes%vprec, diag)
2408 if (handles%id_total_vprec > 0)
then
2410 call post_data(handles%id_total_vprec, total_transport, diag)
2412 if (handles%id_vprec_ga > 0)
then
2414 call post_data(handles%id_vprec_ga, ave_flux, diag)
2418 if (
associated(fluxes%lrunoff))
then
2419 if (handles%id_lrunoff > 0)
call post_data(handles%id_lrunoff, fluxes%lrunoff, diag)
2420 if (handles%id_total_lrunoff > 0)
then
2422 call post_data(handles%id_total_lrunoff, total_transport, diag)
2426 if (
associated(fluxes%frunoff))
then
2427 if (handles%id_frunoff > 0)
call post_data(handles%id_frunoff, fluxes%frunoff, diag)
2428 if (handles%id_total_frunoff > 0)
then
2430 call post_data(handles%id_total_frunoff, total_transport, diag)
2434 if (
associated(fluxes%seaice_melt))
then
2435 if (handles%id_seaice_melt > 0)
call post_data(handles%id_seaice_melt, fluxes%seaice_melt, diag)
2436 if (handles%id_total_seaice_melt > 0)
then
2438 call post_data(handles%id_total_seaice_melt, total_transport, diag)
2444 if ((handles%id_heat_content_lrunoff > 0) .and.
associated(fluxes%heat_content_lrunoff)) &
2445 call post_data(handles%id_heat_content_lrunoff, fluxes%heat_content_lrunoff, diag)
2446 if ((handles%id_total_heat_content_lrunoff > 0) .and.
associated(fluxes%heat_content_lrunoff))
then
2448 call post_data(handles%id_total_heat_content_lrunoff, total_transport, diag)
2451 if ((handles%id_heat_content_frunoff > 0) .and.
associated(fluxes%heat_content_frunoff)) &
2452 call post_data(handles%id_heat_content_frunoff, fluxes%heat_content_frunoff, diag)
2453 if ((handles%id_total_heat_content_frunoff > 0) .and.
associated(fluxes%heat_content_frunoff))
then
2455 call post_data(handles%id_total_heat_content_frunoff, total_transport, diag)
2458 if ((handles%id_heat_content_lprec > 0) .and.
associated(fluxes%heat_content_lprec)) &
2459 call post_data(handles%id_heat_content_lprec, fluxes%heat_content_lprec, diag)
2460 if ((handles%id_total_heat_content_lprec > 0) .and.
associated(fluxes%heat_content_lprec))
then
2462 call post_data(handles%id_total_heat_content_lprec, total_transport, diag)
2465 if ((handles%id_heat_content_fprec > 0) .and.
associated(fluxes%heat_content_fprec)) &
2466 call post_data(handles%id_heat_content_fprec, fluxes%heat_content_fprec, diag)
2467 if ((handles%id_total_heat_content_fprec > 0) .and.
associated(fluxes%heat_content_fprec))
then
2469 call post_data(handles%id_total_heat_content_fprec, total_transport, diag)
2472 if ((handles%id_heat_content_icemelt > 0) .and.
associated(fluxes%heat_content_icemelt)) &
2473 call post_data(handles%id_heat_content_icemelt, fluxes%heat_content_icemelt, diag)
2474 if ((handles%id_total_heat_content_icemelt > 0) .and.
associated(fluxes%heat_content_icemelt))
then
2476 call post_data(handles%id_total_heat_content_icemelt, total_transport, diag)
2479 if ((handles%id_heat_content_vprec > 0) .and.
associated(fluxes%heat_content_vprec)) &
2480 call post_data(handles%id_heat_content_vprec, fluxes%heat_content_vprec, diag)
2481 if ((handles%id_total_heat_content_vprec > 0) .and.
associated(fluxes%heat_content_vprec))
then
2483 call post_data(handles%id_total_heat_content_vprec, total_transport, diag)
2486 if ((handles%id_heat_content_cond > 0) .and.
associated(fluxes%heat_content_cond)) &
2487 call post_data(handles%id_heat_content_cond, fluxes%heat_content_cond, diag)
2488 if ((handles%id_total_heat_content_cond > 0) .and.
associated(fluxes%heat_content_cond))
then
2490 call post_data(handles%id_total_heat_content_cond, total_transport, diag)
2493 if ((handles%id_heat_content_massout > 0) .and.
associated(fluxes%heat_content_massout)) &
2494 call post_data(handles%id_heat_content_massout, fluxes%heat_content_massout, diag)
2495 if ((handles%id_total_heat_content_massout > 0) .and.
associated(fluxes%heat_content_massout))
then
2497 call post_data(handles%id_total_heat_content_massout, total_transport, diag)
2500 if ((handles%id_heat_content_massin > 0) .and.
associated(fluxes%heat_content_massin)) &
2501 call post_data(handles%id_heat_content_massin, fluxes%heat_content_massin, diag)
2502 if ((handles%id_total_heat_content_massin > 0) .and.
associated(fluxes%heat_content_massin))
then
2504 call post_data(handles%id_total_heat_content_massin, total_transport, diag)
2507 if (handles%id_net_heat_coupler > 0 .or. handles%id_total_net_heat_coupler > 0 .or. &
2508 handles%id_net_heat_coupler_ga > 0. )
then
2509 do j=js,je ;
do i=is,ie
2511 if (
associated(fluxes%LW)) res(i,j) = res(i,j) + fluxes%LW(i,j)
2512 if (
associated(fluxes%latent)) res(i,j) = res(i,j) + fluxes%latent(i,j)
2513 if (
associated(fluxes%sens)) res(i,j) = res(i,j) + fluxes%sens(i,j)
2514 if (
associated(fluxes%SW)) res(i,j) = res(i,j) + fluxes%SW(i,j)
2515 if (
associated(fluxes%seaice_melt_heat)) res(i,j) = res(i,j) + fluxes%seaice_melt_heat(i,j)
2517 if (handles%id_net_heat_coupler > 0)
call post_data(handles%id_net_heat_coupler, res, diag)
2518 if (handles%id_total_net_heat_coupler > 0)
then
2520 call post_data(handles%id_total_net_heat_coupler, total_transport, diag)
2522 if (handles%id_net_heat_coupler_ga > 0)
then
2524 call post_data(handles%id_net_heat_coupler_ga, ave_flux, diag)
2528 if (handles%id_net_heat_surface > 0 .or. handles%id_total_net_heat_surface > 0 .or. &
2529 handles%id_net_heat_surface_ga > 0. )
then
2530 do j=js,je ;
do i=is,ie
2532 if (
associated(fluxes%LW)) res(i,j) = res(i,j) + fluxes%LW(i,j)
2533 if (
associated(fluxes%latent)) res(i,j) = res(i,j) + fluxes%latent(i,j)
2534 if (
associated(fluxes%sens)) res(i,j) = res(i,j) + fluxes%sens(i,j)
2535 if (
associated(fluxes%SW)) res(i,j) = res(i,j) + fluxes%SW(i,j)
2536 if (
associated(fluxes%seaice_melt_heat)) res(i,j) = res(i,j) + fluxes%seaice_melt_heat(i,j)
2537 if (
associated(sfc_state%frazil)) res(i,j) = res(i,j) + sfc_state%frazil(i,j) * i_dt
2541 if (
associated(fluxes%heat_content_lrunoff)) &
2542 res(i,j) = res(i,j) + rz_t_conversion*fluxes%heat_content_lrunoff(i,j)
2543 if (
associated(fluxes%heat_content_frunoff)) &
2544 res(i,j) = res(i,j) + rz_t_conversion*fluxes%heat_content_frunoff(i,j)
2545 if (
associated(fluxes%heat_content_lprec)) &
2546 res(i,j) = res(i,j) + rz_t_conversion*fluxes%heat_content_lprec(i,j)
2547 if (
associated(fluxes%heat_content_fprec)) &
2548 res(i,j) = res(i,j) + rz_t_conversion*fluxes%heat_content_fprec(i,j)
2549 if (
associated(fluxes%heat_content_icemelt)) &
2550 res(i,j) = res(i,j) + rz_t_conversion*fluxes%heat_content_icemelt(i,j)
2551 if (
associated(fluxes%heat_content_vprec)) &
2552 res(i,j) = res(i,j) + rz_t_conversion*fluxes%heat_content_vprec(i,j)
2553 if (
associated(fluxes%heat_content_cond)) &
2554 res(i,j) = res(i,j) + rz_t_conversion*fluxes%heat_content_cond(i,j)
2555 if (
associated(fluxes%heat_content_massout)) &
2556 res(i,j) = res(i,j) + rz_t_conversion*fluxes%heat_content_massout(i,j)
2558 if (
associated(fluxes%heat_added)) res(i,j) = res(i,j) + fluxes%heat_added(i,j)
2560 if (handles%id_net_heat_surface > 0)
call post_data(handles%id_net_heat_surface, res, diag)
2562 if (handles%id_total_net_heat_surface > 0)
then
2564 call post_data(handles%id_total_net_heat_surface, total_transport, diag)
2566 if (handles%id_net_heat_surface_ga > 0)
then
2568 call post_data(handles%id_net_heat_surface_ga, ave_flux, diag)
2572 if (handles%id_heat_content_surfwater > 0 .or. handles%id_total_heat_content_surfwater > 0)
then
2573 do j=js,je ;
do i=is,ie
2578 if (
associated(fluxes%heat_content_lrunoff)) res(i,j) = res(i,j) + fluxes%heat_content_lrunoff(i,j)
2579 if (
associated(fluxes%heat_content_frunoff)) res(i,j) = res(i,j) + fluxes%heat_content_frunoff(i,j)
2580 if (
associated(fluxes%heat_content_lprec)) res(i,j) = res(i,j) + fluxes%heat_content_lprec(i,j)
2581 if (
associated(fluxes%heat_content_icemelt)) res(i,j) = res(i,j) + fluxes%heat_content_icemelt(i,j)
2582 if (
associated(fluxes%heat_content_fprec)) res(i,j) = res(i,j) + fluxes%heat_content_fprec(i,j)
2583 if (
associated(fluxes%heat_content_vprec)) res(i,j) = res(i,j) + fluxes%heat_content_vprec(i,j)
2584 if (
associated(fluxes%heat_content_cond)) res(i,j) = res(i,j) + fluxes%heat_content_cond(i,j)
2585 if (
associated(fluxes%heat_content_massout)) res(i,j) = res(i,j) + fluxes%heat_content_massout(i,j)
2588 if (handles%id_heat_content_surfwater > 0)
call post_data(handles%id_heat_content_surfwater, res, diag)
2589 if (handles%id_total_heat_content_surfwater > 0)
then
2591 call post_data(handles%id_total_heat_content_surfwater, total_transport, diag)
2596 if (handles%id_hfrunoffds > 0)
then
2597 do j=js,je ;
do i=is,ie
2599 if (
associated(fluxes%heat_content_lrunoff)) res(i,j) = res(i,j) + fluxes%heat_content_lrunoff(i,j)
2600 if (
associated(fluxes%heat_content_frunoff)) res(i,j) = res(i,j) + fluxes%heat_content_frunoff(i,j)
2602 call post_data(handles%id_hfrunoffds, res, diag)
2606 if (handles%id_hfrainds > 0)
then
2607 do j=js,je ;
do i=is,ie
2609 if (
associated(fluxes%heat_content_lprec)) res(i,j) = res(i,j) + fluxes%heat_content_lprec(i,j)
2610 if (
associated(fluxes%heat_content_fprec)) res(i,j) = res(i,j) + fluxes%heat_content_fprec(i,j)
2611 if (
associated(fluxes%heat_content_cond)) res(i,j) = res(i,j) + fluxes%heat_content_cond(i,j)
2613 call post_data(handles%id_hfrainds, res, diag)
2616 if ((handles%id_LwLatSens > 0) .and.
associated(fluxes%lw) .and. &
2617 associated(fluxes%latent) .and.
associated(fluxes%sens))
then
2618 do j=js,je ;
do i=is,ie
2619 res(i,j) = (fluxes%lw(i,j) + fluxes%latent(i,j)) + fluxes%sens(i,j)
2621 call post_data(handles%id_LwLatSens, res, diag)
2624 if ((handles%id_total_LwLatSens > 0) .and.
associated(fluxes%lw) .and. &
2625 associated(fluxes%latent) .and.
associated(fluxes%sens))
then
2626 do j=js,je ;
do i=is,ie
2627 res(i,j) = (fluxes%lw(i,j) + fluxes%latent(i,j)) + fluxes%sens(i,j)
2630 call post_data(handles%id_total_LwLatSens, total_transport, diag)
2633 if ((handles%id_LwLatSens_ga > 0) .and.
associated(fluxes%lw) .and. &
2634 associated(fluxes%latent) .and.
associated(fluxes%sens))
then
2635 do j=js,je ;
do i=is,ie
2636 res(i,j) = (fluxes%lw(i,j) + fluxes%latent(i,j)) + fluxes%sens(i,j)
2639 call post_data(handles%id_LwLatSens_ga, ave_flux, diag)
2642 if ((handles%id_sw > 0) .and.
associated(fluxes%sw))
then
2643 call post_data(handles%id_sw, fluxes%sw, diag)
2645 if ((handles%id_sw_vis > 0) .and.
associated(fluxes%sw_vis_dir) .and. &
2646 associated(fluxes%sw_vis_dif))
then
2647 call post_data(handles%id_sw_vis, fluxes%sw_vis_dir+fluxes%sw_vis_dif, diag)
2649 if ((handles%id_sw_nir > 0) .and.
associated(fluxes%sw_nir_dir) .and. &
2650 associated(fluxes%sw_nir_dif))
then
2651 call post_data(handles%id_sw_nir, fluxes%sw_nir_dir+fluxes%sw_nir_dif, diag)
2653 if ((handles%id_total_sw > 0) .and.
associated(fluxes%sw))
then
2655 call post_data(handles%id_total_sw, total_transport, diag)
2657 if ((handles%id_sw_ga > 0) .and.
associated(fluxes%sw))
then
2659 call post_data(handles%id_sw_ga, ave_flux, diag)
2662 if ((handles%id_lw > 0) .and.
associated(fluxes%lw))
then
2663 call post_data(handles%id_lw, fluxes%lw, diag)
2665 if ((handles%id_total_lw > 0) .and.
associated(fluxes%lw))
then
2667 call post_data(handles%id_total_lw, total_transport, diag)
2669 if ((handles%id_lw_ga > 0) .and.
associated(fluxes%lw))
then
2671 call post_data(handles%id_lw_ga, ave_flux, diag)
2674 if ((handles%id_lat > 0) .and.
associated(fluxes%latent))
then
2675 call post_data(handles%id_lat, fluxes%latent, diag)
2677 if ((handles%id_total_lat > 0) .and.
associated(fluxes%latent))
then
2679 call post_data(handles%id_total_lat, total_transport, diag)
2681 if ((handles%id_lat_ga > 0) .and.
associated(fluxes%latent))
then
2683 call post_data(handles%id_lat_ga, ave_flux, diag)
2686 if ((handles%id_lat_evap > 0) .and.
associated(fluxes%latent_evap_diag))
then
2687 call post_data(handles%id_lat_evap, fluxes%latent_evap_diag, diag)
2689 if ((handles%id_total_lat_evap > 0) .and.
associated(fluxes%latent_evap_diag))
then
2691 call post_data(handles%id_total_lat_evap, total_transport, diag)
2694 if ((handles%id_lat_fprec > 0) .and.
associated(fluxes%latent_fprec_diag))
then
2695 call post_data(handles%id_lat_fprec, fluxes%latent_fprec_diag, diag)
2697 if ((handles%id_total_lat_fprec > 0) .and.
associated(fluxes%latent_fprec_diag))
then
2699 call post_data(handles%id_total_lat_fprec, total_transport, diag)
2702 if ((handles%id_lat_frunoff > 0) .and.
associated(fluxes%latent_frunoff_diag))
then
2703 call post_data(handles%id_lat_frunoff, fluxes%latent_frunoff_diag, diag)
2705 if (handles%id_total_lat_frunoff > 0 .and.
associated(fluxes%latent_frunoff_diag))
then
2707 call post_data(handles%id_total_lat_frunoff, total_transport, diag)
2710 if ((handles%id_sens > 0) .and.
associated(fluxes%sens))
then
2711 call post_data(handles%id_sens, fluxes%sens, diag)
2714 if ((handles%id_seaice_melt_heat > 0) .and.
associated(fluxes%seaice_melt_heat))
then
2715 call post_data(handles%id_seaice_melt_heat, fluxes%seaice_melt_heat, diag)
2718 if ((handles%id_total_seaice_melt_heat > 0) .and.
associated(fluxes%seaice_melt_heat))
then
2720 call post_data(handles%id_total_seaice_melt_heat, total_transport, diag)
2723 if ((handles%id_total_sens > 0) .and.
associated(fluxes%sens))
then
2725 call post_data(handles%id_total_sens, total_transport, diag)
2727 if ((handles%id_sens_ga > 0) .and.
associated(fluxes%sens))
then
2729 call post_data(handles%id_sens_ga, ave_flux, diag)
2732 if ((handles%id_heat_added > 0) .and.
associated(fluxes%heat_added))
then
2733 call post_data(handles%id_heat_added, fluxes%heat_added, diag)
2736 if ((handles%id_total_heat_added > 0) .and.
associated(fluxes%heat_added))
then
2738 call post_data(handles%id_total_heat_added, total_transport, diag)
2744 if ((handles%id_saltflux > 0) .and.
associated(fluxes%salt_flux)) &
2745 call post_data(handles%id_saltflux, fluxes%salt_flux, diag)
2746 if ((handles%id_total_saltflux > 0) .and.
associated(fluxes%salt_flux))
then
2748 call post_data(handles%id_total_saltflux, total_transport, diag)
2751 if ((handles%id_saltFluxAdded > 0) .and.
associated(fluxes%salt_flux_added)) &
2752 call post_data(handles%id_saltFluxAdded, fluxes%salt_flux_added, diag)
2753 if ((handles%id_total_saltFluxAdded > 0) .and.
associated(fluxes%salt_flux_added))
then
2754 total_transport = ppt2mks*
global_area_integral(fluxes%salt_flux_added, g, scale=rz_t_conversion)
2755 call post_data(handles%id_total_saltFluxAdded, total_transport, diag)
2758 if (handles%id_saltFluxIn > 0 .and.
associated(fluxes%salt_flux_in)) &
2759 call post_data(handles%id_saltFluxIn, fluxes%salt_flux_in, diag)
2760 if ((handles%id_total_saltFluxIn > 0) .and.
associated(fluxes%salt_flux_in))
then
2762 call post_data(handles%id_total_saltFluxIn, total_transport, diag)
2765 if (handles%id_saltFluxGlobalAdj > 0) &
2766 call post_data(handles%id_saltFluxGlobalAdj, fluxes%saltFluxGlobalAdj, diag)
2767 if (handles%id_vPrecGlobalAdj > 0) &
2768 call post_data(handles%id_vPrecGlobalAdj, fluxes%vPrecGlobalAdj, diag)
2769 if (handles%id_netFWGlobalAdj > 0) &
2770 call post_data(handles%id_netFWGlobalAdj, fluxes%netFWGlobalAdj, diag)
2771 if (handles%id_saltFluxGlobalScl > 0) &
2772 call post_data(handles%id_saltFluxGlobalScl, fluxes%saltFluxGlobalScl, diag)
2773 if (handles%id_vPrecGlobalScl > 0) &
2774 call post_data(handles%id_vPrecGlobalScl, fluxes%vPrecGlobalScl, diag)
2775 if (handles%id_netFWGlobalScl > 0) &
2776 call post_data(handles%id_netFWGlobalScl, fluxes%netFWGlobalScl, diag)
2781 if ((handles%id_psurf > 0) .and.
associated(fluxes%p_surf)) &
2782 call post_data(handles%id_psurf, fluxes%p_surf, diag)
2784 if ((handles%id_TKE_tidal > 0) .and.
associated(fluxes%TKE_tidal)) &
2785 call post_data(handles%id_TKE_tidal, fluxes%TKE_tidal, diag)
2787 if ((handles%id_buoy > 0) .and.
associated(fluxes%buoy)) &
2788 call post_data(handles%id_buoy, fluxes%buoy, diag)
2790 if ((handles%id_ustar > 0) .and.
associated(fluxes%ustar)) &
2791 call post_data(handles%id_ustar, fluxes%ustar, diag)
2793 if ((handles%id_ustar_berg > 0) .and.
associated(fluxes%ustar_berg)) &
2794 call post_data(handles%id_ustar_berg, fluxes%ustar_berg, diag)
2796 if ((handles%id_frac_ice_cover > 0) .and.
associated(fluxes%frac_shelf_h)) &
2797 call post_data(handles%id_frac_ice_cover, fluxes%frac_shelf_h, diag)
2799 if ((handles%id_ustar_ice_cover > 0) .and.
associated(fluxes%ustar_shelf)) &
2800 call post_data(handles%id_ustar_ice_cover, fluxes%ustar_shelf, diag)
2805 call cpu_clock_end(handles%id_clock_forcing)
2812 type(
forcing),
intent(inout) :: fluxes
2813 logical,
optional,
intent(in) :: water
2814 logical,
optional,
intent(in) :: heat
2815 logical,
optional,
intent(in) :: ustar
2816 logical,
optional,
intent(in) :: press
2817 logical,
optional,
intent(in) :: shelf
2818 logical,
optional,
intent(in) :: iceberg
2819 logical,
optional,
intent(in) :: salt
2822 integer :: isd, ied, jsd, jed, isdb, iedb, jsdb, jedb
2823 logical :: heat_water
2825 isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
2826 isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
2828 call myalloc(fluxes%ustar,isd,ied,jsd,jed, ustar)
2829 call myalloc(fluxes%ustar_gustless,isd,ied,jsd,jed, ustar)
2831 call myalloc(fluxes%evap,isd,ied,jsd,jed, water)
2832 call myalloc(fluxes%lprec,isd,ied,jsd,jed, water)
2833 call myalloc(fluxes%fprec,isd,ied,jsd,jed, water)
2834 call myalloc(fluxes%vprec,isd,ied,jsd,jed, water)
2835 call myalloc(fluxes%lrunoff,isd,ied,jsd,jed, water)
2836 call myalloc(fluxes%frunoff,isd,ied,jsd,jed, water)
2837 call myalloc(fluxes%seaice_melt,isd,ied,jsd,jed, water)
2838 call myalloc(fluxes%netMassOut,isd,ied,jsd,jed, water)
2839 call myalloc(fluxes%netMassIn,isd,ied,jsd,jed, water)
2840 call myalloc(fluxes%netSalt,isd,ied,jsd,jed, water)
2841 call myalloc(fluxes%seaice_melt_heat,isd,ied,jsd,jed, heat)
2842 call myalloc(fluxes%sw,isd,ied,jsd,jed, heat)
2843 call myalloc(fluxes%lw,isd,ied,jsd,jed, heat)
2844 call myalloc(fluxes%latent,isd,ied,jsd,jed, heat)
2845 call myalloc(fluxes%sens,isd,ied,jsd,jed, heat)
2846 call myalloc(fluxes%latent_evap_diag,isd,ied,jsd,jed, heat)
2847 call myalloc(fluxes%latent_fprec_diag,isd,ied,jsd,jed, heat)
2848 call myalloc(fluxes%latent_frunoff_diag,isd,ied,jsd,jed, heat)
2850 call myalloc(fluxes%salt_flux,isd,ied,jsd,jed, salt)
2852 if (
present(heat) .and.
present(water))
then ;
if (heat .and. water)
then
2853 call myalloc(fluxes%heat_content_cond,isd,ied,jsd,jed, .true.)
2854 call myalloc(fluxes%heat_content_icemelt,isd,ied,jsd,jed, .true.)
2855 call myalloc(fluxes%heat_content_lprec,isd,ied,jsd,jed, .true.)
2856 call myalloc(fluxes%heat_content_fprec,isd,ied,jsd,jed, .true.)
2857 call myalloc(fluxes%heat_content_vprec,isd,ied,jsd,jed, .true.)
2858 call myalloc(fluxes%heat_content_lrunoff,isd,ied,jsd,jed, .true.)
2859 call myalloc(fluxes%heat_content_frunoff,isd,ied,jsd,jed, .true.)
2860 call myalloc(fluxes%heat_content_massout,isd,ied,jsd,jed, .true.)
2861 call myalloc(fluxes%heat_content_massin,isd,ied,jsd,jed, .true.)
2864 call myalloc(fluxes%p_surf,isd,ied,jsd,jed, press)
2866 call myalloc(fluxes%frac_shelf_h,isd,ied,jsd,jed, shelf)
2867 call myalloc(fluxes%ustar_shelf,isd,ied,jsd,jed, shelf)
2868 call myalloc(fluxes%iceshelf_melt,isd,ied,jsd,jed, shelf)
2871 call myalloc(fluxes%ustar_berg,isd,ied,jsd,jed, iceberg)
2872 call myalloc(fluxes%area_berg,isd,ied,jsd,jed, iceberg)
2873 call myalloc(fluxes%mass_berg,isd,ied,jsd,jed, iceberg)
2882 logical,
optional,
intent(in) :: stress
2883 logical,
optional,
intent(in) :: ustar
2884 logical,
optional,
intent(in) :: shelf
2885 logical,
optional,
intent(in) :: press
2886 logical,
optional,
intent(in) :: iceberg
2889 integer :: isd, ied, jsd, jed, isdb, iedb, jsdb, jedb
2890 logical :: heat_water
2892 isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
2893 isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
2895 call myalloc(forces%taux,isdb,iedb,jsd,jed, stress)
2896 call myalloc(forces%tauy,isd,ied,jsdb,jedb, stress)
2898 call myalloc(forces%ustar,isd,ied,jsd,jed, ustar)
2900 call myalloc(forces%p_surf,isd,ied,jsd,jed, press)
2901 call myalloc(forces%p_surf_full,isd,ied,jsd,jed, press)
2902 call myalloc(forces%net_mass_src,isd,ied,jsd,jed, press)
2904 call myalloc(forces%rigidity_ice_u,isdb,iedb,jsd,jed, shelf)
2905 call myalloc(forces%rigidity_ice_v,isd,ied,jsdb,jedb, shelf)
2906 call myalloc(forces%frac_shelf_u,isdb,iedb,jsd,jed, shelf)
2907 call myalloc(forces%frac_shelf_v,isd,ied,jsdb,jedb, shelf)
2910 call myalloc(forces%area_berg,isd,ied,jsd,jed, iceberg)
2911 call myalloc(forces%mass_berg,isd,ied,jsd,jed, iceberg)
2916 subroutine myalloc(array, is, ie, js, je, flag)
2917 real,
dimension(:,:),
pointer :: array
2918 integer,
intent(in) :: is
2919 integer,
intent(in) :: ie
2920 integer,
intent(in) :: js
2921 integer,
intent(in) :: je
2922 logical,
optional,
intent(in) :: flag
2924 if (
present(flag))
then ;
if (flag)
then ;
if (.not.
associated(array))
then
2925 allocate(array(is:ie,js:je)) ; array(is:ie,js:je) = 0.0
2926 endif ;
endif ;
endif
2933 if (
associated(fluxes%ustar))
deallocate(fluxes%ustar)
2934 if (
associated(fluxes%ustar_gustless))
deallocate(fluxes%ustar_gustless)
2935 if (
associated(fluxes%buoy))
deallocate(fluxes%buoy)
2936 if (
associated(fluxes%sw))
deallocate(fluxes%sw)
2937 if (
associated(fluxes%seaice_melt_heat))
deallocate(fluxes%seaice_melt_heat)
2938 if (
associated(fluxes%sw_vis_dir))
deallocate(fluxes%sw_vis_dir)
2939 if (
associated(fluxes%sw_vis_dif))
deallocate(fluxes%sw_vis_dif)
2940 if (
associated(fluxes%sw_nir_dir))
deallocate(fluxes%sw_nir_dir)
2941 if (
associated(fluxes%sw_nir_dif))
deallocate(fluxes%sw_nir_dif)
2942 if (
associated(fluxes%lw))
deallocate(fluxes%lw)
2943 if (
associated(fluxes%latent))
deallocate(fluxes%latent)
2944 if (
associated(fluxes%latent_evap_diag))
deallocate(fluxes%latent_evap_diag)
2945 if (
associated(fluxes%latent_fprec_diag))
deallocate(fluxes%latent_fprec_diag)
2946 if (
associated(fluxes%latent_frunoff_diag))
deallocate(fluxes%latent_frunoff_diag)
2947 if (
associated(fluxes%sens))
deallocate(fluxes%sens)
2948 if (
associated(fluxes%heat_added))
deallocate(fluxes%heat_added)
2949 if (
associated(fluxes%heat_content_lrunoff))
deallocate(fluxes%heat_content_lrunoff)
2950 if (
associated(fluxes%heat_content_frunoff))
deallocate(fluxes%heat_content_frunoff)
2951 if (
associated(fluxes%heat_content_icemelt))
deallocate(fluxes%heat_content_icemelt)
2952 if (
associated(fluxes%heat_content_lprec))
deallocate(fluxes%heat_content_lprec)
2953 if (
associated(fluxes%heat_content_fprec))
deallocate(fluxes%heat_content_fprec)
2954 if (
associated(fluxes%heat_content_cond))
deallocate(fluxes%heat_content_cond)
2955 if (
associated(fluxes%heat_content_massout))
deallocate(fluxes%heat_content_massout)
2956 if (
associated(fluxes%heat_content_massin))
deallocate(fluxes%heat_content_massin)
2957 if (
associated(fluxes%evap))
deallocate(fluxes%evap)
2958 if (
associated(fluxes%lprec))
deallocate(fluxes%lprec)
2959 if (
associated(fluxes%fprec))
deallocate(fluxes%fprec)
2960 if (
associated(fluxes%vprec))
deallocate(fluxes%vprec)
2961 if (
associated(fluxes%lrunoff))
deallocate(fluxes%lrunoff)
2962 if (
associated(fluxes%frunoff))
deallocate(fluxes%frunoff)
2963 if (
associated(fluxes%seaice_melt))
deallocate(fluxes%seaice_melt)
2964 if (
associated(fluxes%salt_flux))
deallocate(fluxes%salt_flux)
2965 if (
associated(fluxes%p_surf_full))
deallocate(fluxes%p_surf_full)
2966 if (
associated(fluxes%p_surf))
deallocate(fluxes%p_surf)
2967 if (
associated(fluxes%TKE_tidal))
deallocate(fluxes%TKE_tidal)
2968 if (
associated(fluxes%ustar_tidal))
deallocate(fluxes%ustar_tidal)
2969 if (
associated(fluxes%ustar_shelf))
deallocate(fluxes%ustar_shelf)
2970 if (
associated(fluxes%iceshelf_melt))
deallocate(fluxes%iceshelf_melt)
2971 if (
associated(fluxes%frac_shelf_h))
deallocate(fluxes%frac_shelf_h)
2972 if (
associated(fluxes%ustar_berg))
deallocate(fluxes%ustar_berg)
2973 if (
associated(fluxes%area_berg))
deallocate(fluxes%area_berg)
2974 if (
associated(fluxes%mass_berg))
deallocate(fluxes%mass_berg)
2976 call coupler_type_destructor(fluxes%tr_fluxes)
2985 if (
associated(forces%taux))
deallocate(forces%taux)
2986 if (
associated(forces%tauy))
deallocate(forces%tauy)
2987 if (
associated(forces%ustar))
deallocate(forces%ustar)
2988 if (
associated(forces%p_surf))
deallocate(forces%p_surf)
2989 if (
associated(forces%p_surf_full))
deallocate(forces%p_surf_full)
2990 if (
associated(forces%net_mass_src))
deallocate(forces%net_mass_src)
2991 if (
associated(forces%rigidity_ice_u))
deallocate(forces%rigidity_ice_u)
2992 if (
associated(forces%rigidity_ice_v))
deallocate(forces%rigidity_ice_v)
2993 if (
associated(forces%frac_shelf_u))
deallocate(forces%frac_shelf_u)
2994 if (
associated(forces%frac_shelf_v))
deallocate(forces%frac_shelf_v)
2995 if (
associated(forces%area_berg))
deallocate(forces%area_berg)
2996 if (
associated(forces%mass_berg))
deallocate(forces%mass_berg)