text:
# Here I build an algorithm for maximum/minimum finder on the go.
# So instead of finding the maximum of the CPR for a constant mesh of data points.
# we adaptively guide the Ic calculation so the maximum up to the desired tolerance is found.
# Assumptions: There was no intent to make this algorithm general, it solves
# the problem for continous and periodic sinusoidal-triangular functions on a -pi:pi interval.

# The algorithm does the following.
#   Compute y values for 5 equally spaced x values (pi/2 separation)
#   Fit using Fourier series and find the correct truncation to avoid overfitting (cross validation)
#   Find numerically the minima (actually it doesn't matter) of this fit
#   Compute the y value at this global minima
#   Add the new data to x and y and fit again until tolerance is met or hard threshold is met
#   returns xmin and ymin

"it computes the minima value of some computationally demanding `generating_function` in a fixed
range (-π, π) using an adaptive algorithm based on numerical optimization of  K-cross validated fits
which guide adaptively the minima value search. `tolerance` sets when the algorithm should stop searching maximums. `max_iter` sets an upper limit for
maximal search. It is important to use an initial discrete xmesh that uses the 2pi periodicity to guide the fit"
function adaptive_max_finder(generating_func::Function, model, hp; tolerance = 1e-4, max_iter = 10, kw...)
    println("adaptive algorithm")
    # xinit = collect(-1:0.1:1)
    # xinit = 2 .* rand(5) .- 1ß
    #xinit = collect(-0.5:.5:π+0.5)+ (2 .* rand(9) .- 1)./20#collect(-1:0.25:1) + (2 .* rand(9) .- 1)./10
    xinit = collect(-π:π/4:π+π/4)#collect(-π:π/8:π+π/4)
    xinit, finit, yinit, fitparams = generate_and_fit(xinit, model, generating_func, hp; kw...)
        # cpr_plot(xinit, yinit, model, fitparams)
    bounds = (minimum(xinit), maximum(xinit))
    xmax = minima_search(model, fitparams, bounds)
        # fmax, ymax = generating_func(xinit, finit, xmax; kw...)
    xmax0 = 10.0
    count = 0
    x = xinit
    f = finit
    y = yinit
    while abs(xmax0-xmax) > tolerance
        # println("error:", abs(xmax0-xmax))
        err_old = abs(xmax0-xmax)
        xmax0 = xmax
        x, f, y, ymax, fitparams = generate_and_fit(xmax, x, f, y, model, generating_func, hp; kw...)
            # cpr_plot(x, y, model, fitparams)
        bounds = (minimum(x), maximum(x))
        xmax = minima_search(model, fitparams, bounds)
        err_new = abs(xmax0-xmax)
        err_new > err_old ? break : nothing
        count += 1

        if count > 20
            println("convergence was not reached")
            break
        else nothing end
        return ymax
        # abs(xmax0-xmax) > tolerance ? nothing : return ymax
    end

end
"generate a new data point using a given `model` and an `x` value"
function generate_and_fit(x::Array, model, generating_func, hp; kw...)
    x, f, y = generating_func(x, hp; kw...)
    # xmean = [x[i]/2+x[i+1]/2 for i in 1:length(x)-1]
    fitparams = msqerror(x, y, model)
    return x, f, y, fitparams
end

function generate_and_fit(x::Number, xlist, flist, ylist, model, generating_func, hp; kw...)
    f, y = generating_func(xlist, flist, x, hp; kw...)
    push!(xlist, x)
    push!(flist, f)
    push!(ylist, y)
    fitparams = msqerror(xlist, ylist, model)
    return xlist, flist, ylist, y, fitparams
end

"numerical search of the global minima"
function minima_search(model, fitparams, bounds)
    opt = optimize(x -> model(x, fitparams), bounds[1], bounds[2])
    return Optim.minimizer(opt)
end

"This function computes the optimal fitting parameters for a given dataset consisting in a
truncated up to `custom_max_order` Fourier expansion, it performs a K-cross validation
analysis as well, to prevent overfitting (ML algorithm)"
function msqerror(xdata, ydata, model; K = 1, custom_max_order = missing, train_tol = 1e-3)
    custom_max_order === missing ? max_order = Int64(floor(length(xdata)/2)) :
        max_order = custom_max_order
    # data partition controlled by data availability
    l = length(xdata)
    part = Int64(floor(l/K))
    println("partition size: ", part)
    fitted_params = []
    errorstest = zeros(max_order)
    errorstrain = zeros(max_order)
    perm_ind = shuffle(1:l) #reshufle data
    xshuffled = xdata[perm_ind]
    yshuffled = ydata[perm_ind]
    l ≥ 2K ? nothing : K = 1   #criteria for cross validation(need at least the double of points)
    for i in 1:K #average over different K partitions
        errors_train = []
        errors_test = []
        xnew = vcat(xshuffled[Int64(part*(i-1)+1):end], xshuffled[1:Int64(part*(i-1))])
        ynew = vcat(yshuffled[Int64(part*(i-1)+1):end], yshuffled[1:Int64(part*(i-1))])
        xtrain = xnew[1:part]
        xtest = xnew[part:end]
        ytrain = ynew[1:part]
        ytest = ynew[part:end]
        # compute the fittings for each model and the train and test data set errors
        for order in 1:max_order
            p0 = rand(Int64(2*order))
            fit = curve_fit(model, xtrain, ytrain, p0)
            push!(errors_train, sum(x -> x, (ytrain - model(xtrain, fit.param)).^2))
            push!(errors_test,sum(x -> x, (ytest - model(xtest, fit.param)).^2))
        end
        errorstest += errors_test
        errorstrain += errors_train
    end
    good_training_indices =  findall(x-> x < train_tol, errorstrain)
    while good_training_indices == []
        train_tol *= 2
        @info "tolerance was increased to $(train_tol)"
        good_training_indices = findall(x-> x < train_tol, errorstrain)
    end
    optimised_order = good_training_indices[findmin(ifelse(K > 1, errorstest, errorstrain)[good_training_indices])[2]]
    # now we fit using the optimal order and all data points
    fit = curve_fit(model, xdata, ydata, rand(Int64(2*optimised_order)))
    println("optimised order: ", optimised_order)
    return fit.param
end

fourier_model(x, p) = sum(i ->
    sin_cos_model(x, p[Int64((i-1)*2 +1):Int64(i*2)], i-1), 1:Int64(length(p)/2))

# gen_model(x) = sin_model(x, [1.0 1. 2.0]) + 0.1*randn(length(x))
# gen_model2(x) = 0.5*triang(x, [1.0 1.]) + 0*randn(length(x)) + sin_model(x, [1, 1, .3])
@. sin_cos_model(x, p, i) = p[1]*sin(i*x) + p[2]*cos(i*x)
@. sin_model(x, p) = p[1]*sin(p[2]x+p[3])
@. cos_model(x, p) = p[1]*cos(x/2+p[3])
@. triang(x, p) = p[1]*sawtoothwave(x+p[2])

# using CairoMakie
# function cpr_plot(phase, I, model, fitted_params)
#     fig = Figure(resolution = (600,400), font = "CMU Serif")
#     ax = Axis(fig, xlabel = "phi (pi)", ylabel = "I(phi)", ylabelsize = 22,
#         xlabelsize= 22, xgridstyle=:dash, ygridstyle=:dash, xtickalign = 1,
#         xticksize=10, ytickalign=1, yticksize=10,  xlabelpadding = -10)
#     phase_range = range(minimum(phase), stop = maximum(phase), length = 100)
#     I_range = model(phase_range, fitted_params)
#     scatter!(ax, phase, I, color = :red, markersize = 5)
#     lines!(ax,  phase_range, I_range, color = :black)
#     fig[1,1] = ax
#     #save("filename.png", fig, px_per_unit = 2.0)
#     display(fig)
# end

## Usage
# xdata = range(-π, stop=π, length=20)
# ydata = sin_model(xdata, [1.0 1. 2.0]) + 0.1*randn(length(xdata))
# how to set the order of the expansion so there is no overfitting
# order_exp = 10 # do cross validation later to improve it.
# p0 = rand(Int64(2*order_exp)) #rand(6)

# fit = curve_fit(mod, xdata, ydata, p0)
# cpr_plot(xdata, ydata, mod, fit.param)

# msqerror(xdata, ydata, fourier_model)

#adaptive_max_finder(gen_model2, fourier_model)


text:
###############################################################################
#
# Calculate empirical gsw from the equations
#
###############################################################################
"""
    stomatal_conductance(
                model::EmpiricalStomatalModel{FT},
                leaf::Leaf{FT},
                envir::AirLayer{FT},
                β::FT
    ) where {FT<:AbstractFloat}
    stomatal_conductance(
                model::EmpiricalStomatalModel{FT},
                canopyi::CanopyLayer{FT},
                envir::AirLayer{FT},
                β::FT
    ) where {FT<:AbstractFloat}
    stomatal_conductance(
                model::EmpiricalStomatalModel{FT},
                canopyi::CanopyLayer{FT},
                envir::AirLayer{FT},
                β::FT,
                ind::Int
    ) where {FT<:AbstractFloat}

Steady state gsw from empirical approach given
- `model` [`EmpiricalStomatalModel`](@ref) type empirical model parameter set
- `leaf` [`Leaf`] type struct
- `canopyi` [`CanopyLayer`](@ref) type struct
- `envir` [`AirLayer`] type struct
- `β` Correction factor over the g1 part of an empirical model
- `ind` Nth leaf in the canopy layer
"""
function stomatal_conductance(
            model::ESMBallBerry{FT},
            leaf::Leaf{FT},
            envir::AirLayer{FT},
            β::FT
) where {FT<:AbstractFloat}
    @unpack g0, g1    = model;
    @unpack P_AIR, rh = envir;

    return g0 + g1 * rh * P_AIR * FT(1e-6) * β * leaf.PSM.a_net / leaf.p_CO₂_s
end

function stomatal_conductance(
            model::ESMBallBerry{FT},
            canopyi::CanopyLayer{FT},
            envir::AirLayer{FT},
            β::FT
) where {FT<:AbstractFloat}
    @unpack g0, g1    = model;
    @unpack P_AIR, rh = envir;

    return g0 .+ g1 * rh * P_AIR * FT(1e-6) * β .* canopyi.ps.PSM.a_net / canopyi.ps.p_CO₂_s
end

function stomatal_conductance(
            model::ESMBallBerry{FT},
            canopyi::CanopyLayer{FT},
            envir::AirLayer{FT},
            β::FT,
            ind::Int
) where {FT<:AbstractFloat}
    @unpack g0, g1    = model;
    @unpack An, p_s   = canopyi;
    @unpack P_AIR, rh = envir;

    return g0 + g1 * rh * P_AIR * FT(1e-6) * β * An[ind] / p_s[ind]
end

function stomatal_conductance(
            model::ESMGentine{FT},
            leaf::Leaf{FT},
            envir::AirLayer{FT},
            β::FT
) where {FT<:AbstractFloat}
    @unpack g0, g1  = model;
    @unpack P_AIR   = envir;

    return g0 + g1 * P_AIR * FT(1e-6) * β * leaf.PSM.a_net / leaf.p_CO₂_i
end

function stomatal_conductance(
            model::ESMGentine{FT},
            canopyi::CanopyLayer{FT},
            envir::AirLayer{FT},
            β::FT
) where {FT<:AbstractFloat}
    @unpack g0, g1  = model;
    @unpack P_AIR   = envir;

    return g0 .+ g1 * P_AIR * FT(1e-6) * β .* canopyi.ps.PSM.a_net / canopyi.ps.p_CO₂_i
end

function stomatal_conductance(
            model::ESMGentine{FT},
            canopyi::CanopyLayer{FT},
            envir::AirLayer{FT},
            β::FT,
            ind::Int
) where {FT<:AbstractFloat}
    @unpack g0, g1  = model;
    @unpack An, p_i = canopyi;
    @unpack P_AIR   = envir;

    return g0 + g1 * P_AIR * FT(1e-6) * β * An[ind] / p_i[ind]
end

function stomatal_conductance(
            model::ESMLeuning{FT},
            leaf::Leaf{FT},
            envir::AirLayer{FT},
            β::FT
) where {FT<:AbstractFloat}
    @unpack d0, g0, g1 = model;
    @unpack P_AIR, p_H₂O = envir;

    if typeof(leaf.PSM) <: C4VJPModel
        _γ_s = FT(0);
    else
        _γ_s = leaf.PSM.γ_star;
    end;

    return g0 + g1 * P_AIR * FT(1e-6) / (1 + (leaf.p_H₂O_sat - p_H₂O)/d0) * β * leaf.PSM.a_net / (leaf.p_CO₂_s - _γ_s)
end

function stomatal_conductance(
            model::ESMLeuning{FT},
            canopyi::CanopyLayer{FT},
            envir::AirLayer{FT},
            β::FT
) where {FT<:AbstractFloat}
    @unpack d0, g0, g1     = model;
    @unpack P_AIR, p_H₂O   = envir;

    if typeof(canopyi.ps.PSM) <: C4VJPModel
        _γ_s = FT(0);
    else
        _γ_s = canopyi.ps.PSM.γ_star;
    end;

    return g0 .+ g1 * P_AIR * FT(1e-6) / (1 + (canopyi.ps.p_H₂O_sat - p_H₂O)/d0) * β .* canopyi.ps.PSM.a_net ./ (canopyi.ps.p_CO₂_s .- _γ_s)
end

function stomatal_conductance(
            model::ESMLeuning{FT},
            canopyi::CanopyLayer{FT},
            envir::AirLayer{FT},
            β::FT,
            ind::Int
) where {FT<:AbstractFloat}
    @unpack d0, g0, g1     = model;
    @unpack An, p_s, p_sat = canopyi;
    @unpack P_AIR, p_H₂O   = envir;

    if typeof(canopyi.ps.PSM) <: C4VJPModel
        _γ_s = FT(0);
    else
        _γ_s = canopyi.ps.PSM.γ_star;
    end;

    return g0 + g1 * P_AIR * FT(1e-6) / (1 + (p_sat - p_H₂O)/d0) * β * An[ind] / (p_s[ind] - _γ_s)
end

function stomatal_conductance(
            model::ESMMedlyn{FT},
            leaf::Leaf{FT},
            envir::AirLayer{FT},
            β::FT
) where {FT<:AbstractFloat}
    @unpack g0, g1            = model;
    @unpack p_CO₂, P_AIR, p_H₂O = envir;
    vpd = max(FT(0.001), leaf.p_H₂O_sat - p_H₂O);

    return g0 + P_AIR * FT(1e-6) / p_CO₂ * (1 + g1/sqrt(vpd)) * β * leaf.PSM.a_net
end

function stomatal_conductance(
            model::ESMMedlyn{FT},
            canopyi::CanopyLayer{FT},
            envir::AirLayer{FT},
            β::FT
) where {FT<:AbstractFloat}
    @unpack g0, g1            = model;
    @unpack p_CO₂, P_AIR, p_H₂O = envir;
    vpd = max(FT(0.001), canopyi.ps.p_H₂O_sat - p_H₂O);

    return g0 .+ P_AIR * FT(1e-6) / p_CO₂ * (1 + g1/sqrt(vpd)) * β .* canopyi.ps.PSM.a_net
end

function stomatal_conductance(
            model::ESMMedlyn{FT},
            canopyi::CanopyLayer{FT},
            envir::AirLayer{FT},
            β::FT,
            ind::Int
) where {FT<:AbstractFloat}
    @unpack g0, g1            = model;
    @unpack An, p_sat         = canopyi;
    @unpack p_CO₂, P_AIR, p_H₂O = envir;
    vpd = max(FT(0.001), p_sat - p_H₂O);

    return g0 + P_AIR * FT(1e-6) / p_CO₂ * (1 + g1/sqrt(vpd)) * β * An[ind]
end

ext: jl


text:
@reexport baremodule LibP4est
  using CBinding: 𝐣𝐥

  # Introduce standard integer types
  const size_t  = 𝐣𝐥.Csize_t
  const ssize_t = 𝐣𝐥.Cssize_t
  const int8_t  = 𝐣𝐥.Int8
  const int16_t = 𝐣𝐥.Int16
  const int32_t = 𝐣𝐥.Int32
  const int64_t = 𝐣𝐥.Int64

  # Forward declare standard library types
  𝐣𝐥.@cstruct FILE
  𝐣𝐥.@cstruct va_list

  𝐣𝐥.Base.include((𝐣𝐥.@__MODULE__), 𝐣𝐥.joinpath(𝐣𝐥.dirname(𝐣𝐥.@__DIR__), "deps", "libp4est.jl"))
end


text:
# syntax: proto3
using ProtoBuf
import ProtoBuf.meta
import ._ProtoBuf_Top_.matrix_io

struct __enum_EnumMalosAction <: ProtoEnum
    MALOS_ACTION_NOT_DEFINED::Int32
    START_DETECTION::Int32
    STOP_DETECTION::Int32
    RESTART_DETECTION::Int32
    CAPTURE_PREVIEW::Int32
    START_STREAMING::Int32
    STOP_STREAMING::Int32
    __enum_EnumMalosAction() = new(0,1,2,3,4,5,6)
end #struct __enum_EnumMalosAction
const EnumMalosAction = __enum_EnumMalosAction()

struct __enum_EnumMalosEyeDetectionType <: ProtoEnum
    MALOS_DETECTION_TYPE_NOT_DEFINED::Int32
    STOP::Int32
    FACE::Int32
    FACE_DEMOGRAPHICS::Int32
    FACE_DESCRIPTOR::Int32
    HAND_THUMB_UP::Int32
    HAND_PALM::Int32
    HAND_PINCH::Int32
    HAND_FIST::Int32
    __enum_EnumMalosEyeDetectionType() = new(0,1,20,21,30,40,41,42,43)
end #struct __enum_EnumMalosEyeDetectionType
const EnumMalosEyeDetectionType = __enum_EnumMalosEyeDetectionType()

mutable struct CameraConfig <: ProtoType
    source::AbstractString
    camera_id::Int32
    width::Int32
    height::Int32
    CameraConfig(; kwargs...) = (o=new(); fillunset(o); isempty(kwargs) || ProtoBuf._protobuild(o, kwargs); o)
end #mutable struct CameraConfig
const __fnum_CameraConfig = Int[4,1,2,3]
meta(t::Type{CameraConfig}) = meta(t, ProtoBuf.DEF_REQ, __fnum_CameraConfig, ProtoBuf.DEF_VAL, true, ProtoBuf.DEF_PACK, ProtoBuf.DEF_WTYPES, ProtoBuf.DEF_ONEOFS, ProtoBuf.DEF_ONEOF_NAMES, ProtoBuf.DEF_FIELD_TYPES)

mutable struct DetectionServerConfig <: ProtoType
    detection_server_address::AbstractString
    detection_server_timeout::Int32
    DetectionServerConfig(; kwargs...) = (o=new(); fillunset(o); isempty(kwargs) || ProtoBuf._protobuild(o, kwargs); o)
end #mutable struct DetectionServerConfig

mutable struct FaceConfig <: ProtoType
    models_path::AbstractString
    send_slow::Bool
    send_text_proto::Bool
    FaceConfig(; kwargs...) = (o=new(); fillunset(o); isempty(kwargs) || ProtoBuf._protobuild(o, kwargs); o)
end #mutable struct FaceConfig

mutable struct MalosEyeConfig <: ProtoType
    camera_config::CameraConfig
    face_config::FaceConfig
    detection_server_config::DetectionServerConfig
    object_to_detect::Base.Vector{Int32}
    action::Int32
    models::Base.Vector{matrix_io.common.Entity}
    config::matrix_io.common.Entity
    MalosEyeConfig(; kwargs...) = (o=new(); fillunset(o); isempty(kwargs) || ProtoBuf._protobuild(o, kwargs); o)
end #mutable struct MalosEyeConfig
const __fnum_MalosEyeConfig = Int[1,2,4,21,5,6,7]
const __pack_MalosEyeConfig = Symbol[:object_to_detect]
meta(t::Type{MalosEyeConfig}) = meta(t, ProtoBuf.DEF_REQ, __fnum_MalosEyeConfig, ProtoBuf.DEF_VAL, true, __pack_MalosEyeConfig, ProtoBuf.DEF_WTYPES, ProtoBuf.DEF_ONEOFS, ProtoBuf.DEF_ONEOF_NAMES, ProtoBuf.DEF_FIELD_TYPES)

export EnumMalosAction, EnumMalosEyeDetectionType, MalosEyeConfig, CameraConfig, DetectionServerConfig, FaceConfig


text:
using NaNMath
using Test

@test isnan(NaNMath.log(-10))
@test isnan(NaNMath.log1p(-100))
@test isnan(NaNMath.pow(-1.5,2.3))
@test isnan(NaNMath.pow(-1.5f0,2.3f0))
@test isnan(NaNMath.pow(-1.5,2.3f0))
@test isnan(NaNMath.pow(-1.5f0,2.3))
@test NaNMath.pow(-1.5f0,2.3f0) isa Float32
@test NaNMath.pow(-1.5f0,2.3) isa Float64
@test NaNMath.pow(-1.5,2.3f0) isa Float64
@test NaNMath.pow(-1.5,2.3) isa Float64
@test isnan(NaNMath.sqrt(-5))
@test NaNMath.sqrt(5) == Base.sqrt(5)
@test NaNMath.sum([1., 2., NaN]) == 3.0
@test NaNMath.sum([1. 2.; NaN 1.]) == 4.0
@test isnan(NaNMath.sum([NaN, NaN]))
@test NaNMath.sum(Float64[]) == 0.0
@test NaNMath.sum([1f0, 2f0, NaN32]) === 3.0f0
@test NaNMath.maximum([1., 2., NaN]) == 2.0
@test NaNMath.maximum([1. 2.; NaN 1.]) == 2.0
@test NaNMath.minimum([1., 2., NaN]) == 1.0
@test NaNMath.minimum([1. 2.; NaN 1.]) == 1.0
@test NaNMath.extrema([1., 2., NaN]) == (1.0, 2.0)
@test NaNMath.extrema([2., 1., NaN]) == (1.0, 2.0)
@test NaNMath.extrema([1. 2.; NaN 1.]) == (1.0, 2.0)
@test NaNMath.extrema([2. 1.; 1. NaN]) == (1.0, 2.0)
@test NaNMath.extrema([NaN, -1., NaN]) == (-1.0, -1.0)
@test NaNMath.mean([1., 2., NaN]) == 1.5
@test NaNMath.mean([1. 2.; NaN 3.]) == 2.0
@test NaNMath.var([1., 2., NaN]) == 0.5
@test NaNMath.std([1., 2., NaN]) == 0.7071067811865476

@test NaNMath.median([1.]) == 1.
@test NaNMath.median([1., NaN]) == 1.
@test NaNMath.median([NaN, 1., 3.]) == 2.
@test NaNMath.median([1., 3., 2., NaN]) == 2.
@test NaNMath.median([NaN, 1, 3]) == 2.
@test NaNMath.median([1, 2, NaN]) == 1.5
@test NaNMath.median([1 2; NaN NaN]) == 1.5
@test NaNMath.median([NaN 2; 1 NaN]) == 1.5
@test isnan(NaNMath.median(Float64[]))
@test isnan(NaNMath.median(Float32[]))
@test isnan(NaNMath.median([NaN]))

@test NaNMath.min(1, 2) == 1
@test NaNMath.min(1.0, 2.0) == 1.0
@test NaNMath.min(1, 2.0) == 1.0
@test NaNMath.min(BigFloat(1.0), 2.0) == BigFloat(1.0)
@test NaNMath.min(BigFloat(1.0), BigFloat(2.0)) == BigFloat(1.0)
@test NaNMath.min(NaN, 1) == 1.0
@test NaNMath.min(NaN32, 1) == 1.0f0
@test isnan(NaNMath.min(NaN, NaN))
@test isnan(NaNMath.min(NaN))
@test NaNMath.min(NaN, NaN, 0.0, 1.0) == 0.0

@test NaNMath.max(1, 2) == 2
@test NaNMath.max(1.0, 2.0) == 2.0
@test NaNMath.max(1, 2.0) == 2.0
@test NaNMath.max(BigFloat(1.0), 2.0) == BigFloat(2.0)
@test NaNMath.max(BigFloat(1.0), BigFloat(2.0)) == BigFloat(2.0)
@test NaNMath.max(NaN, 1) == 1.0
@test NaNMath.max(NaN32, 1) == 1.0f0
@test isnan(NaNMath.max(NaN, NaN))
@test isnan(NaNMath.max(NaN))
@test NaNMath.max(NaN, NaN, 0.0, 1.0) == 1.0


text:
using PyCall
using Pandas

sklearn = pyimport("sklearn")
pd = pyimport("pandas")
lsk = pyimport("lale.lib.sklearn")
llale = pyimport("lale.lib.lale")
Project = llale.Project
OneHotEncoder = lsk.OneHotEncoder
ConcatFeatures = llale.ConcatFeatures
Tree = lsk.DecisionTreeClassifier
KNN = lsk.KNeighborsClassifier
RF =  lsk.RandomForestClassifier
PCA = lsk.PCA
Hyperopt = llale.Hyperopt

train_test_split = pyimport("sklearn.model_selection").train_test_split
mean_squared_error = pyimport("sklearn.metrics").mean_squared_error
metrics = pyimport("sklearn.metrics")
accuracy_scorer = metrics.make_scorer(sklearn.metrics.accuracy_score)

df = Pandas.read_csv("credit.csv")
r,c = size(df)
y = df["class"].pyo
X = iloc(df)[1:r,2:c]

train_X, test_X, train_y, test_y = train_test_split(X, y)

prep_to_numbers = (
  (Project(columns=Dict("type"=>"string")) >> OneHotEncoder(handle_unknown = "ignore")) &
  Project(columns=Dict("type"=>"number")) )>> ConcatFeatures
planned_orig = prep_to_numbers >> ( Tree | KNN)

best_estimator = planned_orig.auto_configure(
    train_X, train_y, optimizer=Hyperopt, cv=3, max_evals=10)
println("accuracy = $(accuracy_scorer(best_estimator, test_X, test_y))")

pca_tree_planned = prep_to_numbers >> PCA >> RF
pca_tree_trained = pca_tree_planned.auto_configure(
    train_X, train_y, optimizer=Hyperopt, cv=3, max_evals=10, verbose=true)
println("accuracy = $(accuracy_scorer(pca_tree_trained, test_X, test_y))")


text:
function surface(mtime)

  global SW,LW,LH,SH,E,TAU,Qnet,MOL
  ########################################################################
  # surface
  ########################################################################

  #roughness length (stensrud, table 2.3)
  #Ice 10^-4
  #Grass (mown) 10^-2
  #Long grass, rocky ground 0.05
  #Pasture 0.20
  #Suburban housing 0.6
  #Forests, cities 1–5

  #Wikipedia
  # Mud flats, snow; no vegetation, no obstacles	0.005
  zm=0.005
  zh=0.005
  ze=0.005

  CD=karman^2/log(ZAC[1]/zm)^2
  CH=sqrt(CD)*karman/log(ZAC[1]/zh)
  CE=sqrt(CD)*karman/log(ZAC[1]/ze)

  #kinematic surface flux
  ustar=-sqrt(CD)*abs(UA[1])*UA[1]
  tstar=-sqrt(CH)*abs(UA[1])*(ΘA[1]-TS[end])
  qstar=-sqrt(CE)*abs(UA[1])*(qA[1]-q_v[end])

  L=0;

  for i=1:5
    wθv=tstar*(1+humid_fac*qA[1])+humid_fac*ΘA[1]*qstar
    L=-ΘA[1]*ustar^3/(karman*gravity_mks*wθv)

    #CD=karman^2/(log(ZAC[1]/zm)-holtslag_psim(ZAC[1]/L))
    #CH=CD/(log(ZAC[1]/zh)-holtslag_psis(ZAC[1]/L))
    #CE=CD/(log(ZAC[1]/ze)-holtslag_psis(ZAC[1]/L))

    CD=CD/(1+CD*(log(ZAC[1]/zm)-holtslag_psim(ZAC[1]/L)/karman))
    CH=CH/(1+CH*(log(ZAC[1]/zm)-holtslag_psis(ZAC[1]/L)/karman))
    CE=CE/(1+CE*(log(ZAC[1]/zm)-holtslag_psis(ZAC[1]/L)/karman))

    ustar=-sqrt(CD)*abs(UA[1])*UA[1]
    tstar=-sqrt(CH)*abs(UA[1])*(ΘA[1]-TS[end])
    qstar=-sqrt(CE)*abs(UA[1])*(qA[1]-q_v[end])

  end

  SH=-ρA[1]*CH*cp*abs(UA[1])*(ΘA[1]-TS[end])
  LH=-ρA[1]*CE*Av*abs(UA[1])*(qA[1]-q_v[end])
  TAU=ρA[1]*CD*abs(UA[1])*UA[1];

  MOL=L;
  E=-LH/Av;

  Qnet=SW-LH-SH-LW;

end


text:
import Base: (+), (-), (*), (//), (/), (==), hash

(+)(a::Padic{P}, b::Padic{P}) where {P} = Padic(a.x + b.x)

(-)(a::Padic{P}) where {P} = Padic(-(a.x))

(-)(a::Padic{P}, b::Padic{P}) where {P} = a + (-b)

(*)(a::Padic{P}, b::Padic{P}) where {P} = Padic(a.x * b.x)

(/)(a::Padic{P}, b::Padic{P}) where {P} = Padic(a.x // b.x)

(==)(a::Padic{P}, b::Padic{Q}) where {P,Q} = (P == Q) && (a.x == b.x)

import Base: promote_rule

promote_rule(::Type{Padic{P}}, ::Type{T}) where {P,T<:Integer} = Padic{P}
promote_rule(::Type{Padic{P}}, ::Type{T}) where {P,T<:Rational} = Padic{P}

(//)(a::Padic, b::Number) = a / b
(//)(a::Number, b::Padic) = a / b
(//)(a::Padic, b::Padic) = a / b

hash(a::Padic{P}, h::UInt64) where {P} = hash(a.x, h)
hash(a::Padic{P}) where {P} = hash(a.x)


text:
using Documenter, ImageEdgeDetection

format = Documenter.HTML(edit_link = "master",
                         prettyurls = get(ENV, "CI", nothing) == "true",
                         assets = String[])

makedocs(;
    modules=[ImageEdgeDetection],
    format = format,
    pages=[
        "Home" => "index.md",
        "Function Reference" => "reference.md"
    ],
    sitename="ImageEdgeDetection.jl",
    authors="Dr. Zygmunt L. Szpak"
)

deploydocs(;
    repo="github.com/JuliaImages/ImageEdgeDetection.jl.git",
)


text:
#processes the V scripts output by the models and turns particles into bootstrapped confidence intervals
#needs the ground_truth_V.csv file from ideal_V.jl

using MetaGen
using JSON
using Pipe: @pipe

include("helper_bootstrap_V.jl")

ppath = "../../scratch_work_07_16_21/09_18/shuffle_0/"

num_videos = 50
num_particles = 100
n_possible_objects = 5

online_data = CSV.read(path * "online_V.csv", DataFrame; delim = "&")
#retro_data = CSV.read(path * "retro_V.csv", DataFrame; delim = "&")
#lesioned_data = CSV.read(path * "lesioned_V.csv", DataFrame; delim = "&")
ground_truth_data = CSV.read(path * "ground_truth_V.csv", DataFrame; delim = "&")

#make CI on the average value for each entry in the V matrix
online_ci_averages = confidence_interval(online_data, num_particles, n_possible_objects, 1000, 0.95)

#get MSE for each particle
online_ci_MSE = MSE_and_confidence_interval(online_data, ground_truth_data, num_particles, n_possible_objects, 1000, 0.95)

#merge them all together
online_V_processed = hcat(select(online_data, 1:2*n_possible_objects+1), online_ci_averages, online_ci_MSE, select(online_data, length(names(online_data))))

CSV.write(path * "online_V_processed.csv", online_V_processed)

max_stars_repo_name: zhangir-azerbayev/MetaGen

max_issues_repo_name: zhangir-azerbayev/MetaGen

max_forks_repo_name: zhangir-azerbayev/MetaGen


text:
#=
Copyright (c) 2015, Intel Corporation

Redistribution and use in source and binary forms, with or without
modification, are permitted provided that the following conditions are met:

    * Redistributions of source code must retain the above copyright notice,
      this list of conditions and the following disclaimer.
    * Redistributions in binary form must reproduce the above copyright
      notice, this list of conditions and the following disclaimer in the
      documentation and/or other materials provided with the distribution.
    * Neither the name of Intel Corporation nor the names of its contributors
      may be used to endorse or promote products derived from this software
      without specific prior written permission.

THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE
FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY,
OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
=#

export drop_fixed_dims
"""
Drop indexing expressions along dimension for an array if the compiler
has determined it to be uniform.  `cbdata` should be a Dict{Symbol,
Vector{Bool}} that maps array names to their uniformity.
"""
function drop_fixed_dims(ast, arg_dim_info)
    function walker(node, cbdata::Dict{Symbol, Vector{Bool}}, index, top_level, read)
        if is_ref(node) && haskey(cbdata, node.args[1])
            info = cbdata[node.args[1]]
            count = 1
            # First argument is reference name so we skip it
            idx = Any[node.args[1]]
            for arg in node.args[2:end]
                if isa(arg, Symbol) && contains(string(arg), "_neuron_index")
                    if count > length(info) || !info[count]
                        push!(idx, arg)
                    end
                    count += 1
                else
                    push!(idx, arg)
                end
            end
            node.args = idx
            return node
        end
        ASTWALK_RECURSE
    end
    AstWalk(ast, walker, arg_dim_info)
end


text:
# TODO : créer un module
#
# module Log
# export level, lg
# export lg0, lg1, lg2, lg3, lg4, lg5
# export ln0, ln1, ln2, ln3, ln4, ln5
# global level = 2
# ...
# end # module

# Quelques fonctions d'aide au debog (EN COURS DE MISE AU POINT)
#
# La principale fonction est lg (pour log) qui affiche ou non un message selon le
# niveau de verbosité courant avec ou sans saut de ligne selon la valeur du
# paramètre suffix.
# Dans tous les cas, cette fonction retourne un booléen indiquant si le niveau
# de verbosité est suffisant.
#
# Le message à afficher est formé par la concaténation des arguments
# vals transformés en String (par la méthode string(...)).
#
# Ce fichier dépend du module Args (gestion des arguments du programme)
# en particulier de la variable args[:loglevel] définisant le niveau de
# verbosité souhaité.
#
# Exemple d'utilisation :
#
#   ln3("coucou")
#   => affiche  "coucou" + prefix si loglevel suffisant (ici si >= 3)
#   lg4(".")
#   => affiche un caractère "." sans saut de ligne pour illustrer une progression
#   if lg3()
#     # pré-calcul... # gros pré-calcul juste pour l'affichage qui suit
#     println("Résultat du pré-calcul")
#   end
#   lg3() n'affiche rien par lui-même, mais autorise une série d'instruction si
#         level est suffisant
#   S'écrit également :
#   lg3() && @show mavariable
#
# BUG ET DÉPENDENCES
# Pour l'instant, ces méthodes ne sont utilisables que après avoir appelé
# la méthode Args.parse_commandline() pour accéder à
# Ceci sera à modifier (création d'un module Log paramétré de l'extérieur)
#

# lg0(vals...) = lg(0, vals..., suffix="", doflush=true  )
lg0(vals...; kwargs... ) = lg(0, vals...; kwargs... )
lg1(vals...; kwargs... ) = lg(1, vals...; kwargs... )
lg2(vals...; kwargs... ) = lg(2, vals...; kwargs... )
lg3(vals...; kwargs... ) = lg(3, vals...; kwargs... )
lg4(vals...; kwargs... ) = lg(4, vals...; kwargs... )
lg5(vals...; kwargs... ) = lg(5, vals...; kwargs... )

ln0(vals...; kwargs... ) = lg(0, vals...; kwargs..., suffix="\n" )
ln1(vals...; kwargs... ) = lg(1, vals...; kwargs..., suffix="\n" )
ln2(vals...; kwargs... ) = lg(2, vals...; kwargs..., suffix="\n" )
ln3(vals...; kwargs... ) = lg(3, vals...; kwargs..., suffix="\n" )
ln4(vals...; kwargs... ) = lg(4, vals...; kwargs..., suffix="\n" )
ln5(vals...; kwargs... ) = lg(5, vals...; kwargs..., suffix="\n" )

function lg(minlevel::Int, vals... ;
            prefix::String="", suffix::String="", doflush=true)
    if Args.args[:level] >= minlevel
        if !isempty(vals)
            print(prefix, join(vals, ""), suffix)
            doflush && flush(stdout)
        end
        return true
    else
        return false
    end
end

# POUR INFORMATION : Á DÉPLACER
#
# Voir aussi (capture de l'affichage d'une commande) :
#  https://stackoverflow.com/questions/54599148/in-julia-1-0-how-do-i-get-strings-using-redirect-stdout
#
# Utilisation :
#   redirect_to_files(prefix * ".log", prefix * ".err") do
#       compute(...)
#   end
#
# function redirect_to_files(dofunc, outfile, errfile)
#     open(outfile, "w") do out
#         open(errfile, "w") do err
#             redirect_stdout(out) do
#                 redirect_stderr(err) do
#                     dofunc()
#                 end
#             end
#         end
#     end
# end


text:
struct CPUEvent <: Event
    task::Core.Task
end

function Event(::CPU)
    return NoneEvent()
end

wait(ev::Union{CPUEvent, NoneEvent, MultiEvent}, progress=nothing) = wait(CPU(), ev, progress)
wait(::CPU, ev::NoneEvent, progress=nothing) = nothing

function wait(cpu::CPU, ev::MultiEvent, progress=nothing)
    dependencies = collect(ev.events)
    cpudeps   = filter(d->d isa CPUEvent && d.task !== nothing, dependencies)
    otherdeps = filter(d->!(d isa CPUEvent), dependencies)
    Base.sync_end(map(e->e.task, cpudeps))
    for event in otherdeps
        wait(CPU(), event, progress)
    end
end

function wait(::CPU, ev::CPUEvent, progress=nothing)
    ev.task === nothing && return

    if progress === nothing
        wait(ev.task)
    else
        while !Base.istaskdone(ev.task)
            progress()
        end
    end
end

function async_copy!(::CPU, A, B; dependencies=nothing)
    wait(CPU(), MultiEvent(dependencies), yield)
    copyto!(A, B)
    return CPUEvent(nothing)
end

function (obj::Kernel{CPU})(args...; ndrange=nothing, workgroupsize=nothing, dependencies=nothing)
    if ndrange isa Integer
        ndrange = (ndrange,)
    end
    if workgroupsize isa Integer
        workgroupsize = (workgroupsize, )
    end
    if dependencies isa Event
        dependencies = (dependencies,)
    end

    if KernelAbstractions.workgroupsize(obj) <: DynamicSize && workgroupsize === nothing
        workgroupsize = (1024,) # Vectorization, 4x unrolling, minimal grain size
    end
    iterspace, dynamic = partition(obj, ndrange, workgroupsize)
    # partition checked that the ndrange's agreed
    if KernelAbstractions.ndrange(obj) <: StaticSize
        ndrange = nothing
    end

    t = Threads.@spawn __run(obj, ndrange, iterspace, args, dependencies, Val(dynamic))
    return CPUEvent(t)
end

# Inference barriers
function __run(obj, ndrange, iterspace, args, dependencies, ::Val{dynamic}) where dynamic
    wait(CPU(), MultiEvent(dependencies), yield)
    N = length(iterspace)
    Nthreads = Threads.nthreads()
    if Nthreads == 1
        len, rem = N, 0
    else
        len, rem = divrem(N, Nthreads)
    end
    # not enough iterations for all the threads?
    if len == 0
        Nthreads = N
        len, rem = 1, 0
    end
    if Nthreads == 1
        __thread_run(1, len, rem, obj, ndrange, iterspace, args, Val(dynamic))
    else
        @sync for tid in 1:Nthreads
            Threads.@spawn __thread_run(tid, len, rem, obj, ndrange, iterspace, args, Val(dynamic))
        end
    end
    return nothing
end

function __thread_run(tid, len, rem, obj, ndrange, iterspace, args, ::Val{dynamic}) where dynamic
    # compute this thread's iterations
    f = 1 + ((tid-1) * len)
    l = f + len - 1
    # distribute remaining iterations evenly
    if rem > 0
        if tid <= rem
            f = f + (tid-1)
            l = l + tid
        else
            f = f + rem
            l = l + rem
        end
    end
    # run this thread's iterations
    for i = f:l
        block = @inbounds blocks(iterspace)[i]
        ctx = mkcontext(obj, block, ndrange, iterspace, Val(dynamic))
        Cassette.overdub(ctx, obj.f, args...)
    end
    return nothing
end

Cassette.@context CPUCtx

function mkcontext(kernel::Kernel{CPU}, I, _ndrange, iterspace, ::Val{dynamic}) where dynamic
    metadata = CompilerMetadata{ndrange(kernel), dynamic}(I, _ndrange, iterspace)
    Cassette.disablehooks(CPUCtx(pass = CompilerPass, metadata=metadata))
end

@inline function Cassette.overdub(ctx::CPUCtx, ::typeof(__index_Local_Linear), idx::CartesianIndex)
    indices = workitems(__iterspace(ctx.metadata))
    return @inbounds LinearIndices(indices)[idx]
end

@inline function Cassette.overdub(ctx::CPUCtx, ::typeof(__index_Group_Linear), idx::CartesianIndex)
    indices = blocks(__iterspace(ctx.metadata))
    return @inbounds LinearIndices(indices)[__groupindex(ctx.metadata)]
end

@inline function Cassette.overdub(ctx::CPUCtx, ::typeof(__index_Global_Linear), idx::CartesianIndex)
    I = @inbounds expand(__iterspace(ctx.metadata), __groupindex(ctx.metadata), idx)
    @inbounds LinearIndices(__ndrange(ctx.metadata))[I]
end

@inline function Cassette.overdub(ctx::CPUCtx, ::typeof(__index_Local_Cartesian), idx::CartesianIndex)
    return idx
end

@inline function Cassette.overdub(ctx::CPUCtx, ::typeof(__index_Group_Cartesian), idx::CartesianIndex)
    __groupindex(ctx.metadata)
end

@inline function Cassette.overdub(ctx::CPUCtx, ::typeof(__index_Global_Cartesian), idx::CartesianIndex)
    return @inbounds expand(__iterspace(ctx.metadata), __groupindex(ctx.metadata), idx)
end

@inline function Cassette.overdub(ctx::CPUCtx, ::typeof(__validindex), idx::CartesianIndex)
    # Turns this into a noop for code where we can turn of checkbounds of
    if __dynamic_checkbounds(ctx.metadata)
        I = @inbounds expand(__iterspace(ctx.metadata), __groupindex(ctx.metadata), idx)
        return I in __ndrange(ctx.metadata)
    else
        return true
    end
end

@inline function Cassette.overdub(ctx::CPUCtx, ::typeof(__print), items...)
    __print(items...)
end

generate_overdubs(CPUCtx)

# Don't recurse into these functions
const cpufuns = (:cos, :cospi, :sin, :sinpi, :tan,
          :acos, :asin, :atan,
          :cosh, :sinh, :tanh,
          :acosh, :asinh, :atanh,
          :log, :log10, :log1p, :log2,
          :exp, :exp2, :exp10, :expm1, :ldexp,
          :isfinite, :isinf, :isnan, :signbit,
          :abs,
          :sqrt, :cbrt,
          :ceil, :floor,)
for f in cpufuns
    @eval function Cassette.overdub(ctx::CPUCtx, ::typeof(Base.$f), x::Union{Float32, Float64})
        @Base._inline_meta
        return Base.$f(x)
    end
end

###
# CPU implementation of shared memory
###
@inline function Cassette.overdub(ctx::CPUCtx, ::typeof(SharedMemory), ::Type{T}, ::Val{Dims}, ::Val) where {T, Dims}
    MArray{__size(Dims), T}(undef)
end

###
# CPU implementation of scratch memory
# - private memory for each workitem
# - memory allocated as a MArray with size `Dims + WorkgroupSize`
###
struct ScratchArray{N, D}
    data::D
    ScratchArray{N}(data::D) where {N, D} = new{N, D}(data)
end

@inline function Cassette.overdub(ctx::CPUCtx, ::typeof(Scratchpad), ::Type{T}, ::Val{Dims}) where {T, Dims}
    return ScratchArray{length(Dims)}(MArray{__size((Dims..., __groupsize(ctx.metadata)...)), T}(undef))
end

Base.@propagate_inbounds function Base.getindex(A::ScratchArray, I...)
    return A.data[I...]
end

Base.@propagate_inbounds function Base.setindex!(A::ScratchArray, val, I...)
    A.data[I...] = val
end


text:
using Test, Random, LSHFunctions

include(joinpath("..", "utils.jl"))

#==================
Tests
==================#
@testset "MinHash tests" begin
    Random.seed!(RANDOM_SEED)

    @testset "Can construct a MinHash hash function" begin
        symbol_collections = [
            collect(1:100),
            Set(1:30),
            ["a", "a", "b", "c"]
        ]
        n_hashes = [1, 8, 4]

        for (symbols, nh) in zip(symbol_collections, n_hashes)
            hashfn = MinHash(nh; symbols=symbols)

            @test isa(hashfn, MinHash{eltype(symbols)})
            @test LSHFunctions.n_hashes(hashfn) == nh
            @test hashtype(hashfn) == UInt32
        end
    end

    @testset "MinHash range is 1:length(symbols)" begin
        symbols = collect(1:10)
        hashfn = MinHash(1000; symbols = symbols)

        # Verify that all hashes computed by hashfn are in the appropriate range
        has_correct_output_types = true
        has_correct_output_ranges = true

        for ii = 1:100
            if !(has_correct_output_types && has_correct_output_ranges)
                break
            end

            n = mod(rand(Int32), length(symbols)) + 1
            x = shuffle(symbols)[1:n]
            hashes = hashfn(x)

            has_correct_output_types  &= isa(hashes, Vector{UInt32})
            has_correct_output_ranges &= all(1 .≤ hashes .≤ 10 - (n - 1))
        end

        @test has_correct_output_types
        @test has_correct_output_ranges
    end

    @testset "Hashes are correctly computed" begin
        symbols = collect(50:100)
        hashfn = MinHash(10; symbols = symbols)

        dataset = shuffle(symbols)[1:10]
        hashes = hashfn(dataset)

        @test eltype(hashes) == hashtype(hashfn)

        @test let hashes_match = true
            for (hash, mapping) in zip(hashes, hashfn.mappings)
                if !hashes_match
                    break
                end

                # Compute MinHash manually
                expected_hash = minimum(mapping[x] for x in dataset)
                hashes_match &= (hash == expected_hash)
            end
            hashes_match
        end
    end

    @testset "Collision probabilities correlated with Jaccard similarity" begin
        # Create three different datasets, such that the first and second datasets
        # have high Jaccard similarity while the first and third have low Jaccard
        # similarity. Test whether we in fact get higher collision rates for the
        # two sets with high similarity than the two sets with low similarity.
        symbols = collect(1:200)
        n_hashes = 100
        hashfn = MinHash(n_hashes; symbols=symbols)

        shuffle!(symbols)
        dataset_1 = symbols[1:100]
        dataset_2 = [symbols[1:75]; symbols[101:125]]
        dataset_3 = [symbols[1:25]; symbols[101:175]]

        hashes_1 = hashfn(dataset_1)
        hashes_2 = hashfn(dataset_2)
        hashes_3 = hashfn(dataset_3)

        @test sum(hashes_1 .== hashes_2) > sum(hashes_1 .== hashes_3)
    end

    @testset "MinHash observed collision frequencies match probabilities" begin
        # In theory, the probability of collision for two datasets should be
        # roughly equal to the Jaccard similarity between those datasets
        symbols = collect(1:200)
        n_hashes = 10_000
        hashfn = MinHash(n_hashes; symbols=symbols)

        shuffle!(symbols)
        dataset_1 = Set(symbols[1:100])
        dataset_2 = Set([symbols[1:75]; symbols[101:125]])
        dataset_3 = Set([symbols[1:25]; symbols[101:175]])

        hashes_1 = hashfn(dataset_1)
        hashes_2 = hashfn(dataset_2)
        hashes_3 = hashfn(dataset_3)

        sim_12 = similarity(hashfn)(Set(dataset_1), Set(dataset_2))
        sim_13 = similarity(hashfn)(Set(dataset_1), Set(dataset_3))
        prob_12 = collision_probability(hashfn, sim_12; n_hashes=1)
        prob_13 = collision_probability(hashfn, sim_13; n_hashes=1)

        mean(x) = sum(x) / length(x)

        @test abs(mean(hashes_1 .== hashes_2) - prob_12) ≤ 0.01
        @test abs(mean(hashes_1 .== hashes_3) - prob_13) ≤ 0.01
    end

    @testset "Can omit symbol set to lazily update hash functions" begin
        # If the users want, they can decide not to pre-specify the symbol
        # set. In that case, the hash functions should be lazily updated.
        hashfn = MinHash(100)

        @test isa(hashfn, MinHash{Any, UInt64})

        symbols = [collect(1:10); "a"; "b"; "c"]
        dataset = shuffle(symbols)[1:3]
        hashes = hashfn(dataset)

        @test isa(hashes, Vector{UInt64})

        # Should still be able to hash new elements, updating the hash
        # functions further.
        new_dataset = ["d"; "e"; "f"]
        new_hashes = hashfn(new_dataset)

        @test sum(new_hashes .== hashes) == 0
    end

    @testset "Achieve collision probabilities with lazily updated hash functions" begin
        symbols = collect(1:1_000)
        dataset_1 = Set(shuffle(symbols)[1:500])
        dataset_2 = Set(shuffle(symbols)[1:500])
        similarity = jaccard(dataset_1, dataset_2)

        n_hashes = 10_000
        hashfn = MinHash(n_hashes; dtype=Int64)
        hashes_1 = hashfn(dataset_1)
        hashes_2 = hashfn(dataset_2)

        collision_proportion = sum(hashes_1 .== hashes_2) / n_hashes
        @test abs(collision_proportion - similarity) ≤ 0.01
    end
end


text:
using MLDataPattern: slidingwindow
import StatsPlots: plot

"""
Calculates the windows used for modeling, the target index, and features useful for plotting.
"""
mutable struct WindowGenerator
    train                       # training windows
    valid                       # validation windows
    h                           # historical steps in window
    f                           # target steps in window
    label_indices               # indices for plotting predictions in a window
    target_idx::Array{Int, 1}   # indices to be predicted
end

"""
Specify `h` historical points, `f` target points. By default, the target points are assumed to follow after all the historical points (`offset = h`).
By setting `offset = 1`, the targets for each historical point will be the next point in time.
"""
function WindowGenerator(h, f, train_df, valid_df, label_columns::Vector{String}; offset=h)
    train = slidingwindow(i -> i+offset:i+offset+f-1, Array(train_df)', h, stride=1)
    valid = slidingwindow(i -> i+offset:i+offset+f-1, Array(valid_df)', h, stride=1)

    label_indices = (offset + 1):(offset + f)

    target_idx = findall(x->x in label_columns, names(train_df))

    return WindowGenerator(train, valid, h, f, label_indices, target_idx)
end

# WindowGenerator with a single label_column.
WindowGenerator(h, f, train_df, valid_df, label_columns::String; offset=h) =
        WindowGenerator(h, f, train_df, valid_df, label_columns=[label_columns]; offset=offset)

# WindowGenerator with multiple label_columns.
WindowGenerator(h, f, train_df, valid_df; label_columns, offset=h) =
        WindowGenerator(h, f, train_df, valid_df, label_columns; offset=offset)

"""
Plots the historical data and the target(s) for prediction.
"""
function plot(wg::WindowGenerator)
    idx = wg.target_idx[1] # Currently, will only plot the first target index but could be redone to plot all label columns.
    sw = wg.valid
    plots = []
    for it in 1:3
        i = rand(1:size(sw,1))
        p = plot(sw[i][1][idx, :], leg=false)
        scatter!(wg.label_indices,sw[i][2][idx,:]')
        it == 3 && xlabel!("Given $(wg.h)h as input, predict $(wg.f)h into the future.")
        push!(plots, p)
    end
    plot(plots..., layout=(3,1))
end

"""
Plots the historical window, the target(s) for prediction, and the model's predictions.
"""
function plot(wg::WindowGenerator, model; set=:valid)
    data = getfield(wg,set)
    i = rand(1:size(data,1))
    plot(1:wg.h, data[i][1][2,:], lab="Inputs", shape=:circ, m=2, leg=:outerright)
    scatter!(wg.label_indices, data[i][2][2,:], lab="Labels", c=:green)
    z = model(unsqueeze(data[i][1],3))[:]
    scatter!(wg.label_indices, z, lab="Predictions", shape=:star5, m=6, c=:orange)
end

max_stars_repo_name: adinhobl/model-zoo

max_issues_repo_name: adinhobl/model-zoo

max_forks_repo_name: adinhobl/model-zoo


text:
#==============================================================================#
# AWSCore/test/runtests.jl
#
# Copyright OC Technology Pty Ltd 2014 - All rights reserved
#==============================================================================#

using Test
using Dates
using AWSCore
using SymDict
using Retry
using XMLDict
using HTTP

using AWSCore: service_query

AWSCore.set_debug_level(1)

@testset "AWSCore" begin

aws = aws_config()

@testset "NoAuth" begin
    pub_request1 = Dict{Symbol, Any}(
        :service => "s3",
        :headers => Dict{String, String}("Range" => "bytes=0-0"),
        :content => "",
        :resource => "/invenia-static-website-content/invenia_ca/index.html",
        :url => "https://s3.us-east-1.amazonaws.com/invenia-static-website-content/invenia_ca/index.html",
        :verb => "GET",
        :region => "us-east-1",
        :creds => nothing,
    )
    pub_request2 = Dict{Symbol, Any}(
        :service => "s3",
        :headers => Dict{String, String}("Range" => "bytes=0-0"),
        :content => "",
        :resource => "ryft-public-sample-data/AWS-x86-AMI-queries.json",
        :url => "https://s3.amazonaws.com/ryft-public-sample-data/AWS-x86-AMI-queries.json",
        :verb => "GET",
        :region => "us-east-1",
        :creds => nothing,
    )
    response = nothing
    try
        response = AWSCore.do_request(pub_request1)
    catch e
        println(e)
        @test ecode(e) in ["AccessDenied", "NoSuchEntity"]
        try
            response = AWSCore.do_request(pub_request2)
        catch e
            println(e)
            @test ecode(e) in ["AccessDenied", "NoSuchEntity"]
        end
    end
    @test response == "<" || response == UInt8['[']
end
@testset "Load Credentials" begin
    user = aws_user_arn(aws)

    @test occursin(r"^arn:aws:iam::[0-9]+:[^:]+$", user)

    println("Authenticated as: $user")

    aws[:region] = "us-east-1"

    println("Testing exceptions...")
    try
        AWSCore.Services.iam("GetFoo", Dict("ContentType" => "JSON"))
        @test false
    catch e
        println(e)
        @test ecode(e) == "InvalidAction"
    end

    try
        AWSCore.Services.iam("GetUser", Dict("UserName" => "notauser",
                                           "ContentType" => "JSON"))
        @test false
    catch e
        println(e)
        @test ecode(e) in ["AccessDenied", "NoSuchEntity"]
    end

    try
        AWSCore.Services.iam("GetUser", Dict("UserName" => "@#!%%!",
                                           "ContentType" => "JSON"))
        @test false
    catch e
        println(e)
        @test ecode(e) == "ValidationError"
    end

    try
        AWSCore.Services.iam("CreateUser", Dict("UserName" => "root",
                                              "ContentType" => "JSON"))
        @test false
    catch e
        println(e)
        @test ecode(e) in ["AccessDenied", "EntityAlreadyExists"]
    end
end

@testset "AWSCredentials" begin
    @testset "Defaults" begin
        creds = AWSCredentials("access_key_id" ,"secret_key")
        @test creds.token == ""
        @test creds.user_arn == ""
        @test creds.account_number == ""
        @test creds.expiry == typemax(DateTime)
        @test creds.renew == nothing
    end

    @testset "Renewal" begin
        # Credentials shouldn't throw an error if no renew function is supplied
        creds = AWSCredentials("access_key_id", "secret_key", renew = nothing)
        newcreds = check_credentials(creds, force_refresh = true)
        # Creds should remain unchanged if no renew function exists
        @test creds === newcreds
        @test creds.access_key_id == "access_key_id"
        @test creds.secret_key == "secret_key"
        @test creds.renew == nothing

        # Creds should error if the renew function returns nothing
        creds = AWSCredentials("access_key_id", "secret_key", renew = () -> nothing)
        @test_throws ErrorException check_credentials(creds, force_refresh = true)
        # Creds should remain unchanged
        @test creds.access_key_id == "access_key_id"
        @test creds.secret_key == "secret_key"

        # Creds should take on value of a returned AWSCredentials except renew function
        function gen_credentials()
            i = 0
            () -> (i += 1; AWSCredentials("NEW_ID_$i", "NEW_KEY_$i"))
        end

        creds = AWSCredentials(
            "access_key_id",
            "secret_key",
            renew = gen_credentials(),
            expiry = now(UTC),
        )

        @test creds.renew !== nothing
        renewed = creds.renew()

        @test creds.access_key_id == "access_key_id"
        @test creds.secret_key == "secret_key"
        @test creds.expiry <= now(UTC)
        @test AWSCore.will_expire(creds)

        @test renewed.access_key_id === "NEW_ID_1"
        @test renewed.secret_key == "NEW_KEY_1"
        @test renewed.renew === nothing
        @test renewed.expiry == typemax(DateTime)
        @test !AWSCore.will_expire(renewed)
        renew = creds.renew

        # Check renewal on time out
        newcreds = check_credentials(creds, force_refresh = false)
        @test creds === newcreds
        @test creds.access_key_id == "NEW_ID_2"
        @test creds.secret_key == "NEW_KEY_2"
        @test creds.renew !== nothing
        @test creds.renew === renew
        @test creds.expiry == typemax(DateTime)
        @test !AWSCore.will_expire(creds)

        # Check renewal doesn't happen if not forced or timed out
        newcreds = check_credentials(creds, force_refresh = false)
        @test creds === newcreds
        @test creds.access_key_id == "NEW_ID_2"
        @test creds.secret_key == "NEW_KEY_2"
        @test creds.renew !== nothing
        @test creds.renew === renew
        @test creds.expiry == typemax(DateTime)

        # Check forced renewal works
        newcreds = check_credentials(creds, force_refresh = true)
        @test creds === newcreds
        @test creds.access_key_id == "NEW_ID_3"
        @test creds.secret_key == "NEW_KEY_3"
        @test creds.renew !== nothing
        @test creds.renew === renew
        @test creds.expiry == typemax(DateTime)
    end

    mktemp() do config_file, config_io
        write(config_io, """[profile test]
                output = json
                region = us-east-1

                [profile test:dev]
                source_profile = test
                role_arn = arn:aws:iam::123456789000:role/Dev

                [profile test:sub-dev]
                source_profile = test:dev
                role_arn = arn:aws:iam::123456789000:role/SubDev

                [profile test2]
                aws_access_key_id = WRONG_ACCESS_ID
                aws_secret_access_key = WRONG_ACCESS_KEY
                output = json
                region = us-east-1

                [profile test3]
                source_profile = test:dev
                role_arn = arn:aws:iam::123456789000:role/test3

                [profile test4]
                source_profile = test:dev
                role_arn = arn:aws:iam::123456789000:role/test3
                aws_access_key_id = RIGHT_ACCESS_ID4
                aws_secret_access_key = RIGHT_ACCESS_KEY4
            """)
        close(config_io)

        mktemp() do creds_file, creds_io
            write(creds_io, """[test]
                aws_access_key_id = TEST_ACCESS_ID
                aws_secret_access_key = TEST_ACCESS_KEY

                [test2]
                aws_access_key_id = RIGHT_ACCESS_ID2
                aws_secret_access_key = RIGHT_ACCESS_KEY2

                [test3]
                aws_access_key_id = RIGHT_ACCESS_ID3
                aws_secret_access_key = RIGHT_ACCESS_KEY3
                """)
            close(creds_io)

            withenv(
                "AWS_SHARED_CREDENTIALS_FILE" => creds_file,
                "AWS_CONFIG_FILE" => config_file,
                "AWS_DEFAULT_PROFILE" => "test",
                "AWS_ACCESS_KEY_ID" => nothing
                ) do

                @testset "Loading" begin
                    # Check credentials load
                    config = aws_config()
                    creds = config[:creds]
                    @test creds isa AWSCredentials

                    @test creds.access_key_id == "TEST_ACCESS_ID"
                    @test creds.secret_key == "TEST_ACCESS_KEY"
                    @test creds.renew !== nothing

                    # Check credential file takes precedence over config
                    ENV["AWS_DEFAULT_PROFILE"] = "test2"
                    config = aws_config()
                    creds = config[:creds]

                    @test creds.access_key_id == "RIGHT_ACCESS_ID2"
                    @test creds.secret_key == "RIGHT_ACCESS_KEY2"

                    # Check credentials take precedence over role
                    ENV["AWS_DEFAULT_PROFILE"] = "test3"
                    config = aws_config()
                    creds = config[:creds]

                    @test creds.access_key_id == "RIGHT_ACCESS_ID3"
                    @test creds.secret_key == "RIGHT_ACCESS_KEY3"

                    ENV["AWS_DEFAULT_PROFILE"] = "test4"
                    config = aws_config()
                    creds = config[:creds]

                    @test creds.access_key_id == "RIGHT_ACCESS_ID4"
                    @test creds.secret_key == "RIGHT_ACCESS_KEY4"
                end

                @testset "Refresh" begin
                    ENV["AWS_DEFAULT_PROFILE"] = "test"
                    # Check credentials refresh on timeout
                    config = aws_config()
                    creds = config[:creds]
                    creds.access_key_id = "EXPIRED_ACCESS_ID"
                    creds.secret_key = "EXPIRED_ACCESS_KEY"
                    creds.expiry = now(UTC)

                    @test creds.renew !== nothing
                    renew = creds.renew
                    @test renew() isa AWSCredentials

                    creds = check_credentials(config[:creds])

                    @test creds.access_key_id == "TEST_ACCESS_ID"
                    @test creds.secret_key == "TEST_ACCESS_KEY"
                    @test creds.expiry > now(UTC)

                    # Check renew function remains unchanged
                    @test creds.renew !== nothing
                    @test creds.renew === renew

                    # Check force_refresh
                    creds.access_key_id = "WRONG_ACCESS_KEY"
                    creds = check_credentials(creds, force_refresh = true)
                    @test creds.access_key_id == "TEST_ACCESS_ID"
                end

                @testset "Profile" begin
                    # Check profile kwarg
                    ENV["AWS_DEFAULT_PROFILE"] = "test"
                    creds = AWSCredentials(profile="test2")
                    @test creds.access_key_id == "RIGHT_ACCESS_ID2"
                    @test creds.secret_key == "RIGHT_ACCESS_KEY2"

                    config = aws_config(profile="test2")
                    creds = config[:creds]
                    @test creds.access_key_id == "RIGHT_ACCESS_ID2"
                    @test creds.secret_key == "RIGHT_ACCESS_KEY2"

                    # Check profile persists on renewal
                    creds.access_key_id = "WRONG_ACCESS_ID2"
                    creds.secret_key = "WRONG_ACCESS_KEY2"
                    creds = check_credentials(creds, force_refresh=true)

                    @test creds.access_key_id == "RIGHT_ACCESS_ID2"
                    @test creds.secret_key == "RIGHT_ACCESS_KEY2"
                end

                @testset "Assume Role" begin
                    # Check we try to assume a role
                    ENV["AWS_DEFAULT_PROFILE"] = "test:dev"

                    try
                        aws_config()
                        @test false
                    catch e
                        @test e isa AWSException
                        @test ecode(e) == "InvalidClientTokenId"
                    end

                    # Check we try to assume a role
                    ENV["AWS_DEFAULT_PROFILE"] = "test:sub-dev"
                    let oldout = stdout
                        r,w = redirect_stdout()
                        try
                            aws_config()
                            @test false
                        catch e
                            @test e isa AWSException
                            @test ecode(e) == "InvalidClientTokenId"
                        end
                        redirect_stdout(oldout)
                        close(w)
                        output = String(read(r))
                        occursin("Assuming \"test:dev\"", output)
                        occursin("Assuming \"test\"", output)
                        close(r)
                    end
                end
            end
        end
    end
end

@testset "XML Parsing" begin
    XML(x)=parse_xml(x)

    xml = """
    <CreateQueueResponse>
        <CreateQueueResult>
            <QueueUrl>
                http://queue.amazonaws.com/123456789012/testQueue
            </QueueUrl>
        </CreateQueueResult>
        <ResponseMetadata>
            <RequestId>
                7a62c49f-347e-4fc4-9331-6e8e7a96aa73
            </RequestId>
        </ResponseMetadata>
    </CreateQueueResponse>
    """

    @assert XML(xml)["CreateQueueResult"]["QueueUrl"] ==
          "http://queue.amazonaws.com/123456789012/testQueue"

    xml = """
    <GetUserResponse xmlns="https://iam.amazonaws.com/doc/2010-05-08/">
      <GetUserResult>
        <User>
          <PasswordLastUsed>2015-12-23T22:45:36Z</PasswordLastUsed>
          <Arn>arn:aws:iam::012541411202:root</Arn>
          <UserId>012541411202</UserId>
          <CreateDate>2015-09-15T01:07:23Z</CreateDate>
        </User>
      </GetUserResult>
      <ResponseMetadata>
        <RequestId>837446c9-abaf-11e5-9f63-65ae4344bd73</RequestId>
      </ResponseMetadata>
    </GetUserResponse>
    """

    @test XML(xml)["GetUserResult"]["User"]["Arn"] == "arn:aws:iam::012541411202:root"

    xml = """
    <GetQueueAttributesResponse>
      <GetQueueAttributesResult>
        <Attribute>
          <Name>ReceiveMessageWaitTimeSeconds</Name>
          <Value>2</Value>
        </Attribute>
        <Attribute>
          <Name>VisibilityTimeout</Name>
          <Value>30</Value>
        </Attribute>
        <Attribute>
          <Name>ApproximateNumberOfMessages</Name>
          <Value>0</Value>
        </Attribute>
        <Attribute>
          <Name>ApproximateNumberOfMessagesNotVisible</Name>
          <Value>0</Value>
        </Attribute>
        <Attribute>
          <Name>CreatedTimestamp</Name>
          <Value>1286771522</Value>
        </Attribute>
        <Attribute>
          <Name>LastModifiedTimestamp</Name>
          <Value>1286771522</Value>
        </Attribute>
        <Attribute>
          <Name>QueueArn</Name>
          <Value>arn:aws:sqs:us-east-1:123456789012:qfoo</Value>
        </Attribute>
        <Attribute>
          <Name>MaximumMessageSize</Name>
          <Value>8192</Value>
        </Attribute>
        <Attribute>
          <Name>MessageRetentionPeriod</Name>
          <Value>345600</Value>
        </Attribute>
      </GetQueueAttributesResult>
      <ResponseMetadata>
        <RequestId>1ea71be5-b5a2-4f9d-b85a-945d8d08cd0b</RequestId>
      </ResponseMetadata>
    </GetQueueAttributesResponse>
    """

    d = Dict(a["Name"] => a["Value"] for a in XML(xml)["GetQueueAttributesResult"]["Attribute"])

    @test d["MessageRetentionPeriod"] == "345600"
    @test d["CreatedTimestamp"] == "1286771522"

    xml = """
    <?xml version="1.0" encoding="UTF-8"?>
    <ListAllMyBucketsResult xmlns="http://s3.amazonaws.com/doc/2006-03-01">
      <Owner>
        <ID>bcaf1ffd86f461ca5fb16fd081034f</ID>
        <DisplayName>webfile</DisplayName>
      </Owner>
      <Buckets>
        <Bucket>
          <Name>quotes</Name>
          <CreationDate>2006-02-03T16:45:09.000Z</CreationDate>
        </Bucket>
        <Bucket>
          <Name>samples</Name>
          <CreationDate>2006-02-03T16:41:58.000Z</CreationDate>
        </Bucket>
      </Buckets>
    </ListAllMyBucketsResult>
    """

    @test map(b->b["Name"], XML(xml)["Buckets"]["Bucket"]) == ["quotes", "samples"]

    xml = """
    <ListDomainsResponse>
      <ListDomainsResult>
        <DomainName>Domain1</DomainName>
        <DomainName>Domain2</DomainName>
        <NextToken>TWV0ZXJpbmdUZXN0RG9tYWluMS0yMDA3MDYwMTE2NTY=</NextToken>
      </ListDomainsResult>
      <ResponseMetadata>
        <RequestId>eb13162f-1b95-4511-8b12-489b86acfd28</RequestId>
        <BoxUsage>0.0000219907</BoxUsage>
      </ResponseMetadata>
    </ListDomainsResponse>
    """

    @test XML(xml)["ListDomainsResult"]["DomainName"] == ["Domain1", "Domain2"]
end

@testset "AWS Signature Version 4" begin
    function aws4_request_headers_test()

        r = @SymDict(
            creds         = AWSCredentials(
                                "AKIDEXAMPLE",
                                "wJalrXUtnFEMI/K7MDENG+bPxRfiCYEXAMPLEKEY"
                            ),
            region        = "us-east-1",
            verb          = "POST",
            service       = "iam",
            url           = "http://iam.amazonaws.com/",
            content       = "Action=ListUsers&Version=2010-05-08",
            headers       = Dict(
                                "Content-Type" =>
                                "application/x-www-form-urlencoded; charset=utf-8",
                                "Host" => "iam.amazonaws.com"
                            )
        )

        AWSCore.sign!(r, DateTime("2011-09-09T23:36:00"))

        h = r[:headers]
        out = join(["$k: $(h[k])\n" for k in sort(collect(keys(h)))])

        expected = (
            "Authorization: AWS4-HMAC-SHA256 " *
            "Credential=AKIDEXAMPLE/20110909/us-east-1/iam/aws4_request, " *
            "SignedHeaders=content-md5;content-type;host;" *
            "x-amz-content-sha256;x-amz-date, " *
            "Signature=1a6db936024345449ef4507f890c5161" *
                     "bbfa2ff2490866653bb8b58b7ba1554a\n" *
            "Content-MD5: r2d9jRneykOuUqFWSFXKCg==\n" *
            "Content-Type: application/x-www-form-urlencoded; " *
                         "charset=utf-8\n" *
            "Host: iam.amazonaws.com\n" *
            "x-amz-content-sha256: b6359072c78d70ebee1e81adcbab4f01" *
                                 "bf2c23245fa365ef83fe8f1f955085e2\n" *
            "x-amz-date: 20110909T233600Z\n")

        @test out == expected
    end

    aws4_request_headers_test()
end

@testset "ARN" begin
    @test arn(aws,"s3","foo/bar") == "arn:aws:s3:::foo/bar"
    @test arn(aws,"s3","foo")     == "arn:aws:s3:::foo"
    @test arn(aws,"sqs", "au-test-queue", "ap-southeast-2", "1234") ==
          "arn:aws:sqs:ap-southeast-2:1234:au-test-queue"

    @test arn(aws,"sns","*","*",1234) == "arn:aws:sns:*:1234:*"
    @test arn(aws,"iam","role/foo-role", "", 1234) ==
          "arn:aws:iam::1234:role/foo-role"
end

@testset "Misc" begin
    @test HTTP.escapepath("invocations/function:f:PROD") ==
                          "invocations/function%3Af%3APROD"
end

@testset "Exception" begin
    code = "InvalidSignatureException"
    message = "Signature expired: ..."
    body = """
        {
           "__type": "$code",
           "message": "$message"
        }
        """
    resp = HTTP.Messages.Response(400, body)
    HTTP.setheader(resp.headers, "Content-Type" => "application/x-amz-json-1.1")

    ex = AWSException(HTTP.StatusError(400, resp))

    @test ex.code == code
    @test ex.message == message
end

instance_type = get(ENV, "AWSCORE_INSTANCE_TYPE", "")
if instance_type == "EC2"
    @testset "EC2" begin
        @test_nowarn AWSCore.ec2_metadata("instance-id")
        @test startswith(AWSCore.ec2_metadata("instance-id"), "i-")

        @test AWSCore.localhost_maybe_ec2()
        @test AWSCore.localhost_is_ec2()
        @test_nowarn AWSCore.ec2_instance_credentials()
        ec2_creds = AWSCore.ec2_instance_credentials()
        @test ec2_creds !== nothing

        default_creds = AWSCredentials()
        @test default_creds.access_key_id == ec2_creds.access_key_id
        @test default_creds.secret_key == ec2_creds.secret_key
    end
elseif instance_type == "ECS"
    @testset "ECS" begin
        @test_nowarn AWSCore.ecs_instance_credentials()
        ecs_creds = AWSCore.ecs_instance_credentials()
        @test ecs_creds !== nothing

        default_creds = AWSCredentials()
        @test default_creds.access_key_id == ecs_creds.access_key_id
        @test default_creds.secret_key == ecs_creds.secret_key
    end
end

end # testset "AWSCore"

#==============================================================================#
# End of file.
#==============================================================================#


text:
export DomainTrans

"""
Transform the domain of a `ParametricFunction`.

Evaluation at `x` with parameters `θ` is equivalent to
`evaluate(inner, θ, transformation(x))`.
"""
immutable DomainTrans{T, S} <: ParametricFamily
    transformation::T
    inner::S
end

function DomainTrans(dom::AbstractInterval, fam::ParametricFamily)
    DomainTrans(bridge(dom, domain(fam)), fam)
end

function DomainTrans(dom::AbstractInterval, transformation::UnivariateTransformation,
                     fam::ParametricFamily)
    DomainTrans(bridge(dom, transformation, domain(fam)), fam)
end

@forward DomainTrans.inner degf

domain(p::DomainTrans) = domain_in_image(p.transformation, domain(p.inner))

points(p::DomainTrans) = inv(p.transformation).(points(p.inner))

function basis!{T}(p::DomainTrans, x::T, b::AbstractVector{T})
    basis!(p.inner, p.transformation(x), b)
end

evaluate(p::DomainTrans, θ, x) = evaluate(p.inner, θ, p.transformation(x))

function evaluate(p::DomainTrans, θ, x::Partial)
    y, dy = p.transformation(x.x, DERIV)
    evaluate(p.inner, θ, Partial(y))*dy
end

function evaluate(p::DomainTrans, θ, x::ValuePartial)
    y, dy = p.transformation(x.x, DERIV)
    z, dz = evaluate(p.inner, θ, ValuePartial(y))
    z, dy*dz
end

fit!(p::DomainTrans, ys::AbstractVector, θ) = fit!(p.inner, ys, θ)


text:
module Surrogates

using LinearAlgebra
using DataFrames
using GLM
using Distributions
using Sobol
using LatinHypercubeSampling
using XGBoost
using LIBSVM
using Flux

abstract type AbstractSurrogate <: Function end
include("Radials.jl")
include("Kriging.jl")
include("Sampling.jl")
include("Optimization.jl")
include("Lobachesky.jl")
include("LinearSurrogate.jl")
include("RandomForestSurrogate.jl")
include("SVMSurrogate.jl")
include("NeuralSurrogate.jl")

export AbstractSurrogate, SamplingAlgorithm
export Kriging, RadialBasis, add_point!, current_estimate, std_error_at_point
export sample, GridSample, UniformSample, SobolSample, LatinHypercubeSample, LowDiscrepancySample
export SRBF,LCBS,EI,DYCORS,surrogate_optimize
export LobacheskySurrogate, lobachesky_integral, lobachesky_integrate_dimension
export LinearSurrogate
export RandomForestSurrogate
export SVMSurrogate
export NeuralSurrogate

end


text:
module LinearSolvers

	abstract type AbstractSolver end
	abstract type AbstractDirectSolver <: AbstractSolver end
	function factorLinearSystem! end
	export AbstractSolver, AbstractDirectSolver, factorLinearSystem!

	using KrylovMethods
	using Distributed
	using SparseArrays
	using LinearAlgebra
	using Pkg
	# check if ParSPMatVec is available
	global hasParSpMatVec = false
	try
		using ParSpMatVec;
		global hasParSpMatVec = ParSpMatVec.isBuilt();
	catch
	end

	export solveLinearSystem!,solveLinearSystem

	solveLinearSystem(A,B,param::AbstractSolver,doTranspose::Int=0) = solveLinearSystem!(A,B,zeros(eltype(B),size(B)),param,doTranspose)

	import Distributed.clear!
	function clear!(M::AbstractSolver)
		M.Ainv = []
	end

	include("iterativeWrapper.jl")
	include("blockIterativeWrapper.jl")
	include("juliaWrapper.jl")

	export clear!

end # module LinearSolvers


text:
import Base: cat, vcat, hcat

# cat(A...; dims): Concatenate the input arrays along the specified
# dimensions in the iterable dims. For each dimension not in dims, all
# input arrays should have the same size, which will also be the size
# of the output array along that dimension. For each dimension in
# dims, the size of the output array is the sum of the sizes of the
# input arrays along that dimension. If dims is a single number, the
# different arrays are tightly stacked along that dimension. If dims
# is an iterable containing several dimensions, this allows one to
# construct block diagonal matrices and their higher-dimensional
# analogues by simultaneously increasing several dimensions for every
# new input array and putting zero blocks elsewhere. For example,
# cat(matrices...; dims=(1,2)) builds a block diagonal matrix, i.e. a
# block matrix with matrices[1], matrices[2], ... as diagonal blocks
# and matching zero blocks away from the diagonal.

# cat can take 0 or more arguments of any type.  In order to catch the
# cases where at least one arg is Value, we need to override the generic
# Base.cat(d,x...). This will call the recording cat_r if a Value is
# found, and the original _cat from the Base definition if not.

const NA = Union{Number,AbstractArray}
const NAR = Union{Number,AbstractArray,Value}

# Base has cat(x...; dims) defined, first specialize this:
cat(X::NA...; dims)=Base._cat(dims, X...)

# Then define the method that catches at least one Value:
cat(X::NAR...; dims)=forw(cat,X...; dims=dims)
back(::typeof(cat),::Type{Arg{N}},y1,y,x...; dims) where {N}=uncat(y1,N,dims,x...)

# In Julia6+ dims can be Val{N} which breaks uncat:
uncat(y1,n,dims::Val{N},x...) where {N}=uncat(y1,n,N,x...)
uncat1(x2,y1,n,dims::Val{N},x...) where {N}=uncat1(x2,y1,n,N,x...)
# This resolves ambiguity with the @primitive version:
uncat(y1::Value,n,dims::Val{N},x...) where {N}=uncat(y1,n,N,x...)
uncat1(x2::Value,y1,n,dims::Val{N},x...) where {N}=uncat1(x2,y1,n,N,x...)

# For the gradient, we need to extract the n'th block from dy which
# has the same shape as y=cat(x...;dims).  Note that the inputs x[i]
# may have fewer dimensions than dy, or even be scalars.  In those
# cases Julia assumes the missing dimensions are 1.  We need to
# reshape dx to the same size as x.

function uncat(y1,n,dims,x...)
    idx = []
    @inbounds for d=1:ndims(y1)
        if in(d,dims)
            pos = 0
            @inbounds for j=1:n-1
                pos += size(x[j],d)
            end
            push!(idx,(pos+1):(pos+size(x[n],d)))
        else
            push!(idx,1:size(y1,d))
        end
    end
    x1 = y1[idx...]
    if isa(value(x[n]),Number)
        length(x1)==1 || error("uncat mismatch")
        x1 = x1[1]
    else
        x1 = reshape(x1, size(x[n]))
    end
    return x1
end

function uncat1(x2,y1,n,dims,x...)
    idx = []
    @inbounds for d=1:ndims(y1)
        if in(d,dims)
            pos = 0
            @inbounds for j=1:n-1
                pos += size(x[j],d)
            end
            push!(idx,(pos+1):(pos+size(x[n],d)))
        else
            push!(idx,1:size(y1,d))
        end
    end
    y2 = zero(y1)
    y2[idx...] = x2
    return y2
end

@primitive  uncat(y1,n,dims,x...),x2  uncat1(x2,y1,n,dims,x...)
@primitive  uncat1(x2,y1,n,dims,x...),y3  uncat(y3,n,dims,x...)

# Here is a graphic that may explain the variable name choice where xi
# stands for the i'th order gradient:
#
# x  → cat    → y
#               ↓
# x1 ← uncat  ← y1
# ↓
# x2 → uncat1 → y2
#               ↓
# x3 ← uncat  ← y3

# cat1(x...)=cat(x...; dims=Val(1))
# cat2(x...)=cat(x...; dims=Val(2))
# addtestN(:cat1, 1., 2.)
# addtestN(:cat1, 1., [2.,3.])
# addtestN(:cat1, [1.,2.], 3.)
# addtestN(:cat1, [1.,2.], [3.,4.])
# addtestN(:cat1, [1. 2.], [3. 4.])
# addtestN(:cat2, 1., 2.)
# addtestN(:cat2, 1., [2. 3.])
# addtestN(:cat2, [1. 2.], 3.)
# addtestN(:cat2, [1.,2.], [3.,4.])
# addtestN(:cat2, [1. 2.], [3. 4.])

# vcat,hcat: should be defined in terms of cat. However base has some
# generic methods that prevent the cat call when all arguments are
# boxed.  This should fix it:

vcat(x::Value...) = cat(x...; dims=Val(1))
hcat(x::Value...) = cat(x...; dims=Val(2))

# vcat(a::Value,b::Value,c...)=cat(1,a,b,c...)
# vcat(a,b::Value,c...)=cat(1,a,b,c...)
# vcat(a::Value,b...)=cat(1,a,b...)
# hcat(a::Value,b::Value,c...)=cat(2,a,b,c...)
# hcat(a,b::Value,c...)=cat(2,a,b,c...)
# hcat(a::Value,b...)=cat(2,a,b...)

# The cat implementation above is slow when called with a lot of
# arguments (100s).  In our parser, the first epoch is a lot slower
# than the other epochs.  This is probably because there is a lot of
# compilation for different argument type patterns (where Values
# appear).  The following version tries to avoid all compilation at
# runtime.  It takes a bunch of arrays (not numbers), could be
# different shapes and sizes, and concatenates them as if they were
# all vectors.  It thus avoids a lot of shape related calculation as
# well.  It works with KnetArrays but slow, probably due to the copy
# kernel call overhead.  To make it fast with KnetArrays we need a
# single GPU kernel call which does all the copying.

"""
    cat1d(args...)

Return `vcat(vec.(args)...)` but possibly more efficiently. Can be used to concatenate the
contents of arrays with different shapes and sizes.
"""
function cat1d(args...)
    @inbounds for arg in args
        if isa(arg,Value)
            return forw(_cat1d, args...)
        end
    end
    return _cat1d(args...)
end

function _cat1d(args...)
    totlen = 0
    @inbounds for arg in args
        totlen += length(arg)
    end
    result = similar(args[1], totlen)
    offset1 = 1
    offset2 = 0
    @inbounds for arg in args
        offset2 += length(arg)
        result[offset1:offset2] = arg
        offset1 = offset2+1
    end
    return result
end

function back(::typeof(_cat1d),::Type{Arg{N}},y1,y,x...) where {N}
    offset2 = 0
    @inbounds for i=1:N
        offset2 += length(x[i])
    end
    offset1 = offset2 - length(x[N]) + 1
    reshape(y1[offset1:offset2], size(x[N]))
end

export cat1d

# addtestN(:cat1d, [1.,2.], [3.,4.])
# addtestN(:cat1d, [1. 2.], [3. 4.])
# addtestN(:cat1d, [1.,2.], [3. 4.])
# addtestN(:cat1d, [1. 2.], [3.,4.])


text:
const ALPHABET = [Char(x) for x in vcat(collect.([42:43, 48:57, 65:90, 97:122])...)]
const INV_ALPHABET = Dict(c => i for (i,c) in enumerate(ALPHABET))

_segment_width(l::Integer) = ceil(Int, log2(l+1))
encode(i::Integer, l::Integer) = string(i, base=2, pad=_segment_width(l))
function decode(c::AbstractString, l::Integer)
    k = max(1, min(_segment_width(l), length(c)))
    parse(Int, c[1:k], base=2), c[k+1:end]
end

function stringify(c::AbstractString)
    if length(c) % 6 != 0
        c = c * '0' ^ mod(6-length(c), 6)
    end
    join(ALPHABET[parse(Int, x, base=2) + 1] for x in [c[i:i+5] for i in 1:6:(length(c)-1)])
end

function destringify(c::AbstractString)
    join(string(INV_ALPHABET[x] - 1, base=2, pad=6) for x in c)
end

walk(n::T, c) where T = _walk(n, destringify(c))

_walk(n::T, c) where T = _walk(NodeType(T), n, c)

function _walk(::LeafNode, n, c::AbstractString)
    !isempty(c) || return n
    i, nc = decode(c, length(c))
    if i == 0 && Set(nc) ⊆ ['0']
        return n
    else
        error("Invalid index!")
    end
end

function _walk(::InnerNode, n, c::AbstractString)
    !isempty(c) || return n
    i, nc = decode(c, nprintchildren(n))
    0 <= i <= nprintchildren(n) || error("Invalid index!")
    if i == 0
        if Set(nc) ⊆ ['0']
            return n
        else
            error("Invalid index!")
        end
    end
    _walk(_ith_child(printchildren(n), i), nc)
end

encode_traversal(t, idxs::Integer...) = stringify(_encode_traversal(t, idxs...))

function _encode_traversal(t, idxs...)
    !isempty(idxs) || return ""
    n = _ith_child(printchildren(t), idxs[1])
    return encode(idxs[1], nchildren(t)) * _encode_traversal(n, idxs[2:end]...)
end

pred_traversal(n::T, p::Function, s::String="") where T = _pred_traversal(NodeType(T), n, p, s)

_pred_traversal(::LeafNode, n, p, s="") = p(n) ? [stringify(s)] : String[]
function _pred_traversal(::InnerNode, n, p, s="")
    d = printchildren(n)
    l = length(d)
    z = Vector{String}[pred_traversal(_ith_child(d, i), p, s * encode(i, l)) for i in 1:l]
    res = isempty(z) ? String[] : reduce(vcat, z)
    p(n) ? vcat(stringify(s), res) : res
end

list_traversal(n) = pred_traversal(n, t -> true)
find_traversal(n, x) = pred_traversal(n, t -> x === t)


text:
module Day20
__revise_mode__ = :eval

using Advent: save_input, eachlinegroup, read_charmatrix

function read_input(filename)
    tiles = Dict{Int,Matrix{Char}}()
    open(filename) do io
        buf = []
        num = nothing
        line = nothing
        while !eof(io)
            line = readline(io)
            line == "" && (tiles[num] = mapreduce(collect, hcat, buf); num = nothing; empty!(buf); continue)
            num == nothing && (num = parse(Int, line[6:9]); continue)
            push!(buf, line)
        end
    end
    @info tiles
    tiles
end

function get_borders(mat)
    Set((mat[1:end, 1], mat[1, 1:end], mat[end, 1:end], mat[1:end, end]))
end

function solve_a(filename)
    tiles = read_input(filename)
    borders = Dict(k => get_borders(v) for (k, v) in tiles)
    counts = Dict(
        k => count(tile) do t
            any(borders) do b
                (t in b.second || reverse(t) in b.second) && k != b.first
            end
        end
        for (k, tile) in borders
    )
    prod(findall(==(2), counts))
end

remove_borders(mat) = mat[2:end-1, 2:end-1]

function solve_b(filename)
    # not solved
    tiles = read_input(filename)
    borders = Dict(k => get_borders(v) for (k, v) in tiles)
    counts = Dict(
        k => count(tile) do t
            any(borders) do b
                (t in b.second || reverse(t) in b.second) && k != b.first
            end
        end
        for (k, tile) in borders
    )
    prod(findall(==(2), counts))
end

solve(filename = save_input(2020, 20)) = solve_a(filename), solve_b(filename)

end

max_stars_repo_name: AnAverageHuman/competitive

max_issues_repo_name: AnAverageHuman/competitive

max_forks_repo_name: AnAverageHuman/competitive


text:
## GETTING AND STTING ITERATION PARAMETER

"""
    rget(model, field::Union{Symbol,Expr})
    rget(mach, field::Union{Symbol,Expr})

Returns the value of the possibly nested parameter, `field`, from
`model` or `mach.model`.

"""
rget(model, field) =
    MLJBase.recursive_getproperty(model, field)
rget(mach::Machine, field) = rget(mach.model, field)

"""
    rset!(model, field::Union{Symbol,Expr}, i)
    rset!(mach, field::Union{Symbol,Expr}, i)

Change the value of the possibly nested parameter of `model` or `mach.model`
to `i`.

"""
function rset!(model, field, i)
    MLJBase.recursive_setproperty!(model, field, i)
    return nothing
end
rset!(mach::Machine, field, i) = rset!(mach.model, field, i)


text:
#########
# Types #
#########

const NULL_INDEX = typemin(Int)
const NULL_TAPE = InstructionTape()

# TrackedReal #
#-------------#

#=
A `TrackedReal` stores a value and a reference back to the original `TrackedArray` which
provided the value.

When performing a forward pass through a previously-recorded `InstructionTape`, any encountered
`TrackedReal` instances which are direct descendents of their origin array must re-validate
themselves by re-retrieving their value from their origin via the given index. A similar
strategy is taken during the reverse pass, when derivs are updated; they are always
re-validated with the origin before and after "local" updates to the TrackedReal. The big
benefit of this strategy is that scalar `getindex` operations don't need to be explicitly
recorded to the tape.

Note that we don't have to worry about the origin's values being invalidated during `InstructionTape`
execution, since the `TrackedArray` type is immutable.

Also note that it's possible to instantiate an origin-less `TrackedReal`. This will most
often happen as a result of conversion, or via a numerical constructor (e.g. `zero`).
In that case, we leave the `origin` field uninitialized, and set the `index` field to
`NULL_INDEX`.

There are two alternative approaches we could've taken for dealing with the potential
missing-ness of the origin, and it's worth explaining why we don't take these approaches.

The first alternative is to encode the lack of an origin as a parameter in the
`TrackedReal` type. Indeed, this is exactly what we do if the origin of a `TrackedReal` is a
constructor like `zero`, setting `O` to `Nothing`. However, this strategy is not a viable
solution for origin-less `TrackedReal`s that result from conversions. For example, take the
expression `y = convert(::Type{T<:TrackedReal}, x::Real)`. `y` will not have an origin, but
`isa(y, T)` MUST hold true regardless. This could break in a world where `T` forcibly
specifies an origin of type `O != Nothing`.

The second alternative is to use `Nullable` + `isnull`. This solution would work, but it
overcomplicates the API and would incur unneccesary pointer loads during re-validation,
making the below implementation preferable.
=#

mutable struct TrackedReal{V<:Real,D<:Real,O} <: Real
    value::V
    deriv::D
    tape::InstructionTape
    index::Int
    origin::O
    TrackedReal{V,D,O}(value, deriv, tape, index, origin) where {V,D,O} = new{V,D,O}(value, deriv, tape, index, origin)
    TrackedReal{V,D,O}(value, deriv, tape) where {V,D,O} = new{V,D,O}(value, deriv, tape, NULL_INDEX)
    TrackedReal{V,D,O}(value, deriv) where {V,D,O} = new{V,D,O}(value, deriv, NULL_TAPE, NULL_INDEX)
    TrackedReal{V,D,O}(value) where {V,D,O} = new{V,D,O}(value, zero(D), NULL_TAPE, NULL_INDEX)
end

TrackedReal(v::V, a::D, tp::InstructionTape, i::Int, o::O) where {V,D,O} = TrackedReal{V,D,O}(v, a, tp, i, o)

TrackedReal(v::V, a::D, tp::InstructionTape = NULL_TAPE) where {V,D} = TrackedReal{V,D,Nothing}(v, a, tp)

# we define these special cases so that the "constructor <--> convert" pun holds for `TrackedReal`
# this is Jarett's favorite piece of code. A true work of art.
@inline TrackedReal{V,D,O}(x::TrackedReal) where {V,D,O} = convert(TrackedReal{V,D,O}, x)

# TrackedArray #
#--------------#

struct TrackedArray{V,D,N,VA,DA} <: AbstractArray{TrackedReal{V,D,TrackedArray{V,D,N,VA,DA}},N}
    value::VA
    deriv::DA
    tape::InstructionTape
    function TrackedArray{V,D,N,VA,DA}(value::AbstractArray{V,N},
                                       deriv::AbstractArray{D,N},
                                       tape::InstructionTape) where {V,D,N,VA,DA}
        @assert IndexStyle(value) === IndexLinear()
        @assert size(value) === size(deriv)
        return new{V,D,N,VA,DA}(value, deriv, tape)
    end
end

function TrackedArray(value::AbstractArray{V,N},
                      deriv::AbstractArray{D,N},
                      tape::InstructionTape) where {V,D,N}
    return TrackedArray{V,D,N,typeof(value),typeof(deriv)}(value, deriv, tape)
end

const TrackedVector{V,D} = TrackedArray{V,D,1}
const TrackedMatrix{V,D} = TrackedArray{V,D,2}
const TrackedVecOrMat{V,D} = Union{TrackedVector{V,D}, TrackedMatrix{V,D}}

###########
# getters #
###########

istracked(x) = false
istracked(::TrackedReal) = true
istracked(::TrackedArray) = true
istracked(::AbstractArray{T}) where {T} = T <: TrackedReal || !(isconcretetype(T))

@inline value(x) = x
@inline value(x::AbstractArray) = istracked(x) ? map(value, x) : x
@inline value(t::TrackedReal) = t.value
@inline value(t::TrackedArray) = t.value

@inline deriv(t::TrackedArray) = t.deriv
@inline deriv(t::TrackedReal) =  t.deriv

@inline valtype(::TrackedReal{V}) where {V} = V
@inline valtype(::Type{TrackedReal{V,D,O}}) where {V,D,O} = V
@inline valtype(::TrackedArray{V}) where {V} = V
@inline valtype(::Type{TrackedArray{V,D,N,VA,DA}}) where {V,D,VA,DA,N} = V

@inline derivtype(::TrackedReal{V,D}) where {V,D} = D
@inline derivtype(::Type{TrackedReal{V,D,O}}) where {V,D,O} = D
@inline derivtype(t::TrackedArray{V,D}) where {V,D} = D
@inline derivtype(::Type{TrackedArray{V,D,N,VA,DA}}) where {V,D,VA,DA,N} = D

@inline origintype(::TrackedReal{V,D,O}) where {V,D,O} = O
@inline origintype(::Type{TrackedReal{V,D,O}}) where {V,D,O} = O

@inline hasorigin(x::Real) = false
@inline hasorigin(t::TrackedReal) = t.index !== NULL_INDEX

@inline hastape(x) = false
@inline hastape(t::TrackedArray) = tape(t) !== NULL_TAPE
@inline hastape(t::TrackedReal) = tape(t) !== NULL_TAPE
@inline hastape(x::AbstractArray) = tape(x) !== NULL_TAPE

@inline tape(x) = NULL_TAPE
@inline tape(t::TrackedArray) = t.tape
@inline tape(t::TrackedReal) = t.tape

function tape(x::AbstractArray)
    if istracked(x)
        for i in x
            hastape(i) && return tape(i)
        end
    end
    return NULL_TAPE
end

function tape(ts...)
    for t in ts
        hastape(t) && return tape(t)
    end
    return NULL_TAPE
end

###########
# setters #
###########

@inline value!(t::TrackedReal, v::Real) = (t.value = v; nothing)
@inline value!(t::TrackedArray, v::AbstractArray) = (copyto!(value(t), v); nothing)

function value!(t::NTuple{N,Any}, v::NTuple{N,Any}) where N
    for i in eachindex(t)
        value!(t[i], v[i])
    end
    return nothing
end

@inline deriv!(t::TrackedReal, v::Real) = (t.deriv = v; nothing)
@inline deriv!(t::TrackedArray, v::AbstractArray) = (copyto!(deriv(t), v); nothing)

function deriv!(t::NTuple{N,Any}, v::NTuple{N,Any}) where N
    for i in eachindex(t)
        deriv!(t[i], v[i])
    end
    return nothing
end

# pulling values from origin #
#----------------------------#

pull_value!(x) = nothing
pull_value!(t::TrackedArray) = nothing
pull_value!(t::TrackedReal) = (hasorigin(t) && value!(t, value(t.origin)[t.index]); nothing)
pull_value!(x::AbstractArray) = (istracked(x) && foreach(pull_value!, x); nothing)

# pulling derivs from origin #
#----------------------------#

pull_deriv!(x) = nothing
pull_deriv!(t::TrackedArray) = nothing
pull_deriv!(t::TrackedReal) = (hasorigin(t) && deriv!(t, deriv(t.origin)[t.index]); nothing)
pull_deriv!(x::AbstractArray) = (istracked(x) && foreach(pull_deriv!, x); nothing)

# pushing derivs back to origin #
#-------------------------------#

push_deriv!(x) = nothing
push_deriv!(t::TrackedArray) = nothing
push_deriv!(t::TrackedReal) = (hasorigin(t) && (t.origin.deriv[t.index] = deriv(t)); nothing)
push_deriv!(x::AbstractArray) = (istracked(x) && foreach(push_deriv!, x); nothing)

# seed/unseed #
#-------------#

seed!(x) = nothing
seed!(t::TrackedReal) = (t.deriv = one(derivtype(t)); push_deriv!(t); nothing)
seed!(t::TrackedArray, i) = (t.deriv[i] = one(derivtype(t)); nothing)
seed!(x::AbstractArray, i) = seed!(x[i])

unseed!(x) = nothing
unseed!(t::TrackedReal) = (t.deriv = zero(derivtype(t)); push_deriv!(t); nothing)
unseed!(t::TrackedArray) = (fill!(deriv(t), zero(derivtype(t))); nothing)
unseed!(x::AbstractArray) = (istracked(x) && foreach(unseed!, x); nothing)
unseed!(t::Tuple) = foreach(unseed!, t)
unseed!(t::TrackedArray, i) = (t.deriv[i] = zero(derivtype(t)); nothing)
unseed!(x::AbstractArray, i) = unseed!(x[i])

#########################
# capture (see tape.jl) #
#########################

# This is type unstable, but that shouldn't be too much of a problem as it's only used below
# the function barrier of `record!`, which just pushes `t` to the tape and always returns
# `nothing`. It enables us to not waste time updating constants during
# `forward_pass!`/`reverse_pass!`.
capture(t::TrackedReal) = ifelse(hastape(t), t, value(t))
capture(t::TrackedArray) = t
capture(t::AbstractArray) = istracked(t) ?  map!(capture, similar(t), t) : copy(t)

########################
# Conversion/Promotion #
########################

_convert(::Type{R}, t::TrackedReal) where {R <: Real} = R(value(t))
_convert(::Type{R}, t::Real) where {R <: TrackedReal} = R(t)
_convert(::Type{R}, t::Real) where {R <: Real} = R(t)
_convert(::Type{R}, t::TrackedReal) where {R <: TrackedReal} = R(t)

# recording a instruction for this preserves the line of references back to the origin's deriv
function Base.convert(::Type{T1}, t::T2) where {T1<:TrackedReal,T2<:TrackedReal}
    V1, D1, O1 = valtype(T1), derivtype(T1), origintype(T1)
    tp = tape(t)
    out = TrackedReal{V1,D1,O1}(_convert(V1, value(t)), _convert(D1, deriv(t)), tp)
    record!(tp, SpecialInstruction, convert, t, out)
    return out
end

@noinline function special_reverse_exec!(instruction::SpecialInstruction{typeof(convert)})
    output = instruction.output
    increment_deriv!(instruction.input, deriv(output))
    unseed!(output)
    return nothing
end

@noinline function special_forward_exec!(instruction::SpecialInstruction{typeof(convert)})
    input = instruction.input
    pull_value!(input)
    value!(instruction.output, value(input))
    return nothing
end

Base.convert(::Type{Real}, t::T) where {T<:TrackedReal} = t
function Base.convert(::Type{R}, t::T) where {R<:Real,T<:TrackedReal}
    throw(ArgumentError("Converting an instance of $T to $R is not defined. Please use `ReverseDiff.value` instead."))
end
function Base.convert(::Type{T}, x::R) where {T<:TrackedReal, R<:Real}
    return TrackedReal{valtype(T), derivtype(T), origintype(T)}(
        convert(valtype(T), value(x)),
    )
end

Base.convert(::Type{T}, t::T) where {T<:TrackedReal} = t
Base.convert(::Type{T}, t::T) where {T<:TrackedArray} = t

for R in REAL_TYPES
    @eval Base.promote_rule(::Type{$R}, ::Type{TrackedReal{V,D,O}}) where {V,D,O} = TrackedReal{promote_type($R,V),D,O}
end

Base.promote_rule(::Type{R}, ::Type{TrackedReal{V,D,O}}) where {R<:Real,V,D,O} = TrackedReal{promote_type(R,V),D,O}
Base.promote_rule(::Type{TrackedReal{V1,D1,O1}}, ::Type{TrackedReal{V2,D2,O2}}) where {V1,V2,D1,D2,O1,O2} = TrackedReal{promote_type(V1,V2),promote_type(D1,D2),Nothing}

###########################
# AbstractArray Interface #
###########################

colon2range(s, i) = i
colon2range(s, ::Colon) = s

function index_iterable(shape::NTuple{N,Any}, i::NTuple{M,Any}) where {N,M}
    if N < M
        return index_iterable(shape, ntuple(n -> i[n], Val(N)))
    elseif M < N && isa(last(i), Colon)
        return index_iterable(shape, ntuple(n -> (n > M ? Colon() : i[n]), Val(N)))
    else
        return Base.Iterators.product(map(colon2range, shape[1:M], i)...)
    end
end

for T in (:AbstractRange, :Colon, :(Union{Colon,AbstractRange}))
    @eval Base.@propagate_inbounds function Base.getindex(t::TrackedArray, i1::$(T), is::$(T)...)
        tp = tape(t)
        out = TrackedArray(value(t)[i1, is...], deriv(t)[i1, is...], tp)
        idx = index_iterable(axes(t), (i1, is...))
        record!(tp, SpecialInstruction, getindex, (t, idx), out)
        return out
    end
end

@noinline function special_reverse_exec!(instruction::SpecialInstruction{typeof(getindex)})
    input, indices = instruction.input
    output = instruction.output
    input_deriv = deriv(input)
    output_deriv = deriv(output)
    i = 0
    for idx in indices
        input_deriv[CartesianIndex(idx)] += output_deriv[i += 1]
    end
    unseed!(output)
    return nothing
end

@noinline function special_forward_exec!(instruction::SpecialInstruction{typeof(getindex)})
    input, indices = instruction.input
    input_value = value(input)
    output_value = value(instruction.output)
    i = 0
    for idx in indices
        output_value[i += 1] = input_value[CartesianIndex(idx)]
    end
    return nothing
end

Base.@propagate_inbounds function Base.getindex(t::TrackedArray, inds::AbstractArray{<:CartesianIndex})
    tp = tape(t)
    out = TrackedArray(value(t)[inds], deriv(t)[inds], tp)
    record!(tp, SpecialInstruction, getindex, (t, inds), out)
    return out
end
Base.@propagate_inbounds function Base.getindex(t::TrackedArray, i1::Integer, is::Integer...)
    ind = LinearIndices(t)[i1, is...]
    return TrackedReal(value(t)[i1, is...], deriv(t)[i1, is...], tape(t), ind, t)
end
Base.@propagate_inbounds function Base.getindex(t::TrackedArray, _inds1::Union{Integer, Colon, AbstractArray{<:Integer}}, _inds2::Union{Integer, Colon, AbstractArray{<:Integer}}...)
    inds1 = _inds1 isa Colon ? axes(t, 1) : _inds1
    inds2 = ntuple(Val(length(_inds2))) do i
        _inds2[i] isa Colon && return axes(t, i+1)
        return _inds2[i]
    end
    tp = tape(t)
    out = TrackedArray(value(t)[inds1, inds2...], deriv(t)[inds1, inds2...], tp)
    record!(tp, SpecialInstruction, (getindex, Val(:generic)), (t, (inds1, inds2...)), out)
    return out
end
@noinline function special_reverse_exec!(instruction::SpecialInstruction{<:Tuple{typeof(getindex), Val{:generic}}})
    input, inds = instruction.input
    output = instruction.output
    cinds = CartesianIndices(map(i -> 1:length(i), inds))
    input_deriv = deriv(input)
    output_deriv = deriv(output)
    i = 0
    for _idx in cinds
        idx = CartesianIndex(map(getindex, inds, Tuple(_idx)))
        input_deriv[idx] += output_deriv[i += 1]
    end
    unseed!(output)
    return nothing
end
@noinline function special_forward_exec!(instruction::SpecialInstruction{<:Tuple{typeof(getindex), Val{:generic}}})
    input, inds = instruction.input
    input_value = value(input)
    output_value = value(instruction.output)
    cinds = CartesianIndices(map(i -> 1:length(i), inds))
    i = 0
    for cind in cinds
        idx = CartesianIndex(map(getindex, inds, Tuple(cind)))
        output_value[i += 1] = input_value[idx]
    end
    return nothing
end

Base.setindex!(t::TrackedArray, args...) = error("TrackedArrays do not support setindex!")

Base.IndexStyle(::TrackedArray) = IndexLinear()

Base.size(t::TrackedArray) = size(value(t))

Base.copy(t::TrackedArray) = t

Base.similar(t::TrackedArray, args::Union{Integer, AbstractUnitRange}...) = similar(value(t), eltype(t), args...)

Base.similar(t::TrackedArray, T::Type, args::Union{Integer, AbstractUnitRange}...) = similar(value(t), T, args...)

reshape_body = :(TrackedArray(reshape(value(t), dims), reshape(deriv(t), dims), tape(t)))
@eval Base.reshape(t::TrackedArray, dims::Val{N}) where {N} = $reshape_body
@eval Base.reshape(t::TrackedArray, dims::Tuple{Vararg{Int,N}}) where {N} = $reshape_body
@eval Base.reshape(t::TrackedArray, dims::Int64...) = $reshape_body
@eval Base.reshape(t::TrackedArray, dims::AbstractUnitRange...) = $reshape_body
@eval Base.reshape(t::TrackedArray, dims::Colon...) = $reshape_body
@eval Base.reshape(t::TrackedArray, dims::Union{AbstractUnitRange,Int64,Colon}...) = $reshape_body

####################
# `Real` Interface #
####################

Base.hash(t::TrackedReal) = hash(value(t))
Base.hash(t::TrackedReal, hsh::UInt64) = hash(value(t), hsh)

Base.deepcopy(t::T) where {T<:TrackedReal} = t
Base.copy(t::T) where {T<:TrackedReal} = t

function Base.float(t::TrackedReal{V,D,O}) where {V,D,O}
    v = float(value(t))
    return TrackedReal{typeof(v),D,O}(v)
end

Base.float(t::TrackedReal{V}) where {V<:AbstractFloat} = t

Base.one(::Type{TrackedReal{V,D,O}}) where {V,D,O} = TrackedReal{V,D,O}(one(V))
Base.zero(::Type{TrackedReal{V,D,O}}) where {V,D,O} = TrackedReal{V,D,O}(zero(V))

Base.rand(::Type{TrackedReal{V,D,O}}) where {V,D,O} = TrackedReal{V,D,O}(rand(V))
Base.rand(rng::Random.AbstractRNG, ::Type{TrackedReal{V,D,O}}) where {V,D,O} = TrackedReal{V,D,O}(rand(rng, V))

Base.eps(t::TrackedReal) = eps(value(t))
Base.eps(::Type{T}) where {T<:TrackedReal} = eps(valtype(T))

Base.floor(t::TrackedReal) = floor(value(t))
Base.floor(::Type{R}, t::TrackedReal) where {R<:Real} = floor(R, value(t))

Base.ceil(t::TrackedReal) = ceil(value(t))
Base.ceil(::Type{R}, t::TrackedReal) where {R<:Real} = ceil(R, value(t))

Base.fld(a::TrackedReal, b::TrackedReal) = fld(value(a), value(b))

Base.cld(a::TrackedReal, b::TrackedReal) = cld(value(a), value(b))

if VERSION ≥ v"1.4"
    Base.div(x::TrackedReal, y::TrackedReal, r::RoundingMode) = div(value(x), value(y), r)
else
    Base.div(x::TrackedReal, y::TrackedReal) = div(value(x), value(y))
end

Base.trunc(t::TrackedReal) = trunc(value(t))
Base.trunc(::Type{R}, t::TrackedReal) where {R<:Real} = trunc(R, value(t))

Base.round(t::TrackedReal) = round(value(t))
Base.round(::Type{R}, t::TrackedReal) where {R<:Real} = round(R, value(t))

Base.oneunit(t::TrackedReal) = one(t)
Base.oneunit(::Type{T}) where {T<:TrackedReal} = one(T)

Base.rtoldefault(::Type{T}) where {T<:TrackedReal} = sqrt(eps(T))

################
# track/track! #
################

track(x::Real, tp::InstructionTape = InstructionTape()) = track(x, typeof(x), tp)

track(x::AbstractArray, tp::InstructionTape = InstructionTape()) = track(x, eltype(x), tp)

track(x::Real, ::Type{D}, tp::InstructionTape = InstructionTape()) where {D} = TrackedReal(x, zero(D), tp)

track(x::AbstractArray, ::Type{D}, tp::InstructionTape = InstructionTape()) where {D} = TrackedArray(x, fill!(similar(x, D), zero(D)), tp)

track!(t::TrackedArray, x::AbstractArray) = (value!(t, x); unseed!(t); t)

track!(t::TrackedReal, x::Real) = (value!(t, x); unseed!(t); t)

function track!(t::AbstractArray{TrackedReal{D,D,Nothing}}, x::AbstractArray, tp::InstructionTape) where D
    for i in eachindex(t)
        t[i] = track(x[i], D, tp)
    end
    return t
end

###################
# Pretty Printing #
###################

idstr(x) = string(objectid(x), base=62)[1:3]

function Base.show(io::IO, t::TrackedReal)
    tape_id = hastape(t) ? idstr(t.tape) : "---"
    origin_id = hasorigin(t) ? "$(t.index), $(idstr(t.origin))" : "---"
    id = idstr(t)
    print(io, "TrackedReal<$(id)>($(value(t)), $(deriv(t)), $(tape_id), $(origin_id))")
end

max_stars_repo_name: jrevels/ReverseDiffPrototype

max_issues_repo_name: jrevels/ReverseDiffPrototype

max_forks_repo_name: jrevels/ReverseDiffPrototype


text:
module HarwellRutherfordBoeing

using Printf
using SparseArrays

export HarwellBoeingMatrix, RutherfordBoeingData

include("hrb_utils.jl")
include("harwell_boeing.jl")
include("rutherford_boeing.jl")

end