GNOM

A Global Neodymium Ocean Model.

This repository holds the code and data for the global steady-state model of the marine neodymium (Nd) cycle as described in Pasquier, Hines, et al. (2022).

This ReadMe serves as documentation for running the GNOM model in Julia. Thanks to Julia's excellent built-in package manager and the DataDeps.jl package, no setup should be required apart from installing Julia, downloading the GEOTRACES dataset, and activating the GNOM environment.

Check the table of contents (there should be a button on the top left) to navigate this file efficiently.

GNOM v2 does not really improve on GNOM v1 except for compatibility with Julia v1.9. Incrementing the version number seemed appropriate since the GNOM v2 update fixes many breaking changes from updated dependencies.

Installation

1. Download and install the latest stable version of Julia from julialang.org (v1.9 as of 17 May 2023).

2. Download or clone this repository on your local machine.

3. (Optional) If you want to use the observational dataset used in GNOM v1.0, you must manually download the neodymium data from the GEOTRACES Intermediate Data Product 2017 (IDP17). (This is the only piece of data that cannot be programmatically downloaded for this project.) Head over to www.geotraces.org to download the IDP21 discrete sample dataset and, following the GEOTRACES.jl recommendation, you should preferably save the NetCDF file locally at

   $HOME/Data/GEOTRACES/GEOTRACES_IDP2021/discrete_sample_data/netcdf/GEOTRACES_IDP2021_Discrete_Sample_Data.nc # name may not be accurate

   where $HOME is your "home" directory. All other data required for running the GNOM and its optimization (including ocean-circulation models, non-IDP17 Nd data, World Ocean Atlas data, and so on) are downloaded programmatically during the model setup.

Single model run

Run the optimized GNOM v1.0

The GNOM v1.0 output was generated by optimizing a bunch of parameters (see Pasquier, Hines, et al. (2022) for parameter values and details). You can reproduce the same output in just a few seconds on your laptop by following these simple 3 steps:

1. Go to (cd) to the GNOM folder.

2. Start Julia

3. In Julia, type

   julia> include("src/Nd_model/single_run.jl")

This sets the model up with the optimized parameters and runs a single simulation of the GNOM, which should take a few seconds to a few minutes, depending if it's the first time that you are running the file. (Unless there is an ERROR thrown out, you can safely ignore the warnings and information printed out.)

The single run output

The src/Nd_model/single_run.jl file creates two variables, DNd and εNd, which contain the vectors for the Nd concentration (in mol/m⁻³) and ε_Nd (unitless). (It also saves the parameters used in a newly created folder in output/single_runs/ for safekeeping.)

Because the εNd vector is the last computed variable, the output should end with something like

200160-element Vector{Float64}:
 -0.0006020528949600701
 -0.000595097968841718
  ⋮
 -0.0013770493809266426
 -0.001364787766450326

although the actual numerical values might be different.

These are the ε_Nd values that you can convert to parts per ten thousand (‱) by typing, e.g.,

julia> εNd .|> per10000
200160-element Vector{Quantity{Float64, NoDims, Unitful.FreeUnits{(‱,), NoDims, nothing}}}:
  -6.020528949600701 ‱
 -5.9509796884171795 ‱
                     ⋮
 -13.770493809266426 ‱
  -13.64787766450326 ‱

where per10000 is the ‱ unit.

Similarly, you can convert DNd to pM by doing:

julia> DNd * mol/m^3 .|> pM
200160-element Vector{Quantity{Float64, 𝐍 𝐋⁻³, Unitful.FreeUnits{(pM,), 𝐍 𝐋⁻³, nothing}}}:
  1.450921127576513 pM
 1.4497239997396394 pM
                     ⋮
 38.646836627631856 pM
 34.017114010235964 pM

where here we additionally had to apply the default unit (mol m⁻³) before converting. (The units used here, e.g., mol and pM, have been preloaded from Unitful.jl.)

Exploring different parameters values

You can edit the parameter values in src/Nd_model/single_run.jl and run it again to simulate the [Nd] and ε_Nd fields for different parameter values.

If you encounter any errors, please share what you did and the stacktrace in a GitHub issue and we will try to troubleshoot the issue with you as fast as possible.

Plotting

There many packages for plotting in Julia. Here we just document two of them: Plots.jl and Makie.jl.

with Plots.jl

The most popular plotting package is Plots.jl and AIBECS.jl provides recipes for it. To use the Plots.jl package, simply type

julia> using Plots

The vectors of tracers, such as DNd and εNd, can be easily plotted with the AIBECS recipes. For example, after running the single run above, you can use

julia> plotverticalmean(εNd .|> per10000, grd, c=:balance)

which will output something like

For more plot types, see the AIBECS documentation, which contains example tutorials and how-to guides with simple plots.

with Makie.jl

Reproducing figures from the GNOM paper

After your single run, you can reproduce the figures of the GNOM paper by running the scripts located in src/plots/GMDpaper/. For example, to reproduce the ε_Nd transects (Fig. 10 in Pasquier, Hines, et al. (2022)), you can type

julia> include("src/plots/GMDpaper/transects_eNd.jl")

and it will produce a figure similar to

Note that Makie.jl comes with multiple backends, of which GNOM includes GLMakie (that I would use for interactive plot windows) and CairoMakie (that I would use for PDFs). The Makie backend is set to GLMakie in src/Nd_model/single_run.jl by setting the variable use_GLMakie = true. If you want to generate PDFs, use CairoMakie by editing the line to use_GLMakie = false. (If you switch between GLMakie and CairoMakie, you must rerun the load.jl file, which you can do by simply setting reload = true before running any of the figure scripts.)

Why Makie and how to plot without the Plots.jl recipes?

The figures in Pasquier, Hines, et al. (2022) were created with Makie.jl because it provides finer control to create detailed publication-quality PDFs. Although AIBECS.jl does not provide recipes for Makie.jl, there are a number of AIBECS functions to rearrange the 1D column vectors into 3D and take slices/averages over given regions/depths. For example, in line with the Plots.jl recipe used above for the vertical mean of ε_Nd, you can get the 2D array of the vertical mean with

julia> verticalmean(εNd .|> per10000, grd)
91×180 Matrix{Quantity{Float64, NoDims, Unitful.FreeUnits{(‱,), NoDims, nothing}}}:
      NaN ‱      NaN ‱  …       NaN ‱       NaN ‱
      NaN ‱      NaN ‱          NaN ‱       NaN ‱
          ⋮             ⋱
 -11.7113 ‱  -11.709 ‱     -11.7166 ‱  -11.7138 ‱
      NaN ‱      NaN ‱  …       NaN ‱       NaN ‱

(where NaNs represent land/dry model grid cells.)

Such functions are extensively used by the plotting scripts in src/plots/GMDpaper/, which you can directly use to reproduce the plots in Pasquier, Hines, et al. (2022).

Most of these plotting functions have been updated with GNOM v2 for compatibility with Julia v1.9 and its improved latency for first plots.

Optimization

To run the optimization, you can type

include("src/Nd_model/setup_and_optimization.jl")

Beware that this will take a few hours on a modern laptop!

This optimization script will randomize the initial value for the parameter values by taking a random sample from the prior distributions determined by the initial guess and the range defined in the parameter type (the Params struct in the model_setup.jl file).

Si model

When running the Nd-cycling single run, the optimized Si-cycling fields required for opal scavenging are automatically downloaded from FigShare. However, if you want, you can also edit the Si-cycle code and re-optimize it.

To run the Si-model optimization on your laptop, call

include("src/Si_model/setup_and_optimization.jl")

Just like for the Nd-cycle model, you can modify the Si-cycle model in src/Si_model/model_setup.jl

Running on a SLURM cluster

The GNOM v1.0 outputs for Pasquier, Hines, et al. (2022) were generated using a SLURM cluster, for which batch scripts were written in this repository.

If you have access to a cluster with the SLURM workload manager, you should be able to run the optimizations by

1. Cloning this repository on your cluster,

2. Have Julia installed,

3. And after cding to the GNOM directory, typing

   sbatch src/slurm/optimize_Nd.sh

This SLURM batch file will request 1 node with 64GB for 20 hours.

Similarly, there is also a batch script for optimizing the Si-cycle model, src/slurm/optimize_Nd.sh.

Citation

To cite the GNOM v1 model, please cite Pasquier, Hines, et al. (2022). You can use the journal's citation tools, directly grab the contents of our CITATION.bib file or simply use:

Pasquier, B., Hines, S. K. V., Liang, H., Wu, Y., Goldstein, S. L., and John, S. G.: GNOM v1.0: an optimized steady-state model of the modern marine neodymium cycle, Geosci. Model Dev., 15, 4625–4656, https://doi.org/10.5194/gmd-15-4625-2022, 2022.

TODO

• Update observations scripts. Currently doiwnloads post IDP17 data, but should be using IDP21 at the same time, so there will be some (many?) duplicates. Also there has been no QC on the IDP21 from the GNOM authors.
• Check that the optimization routines still work with GNOM v2 (updated for Julia v1.9)
• Improve the conecptual model.