MARBL_forcing_mod.F90

1! This file is part of MOM6, the Modular Ocean Model version 6.
2! See the LICENSE file for licensing information.
3! SPDX-License-Identifier: Apache-2.0
4
5!> This module provides a common datatype to provide forcing for MARBL tracers
6!! regardless of driver
8
9!! This module exists to house code used by multiple drivers in config_src/
10!! for passing forcing fields to MARBL
11!! (This comment can go in the wiki on the NCAR fork?)
12
13use mom_diag_mediator, only : safe_alloc_ptr, diag_ctrl, register_diag_field, post_data
14use mom_time_manager, only : time_type
15use mom_error_handler, only : mom_error, warning, fatal
16use mom_file_parser, only : get_param, log_param, param_file_type
17use mom_grid, only : ocean_grid_type
19use mom_interpolate, only : external_field, init_external_field, time_interp_external
20use mom_io, only : slasher
21use marbl_constants_mod, only : molw_fe
22use mom_forcing_type, only : forcing
23
24implicit none ; private
25
26#include <MOM_memory.h>
27
28public :: marbl_forcing_init
30
31!> Data type used to store diagnostic index returned from register_diag_field()
32!! For the forcing fields that can be written via post_data()
33type, private :: marbl_forcing_diag_ids
34 integer :: atm_fine_dust !< Atmospheric fine dust component of dust_flux
35 integer :: atm_coarse_dust !< Atmospheric coarse dust component of dust_flux
36 integer :: atm_bc !< Atmospheric black carbon component of iron_flux
37 integer :: ice_dust !< Sea-ice dust component of dust_flux
38 integer :: ice_bc !< Sea-ice black carbon component of iron_flux
40
41!> Control structure for this module
42type, public :: marbl_forcing_cs ; private
43 type(diag_ctrl), pointer :: diag => null() !< A structure that is used to
44 !! regulate the timing of diagnostic output.
45
46 real :: dust_ratio_thres !< coarse/fine dust ratio threshold [1]
47 real :: dust_ratio_to_fe_bioavail_frac !< ratio of dust to iron bioavailability fraction [1]
48 real :: fe_bioavail_frac_offset !< offset for iron bioavailability fraction [1]
49 real :: atm_fe_to_bc_ratio !< atmospheric iron to black carbon ratio [1]
50 real :: atm_bc_fe_bioavail_frac !< atmospheric black carbon to iron bioavailablity fraction ratio [1]
51 real :: seaice_fe_to_bc_ratio !< sea-ice iron to black carbon ratio [1]
52 real :: seaice_bc_fe_bioavail_frac !< sea-ice black carbon to iron bioavailablity fraction ratio [1]
53 real :: iron_frac_in_atm_fine_dust !< Fraction of fine dust from the atmosphere that is iron [1]
54 real :: iron_frac_in_atm_coarse_dust !< Fraction of coarse dust from the atmosphere that is iron [1]
55 real :: iron_frac_in_seaice_dust !< Fraction of dust from the sea ice that is iron [1]
56 real :: atm_co2_const !< atmospheric CO2 (if specifying a constant value) [ppm]
57 real :: atm_alt_co2_const !< alternate atmospheric CO2 for _ALT_CO2 tracers
58 !! (if specifying a constant value) [ppm]
59
60 type(marbl_forcing_diag_ids) :: diag_ids !< used for registering and posting some MARBL forcing fields as diagnostics
61
62 logical :: use_marbl_tracers !< most functions can return immediately
63 !! MARBL tracers are turned off
64 integer :: atm_co2_iopt !< Integer version of atm_co2_opt, which determines source of atm_co2
65 integer :: atm_alt_co2_iopt !< Integer version of atm_alt_co2_opt, which determines source of atm_alt_co2
66
67end type marbl_forcing_cs
68
69! Module parameters
70integer, parameter :: atm_co2_constant_iopt = 0 !< module parameter denoting atm_co2_opt = 'constant'
71integer, parameter :: atm_co2_prognostic_iopt = 1 !< module parameter denoting atm_co2_opt = 'diagnostic'
72integer, parameter :: atm_co2_diagnostic_iopt = 2 !< module parameter denoting atm_co2_opt = 'prognostic'
73
74contains
75
76 subroutine marbl_forcing_init(G, US, param_file, diag, day, inputdir, use_marbl, CS)
77 type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
78 type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
79 type(param_file_type), intent(in) :: param_file !< A structure to parse for run-time parameters
80 type(diag_ctrl), target, intent(in) :: diag !< Structure used to regulate diagnostic output.
81 type(time_type), target, intent(in) :: day !< Time of the start of the run.
82 character(len=*), intent(in) :: inputdir !< Directory containing input files
83 logical, intent(in) :: use_marbl !< Is MARBL tracer package active?
84 type(marbl_forcing_cs), pointer, intent(inout) :: cs !< A pointer that is set to point to control
85 !! structure for MARBL forcing
86
87 character(len=40) :: mdl = "MARBL_forcing_mod" ! This module's name.
88 character(len=15) :: atm_co2_opt
89 character(len=200) :: err_message
90
91 if (associated(cs)) then
92 call mom_error(warning, "marbl_forcing_init called with an associated control structure.")
93 return
94 endif
95
96 allocate(cs)
97 cs%diag => diag
98
99 cs%use_marbl_tracers = .true.
100 if (.not. use_marbl) then
101 cs%use_marbl_tracers = .false.
102 return
103 endif
104
105 call get_param(param_file, mdl, "DUST_RATIO_THRES", cs%dust_ratio_thres, &
106 "coarse/fine dust ratio threshold", units="1", default=69.00594)
107 call get_param(param_file, mdl, "DUST_RATIO_TO_FE_BIOAVAIL_FRAC", cs%dust_ratio_to_fe_bioavail_frac, &
108 "ratio of dust to iron bioavailability fraction", units="1", default=1./366.314)
109 call get_param(param_file, mdl, "FE_BIOAVAIL_FRAC_OFFSET", cs%fe_bioavail_frac_offset, &
110 "offset for iron bioavailability fraction", units="1", default=0.0146756)
111 call get_param(param_file, mdl, "ATM_FE_TO_BC_RATIO", cs%atm_fe_to_bc_ratio, &
112 "atmospheric iron to black carbon ratio", units="1", default=1.)
113 call get_param(param_file, mdl, "ATM_BC_FE_BIOAVAIL_FRAC", cs%atm_bc_fe_bioavail_frac, &
114 "atmospheric black carbon to iron bioavailablity fraction ratio", units="1", default=0.06)
115 call get_param(param_file, mdl, "SEAICE_FE_TO_BC_RATIO", cs%seaice_fe_to_bc_ratio, &
116 "sea-ice iron to black carbon ratio", units="1", default=1.)
117 call get_param(param_file, mdl, "SEAICE_BC_FE_BIOAVAIL_FRAC", cs%seaice_bc_fe_bioavail_frac, &
118 "sea-ice black carbon to iron bioavailablity fraction ratio", units="1", default=0.06)
119 call get_param(param_file, mdl, "IRON_FRAC_IN_ATM_FINE_DUST", cs%iron_frac_in_atm_fine_dust, &
120 "Fraction of fine dust from the atmosphere that is iron", units="1", default=0.035)
121 call get_param(param_file, mdl, "IRON_FRAC_IN_ATM_COARSE_DUST", cs%iron_frac_in_atm_coarse_dust, &
122 "Fraction of coarse dust from the atmosphere that is iron", units="1", default=0.035)
123 call get_param(param_file, mdl, "IRON_FRAC_IN_SEAICE_DUST", cs%iron_frac_in_seaice_dust, &
124 "Fraction of dust from sea ice that is iron", units="1", default=0.035)
125 call get_param(param_file, mdl, "ATM_CO2_OPT", atm_co2_opt, &
126 "Source of atmospheric CO2 [constant, diagnostic, or prognostic]", &
127 default="constant")
128 select case (trim(atm_co2_opt))
129 case("prognostic")
130 cs%atm_co2_iopt = atm_co2_prognostic_iopt
131 case("diagnostic")
132 cs%atm_co2_iopt = atm_co2_diagnostic_iopt
133 case("constant")
134 cs%atm_co2_iopt = atm_co2_constant_iopt
135 case DEFAULT
136 write(err_message, "(3A)") "'", trim(atm_co2_opt), "' is not a valid ATM_CO2_OPT value"
137 call mom_error(fatal, err_message)
138 end select
139 if (cs%atm_co2_iopt == atm_co2_constant_iopt) then
140 call get_param(param_file, mdl, "ATM_CO2_CONST", cs%atm_co2_const, &
141 "Value to send to MARBL as xco2", &
142 default=284.317, units="ppm")
143 endif
144 call get_param(param_file, mdl, "ATM_ALT_CO2_OPT", atm_co2_opt, &
145 "Source of alternate atmospheric CO2 [constant, diagnostic, or prognostic]", &
146 default="constant")
147 select case (trim(atm_co2_opt))
148 case("prognostic")
149 cs%atm_alt_co2_iopt = atm_co2_prognostic_iopt
150 case("diagnostic")
151 cs%atm_alt_co2_iopt = atm_co2_diagnostic_iopt
152 case("constant")
153 cs%atm_alt_co2_iopt = atm_co2_constant_iopt
154 case DEFAULT
155 write(err_message, "(3A)") "'", trim(atm_co2_opt), "' is not a valid ATM_ALT_CO2_OPT value"
156 call mom_error(fatal, err_message)
157 end select
158 if (cs%atm_alt_co2_iopt == atm_co2_constant_iopt) then
159 call get_param(param_file, mdl, "ATM_ALT_CO2_CONST", cs%atm_alt_co2_const, &
160 "Value to send to MARBL as xco2_alt_co2", &
161 default=284.317, units="ppm")
162 endif
163
164 ! Register diagnostic fields for outputing forcing values
165 ! These fields are posted from convert_driver_fields_to_forcings(), and they are received
166 ! in physical units so no conversion is necessary here.
167 cs%diag_ids%atm_fine_dust = register_diag_field("ocean_model", "ATM_FINE_DUST_FLUX_CPL", &
168 cs%diag%axesT1, & ! T=> tracer grid? 1 => no vertical grid
169 day, "ATM_FINE_DUST_FLUX from cpl", "kg/m^2/s")
170 cs%diag_ids%atm_coarse_dust = register_diag_field("ocean_model", "ATM_COARSE_DUST_FLUX_CPL", &
171 cs%diag%axesT1, & ! T=> tracer grid? 1 => no vertical grid
172 day, "ATM_COARSE_DUST_FLUX from cpl", "kg/m^2/s")
173 cs%diag_ids%atm_bc = register_diag_field("ocean_model", "ATM_BLACK_CARBON_FLUX_CPL", &
174 cs%diag%axesT1, & ! T=> tracer grid? 1 => no vertical grid
175 day, "ATM_BLACK_CARBON_FLUX from cpl", "kg/m^2/s")
176
177 cs%diag_ids%ice_dust = register_diag_field("ocean_model", "SEAICE_DUST_FLUX_CPL", &
178 cs%diag%axesT1, & ! T=> tracer grid? 1 => no vertical grid
179 day, "SEAICE_DUST_FLUX from cpl", "kg/m^2/s")
180 cs%diag_ids%ice_bc = register_diag_field("ocean_model", "SEAICE_BLACK_CARBON_FLUX_CPL", &
181 cs%diag%axesT1, & ! T=> tracer grid? 1 => no vertical grid
182 day, "SEAICE_BLACK_CARBON_FLUX from cpl", "kg/m^2/s")
183
184 end subroutine marbl_forcing_init
185
186 ! Note: ice fraction and u10_sqr are handled in mom_surface_forcing because of CFCs
187 subroutine convert_driver_fields_to_forcings(atm_fine_dust_flux, atm_coarse_dust_flux, &
188 seaice_dust_flux, atm_bc_flux, seaice_bc_flux, &
189 nhx_dep, noy_dep, atm_co2_prog, atm_co2_diag, &
190 afracr, swnet_afracr, ifrac_n, &
191 swpen_ifrac_n, Time, G, US, i0, j0, fluxes, CS)
192
193 real, dimension(:,:), pointer, intent(in) :: atm_fine_dust_flux !< atmosphere fine dust flux from IOB
194 !! [kg m-2 s-1]
195 real, dimension(:,:), pointer, intent(in) :: atm_coarse_dust_flux !< atmosphere coarse dust flux from IOB
196 !! [kg m-2 s-1]
197 real, dimension(:,:), pointer, intent(in) :: seaice_dust_flux !< sea ice dust flux from IOB [kg m-2 s-1]
198 real, dimension(:,:), pointer, intent(in) :: atm_bc_flux !< atmosphere black carbon flux from IOB
199 !! [kg m-2 s-1]
200 real, dimension(:,:), pointer, intent(in) :: seaice_bc_flux !< sea ice black carbon flux from IOB
201 !! [kg m-2 s-1]
202 real, dimension(:,:), pointer, intent(in) :: afracr !< open ocean fraction [1]
203 real, dimension(:,:), pointer, intent(in) :: nhx_dep !< NHx flux from atmosphere [kg m-2 s-1]
204 real, dimension(:,:), pointer, intent(in) :: noy_dep !< NOy flux from atmosphere [kg m-2 s-1]
205 real, dimension(:,:), pointer, intent(in) :: atm_co2_prog !< Prognostic atmospheric CO2 concentration
206 !! [ppm]
207 real, dimension(:,:), pointer, intent(in) :: atm_co2_diag !< Diagnostic atmospheric CO2 concentration
208 !! [ppm]
209 real, dimension(:,:), pointer, intent(in) :: swnet_afracr !< shortwave flux * open ocean fraction
210 !! [W m-2]
211 real, dimension(:,:,:), pointer, intent(in) :: ifrac_n !< per-category ice fraction [1]
212 real, dimension(:,:,:), pointer, intent(in) :: swpen_ifrac_n !< per-category shortwave flux * ice fraction
213 !! [W m-2]
214 type(time_type), intent(in) :: time !< The time of the fluxes, used for
215 !! interpolating the salinity to the
216 !! right time, when it is being
217 !! restored.
218 type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
219 type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
220 integer, intent(in) :: i0 !< i index offset
221 integer, intent(in) :: j0 !< j index offset
222 type(forcing), intent(inout) :: fluxes !< MARBL-specific forcing fields
223 type(marbl_forcing_cs), pointer, intent(inout) :: cs !< A pointer that is set to point to
224 !! control structure for MARBL forcing
225
226 integer :: i, j, is, ie, js, je, m
227 real :: atm_fe_bioavail_frac !< Fraction of iron from the atmosphere available for biological uptake [1]
228 real :: seaice_fe_bioavail_frac !< Fraction of iron from sea ice available for biological uptake [1]
229 ! Note: following two conversion factors are used to both convert from km m-2 s-1 -> mmol m-2 s-1
230 !! AND cast in MOM6's unique dimensional consistency scaling system [conc Z T-1]
231 real :: iron_flux_conversion !< Factor to convert iron flux from kg m-2 s-1 -> mmol m-3 (m s-1)
232 !! [s m2 kg-1 conc Z T-1 ~> mmol kg-1]
233 real :: ndep_conversion !< Factor to convert nitrogen deposition from kg m-2 s-1 -> mmol m-3 (m s-1)
234 !! [s m2 kg-1 conc Z T-1 ~> mmol kg-1]
235
236 if (.not. cs%use_marbl_tracers) return
237
238 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
239 ndep_conversion = (1.e6/14.) * (us%m_to_Z * us%T_to_s)
240 iron_flux_conversion = (1.e6 / molw_fe) * (us%m_to_Z * us%T_to_s)
241
242 ! Post fields from coupler to diagnostics
243 ! TODO: units from diag register are incorrect; we should be converting these in the cap, I think
244 if (cs%diag_ids%atm_fine_dust > 0) &
245 call post_data(cs%diag_ids%atm_fine_dust, atm_fine_dust_flux(is-i0:ie-i0,js-j0:je-j0), &
246 cs%diag, mask=g%mask2dT(is:ie,js:je))
247 if (cs%diag_ids%atm_coarse_dust > 0) &
248 call post_data(cs%diag_ids%atm_coarse_dust, atm_coarse_dust_flux(is-i0:ie-i0,js-j0:je-j0), &
249 cs%diag, mask=g%mask2dT(is:ie,js:je))
250 if (cs%diag_ids%atm_bc > 0) &
251 call post_data(cs%diag_ids%atm_bc, atm_bc_flux(is-i0:ie-i0,js-j0:je-j0), cs%diag, &
252 mask=g%mask2dT(is:ie,js:je))
253 if (cs%diag_ids%ice_dust > 0) &
254 call post_data(cs%diag_ids%ice_dust, seaice_dust_flux(is-i0:ie-i0,js-j0:je-j0), cs%diag, &
255 mask=g%mask2dT(is:ie,js:je))
256 if (cs%diag_ids%ice_bc > 0) &
257 call post_data(cs%diag_ids%ice_bc, seaice_bc_flux(is-i0:ie-i0,js-j0:je-j0), cs%diag, &
258 mask=g%mask2dT(is:ie,js:je))
259
260 do j=js,je ; do i=is,ie
261 ! Nitrogen Deposition
262 fluxes%nhx_dep(i,j) = (g%mask2dT(i,j) * ndep_conversion) * nhx_dep(i-i0,j-j0)
263 fluxes%noy_dep(i,j) = (g%mask2dT(i,j) * ndep_conversion) * noy_dep(i-i0,j-j0)
264 enddo ; enddo
265
266 ! Atmospheric CO2
267 select case (cs%atm_co2_iopt)
268 case (atm_co2_prognostic_iopt)
269 if (associated(atm_co2_prog)) then
270 do j=js,je ; do i=is,ie
271 fluxes%atm_co2(i,j) = g%mask2dT(i,j) * atm_co2_prog(i-i0,j-j0)
272 enddo ; enddo
273 else
274 call mom_error(fatal, &
275 "ATM_CO2_OPT = 'prognostic' but atmosphere is not providing this field")
276 endif
277 case (atm_co2_diagnostic_iopt)
278 if (associated(atm_co2_diag)) then
279 do j=js,je ; do i=is,ie
280 fluxes%atm_co2(i,j) = g%mask2dT(i,j) * atm_co2_diag(i-i0,j-j0)
281 enddo ; enddo
282 else
283 call mom_error(fatal, &
284 "ATM_CO2_OPT = 'diagnostic' but atmosphere is not providing this field")
285 endif
286 case (atm_co2_constant_iopt)
287 do j=js,je ; do i=is,ie
288 fluxes%atm_co2(i,j) = g%mask2dT(i,j) * cs%atm_co2_const
289 enddo ; enddo
290 end select
291
292 ! Alternate Atmospheric CO2
293 select case (cs%atm_alt_co2_iopt)
294 case (atm_co2_prognostic_iopt)
295 if (associated(atm_co2_prog)) then
296 do j=js,je ; do i=is,ie
297 fluxes%atm_alt_co2(i,j) = g%mask2dT(i,j) * atm_co2_prog(i-i0,j-j0)
298 enddo ; enddo
299 else
300 call mom_error(fatal, &
301 "ATM_ALT_CO2_OPT = 'prognostic' but atmosphere is not providing this field")
302 endif
303 case (atm_co2_diagnostic_iopt)
304 if (associated(atm_co2_diag)) then
305 do j=js,je ; do i=is,ie
306 fluxes%atm_alt_co2(i,j) = g%mask2dT(i,j) * atm_co2_diag(i-i0,j-j0)
307 enddo ; enddo
308 else
309 call mom_error(fatal, &
310 "ATM_ALT_CO2_OPT = 'diagnostic' but atmosphere is not providing this field")
311 endif
312 case (atm_co2_constant_iopt)
313 do j=js,je ; do i=is,ie
314 fluxes%atm_alt_co2(i,j) = g%mask2dT(i,j) * cs%atm_co2_const
315 enddo ; enddo
316 end select
317
318 ! Dust flux
319 if (associated(atm_fine_dust_flux)) then
320 do j=js,je ; do i=is,ie
321 fluxes%dust_flux(i,j) = (us%kg_m2s_to_RZ_T * g%mask2dT(i,j)) * &
322 ((atm_fine_dust_flux(i-i0,j-j0) + atm_coarse_dust_flux(i-i0,j-j0)) + &
323 seaice_dust_flux(i-i0,j-j0))
324 enddo ; enddo
325 endif
326
327 if (associated(atm_bc_flux)) then
328 do j=js,je ; do i=is,ie
329 ! TODO: abort if atm_fine_dust_flux and atm_coarse_dust_flux are not associated?
330 ! Contribution of atmospheric dust to iron flux
331 if (atm_coarse_dust_flux(i-i0,j-j0) < &
332 cs%dust_ratio_thres * atm_fine_dust_flux(i-i0,j-j0)) then
333 atm_fe_bioavail_frac = cs%fe_bioavail_frac_offset + cs%dust_ratio_to_fe_bioavail_frac * &
334 (cs%dust_ratio_thres - atm_coarse_dust_flux(i-i0,j-j0) / atm_fine_dust_flux(i-i0,j-j0))
335 else
336 atm_fe_bioavail_frac = cs%fe_bioavail_frac_offset
337 endif
338
339 ! Contribution of atmospheric dust to iron flux
340 fluxes%iron_flux(i,j) = (atm_fe_bioavail_frac * &
341 (cs%iron_frac_in_atm_fine_dust * atm_fine_dust_flux(i-i0,j-j0) + &
342 cs%iron_frac_in_atm_coarse_dust * atm_coarse_dust_flux(i-i0,j-j0)))
343
344 ! Contribution of atmospheric black carbon to iron flux
345 fluxes%iron_flux(i,j) = fluxes%iron_flux(i,j) + (cs%atm_bc_fe_bioavail_frac * &
346 (cs%atm_fe_to_bc_ratio * atm_bc_flux(i-i0,j-j0)))
347
348 seaice_fe_bioavail_frac = atm_fe_bioavail_frac
349 ! Contribution of seaice dust to iron flux
350 fluxes%iron_flux(i,j) = fluxes%iron_flux(i,j) + (seaice_fe_bioavail_frac * &
351 (cs%iron_frac_in_seaice_dust * seaice_dust_flux(i-i0,j-j0)))
352
353 ! Contribution of seaice black carbon to iron flux
354 fluxes%iron_flux(i,j) = fluxes%iron_flux(i,j) + (cs%seaice_bc_fe_bioavail_frac * &
355 (cs%seaice_fe_to_bc_ratio * seaice_bc_flux(i-i0,j-j0)))
356
357 ! Unit conversion (kg m-2 s-1 -> conc Z T-1)
358 fluxes%iron_flux(i,j) = (g%mask2dT(i,j) * iron_flux_conversion) * fluxes%iron_flux(i,j)
359
360 enddo ; enddo
361 endif
362
363 ! Per ice-category forcings
364 ! If the cap receives per-category fields, memory should be allocated in fluxes
365 if (associated(ifrac_n)) then
366 do j=js,je ; do i=is,ie
367 fluxes%fracr_cat(i,j,1) = min(1., afracr(i-i0,j-j0))
368 fluxes%qsw_cat(i,j,1) = swnet_afracr(i-i0,j-j0)
369 do m=1,size(ifrac_n, 3)
370 fluxes%fracr_cat(i,j,m+1) = min(1., ifrac_n(i-i0,j-j0,m))
371 fluxes%qsw_cat(i,j,m+1) = swpen_ifrac_n(i-i0,j-j0,m)
372 enddo
373 where (fluxes%fracr_cat(i,j,:) > 0.)
374 fluxes%qsw_cat(i,j,:) = fluxes%qsw_cat(i,j,:) / fluxes%fracr_cat(i,j,:)
375 elsewhere
376 fluxes%fracr_cat(i,j,:) = 0.
377 fluxes%qsw_cat(i,j,:) = 0.
378 endwhere
379 fluxes%fracr_cat(i,j,:) = g%mask2dT(i,j) * fluxes%fracr_cat(i,j,:)
380 fluxes%qsw_cat(i,j,:) = (us%W_m2_to_QRZ_T * g%mask2dT(i,j)) * fluxes%qsw_cat(i,j,:)
381 enddo ; enddo
382 endif
383
385
386end module marbl_forcing_mod