API Reference

Base

IMASdd.infoFunction
info(uloc::AbstractString, extras::Bool=true)

Return information of a node in the IMAS data structure, possibly including extra structures

source
function info(ids_type::Type, field::Symbol)

Return information of a filed of an IDS

source
IMASdd.unitsFunction
units(uloc::String)
source
units(@nospecialize(ids::IDS), field::Symbol)

Return string with units for a given IDS field

source
IMASdd.coordinatesFunction
coordinates(ids::IDS, field::Symbol; override_coord_leaves::Union{Nothing,Vector{<:Union{Nothing,Symbol}}}=nothing)

Return a vector of Coordinate with the .ids and .field filled to point at the coordinate entries in the dd

if field === :- then there's no coordinate

Use override_coord_leaves to override fetching coordinates of a given field

NOTE: getproperty(coords[X]) value is nothing when the data does not have a coordinate

getproperty(coords[X]) Coordinate value is `missing` if the coordinate is missing in the data structure
source
IMASdd.time_coordinate_indexFunction
time_coordinate_index(@nospecialize(ids::IDS), field::Symbol; error_if_not_time_dependent::Bool)

Return index of time coordinate

If error_if_not_time_dependent == false it will return 0 for arrays that are not time dependent

source
IMASdd.time_coordinateFunction
time_coordinate(@nospecialize(ids::IDS), field::Symbol)

Return a vector of Coordinate with the .ids and .field filled to point at the time coordinate of the field

source
IMASdd.access_logConstant
IMASdd.access_log

IMASdd.access_log.enable = true / false

@show IMASdd.access_log

empty!(IMASdd.access_log) # to reset

Track access to the data dictionary

source
Base.getpropertyFunction
Base.getproperty(ids::Union{IDSraw, IDSvectorRawElement}, field::Symbol)

No processing for IDSraw and IDSvectorRawElement

source
getproperty(ids::IDS, field::Symbol; to_cocos::Int=user_cocos)

Return IDS value for requested field

source
getproperty(ids::IDS, field::Symbol, @nospecialize(default::Any); to_cocos::Int=user_cocos)

Return IDS value for requested field or default if field is missing

NOTE: This is useful because accessing a missing field in an IDS would raise an error

source
Base.isemptyFunction
isempty(@nospecialize(ids::IDSvector))

returns true if IDSvector is empty

source
isempty(@nospecialize(ids::IDS); include_expr::Bool=false, eval_expr::Bool=false)

Returns true if none of the IDS fields downstream have data (or expressions)

NOTE: By default it does not include nor evaluate expressions

source
isempty(@nospecialize(ids::IDS), field::Symbol; include_expr::Bool=false, eval_expr::Bool=false)

Returns true if the ids field has no data (or expression)

NOTE: By default it does not include nor evaluate expressions

source
Base.setproperty!Function
Base.setproperty!(@nospecialize(ids::IDS), field::Symbol, value; skip_non_coordinates::Bool=false, error_on_missing_coordinates::Bool=true)
source
Base.setproperty!(@nospecialize(ids::IDS), field::Symbol, value::AbstractArray{<:IDS}; skip_non_coordinates::Bool=false, error_on_missing_coordinates::Bool=true)

Handle setproperty of entire vectors of IDS structures at once (ids.field is of type IDSvector)

source
Base.setproperty!(@nospecialize(ids::IDS), field::Symbol, value::AbstractArray; skip_non_coordinates::Bool=false, error_on_missing_coordinates::Bool=true)

Ensures coordinates are set before the data that depends on those coordinates.

If skip_non_coordinates is set, then fields that are not coordinates will be silently skipped.

source
IMASdd.indexFunction
index(@nospecialize(ids::IDSvectorElement))

Returns index of the IDSvectorElement in the parent IDSvector

source
IMASdd.keys_no_missingFunction
keys_no_missing(@nospecialize(ids::IDS); include_expr::Bool=true, eval_expr::Bool=false)

Returns generator of fields with data in a IDS

NOTE: By default it includes expressions, but does not evaluate them. It assumes that a IDStop without data will also have no valid expressions.

source
Base.resize!Function
Base.resize!(@nospecialize(ids::IDSvector{T}), condition::Pair{String}, conditions::Pair{String}...; wipe::Bool=true, error_multiple_matches::Bool=true) where {T<:IDSvectorElement}

Resize if a set of conditions are not met.

If wipe=true and an entry matching the condition is found, then the content of the matching IDS is emptied.

Either way, the IDS is populated with the conditions.

NOTE: error_multiple_matches will delete all extra entries matching the conditions.

Returns the selected IDS

source
resize!(
    @nospecialize(ids::IDSvector{T}),
    identifier_name::Symbol,
    conditions::Pair{String}...;
    wipe::Bool=true,
    error_multiple_matches::Bool=true
)::T where {T<:IDSvectorElement}

Resize ids if identifier_name is not found based on index of index_2_name(ids) and a set of conditions are not met.

If wipe=true and an entry matching the condition is found, then the content of the matching IDS is emptied.

Either way, the IDS is populated with the conditions.

NOTE: error_multiple_matches will delete all extra entries matching the conditions.

Returns the selected IDS

source
resize!(@nospecialize(ids::IDSvector{T}); wipe::Bool=true) where {T<:IDSvectorTimeElement}

Resize time dependent array at global_time

source
resize!(@nospecialize(ids::IDSvector{T}), time0::Float64; wipe::Bool=true) where {T<:IDSvectorTimeElement}

Resize time dependent array based on time

source
Base.deleteat!Function
Base.deleteat!(@nospecialize(ids::T), condition::Pair{String}, conditions::Pair{String}...) where {T<:IDSvector}

If an entry matching the condition is found, then the content of the matching IDS is emptied

source
deleteat!(@nospecialize(ids::T), identifier_name::Symbol, conditions::Pair{String}...)::T where {T<:IDSvector}

Deletes all entries that match based on index of index_2_name(ids)

source
Base.ismissingFunction
Base.ismissing(@nospecialize(ids::IDS), field::Symbol)

returns true/false if field is missing in IDS

source
Base.diffFunction
Base.diff(
    @nospecialize(ids1::T),
    @nospecialize(ids2::T);
    tol::Float64=1E-2,
    recursive::Bool=true,
    verbose::Bool=false) where {T<:IDS}

Compares two IDSs and returns dictionary with differences

NOTE: This function does not evaluate expressions (use freeze() on the IDSs to compare values instead of functions)

source
IMASdd.top_idsFunction
top_ids(@nospecialize(ids::Union{IDS,IDSvector}))

Return top-level IDS in the hierarchy and nothing if top level is not a top-level IDS

source
IMASdd.top_ddFunction
top_dd(@nospecialize(ids::Union{IDS,IDSvector}))

Return top-level dd in the hierarchy, and nothing if top level is not dd

source
Base.parentFunction
parent(@nospecialize(ids::Union{IDS,IDSvector}); IDS_is_absolute_top::Bool=false, error_parent_of_nothing::Bool=true)

Return parent IDS/IDSvector in the hierarchy

If IDS_is_absolute_top=true then returns nothing instead of dd()

If error_parent_of_nothing=true then asking parent(nothing) will just return nothing

source
IMASdd.gotoFunction
goto(@nospecialize(ids::Union{IDS,IDSvector}), loc_fs::String)

Reach location in a given IDS

NOTE: loc_fs is the path expressed in fs format

source
goto(@nospecialize(ids::Union{IDS,IDSvector}), path::Union{AbstractVector,Tuple})

Reach location in a given IDS

source
IMASdd.leavesFunction
leaves(@nospecialize(ids::IDS))

Returns iterator with (filled) leaves in the IDS

source
IMASdd.filled_ids_fieldsFunction
filled_ids_fields(@nospecialize(ids::IDS); eval_expr::Bool=false)

Returns a vector with tuples pointing to all the (ids, field) that have data downstream

source
IMASdd.pathsFunction
paths(@nospecialize(ids::IDS); eval_expr::Bool=false)

Returns the locations in the IDS that have data downstream

source
IMASdd.selective_copy!Function
selective_copy!(@nospecialize(h_in::IDS), @nospecialize(h_out::IDS), path::Vector{<:AbstractString}, time0::Float64)

Copies the content of a path from one IDS to another (if the path exists) at a given time0

NOTE:

  • the path is a i2p(ulocation)
  • if time0 is NaN then all times are retained
source
IMASdd.selective_delete!Function
selective_delete!(@nospecialize(h_in::IDS), path::Vector{<:AbstractString})

Deletes a path from one IDS

NOTE:

  • the path is a i2p(ulocation)
source
IMASdd.@findallMacro
@findall ids :symbol
@findall ids r"Regular Expression"
@findall [ids1, ids2] [:sybmol1, :symbol2]
@findall [ids1, ids2] r"Regular Expression"

Searches for specified fields within single/multiple IDS objects, while capturing their names into IDSFieldFinder.root_name

See also findall(root_ids::Union{IDS,IDSvector}, target::Union{Symbol,AbstractArray{Symbol},Regex}=r""; kwargs) which this macro calls after expansion.

Arguments

  • `root_ids:: Root IDS objects to search.
  • `target_fields::Union{Symbol, AbstractArray{Symbol}, Regex}: Fields to search for, specified by a single symbol, array of symbols, or regular expression.

Returns

  • Vector{IDS_Field_Finder}: A vector of IDS_Field_Finder structures, each containing details on a located field such as parent IDS, root IDS, field type, and full field path.

Example

julia> @findall [dd1, dd2] [:psi, :j_tor]
julia> @findall [dd1, dd2] r"psi"

julia> eqt = dd.equilibrium.time_slice[]

julia> @findall eqt :psi
julia> @findall eqt r"global.*psi"
source
Base.findallFunction
findall(ids::Union{IDS, IDSvector}, target_fields::Union{Symbol,AbstractArray{Symbol},Regex}=r""; include_subfields::Bool=true)

Searches for specified fields within IDS objects, supporting nested field exploration and customizable filtering.

Arguments

  • root_ids::Union{IDS, IDSvector}: Root IDS objects to search.
  • target_fields::Union{Symbol, AbstractArray{Symbol}, Regex} = r"": Fields to search for, specified by a single symbol, array of symbols, or regular expression.
  • include_subfields::Bool = true: If true, retrieves nested fields below the target field when found; if false, stops at the matching field.

Returns

  • Vector{IDS_Field_Finder}: A vector of IDS_Field_Finder structures, each containing details on a located field such as parent IDS, root IDS, field type, and full field path.

Example

julia> findall(dd.equilibrium.time_slice[].global_quantities) # By default, it searches everything under given IDS
julia> findall(dd.equilibrium.time_slice[].global_quantities, r"") # Same behavior (Default)

# Find fields matching a single symbol within a IDS structure
julia> IFF = findall(dd.equilibrium.time_slice, :psi)

# Search for multiple symbols within multiple root IDS objects
julia> IFF = findall([dd.equilibrium, dd.core_profiles], [:psi, :j_tor])

# Use regular expressions for flexible and powerful search patterns
julia> IFF = findall(dd, r"prof.*1d.*psi")

# Control subfield inclusion using the `include_subfields` keyword
julia> IFF = findall(dd, r"prof.*2d"; include_subfields=false)
julia> IFF = findall(dd, r"prof.*2d"; include_subfields=true) # Default behavior

# Default show for IFF (IDF_Field_Finder) structure
julia> IFF

# Retrieve actual values of found IDS objects (lazy evaluation)
julia> IFF[1].value
julia> IFF[end].value
source
findall(identifier_name::Symbol, @nospecialize(ids::IDSvector))

Return items from IDSvector based on index of index_2_name(ids)

source

IO

IMASdd.file2imasFunction
file2imas(filename::AbstractString; kw...)

Load IDS from a file that can be in different formats .json or .h5 both ITER tensorized (h5i) or OMAS hierarchical (hdf)

source
IMASdd.is_h5iFunction
is_h5i(filename::AbstractString)

Returns true if a file is in ITER tensorized (h5i) format

source
IMASdd.dict2imasFunction
dict2imas(dct::AbstractDict, @nospecialize(ids::IDS); show_warnings::Bool=true)

Populate IMAS data structure ids based on data contained in Julia dictionary dct.

source
IMASdd.imas2dictFunction
imas2dict(ids::Union{IDS,IDSvector}; freeze::Bool=false, strict::Bool=false)

Populate Julia structure of dictionaries and vectors with data from IMAS data structure ids

source
IMASdd.json2imasFunction
json2imas(filename::AbstractString, @nospecialize(ids::IDS)=dd(); error_on_missing_coordinates::Bool=true, show_warnings::Bool=true)

Load the IMAS data structure from a JSON file with given filename

source
IMASdd.jstr2imasFunction
jstr2imas(json_string::String, @nospecialize(ids::IDS)=dd(); error_on_missing_coordinates::Bool=true, show_warnings::Bool=true)

Load the IMAS data structure from a JSON string

source
IMASdd.imas2jsonFunction
imas2json(@nospecialize(ids::Union{IDS,IDSvector}), filename::AbstractString; freeze::Bool=false, strict::Bool=false, indent::Int=0, kw...)

Save the IMAS data structure to a JSON file with given filename.

Arguments

  • freeze evaluates expressions
  • strict dumps fields that are strictly in ITER IMAS only
  • kw... arguments are passed to the JSON.print function
source
IMASdd.hdf2imasFunction
hdf2imas(filename::AbstractString, target_path::AbstractString; error_on_missing_coordinates::Bool=true, verbose::Bool=false, kw...)

Load an object from an HDF5 file using a simple entry point. Given a file and an internal target path (as a string), the function returns either the dataset’s value or an IMAS ids structure constructed from a group.

If the object at target_path is a group and has a "concretetype" attribute, that type is evaluated and instantiated; otherwise, a default (dd()) is used. Coordinate data is processed based on `erroronmissingcoordinates`.

Arguments

  • filename: Path to the HDF5 file
  • target_path: Internal HDF5 path to the desired dataset or group
  • error_on_missing_coordinates (default: true): Enforce coordinate checks
  • verbose (default: false): Enable verbose logging
  • kw...: Additional keyword arguments passed to HDF5.h5open

Returns

The value of the dataset or the constructed IMAS ids.

source
hdf2imas(filename::AbstractString; error_on_missing_coordinates::Bool=true, kw...)

Load data from a HDF5 file generated by OMAS Python platform (ie. hierarchical HDF5)

kw... arguments are passed to the HDF5.h5open function

source
IMASdd.hdf2dict!Function
hdf2dict!(gparent::Union{HDF5.File,HDF5.Group}, ids::AbstractDict)

Load data from a HDF5 file into a dictionary

source
IMASdd.imas2hdfFunction
imas2hdf(@nospecialize(ids::Union{IDS,IDSvector}), filename::AbstractString;
         mode::String="w", target_group::String="/", overwrite::Bool=false,
         freeze::Bool=false, strict::Bool=false, desc::String="", kw...)

Save an IMAS data structure to an OMAS HDF5 file.

Arguments:

  • filename: HDF5 file path
  • mode: File open mode ("w", "a", or "r+"); "a" is converted to "r+"
  • target_group: Group where data will be stored (default is "/")
  • overwrite: If true, overwrite the target group if it exists
  • show_warnings: If true, display warn messages
  • freeze evaluates expressions
  • strict dumps fields that are strictly in ITER IMAS only
  • desc: description of additional information (e.g., Shot number)
  • compress: compression level, an integer between 0 (no compression) and 9 (highest)
  • kw...: Options passed to the internal dispatch

Returns:

The result of imas2hdf(ids, gparent; freeze, strict, desc).

source
IMASdd.h5i2imasFunction
h5i2imas(filename::AbstractString, @nospecialize(ids::IDS)=dd(); kw...)

Load data from a HDF5 file generated by IMAS platform (ie. tensorized HDF5)

kw... arguments are passed to the HDF5.h5open function

source
IMASdd.imas2h5iFunction
imas2h5i(
    @nospecialize(ids::Union{IDS,IDSvector}),
    filename::AbstractString;
    freeze::Bool=false,
    strict::Bool=false,
    run::Int=0,
    shot::Int=0,
    hdf5_backend_version::String="1.0",
    kw...
)

Save data to a HDF5 file generated by IMAS platform (ie. tensorized HDF5)

Arguments

  • kw... arguments are passed to the HDF5.h5open function
  • freeze evaluates expressions
  • strict dumps fields that are strictly in ITER IMAS only
  • run, shot, hdf5_backend_version arguments are used to set the HDF5 attributes
source
IMASdd.h5mergeFunction
h5merge(output_file::AbstractString, keys_files::Union{AbstractDict{<:AbstractString,<:AbstractString},AbstractVector{<:Pair{<:AbstractString,<:AbstractString}};
      mode::AbstractString="a", skip_existing_entries::Bool=false,
      h5_group_search_depth::Integer=0, h5_strip_group_prefix::Bool=false,
      verbose::Bool=false)

Merges multiple files into a single HDF5 output file.

Arguments

  • output_file: Path to the HDF5 file where data will be merged.

  • keys_files: A vector or dictionary mapping target group names to input filenames.

  • mode: "w" to create a new file or "a" to append to an existing one.

  • skip_existing_entries: If true, groups already present in the output file are not overwritten.

  • h5_group_search_depth: For input HDF5 files, the depth at which to collect group paths.

    • 0 means use the root ("/").
    • 1 means collect immediate children of the root.
    • Higher values collect groups deeper in the hierarchy.
  • h5_strip_group_prefix: If true, the target group name (the key from keys_files) is omitted from the output HDF5 path. For example, if an input file contains a group path "/level1/level2" and the key is "parent", then:

    • With h5_strip_group_prefix = false, the output path becomes "/parent/level1/level2".
    • With h5_strip_group_prefix = true, the output path becomes "/level1/level2".
  • verbose: If true, additional logging information is printed.

Behavior

  • For input files with the .h5 extension, the function opens the file and collects group paths up to h5_group_search_depth. Each collected path is modified by stripping the first N components (using "/" as the delimiter) according to the flag h5_strip_group_prefix (if true, the parent key is omitted). Then, the corresponding objects are copied into the output file.
  • For other file types (e.g., JSON, YAML, text, markdown), the file is read and stored as text or raw binary data.
  • The function records attributes for each copied group that indicate the original file paths.

Returns a vector of group names (as strings) that were processed.

source
h5merge(
    output_file::AbstractString,
    directories::AbstractVector{<:AbstractString};
    include_base_dir::Bool=true,
    cleanup::Bool=false,
    kwargs...)

Add all files in multiple directories (and their subdirectories) to an HDF5 output_file

source
h5merge(
    output_file::AbstractString,
    directory::AbstractString;
    mode::AbstractString="a",
    skip_existing_entries::Bool=false,
    follow_symlinks::Bool=false,
    verbose::Bool=false,
    include_base_dir::Bool=false,
    pattern::Union{Regex,Nothing}=nothing,
    kwargs...
    )

Add all files in a directory (and subdirectories) to an HDF5 output_file

source
IMASdd.read_combined_h5Function
read_combined_h5(filename::AbstractString; show_warnings::Bool=true, error_on_missing_coordinates::Bool=true, pattern::Regex=r"", kw...)

Iteratively traverse an HDF5 file from the root ("/") using a stack.

Arguments

  • filename: Path to the combined HDF5 file
  • error_on_missing_coordinates (default true): Enforce coordinate checks during dispatch
  • pattern (default r""): A regex used to filter which paths are processed
  • kw...: Additional keyword arguments passed to hdf2imas

Returns

  • Dict{String,Any}: loaded data with keys (path as string)
source

Time

IMASdd.global_timeFunction
global_time(ids::Union{IDS,IDSvector})

Get the dd.global_time of a given IDS

If top-level dd cannot be reached then returns Inf

source
global_time(ids::Union{IDS,IDSvector}, time0::Float64)

Set the dd.global_time of a given IDS

source
IMASdd.set_time_arrayFunction
set_time_array(@nospecialize(ids::IDS), field::Symbol, value)

Set value of a time-dependent array at the dd.global_time

source
set_time_array(@nospecialize(ids::IDS), field::Symbol, time0::Float64, value)

Set value of a time-dependent array at time0

NOTE: updates the closest causal element of an array

source
IMASdd.get_time_arrayFunction
get_time_array(@nospecialize(ids::IDS{T}), field::Symbol, scheme::Symbol=:constant) where {T<:Real}

Get data from a time-dependent array at the dd.global_time

source
get_time_array(@nospecialize(ids::IDS), field::Symbol, time0::Float64, scheme::Symbol=:constant)

Get data from time dependent array

NOTE: logic for @ddtime array handling:

  • interpolation (i) scheme between array bounds
  • constant (c) extrapolation within bounds of time array
  • error (e) when time0 is before minimum(time)

For example:

time:   -oooo-
data:   -o-o--
ddtime: eiiicc
source
IMASdd.@ddtimeMacro
@ddtime( X.Y )

Get data from time dependent array. Equivalent to:

get_time_array(X, :Y)

and

@ddtime( X.Y = V)

Set data in a time dependent array. Equivalent to:

set_time_array(X, :Y, V)
source
IMASdd.last_timeFunction
last_time(dd::DD)

Returns the last time referenced in all the IDSs dd.XXX.time vectors (including dd.global_time)

source
IMASdd.last_global_timeFunction
last_global_time(dd::DD)

Returns the last time referenced in all the IDSs dd.XXX.time vectors (including dd.global_time)

source
IMASdd.new_timeslice!Function
new_timeslice!(@nospecialize(ids::IDS), time0::Float64=global_time(ids))

Recursively appends a deepcopy at time time0 of the last time-slice of all time-dependent array structures under a given ids

source
new_timeslice!(@nospecialize(ids::IDS{T}), times::AbstractVector{Float64}) where {T<:Real}

Extend IDSvector{<:IDSvectorTimeElement} and time dependent data arrays with times

source
Missing docstring.

Missing docstring for resize!. Check Documenter's build log for details.

IMASdd.retime!Function
retime!(@nospecialize(ids::IDS), time0::Float64=global_time(ids))

Recursively change the time of the last time-slices or last time-depedent vector elements in a IDS

source
IMASdd.get_timesliceFunction
get_timeslice(@nospecialize(ids::IDS), time0::Float64=global_time(ids), scheme::Symbol=:constant; slice_pulse_schedule::Bool=true)

Returns data at the given time0 (by default at the global_time)

Data is selected from time dependent arrays of structures using closest causal time point.

Data is selected from time dependent arrays using these possible schemes [:constant, :linear, :quadratic, :cubic, :pchip, :lagrange]

source
get_timeslice(el_type::Type{Z}, @nospecialize(ids::IDS), time0::Float64=global_time(ids), scheme::Symbol=:constant; slice_pulse_schedule::Bool=false) where {Z<:Real}

gettimeslice that retuns IDS of type `eltype`

source
IMASdd.copy_timeslice!Function
copy_timeslice!(
    @nospecialize(ids0::IDS{T1}),
    @nospecialize(ids::IDS{T2}),
    time0::Float64,
    scheme::Symbol;
    slice_pulse_schedule::Bool) where {T1<:Real,T2<:Real}

Copy data at a given time from ids to ids0

source
IMASdd.trim_time!Function
trim_time!(@nospecialize(ids::IDS); trim_pulse_schedule::Bool=false)

Recursively remove all time dependent data tha occurs after global_time

source
trim_time!(@nospecialize(ids::IDS), time_range::Tuple{Float64,Float64}; trim_pulse_schedule::Bool=false)

Recursively remove all time dependent data tha occurs outside of time_range

source
IMASdd.time_dependent_leavesFunction
time_dependent_leaves(ids::IDS{T}) where {T<:Real}

Returns Dict{String,Vector{IMASnodeRepr{T}}} mapping time coordinate locations to vectors of data fields that use that time coordinate.

NOTE: Excludes :time fields and error fields ending in "_σ"

source
IMASdd.time_groupsFunction
time_groups(ids::IDS{T}; min_channels::Int) where {T<:Real}

Groups identical time vectors and optionally filters by minimum group size.

Returns Vector{Vector{IMASnodeRepr{T}}} containing groups of time fields that share identical time arrays, keeping only groups with at least min_channels members.

source

Math

IMASdd.interp1dFunction
interp1d(x, y, scheme::Symbol=:linear)

One dimensional curve interpolations with scheme [:constant, :linear, :quadratic, :cubic, :pchip, :lagrange]

NOTE: this interpolation method will extrapolate

source
IMASdd.extrap1dFunction
extrap1d(itp::DataInterpolations.AbstractInterpolation; first=:extrapolate, last=:extrapolate) where {T<:Real}

first and last can be [:extrapolate, :constant, :nan, --value--] affect how the extrapolation is done at the either end of the array

source
IMASdd.gradientFunction
gradient(coord::AbstractVector{C}, arr::AbstractVector{A}; method::Symbol=:second_order) where {C<:Real, A<:Real}

The finite difference gradient. The returned gradient has the same shape as the input array.

method of the gradient can be one of [:backward, :central, :forward, :secondorder, :thirdorder]

For :central the gradient is computed using second order accurate central differences in the interior points and first order accurate one-sides (forward or backward) differences at the boundaries.

For :second_order the gradient is computed using second order accurate central differences in the interior points, and 2nd order differences at the boundaries.

For :third_order the gradient is computed from the cubic spline passing through the points

source
gradient(coord1::AbstractVector, coord2::AbstractVector, mat::Matrix, dim::Int; method::Symbol=:second_order)

Finite difference method of the gradient: [:backward, :central, :forward, :secondorder, :thirdorder]

Can be applied to either the first (dim=1) or second (dim=2) dimension

source
gradient(coord1::AbstractVector, coord2::AbstractVector, mat::Matrix; method::Symbol=:second_order)

Finite difference method of the gradient: [:backward, :central, :forward, :secondorder, :thirdorder]

Computes the gradient in both dimensions

source