MOM6
MOM_marine_ice.F90
Go to the documentation of this file.
1 !> Routines incorporating the effects of marine ice (sea-ice and icebergs) into
2 !! the ocean model dynamics and thermodynamics.
4 
5 ! This file is part of MOM6. See LICENSE.md for the license.
6 
7 use mom_constants, only : hlf
9 use mom_domains, only : pass_var, pass_vector, agrid, bgrid_ne, cgrid_ne
10 use mom_domains, only : to_all, omit_corners
11 use mom_error_handler, only : mom_error, fatal, warning
12 use mom_error_handler, only : mom_error, mom_mesg, fatal, warning
16 use mom_grid, only : ocean_grid_type
17 use mom_time_manager, only : time_type
19 use mom_variables, only : surface
20 
21 implicit none ; private
22 
23 #include <MOM_memory.h>
24 
26 
27 !> Control structure for MOM_marine_ice
28 type, public :: marine_ice_cs ; private
29  real :: kv_iceberg !< The viscosity of the icebergs [m2 s-1] (for ice rigidity)
30  real :: berg_area_threshold !< Fraction of grid cell which iceberg must occupy
31  !! so that fluxes below are set to zero. (0.5 is a
32  !! good value to use.) Not applied for negative values.
33  real :: latent_heat_fusion !< Latent heat of fusion [J kg-1]
34  real :: density_iceberg !< A typical density of icebergs [kg m-3] (for ice rigidity)
35 
36  type(time_type), pointer :: time !< A pointer to the ocean model's clock.
37  type(diag_ctrl), pointer :: diag !< A structure that is used to regulate the timing of diagnostic output.
38 end type marine_ice_cs
39 
40 contains
41 
42 !> add_berg_flux_to_shelf adds rigidity and ice-area coverage due to icebergs
43 !! to the forces type fields, and adds ice-areal coverage and modifies various
44 !! thermodynamic fluxes due to the presence of icebergs.
45 subroutine iceberg_forces(G, forces, use_ice_shelf, sfc_state, &
46  time_step, CS)
47  type(ocean_grid_type), intent(inout) :: g !< The ocean's grid structure
48  type(mech_forcing), intent(inout) :: forces !< A structure with the driving mechanical forces
49  type(surface), intent(inout) :: sfc_state !< A structure containing fields that
50  !! describe the surface state of the ocean.
51  logical, intent(in) :: use_ice_shelf !< If true, this configuration uses ice shelves.
52  real, intent(in) :: time_step !< The coupling time step [s].
53  type(marine_ice_cs), pointer :: cs !< Pointer to the control structure for MOM_marine_ice
54 
55  real :: kv_rho_ice ! The viscosity of ice divided by its density [m5 kg-1 s-1].
56  integer :: i, j, is, ie, js, je
57  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
58  !This routine adds iceberg data to the ice shelf data (if ice shelf is used)
59  !which can then be used to change the top of ocean boundary condition used in
60  !the ocean model. This routine is taken from the add_shelf_flux subroutine
61  !within the ice shelf model.
62 
63  if (.not.associated(cs)) return
64 
65  if (.not.(associated(forces%area_berg) .and. associated(forces%mass_berg) ) ) return
66 
67  if (.not.(associated(forces%frac_shelf_u) .and. associated(forces%frac_shelf_v) .and. &
68  associated(forces%rigidity_ice_u) .and. associated(forces%rigidity_ice_v)) ) return
69 
70  ! This section sets or augments the values of fields in forces.
71  if (.not. use_ice_shelf) then
72  forces%frac_shelf_u(:,:) = 0.0 ; forces%frac_shelf_v(:,:) = 0.0
73  endif
74  if (.not. forces%accumulate_rigidity) then
75  forces%rigidity_ice_u(:,:) = 0.0 ; forces%rigidity_ice_v(:,:) = 0.0
76  endif
77 
78  call pass_var(forces%area_berg, g%domain, to_all+omit_corners, halo=1, complete=.false.)
79  call pass_var(forces%mass_berg, g%domain, to_all+omit_corners, halo=1, complete=.true.)
80  kv_rho_ice = cs%kv_iceberg / cs%density_iceberg
81  do j=js,je ; do i=is-1,ie
82  if ((g%areaT(i,j) + g%areaT(i+1,j) > 0.0)) & ! .and. (G%dxdy_u(I,j) > 0.0)) &
83  forces%frac_shelf_u(i,j) = forces%frac_shelf_u(i,j) + &
84  (((forces%area_berg(i,j)*g%US%L_to_m**2*g%areaT(i,j)) + &
85  (forces%area_berg(i+1,j)*g%US%L_to_m**2*g%areaT(i+1,j))) / &
86  (g%US%L_to_m**2*g%areaT(i,j) + g%US%L_to_m**2*g%areaT(i+1,j)) )
87  forces%rigidity_ice_u(i,j) = forces%rigidity_ice_u(i,j) + kv_rho_ice * &
88  min(forces%mass_berg(i,j), forces%mass_berg(i+1,j))
89  enddo ; enddo
90  do j=js-1,je ; do i=is,ie
91  if ((g%areaT(i,j) + g%areaT(i,j+1) > 0.0)) & ! .and. (G%dxdy_v(i,J) > 0.0)) &
92  forces%frac_shelf_v(i,j) = forces%frac_shelf_v(i,j) + &
93  (((forces%area_berg(i,j)*g%US%L_to_m**2*g%areaT(i,j)) + &
94  (forces%area_berg(i,j+1)*g%US%L_to_m**2*g%areaT(i,j+1))) / &
95  (g%US%L_to_m**2*g%areaT(i,j) + g%US%L_to_m**2*g%areaT(i,j+1)) )
96  forces%rigidity_ice_v(i,j) = forces%rigidity_ice_v(i,j) + kv_rho_ice * &
97  min(forces%mass_berg(i,j), forces%mass_berg(i,j+1))
98  enddo ; enddo
99  !### This halo update may be unnecessary. Test it. -RWH
100  call pass_vector(forces%frac_shelf_u, forces%frac_shelf_v, g%domain, to_all, cgrid_ne)
101 
102 end subroutine iceberg_forces
103 
104 !> iceberg_fluxes adds ice-area-coverage and modifies various
105 !! thermodynamic fluxes due to the presence of icebergs.
106 subroutine iceberg_fluxes(G, US, fluxes, use_ice_shelf, sfc_state, &
107  time_step, CS)
108  type(ocean_grid_type), intent(inout) :: g !< The ocean's grid structure
109  type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
110  type(forcing), intent(inout) :: fluxes !< A structure with pointers to themodynamic,
111  !! tracer and mass exchange forcing fields
112  type(surface), intent(inout) :: sfc_state !< A structure containing fields that
113  !! describe the surface state of the ocean.
114  logical, intent(in) :: use_ice_shelf !< If true, this configuration uses ice shelves.
115  real, intent(in) :: time_step !< The coupling time step [s].
116  type(marine_ice_cs), pointer :: cs !< Pointer to the control structure for MOM_marine_ice
117 
118  real :: fraz ! refreezing rate [R Z T-1 ~> kg m-2 s-1]
119  real :: i_dt_lhf ! The inverse of the timestep times the latent heat of fusion [kg J-1 T-1 ~> kg J-1 s-1].
120  integer :: i, j, is, ie, js, je, isd, ied, jsd, jed
121  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
122  isd = g%isd ; jsd = g%jsd ; ied = g%ied ; jed = g%jed
123  !This routine adds iceberg data to the ice shelf data (if ice shelf is used)
124  !which can then be used to change the top of ocean boundary condition used in
125  !the ocean model. This routine is taken from the add_shelf_flux subroutine
126  !within the ice shelf model.
127 
128  if (.not.associated(cs)) return
129  if (.not.(associated(fluxes%area_berg) .and. associated(fluxes%ustar_berg) .and. &
130  associated(fluxes%mass_berg) ) ) return
131  if (.not.(associated(fluxes%frac_shelf_h) .and. associated(fluxes%ustar_shelf)) ) return
132 
133 
134  if (.not.(associated(fluxes%area_berg) .and. associated(fluxes%ustar_berg) .and. &
135  associated(fluxes%mass_berg) ) ) return
136  if (.not. use_ice_shelf) then
137  fluxes%frac_shelf_h(:,:) = 0.
138  fluxes%ustar_shelf(:,:) = 0.
139  endif
140  do j=jsd,jed ; do i=isd,ied ; if (g%areaT(i,j) > 0.0) then
141  fluxes%frac_shelf_h(i,j) = fluxes%frac_shelf_h(i,j) + fluxes%area_berg(i,j)
142  fluxes%ustar_shelf(i,j) = fluxes%ustar_shelf(i,j) + fluxes%ustar_berg(i,j)
143  endif ; enddo ; enddo
144 
145  !Zero'ing out other fluxes under the tabular icebergs
146  if (cs%berg_area_threshold >= 0.) then
147  i_dt_lhf = 1.0 / (us%s_to_T*time_step * cs%latent_heat_fusion)
148  do j=jsd,jed ; do i=isd,ied
149  if (fluxes%frac_shelf_h(i,j) > cs%berg_area_threshold) then
150  ! Only applying for ice shelf covering most of cell.
151 
152  if (associated(fluxes%sw)) fluxes%sw(i,j) = 0.0
153  if (associated(fluxes%lw)) fluxes%lw(i,j) = 0.0
154  if (associated(fluxes%latent)) fluxes%latent(i,j) = 0.0
155  if (associated(fluxes%evap)) fluxes%evap(i,j) = 0.0
156 
157  ! Add frazil formation diagnosed by the ocean model [J m-2] in the
158  ! form of surface layer evaporation [R Z T-1 ~> kg m-2 s-1]. Update lprec in the
159  ! control structure for diagnostic purposes.
160 
161  if (associated(sfc_state%frazil)) then
162  fraz = us%kg_m3_to_R*us%m_to_Z*sfc_state%frazil(i,j) * i_dt_lhf
163  if (associated(fluxes%evap)) &
164  fluxes%evap(i,j) = fluxes%evap(i,j) - fraz
165  ! fluxes%lprec(i,j) = fluxes%lprec(i,j) - fraz
166  sfc_state%frazil(i,j) = 0.0
167  endif
168 
169  !Alon: Should these be set to zero too?
170  if (associated(fluxes%sens)) fluxes%sens(i,j) = 0.0
171  if (associated(fluxes%salt_flux)) fluxes%salt_flux(i,j) = 0.0
172  if (associated(fluxes%lprec)) fluxes%lprec(i,j) = 0.0
173  endif
174  enddo ; enddo
175  endif
176 
177 end subroutine iceberg_fluxes
178 
179 !> Initialize control structure for MOM_marine_ice
180 subroutine marine_ice_init(Time, G, param_file, diag, CS)
181  type(time_type), target, intent(in) :: time !< Current model time
182  type(ocean_grid_type), intent(in) :: g !< Ocean grid structure
183  type(param_file_type), intent(in) :: param_file !< Runtime parameter handles
184  type(diag_ctrl), target, intent(inout) :: diag !< Diagnostics control structure
185  type(marine_ice_cs), pointer :: cs !< Pointer to the control structure for MOM_marine_ice
186 ! This include declares and sets the variable "version".
187 #include "version_variable.h"
188  character(len=40) :: mdl = "MOM_marine_ice" ! This module's name.
189 
190  if (associated(cs)) then
191  call mom_error(warning, "marine_ice_init called with an associated control structure.")
192  return
193  else ; allocate(cs) ; endif
194 
195  ! Write all relevant parameters to the model log.
196  call log_version(mdl, version)
197 
198  call get_param(param_file, mdl, "KV_ICEBERG", cs%kv_iceberg, &
199  "The viscosity of the icebergs", units="m2 s-1",default=1.0e10)
200  call get_param(param_file, mdl, "DENSITY_ICEBERGS", cs%density_iceberg, &
201  "A typical density of icebergs.", units="kg m-3", default=917.0)
202  call get_param(param_file, mdl, "LATENT_HEAT_FUSION", cs%latent_heat_fusion, &
203  "The latent heat of fusion.", units="J/kg", default=hlf)
204  call get_param(param_file, mdl, "BERG_AREA_THRESHOLD", cs%berg_area_threshold, &
205  "Fraction of grid cell which iceberg must occupy, so that fluxes "//&
206  "below berg are set to zero. Not applied for negative "//&
207  "values.", units="non-dim", default=-1.0)
208 
209 end subroutine marine_ice_init
210 
211 end module mom_marine_ice
mom_time_manager
Wraps the FMS time manager functions.
Definition: MOM_time_manager.F90:2
mom_forcing_type::allocate_forcing_type
subroutine, public allocate_forcing_type(G, fluxes, water, heat, ustar, press, shelf, iceberg, salt)
Conditionally allocate fields within the forcing type.
Definition: MOM_forcing_type.F90:2811
mom_forcing_type::mech_forcing
Structure that contains pointers to the mechanical forcing at the surface used to drive the liquid oc...
Definition: MOM_forcing_type.F90:187
mom_marine_ice::iceberg_forces
subroutine, public iceberg_forces(G, forces, use_ice_shelf, sfc_state, time_step, CS)
add_berg_flux_to_shelf adds rigidity and ice-area coverage due to icebergs to the forces type fields,...
Definition: MOM_marine_ice.F90:47
mom_variables::surface
Pointers to various fields which may be used describe the surface state of MOM, and which will be ret...
Definition: MOM_variables.F90:38
mom_diag_mediator::query_averaging_enabled
logical function, public query_averaging_enabled(diag_cs, time_int, time_end)
Call this subroutine to determine whether the averaging is currently enabled. .true....
Definition: MOM_diag_mediator.F90:1850
mom_marine_ice
Routines incorporating the effects of marine ice (sea-ice and icebergs) into the ocean model dynamics...
Definition: MOM_marine_ice.F90:3
mom_file_parser::log_version
An overloaded interface to log version information about modules.
Definition: MOM_file_parser.F90:109
mom_constants
Provides a few physical constants.
Definition: MOM_constants.F90:2
mom_diag_mediator
The subroutines here provide convenient wrappers to the fms diag_manager interfaces with additional d...
Definition: MOM_diag_mediator.F90:3
mom_error_handler::mom_mesg
subroutine, public mom_mesg(message, verb, all_print)
This provides a convenient interface for writing an informative comment.
Definition: MOM_error_handler.F90:53
mom_domains::pass_var
Do a halo update on an array.
Definition: MOM_domains.F90:49
mom_file_parser::param_file_type
A structure that can be parsed to read and document run-time parameters.
Definition: MOM_file_parser.F90:54
mom_file_parser::get_param
An overloaded interface to read and log the values of various types of parameters.
Definition: MOM_file_parser.F90:102
mom_unit_scaling::unit_scale_type
Describes various unit conversion factors.
Definition: MOM_unit_scaling.F90:14
mom_diag_mediator::post_data
Make a diagnostic available for averaging or output.
Definition: MOM_diag_mediator.F90:70
mom_domains::pass_vector
Do a halo update on a pair of arrays representing the two components of a vector.
Definition: MOM_domains.F90:54
mom_forcing_type
This module implements boundary forcing for MOM6.
Definition: MOM_forcing_type.F90:2
mom_domains::to_all
integer, parameter, public to_all
A flag for passing in all directions.
Definition: MOM_domains.F90:132
mom_domains
Describes the decomposed MOM domain and has routines for communications across PEs.
Definition: MOM_domains.F90:2
mom_variables
Provides transparent structures with groups of MOM6 variables and supporting routines.
Definition: MOM_variables.F90:2
mom_file_parser
The MOM6 facility to parse input files for runtime parameters.
Definition: MOM_file_parser.F90:2
mom_grid
Provides the ocean grid type.
Definition: MOM_grid.F90:2
mom_unit_scaling
Provides a transparent unit rescaling type to facilitate dimensional consistency testing.
Definition: MOM_unit_scaling.F90:2
mom_marine_ice::iceberg_fluxes
subroutine, public iceberg_fluxes(G, US, fluxes, use_ice_shelf, sfc_state, time_step, CS)
iceberg_fluxes adds ice-area-coverage and modifies various thermodynamic fluxes due to the presence o...
Definition: MOM_marine_ice.F90:108
mom_forcing_type::forcing
Structure that contains pointers to the boundary forcing used to drive the liquid ocean simulated by ...
Definition: MOM_forcing_type.F90:50
mom_error_handler::mom_error
subroutine, public mom_error(level, message, all_print)
This provides a convenient interface for writing an mpp_error message with run-time filter based on a...
Definition: MOM_error_handler.F90:72
mom_marine_ice::marine_ice_init
subroutine, public marine_ice_init(Time, G, param_file, diag, CS)
Initialize control structure for MOM_marine_ice.
Definition: MOM_marine_ice.F90:181
mom_marine_ice::marine_ice_cs
Control structure for MOM_marine_ice.
Definition: MOM_marine_ice.F90:28
mom_error_handler
Routines for error handling and I/O management.
Definition: MOM_error_handler.F90:2
mom_grid::ocean_grid_type
Ocean grid type. See mom_grid for details.
Definition: MOM_grid.F90:26
mom_diag_mediator::diag_ctrl
The following data type a list of diagnostic fields an their variants, as well as variables that cont...
Definition: MOM_diag_mediator.F90:239