text:
function minerva_scanner(records::Dict{Int,Any}, scanFunction::Function; logger::Union{Nothing,SimpleLogger} = nothing)::VLResult

    # initialize -
    scanned_records_dictionary = Dict{Int64,Array{MinervaToken,1}}()

    try

        # connect to the database -
        path_to_database_file = joinpath(_PATH_TO_DATABASE,"Database.db")
        CFGM_DB_CONNECTION = SQLite.DB(path_to_database_file)

        # load the synonym table -
        data_op_result = load_synonym_table(CFGM_DB_CONNECTION; logger=logger)
        synonym_dictionary = check(data_op_result; logger=logger)

        # so for each record in my records array, lets scan -
        for (key, record) in records
            scan_result = scanFunction(record, synonym_dictionary; logger=logger)
            scanned_record = check(scan_result; logger=logger)
            scanned_records_dictionary[key] = scanned_record
        end

        # return -
        return VLResult(scanned_records_dictionary)
    catch error
        return VLResult(error)
    end

end

function minerva_bound_type_assignment_scan_function(record::Array{String,1}, synonym_dictionary::Union{Nothing,Dict{String,String}} = nothing;
    logger::Union{Nothing,SimpleLogger})::VLResult

    # initialize -
    canonical_record_array = Array{String,1}()
    canonical_token_array = Array{MinervaToken,1}()
    token_type_dictionary = Dict{String,MinervaToken}()

    # hardcode the token type dictionary -
    token_type_dictionary[string(hash("("))] = MinervaToken("(",LPAREN())
    token_type_dictionary[string(hash(")"))] = MinervaToken(")",RPAREN())
    token_type_dictionary[string(hash("SOURCE"))] = MinervaToken("source",SOURCE())
    token_type_dictionary[string(hash("SINK"))] = MinervaToken("sink",SINK())
    token_type_dictionary[string(hash("BOUND"))] = MinervaToken("bound",BOUND())
    token_type_dictionary[string(hash("UNBOUND"))] = MinervaToken("unbound",UNBOUND())
    token_type_dictionary[string(hash("source"))] = MinervaToken("source",SOURCE())
    token_type_dictionary[string(hash("sink"))] = MinervaToken("sink",SINK())
    token_type_dictionary[string(hash("bound"))] = MinervaToken("bound",BOUND())
    token_type_dictionary[string(hash("unbound"))] = MinervaToken("unbound",UNBOUND())
    token_type_dictionary[string(hash("is"))] = MinervaToken("is",IS())
    token_type_dictionary[string(hash("a"))] = MinervaToken("a",A())
    token_type_dictionary[string(hash("type"))] = MinervaToken("TYPE",TYPE())

    try

        # First = remove all synonyms and reduce the sentence to canonical form -
        tmp_record_array = reverse(copy(record))
        while (isempty(tmp_record_array) == false)

            # check: pop the component, check against the syn data table -
            record_component = pop!(tmp_record_array)

            # do we have this component as a key?
            if (haskey(synonym_dictionary,record_component) == true)
                # we have a synonym -
                new_record_component = synonym_dictionary[record_component]
                push!(canonical_record_array, new_record_component)
            else
                push!(canonical_record_array,record_component)
            end
        end

        # ok, so now we have a canonical forms - lets tokenize the cononical elements -
        # first, lets string match -
        reverse_canonical_record_array_copy = reverse(copy(canonical_record_array))
        while (isempty(reverse_canonical_record_array_copy) == false)

            # item -
            item_to_classify = pop!(reverse_canonical_record_array_copy)

            # check: do we have this item in the token dictionary -
            test_key = string(hash(item_to_classify))
            if (haskey(token_type_dictionary, test_key) == true)
                token_item = token_type_dictionary[test_key]
                push!(canonical_token_array, token_item)
            else
                minerva_token = MinervaToken(item_to_classify, BIOLOGICAL_SYMBOL())
                push!(canonical_token_array,minerva_token)
            end
        end

        # finally - add a SEMICOLON -
        push!(canonical_token_array, MinervaToken(";", SEMICOLON()))

        # return -
        return VLResult(canonical_token_array)

    catch error
        return VLResult(error)
    end
end

function minerva_biological_type_assignment_scan_function(record::Array{String,1}, synonym_dictionary::Union{Nothing,Dict{String,String}} = nothing;
    logger::Union{Nothing,SimpleLogger})::VLResult

    # initialize -
    canonical_record_array = Array{String,1}()
    canonical_token_array = Array{MinervaToken,1}()
    token_type_dictionary = Dict{String,MinervaToken}()

    # hardcode the token type dictionary -
    token_type_dictionary[string(hash("("))] = MinervaToken("(",LPAREN())
    token_type_dictionary[string(hash(")"))] = MinervaToken(")",RPAREN())
    token_type_dictionary[string(hash("SOURCE"))] = MinervaToken("source",SOURCE())
    token_type_dictionary[string(hash("SINK"))] = MinervaToken("sink",SINK())
    token_type_dictionary[string(hash("BOUND"))] = MinervaToken("bound",BOUND())
    token_type_dictionary[string(hash("UNBOUND"))] = MinervaToken("unbound",UNBOUND())
    token_type_dictionary[string(hash("is"))] = MinervaToken("is",IS())
    token_type_dictionary[string(hash("a"))] = MinervaToken("a",A())
    token_type_dictionary[string(hash("type"))] = MinervaToken("TYPE",TYPE())
    token_type_dictionary[string(hash("GENE"))] = MinervaToken("GENE",GENE_TYPE_SYMBOL())
    token_type_dictionary[string(hash("MESSENGER-RNA"))] = MinervaToken("MESSENGER-RNA",mRNA_TYPE_SYMBOL())
    token_type_dictionary[string(hash("TRANSFER-RNA"))] = MinervaToken("TRANSFER-RNA",tRNA_TYPE_SYMBOL())
    token_type_dictionary[string(hash("REGULATORY-RNA"))] = MinervaToken("REGULATORY-RNA",regRNA_TYPE_SYMBOL())
    token_type_dictionary[string(hash("PROTEIN"))] = MinervaToken("PROTEIN",PROTEIN_TYPE_SYMBOL())
    token_type_dictionary[string(hash("METABOLITE"))] = MinervaToken("METABOLITE",METABOLITE_TYPE_SYMBOL())
    token_type_dictionary[string(hash("RNA_POLYMERASE_II_SYMBOL"))] = MinervaToken("RNAP_SYMBOL",RNA_POLYMERASE_II_SYMBOL())
    token_type_dictionary[string(hash("RIBOSOME_SYMBOL"))] = MinervaToken("RIBOSOME_SYMBOL",RIBOSOME_SYMBOL())

    try

        # First = remove all synonyms and reduce the sentence to canonical form -
        tmp_record_array = reverse(copy(record))
        while (isempty(tmp_record_array) == false)

            # check: pop the component, check against the syn data table -
            record_component = pop!(tmp_record_array)

            # do we have this component as a key?
            if (haskey(synonym_dictionary,record_component) == true)
                # we have a synonym -
                new_record_component = synonym_dictionary[record_component]
                push!(canonical_record_array, new_record_component)
            else
                push!(canonical_record_array,record_component)
            end
        end

        # ok, so now we have a canonical forms - lets tokenize the cononical elements -
        # first, lets string match -
        reverse_canonical_record_array_copy = reverse(copy(canonical_record_array))
        while (isempty(reverse_canonical_record_array_copy) == false)

            # item -
            item_to_classify = pop!(reverse_canonical_record_array_copy)

            # check: do we have this item in the token dictionary -
            test_key = string(hash(item_to_classify))
            if (haskey(token_type_dictionary, test_key) == true)
                token_item = token_type_dictionary[test_key]
                push!(canonical_token_array, token_item)
            else
                minerva_token = MinervaToken(item_to_classify, BIOLOGICAL_TYPE_PREFIX())
                push!(canonical_token_array,minerva_token)
            end
        end

        # finally - add a SEMICOLON -
        push!(canonical_token_array, MinervaToken(";", SEMICOLON()))

        # return -
        return VLResult(canonical_token_array)
    catch error
        return VLResult(error)
    end

end

function minerva_grn_scan_function(record::Array{String,1}, synonym_dictionary::Union{Nothing,Dict{String,String}} = nothing;
    logger::Union{Nothing,SimpleLogger})::VLResult

    # initialize -
    canonical_record_array = Array{String,1}()
    canonical_token_array = Array{MinervaToken,1}()
    token_type_dictionary = Dict{String,MinervaToken}()

    # stuff -
    token_type_dictionary[string(hash("("))] = MinervaToken("(",LPAREN())
    token_type_dictionary[string(hash(")"))] = MinervaToken(")",RPAREN())
    token_type_dictionary[string(hash("transcription"))] = MinervaToken("transcription",TRANSCRIPTION())
    token_type_dictionary[string(hash("translation"))] = MinervaToken("translation",TRANSLATION())
    token_type_dictionary[string(hash("catalyze"))] = MinervaToken("catalyze",CATALYZE())
    token_type_dictionary[string(hash("bind"))] = MinervaToken("bind",BIND())
    token_type_dictionary[string(hash("positive"))] = MinervaToken("positive",POSITIVE())
    token_type_dictionary[string(hash("negative"))] = MinervaToken("negative",NEGATIVE())
    token_type_dictionary[string(hash("phosphorylate"))] = MinervaToken("phosphorylate",PHOSPHORYLATE())
    token_type_dictionary[string(hash("dephosphorylate"))] = MinervaToken("dephosphorylate",DEPHOSPHORYLATE())
    token_type_dictionary[string(hash("complex"))] = MinervaToken("complex",COMPLEX())
    token_type_dictionary[string(hash("form"))] = MinervaToken("form",FORM())
    token_type_dictionary[string(hash("and"))] = MinervaToken("and",AND())
    token_type_dictionary[string(hash(","))] = MinervaToken(",",AND())
    token_type_dictionary[string(hash("|"))] = MinervaToken("|",OR())
    token_type_dictionary[string(hash("or"))] = MinervaToken("or",OR())
    token_type_dictionary[string(hash(" "))] = MinervaToken(" ",SPACE())
    token_type_dictionary[string(hash("to"))] = MinervaToken("to",TO())
    token_type_dictionary[string(hash("->"))] = MinervaToken("to",TO())
    token_type_dictionary[string(hash("are"))] = MinervaToken("are",ARE())
    token_type_dictionary[string(hash("="))] = MinervaToken("=",ARE())
    token_type_dictionary[string(hash("is"))] = MinervaToken("is",IS())
    token_type_dictionary[string(hash("a"))] = MinervaToken("a",A())
    token_type_dictionary[string(hash("type"))] = MinervaToken("TYPE",TYPE())
    token_type_dictionary[string(hash("of"))] = MinervaToken("of",OF())
    token_type_dictionary[string(hash("the"))] = MinervaToken("the",THE())
    token_type_dictionary[string(hash("at"))] = MinervaToken("at",AT())
    token_type_dictionary[string(hash("in"))] = MinervaToken("in",IN())
    token_type_dictionary[string(hash("promoter"))] = MinervaToken("promoter",PROMOTER())
    token_type_dictionary[string(hash("gene"))] = MinervaToken("gene",GENE())
    token_type_dictionary[string(hash("RNA_POLYMERASE_II_SYMBOL"))] = MinervaToken("RNAP_SYMBOL",RNA_POLYMERASE_II_SYMBOL())
    token_type_dictionary[string(hash("RIBOSOME_SYMBOL"))] = MinervaToken("RIBOSOME_SYMBOL",RIBOSOME_SYMBOL())

    # token_type_dictionary[string(hash("GENE_TYPE_SYMBOL"))] = MinervaToken("GENE_SYMBOL",GENE_SYMBOL())
    # token_type_dictionary[string(hash("mRNA_TYPE_SYMBOL"))] = MinervaToken("mRNA_SYMBOL",mRNA_SYMBOL())
    # token_type_dictionary[string(hash("tRNA_TYPE_SYMBOL"))] = MinervaToken("tRNA_SYMBOL",tRNA_SYMBOL())
    # token_type_dictionary[string(hash("regRNA_TYPE_SYMBOL"))] = MinervaToken("regRNA_SYMBOL",regRNA_SYMBOL())
    # token_type_dictionary[string(hash("PROTEIN_TYPE_SYMBOL"))] = MinervaToken("PROTEIN_SYMBOL",PROTEIN_SYMBOL())
    # token_type_dictionary[string(hash("METABOLITE_TYPE_SYMBOL"))] = MinervaToken("METABOLITE_SYMBOL",METABOLITE_SYMBOL())

    try

        # First = remove all synonyms and reduce the sentence to canonical form -
        tmp_record_array = reverse(copy(record))
        while (isempty(tmp_record_array) == false)

            # check: pop the component, check against the syn data table -
            record_component = pop!(tmp_record_array)

            # do we have this component as a key?
            if (haskey(synonym_dictionary,record_component) == true)
                # we have a synonym -
                new_record_component = synonym_dictionary[record_component]
                push!(canonical_record_array, new_record_component)
            else
                push!(canonical_record_array,record_component)
            end
        end

        # ok, so now we have a canonical forms - lets tokenize the cononical elements -
        # first, lets string match -
        reverse_canonical_record_array_copy = reverse(copy(canonical_record_array))
        while (isempty(reverse_canonical_record_array_copy) == false)

            # item -
            item_to_classify = pop!(reverse_canonical_record_array_copy)

            # check: do we have this item in the token dictionary -
            test_key = string(hash(item_to_classify))
            if (haskey(token_type_dictionary, test_key) == true)
                token_item = token_type_dictionary[test_key]
                push!(canonical_token_array, token_item)
            else
                minerva_token = MinervaToken(item_to_classify, BIOLOGICAL_SYMBOL())
                push!(canonical_token_array,minerva_token)
            end
        end

        # finally - add a SEMICOLON -
        push!(canonical_token_array, MinervaToken(";", SEMICOLON()))

        # return -
        return VLResult(canonical_token_array)
    catch error
        return VLResult(error)
    end
end

max_stars_repo_name: varnerlab/CellFreeModelGenerationKit

max_issues_repo_name: varnerlab/CellFreeModelGenerationKit

max_forks_repo_name: varnerlab/CellFreeModelGenerationKit


text:
## Real-to-real transforms (requires FFTW.jl)

import FFTW: kind2string

const R2R_SUPPORTED_KINDS = (
    FFTW.DHT,
    FFTW.REDFT00,
    FFTW.REDFT01,
    FFTW.REDFT10,
    FFTW.REDFT11,
    FFTW.RODFT00,
    FFTW.RODFT01,
    FFTW.RODFT10,
    FFTW.RODFT11,
)

"""
    R2R(kind)

Real-to-real transform of type `kind`.

The possible values of `kind` are those described in the
[`FFTW.r2r`](https://juliamath.github.io/FFTW.jl/stable/fft/#FFTW.r2r)
docs and the [`FFTW`](http://www.fftw.org/doc/) manual:

- [discrete cosine transforms](http://www.fftw.org/doc/1d-Real_002deven-DFTs-_0028DCTs_0029.html#g_t1d-Real_002deven-DFTs-_0028DCTs_0029):
  `FFTW.REDFT00`, `FFTW.REDFT01`, `FFTW.REDFFT10`, `FFTW.REDFFT11`

- [discrete sine transforms](http://www.fftw.org/doc/1d-Real_002dodd-DFTs-_0028DSTs_0029.html#g_t1d-Real_002dodd-DFTs-_0028DSTs_0029):
  `FFTW.RODFT00`, `FFTW.RODFT01`, `FFTW.RODFFT10`, `FFTW.RODFFT11`

- [discrete Hartley transform](http://www.fftw.org/doc/1d-Discrete-Hartley-Transforms-_0028DHTs_0029.html#g_t1d-Discrete-Hartley-Transforms-_0028DHTs_0029):
  `FFTW.DHT`

Note: [half-complex format
DFTs](http://www.fftw.org/doc/The-Halfcomplex_002dformat-DFT.html#The-Halfcomplex_002dformat-DFT)
(`FFTW.R2HC`, `FFTW.HC2R`) are not currently supported.

"""
struct R2R{kind} <: AbstractTransform
    function R2R{kind}() where kind
        if kind ∉ R2R_SUPPORTED_KINDS
            throw(ArgumentError(
                "unsupported r2r transform kind: $(kind2string(kind))"))
        end
        new()
    end
end

"""
    R2R!(kind)

In-place version of [`R2R`](@ref).

See also [`FFTW.r2r!`](https://juliamath.github.io/FFTW.jl/stable/fft/#FFTW.r2r!).
"""
struct R2R!{kind} <: AbstractTransform
    function R2R!{kind}() where kind
        R2R{kind}()  # executes verification code above
        new()
    end
end

@inline R2R(kind) = R2R{kind}()
@inline R2R!(kind) = R2R!{kind}()

const AnyR2R{kind} = Union{R2R{kind}, R2R!{kind}} where {kind}

is_inplace(::R2R) = false
is_inplace(::R2R!) = true

Base.show(io::IO, tr::R2R) = print(io, "R2R{", kind2string(kind(tr)), "}")
Base.show(io::IO, tr::R2R!) = print(io, "R2R!{", kind2string(kind(tr)), "}")

"""
    kind(transform::R2R)

Get `kind` of real-to-real transform.
"""
kind(::AnyR2R{K}) where {K} = K

length_output(::AnyR2R, length_in::Integer) = length_in
eltype_input(::AnyR2R, ::Type) = Nothing  # both real and complex inputs are accepted
eltype_output(::AnyR2R, ::Type{T}) where {T} = T

# NOTE: plan_r2r is type-unstable!!
# More precisely, plan_r2r returns a FFTW.r2rFFTWPlan{T, K, inplace, N},
# whose second parameter `K` is "a tuple of the transform kinds along each
# dimension". That is, `K` is a tuple of the form `(kind1, kind2, ..., kindN)`
# with the same length as `dims`.
#
# Since we have static information on the kind (it's a parameter of the R2R
# type), we try to work around the issue by typing the return type. This
# function will be type-stable if `dims` has a static length, i.e. if
# `length(dims)` is known by the compiler. This will be the case if `dims` is a
# tuple or a scalar value (e.g. `(1, 3)` or `1`), but not if it is a range (e.g.
# `2:3`).
function plan(transform::AnyR2R, A, dims; kwargs...)
    kd = kind(transform)
    K = ntuple(_ -> kd, length(dims))
    R = FFTW.r2rFFTWPlan{T,K} where {T}  # try to guess the return type
    _plan_r2r(transform, A, kd, dims; kwargs...) :: R
end

_plan_r2r(::R2R, args...; kwargs...) = FFTW.plan_r2r(args...; kwargs...)
_plan_r2r(::R2R!, args...; kwargs...) = FFTW.plan_r2r!(args...; kwargs...)

# Scale factors for r2r transforms.
scale_factor(::AnyR2R, A, dims) = _prod_dims(2 .* size(A), dims)
scale_factor(::AnyR2R{FFTW.REDFT00}, A, dims) =
    _prod_dims(2 .* (size(A) .- 1), dims)
scale_factor(::AnyR2R{FFTW.RODFT00}, A, dims) =
    _prod_dims(2 .* (size(A) .+ 1), dims)
scale_factor(::AnyR2R{FFTW.DHT}, A, dims) = _prod_dims(A, dims)

for T in (:R2R, :R2R!)
    @eval begin
        # From FFTW docs (4.8.3 1d Real-even DFTs (DCTs)):
        #   The unnormalized inverse of REDFT00 is REDFT00, of REDFT10 is REDFT01 and
        #   vice versa, and of REDFT11 is REDFT11.
        #   Each unnormalized inverse results in the original array multiplied by N,
        #   where N is the logical DFT size. For REDFT00, N=2(n-1) (note that n=1 is not
        #   defined); otherwise, N=2n.
        binv(::$T{FFTW.REDFT00}) = $T{FFTW.REDFT00}()
        binv(::$T{FFTW.REDFT01}) = $T{FFTW.REDFT10}()
        binv(::$T{FFTW.REDFT10}) = $T{FFTW.REDFT01}()
        binv(::$T{FFTW.REDFT11}) = $T{FFTW.REDFT11}()

        # From FFTW docs (4.8.4 1d Real-odd DFTs (DSTs)):
        #    The unnormalized inverse of RODFT00 is RODFT00, of RODFT10 is RODFT01 and
        #    vice versa, and of RODFT11 is RODFT11.
        #    Each unnormalized inverse results in the original array multiplied by N,
        #    where N is the logical DFT size. For RODFT00, N=2(n+1); otherwise, N=2n.
        binv(::$T{FFTW.RODFT00}) = $T{FFTW.RODFT00}()
        binv(::$T{FFTW.RODFT01}) = $T{FFTW.RODFT10}()
        binv(::$T{FFTW.RODFT10}) = $T{FFTW.RODFT01}()
        binv(::$T{FFTW.RODFT11}) = $T{FFTW.RODFT11}()

        # From FFTW docs (4.8.5 1d Discrete Hartley Transforms (DHTs)):
        #    [...] applying the transform twice (the DHT is its own inverse) will
        #    multiply the input by n.
        binv(::$T{FFTW.DHT}) = $T{FFTW.DHT}()
    end
end

expand_dims(::F, ::Val{N}) where {F <: AnyR2R, N} =
    N === 0 ? () : (F(), expand_dims(F(), Val(N - 1))...)

ext: jl


text:
abstract type PowerNetworkMatrix{T} <: AbstractArray{T, 2} end

#  The container code for PowerNetworkMatrix is based in JuMP's Container in order to
#  remove the limitations of AxisArrays and the doubts about long term maintenance
#  https://github.com/JuliaOpt/JuMP.jl/blob/master/src/Containers/DenseAxisArray.jl
#  JuMP'sCopyright 2017, Iain Dunning, Joey Huchette, Miles Lubin, and contributors
#  This Source Code Form is subject to the terms of the Mozilla Public
#  License, v. 2.0.

function get_bus_indices(branch, bus_lookup)
    fr_b = bus_lookup[get_number(get_from(get_arc(branch)))]
    to_b = bus_lookup[get_number(get_to(get_arc(branch)))]
    return fr_b, to_b
end

function _make_ax_ref(buses::Vector{Bus})
    return _make_ax_ref(get_number.(buses))
end

function _make_ax_ref(ax::Vector)
    ref = Dict{eltype(ax), Int}()
    for (ix, el) in enumerate(ax)
        if haskey(ref, el)
            @error("Repeated index element $el. Index sets must have unique elements.")
        end
        ref[el] = ix
    end
    return ref
end

# AbstractArray interface
Base.isempty(A::PowerNetworkMatrix) = isempty(A.data)
Base.size(A::PowerNetworkMatrix) = size(A.data)
Base.LinearIndices(A::PowerNetworkMatrix) =
    error("PowerSystems.$(typeof(A)) does not support this operation.")
Base.axes(A::PowerNetworkMatrix) = A.axes
Base.CartesianIndices(A::PowerNetworkMatrix) =
    error("PowerSystems.$(typeof(A)) does not support this operation.")

############
# Indexing #
############

Base.eachindex(A::PowerNetworkMatrix) = CartesianIndices(size(A.data))

function lookup_index(i, lookup::Dict)
    return isa(i, Colon) ? Colon() : lookup[i]
end

function lookup_index(i::Component, lookup::Dict)
    return isa(i, Colon) ? Colon() : lookup[Base.to_index(i)]
end

# Lisp-y tuple recursion trick to handle indexing in a nice type-
# stable way. The idea here is that `_to_index_tuple(idx, lookup)`
# performs a lookup on the first element of `idx` and `lookup`,
# then recurses using the remaining elements of both tuples.
# The compiler knows the lengths and types of each tuple, so
# all of the types are inferable.
function _to_index_tuple(idx::Tuple, lookup::Tuple)
    tuple(
        lookup_index(first(idx), first(lookup)),
        _to_index_tuple(Base.tail(idx), Base.tail(lookup))...,
    )
end

# Handle the base case when we have more indices than lookups:
function _to_index_tuple(idx::NTuple{N}, ::NTuple{0}) where {N}
    ntuple(k -> begin
        i = idx[k]
        (i == 1) ? 1 : error("invalid index $i")
    end, Val(N))
end

# Handle the base case when we have fewer indices than lookups:
_to_index_tuple(idx::NTuple{0}, lookup::Tuple) = ()

# Resolve ambiguity with the above two base cases
_to_index_tuple(idx::NTuple{0}, lookup::NTuple{0}) = ()
to_index(A::PowerNetworkMatrix, idx...) = _to_index_tuple(idx, A.lookup)

# Doing `Colon() in idx` is relatively slow because it involves
# a non-unrolled loop through the `idx` tuple which may be of
# varying element type. Another lisp-y recursion trick fixes that
has_colon(idx::Tuple{}) = false
has_colon(idx::Tuple) = isa(first(idx), Colon) || has_colon(Base.tail(idx))

# TODO: better error (or just handle correctly) when user tries to index with a range like a:b
# The only kind of slicing we support is dropping a dimension with colons
function Base.getindex(A::PowerNetworkMatrix, idx...)
    #=
    if has_colon(idx)
        PTDF(A.data[to_index(A,idx...)...], (ax for (i,ax) in enumerate(A.axes) if idx[i] == Colon())...)
    else
        return A.data[to_index(A,idx...)...]
    end
    =#
    return A.data[to_index(A, idx...)...]
end
Base.getindex(A::PowerNetworkMatrix, idx::CartesianIndex) = A.data[idx]
Base.setindex!(A::PowerNetworkMatrix, v, idx...) = A.data[to_index(A, idx...)...] = v
Base.setindex!(A::PowerNetworkMatrix, v, idx::CartesianIndex) = A.data[idx] = v

Base.IndexStyle(::Type{PowerNetworkMatrix}) = IndexAnyCartesian()

########
# Keys #
########

struct PowerNetworkMatrixKey{T <: Tuple}
    I::T
end
Base.getindex(k::PowerNetworkMatrixKey, args...) = getindex(k.I, args...)

struct PowerNetworkMatrixKeys{T <: Tuple}
    product_iter::Base.Iterators.ProductIterator{T}
end
Base.length(iter::PowerNetworkMatrixKeys) = length(iter.product_iter)
function Base.eltype(iter::PowerNetworkMatrixKeys)
    return PowerNetworkMatrixKey{eltype(iter.product_iter)}
end
function Base.iterate(iter::PowerNetworkMatrixKeys)
    next = iterate(iter.product_iter)
    return next === nothing ? nothing : (PowerNetworkMatrixKey(next[1]), next[2])
end
function Base.iterate(iter::PowerNetworkMatrixKeys, state)
    next = iterate(iter.product_iter, state)
    return next === nothing ? nothing : (PowerNetworkMatrixKey(next[1]), next[2])
end
function Base.keys(a::PowerNetworkMatrix)
    return PowerNetworkMatrixKeys(Base.Iterators.product(a.axes...))
end
Base.getindex(a::PowerNetworkMatrix, k::PowerNetworkMatrixKey) = a[k.I...]

########
# Show #
########

# Adapted printing from JuMP's implementation of the Julia's show.jl
# used in PowerNetworkMatrixs

# Copyright (c) 2009-2016: Jeff Bezanson, Stefan Karpinski, Viral B. Shah,
# and other contributors:
#
# https://github.com/JuliaLang/julia/contributors
#
# Permission is hereby granted, free of charge, to any person obtaining
# a copy of this software and associated documentation files (the
# "Software"), to deal in the Software without restriction, including
# without limitation the rights to use, copy, modify, merge, publish,
# distribute, sublicense, and/or sell copies of the Software, and to
# permit persons to whom the Software is furnished to do so, subject to
# the following conditions:
#
# The above copyright notice and this permission notice shall be
# included in all copies or substantial portions of the Software.

function Base.summary(io::IO, A::PowerNetworkMatrix)
    _summary(io, A)
    for (k, ax) in enumerate(A.axes)
        print(io, "    Dimension $k, ")
        show(IOContext(io, :limit => true), ax)
        println(io)
    end
    print(io, "And data, a ", size(A.data))
end
_summary(io::IO, A::PowerNetworkMatrix) = println(io, "PowerNetworkMatrix")

function Base.summary(
    io::IOContext{Base.GenericIOBuffer{Array{UInt8, 1}}},
    ::PowerNetworkMatrix,
)
    println(io, "PowerNetworkMatrix")
end

function Base.summary(A::PowerNetworkMatrix)
    io = IOBuffer()
    Base.summary(io, A)
    String(take!(io))
end

if isdefined(Base, :print_array) # 0.7 and later
    Base.print_array(io::IO, X::PowerNetworkMatrix) = Base.print_matrix(io, X.data)
end

# n-dimensional arrays
function Base.show_nd(
    io::IO,
    a::PowerNetworkMatrix,
    print_matrix::Function,
    label_slices::Bool,
)
    limit::Bool = get(io, :limit, false)
    if isempty(a)
        return
    end
    tailinds = Base.tail(Base.tail(axes(a.data)))
    nd = ndims(a) - 2
    for I in CartesianIndices(tailinds)
        idxs = I.I
        if limit
            for i in 1:nd
                ii = idxs[i]
                ind = tailinds[i]
                if length(ind) > 10
                    if ii == ind[4] && all(d -> idxs[d] == first(tailinds[d]), 1:(i - 1))
                        for j in (i + 1):nd
                            szj = size(a.data, j + 2)
                            indj = tailinds[j]
                            if szj > 10 && first(indj) + 2 < idxs[j] <= last(indj) - 3
                                @goto skip
                            end
                        end
                        #println(io, idxs)
                        print(io, "...\n\n")
                        @goto skip
                    end
                    if ind[3] < ii <= ind[end - 3]
                        @goto skip
                    end
                end
            end
        end
        if label_slices
            print(io, "[:, :, ")
            for i in 1:(nd - 1)
                show(io, a.axes[i + 2][idxs[i]])
                print(io, ", ")
            end
            show(io, a.axes[end][idxs[end]])
            println(io, "] =")
        end
        slice = view(a.data, axes(a.data, 1), axes(a.data, 2), idxs...)
        Base.print_matrix(io, slice)
        print(io, idxs == map(last, tailinds) ? "" : "\n\n")
        @label skip
    end
end

function Base.show(io::IO, array::PowerNetworkMatrix)
    summary(io, array)
    isempty(array) && return
    println(io, ":")
    Base.print_array(io, array)
end

Base.to_index(b::Bus) = get_number(b)
Base.to_index(b::T) where {T <: ACBranch} = get_name(b)
Base.to_index(ix::Component...) = to_index.(ix)

"""returns the raw array data of the `PowerNetworkMatrix`"""
get_data(mat::PowerNetworkMatrix) = mat.data

"""
    returns the lookup tuple of the `PowerNetworkMatrix`. The first entry corresponds
    to the first dimension and the second entry corresponds to the second dimension. For
    instance in Ybus the first dimension is buses and second dimension is buses too, and in
    PTDF the first dimension is branches and the second dimension is buses
"""
get_lookup(mat::PowerNetworkMatrix) = mat.lookup


text:
# Undergraduate Student: Arturo Burgos
# Professor: João Rodrigo Andrade
# Federal University of Uberlândia - UFU, Fluid Mechanics Laboratory - MFLab, Block 5P, Uberlândia, MG, Brazil

using Printf # Just the same as import a module
using Statistics
using CPUTime

# A = rand(1:50,3,3)
# B = rand(1:50,3,3)

A = [1 2; 3 4]
println("\n\n\nThe A Matrix is: \n\n",A)

B = [1 2; 3 4]
println("\n\nThe B Matrix is: \n\n",B)

R = zeros(2,2) # Here I create a zero matrix just to get the multiply of A and B

D = A*B # This is the vectorized way for this problem, like MATLAB
C = A.*B # This multiplies point by point

println("\n\nThe D Matrix is: \n\n",D)
println("\n\nThe C Matrix is: \n\n",C)

# println("\n\n",D[:,1])

function multivec(mtx1,mtx2,mtx3)

    mtx3 = mtx1*mtx2

    return mtx3 # Returns a value to the variable
end

R = multivec(A,B,R);
println("\n\nThe R Matrix is: \n\n",R)
println("\n\n")

@time @CPUtime multivec(A,B,R)

max_stars_repo_name: bangyen/Comparing-Languages

max_issues_repo_name: bangyen/Comparing-Languages

max_forks_repo_name: bangyen/Comparing-Languages


text:
module Shapes

import Adapt

import Base: getindex, length, size
using Base: @pure, @propagate_inbounds, @_inline_meta

import Random
using Random: AbstractRNG

using Requires: @init, @require

import UnsafeArrays

import StaticArrays: Size, Length, similar_type, get
using StaticArrays: tuple_prod,
                    tuple_length,
                    tuple_minimum,
                    promote_tuple_eltype,
                    size_to_tuple,
                    SArray,
                    MArray,
                    StaticArray,
                    SOneTo

export AbstractShape,
       checkaxes,
       checksize,
       Shape,
       Size,
       Length,
       ScalarShape,
       VectorShape,
       MatrixShape,
       MultiShape,
       concrete_eltype,

       ShapedView,
       allocate

include("util.jl")
include("traits.jl")
include("ShapedView.jl")

getdims(shape::AbstractShape, dims) = getdims(Size(shape), dims)
getdims(s::Size, dims) = ((Tuple(s)..., dims...))
getdims(::Size{()}, dims) = dims

Base.Array{T}(::UndefInitializer, shape::AbstractShape, dims::Dims) where {T} =
    Array{T}(undef, getdims(shape, dims))
Base.Array(::UndefInitializer, shape::AbstractShape, dims::Dims) =
    Array{concrete_eltype(shape)}(undef, shape, dims)

Base.Array{T}(::UndefInitializer, shape::AbstractShape, dims::Integer...) where {T} =
    Array{T}(undef, shape, dims)
Base.Array(::UndefInitializer, shape::AbstractShape, dims::Integer...) =
    Array(undef, shape, dims)

Base.Matrix{T}(::UndefInitializer, shape::AbstractShape, dim::Integer) where {T} =
    Matrix{T}(undef, getdims(shape, (dim,)))
Base.Matrix(::UndefInitializer, shape::AbstractShape, dim::Integer) =
    Matrix{concrete_eltype(shape)}(undef, shape, dim)

Base.Vector{T}(::UndefInitializer, shape::AbstractShape) where {T} =
    Vector{T}(undef, length(shape))
Base.Vector(::UndefInitializer, shape::AbstractShape) =
    Vector{concrete_eltype(shape)}(undef, shape)

Base.zeros(T, shape::AbstractShape, dims::Dims) = zeros(T, getdims(shape, dims))
Base.zeros(shape::AbstractShape, dims::Dims) = zeros(concrete_eltype(shape), shape, dims)

Base.zeros(T, shape::AbstractShape, dims::Integer...) = zeros(T, shape, dims)
Base.zeros(shape::AbstractShape, dims::Integer...) = zeros(shape, dims)

allocate(shape::AbstractShape, dims::Dims) = Array(undef, shape, dims)
allocate(shape::AbstractShape, dims::Integer...) = allocate(shape, dims)

allocate(T::DataType, dims::Dims) = Array{T}(undef, dims)
allocate(T::DataType, dims::Integer...) = allocate(T, dims)

allocate(T::DataType, shape::AbstractShape, dims::Dims) = allocate(T, size(shape)..., dims...)
allocate(T::DataType, shape::AbstractShape, dims::Integer...) = allocate(T, shape, dims)

function Random.rand(rng::AbstractRNG, shape::AbstractShape, dims::Dims)
    _rand(rng, concrete_eltype(shape), getdims(shape, dims))
end
Random.rand(rng::AbstractRNG, shape::AbstractShape, dims::Integer...) = rand(rng, shape, dims)

Random.rand(shape::AbstractShape, dims::Dims) = _rand(concrete_eltype(shape), getdims(shape, dims))
Random.rand(shape::AbstractShape, dims::Integer...) = rand(shape, dims)

_rand(rng, T, ::Dims{0}) = rand(rng, T)
_rand(rng, T, dims::Dims) = rand(rng, T, dims)
_rand(T, ::Dims{0}) = rand(T)
_rand(T, dims::Dims) = rand(T, dims)

# StaticArrays support
similar_type(::Type{SH}) where {SH<:AbstractShape} =
    SArray{tuple_to_size(Tuple(Size(SH))),concrete_eltype(SH),ndims(SH),length(SH)}
similar_type(::SH) where {SH<:AbstractShape} = similar_type(SH)

similar_type(::Type{SA}, ::Type{SH}) where {SA<:StaticArray,SH<:AbstractShape} =
    similar_type(SA, concrete_eltype(SH), Size(SH))
similar_type(::SA, ::Type{SH}) where {SA<:StaticArray,SH<:AbstractShape} =
    similar_type(SA, concrete_eltype(SH), Size(SH))
similar_type(::Type{SA}, ::SH) where {SA<:StaticArray,SH<:AbstractShape} =
    similar_type(SA, concrete_eltype(SH), Size(SH))
similar_type(::SA, ::SH) where {SA<:StaticArray,SH<:AbstractShape} =
    similar_type(SA, concrete_eltype(SH), Size(SH))

@init @require ElasticArrays="fdbdab4c-e67f-52f5-8c3f-e7b388dad3d4" begin
    import .ElasticArrays
    ElasticArrays.ElasticArray(::UndefInitializer, shape::AbstractShape, dims::Dims) = ElasticArrays.ElasticArray{eltype(shape)}(undef, getdims(shape, dims)...)
    ElasticArrays.ElasticArray(::UndefInitializer, shape::AbstractShape, dims::Integer...) = ElasticArrays.ElasticArray(undef, shape, dims)
end

end # module


text:
module CloudGraphs_RegisteredTypes

include("RegisterTypes.jl")

end


text:
# This file is a part of Julia. License is MIT: https://julialang.org/license

using Test, Base.BinaryPlatforms, Base.BinaryPlatforms.CPUID

@testset "CPUID" begin
    @test CPUID.cpu_isa() isa CPUID.ISA

    get_x86_64(n) = (CPUID.ISAs_by_family["x86_64"][n].second)
    @test get_x86_64(2) <  get_x86_64(4)
    @test get_x86_64(5) <= get_x86_64(5)
    @test get_x86_64(3) >= get_x86_64(3)
    @test get_x86_64(7) >= get_x86_64(1)
    @test sort([get_x86_64(6), get_x86_64(4), get_x86_64(2), get_x86_64(4)]) ==
        [get_x86_64(2), get_x86_64(4), get_x86_64(4), get_x86_64(6)]
end

# Helper constructor to create a Platform with `validate_strict` set to `true`.
P(args...; kwargs...) = Platform(args...; validate_strict=true, kwargs...)

# Ensure the platform type constructors are well behaved
@testset "Platform validation" begin
    for os_name in ("Linux", "linux", "FREEBSD", "wiNDows", "MacOS")
        p = P("x86_64", os_name)
        @test isa(p, Platform)
        @test os(p) == lowercase(os_name)
    end
    for arch_name in ("x86_64", "ARMv7l", "armv6l", "I686", "pOWerpc64le")
        p = P(arch_name, "linux")
        @test isa(p, Platform)
        @test arch(p) == lowercase(arch_name)
    end
    # Test some normalization
    @test arch(P("amd64", "freebsd")) == "x86_64"
    @test arch(P("ARM", "linux")) == "armv7l"
    @test arch(P("ppc64le", "linux")) == "powerpc64le"

    # test some error cases
    @test_throws ArgumentError P("x64_86", "linux")
    @test_throws ArgumentError P("x86_64", "not_an_os")
    @test_throws ArgumentError P("x86_64", "linux"; libc="crazy_libc")
    @test_throws ArgumentError P("x86_64", "linux"; libc="glibc", call_abi="crazy_abi")
    @test_throws ArgumentError P("x86_64", "linux"; libc="glibc", call_abi="eabihf")
    @test_throws ArgumentError P("i686", "linux"; libc="musl", call_abi="eabi")
    @test_throws ArgumentError P("arm", "linux"; call_abi="")
    @test_throws ArgumentError P("armv7l", "linux"; call_abi="kekeke")
    @test_throws ArgumentError P("armv6l", "linux"; call_abi="kekeke")
    @test_throws ArgumentError P("armv6l", "linux"; libgfortran_version="lel")
    @test_throws ArgumentError P("x86_64", "linux"; cxxstring_abi="lel")
    @test_throws ArgumentError P("x86_64", "windows"; libstdcxx_version="lel")
    @test_throws ArgumentError P("i686", "macos")
    @test_throws ArgumentError P("x86_64", "macos"; libc="glibc")
    @test_throws ArgumentError P("x86_64", "macos"; call_abi="eabihf")
    @test_throws ArgumentError P("powerpc64le", "windows")
    @test_throws ArgumentError P("x86_64", "windows"; call_abi="eabihf")
    @test_throws ArgumentError P("x86_64", "freebsd"; libc="crazy_libc")
    @test_throws ArgumentError P("x86_64", "freebsd"; call_abi="crazy_abi")
    @test_throws ArgumentError P("x86_64", "freebsd"; call_abi="eabihf")
    @test_throws ArgumentError P("x86_64", "linux"; arch="i686")
    @test_throws ArgumentError P("x86_64", "linux"; ARCH="i686")
    @test_throws ArgumentError P("x86_64", "linux"; os="windows")
end

@testset "Platform properties" begin
    # Test that `platform_name()` works
    platforms = ("Linux", "macOS", "Windows", "FreeBSD")
    for platform in platforms
        @test platform_name(P("x86_64", platform)) == platform
    end

    # Test `arch()`
    arch_names = ("x86_64", "i686", "powerpc64le", "armv7l", "armv6l", "aarch64")
    for arch_name in arch_names
        @test arch(P(arch_name, "linux")) == arch_name
    end

    # Test that if we aren't using strict validation, we can actually use new names too:
    @test arch(Platform("jpu", "linux")) == "jpu"

    # platform_dlext()
    @test platform_dlext(P("x86_64", "linux")) == "so"
    @test platform_dlext(P("armv7l", "windows")) == "dll"
    @test platform_dlext(P("x86_64", "freebsd")) == "so"
    @test platform_dlext(P("aarch64", "macos")) == "dylib"
    @test platform_dlext() == platform_dlext(HostPlatform())

    # wordsize()
    @test wordsize(P("i686", "linux")) == wordsize(P("armv7l", "windows")) == 32
    @test wordsize(P("aarch64", "macos")) == wordsize(P("x86_64", "freebsd")) == 64
    @test wordsize(P("x86_64", "windows")) == wordsize(P("powerpc64le", "linux")) == 64

    # call_abi()
    for platform in platforms
        @test call_abi(P("x86_64", platform)) === nothing
    end
    @test call_abi(P("armv7l", "linux")) == "eabihf"
    @test call_abi(P("armv7l", "linux"; call_abi="eabihf")) == "eabihf"
    @test call_abi(P("armv7l", "linux"; call_abi="eabi")) == "eabi"
    @test call_abi(P("armv6l", "linux")) == "eabihf"
    # Test that we can at least set an `eabi` call ABI, not that Julia actually supports it...
    @test call_abi(P("armv7l", "linux"; call_abi="eabi")) == "eabi"

    # Test some different OS's and libc/call ABIs
    @test triplet(P("i686", "Windows")) == "i686-w64-mingw32"
    @test triplet(P("x86_64", "linux"; libc="musl")) == "x86_64-linux-musl"
    @test triplet(P("armv7l", "linux"; libc="musl")) == "armv7l-linux-musleabihf"
    @test triplet(P("armv6l", "linux"; libc="musl", call_abi="eabihf")) == "armv6l-linux-musleabihf"
    @test triplet(P("armv6l", "linux"; call_abi="eabi")) == "armv6l-linux-gnueabi"
    @test triplet(P("x86_64", "linux")) == "x86_64-linux-gnu"
    @test triplet(P("armv6l", "linux")) == "armv6l-linux-gnueabihf"
    @test triplet(P("x86_64", "macos")) == "x86_64-apple-darwin"
    @test triplet(P("x86_64", "macos"; os_version=v"16")) == "x86_64-apple-darwin16"
    @test triplet(P("x86_64", "freebsd")) == "x86_64-unknown-freebsd"
    @test triplet(P("i686", "freebsd")) == "i686-unknown-freebsd"

    # Now test libgfortran/cxxstring ABIs
    @test triplet(P("x86_64", "linux"; libgfortran_version=v"3", cxxstring_abi="cxx11")) == "x86_64-linux-gnu-libgfortran3-cxx11"
    @test triplet(P("armv7l", "linux"; libc="musl", cxxstring_abi="cxx03")) == "armv7l-linux-musleabihf-cxx03"
    if !isnothing(detect_libgfortran_version())
        # When `libgfortran` can be detected at runtime, make sure
        # `HostPlatform` has the appropriate key.
        @test tags(HostPlatform())["libgfortran_version"] == string(detect_libgfortran_version())
    end

    # Test tags()
    t = tags(P("x86_64", "linux"))
    @test all(haskey.(Ref(t), ("arch", "os", "libc")))
    @test haskey(tags(P("x86_64", "linux"; customtag="foo")), "customtag")
    @test tags(HostPlatform())["julia_version"] == string(VERSION.major, ".", VERSION.minor, ".", VERSION.patch)

    # Test that we can modify tags at will using the dict-like interface:
    p = P("x86_64", "linux")
    p["foo"] = "bar"
    @test tags(p)["foo"] == "bar"
    @test p["foo"] == "bar"
    @test p["os"] == "linux"
    p["os"] = "JuliaOS"
    @test p["os"] == "juliaos"

    # Test that trying to set illegal tags fails
    @test_throws ArgumentError p["os"] = "a+b"

    # Test that our `hash()` is stable
    @test hash(HostPlatform()) == hash(HostPlatform())

    # Test that round-tripping through `triplet` for a does not
    # maintain equality, as we end up losing the `compare_strategies`:
    p = Platform("x86_64", "linux"; cuda = v"11")
    Base.BinaryPlatforms.set_compare_strategy!(p, "cuda", Base.BinaryPlatforms.compare_version_cap)
    q = parse(Platform, triplet(p))
    @test q != p
end

@testset "Triplet parsing" begin
    # Make sure the Platform() with explicit triplet works
    R(str) = parse(Platform, str; validate_strict=true)
    @test R("x86_64-linux-gnu") == P("x86_64", "linux")
    @test R("x86_64-linux-musl") == P("x86_64", "linux"; libc="musl")
    @test R("i686-unknown-linux-gnu") == P("i686", "linux")
    @test R("x86_64-apple-darwin") == P("x86_64", "macos")
    @test R("x86_64-apple-darwin14") == P("x86_64", "macos"; os_version="14")
    @test R("x86_64-apple-darwin17.0.0") == P("x86_64", "macos"; os_version="17")
    @test R("armv7l-pc-linux-gnueabihf") == P("armv7l", "linux")
    @test R("armv7l-linux-musleabihf") == P("armv7l", "linux"; libc="musl")
    @test R("armv6l-linux-gnueabi") == P("armv6l", "linux"; call_abi="eabi")
    # Test that the short name "arm" goes to `armv7l`
    @test R("arm-linux-gnueabihf") == P("armv7l", "linux")
    @test R("aarch64-unknown-linux-gnu") == P("aarch64", "linux")
    @test R("powerpc64le-linux-gnu") == P("powerpc64le", "linux")
    @test R("ppc64le-linux-gnu") == P("powerpc64le", "linux")
    @test R("x86_64-w64-mingw32") == P("x86_64", "windows")
    @test R("i686-w64-mingw32") == P("i686", "windows")

    # FreeBSD has lots of arch names that don't match elsewhere
    @test R("x86_64-unknown-freebsd11.1") == P("x86_64", "freebsd"; os_version=v"11.1")
    @test R("i686-unknown-freebsd11.1") == P("i686", "freebsd"; os_version=v"11.1")
    @test R("amd64-unknown-freebsd12.0") == P("x86_64", "freebsd"; os_version=v"12.0")
    @test R("i386-unknown-freebsd10.3") == P("i686", "freebsd"; os_version=v"10.3")
    @test R("aarch64-apple-darwin18.7") == P("aarch64", "macos"; os_version=v"18.7")
    @test R("arm64-apple-darwin20") == P("aarch64", "macos"; os_version=v"20")

    # Test inclusion of ABI stuff, both old-style and new-style
    @test R("x86_64-linux-gnu-gcc7") == P("x86_64", "linux"; libgfortran_version=v"4")
    @test R("x86_64-linux-gnu-gcc4-cxx11") == P("x86_64", "linux"; libgfortran_version=v"3", cxxstring_abi="cxx11")
    @test R("x86_64-linux-gnu-cxx11") == P("x86_64", "linux"; cxxstring_abi="cxx11")
    @test R("x86_64-linux-gnu-libgfortran3-cxx03") == P("x86_64", "linux"; libgfortran_version=v"3", cxxstring_abi="cxx03")
    @test R("x86_64-linux-gnu-libstdcxx26") ==  P("x86_64", "linux"; libstdcxx_version=v"3.4.26")

    @test_throws ArgumentError R("totally FREEFORM text!!1!!!1!")
    @test_throws ArgumentError R("invalid-triplet-here")
    @test parse(Platform, "aarch64-linux-gnueabihf") == Platform("aarch64", "linux"; call_abi="eabihf")
    @test_throws ArgumentError R("aarch64-linux-gnueabihf")
    @test_throws ArgumentError R("x86_64-w32-mingw64")

    # Test extended attributes
    @test R("x86_64-linux-gnu-march+avx2") == P("x86_64", "linux"; march="avx2")
    @test R("x86_64-linux-gnu-march+x86_64-cuda+10.1") == P("x86_64", "linux"; march="x86_64", cuda="10.1")

    # Round-trip our little homie through `triplet()`, with some bending
    # of the rules for MacOS and FreeBSD, who have incomplete `os_version`
    # numbers embedded within their triplets.
    p = Platform("x86_64", "linux")
    @test parse(Platform, triplet(p)) == p

    # Also test round-tripping through `repr()`:
    p = Platform("aarch64", "macos"; os_version=v"20", march="armv8_4_crypto_sve")
    @test eval(Meta.parse(repr(p))) == p
end

@testset "platforms_match()" begin
    # Just do a quick combinatorial sweep for completeness' sake for platform matching
    linux = P("x86_64", "linux")
    for libgfortran_version in (nothing, v"3", v"5"),
        libstdcxx_version in (nothing, v"3.4.18", v"3.4.26"),
        cxxstring_abi in (nothing, :cxx03, :cxx11)

        p = P("x86_64", "linux"; libgfortran_version, libstdcxx_version, cxxstring_abi)
        @test platforms_match(linux, p)
        @test platforms_match(p, linux)

        # Also test auto-string-parsing
        @test platforms_match(triplet(linux), p)
        @test platforms_match(linux, triplet(p))
    end

    # Test that Julia version is matched only on major.minor by default
    @test platforms_match(P("x86_64", "linux"; julia_version=v"1.5.0"),
                          P("x86_64", "linux"; julia_version=v"1.5.1"))
    @test !platforms_match(P("x86_64", "linux"; julia_version=v"1.5.0"),
                           P("x86_64", "linux"; julia_version=v"1.6.0"))

    # Ensure many of these things do NOT match
    @test !platforms_match(linux, P("i686", "linux"))
    @test !platforms_match(linux, P("x86_64", "windows"))
    @test !platforms_match(linux, P("x86_64", "macos"))

    # Make some explicitly non-matching compiler ABI platforms
    host = P("x86_64", "linux"; libgfortran_version=v"5", cxxstring_abi="cxx11")
    for arch in ("x86_64", "i686", "aarch64", "armv6l", "armv7l", "powerpc64le"),
        kwargs in ((:libgfortran_version => v"3",), (:cxxstring_abi => "cxx03",),
                   (:libgfortran_version => v"4", :cxxstring_abi => "cxx11"),
                   (:libgfortran_version => v"3", :cxxstring_abi => "cxx03"))
        a = P(arch, "linux"; libgfortran_version=v"5", cxxstring_abi="cxx11")
        b = P(arch, "linux"; kwargs...)
        @test !platforms_match(a, b)
    end

    # Test version bounds with HostPlatform()
    host = HostPlatform(P("x86_64", "macos"; os_version="14", libstdcxx_version=v"3.4.26"))
    @test platforms_match(host, P("x86_64", "macos"))
    @test platforms_match(host, P("x86_64", "macos"; os_version="14"))
    @test platforms_match(host, P("x86_64", "macos"; os_version="13"))
    @test !platforms_match(host, P("x86_64", "macos"; os_version="15"))
    @test platforms_match(host, P("x86_64", "macos"; libstdcxx_version="3.4.18"))
    @test platforms_match(host, P("x86_64", "macos"; os_version=v"10", libstdcxx_version="3.4.18"))
    @test !platforms_match(host, P("x86_64", "macos"; os_version=v"10", libstdcxx_version="3.4.27"))
    @test !platforms_match(host, P("x86_64", "macos"; os_version=v"14", libstdcxx_version=v"4"))
end

@testset "DL name/version parsing" begin
    # Make sure our version parsing code is working
    @test parse_dl_name_version("libgfortran.so", "linux") == ("libgfortran", nothing)
    @test parse_dl_name_version("libgfortran.so.3", "linux") == ("libgfortran", v"3")
    @test parse_dl_name_version("libgfortran.so.3.4", "linux") == ("libgfortran", v"3.4")
    @test_throws ArgumentError parse_dl_name_version("libgfortran.so.3.4a", "linux")
    @test_throws ArgumentError parse_dl_name_version("libgfortran", "linux")
    @test_throws ArgumentError parse_dl_name_version("libgfortranso", "linux")
    @test parse_dl_name_version("libgfortran.so", "freebsd") == ("libgfortran", nothing)
    @test parse_dl_name_version("libgfortran.so.3", "freebsd") == ("libgfortran", v"3")
    @test parse_dl_name_version("libgfortran.so.3.4", "freebsd") == ("libgfortran", v"3.4")
    @test_throws ArgumentError parse_dl_name_version("libgfortran.so.3.4a", "freebsd")
    @test_throws ArgumentError parse_dl_name_version("libgfortran", "freebsd")
    @test_throws ArgumentError parse_dl_name_version("libgfortranso", "freebsd")
    @test parse_dl_name_version("libgfortran.dylib", "macos") == ("libgfortran", nothing)
    @test parse_dl_name_version("libgfortran.3.dylib", "macos") == ("libgfortran", v"3")
    @test parse_dl_name_version("libgfortran.3.4.dylib", "macos") == ("libgfortran", v"3.4")
    @test parse_dl_name_version("libgfortran.3.4a.dylib", "macos") == ("libgfortran.3.4a", nothing)
    @test_throws ArgumentError parse_dl_name_version("libgfortran", "macos")
    @test_throws ArgumentError parse_dl_name_version("libgfortrandylib", "macos")
    @test parse_dl_name_version("libgfortran.dll", "windows") == ("libgfortran", nothing)
    @test parse_dl_name_version("libgfortran-3.dll", "windows") == ("libgfortran", v"3")
    @test parse_dl_name_version("libgfortran-3.4.dll", "windows") == ("libgfortran", v"3.4")
    @test parse_dl_name_version("libgfortran-3.4a.dll", "windows") == ("libgfortran-3.4a", nothing)
    @test_throws ArgumentError parse_dl_name_version("libgfortran", "windows")
    @test_throws ArgumentError parse_dl_name_version("libgfortrandll", "windows")
end

@testset "Sys.is* overloading" begin
    # Test that we can indeed ask if something is linux or windows, etc...
    @test Sys.islinux(P("aarch64", "linux"))
    @test !Sys.islinux(P("x86_64", "windows"))
    @test Sys.iswindows(P("i686", "windows"))
    @test !Sys.iswindows(P("powerpc64le", "linux"))
    @test Sys.isapple(P("x86_64", "macos"))
    @test !Sys.isapple(P("armv7l", "windows"))
    @test Sys.isbsd(P("aarch64", "macos"))
    @test Sys.isbsd(P("x86_64", "freebsd"))
    @test !Sys.isbsd(P("x86_64", "linux"; libc="musl"))
end

@testset "Compiler ABI detection" begin
    # Let's check and ensure that we can autodetect the currently-running Julia process
    @test detect_libgfortran_version() !== nothing

    # We run these to get coverage, but we can't test anything, because we could be built
    # with `clang`, which wouldn't have any `libstdc++` constraints at all
    detect_libstdcxx_version()
    detect_cxxstring_abi()
end

@testset "select_platform" begin
    platforms = Dict(
        # Typical binning test
        P("x86_64", "linux"; libgfortran_version=v"3") => "linux4",
        P("x86_64", "linux"; libgfortran_version=v"4") => "linux7",
        P("x86_64", "linux"; libgfortran_version=v"5") => "linux8",

        # Ambiguity test
        P("aarch64", "linux"; libgfortran_version=v"3") => "linux4",
        P("aarch64", "linux"; libgfortran_version=v"3", libstdcxx_version=v"3.4.18") => "linux5",

        # OS test
        P("x86_64", "macos"; libgfortran_version=v"3") => "mac4",
        P("x86_64", "windows"; cxxstring_abi=:cxx11) => "win",
    )

    @test select_platform(platforms, P("x86_64", "linux")) == "linux8"
    @test select_platform(platforms, P("x86_64", "linux"; libgfortran_version=v"4")) == "linux7"

    # Ambiguity test
    @test select_platform(platforms, P("aarch64", "linux")) == "linux5"
    @test select_platform(platforms, P("aarch64", "linux"; libgfortran_version=v"3")) == "linux5"
    @test select_platform(platforms, P("aarch64", "linux"; libgfortran_version=v"4")) === nothing

    @test select_platform(platforms, P("x86_64", "macos")) == "mac4"
    @test select_platform(platforms, P("x86_64", "macos"; libgfortran_version=v"4")) === nothing

    @test select_platform(platforms, P("x86_64", "windows"; cxxstring_abi="cxx11")) == "win"
    @test select_platform(platforms, P("x86_64", "windows"; cxxstring_abi="cxx03")) === nothing

    # Sorry, Alex. ;)
    @test select_platform(platforms, P("x86_64", "freebsd")) === nothing
end

@testset "Custom comparators" begin
    # We're going to define here some custom comparators for Platform objects to ensure they work.
    # First, a symmetric one, which doesn't care which `Platform` object requested this comparison:
    function matches_oddness(a::String, b::String, a_requested::Bool, b_requested::Bool)
        return (parse(Int, a) % 2) == (parse(Int, b) % 2)
    end

    comp_strat = Dict("vally" => matches_oddness)

    # First, test that these two do not match, because it's using equality to test the `vally` tag
    a = Platform("x86_64", "linux"; vally="2")
    b = Platform("x86_64", "linux"; vally="4")
    @test !platforms_match(a, b)

    # Now, test that setting one or both `Platform`'s to use the `matches_oddness()` comparator works:
    ac = Platform("x86_64", "linux"; vally="2", compare_strategies=comp_strat)
    bc = Platform("x86_64", "linux"; vally="4", compare_strategies=comp_strat)
    @test platforms_match(ac, b)
    @test platforms_match(a, bc)
    @test platforms_match(ac, bc)

    # Test that even with the comparison strat, we don't match if they're not both even:
    bfc = Platform("x86_64", "linux"; vally="3", compare_strategies=comp_strat)
    @test !platforms_match(ac, bfc)

    # Next, an asymmetric comparison strategy.  We'll create a "less than or equal to" constraint
    # that uses the `{a,b}_requested` parameters to determine which number represents the limit.
    function less_than_constraint(a::String, b::String, a_requested::Bool, b_requested::Bool)
        a = parse(Int, a)
        b = parse(Int, b)
        if a_requested && !b_requested
            return b < a
        end
        if b_requested && !a_requested
            return a < b
        end
        # If two constraints have been requested, return true if they are the same constraint.
        return a == b
    end

    comp_strat = Dict("vally" => less_than_constraint)
    a = Platform("x86_64", "linux"; vally="2")
    b = Platform("x86_64", "linux"; vally="4")

    ac = Platform("x86_64", "linux"; vally="2", compare_strategies=comp_strat)
    bc = Platform("x86_64", "linux"; vally="4", compare_strategies=comp_strat)

    # Vanilla comparison doesn't work
    @test !platforms_match(a, b)

    # a and bc match, but not ac and b.  Also test reciprocity.
    @test platforms_match(a, bc)
    @test platforms_match(bc, a)
    @test !platforms_match(ac, b)
    @test !platforms_match(b, ac)

    # ac and bc do not match, but ac and ac do
    @test !platforms_match(ac, bc)
    @test platforms_match(ac, ac)
    @test platforms_match(bc, bc)
end

max_stars_repo_name: jonas-schulze/julia

max_issues_repo_name: jonas-schulze/julia

max_forks_repo_name: jonas-schulze/julia


text:
function init_cfg()
    cfg = Dict{String, String}()

    cfg["resistance_file"] = "(browse to resistance file)"
    cfg["resistance_file_is_conductance"] = "false"
    cfg["source_file"] = "(browse to source file)"
    cfg["project_name"] = "(filename prefix for project)"

    cfg["parallelize"] = "true"
    cfg["parallel_batch_size"] = "10"

    cfg["block_size"] = "1"
    cfg["radius"] = "(choose source radius)"
    cfg["buffer"] = "0"
    cfg["source_threshold"] = "0"
    cfg["source_from_resistance"] = "false"
    cfg["r_cutoff"] = "0.0"
    cfg["precision"] = "double"

    cfg["calc_flow_potential"] = "false"
    cfg["calc_normalized_current"] = "false"

    cfg["correct_artifacts"] = "true"

    cfg["write_raw_currmap"] = "true"
    cfg["write_as_tif"] = "true"
    cfg["mask_nodata"] = "true"

    cfg["conditional"] = "false"
    cfg["n_conditions"] = "1"
    cfg["compare_to_future"] = "none"
    cfg["condition1_file"] = "(browse to condition 1 file)"
    cfg["condition2_file"] = "(browse to condition 2 file)"
    cfg["condition1_future_file"] = "(browse to future condition 1 file)"
    cfg["condition2_future_file"] = "(browse to future condition 2 file)"
    cfg["comparison1"] = "within"
    cfg["comparison2"] = "within"
    cfg["condition1_lower"] = "0"
    cfg["condition2_lower"] = "0"
    cfg["condition1_upper"] = "0"
    cfg["condition2_upper"] = "0"

    cfg["allow_different_projections"] = "false"

    cfg
end

function update_cfg!(cfg, cfg_new)
    for (key,val) in cfg_new
        cfg[key] = val
    end
end

function parse_cfg(path::String)
    cf = Dict{String, String}()
    f = open(path, "r")
    for i in eachline(f, keep = true)
        if first(i) == '['
            continue
        end
        idx = something(findfirst(isequal('='), i), 0)
        var = rstrip(i[1:idx - 1])
        val = strip(i[idx + 1:end])
        cf[var] = val
    end
    close(f)
    cf
end

function init_csdict(cfg)
    a = Dict{String, String}()

    a["ground_file_is_resistances"] = "True"
    a["use_direct_grounds"] = "False"

    a["output_file"] = "temp"
    a["write_cum_cur_map_only"] =  "False"

    a["scenario"] = "Advanced"

    a
end


text:
function _get_numsteps(Tsample, δ, ζ, ::DeterministicSwitching)
    αlow = Tsample / δ
    NLOW = ceil(Int, αlow)
    return NLOW, NLOW
end

function _get_numsteps(Tsample, δ, ζ, ::NonDeterministicSwitching)
    ζ⁻ = inf(ζ)
    αlow = (Tsample + ζ⁻)/δ
    NLOW = ceil(Int, αlow)

    ζ⁺ = sup(ζ)
    αhigh = (Tsample + ζ⁺)/δ
    NHIGH = ceil(Int, αhigh)

    return NLOW, NHIGH
end

function _get_max_jumps(tspan, Tsample, δ, ζ, ::DeterministicSwitching)
    NLOW = ceil(Tsample / δ)
    T = tend(tspan)
    aux = (T - (NLOW - 1) * δ) / (NLOW * δ)
    max_jumps = ceil(Int, aux)
    @assert max_jumps >= 0 "inconsistent choice of parameters, `max_jumps` should " *
                           "be positive but it is $max_jumps"

    return max_jumps
end

function _get_max_jumps(tspan, Tsample, δ, ζ, ::NonDeterministicSwitching)
    ζ⁻ = inf(ζ)
    αlow = (Tsample + ζ⁻)/δ
    NLOW = ceil(Int, αlow)
    T = tend(tspan)
    aux = (T - (NLOW - 1) * δ) / (NLOW * δ)
    max_jumps = ceil(Int, aux)
    @assert max_jumps >= 0 "inconsistent choice of parameters, `max_jumps` should " *
                           "be positive but it is $max_jumps"

    return max_jumps
end

const no_tspan = emptyinterval()

function solve(ivp::IVP{<:HACLD1}, args...; kwargs...)

    # preliminary checks
    _check_dim(ivp)

    # get time span (or the emptyset if NSTEPS was specified)
    tsp = _get_tspan(args...; kwargs...)

    # get the continuous post or find a default one
    alg = _get_cpost(ivp, args...; kwargs...)
    if alg == nothing
        alg = _default_cpost(ivp, tsp; kwargs...)
    end

    X0 = initial_state(ivp)
    ha = system(ivp)
    @unpack sys, rmap, Tsample, ζ, switching = ha
    t0 = isempty(tsp) ? 0.0 : tstart(tsp)
    δ = step_size(alg)

    # get maximum number of jumps
    if haskey(kwargs, :max_jumps)
        max_jumps = kwargs[:max_jumps]
    else
        tsp ≠ emptyinterval() || throw(ArgumentError("either the time span `tspan` or the maximum number " *
                                                       "of jumps `max_jumps` should be specified"))
        max_jumps = _get_max_jumps(tsp, Tsample, δ, ζ, switching)
    end

    # number of steps for the continuous post
    prob = IVP(sys, X0)
    NLOW, NHIGH = _get_numsteps(Tsample, δ, ζ, switching)
    ζint = jitter(ha)
    ζ⁻ = inf(ζint)
    ζ⁺ = sup(ζint)
    @assert NLOW > 0 throw(ArgumentError("inconsistent choice of parameters"))

    # solve first flowpipe
    sol = post(alg, prob, no_tspan; NSTEPS=NHIGH)

    if max_jumps == 0
        return ReachSolution(Flowpipe(sol[1:NLOW-1]), alg)
    end

    # preallocate output vector of flowpipes
    N = numtype(alg)
    RT = rsetrep(alg)
    SRT = SubArray{RT, 1, Vector{RT}, Tuple{UnitRange{Int}}, true}
    FT = Flowpipe{N, RT, SRT}
    out = Vector{ShiftedFlowpipe{FT, N}}()
    sizehint!(out, max_jumps+1)

    push!(out, ShiftedFlowpipe(Flowpipe(view(array(sol), 1:NHIGH)), t0))

    # prepare successor for next jump
    Xend = _transition_successors(sol, NLOW, NHIGH, switching)
    t0 += Tsample + ζ⁻

    # adjust integer bounds for subsequent jumps
    NHIGH += ceil(Int, abs(ζ⁻) / δ)

    @inbounds for k in 2:max_jumps+1

        # apply reset map and instantiate new initial-value problem
        prob = IVP(sys, rmap(Xend))

        # solve next chunk
        sol = post(alg, prob, no_tspan; NSTEPS=NHIGH)

        # store flowpipe until first intersection with the guard
        aux = view(array(sol), 1:NHIGH)

        push!(out, ShiftedFlowpipe(Flowpipe(aux), t0))

        # get successors after discrete jump from Tsample + ζ⁻ .. Tsample + ζ⁺
        Xend = _transition_successors(sol, NLOW, NHIGH, switching)

        # adjust initial time for next jump
        t0 += Tsample
    end

    return ReachSolution(HybridFlowpipe(out), alg)
end

# deterministic switching
@inline function _transition_successors(sol, NLOW, NHIGH, ::DeterministicSwitching)
    # in this scenario, NLOW == NHIGH most of the time
    # sol[NLOW] |> set
    # however, it may happen that NHIGH = NLOW + 1 if Tsample is a multiple of δ
    if NLOW == NHIGH
        # TODO: can we use dispach and remove this branch?
        sol[NLOW] |> set
    else
        view(array(sol), NLOW:NHIGH) |> Flowpipe |> convexify |> set
    end
end

# non-deterministic switching
@inline function _transition_successors(sol, NLOW, NHIGH, ::NonDeterministicSwitching)
    view(array(sol), NLOW:NHIGH) |> Flowpipe |> convexify |> set
end


text:
struct VroomSolver
    exploration_level::Int
    nb_threads::Int
    with_tw::Bool
end

function VroomSolver(;exploration_level::Int = 1, nb_threads::Int = 1,
                     with_tw::Bool = false)
    return VroomSolver(exploration_level, nb_threads, with_tw)
end

const VROOM_ERROR_INPUT_INSTANTIATION = 1
const VROOM_ERROR_ADD_MATRIX = 2
const VROOM_ERROR_ADD_VEHICLE = 3
const VROOM_ERROR_ADD_JOB = 4
const VROOM_ERROR_SOLVE = 5

macro vroom_ccall(func, args...)
    args = map(esc,args)
    f = "vroom_$(func)"
    quote
        ccall(($f,$(string(dirname(@__FILE__),"/c_interface/libvroominterface"))), $(args...))
    end
end

abstract type VroomException <: Exception end
struct VroomInputException   <: VroomException end
struct VroomMatrixException  <: VroomException end
struct VroomVehicleException <: VroomException
    vehicle_idx::Int
end
struct VroomJobException     <: VroomException
    job_idx::Int
end
struct VroomSolveException   <: VroomException end

function Base.showerror(io::IO, e::VroomInputException)
    print(io, "Vroom failed to instantiate the input")
end
function Base.showerror(io::IO, e::VroomMatrixException)
    print(io, "Vroom failed to load the matrix in the input")
end
function Base.showerror(io::IO, e::VroomJobException)
    print(io, "Vroom failed to load the job $(e.job_idx) in the input")
end
function Base.showerror(io::IO, e::VroomVehicleException)
    print(io, "Vroom failed to load the vehicle $(e.vehicle_idx) in the input")
end
function Base.showerror(io::IO, e::VroomSolveException)
    print(io, "Vroom failed to solve the input")
end

########################################################
############## vroom interface functions ###############
########################################################
function vroom_create_input(data::RvrpInstance,
                            computed_data::RvrpComputedData, with_tw::Bool)

    vroomdata = Ref{Ptr{Cvoid}}()
    status = @vroom_ccall input Cint (Bool, Ref{Ptr{Cvoid}},) false vroomdata
    if status == VROOM_ERROR_INPUT_INSTANTIATION
        throw(VroomInputException())
    end
    @assert status == 0

    # only one matrix is supported
    mat = data.travel_specifications[1].travel_distance_matrix
    # Distance matrix: vroom does not support floats
    mat_size = size(mat, 1)
    array::Vector{Int32} = reshape(mat, (mat_size*mat_size))
    @vroom_ccall add_matrix Cvoid (Ptr{Cvoid}, Cint, Ptr{Cint},
                 ) vroomdata[] Int32(mat_size) array

    vroom_add_jobs(vroomdata[], computed_data, data.requests, data.location_groups,
             data.product_specification_classes)
    vroom_add_vehicles(vroomdata[], with_tw, computed_data, data.vehicle_categories,
                 data.location_groups, data.vehicle_sets)
    return vroomdata[]
end

function vroom_add_jobs(data_ptr::Ptr{Cvoid},
                  computed_data::RvrpComputedData,
                  requests::Vector{Request},
                  lgs::Vector{LocationGroup},
                  product_specification_classes::Vector{ProductSpecificationClass})
    for r_idx in 1:length(requests)
        req = requests[r_idx]
        nb_tw = length(req.pickup_time_windows)
        tw_starts::Vector{Cint} = [
            req.pickup_time_windows[i].soft_range.lb for i in 1:nb_tw
        ]
        tw_ends::Vector{Cint} = [
            req.pickup_time_windows[i].soft_range.ub for i in 1:nb_tw
        ]
        lg_idx = computed_data.location_group_id_2_index[req.pickup_location_group_id]
        l_idx = computed_data.location_id_2_index[lgs[lg_idx].location_ids[1]]
        d = req.pickup_service_time
        c = get_capacity_consumptions(req, product_specification_classes,
                                           computed_data)

        status = @vroom_ccall(add_job, Cint, (Ptr{Cvoid}, Cint, Cint, Cint,
                Cint, Cint, Ptr{Cint}, Ptr{Cint}), data_ptr, UInt64(r_idx-1),
                UInt16(l_idx-1), UInt32(d), Int64(first(c)[2]), nb_tw,
                tw_starts, tw_ends)
        if status == VROOM_ERROR_ADD_JOB
            throw(VroomJobException(r_idx))
        end
        @assert status == 0
    end
end

function vroom_add_vehicles(data_ptr::Ptr{Cvoid}, with_tw::Bool,
                      computed_data::RvrpComputedData,
                      categories::Vector{VehicleCategory},
                      lgs::Vector{LocationGroup},
                      vehicle_sets::Vector{HomogeneousVehicleSet})
    v_idx = 1
    for v_set in vehicle_sets
        cat_idx = computed_data.vehicle_category_id_2_index[v_set.vehicle_category_id]
        caps = [v for (k,v) in categories[cat_idx].capacities.of_vehicle]

        lg_idx = computed_data.location_group_id_2_index[v_set.departure_location_group_id]
        s_idx = computed_data.location_id_2_index[lgs[lg_idx].location_ids[1]]
        lg_idx = computed_data.location_group_id_2_index[v_set.arrival_location_group_id]
        e_idx = computed_data.location_id_2_index[lgs[lg_idx].location_ids[1]]

        tws = v_set.work_periods[1].soft_range.lb
        twe = v_set.work_periods[1].soft_range.ub
        for i in 1:v_set.nb_of_vehicles_range.soft_range.ub
            status = @vroom_ccall(add_vehicle, Cint, (Ptr{Cvoid}, Cint, Cint,
                    Cint, Cint, Bool, Cint, Cint), data_ptr, UInt64(v_idx-1),
                    Int64(caps[1]), UInt16(s_idx-1), UInt16(e_idx-1), with_tw,
                    UInt32(tws), UInt32(twe))
            if status == VROOM_ERROR_ADD_VEHICLE
                throw(VroomVehicleException(v_idx))
            end
            @assert status == 0
            v_idx += 1
        end
    end
end

function vroom_solve(vroom_input::Ptr{Cvoid}, exploration_level::Int,
                     nb_thread::Int)
    sol = Ref{Ptr{Cvoid}}()
    status = @vroom_ccall(solve, Cint, (Ptr{Cvoid}, Cint, Cint, Ref{Ptr{Cvoid}}),
            vroom_input, Cint(exploration_level), Cint(nb_thread), sol)
    if status == VROOM_ERROR_SOLVE
        throw(VroomSolveException())
    end
    @assert status == 0

    return sol[]
end

function vroom_build_routes(data::RvrpInstance, solver::VroomSolver,
                      computed_data::RvrpComputedData,
                      sol_ptr::Ptr{Cvoid})
    vehicle_index_to_set_index = Dict{Int,Int}()
    vroom_idx = 0 # vehicle index in vroom
    for v_set_idx in 1:length(data.vehicle_sets)
        for i in 1:data.vehicle_sets[v_set_idx].nb_of_vehicles_range.soft_range.ub
            vehicle_index_to_set_index[vroom_idx] = v_set_idx
            vroom_idx += 1
        end
    end
    routes = Route[]
    nb_routes = @vroom_ccall get_nb_routes Cint (Ptr{Cvoid},) sol_ptr
    for r_idx in 1:nb_routes
        r_id = string("route_", r_idx)
        vroom_v_idx = @vroom_ccall get_route_vehicle_id Cint (Ptr{Cvoid}, Cint, ) sol_ptr Cint(r_idx-1)
        v_set_idx = vehicle_index_to_set_index[vroom_v_idx]
        nb_actions = @vroom_ccall get_route_nb_actions Cint (Ptr{Cvoid}, Cint, ) sol_ptr Cint(r_idx-1)
        sequence = Action[]
        for action_idx in 1:nb_actions
            if action_idx == 1 || action_idx == nb_actions
                req_id = ""; action_type = 0
                loc_id = data.location_groups[computed_data.location_group_id_2_index[data.vehicle_sets[v_set_idx].departure_location_group_id]].location_ids[1]
            else
                vroom_job_idx = @vroom_ccall get_step_job_id Cint (Ptr{Cvoid}, Cint, Cint,) sol_ptr Cint(r_idx-1) Cint(action_idx-1)
                req_id = data.requests[vroom_job_idx+1].id; action_type = 1
                loc_id = data.location_groups[computed_data.location_group_id_2_index[data.requests[vroom_job_idx+1].pickup_location_group_id]].location_ids[1]
            end
            scheduled_start_time = @vroom_ccall get_step_arrival Cint (Ptr{Cvoid}, Cint, Cint,) sol_ptr Cint(r_idx-1) Cint(action_idx-1)
            push!(sequence, Action(string("action_", action_idx), loc_id,
                                        action_type, req_id,
                                        scheduled_start_time))
        end
        push!(routes, Route(
            r_id, data.vehicle_sets[v_set_idx].id, sequence, 0
        ))
    end
    return routes
end

function vroom_transform_solution(data::RvrpInstance,
        solver::VroomSolver,computed_data::RvrpComputedData,
        sol_ptr::Ptr{Cvoid})

    id = string(data.id, "_vroom_SOL_", rand(1:10000))
    problem_id = data.id
    cost = @vroom_ccall get_sol_cost Cint (Ptr{Cvoid},) sol_ptr
    routes = vroom_build_routes(data, solver, computed_data, sol_ptr)
    unassigned = String[]
    nb_unassigned = @vroom_ccall get_nb_unassigned Cint (Ptr{Cvoid},) sol_ptr
    for i in 1:nb_unassigned
        vroom_idx = @vroom_ccall get_unassigned_job_id Cint (Ptr{Cvoid}, Cint,) sol_ptr Cint(i-1)
        unassigned_id = data.requests[vroom_idx+1].id
        push!(unassigned, unassigned_id)
    end
    return RvrpSolution(id, problem_id, cost, routes, unassigned)
end

########################################################
############## RvrpSolver functions ####################
########################################################

function solve(data::RvrpInstance, solver::VroomSolver)
    if !check_initial_feasibility(data)
        id = string(data.id, "_vroom_SOL_", rand(1:10000))
        return build_empty_sol(data, id)
    end
    computed_data = build_computed_data(data)
    # check_data(data, solver)
    vroom_ptr = vroom_create_input(data, computed_data, solver.with_tw)
    sol_ptr = vroom_solve(vroom_ptr, solver.exploration_level,
                          solver.nb_threads)
    return vroom_transform_solution(data, solver, computed_data, sol_ptr)
end

function supported_features(::Type{VroomSolver})
    pickonly_features = BitSet()

    # Product features
    union!(pickonly_features, HAS_PRODUCT_CAPACITY_CONSUMPTIONS)

    # Request based features
    union!(pickonly_features, HAS_PICKUPONLY_REQUESTS)
    union!(pickonly_features, HAS_PICKUP_TIME_WINDOWS)
    union!(pickonly_features, HAS_MULTIPLE_PICKUP_TIME_WINDOWS)

    # VehicleCategory based features
    union!(pickonly_features, HAS_VEHICLE_CAPACITIES)

    # HomogeneousVehicleSet based features
    union!(pickonly_features, HAS_MAX_NB_VEHICLES)
    union!(pickonly_features, HAS_WORKING_TIME_WINDOW)
    union!(pickonly_features, HAS_TRAVEL_DISTANCE_UNIT_COST)
    union!(pickonly_features, HAS_ARRIVAL_DIFFERENT_FROM_DEPARTURE)

    # Instance based features
    union!(pickonly_features, HAS_MULTIPLE_VEHICLE_CATEGORIES)
    union!(pickonly_features, HAS_MULTIPLE_VEHICLE_SETS)

    #TODO support these features with a transformation
    deliveronly_features = BitSet()

    # Product features
    union!(pickonly_features, HAS_PRODUCT_CAPACITY_CONSUMPTIONS)

    # Request based features
    union!(deliveronly_features, HAS_DELIVERYONLY_REQUESTS)
    union!(deliveronly_features, HAS_DELIVERY_TIME_WINDOWS)
    union!(deliveronly_features, HAS_MULTIPLE_DELIVERY_TIME_WINDOWS)

    # VehicleCategory based features
    union!(deliveronly_features, HAS_VEHICLE_CAPACITIES)

    # HomogeneousVehicleSet based features
    union!(deliveronly_features, HAS_MAX_NB_VEHICLES)
    union!(deliveronly_features, HAS_WORKING_TIME_WINDOW)
    union!(deliveronly_features, HAS_TRAVEL_DISTANCE_UNIT_COST)
    union!(deliveronly_features, HAS_ARRIVAL_DIFFERENT_FROM_DEPARTURE)

    # Instance based features
    union!(deliveronly_features, HAS_MULTIPLE_VEHICLE_CATEGORIES)
    union!(deliveronly_features, HAS_MULTIPLE_VEHICLE_SETS)

    return [pickonly_features, deliveronly_features]
end


text:
import OrdinaryDiffEq

@testset "test_timeresp" begin

A = [0 1; 0 0]
B = [0;1]
C = [1 0]
sys = ss(A,B,C,0)
Q = I # eye(2)
R = I # eye(1)
L = lqr(sys,Q,R)

u(x,i) = -L*x # Form control law
t=0:0.1:50
x0 = [1.,0]
y, t, x, uout = lsim(sys,u,t,x0=x0) # Continuous time

th = 1e-6
@test sum(abs.(x[end,:])) < th

y, t, x, uout = lsim(c2d(sys,0.1)[1],u,t,x0=x0) # Discrete time
@test sum(abs.(x[end,:])) < th

#Do a manual simulation with uout
ym, tm, xm = lsim(sys, uout, t, x0=x0)
@test y ≈ ym
@test x ≈ xm

# Now compare to closed loop
# Discretization is needed before feedback
sysd = c2d(sys, 0.1)[1]
# Create the closed loop system
sysdfb = ss(sysd.A-sysd.B*L, sysd.B, sysd.C, sysd.D, 0.1)
#Simulate without input
yd, td, xd = lsim(sysdfb, zeros(501), t, x0=x0)
@test y ≈ yd
@test x ≈ xd

# Error for nonuniformly spaced vector
@test_throws ErrorException lsim(sys, [1, 2, 3, 4], [0, 1, 1, 2])

# Error for mismatch of sample rates of discrete-system and signal
@test_throws ErrorException lsim(sysd, [1, 2, 3, 4], 0:0.2:0.6)

# Test for problem with broadcast dimensions
@test lsim(sys, zeros(5), 0:0.2:0.8)[1][:] == zeros(5)

# lsim for Int system with Float64 input (regression test for #264)
@test lsim(ss(1,1,1,1,1), ones(5), 0:4)[1][:] == 1:5

# Various combinations of BigFloat
@test lsim(ss(big.(1),1,1,1,1), ones(5), 0:4)[1][:] == 1:5
@test lsim(ss(big.(1),1,1,1,1), big.(ones(5)), 0:4)[1][:] == 1:5
@test lsim(ss(1,1,1,1,1), big.(ones(5)), 0:4)[1][:] == 1:5

# Tests for HeteroStateSpace
@test lsim(HeteroStateSpace(big.(1.0),1,1,1,1), ones(5), 0:4)[1][:] == 1:5

# lsim for discrete-time complex-coefficient systems

# Complex system, real input signal
G1 = ss(1.0im, 1)
@test lsim(G1, [1, 2], 0:1)[1][:] == 1.0im*[1, 2]
@test impulse(G1, 4)[1][:] == [1.0im; zeros(4)]
@test step(G1, 4)[1][:] == fill(1.0im, 5)

G2 = ss(1.0im, 1, 1, 0, 1)
@test lsim(G2, ones(3), 0:2)[1][:] == [0.0, 1, 1 + im]
@test impulse(G2, 4)[1][:] == [0.0, 1, im, -1, -im]
@test step(G2, 4)[1][:] == [0.0, 1, 1+im, 1+im-1, 0]

# Real system, complex input signal
@test lsim(ss(1, 1), [1.0, 1.0im], 0:1)[1][:] == [1.0, 1.0im]
@test lsim(ss(1.0, 1), [1.0, 1.0im], 0:1)[1][:] == [1.0, 1.0im]
@test lsim(ss(1, 1, 1, 0, 1), [1.0, 1im, 1], 0:2)[1][:] == [0.0, 1, 1 + im]

# Complex system, complex input signal
@test lsim(ss(1.0im, 1), [1im, 2], 0:1)[1][:] == [-1, 2im]
@test lsim(ss(1.0im, 1, 1, 0, 1), [1.0, 1im, 1], 0:2)[1][:] == [0.0, 1, 2im]

# Test that the discrete lsim accepts u function that returns scalar
L = lqr(sysd,Q,R)
u(x,i) = -L*x
yd, td, xd = lsim(sysd, u, t, x0=x0)
@test norm(y - yd)/norm(y) < 0.05 # Since the cost matrices are not discretized, these will differ a bit
@test norm(x - xd)/norm(x) < 0.05

# Test lsim with default time vector
uv = randn(length(t))
y,t = lsim(c2d(sys,0.1)[1],uv,t,x0=x0)
yd,td = lsim(c2d(sys,0.1)[1],uv,x0=x0)
@test yd == y
@test td == t

#Test step and impulse Continuous
t0 = 0:0.05:2
systf = [tf(1,[1,1]) 0; 0 tf([1,-1],[1,2,1])]
sysss = ss([-1 0 0; 0 -2 -1; 0 1 0], [1 0; 0 1; 0 0], [1 0 0; 0 1 -1], 0)
yreal = zeros(length(t0), 2, 2)
xreal = zeros(length(t0), 3, 2)

#Step tf
y, t, x = step(systf, t0)
yreal[:,1,1] = 1 .- exp.(-t)
yreal[:,2,2] = -1 .+ exp.(-t)+2*exp.(-t).*t
@test y ≈ yreal atol=1e-4
#Step ss
y, t, x = step(sysss, t)
@test y ≈ yreal atol=1e-4
xreal[:,1,1] = yreal[:,1,1]
xreal[:,2,2] = exp.(-t).*t
xreal[:,3,2] = exp.(-t).*(-t .- 1) .+ 1
@test x ≈ xreal atol=1e-5

#Impulse tf
y, t, x = impulse(systf, t)
yreal[:,1,1] = exp.(-t)
yreal[:,2,2] = exp.(-t).*(1 .- 2*t)
@test y ≈ yreal atol=1e-2
#Impulse ss
y, t, x = impulse(sysss, t)
@test y ≈ yreal atol=1e-2
xreal[:,1,1] = yreal[:,1,1]
xreal[:,2,2] = -exp.(-t).*t + exp.(-t)
xreal[:,3,2] = exp.(-t).*t
@test x ≈ xreal atol=1e-1

#Test step and impulse Discrete
t0 = 0:0.05:2
systf = [tf(1,[1,1]) 0; 0 tf([1,-1],[1,2,1])]
sysss = ss([-1 0 0; 0 -2 -1; 0 1 0], [1 0; 0 1; 0 0], [1 0 0; 0 1 -1], 0)
yreal = zeros(length(t0), 2, 2)
xreal = zeros(length(t0), 3, 2)

#Step tf
y, t, x = step(systf, t0, method=:zoh)
yreal[:,1,1] = 1 .- exp.(-t)
yreal[:,2,2] = -1 .+ exp.(-t) + 2*exp.(-t).*t
@test y ≈ yreal atol=1e-14
#Step ss
y, t, x = step(sysss, t, method=:zoh)
@test y ≈ yreal atol=1e-13
xreal[:,1,1] = yreal[:,1,1]
xreal[:,2,2] = exp.(-t).*t
xreal[:,3,2] = exp.(-t).*(-t .- 1) .+ 1
@test x ≈ xreal atol=1e-14

#Impulse tf
y, t, x = impulse(systf, t, method=:zoh)
yreal[:,1,1] = exp.(-t)
yreal[:,2,2] = exp.(-t).*(1 .- 2*t)
@test y ≈ yreal atol=1e-14
#Impulse ss
y, t, x = impulse(1.0sysss, t, method=:zoh)
@test y ≈ yreal atol=1e-13
xreal[:,1,1] = yreal[:,1,1]
xreal[:,2,2] = -exp.(-t).*t + exp.(-t)
xreal[:,3,2] = exp.(-t).*t
@test x ≈ xreal atol=1e-13

#Step response of discrete system with specified final time
G = tf([1], [1; zeros(3)], 1)
y, t2, x = step(G, 10)
@test y ≈ [zeros(3); ones(8)] atol=1e-5
@test t2 == 0:1:10 # isapprox is broken for ranges (julia 1.3.1)

#Impulse response of discrete system to final time that is not mulitple of the sample time
G = tf([1], [1; zeros(3)], 0.3)
y, t2, x = step(G, 2)
@test y ≈ [zeros(3); ones(4)] atol=1e-5
@test t2 ≈ 0:0.3:1.8 atol=1e-5

#Make sure t was never changed
@test t0 == t

@testset "Simulators" begin

h              = 0.1
Tf             = 20
t              = 0:h:Tf
P              = ss(tf(1,[2,1])^2)
reference(x,t) = [1.]
s              = Simulator(P, reference)
x0             = [0.,0]
tspan          = (0.0,Tf)
sol            = solve(s, x0, tspan, OrdinaryDiffEq.Tsit5())
@test step(P,t)[3] ≈ reduce(hcat,sol.(t))'
end

end


text:
using PMDPs
using PMDPs.LP
using POMDPs
using BSON, CSV
using DrWatson
using RandomNumbers.Xorshifts
using Random
using MCTS, DiscreteValueIteration
using DataFrames

using POMDPSimulators
using StaticArrays, Distributions # load pp

N = 3
# LP
for expected_res in 50:50:1200
    # pp_params = Dict(pairs((nᵣ=3, c=3, T=10, expected_res=3., res_budget_μ=5.)))
    # pp_params = Dict(pairs((nᵣ=10, c=5, T=100, expected_res=100., res_budget_μ=5.)))
    pp_params = Dict(pairs((nᵣ=10, c=40, T=1000, expected_res=Float64(expected_res), res_budget_μ=5.)))
    name = "linear_problem"
    display("Evaluating $name with $pp_params")
    sname = savename("traces_lp", pp_params,  "bson")
    # data = load(datadir("traces", sname))
    data = PMDPs.load_traces(datadir("traces", sname))

    # PMDPs.process_data(data, PMDPs.flatrate; N=N)
    PMDPs.process_data(data, PMDPs.hindsight; N=N)
    # PMDPs.process_data(data, PMDPs.vi; N=N)
    # PMDPs.process_data(data, PMDPs.mcts; N=N)
end
println("Done.")

for seed in 1:2
    for expected_res in 25:25:600
        pp_params = Dict(pairs((NV=8, NE=20, seed=seed, NP=50, c=10, T=1000, expected_res=Float64(expected_res), res_budget_μ=5.)))
        name = "graph_problem"
        display("Evaluating $name with $pp_params,")
        sname = savename("traces_gp", pp_params,  "bson")
        data = load(datadir("traces", sname))
        data = PMDPs.load_traces(datadir("traces", sname))

        # PMDPs.process_data(data, PMDPs.flatrate; N=N)
        PMDPs.process_data(data, PMDPs.hindsight; N=N)
        # PMDPs.process_data(data, PMDPs.vi; N=N)
        # PMDPs.process_data(data, PMDPs.mcts; N=N)
    end
end


text:
using ForwardDiff
using StructArrays

function zigzag(m::ConditionalModel, T = 1000.0; c=10.0, adapt=false) where {A,B}

    ℓ = Base.Fix1(logdensity, m)

    t = as(m)

    function f(x)
        (θ, logjac) = transform_and_logjac(t, x)
        -ℓ(θ) - logjac
    end

    d = t.dimension

    function partiali()
        ith = zeros(d)
        function (x,i)
            ith[i] = 1
            sa = StructArray{ForwardDiff.Dual{}}((x, ith))
            δ = f(sa).partials[]
            ith[i] = 0
            return δ
        end
    end

    ∇ϕi = partiali()

    # Draw a random starting points and velocity

    t0 = 0.0
    x0 = randn(d)
    θ0 = randn(d)

    pdmp(∇ϕi, t0, x0, θ0, T, c*ones(d), ZigZag(sparse(I(d)), 0*x0); adapt=adapt)

end

m = @model x begin
    α ~ Normal()
    β ~ Normal()
    yhat = α .+ β .* x
    y ~ For(eachindex(x)) do j
        Normal(yhat[j], 2.0)
    end
end

x = randn(3);
truth = [0.61, -0.34, -1.74];

post = m(x=x) | (y=truth,)

trace, final, (num, acc) = @time zigzag(post, c=10)

ts, xs = ZigZagBoomerang.sep(discretize(trace, 0.1))

p = lines(ts, getindex.(xs, 1))
lines!(ts, getindex.(xs, 2), color=:red)

p


text:
export kernel_pca

using mlpack._Internal.io

import mlpack_jll
const kernel_pcaLibrary = mlpack_jll.libmlpack_julia_kernel_pca

# Call the C binding of the mlpack kernel_pca binding.
function kernel_pca_mlpackMain()
  success = ccall((:kernel_pca, kernel_pcaLibrary), Bool, ())
  if !success
    # Throw an exception---false means there was a C++ exception.
    throw(ErrorException("mlpack binding error; see output"))
  end
end

" Internal module to hold utility functions. "
module kernel_pca_internal
  import ..kernel_pcaLibrary

end # module

"""
    kernel_pca(input, kernel; [bandwidth, center, degree, kernel_scale, new_dimensionality, nystroem_method, offset, sampling, verbose])

This program performs Kernel Principal Components Analysis (KPCA) on the
specified dataset with the specified kernel.  This will transform the data onto
the kernel principal components, and optionally reduce the dimensionality by
ignoring the kernel principal components with the smallest eigenvalues.

For the case where a linear kernel is used, this reduces to regular PCA.

The kernels that are supported are listed below:

 * 'linear': the standard linear dot product (same as normal PCA):
    K(x, y) = x^T y

 * 'gaussian': a Gaussian kernel; requires bandwidth:
    K(x, y) = exp(-(|| x - y || ^ 2) / (2 * (bandwidth ^ 2)))

 * 'polynomial': polynomial kernel; requires offset and degree:
    K(x, y) = (x^T y + offset) ^ degree

 * 'hyptan': hyperbolic tangent kernel; requires scale and offset:
    K(x, y) = tanh(scale * (x^T y) + offset)

 * 'laplacian': Laplacian kernel; requires bandwidth:
    K(x, y) = exp(-(|| x - y ||) / bandwidth)

 * 'epanechnikov': Epanechnikov kernel; requires bandwidth:
    K(x, y) = max(0, 1 - || x - y ||^2 / bandwidth^2)

 * 'cosine': cosine distance:
    K(x, y) = 1 - (x^T y) / (|| x || * || y ||)

The parameters for each of the kernels should be specified with the options
`bandwidth`, `kernel_scale`, `offset`, or `degree` (or a combination of those
parameters).

Optionally, the Nystroem method ("Using the Nystroem method to speed up kernel
machines", 2001) can be used to calculate the kernel matrix by specifying the
`nystroem_method` parameter. This approach works by using a subset of the data
as basis to reconstruct the kernel matrix; to specify the sampling scheme, the
`sampling` parameter is used.  The sampling scheme for the Nystroem method can
be chosen from the following list: 'kmeans', 'random', 'ordered'.

For example, the following command will perform KPCA on the dataset `input`
using the Gaussian kernel, and saving the transformed data to `transformed`:

```julia
julia> using CSV
julia> input = CSV.read("input.csv")
julia> transformed = kernel_pca(input, "gaussian")
```

# Arguments

 - `input::Array{Float64, 2}`: Input dataset to perform KPCA on.
 - `kernel::String`: The kernel to use; see the above documentation for
      the list of usable kernels.
 - `bandwidth::Float64`: Bandwidth, for 'gaussian' and 'laplacian'
      kernels.  Default value `1`.

 - `center::Bool`: If set, the transformed data will be centered about the
      origin.  Default value `false`.

 - `degree::Float64`: Degree of polynomial, for 'polynomial' kernel.
      Default value `1`.

 - `kernel_scale::Float64`: Scale, for 'hyptan' kernel.  Default value
      `1`.

 - `new_dimensionality::Int`: If not 0, reduce the dimensionality of the
      output dataset by ignoring the dimensions with the smallest eigenvalues.
      Default value `0`.

 - `nystroem_method::Bool`: If set, the Nystroem method will be used.
      Default value `false`.

 - `offset::Float64`: Offset, for 'hyptan' and 'polynomial' kernels.
      Default value `0`.

 - `sampling::String`: Sampling scheme to use for the Nystroem method:
      'kmeans', 'random', 'ordered'  Default value `kmeans`.

 - `verbose::Bool`: Display informational messages and the full list of
      parameters and timers at the end of execution.  Default value `false`.

# Return values

 - `output::Array{Float64, 2}`: Matrix to save modified dataset to.

"""
function kernel_pca(input,
                    kernel::String;
                    bandwidth::Union{Float64, Missing} = missing,
                    center::Union{Bool, Missing} = missing,
                    degree::Union{Float64, Missing} = missing,
                    kernel_scale::Union{Float64, Missing} = missing,
                    new_dimensionality::Union{Int, Missing} = missing,
                    nystroem_method::Union{Bool, Missing} = missing,
                    offset::Union{Float64, Missing} = missing,
                    sampling::Union{String, Missing} = missing,
                    verbose::Union{Bool, Missing} = missing,
                    points_are_rows::Bool = true)
  # Force the symbols to load.
  ccall((:loadSymbols, kernel_pcaLibrary), Nothing, ());

  IORestoreSettings("Kernel Principal Components Analysis")

  # Process each input argument before calling mlpackMain().
  IOSetParamMat("input", input, points_are_rows)
  IOSetParam("kernel", kernel)
  if !ismissing(bandwidth)
    IOSetParam("bandwidth", convert(Float64, bandwidth))
  end
  if !ismissing(center)
    IOSetParam("center", convert(Bool, center))
  end
  if !ismissing(degree)
    IOSetParam("degree", convert(Float64, degree))
  end
  if !ismissing(kernel_scale)
    IOSetParam("kernel_scale", convert(Float64, kernel_scale))
  end
  if !ismissing(new_dimensionality)
    IOSetParam("new_dimensionality", convert(Int, new_dimensionality))
  end
  if !ismissing(nystroem_method)
    IOSetParam("nystroem_method", convert(Bool, nystroem_method))
  end
  if !ismissing(offset)
    IOSetParam("offset", convert(Float64, offset))
  end
  if !ismissing(sampling)
    IOSetParam("sampling", convert(String, sampling))
  end
  if verbose !== nothing && verbose === true
    IOEnableVerbose()
  else
    IODisableVerbose()
  end

  IOSetPassed("output")
  # Call the program.
  kernel_pca_mlpackMain()

  return IOGetParamMat("output", points_are_rows)
end


text:
# the following methods agglomerate pattern-repl pairs.
const TFC = Union{TF,Char}

# if pattern a is found, pattern b will not be used
is_shadowing(a::Char, b::Char) = a == b
is_shadowing(a::Char, b::String) = length(b) == 1 && a == b[1]
is_shadowing(a::String, b::Char) = length(a) == 1 && a[1] == b
is_shadowing(a::String, b::String) = a == b
function is_shadowing(a::Regex, b::Char)
    m = match(a, string(b))
    m != nothing && m.match == string(b)
end
function is_shadowing(a::Regex, b::String)
    m = match(a, string(b))
    m != nothing && m.match == b
end
is_shadowing(a::Callable, b::Char) = a(b)
is_shadowing(a::Callable, b::String) = length(b) == 1 && a(b[1])
is_shadowing(a::TFC, b::TFC) = a == b

function remove_all_shadowed!(pairs::Vector{<:Pair{<:TFC,<:TS}})
    for i = length(pairs):-1:1
        if any(is_shadowing.(first.(pairs[1:i-1]), first(pairs[i])))
            deleteat!(pairs, i)
        end
    end
    pairs
end

isless_pattern(a::Pair, b::Pair) = isless_pattern(first(a), first(b))

isless_pattern(a::String, b::String) = length(b) < length(a)
isless_pattern(a::String, ::Char) = 1 < length(a)
isless_pattern(::Char, b::String) = length(b) < 1
isless_pattern(a::Union{Char,String}, ::TFC) = true
isless_pattern(a::TFC, b::String) = false
isless_pattern(a::TFC, b::TFC) = false

function sortpatterns!(pairs::AbstractVector)
    sort!(pairs, lt=isless_pattern, alg=Base.DEFAULT_STABLE)
end

# if regex is of the form "aaa(|bbb)*" replace it with array of strings
function expand_pattern(re::Regex)
    s = re.pattern
    words = Union{String,Regex}[]
    io = IOBuffer()
    k = 1
    i = 1
    while !done(s, i)
        c, i = next(s, i)
        if c == '\\'
            done(s, i) && error("regex ends with \\")
            c, i = next(s, i)
            write(io, c)
        elseif c == '|'
            push!(words, String(take!(io)))
            k = i
        elseif c in "[|+*^\\\$"
            push!(words, Regex(s[k:end]))
            k = 0
            break
        else
            write(io, c)
        end
    end
    if k > 0
        push!(words, String(take!(io)))
    end
    n = length(words) - (k == 0)
    sortpatterns!(view(words, 1:n))
    unique(words)
end

expand_pattern(s::AbstractString) = [s]
expand_pattern(f::Base.EqualTo{Char}) = [f.x]
expand_pattern(f::Base.OccursIn) = unique([c for c in f.x])
expand_pattern(x) = [x]

function expand_pairs(pairs::AbstractVector{<:Pair})
    res = Pair{<:TFC,<:TS}[]
    for p in pairs
        append!(res, Pair.(expand_pattern(p.first), p.second))
    end
    remove_all_shadowed!(res)
    sortpatterns!(res)
    unique(res)
end

# Condense strings and characters to a single regex.
# The repl is a mapping of patterns to replacement strings.
# it is assumed, that the pairs are sorted with the
# longer strings before the 1-character strings before
# before regexes, functions, and empty string.
#
function condense_pairs(pairs::AbstractVector{<:Pair})
    n = findlast(p->(p.first isa Union{Char,String}), pairs)
    ( n == nothing || n <= 1 ) && return pairs
    ps = view(pairs, 1:n)
    dict = Dict{String,Any}(Pair.(string.(first.(ps)), last.(ps)))
    pout = Pair[]
    fdict(c) = dict[string(c)]
    if length(ps[1].first) > 1 || length(ps[n].first) < 1
        io = IOBuffer()
        frst = true
        for p in ps
            frst || write(io, '|')
            escape_string(io, string(p.first), "\\\$^[*+()/{}|")
            frst = false
        end
        push!(pout, Pair(Regex(String(take!(io))), fdict))

    elseif length(pairs[n].first) == 1
        push!(pout, Pair(occursin(join(first.(ps))), fdict))
    end
    append!(pout, pairs[n+1:end])
end


text:
"""
Shorthand for `isnothing(optional) ? fallback : optional`
"""
@inline ifnothing(optional, fallback) = isnothing(optional) ? fallback : optional

function Base.foreach(f::Function, contenttype::Type, layout::GridLayout; recursive = true)
    for c in layout.content
        if recursive && c.content isa GridLayout
            foreach(f, contenttype, c.content)
        elseif c.content isa contenttype
            f(c.content)
        end
    end
end

"""
Swaps or rotates the layout positions of the given elements to their neighbor's.
"""
function swap!(layout_elements...)
    gridcontents = gridcontent.(collect(layout_elements))

    # copy relevant fields before gridcontents are mutated
    parents = map(gc -> gc.parent, gridcontents)
    spans = map(gc -> gc.span, gridcontents)
    sides = map(gc -> gc.side, gridcontents)

    for (gc, parent, span, side) in zip(circshift(gridcontents, 1), parents, spans, sides)
        parent[span.rows, span.cols, side] = gc.content
    end
end


text:
type Grammar
  rules::Dict{Symbol, Rule}
end

function show(io::IO,grammar::Grammar)
  println("StringParserPEG.Grammar(Dict{Symbol,StringParserPEG.Rule}(")
  for (sym,rule) in grammar.rules
    println("  $sym => $(string(rule)),")
  end
  println(")")
end

function Grammar(definition::AbstractString)
  (ast,pos,err) = parse(grammargrammar,definition)
  if err!=nothing
    throw(err)
  end
  transform(togrammar,ast)
end

function Grammar(definition::AbstractString, refgrammar::Grammar)
  newgrammar = Grammar(definition)
  return Grammar(merge(refgrammar.rules,newgrammar.rules))
end


text:
function compute_derivatives(obj::AdditionalEquations, arg0::SpacecraftState, arg1::Vector{jdouble})
    return jcall(obj, "computeDerivatives", Vector{jdouble}, (SpacecraftState, Vector{jdouble}), arg0, arg1)
end

function get_name(obj::AdditionalEquations)
    return jcall(obj, "getName", JString, ())
end

function init(obj::AdditionalEquations, arg0::SpacecraftState, arg1::AbsoluteDate)
    return jcall(obj, "init", void, (SpacecraftState, AbsoluteDate), arg0, arg1)
end


text:
using JuMP, EAGO

                     m = Model()

                     EAGO.register_eago_operators!(m)

                     @variable(m, -1 <= x[i=1:4] <= 1)

                     @variable(m, -1.3634297333714107 <= q <= -0.20981772159177886)

                     add_NL_constraint(m, :(sigmoid(-0.25086622252953283 + -0.4158912166699471*sigmoid(0.18862560740304835 + -0.15744636994779526*sigmoid(0.5575683410422836 + 0.3104874388021095*sigmoid(0.8764207742833086 + -0.10444718155116828*$(x[1]) + 0.8852151451027264*$(x[2]) + -0.8060184127084176*$(x[3]) + -0.617200197422199*$(x[4])) + 0.12330243724213741*sigmoid(0.934307327106453 + 0.477402449923221*$(x[1]) + -0.9131747528764835*$(x[2]) + 0.07837215353953164*$(x[3]) + 0.2564992344946111*$(x[4]))) + 0.6894873557418655*sigmoid(0.5818927234969289 + 0.4060546195732435*sigmoid(0.8764207742833086 + -0.10444718155116828*$(x[1]) + 0.8852151451027264*$(x[2]) + -0.8060184127084176*$(x[3]) + -0.617200197422199*$(x[4])) + 0.22410021012329473*sigmoid(0.934307327106453 + 0.477402449923221*$(x[1]) + -0.9131747528764835*$(x[2]) + 0.07837215353953164*$(x[3]) + 0.2564992344946111*$(x[4])))) + 0.15556574286983027*sigmoid(-0.6464546506721769 + -0.3052437329487665*sigmoid(0.5575683410422836 + 0.3104874388021095*sigmoid(0.8764207742833086 + -0.10444718155116828*$(x[1]) + 0.8852151451027264*$(x[2]) + -0.8060184127084176*$(x[3]) + -0.617200197422199*$(x[4])) + 0.12330243724213741*sigmoid(0.934307327106453 + 0.477402449923221*$(x[1]) + -0.9131747528764835*$(x[2]) + 0.07837215353953164*$(x[3]) + 0.2564992344946111*$(x[4]))) + 0.4287563119683666*sigmoid(0.5818927234969289 + 0.4060546195732435*sigmoid(0.8764207742833086 + -0.10444718155116828*$(x[1]) + 0.8852151451027264*$(x[2]) + -0.8060184127084176*$(x[3]) + -0.617200197422199*$(x[4])) + 0.22410021012329473*sigmoid(0.934307327106453 + 0.477402449923221*$(x[1]) + -0.9131747528764835*$(x[2]) + 0.07837215353953164*$(x[3]) + 0.2564992344946111*$(x[4]))))) + sigmoid(-0.9805026843034819 + -0.19393497198753096*sigmoid(0.18862560740304835 + -0.15744636994779526*sigmoid(0.5575683410422836 + 0.3104874388021095*sigmoid(0.8764207742833086 + -0.10444718155116828*$(x[1]) + 0.8852151451027264*$(x[2]) + -0.8060184127084176*$(x[3]) + -0.617200197422199*$(x[4])) + 0.12330243724213741*sigmoid(0.934307327106453 + 0.477402449923221*$(x[1]) + -0.9131747528764835*$(x[2]) + 0.07837215353953164*$(x[3]) + 0.2564992344946111*$(x[4]))) + 0.6894873557418655*sigmoid(0.5818927234969289 + 0.4060546195732435*sigmoid(0.8764207742833086 + -0.10444718155116828*$(x[1]) + 0.8852151451027264*$(x[2]) + -0.8060184127084176*$(x[3]) + -0.617200197422199*$(x[4])) + 0.22410021012329473*sigmoid(0.934307327106453 + 0.477402449923221*$(x[1]) + -0.9131747528764835*$(x[2]) + 0.07837215353953164*$(x[3]) + 0.2564992344946111*$(x[4])))) + -0.01642706539433636*sigmoid(-0.6464546506721769 + -0.3052437329487665*sigmoid(0.5575683410422836 + 0.3104874388021095*sigmoid(0.8764207742833086 + -0.10444718155116828*$(x[1]) + 0.8852151451027264*$(x[2]) + -0.8060184127084176*$(x[3]) + -0.617200197422199*$(x[4])) + 0.12330243724213741*sigmoid(0.934307327106453 + 0.477402449923221*$(x[1]) + -0.9131747528764835*$(x[2]) + 0.07837215353953164*$(x[3]) + 0.2564992344946111*$(x[4]))) + 0.4287563119683666*sigmoid(0.5818927234969289 + 0.4060546195732435*sigmoid(0.8764207742833086 + -0.10444718155116828*$(x[1]) + 0.8852151451027264*$(x[2]) + -0.8060184127084176*$(x[3]) + -0.617200197422199*$(x[4])) + 0.22410021012329473*sigmoid(0.934307327106453 + 0.477402449923221*$(x[1]) + -0.9131747528764835*$(x[2]) + 0.07837215353953164*$(x[3]) + 0.2564992344946111*$(x[4]))))) - $q <= 0.0))

                     @objective(m, Min, q)

                     return m

max_stars_repo_name: PSORLab/RSActivationFunctions

max_issues_repo_name: PSORLab/RSActivationFunctions

max_forks_repo_name: PSORLab/RSActivationFunctions


text:
"""
    agent_step!(agent, model)
    Generalised step function
"""
function agent_step!(agent, model)
    steps = model.step_function
    for func in steps
        func(agent, model)
    end
end

"""
    ABM(agents::Array{<:AbstractAgent}, args...; kwargs...)
    A method for ABM function to accept a list of agent types
"""
function AxBM(agents::Array{<:AbstractAgent}, args...; kwargs...)
    return ABM(Union{agents...}, args...; kwargs...)
end


text:
using ITensors
using ITensorMPOCompression

N=10; #10 sites
NNN=7; #Include up to 7th nearest neighbour interactions
sites = siteinds("S=1/2",N);
H=make_transIsing_MPO(sites,NNN);
is_lower_regular_form(H)==true
pprint(H[2])
orthogonalize!(H)
pprint(H[2])
get_Dw(H)
is_lower_regular_form(H)==true
is_orthogonal(H,left)==true


text:
# %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
# Pass through the 4f Fourier system to get the camera image
# %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

using MAT
cd("C:\\Users\\Xuanwen\\Google Drive\\thesisproj\\PostOBJ")
# %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
# %%%%% Convolution to get the image of microscope  %%%%%
# %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
tic()
STR = matread("bead1um8500_sf20nm.mat")
OTF = STR["otf_bead"]
OTF = convert(Array{Float64,2}, OTF)
STR = matread("PSF8500_sf2um.mat")
PSF = STR["PSF_res"]

display("Convoluting PSF with OTF......")
Img = conv2(OTF,PSF)

toc()

max_stars_repo_name: xwghua/DRiP-public

max_issues_repo_name: xwghua/DRiP-public

max_forks_repo_name: xwghua/DRiP-public