MITgcmTools.jl
Set of tools for running MITgcm, analyzing its output, and/or modifying its inputs. A set of Pluto.jl notebooks, which e.g. run MITgcm
interactively, can be found in the examples/
folder.
Example Guide
To load one of the notebooks using Pluto.jl
:
- open
julia
in terminal window - type the following commands at the
Julia
prompt - in web-browser, open one of the notebooks listed hereafter
cd("examples/")
using Pluto
Pluto.run()
Examples / Running Models
MITgcm_configurations.jl
: explore MITgcm configurations and their parameters.MITgcm_worklow.jl
: build, setup, run, and plot for any standard configuration.run_MITgcm.jl
: a more detailed look into compiling and running the model.
Examples / Analyzing Results
HS94_plotmap.jl
: readhs94.cs-32x32x5
output, interpolate, and plot mapHS94_particles.jl
: compute particle trajectories fromhs94.cs-32x32x5
outputHS94_Makie.jl
: example usingMakie.jl
instead ofPlots.jl
Explore And Run MITgcm
The recommended, simple, method to run the model is via the climate model interface (see docs@ClimateModels.jl for detail). The MITgcm_launch
function can be used to run a MITgcm
configuration after setting up the MITgcm_config
. Using this interface facilitates operations like compiling and setting up a temporary folder to run the model. Key functions, incl. the climate model interface, are documented further down in the docs.
The verification_experiments
function provides a list of the most-standard MITgcm configurations that can all be run either in batch mode or interactively. The MITgcm_path
variable points to where MITgcm is compiled. Interactive / reactive notebooks are found in the examples/
folder (e.g. run_MITgcm.jl
seen just below).
MITgcmTools.testreport
— Functiontestreport(nam::String,ext="")
Run the testreport script for one model config nam
(or "all"), with additional options (optional) speficied in ext
using MITgcmTools
testreport(MITgcm_config(configuration="front_relax"),"-norun")
#testreport(MITgcm_config(configuration="all"),"-norun")
MITgcmTools.verification_experiments
— Functionverification_experiments()
Get list of all most-standard
configurations of MITgcm
and return as an Array of MITgcm_config
using MITgcmTools
exps=verification_experiments()
MITgcmTools.MITgcm_namelist
— TypeMITgcm_namelist(groups,params)
Data structure representing a MITgcm namelist file, such as data.pkg
, which contains
groups :: Array{Symbol,1} = Array{Symbol,1}(undef, 0)
params :: Array{OrderedDict{Symbol,Any},1} = Array{OrderedDict{Symbol,Any},1}(undef, 0)
with model parameters (params
) being organized in groups
as done in the files.
using MITgcmTools
fil=joinpath(MITgcm_path,"verification","advect_xy","run","data")
nml=read_namelist(fil)
MITgcm_namelist(nml.groups,nml.params)
MITgcm_namelist(groups=nml.groups,params=nml.params)
MITgcm_namelist(groups=nml.groups)
ClimateModels / MITgcm interface
MITgcmTools.MITgcm_config
— TypeMITgcm_config()
Concrete type of AbstractModelConfig
for MITgcm
which contains
model :: String = "MITgcm"
configuration :: String = ""
options :: Array{String,1} = Array{String,1}(undef, 0)
inputs :: Array{String,1} = Array{String,1}(undef, 0)
outputs :: Array{String,1} = Array{String,1}(undef, 0)
status :: Array{String,1} = Array{String,1}(undef, 0)
channel :: Channel{Any} = Channel{Any}(10)
folder :: String = tempdir()
ID :: UUID = UUIDs.uuid4()
and with defaults that can be constructed as follows for example
using MITgcmTools
tmp=MITgcm_config()
exps=verification_experiments()
exps[end]
(part of the climate model interface as specialized for MITgcm
)
Missing docstring for build
. Check Documenter's build log for details.
ClimateModels.compile
— Functioncompile(config::MITgcm_config)
Compile the model using make
in build/
that has already been setup.
(part of the climate model interface as specialized for MITgcm
)
ClimateModels.setup
— Functionsetup(config::MITgcm_config)
Create a run/
folder and link everything there as needed to be ready to run model as normally done for most-standard MITgcm configurations (incl. prepare_run
and mitgcmuv
). Call init_git_log(config)
to setup git tracker and put!(config.channel,MITgcm_launch)
to be executed via launch(config)
later.
(part of the climate model interface as specialized for MITgcm
)
MITgcmTools.MITgcm_launch
— FunctionMITgcm_launch(config::MITgcm_config)
Go to run/
folder and effectively call mitgcmuv > output.txt
(part of the climate model interface as specialized for MITgcm
)
ClimateModels.clean
— Functionclean(config::MITgcm_config)
Cancel any remaining task (config.channel), clean up the build directory (via testreport), and clean the run directory (via rm).
(part of the climate model interface as specialized for MITgcm
)
Reading MITgcm outputs
MITgcmTools.read_mdsio
— Functionread_mdsio(datafile)
Read a MITgcm
mdsio file (".data" binary + ".meta" text pair), and return as an Array
p="./hs94.cs-32x32x5/run/"
x=read_mdsio(p*"surfDiag.0000000020.002.001.data")
y=read_mdsio(p*"pickup.ckptA.002.001.data")
z=read_mdsio(p*"T.0000000000.002.001.data")
read_mdsio(pth::String,fil::String)
Read a MITgcm
's MDSIO files (".data" binary + ".meta" text pair), combine, and return as an Array
p="./hs94.cs-32x32x5/run/"
x=read_mdsio(p,"surfDiag.0000000020")
y=read_mdsio(p,"pickup.ckptA")
z=read_mdsio(p,"T.0000000000")
MITgcmTools.read_meta
— Functionread_meta(metafile)
Read a MITgcm
metadata file, parse it, and return as a NamedTuple
p="./hs94.cs-32x32x5/run/"
meta=read_meta(p*"surfDiag.0000000020.002.001.meta")
pairs(meta)
meta.dimList
read_meta(pth::String,fil::String)
Read a MITgcm
metadata files, parse them, and return as an array of NamedTuple
p="./hs94.cs-32x32x5/run/"
meta=read_meta(p,"surfDiag.0000000020")
pairs(meta[end])
[meta[i].dimList for i in 1:length(meta)]
MITgcmTools.read_namelist
— Functionread_namelist(fil)
Read a MITgcm
namelist file, parse it, and return as a NamedTuple
using MITgcmTools
testreport("advect_xy")
fil=joinpath(MITgcm_path,"verification","advect_xy","run","data")
namelist=read_namelist(fil)
MITgcmTools.write_namelist
— Functionwrite_namelist(fil)
Save a MITgcm
namelist file. In the example below, one is read from file, modified, and then saved to a new file using write_namelist.
using MITgcmTools
fil=joinpath(MITgcm_path,"verification","advect_xy","run","data")
nml=read_namelist(fil)
write_namelist(fil*"_new",namelist)
or
nml=read(fil,MITgcm_namelist())
write(fil*"_new",nml)
MITgcmTools.read_available_diagnostics
— Functionread_available_diagnostics(fldname::String; filename="available_diagnostics.log")
Get the information for a particular variable fldname
from the available_diagnostics.log
text file generated by MITgcm
.
MITgcmTools.read_bin
— Functionread_bin(fil::String,kt::Union{Int,Missing},kk::Union{Int,Missing},prec::DataType,mygrid::gcmgrid)
Read model output from binary file and convert to MeshArray. Other methods:
read_bin(fil::String,prec::DataType,mygrid::gcmgrid)
read_bin(fil::String,mygrid::gcmgrid)
MITgcmTools.read_flt
— Functionread_flt(dirIn::String,prec::DataType)
Read displacements from MITgcm/pkg/flt output file into a DataFrame.
MITgcmTools.read_nctiles
— Functionread_nctiles(fileName,fldName,mygrid)
Read model output from NCTiles file and convert to MeshArray instance.
mygrid=GridSpec("LatLonCap")
fileName="nctiles_grid/GRID"
Depth=read_nctiles(fileName,"Depth",mygrid)
hFacC=read_nctiles(fileName,"hFacC",mygrid)
hFacC=read_nctiles(fileName,"hFacC",mygrid,I=(:,:,1))
Format conversions
MITgcmTools.findtiles
— Functionfindtiles(ni::Int,nj::Int,mygrid::gcmgrid)
findtiles(ni::Int,nj::Int,grid::String="LatLonCap",GridParentDir="../inputs/GRID_LLC90/")
Return a MeshArray
map of tile indices, mytiles["tileNo"]
, for tile size ni,nj
and extract grid variables accordingly.
MITgcmTools.cube2compact
— Functioncube2compact(x::Array)
Reshape from e.g. size (192, 32, 5) in cube format to (32, 192, 5) in compact format.
MITgcmTools.compact2cube
— Functioncompact2cube(x::Array)
Reshape from e.g. size (32, 192, 5) in cube format to (192, 32, 5) in compact format.
MITgcmTools.convert2array
— Functionconvert2array(fld::MeshArray)
Convert MeshArray to Array (or vice versa otherwise)
MITgcmTools.convert2gcmfaces
— Functionconvert2gcmfaces(fld::MeshArray)
Convert mitgcm output to MeshArray (or vice versa otherwise)
Formulae etc
MITgcmTools.SeaWaterDensity
— FunctionSeaWaterDensity(Θ,Σ,Π,Π0)
Compute potential density (ρP), in situ density (ρI), and density referenced to PREF (Π0 in decibars) from potential temperature (Θ in °C), salinity (Σ in psu) and pressure (Π in decibars) according to the UNESCO / Jackett & McDougall 1994 equation of state.
Credits: code based on a Matlab implementation by B. Ferron Reference: https://www.jodc.go.jp/info/iocdoc/UNESCOtech/059832eb.pdf Check value: ρI = 1041.83267kg/m^3
for Θ=3°Celcius
, Σ=35psu
, Π=3000dbar
(ρP,ρI,ρR) = SeaWaterDensity(3.,35.5,3000.)
isapprox(ρI,1041.83267, rtol=1e-6)
MITgcmTools.MixedLayerDepth
— FunctionMixedLayerDepth(Θ,Σ,Δ,mthd)
Compute mixed layer depth from potential temperature (Θ in °C), salinity (Σ in psu) and depth (Δ in method) according to various formulas (mthd == "BM", "Suga", "Kara"). Inputs must be dense vectors without any missing value (or NaN, etc).
D=collect(0.0:1.0:500.0); tmp=(1.0.-tanh.(5*(-1 .+ 2/D[end]*D)));
T=2.0 .+ 8.0*tmp; S=34.0 .+ 0.5*tmp;
(ρP,ρI,ρR) = SeaWaterDensity(T,S,D);
mld=MixedLayerDepth(T,S,D,"BM"); isapprox(mld,134.0)
using Plots
plot(ρP,-D,w=2,label="Potential Density",ylabel="Depth")
plot!(vec([ρP[1] ρP[end]]),-fill(mld,2),label="Mixed Layer Depth",w=2,c="black",s=:dash)
Index
MITgcmTools.MITgcm_config
MITgcmTools.MITgcm_namelist
ClimateModels.clean
ClimateModels.compile
ClimateModels.setup
MITgcmTools.MITgcm_launch
MITgcmTools.MixedLayerDepth
MITgcmTools.SeaWaterDensity
MITgcmTools.compact2cube
MITgcmTools.convert2array
MITgcmTools.convert2gcmfaces
MITgcmTools.cube2compact
MITgcmTools.findtiles
MITgcmTools.read_available_diagnostics
MITgcmTools.read_bin
MITgcmTools.read_flt
MITgcmTools.read_mdsio
MITgcmTools.read_meta
MITgcmTools.read_namelist
MITgcmTools.read_nctiles
MITgcmTools.testreport
MITgcmTools.verification_experiments
MITgcmTools.write_namelist