Applying an external Si (or P) flux in MARBL

[pdinezio]

I am trying to force a B case with active ocean BGC with a separate Si or P flux in MARBL.

I have been experimenting with dust_flux_source = 'monthly-calendar' and modifying the input file in dust_flux_input%filename. However, this does not let me change Si or P separately. I was browsing the source code and found that the coupler can initialize a variable ext_Si_flux in ecosys_forcing_mod.F90.

I was wondering if anyone had an example of how to read ext_Si_flux and ext_P_flux or conversely how to read a separate Si or P flux into MARBL in addition to the dust flux.

Thank you so much!

Setup:
Running a B1850C4L45BGCRBPRP compset using MARBL with CESM2_1_5 release (git describe: release-cesm2.1.5-0-g7a6c5b0)

dust_flux_source = 'monthly-calendar'
iron_flux_source = 'monthly-calendar'

 

[mlevy]

The ext_Si_flux forcing field is used when running with ladjust_bury_coeff = .true. to tune the various burial coefficients (POC, POP, and bSi burial coefficients are typically set to a constant, but the constant is determined from a long spin-up where the coefficients can vary in time). The forcing in that case just comes from the river runoff, which is handled entirely by POP:

                else if (index == ext_Si_flux_ind) then
                   riv_flux_ind = dsi_riv_flux_ind
                   if (riv_flux_ind > 0) then
                     forcing_field%field_0d(:,:,iblock) = riv_flux_forcing_fields(riv_flux_ind)%field_0d(:,:,iblock)
                   else
                     forcing_field%field_0d(:,:,iblock) = c0
                   end if

Which brings me to your second question - we don't have a mechanism for reading a separate Si or P flux (aside from the runoff fluxes, which probably isn't what you want to modify). I think the best suggestion is to copy marbl_surface_flux_mod.F90 into SourceMods and change the hard-coded parameters used to convert dust to P and Si fluxes:

    !-----------------------------------------------------------------------
    !  Add phosphate and silicate from dust after Krishnamurthy et al. (2010)
    !  factors convert from g/cm2/s to nmol/cm2/s
    !  ( P frac in dust by weight) * ( P solubility) / ( P molecular weight) * (mol->nmol)
    !  (Si frac in dust by weight) * (Si solubility) / (Si molecular weight) * (mol->nmol)
    !-----------------------------------------------------------------------

    surface_fluxes(:, po4_ind) = surface_fluxes(:, po4_ind)   + (dust_flux_in * (0.00105_r8 *  0.15_r8 / 30.974_r8 * 1.0e9_r8))

    surface_fluxes(:, sio3_ind) = surface_fluxes(:, sio3_ind) + (dust_flux_in * (  0.308_r8 * 0.075_r8 / 28.085_r8 * 1.0e9_r8))

 

[pdinezio]

Hi Mike - thanks for the quick reply. Do you have any guidance on how we modify marbl code to read two variables from .dust_flux_input%filename ?

I spent some time understanding ecosys_forcing_mod.F90 - but could not figure out how to initialize another variable like dust_flux_in so we can feed dust for Si and dust for P into the model.

I hope this makes sense. Thanks!

 

[mlevy]

Do you have any guidance on how we modify marbl code to read two variables from .dust_flux_input%filename

This one would be tricky, because it requires modifying both POP and MARBL. Unfortunately I don't see great documentation on the process for the MARBL side, though there is some discussion on what POP would need. On the MARBL side, you would need to do the following:

1. Add a new index to the surface_flux_indexing_type and make sure it gets initialized
2. Update marbl_init_surface_flux_forcings() to include a block of code to initialize your new forcing
3. Modify the code in marbl_surface_flux_mod.F90 that actually uses the code to point to your new index instead of the original dust flux. Pay attention to the unit conversion portion of the block of code I commented in my previous comment; you'll probably still need to convert from g/cm2/s to nmol/cm2/s.
4. Modify POP according to the first link in this comment to read in the new forcing; you can probably copy the dust_flux block and make sure the variable name you use in step [2] matches the variable you add to dust_flux_input%filename

I'm sure I'm missing some important details, but if you give those steps a try and report back, I'm happy to help troubleshoot any errors you run into along the way. Any chance you are putting this run together on Derecho? It'll be easy to help if I have access to your case and run directories, but if you're running this on a university cluster I'm sure we'll figure something out.