the Creative Commons Attribution 4.0 License.
the Creative Commons Attribution 4.0 License.
| 22 Jul 2025
A high-resolution nested model to study the effects of alkalinity additions in Halifax Harbour, a mid-latitude coastal fjord
Abstract. Surface ocean alkalinity enhancement (OAE), through the release of alkaline materials, is an emerging marine carbon dioxide removal technology that could increase the storage of anthropogenic carbon in the ocean. Observations collected during recent and on-going field trials will provide important information on the feasibility and effects of alkalinity additions on carbon cycling and study ecological responses. However, given the scales involved (24/7 continuous addition, meters to hundreds/thousands of kilometers and minutes to months for alkalinity dispersion) observations, even with the use of autonomous platforms, will remain inherently sparse and limited. Alone, they cannot provide a comprehensive quantification of the effects of OAE on the carbonate system, and ultimately of the net air-sea CO2 fluxes. Numerical models, informed and validated by field observations, are therefore essential to OAE deployments and the measurement, reporting, and verification (MRV) of any resulting carbon uptake. They can help guide fieldwork design, including optimal design of measurement monitoring networks, provide forecasts of the ocean state, simulate the effects of alkalinity additions on the seawater carbonate system, and allow one to quantify net CO2 uptake. Here we describe a coupled physical-biogeochemical model that is specifically designed for coastal OAE. The model is an implementation of the Regional Ocean Modelling System (ROMS) in a nested grid configuration with increasing spatial resolution from the Scotian Shelf to Halifax Harbour (coastal fjord, eastern Canada), a current test site for operational alkalinity addition. The biogeochemical model simulates oxygen dynamics, carbonate system processes (including air-sea gas exchange), and feedstock properties (dissolution, sinking). We present a multi-year hindcast validated against the long-term weekly time series available for a long-term monitoring station at the deepest part of Halifax Harbour, as well as alkalinity addition simulations at various locations inside and outside the harbour to show the model’s capabilities for assessing the effects of OAE at this coastal site.
- Preprint
(12244 KB) - Metadata XML
-
Supplement
(16893 KB) - BibTeX
- EndNote
Status: open (until 24 Sep 2025)
-
RC1: 'Comment on egusphere-2025-3361', Anonymous Referee #1, 23 Aug 2025
reply
General comments:
The authors present a nested regional ocean model of Halifax Harbour and part of the Scotian shelf which is validated against measurements. A simple dissolution model is implemented and pulse releases of an alkaline effluent are modelled, consisting of a mix of dissolved and particulate alkalinity. The subsequent changes in alkalinity and DIC (from the induced CO2 uptake) are evaluated and analyzed.
Overall the manuscript is well laid out, focused and easy to follow. The simulations presented establish an important standard of rigor for future OAE deployments in other areas. I recommend publication.
Specific comments
The authors show that alkalinity addition inside a natural enclosed harbour enables a substantial fraction of the theoretically maximal CO2 uptake to occur quickly and within the simulation domain, due to the long residence time and relatively shallow waters. As pointed out in L556-559, this makes MRV much easier both experimentally and from a simulation perspective. Of course the flipside of this is that a confined body of water which does not quickly spread any added ∆TA over large ocean areas will also limit the total sustained alkalinity addition rate in that area, limiting scaling of OAE.
It would be useful to add an estimation of this in the manuscript. For a rough, first pass estimate, perhaps one could assume that the response of ∆TA and ∆DIC are roughly additive and linear with respect to addition rate. Then, for each of the three locations, one could calculate what the maximum addition rate would be which would raise the maximal ∆pH to some acceptable limit (what that limit is is of course arbitrary, but perhaps something conservative like +0.1 or +0.05 units would be illustrative).
Another approach would be perhaps to examine the export rate of alkalinity out of the simulation boundary and try to estimate what sustained alkalinity addition rate (rather than a pulse) could be achieved, again within some ∆pH or ∆TA limit set within the domain.
A discussion of this and the tradeoffs of release locations would be useful to the reader to understand better what sort of scale OAE can achieve.L317 k_{diss}TA_p term:
The treatment of dissolution as an exponential decay process (i.e. dTAp/dt = -k TAp) was surprising at first glance. Usually dissolution of particular matter is treated with a shrinking core model, where the dissolution rate has units of mol cm-2 s-1, the radius of particles shrinks linearly and fully dissolves in a finite amount of time. For a very narrow (as indicated in L335, “a particle size of 12µm”) or uniform distribution of particle sizes I believe an exponential dissolution curve is only a mediocre fit.
I can see that an exponential model could perhaps capture the behaviour of a gaussian or log-normal distribution of particle sizes, but a short discussion of this and a justification of the choice of model here would be helpful.
L317 w_{p}TA_{p} term:
It’s unclear to me how the sinking term is applied. As written it looks like there is an exponential decay, i.e. each time step some fraction of TA_p is lost to sinking from any given simulation grid voxel. What happens to that TA_p ? Does it get added to the cell below, until the bottom cell is reached after which it disappears in to the sediment ? Or does the model assume the sunk particles are removed completely (i.e. they sink out entirely at a rate of W_p*TA_p from anywhere in the column ?). As currently written it seems more like it’s the latter, as there is no term that accounts for sinking particles that arrive from a cell above (i was expecting a second term like +w_p*TA_p^{z=i-1} )
Please clarify how the sinking mechanism is implemented and justify its construction.
The sinking rate is stated as 5.5 m^{-1} later (L337) but that can’t be w_p since the units wouldn’t be right (w_p should have units of inverse time, like k_{diss}). How is w_p calculated from the 5.5m^{-1} ?
L326 The treatment of sediment loss in layer N is a little unclear. It says a term is “added” to ∂∆TA/∂t ? Or does this replace the regular dissolution term in ∂∆TA/∂t (last term in Equation 9) ? It might be clearer here to just rewrite the full Equation 9 (and perhaps Equation 8) in the case of the bottom cell, for clarity.
It’s also confusing to me that the loss of TAp due to sinking/burial is already explicitly treated in equation 8 using w_p and then it’s treated again here with the \theta_{loss} term. Is \theta_{loss} a constant ? Or is it calculated from w_p ?
L424ff The comparison of H2 and H3 is very interesting and suggests perhaps a resolution as high as H3 isn’t necessary. A similar comparison of H1 vs H2 would also be useful if the releases can be reasonably implemented at the coarsest level. Even if the release location would have to be assumed to be wider or poorly matched in terms of exact location, injection of the same amount of alkalinity in the coarsest model could be interesting to determine to what extent the H2 level is required.
L769 It was a surprise to read here that the sediment loss term was set to zero. I feel like this should have been mentioned earlier, perhaps even right when the loss term(s) are introduced in L317ff. Is both wp and \theta_{loss} set to zero or just the latter ? If it’s just the latter, does the model currently just settle all the particles on the floor and let them dissolve from there until completely dissolved ?
Technical corrections:
L120: I assume the conversion factor is 1025 kg m^-3, not 1.025kg m^-3 (remove dot or change dot to comma)
L243 In equation (3), it appears that the parameter “c1” is duplicate as a coefficient to t and as an exponent. Likely it is meant to be c2 instead ?
L325 change to “is added that mimics” or “is added to mimic”
L331 “1.29 ml s-1”, exponentiate the “-1”
L475 In such cases,
Fig.1D consider using a different color scheme for the bathymetry as the scale is different.
Figs. 3, 5,6,7, 10: Is it possible to indicate the release location in these plots with a small black arrow or similar. I know they are shown in Fig 1 D, but it would be very helpful to have that info on each of the other plots too.
Figure 7: It would be nice to add a horizontal dashed line to the two graphs indicating the theoretical maximum uptake (at your CO2 efficiency of 0.89) to get a sense for what fraction of the ultimate uptake occurs within the simulation domains.
Fig S4-S8 The observations of the depth profiles are sparse enough in time that it’s difficult to assess visually how closely the corresponding model predictions match. Perhaps, for each observation time and depth simply make a scatter plot against the corresponding prediction value ? Could be color coded by depth perhaps to see if correlation is better at surface vs depth.
Citation: https://doi.org/10.5194/egusphere-2025-3361-RC1
Viewed
| HTML | XML | Total | Supplement | BibTeX | EndNote | |
|---|---|---|---|---|---|---|
| 346 | 99 | 11 | 456 | 29 | 5 | 9 |
- HTML: 346
- PDF: 99
- XML: 11
- Total: 456
- Supplement: 29
- BibTeX: 5
- EndNote: 9
Viewed (geographical distribution)
| Country | # | Views | % |
|---|
| Total: | 0 |
| HTML: | 0 |
| PDF: | 0 |
| XML: | 0 |
- 1