API Reference
Base
IMASdd.info
— Functioninfo(uloc::AbstractString, extras::Bool=true)
Return information of a node in the IMAS data structure, possibly including extra structures
function info(ids_type::Type, field::Symbol)
Return information of a filed of an IDS
IMASdd.units
— Functionunits(uloc::String)
units(@nospecialize(ids::IDS), field::Symbol)
Return string with units for a given IDS field
IMASdd.coordinates
— Functioncoordinates(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
IMASdd.time_coordinate_index
— Functiontime_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
IMASdd.time_coordinate
— Functiontime_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
IMASdd.access_log
— ConstantIMASdd.access_log
IMASdd.access_log.enable = true / false
@show IMASdd.access_log
empty!(IMASdd.access_log) # to reset
Track access to the data dictionary
Base.getproperty
— FunctionBase.getproperty(ids::Union{IDSraw, IDSvectorRawElement}, field::Symbol)
No processing for IDSraw and IDSvectorRawElement
getproperty(ids::IDS, field::Symbol; to_cocos::Int=user_cocos)
Return IDS value for requested field
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
Base.isempty
— Functionisempty(@nospecialize(ids::IDSvector))
returns true if IDSvector is empty
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
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
IMASdd.isfrozen
— Functionisfrozen(@nospecialize(ids::IDS))
Returns if the ids has been frozen
Base.setproperty!
— FunctionBase.setproperty!(@nospecialize(ids::IDS), field::Symbol, value; skip_non_coordinates::Bool=false, error_on_missing_coordinates::Bool=true)
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)
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.
IMASdd.index
— Functionindex(@nospecialize(ids::IDSvectorElement))
Returns index of the IDSvectorElement in the parent IDSvector
IMASdd.keys_no_missing
— Functionkeys_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.
Base.resize!
— FunctionBase.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
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
resize!(@nospecialize(ids::IDSvector{T}); wipe::Bool=true) where {T<:IDSvectorTimeElement}
Resize time dependent array at global_time
resize!(@nospecialize(ids::IDSvector{T}), time0::Float64; wipe::Bool=true) where {T<:IDSvectorTimeElement}
Resize time dependent array based on time
Base.deleteat!
— FunctionBase.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
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)
Base.ismissing
— FunctionBase.ismissing(@nospecialize(ids::IDS), field::Symbol)
returns true/false if field is missing in IDS
Base.diff
— FunctionBase.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)
IMASdd.top_ids
— Functiontop_ids(@nospecialize(ids::Union{IDS,IDSvector}))
Return top-level IDS in the hierarchy and nothing
if top level is not a top-level IDS
IMASdd.top_dd
— Functiontop_dd(@nospecialize(ids::Union{IDS,IDSvector}))
Return top-level dd
in the hierarchy, and nothing
if top level is not dd
Base.parent
— Functionparent(@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
IMASdd.goto
— Functiongoto(@nospecialize(ids::Union{IDS,IDSvector}), loc_fs::String)
Reach location in a given IDS
NOTE: loc_fs is the path expressed in fs format
goto(@nospecialize(ids::Union{IDS,IDSvector}), path::Union{AbstractVector,Tuple})
Reach location in a given IDS
IMASdd.leaves
— Functionleaves(@nospecialize(ids::IDS))
Returns iterator with (filled) leaves in the IDS
IMASdd.filled_ids_fields
— Functionfilled_ids_fields(@nospecialize(ids::IDS); eval_expr::Bool=false)
Returns a vector with tuples pointing to all the (ids, field) that have data downstream
IMASdd.paths
— Functionpaths(@nospecialize(ids::IDS); eval_expr::Bool=false)
Returns the locations in the IDS that have data downstream
IMASdd.selective_copy!
— Functionselective_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
IMASdd.selective_delete!
— Functionselective_delete!(@nospecialize(h_in::IDS), path::Vector{<:AbstractString})
Deletes a path from one IDS
NOTE:
- the path is a i2p(ulocation)
IMASdd.@findall
— Macro@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 ofIDS_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"
Base.findall
— Functionfindall(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
: Iftrue
, retrieves nested fields below the target field when found; iffalse
, stops at the matching field.
Returns
Vector{IDS_Field_Finder}
: A vector ofIDS_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
findall(identifier_name::Symbol, @nospecialize(ids::IDSvector))
Return items from IDSvector based on index
of index_2_name(ids)
IO
IMASdd.file2imas
— Functionfile2imas(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)
IMASdd.is_h5i
— Functionis_h5i(filename::AbstractString)
Returns true if a file is in ITER tensorized (h5i) format
IMASdd.dict2imas
— Functiondict2imas(dct::AbstractDict, @nospecialize(ids::IDS); show_warnings::Bool=true)
Populate IMAS data structure ids
based on data contained in Julia dictionary dct
.
IMASdd.imas2dict
— Functionimas2dict(ids::Union{IDS,IDSvector}; freeze::Bool=false, strict::Bool=false)
Populate Julia structure of dictionaries and vectors with data from IMAS data structure ids
IMASdd.json2imas
— Functionjson2imas(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
IMASdd.jstr2imas
— Functionjstr2imas(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
IMASdd.imas2json
— Functionimas2json(@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 expressionsstrict
dumps fields that are strictly in ITER IMAS onlykw...
arguments are passed to theJSON.print
function
IMASdd.hdf2imas
— Functionhdf2imas(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 filetarget_path
: Internal HDF5 path to the desired dataset or grouperror_on_missing_coordinates
(default:true
): Enforce coordinate checksverbose
(default:false
): Enable verbose loggingkw...
: Additional keyword arguments passed toHDF5.h5open
Returns
The value of the dataset or the constructed IMAS ids.
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
IMASdd.hdf2dict!
— Functionhdf2dict!(gparent::Union{HDF5.File,HDF5.Group}, ids::AbstractDict)
Load data from a HDF5 file into a dictionary
IMASdd.imas2hdf
— Functionimas2hdf(@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 pathmode
: 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 existsshow_warnings
: If true, display warn messagesfreeze
evaluates expressionsstrict
dumps fields that are strictly in ITER IMAS onlydesc
: 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)
.
IMASdd.h5i2imas
— Functionh5i2imas(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
IMASdd.imas2h5i
— Functionimas2h5i(
@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 theHDF5.h5open
functionfreeze
evaluates expressionsstrict
dumps fields that are strictly in ITER IMAS onlyrun
,shot
,hdf5_backend_version
arguments are used to set the HDF5 attributes
IMASdd.h5merge
— Functionh5merge(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
: Iftrue
, 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
: Iftrue
, the target group name (the key fromkeys_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"
.
- With
verbose
: Iftrue
, additional logging information is printed.
Behavior
- For input files with the
.h5
extension, the function opens the file and collects group paths up toh5_group_search_depth
. Each collected path is modified by stripping the first N components (using "/" as the delimiter) according to the flagh5_strip_group_prefix
(iftrue
, 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.
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
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
IMASdd.read_combined_h5
— Functionread_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 fileerror_on_missing_coordinates
(defaulttrue
): Enforce coordinate checks during dispatchpattern
(defaultr""
): A regex used to filter which paths are processedkw...
: Additional keyword arguments passed tohdf2imas
Returns
Dict{String,Any}
: loaded data with keys (path as string)
Time
IMASdd.global_time
— Functionglobal_time(ids::Union{IDS,IDSvector})
Get the dd.global_time of a given IDS
If top-level dd cannot be reached then returns Inf
global_time(ids::Union{IDS,IDSvector}, time0::Float64)
Set the dd.global_time of a given IDS
IMASdd.set_time_array
— Functionset_time_array(@nospecialize(ids::IDS), field::Symbol, value)
Set value of a time-dependent array at the dd.global_time
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
IMASdd.get_time_array
— Functionget_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
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
IMASdd.@ddtime
— Macro@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)
IMASdd.last_time
— Functionlast_time(dd::DD)
Returns the last time referenced in all the IDSs dd.XXX.time
vectors (including dd.global_time
)
IMASdd.last_global_time
— Functionlast_global_time(dd::DD)
Returns the last time referenced in all the IDSs dd.XXX.time
vectors (including dd.global_time
)
IMASdd.new_timeslice!
— Functionnew_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
new_timeslice!(@nospecialize(ids::IDS{T}), times::AbstractVector{Float64}) where {T<:Real}
Extend IDSvector{<:IDSvectorTimeElement} and time dependent data arrays with times
IMASdd.retime!
— Functionretime!(@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
IMASdd.get_timeslice
— Functionget_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]
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`
IMASdd.copy_timeslice!
— Functioncopy_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
IMASdd.trim_time!
— Functiontrim_time!(@nospecialize(ids::IDS); trim_pulse_schedule::Bool=false)
Recursively remove all time dependent data tha occurs after global_time
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
IMASdd.time_dependent_leaves
— Functiontime_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 "_σ"
IMASdd.time_groups
— Functiontime_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.
Math
IMASdd.interp1d
— Functioninterp1d(x, y, scheme::Symbol=:linear)
One dimensional curve interpolations with scheme [:constant, :linear, :quadratic, :cubic, :pchip, :lagrange]
NOTE: this interpolation method will extrapolate
IMASdd.extrap1d
— Functionextrap1d(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
IMASdd.gradient
— Functiongradient(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
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
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
IMASdd.nanmaximum
— Functionnanmaximum(a::AbstractArray)
Maximum ignoring NaNs in an array
IMASdd.nanminimum
— Functionnanminimum(a::AbstractArray)
Minimum ignoring NaNs in an array