text:
"""
 M6 function
 Quintic spline
               |  (3-q)^5-6(2-q)^5+15(1-q)^5  for 0 <= q < 1
 M6(x) = 1/120 |  (3-q)^5-6(2-q)^5            for 1 <= q < 2
               |  (3-q)^5                     for 2 <= q < 3
               |  0                           for q >= 3
"""

@inline function f_m6( q :: Float64 )

    if ( q < 1.0 )
        f_m6 = (3-q)^5-6*(2-q)^5+15*(1-q)^5
    elseif ( q >= 1.0 && q < 2.0 )
        f_m6 = (3-q)^5-6*(2-q)^5
    elseif ( q >= 2 && q < 3 )
        f_m6 = (3-q)^5
    else
	f_m6 = 0.0
    end

    f_m6 / 120

end

export compute_rho_m6!

function compute_rho_m6!( fields    :: MeshFields,
                          particles :: Particles,
                          xt        :: Array{ComplexF64,3},
                          ua        :: UA )

    fill!(fields.ρ , 0.0)
    nx, ny = fields.mesh.nx, fields.mesh.ny
    dx, dy = fields.mesh.dx, fields.mesh.dy

    xmin, xmax = fields.mesh.xmin, fields.mesh.xmax
    ymin, ymax = fields.mesh.ymin, fields.mesh.ymax
    dimx = xmax - xmin
    dimy = ymax - ymin

    for k = 1:particles.nbpart

        t = particles.t[k]

        mul!(ua.ftau, ua.ptau, view(xt,:,1,k))

        for n = 1:ua.ntau
            ua.ftau[n] *= exp(1im*ua.ltau[n]*t/ua.ε)/ua.ntau
        end

        xt1 = real(sum(ua.ftau))

        mul!(ua.ftau, ua.ptau, view(xt,:,2,k))

        for n = 1:ua.ntau
            ua.ftau[n] *= exp.(1im*ua.ltau[n]*t/ua.ε)/ua.ntau
        end

        xt2 = real(sum(ua.ftau))

        xt1 = mod( xt1 - xmin , dimx)
        xt2 = mod( xt2 - ymin , dimy)

        px = xt1/dx
        py = xt2/dy

        particles.x[1,k] = xt1 + xmin
        particles.x[2,k] = xt2 + ymin

        i  = trunc(Int32, px)
        j  = trunc(Int32, py)
        dpx = px - i
        dpy = py - j

        weight = particles.w

        im3 = mod(i-3,nx)+1
        im2 = mod(i-2,nx)+1
        im1 = mod(i-1,nx)+1
        ip1 = mod(i+1,nx)+1
        ip2 = mod(i+2,nx)+1
        ip3 = mod(i+3,nx)+1
        jm3 = mod(j-3,ny)+1
        jm2 = mod(j-2,ny)+1
        jm1 = mod(j-1,ny)+1
        jp1 = mod(j+1,ny)+1
        jp2 = mod(j+2,ny)+1
        jp3 = mod(j+3,ny)+1

        i = i+1
        j = j+1

        cm3x = f_m6(3+dpx)
        cp3x = f_m6(3-dpx)
        cm2x = f_m6(2+dpx)
        cp2x = f_m6(2-dpx)
        cm1x = f_m6(1+dpx)
        cp1x = f_m6(1-dpx)
        cx   = f_m6(dpx)
        cy   = f_m6(dpy)
        cp1y = f_m6(1-dpy)
        cm1y = f_m6(1+dpy)
        cp2y = f_m6(2-dpy)
        cm2y = f_m6(2+dpy)
        cp3y = f_m6(3-dpy)
        cm3y = f_m6(3+dpy)

	    fields.ρ[im3,jm3] += cm3x * cm3y * weight
        fields.ρ[im3,jm2] += cm3x * cm2y * weight
        fields.ρ[im3,jm1] += cm3x * cm1y * weight
        fields.ρ[im3,j  ] += cm3x * cy   * weight
        fields.ρ[im3,jp1] += cm3x * cp1y * weight
        fields.ρ[im3,jp2] += cm3x * cp2y * weight
        fields.ρ[im3,jp3] += cm3x * cp3y * weight

	    fields.ρ[im2,jm3] += cm2x * cm3y * weight
        fields.ρ[im2,jm2] += cm2x * cm2y * weight
        fields.ρ[im2,jm1] += cm2x * cm1y * weight
        fields.ρ[im2,j  ] += cm2x * cy   * weight
        fields.ρ[im2,jp1] += cm2x * cp1y * weight
        fields.ρ[im2,jp2] += cm2x * cp2y * weight
        fields.ρ[im2,jp3] += cm2x * cp3y * weight

	    fields.ρ[im1,jm3] += cm1x * cm3y * weight
        fields.ρ[im1,jm2] += cm1x * cm2y * weight
        fields.ρ[im1,jm1] += cm1x * cm1y * weight
        fields.ρ[im1,j  ] += cm1x * cy   * weight
        fields.ρ[im1,jp1] += cm1x * cp1y * weight
        fields.ρ[im1,jp2] += cm1x * cp2y * weight
        fields.ρ[im1,jp3] += cm1x * cp3y * weight

        fields.ρ[i  ,jm3] += cx   * cm3y * weight
        fields.ρ[i  ,jm2] += cx   * cm2y * weight
        fields.ρ[i  ,jm1] += cx   * cm1y * weight
        fields.ρ[i  ,j  ] += cx   * cy   * weight
        fields.ρ[i  ,jp1] += cx   * cp1y * weight
        fields.ρ[i  ,jp2] += cx   * cp2y * weight
        fields.ρ[i  ,jp3] += cx   * cp3y * weight

        fields.ρ[ip1,jm3] += cp1x * cm3y * weight
        fields.ρ[ip1,jm2] += cp1x * cm2y * weight
        fields.ρ[ip1,jm1] += cp1x * cm1y * weight
        fields.ρ[ip1,j  ] += cp1x * cy   * weight
        fields.ρ[ip1,jp1] += cp1x * cp1y * weight
        fields.ρ[ip1,jp2] += cp1x * cp2y * weight
        fields.ρ[ip1,jp3] += cp1x * cp3y * weight

        fields.ρ[ip2,jm3] += cp2x * cm3y * weight
        fields.ρ[ip2,jm2] += cp2x * cm2y * weight
        fields.ρ[ip2,jm1] += cp2x * cm1y * weight
        fields.ρ[ip2,j  ] += cp2x * cy   * weight
        fields.ρ[ip2,jp1] += cp2x * cp1y * weight
        fields.ρ[ip2,jp2] += cp2x * cp2y * weight
        fields.ρ[ip2,jp3] += cp2x * cp3y * weight

        fields.ρ[ip3,jm3] += cp3x * cm3y * weight
        fields.ρ[ip3,jm2] += cp3x * cm2y * weight
        fields.ρ[ip3,jm1] += cp3x * cm1y * weight
        fields.ρ[ip3,j  ] += cp3x * cy   * weight
        fields.ρ[ip3,jp1] += cp3x * cp1y * weight
        fields.ρ[ip3,jp2] += cp3x * cp2y * weight
        fields.ρ[ip3,jp3] += cp3x * cp3y * weight

    end

    fields.ρ[1:nx,ny+1] .= fields.ρ[1:nx,1]
    fields.ρ[nx+1,1:ny] .= fields.ρ[1,1:ny]
    fields.ρ[nx+1,ny+1]  = fields.ρ[1,1]

    fields.ρ ./= (dx*dy)

    rho_total = sum(view(fields.ρ,1:nx,1:ny)) * dx * dy

    fields.ρ .-= rho_total/dimx/dimy

end

function compute_rho_m6!( fields  :: MeshFields, particles :: Particles)

    fill!(fields.ρ , 0.0)
    nx = fields.mesh.nx
    ny = fields.mesh.ny

    dx = fields.mesh.dx
    dy = fields.mesh.dy

    xmin, xmax = fields.mesh.xmax, fields.mesh.xmin
    ymin, ymax = fields.mesh.ymax, fields.mesh.ymin
    dimx = fields.mesh.xmax - fields.mesh.xmin
    dimy = fields.mesh.ymax - fields.mesh.ymin

    for k = 1:particles.nbpart

        xt1 = mod( particles.x[1,k] - xmin , dimx)
        xt2 = mod( particles.x[2,k] - ymin , dimy)

        particles.x[1,k] = xt1 + xmin
        particles.x[2,k] = xt2 + ymin

        px = xt1/dx
        py = xt2/dy

        i = trunc(Int32, px)
        j = trunc(Int32, py)

        dpx = px - i
        dpy = py - j

        weight = particles.w

        im3 = mod(i-3,nx)+1
        im2 = mod(i-2,nx)+1
        im1 = mod(i-1,nx)+1
        ip1 = mod(i+1,nx)+1
        ip2 = mod(i+2,nx)+1
        ip3 = mod(i+3,nx)+1
        jm3 = mod(j-3,ny)+1
        jm2 = mod(j-2,ny)+1
        jm1 = mod(j-1,ny)+1
        jp1 = mod(j+1,ny)+1
        jp2 = mod(j+2,ny)+1
        jp3 = mod(j+3,ny)+1

        i = i+1
        j = j+1

        cm3x = f_m6(3+dpx)
        cp3x = f_m6(3-dpx)
        cm2x = f_m6(2+dpx)
        cp2x = f_m6(2-dpx)
        cm1x = f_m6(1+dpx)
        cp1x = f_m6(1-dpx)
        cx   = f_m6(dpx)
        cy   = f_m6(dpy)
        cp1y = f_m6(1-dpy)
        cm1y = f_m6(1+dpy)
        cp2y = f_m6(2-dpy)
        cm2y = f_m6(2+dpy)
        cp3y = f_m6(3-dpy)
        cm3y = f_m6(3+dpy)

	    fields.ρ[im3,jm3] += cm3x * cm3y * weight
        fields.ρ[im3,jm2] += cm3x * cm2y * weight
        fields.ρ[im3,jm1] += cm3x * cm1y * weight
        fields.ρ[im3,j  ] += cm3x * cy   * weight
        fields.ρ[im3,jp1] += cm3x * cp1y * weight
        fields.ρ[im3,jp2] += cm3x * cp2y * weight
        fields.ρ[im3,jp3] += cm3x * cp3y * weight

	    fields.ρ[im2,jm3] += cm2x * cm3y * weight
        fields.ρ[im2,jm2] += cm2x * cm2y * weight
        fields.ρ[im2,jm1] += cm2x * cm1y * weight
        fields.ρ[im2,j  ] += cm2x * cy   * weight
        fields.ρ[im2,jp1] += cm2x * cp1y * weight
        fields.ρ[im2,jp2] += cm2x * cp2y * weight
        fields.ρ[im2,jp3] += cm2x * cp3y * weight

	    fields.ρ[im1,jm3] += cm1x * cm3y * weight
        fields.ρ[im1,jm2] += cm1x * cm2y * weight
        fields.ρ[im1,jm1] += cm1x * cm1y * weight
        fields.ρ[im1,j  ] += cm1x * cy   * weight
        fields.ρ[im1,jp1] += cm1x * cp1y * weight
        fields.ρ[im1,jp2] += cm1x * cp2y * weight
        fields.ρ[im1,jp3] += cm1x * cp3y * weight

        fields.ρ[i  ,jm3] += cx   * cm3y * weight
        fields.ρ[i  ,jm2] += cx   * cm2y * weight
        fields.ρ[i  ,jm1] += cx   * cm1y * weight
        fields.ρ[i  ,j  ] += cx   * cy   * weight
        fields.ρ[i  ,jp1] += cx   * cp1y * weight
        fields.ρ[i  ,jp2] += cx   * cp2y * weight
        fields.ρ[i  ,jp3] += cx   * cp3y * weight

        fields.ρ[ip1,jm3] += cp1x * cm3y * weight
        fields.ρ[ip1,jm2] += cp1x * cm2y * weight
        fields.ρ[ip1,jm1] += cp1x * cm1y * weight
        fields.ρ[ip1,j  ] += cp1x * cy   * weight
        fields.ρ[ip1,jp1] += cp1x * cp1y * weight
        fields.ρ[ip1,jp2] += cp1x * cp2y * weight
        fields.ρ[ip1,jp3] += cp1x * cp3y * weight

        fields.ρ[ip2,jm3] += cp2x * cm3y * weight
        fields.ρ[ip2,jm2] += cp2x * cm2y * weight
        fields.ρ[ip2,jm1] += cp2x * cm1y * weight
        fields.ρ[ip2,j  ] += cp2x * cy   * weight
        fields.ρ[ip2,jp1] += cp2x * cp1y * weight
        fields.ρ[ip2,jp2] += cp2x * cp2y * weight
        fields.ρ[ip2,jp3] += cp2x * cp3y * weight

        fields.ρ[ip3,jm3] += cp3x * cm3y * weight
        fields.ρ[ip3,jm2] += cp3x * cm2y * weight
        fields.ρ[ip3,jm1] += cp3x * cm1y * weight
        fields.ρ[ip3,j  ] += cp3x * cy   * weight
        fields.ρ[ip3,jp1] += cp3x * cp1y * weight
        fields.ρ[ip3,jp2] += cp3x * cp2y * weight
        fields.ρ[ip3,jp3] += cp3x * cp3y * weight

    end

    fields.ρ[1:nx,ny+1] .= fields.ρ[1:nx,1]
    fields.ρ[nx+1,1:ny] .= fields.ρ[1,1:ny]
    fields.ρ[nx+1,ny+1]  = fields.ρ[1,1]

    fields.ρ ./= (dx*dy)

    rho_total = sum(fields.ρ[1:nx,1:ny]) * dx * dy

    println( " rho_total = $rho_total ")

    fields.ρ .-= rho_total/dimx/dimy

end


text:
using Nemo, FastArithmetic, BenchmarkTools

const k=ResidueRing(ZZ,5) # ,u=FiniteField(5,1,"u")
const K,t=PolynomialRing(k,"t")

Base.size(k::Nemo.GenResRing{Nemo.fmpz}) = Int(modulus(k))
Base.size{T<:Nemo.FinField}(k::T) = Int(order(k))

function create(n::Int64,k::Nemo.Ring)
  A=Array(k,n)
  c=size(k)-1
  for i in 1:n
    A[i]=k(rand(0:c))
  end
  return A
end

function createPol(n::Int64,K::Nemo.Ring)
  k=base_ring(K)
  return K(create(n,k))
end

function createIrrPol(n::Int64,K::Nemo.Ring)
  P=createPol(n,K)
  while degree(P)!=n || length(factor(P))!=1
    P=createPol(n,K)
  end
  return P
end

function create2d(m::Int64,n::Int64,k::Nemo.Ring)
  A=Array(k,n,m)
  c=size(k)-1
  for i in 1:n, j in 1:m
    A[i,j]=k(rand(0:c))
  end
  return A
end

ext: jl

max_stars_repo_name: artichaut/fastarithmetic

max_issues_repo_name: artichaut/fastarithmetic

max_forks_repo_name: artichaut/fastarithmetic


text:
using HTTP
using JSON

"""
    Subreddit(sub::AbstractString; feed::AbstractString="new")

Initialize a new subreddit source.

# Arguments
* `sub::AbstractString`: Subreddit to monitor.

# Keyword Arguments
* `feed::AbstractString="new"`: Retrieve new posts or hot posts.
"""
struct Subreddit <: Source
    sub::AbstractString  # Subreddit name.
    feed::AbstractString  # Feed to refresh.

    function Subreddit(sub::AbstractString; feed::AbstractString="new")
        sub = lowercase(sub)
        feed = lowercase(feed)
        if !in(feed, ["new", "hot"])
            error("Invalid feed.")
        end
        new(sub, feed)
    end
end

function listen(source::Subreddit, channel::Channel{Message})
    uri = HTTP.URI("http://reddit.com/r/$(source.sub)/$(source.feed).json?limit=1")
    id = ""
    while true
        try
            data = JSON.parse(HTTP.request(uri).body)["data"]["children"][end]["data"]
            if data["id"] != id
                id = data["id"]
                title = data["title"]
                link = "https://redd.it/$id"
                put!(channel, SubredditMessage(title, link))
            end
        catch
            warn("Parsing response failed.")
        end
        sleep(60)
    end
end

"""
    SubredditMessage(title::AbstractString, link::AbstractString) -> SubredditMessage

Create a new subreddit message.

# Arguments
* `title::AbstractString`: Title of the post.
* `link::AbstractString`: Shortened URL.
"""
struct SubredditMessage <: Message
    title::AbstractString  # Title of the post.
    link::AbstractString  # Shortened URL.

    SubredditMessage(title, link) = new(title, link)
end


text:
#------------------------------
# Linear ODE system with two unkown parameters
#------------------------------

#------------------------------
# Setup
#------------------------------

using Pkg;

# all paths are relative to the repository main folder

Pkg.activate(".")
Pkg.instantiate()
Pkg.status()

using Distributions
using Random
using Flux
using DiffEqFlux
using OrdinaryDiffEq
using Plots
using LaTeXStrings
gr()

include("../src/simulation.jl")
include("../src/model.jl")
include("../src/plotting.jl")

#------------------------------
# Define and obtain ground-truth developments
#------------------------------

# define initial condition
true_u0 = Float32[2, 1]
# define time span on which to solve the ODE
tspan = (0.0f0, 10.0f0)
# define parameters for the two distinct groups
true_odeparams_group1 = Float32[-0.2, 0.00, 0.00, -0.2]
true_odeparams_group2 = Float32[-0.2, 0.00, 0.00, 0.2]

# define corresponding ODE problems for the two groups
prob1 = ODEProblem(linear_2d_system,true_u0,tspan,true_odeparams_group1)
prob2 = ODEProblem(linear_2d_system,true_u0,tspan,true_odeparams_group2)

# solve ODE systems to obtain "true" underlying trajectory in each group
dt=0.1
sol_group1 = solve(prob1, Tsit5(), saveat = dt);
sol_group2 = solve(prob2, Tsit5(), saveat = dt);

#------------------------------
# Simulate data based on true trajectories
#------------------------------

# define numbers of individuals and variables
n = 100
p = 10
q, q_info = 50, 10

# choose simulation setting for baseline variables
baseline_simulation = "groupsonly" # alternatives: "trueparams" or "groupsonly"

# set seed for reproducibility
Random.seed!(12);

# generate time dependent variables
xs, tvals, group1, group2 = generate_xs(n, p, true_u0, sol_group1, sol_group2, maxntps=1, σ_var=0.1, σ_ind=0.1); # default vals: t_start=1.5, t_end=10, maxntps=10, dt=0.1, σ_var=0.1, σ_ind=0.5
# generate baseline variables
if baseline_simulation == "trueparams"
    x_baseline = generate_baseline(n, q, q_info, group1, true_odeparams_group1, true_odeparams_group2); # defaul vals: σ_info=0.1, σ_noise=0.1
elseif baseline_simulation == "groupsonly"
    x_baseline = generate_baseline(n, q, q_info, group1,  σ_info=0.5, σ_noise=0.5); # default vals: σ_info=1, σ_noise=1
else
    error("Please select either 'trueparams' or 'groupsonly' as mode for simulation of baseline variables")
end

# look at data:
data = simdata(xs, x_baseline, tvals, group1, group2);

plot(plot_truesolution(2, data, sol_group1, sol_group2, showdata=false),
    plot_truesolution(1,  data, sol_group1, sol_group2, showdata=false),
    plot_truesolution(2, data, sol_group1, sol_group2, showdata=true),
    plot_truesolution(1,  data, sol_group1, sol_group2, showdata=true),
    layout = (2,2),
    legend = false
)

#------------------------------
# Define and train model
#------------------------------

zdim = nODEparams = 2
m = init_vae(p, q, zdim, nODEparams, prob1)
L = loss_wrapper(m)
ps = getparams(m)
opt = ADAM(0.0005)
trainingdata = zip(xs, x_baseline, tvals);
evalcb() = @show(mean(L.(xs, x_baseline, tvals)))
evalcb_zs() = eval_z_trajectories(xs, x_baseline, tvals, group1, sol_group1, sol_group2, m, dt, swapcolorcoding=true)

for epoch in 1:35
    Flux.train!(L, ps, trainingdata, opt)
    evalcb()
    evalcb_zs()
end

#------------------------------
# Look at results
#------------------------------

# individual plots
p1 = plot_individual_solutions(20,xs, x_baseline, tvals, group1, sol_group1, sol_group2, m, dt,
                                swapcolorcoding=true, showlegend=true)
p2 = plot_individual_solutions(28,xs, x_baseline, tvals, group1, sol_group1, sol_group2, m, dt,
                                swapcolorcoding=true, showlegend=false)
p3 = allindsplot(2, data, m, sol_group1, sol_group2, swapcolorcoding=true)
p4 = allindsplot(1, data, m, sol_group1, sol_group2, swapcolorcoding=true, showlegend=false)

# solutions across the entire dataset
plot(p3,p4, layout = (1,2))


text:
## setup_3rd_party_packages.jl
#
## installs julia packages
## julia version is 0.6.2

## some annotations
## Tom Price December 2018
#
## Currently the following are not used as much as they might be:
# DataFrames
# Nullables

Pkg.init()
dir = Pkg.Dir.path()
info("Pinning repository $dir")
metadata_dir = joinpath(dir, "METADATA")
if Pkg.Dir.isdir(metadata_dir)
    info("Package directory $dir exists.")
    run(`git -C $metadata_dir reset --hard 62c084ec638baee70c4f361c75304583bf4647da`)
end

function install(
    library ::String,
    vers 	::Vararg{VersionNumber,2}
)
	println("installing library: ", library)
	Pkg.add(library, vers...)
end

install("FunctionalData", v"0.1.2", v"0.1.3-")
install("Match", v"0.4.0", v"0.4.1-")
install("Compat", v"0.61.0", v"0.61.1-")
install("FactCheck", v"0.4.3", v"0.4.4-") ## for testing and precompiling
install("GZip", v"0.3.0", v"0.3.1-")
install("Polynomials", v"0.2.2", v"0.2.3-")
install("IterTools", v"0.2.1", v"0.2.2-")
install("Combinatorics", v"0.6.0", v"0.6.1-")
install("SHA", v"0.5.7", v"0.5.8-")
install("URIParser", v"0.3.1", v"0.3.2-")
install("BinDeps", v"0.8.7", v"0.8.8-")
install("Blosc", v"0.4.2", v"0.4.3-")
install("BinaryProvider", v"0.3.0", v"0.3.1-")
install("Reexport", v"0.1.0", v"0.1.1-")
install("NaNMath", v"0.3.1", v"0.3.2-")
install("MathProgBase", v"0.6.0", v"0.6.1-")
install("DataStructures", v"0.7.4", v"0.7.5-")
install("Missings", v"0.2.7", v"0.2.8-")
install("WeakRefStrings", v"0.4.3", v"0.4.4-")
install("NamedTuples", v"4.0.2", v"4.0.3-")
install("Nullables", v"0.0.3", v"0.0.4-")
install("JSON", v"0.16.4", v"0.16.5-")
install("Syslogs", v"0.1.1", v"0.1.2-")
install("Memento", v"0.6.0", v"0.6.1-")
install("CategoricalArrays", v"0.3.6", v"0.3.7-")
install("Calculus", v"0.2.2", v"0.2.3-")
install("DataStreams", v"0.3.4", v"0.3.5-")
install("DiffBase", v"0.2.0", v"0.2.1-")
install("Distances", v"0.6.0", v"0.6.1-")
install("FileIO", v"0.7.0", v"0.7.1-")
install("SpecialFunctions", v"0.3.6", v"0.3.7-")
install("ForwardDiff", v"0.4.2", v"0.4.3-")
install("LegacyStrings", v"0.3.0", v"0.3.1-")
install("HttpCommon", v"0.4.0", v"0.4.1-")
install("HttpParser", v"0.3.1", v"0.3.2-")
install("MbedTLS", v"0.5.8", v"0.5.9-")
install("IniFile", v"0.4.0", v"0.4.1-")
install("HTTP", v"0.6.9", v"0.7.0-")

install("StaticArrays", v"0.7.0", v"0.7.1-")
install("NearestNeighbors", v"0.3.0", v"0.3.1-")
install("ReverseDiffSparse", v"0.7.3", v"0.7.4-")
install("SortingAlgorithms", v"0.2.0", v"0.2.1-")
install("StatsBase", v"0.22.0", v"0.22.1-")

println("All additional packages")

install("DataFrames", v"0.11.0", v"0.11.1-")
install("HDF5", v"0.8.8", v"0.8.9-")
install("JLD", v"0.8.3", v"0.8.4-") ## needed by Ipopt, JuMP, and NLopt
install("Clustering", v"0.9.1", v"0.9.2-")
install("Ipopt", v"0.2.4", v"0.2.9-")
install("JuMP", v"0.17.1", v"0.17.2-")
#install("HttpServer", v"0.2.0", v"0.2.1-")
install("Dierckx", v"0.3.0", v"0.3.1-")
install("DataArrays", v"0.7.0", v"0.7.1-")
## Packages that are no longer used:

# install("Cairo", v"0.5.1", v"0.5.2-")
# install("CodecZlib", v"0.4.2", v"0.4.3-")
# install("Colors", v"0.8.2", v"0.8.3-")
# install("ColorTypes", v"0.6.7", v"0.6.8-")
# install("FixedPointNumbers", v"0.4.6", v"0.4.7-")
# install("Graphics", v"0.2.0", v"0.2.1-")
# install("Gtk", v"0.13.1", v"0.13.2-")
# install("GtkReactive", v"0.4.0", v"0.4.1-")
# install("IntervalSets", v"0.2.0", v"0.2.1-")
# install("MicroLogging", v"0.2.0", v"0.2.0-")
# install("MySQL", v"0.3.0", v"0.3.1-") # remove MySQL dependency
# install("NLopt", v"0.3.6", v"0.3.7-")
# install("ProfileView", v"0.3.0", v"0.3.1-")
# install("Reactive", v"0.6.0", v"0.6.1-")
# install("RoundingIntegers", v"0.0.3", v"0.0.4-")
# install("TranscodingStreams", v"0.5.1", v"0.5.2-")

Pkg.build(
	"DataArrays",
	"DataStructures",
	"SpecialFunctions",
	"MathProgBase",
	"HDF5",
	"Clustering",
	"Nullables",
	"JSON",
	"HTTP",
	"JLD",
	"JuMP",
	"Dierckx",
	"Ipopt",
	"DataFrames")
println("Building: Done")

include("/root/chaipcr/bioinformatics/build/useallpackages.jl")

Pkg.status()
println("Using: Done")

#

max_stars_repo_name: MakerButt/chaipcr

max_issues_repo_name: MakerButt/chaipcr

max_forks_repo_name: MakerButt/chaipcr


text:
### Views

jvec_rgn(x::Vector, first::Int, n::Int) = pointer_to_array(pointer(x, first), n, false)

as_jvec{T}(x::ContiguousView{T,1,Vector{T}}) = pointer_to_array(pointer(x), length(x), false)

function view(x::SparseMatrixCSC, ::Colon, j::Integer)
    1 <= j <= x.n || throw(BoundsError())
    r1 = convert(Int, x.colptr[j])
    r2 = convert(Int, x.colptr[j+1]) - 1
    rgn = r1:r2
    SparseVectorView(x.m, view(x.rowval, rgn), view(x.nzval, rgn))
end

function getcol(x::SparseMatrixCSC, j::Integer)
    1 <= j <= x.n || throw(BoundsError())
    r1 = convert(Int, x.colptr[j])
    r2 = convert(Int, x.colptr[j+1]) - 1
    SparseVector(x.m, x.rowval[r1:r2], x.nzval[r1:r2])
end

function unsafe_colrange{Tv,Ti}(x::SparseMatrixCSC{Tv,Ti}, J::UnitRange)
    jfirst = first(J)
    jlast = last(J)
    (1 <= jfirst <= x.n && jlast <= x.n) || throw(BoundsError())
    r1 = x.colptr[jfirst]
    r2 = x.colptr[jlast+1] - one(r1)
    newcolptr = view(x.colptr, jfirst:jlast+1) - (r1 - one(r1))

    fi = convert(Int, r1)
    nc = convert(Int, r2 - r1) + 1
    SparseMatrixCSC{Tv, Ti}(x.m, length(J), newcolptr,
        jvec_rgn(x.rowval, fi, nc), jvec_rgn(x.nzval, fi, nc))
end


text:
using aoc19
using Test, aoc19

@testset "aoc19.jl" begin
    @test_nowarn aoc19.benchmarkAll(onlyOnce = true)
    @test aoc19.formatTime(0.1) == "100.000 ms"
end

@testset "day1" begin
    @testset "part1" begin
        @test aoc19.day1.fuelRequired(12) == 2
        @test aoc19.day1.fuelRequired(14) == 2
        @test aoc19.day1.fuelRequired(1969) == 654
        @test aoc19.day1.fuelRequired(100756) == 33583
        @test part1(day = 1) == 3373568
    end

    @testset "part2" begin
        @test aoc19.day1.fuelRequired2(14) == 2
        @test aoc19.day1.fuelRequired2(1969) == 966
        @test aoc19.day1.fuelRequired2(100756) == 50346
        @test part2(day = 1) == 5057481
    end
end

@testset "day2" begin
    @testset "part1" begin
        @test part1(day = 2) == 9706670
    end

    @testset "part2" begin
        @test part2(day = 2) == 2552
    end
end

@testset "day3" begin
    ex1 = ["R8" "U5" "L5" "D3" ; "U7" "R6" "D4" "L4"]
    ex2 = ["R75" "D30" "R83" "U83" "L12" "D49" "R71" "U7" "L72" ; "U62" "R66" "U55" "R34" "D71" "R55" "D58" "R83" ""]
    ex3 = ["R98" "U47" "R26" "D63" "R33" "U87" "L62" "D20" "R33" "U53" "R51" ; "U98" "R91" "D20" "R16" "D67" "R40" "U7" "R15" "U6" "R7" ""]
    @testset "part1" begin
        @test aoc19.day3.part1(ex1) == 6
        @test aoc19.day3.part1(ex2) == 159
        @test aoc19.day3.part1(ex3) == 135
        @test part1(day = 3) == 709
    end

    @testset "part2" begin
        @test aoc19.day3.part2(ex1) == 30
        @test aoc19.day3.part2(ex2) == 610
        @test aoc19.day3.part2(ex3) == 410
        @test part2(day = 3) == 13836
    end
end

@testset "day4" begin
    @testset "part1" begin
        @test part1(day = 4) == 2814
    end
    @testset "part2" begin
        @test part2(day = 4) == 1991
    end
end

@testset "day5" begin
    @testset "part1" begin
        @test part1(day = 5) == 10987514
    end

    @testset "part2" begin
        program = [3,21,1008,21,8,20,1005,20,22,107,8,21,20,1006,20,31,
        1106,0,36,98,0,0,1002,21,125,20,4,20,1105,1,46,104,
        999,1105,1,46,1101,1000,1,20,4,20,1105,1,46,98,99]
        @test aoc19.day5.part1(program, -1) == 999
        @test aoc19.day5.part1(program, 7) == 999
        @test aoc19.day5.part1(program, 8) == 1000
        @test aoc19.day5.part1(program, 9) == 1001
        @test part2(day = 5) == 14195011
    end
end

@testset "day6" begin
    @testset "part1" begin
        input = [["COM","B"], ["B","C"], ["C","D"], ["D","E"], ["E","F"], ["B","G"], ["G","H"], ["D","I"], ["E","J"], ["J","K"], ["K","L"]]
        @test aoc19.day6.part1(input) == 42
        @test part1(day = 6) == 245089
    end

    @testset "part2" begin
        input = [["COM","B"], ["B","C"], ["C","D"], ["D","E"], ["E","F"], ["B","G"], ["G","H"], ["D","I"], ["E","J"], ["J","K"], ["K","L"], ["K", "YOU"], ["I", "SAN"]]
        @test aoc19.day6.part2(input) == 4
        @test part2(day = 6) == 511
    end
end

@testset "day7" begin
    @testset "part1" begin
        @test aoc19.day7.part1([3,15,3,16,1002,16,10,16,1,16,15,15,4,15,99,0,0]) == 43210
        @test aoc19.day7.part1([3,23,3,24,1002,24,10,24,1002,23,-1,23,101,5,23,23,1,24,23,23,4,23,99,0,0]) == 54321
        @test aoc19.day7.part1([3,31,3,32,1002,32,10,32,1001,31,-2,31,1007,31,0,33,1002,33,7,33,1,33,31,31,1,32,31,31,4,31,99,0,0,0]) == 65210
        @test part1(day = 7) == 87138
    end

    @testset "part2" begin
        @test aoc19.day7.part2([3,26,1001,26,-4,26,3,27,1002,27,2,27,1,27,26,27,4,27,1001,28,-1,28,1005,28,6,99,0,0,5]) == 139629729
        @test aoc19.day7.part2([3,52,1001,52,-5,52,3,53,1,52,56,54,1007,54,5,55,1005,55,26,1001,54,
        -5,54,1105,1,12,1,53,54,53,1008,54,0,55,1001,55,1,55,2,53,55,53,4,
        53,1001,56,-1,56,1005,56,6,99,0,0,0,0,10]) == 18216
        @test part2(day = 7) == 17279674
    end
end

@testset "day8" begin
    @testset "part1" begin
        @test part1(day = 8) == 2125
    end

    @testset "part2" begin
        @test part2(day = 8) isa aoc19.day8.UnicodePlots.Plot
    end
end

@testset "day10" begin
    @testset "part1" begin
        @test part1(day = 10) == 260
    end

    @testset "part2" begin
        @test part2(day = 10) == 608
    end
end

@testset "day11" begin
    @testset "part1" begin
        @test part1(day = 11) == 2720
    end

    @testset "part2" begin
        @test part2(day = 11) isa aoc19.day11.UnicodePlots.Plot
    end
end

@testset "day12" begin
    @testset "part1" begin
        @test aoc19.day12.part1([[-1, 0, 2], [2, -10, -7], [4, -8, 8], [3, 5, -1]], 10) == 179
        @test aoc19.day12.part1([[-8, -10, 0], [5, 5, 10], [2, -7, 3], [9, -8, -3]], 100) == 1940
        @test part1(day = 12) == 14606
    end

    @testset "part1" begin
        @test aoc19.day12.part2([[-1, 0, 2], [2, -10, -7], [4, -8, 8], [3, 5, -1]]) == 2772
        @test aoc19.day12.part2([[-8, -10, 0], [5, 5, 10], [2, -7, 3], [9, -8, -3]]) == 4686774924
        @test part2(day = 12) == 543673227860472
    end
end

@testset "day13" begin
    @testset "part1" begin
        @test part1(day = 13) == 318
    end

    @testset "part1" begin
        @test part2(day = 13) == 16309
    end
end

@testset "day14" begin
    @testset "part1" begin
        @test part1(day = 14) == 136771
    end

    @testset "part1" begin
        @test part2(day = 14) == 8193614
    end
end


text:
@testset "permutedims" begin

    @testset "Matrix" begin
        @test AxisIndices.permute_axes(rand(4, 2)) === (Base.OneTo(2), Base.OneTo(4))
        axs = @inferred(AxisIndices.permute_axes((Axis(1:4), Axis(1:2))))
        @test axs isa Tuple{Axis,Axis}
        @test length.(axs) == (2, 4)
    end

    @testset "Vector" begin
        @test length.(AxisIndices.permute_axes(rand(4))) == (1, 4)
        axs = AxisIndices.permute_axes((Axis(1:4),))
        @test axs isa Tuple{SimpleAxis,Axis}
        @test length.(axs) == (1, 4)

        v = AxisArray([10, 20, 30], (2:4,))
        @test permutedims(v) == [10 20 30]
        @test keys.(axes(permutedims(v))) == (1:1, 2:4)
        @test permutedims(permutedims(v)) == reshape([10, 20, 30], Val(2))
        @test keys.(axes(permutedims(permutedims(v)))) == (2:4, 1:1)
    end

    @testset "Array" begin
        @test length.(AxisIndices.permute_axes(rand(2, 4, 6), (1, 3, 2))) == (2, 6, 4)

        axs = AxisIndices.permute_axes((Axis(1:2), Axis(1:4), Axis(1:6)), (1, 3, 2))
        @test axs isa Tuple{Axis,Axis,Axis}
        @test length.(axs) == (2, 6, 4)

        a = AxisArray(ones(10, 20, 30, 40), (2:11, 2:21, 2:31, 2:41));
        @test (keys.(axes(permutedims(a, (1, 2, 3, 4)))) ==
            keys.(axes(permutedims(a, 1:4))) ==
            (2:11, 2:21, 2:31, 2:41)
        )
        @test (keys.(axes(permutedims(a, (1, 3, 2, 4)))) == (2:11, 2:31, 2:21, 2:41))
    end
end

@testset "permuteddimsview" begin
    @testset "standard" begin
        a = [1 3; 2 4]
        v = permuteddimsview(a, (1,2))
        @test v == a
        v = permuteddimsview(a, (2,1))
        @test v == a'
        a = rand(3,7,5)
        v = permuteddimsview(a, (2,3,1))
        @test v == permutedims(a, (2,3,1))
    end

    @testset "AxisArray" begin
        a = AxisArray([1 3; 2 4], [:one, :two], ["three", "four"]);
        v = permuteddimsview(a, (1,2))
        @test v == a
        @test keys.(axes(v)) == ([:one, :two], ["three", "four"])
        v = permuteddimsview(a, (2,1))
        @test v == a'
        @test keys.(axes(v)) == (["three", "four"], [:one, :two])
        a = AxisArray(rand(2,3,4), ["a", "b"], [:a, :b, :c], [1,2,3,4])
        v = permuteddimsview(a, (2,3,1))
        @test v == permutedims(a, (2,3,1))
        @test keys.(axes(v)) == ([:a, :b, :c], [1,2,3,4], ["a", "b"])
    end
end

@testset "PermuteDimsArray" begin
    x = AxisArray(ones(2,2))
    y = PermutedDimsArray(x, (2, 1))
    @test axes(y) isa Tuple{SimpleAxis, SimpleAxis}
    @test axes(y, 1) isa SimpleAxis
end


text:
# This file was generated by the Julia Swagger Code Generator
# Do not modify this file directly. Modify the swagger specification instead.

struct FilesApi <: SwaggerApi
    client::Swagger.Client
end

function _swaggerinternal_filesCommentsDelete(_api::FilesApi; file=nothing, id=nothing, _mediaType=nothing)
    _ctx = Swagger.Ctx(_api.client, "POST", Dict{String, Any}, "/files.comments.delete", ["slackAuth"])
    Swagger.set_param(_ctx.form, "file", file)  # type String
    Swagger.set_param(_ctx.form, "id", id)  # type String
    Swagger.set_header_accept(_ctx, ["application/json"])
    Swagger.set_header_content_type(_ctx, (_mediaType === nothing) ? ["application/x-www-form-urlencoded", "application/json"] : [_mediaType])
    return _ctx
end

"""

Deletes an existing comment on a file.
Param: file::String
Param: id::String
Return: Dict{String, Any}
"""
function filesCommentsDelete(_api::FilesApi; file=nothing, id=nothing, _mediaType=nothing)
    _ctx = _swaggerinternal_filesCommentsDelete(_api; file=file, id=id, _mediaType=_mediaType)
    Swagger.exec(_ctx)
end

function filesCommentsDelete(_api::FilesApi, response_stream::Channel; file=nothing, id=nothing, _mediaType=nothing)
    _ctx = _swaggerinternal_filesCommentsDelete(_api; file=file, id=id, _mediaType=_mediaType)
    Swagger.exec(_ctx, response_stream)
end

function _swaggerinternal_filesDelete(_api::FilesApi; file=nothing, _mediaType=nothing)
    _ctx = Swagger.Ctx(_api.client, "POST", Dict{String, Any}, "/files.delete", ["slackAuth"])
    Swagger.set_param(_ctx.form, "file", file)  # type String
    Swagger.set_header_accept(_ctx, ["application/json"])
    Swagger.set_header_content_type(_ctx, (_mediaType === nothing) ? ["application/x-www-form-urlencoded", "application/json"] : [_mediaType])
    return _ctx
end

"""

Deletes a file.
Param: file::String
Return: Dict{String, Any}
"""
function filesDelete(_api::FilesApi; file=nothing, _mediaType=nothing)
    _ctx = _swaggerinternal_filesDelete(_api; file=file, _mediaType=_mediaType)
    Swagger.exec(_ctx)
end

function filesDelete(_api::FilesApi, response_stream::Channel; file=nothing, _mediaType=nothing)
    _ctx = _swaggerinternal_filesDelete(_api; file=file, _mediaType=_mediaType)
    Swagger.exec(_ctx, response_stream)
end

function _swaggerinternal_filesInfo(_api::FilesApi; file=nothing, count=nothing, page=nothing, limit=nothing, cursor=nothing, _mediaType=nothing)
    _ctx = Swagger.Ctx(_api.client, "GET", Dict{String, Any}, "/files.info", ["slackAuth"])
    Swagger.set_param(_ctx.query, "file", file)  # type String
    Swagger.set_param(_ctx.query, "count", count)  # type String
    Swagger.set_param(_ctx.query, "page", page)  # type String
    Swagger.set_param(_ctx.query, "limit", limit)  # type Int32
    Swagger.set_param(_ctx.query, "cursor", cursor)  # type String
    Swagger.set_header_accept(_ctx, ["application/json"])
    Swagger.set_header_content_type(_ctx, (_mediaType === nothing) ? ["application/x-www-form-urlencoded"] : [_mediaType])
    return _ctx
end

"""

Gets information about a file.
Param: file::String
Param: count::String
Param: page::String
Param: limit::Int32
Param: cursor::String
Return: Dict{String, Any}
"""
function filesInfo(_api::FilesApi; file=nothing, count=nothing, page=nothing, limit=nothing, cursor=nothing, _mediaType=nothing)
    _ctx = _swaggerinternal_filesInfo(_api; file=file, count=count, page=page, limit=limit, cursor=cursor, _mediaType=_mediaType)
    Swagger.exec(_ctx)
end

function filesInfo(_api::FilesApi, response_stream::Channel; file=nothing, count=nothing, page=nothing, limit=nothing, cursor=nothing, _mediaType=nothing)
    _ctx = _swaggerinternal_filesInfo(_api; file=file, count=count, page=page, limit=limit, cursor=cursor, _mediaType=_mediaType)
    Swagger.exec(_ctx, response_stream)
end

function _swaggerinternal_filesList(_api::FilesApi; user=nothing, channel=nothing, ts_from=nothing, ts_to=nothing, types=nothing, count=nothing, page=nothing, show_files_hidden_by_limit=nothing, _mediaType=nothing)
    _ctx = Swagger.Ctx(_api.client, "GET", Dict{String, Any}, "/files.list", ["slackAuth"])
    Swagger.set_param(_ctx.query, "user", user)  # type String
    Swagger.set_param(_ctx.query, "channel", channel)  # type String
    Swagger.set_param(_ctx.query, "ts_from", ts_from)  # type Float32
    Swagger.set_param(_ctx.query, "ts_to", ts_to)  # type Float32
    Swagger.set_param(_ctx.query, "types", types)  # type String
    Swagger.set_param(_ctx.query, "count", count)  # type String
    Swagger.set_param(_ctx.query, "page", page)  # type String
    Swagger.set_param(_ctx.query, "show_files_hidden_by_limit", show_files_hidden_by_limit)  # type Bool
    Swagger.set_header_accept(_ctx, ["application/json"])
    Swagger.set_header_content_type(_ctx, (_mediaType === nothing) ? ["application/x-www-form-urlencoded"] : [_mediaType])
    return _ctx
end

"""

List for a team, in a channel, or from a user with applied filters.
Param: user::String
Param: channel::String
Param: ts_from::Float32
Param: ts_to::Float32
Param: types::String
Param: count::String
Param: page::String
Param: show_files_hidden_by_limit::Bool
Return: Dict{String, Any}
"""
function filesList(_api::FilesApi; user=nothing, channel=nothing, ts_from=nothing, ts_to=nothing, types=nothing, count=nothing, page=nothing, show_files_hidden_by_limit=nothing, _mediaType=nothing)
    _ctx = _swaggerinternal_filesList(_api; user=user, channel=channel, ts_from=ts_from, ts_to=ts_to, types=types, count=count, page=page, show_files_hidden_by_limit=show_files_hidden_by_limit, _mediaType=_mediaType)
    Swagger.exec(_ctx)
end

function filesList(_api::FilesApi, response_stream::Channel; user=nothing, channel=nothing, ts_from=nothing, ts_to=nothing, types=nothing, count=nothing, page=nothing, show_files_hidden_by_limit=nothing, _mediaType=nothing)
    _ctx = _swaggerinternal_filesList(_api; user=user, channel=channel, ts_from=ts_from, ts_to=ts_to, types=types, count=count, page=page, show_files_hidden_by_limit=show_files_hidden_by_limit, _mediaType=_mediaType)
    Swagger.exec(_ctx, response_stream)
end

function _swaggerinternal_filesRemoteAdd(_api::FilesApi; external_id=nothing, title=nothing, filetype=nothing, external_url=nothing, preview_image=nothing, indexable_file_contents=nothing, _mediaType=nothing)
    _ctx = Swagger.Ctx(_api.client, "POST", Dict{String, Any}, "/files.remote.add", ["slackAuth"])
    Swagger.set_param(_ctx.form, "external_id", external_id)  # type String
    Swagger.set_param(_ctx.form, "title", title)  # type String
    Swagger.set_param(_ctx.form, "filetype", filetype)  # type String
    Swagger.set_param(_ctx.form, "external_url", external_url)  # type String
    Swagger.set_param(_ctx.form, "preview_image", preview_image)  # type String
    Swagger.set_param(_ctx.form, "indexable_file_contents", indexable_file_contents)  # type String
    Swagger.set_header_accept(_ctx, ["application/json"])
    Swagger.set_header_content_type(_ctx, (_mediaType === nothing) ? ["application/x-www-form-urlencoded"] : [_mediaType])
    return _ctx
end

"""

Adds a file from a remote service
Param: external_id::String
Param: title::String
Param: filetype::String
Param: external_url::String
Param: preview_image::String
Param: indexable_file_contents::String
Return: Dict{String, Any}
"""
function filesRemoteAdd(_api::FilesApi; external_id=nothing, title=nothing, filetype=nothing, external_url=nothing, preview_image=nothing, indexable_file_contents=nothing, _mediaType=nothing)
    _ctx = _swaggerinternal_filesRemoteAdd(_api; external_id=external_id, title=title, filetype=filetype, external_url=external_url, preview_image=preview_image, indexable_file_contents=indexable_file_contents, _mediaType=_mediaType)
    Swagger.exec(_ctx)
end

function filesRemoteAdd(_api::FilesApi, response_stream::Channel; external_id=nothing, title=nothing, filetype=nothing, external_url=nothing, preview_image=nothing, indexable_file_contents=nothing, _mediaType=nothing)
    _ctx = _swaggerinternal_filesRemoteAdd(_api; external_id=external_id, title=title, filetype=filetype, external_url=external_url, preview_image=preview_image, indexable_file_contents=indexable_file_contents, _mediaType=_mediaType)
    Swagger.exec(_ctx, response_stream)
end

function _swaggerinternal_filesRemoteInfo(_api::FilesApi; file=nothing, external_id=nothing, _mediaType=nothing)
    _ctx = Swagger.Ctx(_api.client, "GET", Dict{String, Any}, "/files.remote.info", ["slackAuth"])
    Swagger.set_param(_ctx.query, "file", file)  # type String
    Swagger.set_param(_ctx.query, "external_id", external_id)  # type String
    Swagger.set_header_accept(_ctx, ["application/json"])
    Swagger.set_header_content_type(_ctx, (_mediaType === nothing) ? ["application/x-www-form-urlencoded"] : [_mediaType])
    return _ctx
end

"""

Retrieve information about a remote file added to Slack
Param: file::String
Param: external_id::String
Return: Dict{String, Any}
"""
function filesRemoteInfo(_api::FilesApi; file=nothing, external_id=nothing, _mediaType=nothing)
    _ctx = _swaggerinternal_filesRemoteInfo(_api; file=file, external_id=external_id, _mediaType=_mediaType)
    Swagger.exec(_ctx)
end

function filesRemoteInfo(_api::FilesApi, response_stream::Channel; file=nothing, external_id=nothing, _mediaType=nothing)
    _ctx = _swaggerinternal_filesRemoteInfo(_api; file=file, external_id=external_id, _mediaType=_mediaType)
    Swagger.exec(_ctx, response_stream)
end

function _swaggerinternal_filesRemoteList(_api::FilesApi; channel=nothing, ts_from=nothing, ts_to=nothing, limit=nothing, cursor=nothing, _mediaType=nothing)
    _ctx = Swagger.Ctx(_api.client, "GET", Dict{String, Any}, "/files.remote.list", ["slackAuth"])
    Swagger.set_param(_ctx.query, "channel", channel)  # type String
    Swagger.set_param(_ctx.query, "ts_from", ts_from)  # type Float32
    Swagger.set_param(_ctx.query, "ts_to", ts_to)  # type Float32
    Swagger.set_param(_ctx.query, "limit", limit)  # type Int32
    Swagger.set_param(_ctx.query, "cursor", cursor)  # type String
    Swagger.set_header_accept(_ctx, ["application/json"])
    Swagger.set_header_content_type(_ctx, (_mediaType === nothing) ? ["application/x-www-form-urlencoded"] : [_mediaType])
    return _ctx
end

"""

Retrieve information about a remote file added to Slack
Param: channel::String
Param: ts_from::Float32
Param: ts_to::Float32
Param: limit::Int32
Param: cursor::String
Return: Dict{String, Any}
"""
function filesRemoteList(_api::FilesApi; channel=nothing, ts_from=nothing, ts_to=nothing, limit=nothing, cursor=nothing, _mediaType=nothing)
    _ctx = _swaggerinternal_filesRemoteList(_api; channel=channel, ts_from=ts_from, ts_to=ts_to, limit=limit, cursor=cursor, _mediaType=_mediaType)
    Swagger.exec(_ctx)
end

function filesRemoteList(_api::FilesApi, response_stream::Channel; channel=nothing, ts_from=nothing, ts_to=nothing, limit=nothing, cursor=nothing, _mediaType=nothing)
    _ctx = _swaggerinternal_filesRemoteList(_api; channel=channel, ts_from=ts_from, ts_to=ts_to, limit=limit, cursor=cursor, _mediaType=_mediaType)
    Swagger.exec(_ctx, response_stream)
end

function _swaggerinternal_filesRemoteRemove(_api::FilesApi; file=nothing, external_id=nothing, _mediaType=nothing)
    _ctx = Swagger.Ctx(_api.client, "POST", Dict{String, Any}, "/files.remote.remove", ["slackAuth"])
    Swagger.set_param(_ctx.form, "file", file)  # type String
    Swagger.set_param(_ctx.form, "external_id", external_id)  # type String
    Swagger.set_header_accept(_ctx, ["application/json"])
    Swagger.set_header_content_type(_ctx, (_mediaType === nothing) ? ["application/x-www-form-urlencoded"] : [_mediaType])
    return _ctx
end

"""

Remove a remote file.
Param: file::String
Param: external_id::String
Return: Dict{String, Any}
"""
function filesRemoteRemove(_api::FilesApi; file=nothing, external_id=nothing, _mediaType=nothing)
    _ctx = _swaggerinternal_filesRemoteRemove(_api; file=file, external_id=external_id, _mediaType=_mediaType)
    Swagger.exec(_ctx)
end

function filesRemoteRemove(_api::FilesApi, response_stream::Channel; file=nothing, external_id=nothing, _mediaType=nothing)
    _ctx = _swaggerinternal_filesRemoteRemove(_api; file=file, external_id=external_id, _mediaType=_mediaType)
    Swagger.exec(_ctx, response_stream)
end

function _swaggerinternal_filesRemoteShare(_api::FilesApi; file=nothing, external_id=nothing, channels=nothing, _mediaType=nothing)
    _ctx = Swagger.Ctx(_api.client, "GET", Dict{String, Any}, "/files.remote.share", ["slackAuth"])
    Swagger.set_param(_ctx.query, "file", file)  # type String
    Swagger.set_param(_ctx.query, "external_id", external_id)  # type String
    Swagger.set_param(_ctx.query, "channels", channels)  # type String
    Swagger.set_header_accept(_ctx, ["application/json"])
    Swagger.set_header_content_type(_ctx, (_mediaType === nothing) ? ["application/x-www-form-urlencoded"] : [_mediaType])
    return _ctx
end

"""

Share a remote file into a channel.
Param: file::String
Param: external_id::String
Param: channels::String
Return: Dict{String, Any}
"""
function filesRemoteShare(_api::FilesApi; file=nothing, external_id=nothing, channels=nothing, _mediaType=nothing)
    _ctx = _swaggerinternal_filesRemoteShare(_api; file=file, external_id=external_id, channels=channels, _mediaType=_mediaType)
    Swagger.exec(_ctx)
end

function filesRemoteShare(_api::FilesApi, response_stream::Channel; file=nothing, external_id=nothing, channels=nothing, _mediaType=nothing)
    _ctx = _swaggerinternal_filesRemoteShare(_api; file=file, external_id=external_id, channels=channels, _mediaType=_mediaType)
    Swagger.exec(_ctx, response_stream)
end

function _swaggerinternal_filesRemoteUpdate(_api::FilesApi; file=nothing, external_id=nothing, title=nothing, filetype=nothing, external_url=nothing, preview_image=nothing, indexable_file_contents=nothing, _mediaType=nothing)
    _ctx = Swagger.Ctx(_api.client, "POST", Dict{String, Any}, "/files.remote.update", ["slackAuth"])
    Swagger.set_param(_ctx.form, "file", file)  # type String
    Swagger.set_param(_ctx.form, "external_id", external_id)  # type String
    Swagger.set_param(_ctx.form, "title", title)  # type String
    Swagger.set_param(_ctx.form, "filetype", filetype)  # type String
    Swagger.set_param(_ctx.form, "external_url", external_url)  # type String
    Swagger.set_param(_ctx.form, "preview_image", preview_image)  # type String
    Swagger.set_param(_ctx.form, "indexable_file_contents", indexable_file_contents)  # type String
    Swagger.set_header_accept(_ctx, ["application/json"])
    Swagger.set_header_content_type(_ctx, (_mediaType === nothing) ? ["application/x-www-form-urlencoded"] : [_mediaType])
    return _ctx
end

"""

Updates an existing remote file.
Param: file::String
Param: external_id::String
Param: title::String
Param: filetype::String
Param: external_url::String
Param: preview_image::String
Param: indexable_file_contents::String
Return: Dict{String, Any}
"""
function filesRemoteUpdate(_api::FilesApi; file=nothing, external_id=nothing, title=nothing, filetype=nothing, external_url=nothing, preview_image=nothing, indexable_file_contents=nothing, _mediaType=nothing)
    _ctx = _swaggerinternal_filesRemoteUpdate(_api; file=file, external_id=external_id, title=title, filetype=filetype, external_url=external_url, preview_image=preview_image, indexable_file_contents=indexable_file_contents, _mediaType=_mediaType)
    Swagger.exec(_ctx)
end

function filesRemoteUpdate(_api::FilesApi, response_stream::Channel; file=nothing, external_id=nothing, title=nothing, filetype=nothing, external_url=nothing, preview_image=nothing, indexable_file_contents=nothing, _mediaType=nothing)
    _ctx = _swaggerinternal_filesRemoteUpdate(_api; file=file, external_id=external_id, title=title, filetype=filetype, external_url=external_url, preview_image=preview_image, indexable_file_contents=indexable_file_contents, _mediaType=_mediaType)
    Swagger.exec(_ctx, response_stream)
end

function _swaggerinternal_filesRevokePublicURL(_api::FilesApi; file=nothing, _mediaType=nothing)
    _ctx = Swagger.Ctx(_api.client, "POST", Dict{String, Any}, "/files.revokePublicURL", ["slackAuth"])
    Swagger.set_param(_ctx.form, "file", file)  # type String
    Swagger.set_header_accept(_ctx, ["application/json"])
    Swagger.set_header_content_type(_ctx, (_mediaType === nothing) ? ["application/x-www-form-urlencoded", "application/json"] : [_mediaType])
    return _ctx
end

"""

Revokes public/external sharing access for a file
Param: file::String
Return: Dict{String, Any}
"""
function filesRevokePublicURL(_api::FilesApi; file=nothing, _mediaType=nothing)
    _ctx = _swaggerinternal_filesRevokePublicURL(_api; file=file, _mediaType=_mediaType)
    Swagger.exec(_ctx)
end

function filesRevokePublicURL(_api::FilesApi, response_stream::Channel; file=nothing, _mediaType=nothing)
    _ctx = _swaggerinternal_filesRevokePublicURL(_api; file=file, _mediaType=_mediaType)
    Swagger.exec(_ctx, response_stream)
end

function _swaggerinternal_filesSharedPublicURL(_api::FilesApi; file=nothing, _mediaType=nothing)
    _ctx = Swagger.Ctx(_api.client, "POST", Dict{String, Any}, "/files.sharedPublicURL", ["slackAuth"])
    Swagger.set_param(_ctx.form, "file", file)  # type String
    Swagger.set_header_accept(_ctx, ["application/json"])
    Swagger.set_header_content_type(_ctx, (_mediaType === nothing) ? ["application/x-www-form-urlencoded", "application/json"] : [_mediaType])
    return _ctx
end

"""

Enables a file for public/external sharing.
Param: file::String
Return: Dict{String, Any}
"""
function filesSharedPublicURL(_api::FilesApi; file=nothing, _mediaType=nothing)
    _ctx = _swaggerinternal_filesSharedPublicURL(_api; file=file, _mediaType=_mediaType)
    Swagger.exec(_ctx)
end

function filesSharedPublicURL(_api::FilesApi, response_stream::Channel; file=nothing, _mediaType=nothing)
    _ctx = _swaggerinternal_filesSharedPublicURL(_api; file=file, _mediaType=_mediaType)
    Swagger.exec(_ctx, response_stream)
end

function _swaggerinternal_filesUpload(_api::FilesApi; file=nothing, content=nothing, filetype=nothing, filename=nothing, title=nothing, initial_comment=nothing, channels=nothing, thread_ts=nothing, _mediaType=nothing)
    _ctx = Swagger.Ctx(_api.client, "POST", Dict{String, Any}, "/files.upload", ["slackAuth"])
    Swagger.set_param(_ctx.form, "file", file)  # type String
    Swagger.set_param(_ctx.form, "content", content)  # type String
    Swagger.set_param(_ctx.form, "filetype", filetype)  # type String
    Swagger.set_param(_ctx.form, "filename", filename)  # type String
    Swagger.set_param(_ctx.form, "title", title)  # type String
    Swagger.set_param(_ctx.form, "initial_comment", initial_comment)  # type String
    Swagger.set_param(_ctx.form, "channels", channels)  # type String
    Swagger.set_param(_ctx.form, "thread_ts", thread_ts)  # type Float32
    Swagger.set_header_accept(_ctx, ["application/json"])
    Swagger.set_header_content_type(_ctx, (_mediaType === nothing) ? ["application/x-www-form-urlencoded"] : [_mediaType])
    return _ctx
end

"""

Uploads or creates a file.
Param: file::String
Param: content::String
Param: filetype::String
Param: filename::String
Param: title::String
Param: initial_comment::String
Param: channels::String
Param: thread_ts::Float32
Return: Dict{String, Any}
"""
function filesUpload(_api::FilesApi; file=nothing, content=nothing, filetype=nothing, filename=nothing, title=nothing, initial_comment=nothing, channels=nothing, thread_ts=nothing, _mediaType=nothing)
    _ctx = _swaggerinternal_filesUpload(_api; file=file, content=content, filetype=filetype, filename=filename, title=title, initial_comment=initial_comment, channels=channels, thread_ts=thread_ts, _mediaType=_mediaType)
    Swagger.exec(_ctx)
end

function filesUpload(_api::FilesApi, response_stream::Channel; file=nothing, content=nothing, filetype=nothing, filename=nothing, title=nothing, initial_comment=nothing, channels=nothing, thread_ts=nothing, _mediaType=nothing)
    _ctx = _swaggerinternal_filesUpload(_api; file=file, content=content, filetype=filetype, filename=filename, title=title, initial_comment=initial_comment, channels=channels, thread_ts=thread_ts, _mediaType=_mediaType)
    Swagger.exec(_ctx, response_stream)
end

export filesCommentsDelete, filesDelete, filesInfo, filesList, filesRemoteAdd, filesRemoteInfo, filesRemoteList, filesRemoteRemove, filesRemoteShare, filesRemoteUpdate, filesRevokePublicURL, filesSharedPublicURL, filesUpload


text:
struct FillArray{T,N} <: AbstractArray{T,N}
  value::T
  size::NTuple{N,Int}
end

Base.size(xs::FillArray) = xs.size

Base.getindex(xs::FillArray, ::Int...) = xs.value


text:
function solve_minlp(ROPF, fixing, index_var, nb_max_shunts, typeofconstraint::T) where T<: Type{NoConstraint}
    instance = ROPF.instance_name
    root = pwd()
    # date = Dates.format(now(), "yy_u_dd_HH_MM_SS")
    isdir("knitro_runs") || mkpath("knitro_runs")
    if length(fixing) > 0
        outlog = joinpath(pwd(), "knitro_runs", "BB_$(instance).log")
    else
        outlog = joinpath(pwd(), "knitro_runs", "$(instance).log")
    end
    pb_path = ROPF.output_instance_path
    src_ampl_path = joinpath(pwd(), "src_ampl")

    cd(pb_path)

    instance_dat_file = "$(instance).dat"
    cp(instance_dat_file, "minlp_instance.dat", force=true)

    f = open("fixing.dat", "w")
    all_fixing = true
    for (var, index) in index_var
        value = fixing[index]
        if value != -1
            write(f, "$var    $value \n")
        else
            all_fixing = false
        end
    end
    close(f)
    if length(fixing) == 0
        all_fixing = false
    end

    if all_fixing
        open("phase3.run", "w") do f
            println(f, "include $(joinpath(src_ampl_path, "phase3.run"));")
        end
    else
        open("minlp.run", "w") do f
            println(f, "include $(joinpath(src_ampl_path, "minlp.run"));")
        end
    end

    open("minlp.mod", "w") do f
      println(f, "include $(joinpath(src_ampl_path, "minlp.mod"));")
    end

    try
        if all_fixing
            run(`cmd /c ampl phase3.run '>' $(outlog)`)
        else
            run(`cmd /c ampl minlp.run '>' $(outlog)`)
        end
        # run(`cmd /c ampl real_minlp.run `)
    catch
        @warn("AMPL/Knitro failed, returning.")
    end

    mv("knitro_solution.csv", joinpath(root, "knitro_solution.csv"), force=true)
    # mv("knitro_solution.csv", "D:\\repo\\ROPF.jl\\knitro_optimal_solutions\\$(instance).csv", force=true)
    rm("minlp_instance.dat")
    rm("fixing.dat")
    cd(root)
    return read_obj_outlog(outlog)
end

function read_obj_outlog(outlog)
    objective_value = 0.0
    status = ""
    time = 0.0
    lines = readlines(outlog)
    for line in lines
        splitted_line = split(line, ":")
        if splitted_line[1] == "EXIT"
            status = splitted_line[2]
        end
        splitted_line = split(line, "=")
        if splitted_line[1] == "Final objective value               "
            objective_value = parse(Float64,splitted_line[end])
        elseif splitted_line[1] == "Total program time (secs)           "
            t = split(splitted_line[2])[1]
            time += parse(Float64,t)
        end
    end

    if status == " Locally optimal solution found." || status == " Optimal solution found."
        return objective_value
    else
        return +Inf
    end
end

function solve_minlp(ROPF, fixing, index_var, nb_max_shunts, typeofconstraint::T) where T<: Type{MAXkshuntsConstraint}
    instance = ROPF.instance_name
    root = pwd()
    # date = Dates.format(now(), "yy_u_dd_HH_MM_SS")
    isdir("knitro_runs") || mkpath("knitro_runs")
    if length(fixing) > 0
        outlog = joinpath(pwd(), "knitro_runs", "MAX$(nb_max_shunts)shunts_BB_$(instance).log")
    else
        outlog = joinpath(pwd(), "knitro_runs", "MAX$(nb_max_shunts)shunts_$(instance).log")
    end
    pb_path = ROPF.output_instance_path
    src_ampl_path = joinpath(pwd(), "src_ampl")

    cd(pb_path)

    instance_dat_file = "$(instance).dat"
    cp(instance_dat_file, "minlp_instance.dat", force=true)

    f = open("k.dat", "w")
    write(f, "$nb_max_shunts")
    close(f)

    f = open("fixing.dat", "w")
    all_fixing = true
    for (var, index) in index_var
        value = fixing[index]
        if value != -1
            write(f, "$var    $value \n")
        else
            all_fixing = false
        end
    end
    close(f)
    if length(fixing) == 0
        all_fixing = false
    end

    if all_fixing
        open("phase3_maxk.run", "w") do f
            println(f, "include $(joinpath(src_ampl_path, "phase3_maxk.run"));")
        end
    else
        open("minlp_maxk.run", "w") do f
            println(f, "include $(joinpath(src_ampl_path, "minlp_maxk.run"));")
        end
    end

    open("minlp_maxk.mod", "w") do f
      println(f, "include $(joinpath(src_ampl_path, "minlp_maxk.mod"));")
    end

    try
        if all_fixing
            run(`cmd /c ampl phase3_maxk.run '>' $(outlog)`)
        else
            run(`cmd /c ampl minlp_maxk.run '>' $(outlog)`)
        end
        # run(`cmd /c ampl real_minlp.run `)
    catch
        @warn("AMPL/Knitro failed, returning.")
    end

    mv("knitro_solution.csv", joinpath(root, "knitro_solution.csv"), force=true)
    # mv("knitro_solution.csv", "D:\\repo\\ROPF.jl\\knitro_optimal_solutions\\$(instance).csv", force=true)
    rm("minlp_instance.dat")
    rm("fixing.dat")
    cd(root)
    return read_obj_outlog(outlog)
end

function solve_minlp(ROPF, fixing, index_var, nb_max_shunts, typeofconstraint::T) where T<: Type{MAXkmovesConstraint}
    instance = ROPF.instance_name
    root = pwd()
    # date = Dates.format(now(), "yy_u_dd_HH_MM_SS")
    isdir("knitro_runs") || mkpath("knitro_runs")
    if length(fixing) > 0
        outlog = joinpath(pwd(), "knitro_runs", "MAX$(nb_max_shunts)moves_BB_$(instance).log")
    else
        outlog = joinpath(pwd(), "knitro_runs", "MAX$(nb_max_shunts)moves_$(instance).log")
    end
    pb_path = ROPF.output_instance_path
    src_ampl_path = joinpath(pwd(), "src_ampl")

    cd(pb_path)

    instance_dat_file = "$(instance).dat"
    cp(instance_dat_file, "minlp_instance.dat", force=true)

    f = open("k.dat", "w")
    write(f, "$nb_max_shunts")
    close(f)

    cp(joinpath(root,"initial_state\\$(instance).dat"), "init_pt.dat", force=true)

    f = open("fixing.dat", "w")
    all_fixing = true
    for (var, index) in index_var
        value = fixing[index]
        if value != -1
            write(f, "$var    $value \n")
        else
            all_fixing = false
        end
    end
    close(f)
    if length(fixing) == 0
        all_fixing = false
    end

    if all_fixing
        open("phase3_maxkmoves.run", "w") do f
            println(f, "include $(joinpath(src_ampl_path, "phase3_maxkmoves.run"));")
        end
    else
        open("minlp_maxkmoves.run", "w") do f
            println(f, "include $(joinpath(src_ampl_path, "minlp_maxkmoves.run"));")
        end
    end

    open("minlp_maxkmoves.mod", "w") do f
      println(f, "include $(joinpath(src_ampl_path, "minlp_maxkmoves.mod"));")
    end

    try
        if all_fixing
            run(`cmd /c ampl phase3_maxkmoves.run '>' $(outlog)`)
        else
            run(`cmd /c ampl minlp_maxkmoves.run '>' $(outlog)`)
        end
        # run(`cmd /c ampl real_minlp.run `)
    catch
        @warn("AMPL/Knitro failed, returning.")
    end

    mv("knitro_solution.csv", joinpath(root, "knitro_solution.csv"), force=true)
    # mv("knitro_solution.csv", "D:\\repo\\ROPF.jl\\knitro_optimal_solutions\\$(instance).csv", force=true)
    rm("minlp_instance.dat")
    rm("fixing.dat")
    cd(root)
    return read_obj_outlog(outlog)
end


text:
# This file was generated by the Julia Swagger Code Generator
# Do not modify this file directly. Modify the swagger specification instead.

mutable struct IoK8sApiNetworkingV1NetworkPolicyEgressRule <: SwaggerModel
    ports::Any # spec type: Union{ Nothing, Vector{IoK8sApiNetworkingV1NetworkPolicyPort} } # spec name: ports
    to::Any # spec type: Union{ Nothing, Vector{IoK8sApiNetworkingV1NetworkPolicyPeer} } # spec name: to

    function IoK8sApiNetworkingV1NetworkPolicyEgressRule(;ports=nothing, to=nothing)
        o = new()
        validate_property(IoK8sApiNetworkingV1NetworkPolicyEgressRule, Symbol("ports"), ports)
        setfield!(o, Symbol("ports"), ports)
        validate_property(IoK8sApiNetworkingV1NetworkPolicyEgressRule, Symbol("to"), to)
        setfield!(o, Symbol("to"), to)
        o
    end
end # type IoK8sApiNetworkingV1NetworkPolicyEgressRule

const _property_map_IoK8sApiNetworkingV1NetworkPolicyEgressRule = Dict{Symbol,Symbol}(Symbol("ports")=>Symbol("ports"), Symbol("to")=>Symbol("to"))
const _property_types_IoK8sApiNetworkingV1NetworkPolicyEgressRule = Dict{Symbol,String}(Symbol("ports")=>"Vector{IoK8sApiNetworkingV1NetworkPolicyPort}", Symbol("to")=>"Vector{IoK8sApiNetworkingV1NetworkPolicyPeer}")
Base.propertynames(::Type{ IoK8sApiNetworkingV1NetworkPolicyEgressRule }) = collect(keys(_property_map_IoK8sApiNetworkingV1NetworkPolicyEgressRule))
Swagger.property_type(::Type{ IoK8sApiNetworkingV1NetworkPolicyEgressRule }, name::Symbol) = Union{Nothing,eval(Meta.parse(_property_types_IoK8sApiNetworkingV1NetworkPolicyEgressRule[name]))}
Swagger.field_name(::Type{ IoK8sApiNetworkingV1NetworkPolicyEgressRule }, property_name::Symbol) =  _property_map_IoK8sApiNetworkingV1NetworkPolicyEgressRule[property_name]

function check_required(o::IoK8sApiNetworkingV1NetworkPolicyEgressRule)
    true
end

function validate_property(::Type{ IoK8sApiNetworkingV1NetworkPolicyEgressRule }, name::Symbol, val)
end


text:
using GPUArrays

# CuArrays#207
batched_tr(A::CuArray{T, 3}) where T = batched_tr!(similar(A, size(A, 3)), A)

function batched_tr!(B::CuArray{T, 1}, A::CuArray{T, 3}) where T
    gpu_call(B, (A, B)) do state, A, B
        batch = @linearidx(B)
        s = zero(eltype(B))
        for i in axes(A, 1)
            s += A[i, i, batch]
        end
        B[batch] = s
        return
    end
    return B
end


text:
#=
--------------------------------------------------------------------------------------------------------------------------------
Base: time series
--------------------------------------------------------------------------------------------------------------------------------
=#

"""
    companion_form(Ψ::Array{Float64,2}, Σ::SymMatrix)

Construct the companion form parameters of a VAR(p) with coefficients Ψ and var-cov matrix of the residuals Σ.
"""
function companion_form(Ψ::Array{Float64,2}, Σ::SymMatrix)

    # Dimensions
    n = size(Σ,2);
    p = Int64(size(Ψ,2)/n);
    np_1 = n*(p-1);

    # Companion form VAR(p)
    C = [Ψ; Matrix(I, np_1, np_1) zeros(np_1, n)];
    V = Symmetric([Σ zeros(n, np_1); zeros(np_1, n*p)])::SymMatrix;

    # Return output
    return C, V;
end

"""
    ext_companion_form(Ψ::Array{Float64,2}, Σ::SymMatrix)

Construct the companion form parameters of a VAR(p) with coefficients Ψ and var-cov matrix of the residuals Σ. The companion form is extend with additional n entries.
"""
function ext_companion_form(Ψ::Array{Float64,2}, Σ::SymMatrix)

    # Dimensions
    n = size(Σ,2);
    p = Int64(size(Ψ,2)/n);
    np = n*p;

    # Companion form VAR(p)
    C = [Ψ zeros(n, n); Matrix(I, np, np) zeros(np, n)];
    V = Symmetric([Σ zeros(n, np); zeros(np, np+n)])::SymMatrix;

    # Return output
    return C, V;
end

#=
-------------------------------------------------------------------------------------------------------------------------------
Combinatorics and probability
-------------------------------------------------------------------------------------------------------------------------------
=#

"""
    no_combinations(n::Int64, k::Int64)

Compute the binomial coefficient of `n` observations and `k` groups, for big integers.

# Examples
```jldoctest
julia> no_combinations(1000000,100000)
7.333191945934207610471288280331309569215030711272858517142085449641265002716664e+141178
```
"""
no_combinations(n::Int64, k::Int64) = factorial(big(n))/(factorial(big(k))*factorial(big(n-k)));

"""
    rand_without_replacement(nT::Int64, d::Int64)

Draw `length(P)-d` elements from the positional vector `P` without replacement.
`P` is permanently changed in the process.

rand_without_replacement(n::Int64, T::Int64, d::Int64)

Draw `length(P)-d` elements from the positional vector `P` without replacement.
In the sampling process, no more than n-1 elements are removed for each point in time.
`P` is permanently changed in the process.

# Examples
```jldoctest
julia> rand_without_replacement(20, 5)
15-element Array{Int64,1}:
  1
  2
  3
  5
  7
  8
 10
 11
 13
 14
 16
 17
 18
 19
 20
```
"""
function rand_without_replacement(nT::Int64, d::Int64)

    # Positional vector
    P = collect(1:nT);

    # Draw without replacement d times
    for i=1:d
        deleteat!(P, findall(P.==rand(P)));
    end

    # Return output
    return setdiff(1:nT, P);
end

function rand_without_replacement(n::Int64, T::Int64, d::Int64)

    # Positional vector
    P = collect(1:n*T);

    # Full set of coordinates
    coord = [repeat(1:n, T) kron(1:T, convert(Array{Int64}, ones(n)))];

    # Counter
    coord_counter = convert(Array{Int64}, zeros(T));

    # Loop over d
    for i=1:d

        while true

            # New candidate draw
            draw = rand(P);
            coord_draw = @view coord[draw, :];

            # Accept the draw if all observations are not missing for time t = coord[draw, :][2]
            if coord_counter[coord_draw[2]] < n-1
                coord_counter[coord_draw[2]] += 1;

                # Draw without replacement
                deleteat!(P, findall(P.==draw));
                break;
            end
        end
    end

    # Return output
    return setdiff(1:n*T, P);
end


text:
export NLDS, updatemaxncuts!

# π, σ and ρ do not really make sense alone so only
# their product will T, h, d, e is stored
mutable struct Solution <: SOI.AbstractSolution
    status::Symbol
    objval
    objvalx
    objvalxuray
    x
    xuray # unbouded ray
    θ
    θuray # unbouded ray
    πT
    πh
    σd
    ρe
    function Solution(status::Symbol, objval)
        new(status, objval, nothing, nothing, nothing, nothing, nothing, nothing, nothing, nothing, nothing, nothing)
    end
end

# Feasibility cut
# D = π T
# d = π h + σ d
function SOI.feasibility_cut(sol::Solution)
    (sol.πT, sol.πh + sol.σd)
end

# Optimality cut
# E = π T
# e = π h + ρ e + σ d
function SOI.optimality_cut(sol::Solution)
    (sol.πT, sol.πh + sol.σd + sol.ρe)
end

SOI.getstatus(sol::Solution) = sol.status

SOI.getobjectivevalue(sol::Solution) = sol.objval

SOI.getnodeobjectivevalue(sol::Solution) = sol.objvalx

SOI.getnodevalue(sol::Solution) = sol.x

function SOI.getbellmanvalue(sol::Solution, i)
    sol.θ[i]
end

# Nested L-Shaped Decomposition Subproblem (NLDS)

# Primal
# min c x
#     + θ                         -> AvgCutGenerator
#     + sum θ_i                   -> MultiCutGenerator
#     + λ c_a x_a                 -> parent unbounded
#     h -   Tx_a - Wx ∈ K         -> parent bounded
#     h - λ Tx_a - Wx ∈ K         -> parent unbounded
#     Dx >= d                     (feasibility cuts)
#     Ex + θ   >= e               -> AvgCutGenerator
#     Ex + θ_i >= e               -> MultiCutGenerator
#     x in C
#     λ >= 0

# Note that we do not have θ >= 0

# Dual
# max π h + σ d + ρ e
#     - π Tx_a                    -> parent bounded
#     c - (π W + σ D + ρ E) ∈ C^* -> x dual
#     c_a x_a - π Tx_a >= 0       -> parent unbounded
#     π in K^*
#     1'ρ = 1 (or <= 1 if θ >= 0)
#     σ >= 0
#     ρ >= 0
mutable struct NLDS{S}
    W::AbstractMatrix{S}
    h::AbstractVector{S}
    T::AbstractMatrix{S}
    K
    C
    c::AbstractVector{S}

    # parent solution
    x_a::AbstractVector{S}
    xuray_a::Union{Nothing, AbstractVector{S}}
    objvalxuray_a::Union{Nothing, S}

    childFC::Vector{CutStore{S}}
    childOC::Vector{CutStore{S}}
    localFC::CutStore{S}
    localOC::CutStore{S}
    proba::Vector{Float64}
    θfree::BitSet
    θlb::Vector{Float64}
    childT::Union{Nothing, Vector{AbstractMatrix{S}}}
    cutgen::AbstractOptimalityCutGenerator

    nx::Int
    nθ::Int
    nπ::Int
    πs::Vector{Int}
    # Number of feasibility cuts
    nσ::Int
    # Location of each feasibility cut in the constraint matrix
    σs::Vector{Int}
    # Number of Optimality Cuts for each θ_i
    nρ::Vector{Int}
    # Location of each optimality cut in the constraint matrix
    ρs::Vector{Vector{Int}}

    model
    loaded
    solved
    prevsol::Union{Nothing, Solution}

    newcut::Symbol
    pruningalgo::AbstractCutPruningAlgo
    FCpruner::AbstractCutPruner
    OCpruners::Vector

    function NLDS{S}(W::AbstractMatrix{S}, h::AbstractVector{S}, T::AbstractMatrix{S}, K, C, c::AbstractVector{S}, solver, pruningalgo::AbstractCutPruningAlgo, newcut::Symbol=:AddImmediately) where S
        nx = size(W, 2)
        nθ = 0
        nπ = length(h)
        localFC = CutStore{S}(nx)
        localOC = CutStore{S}(nx)
        FCpruner = CutPruner{nx, S}(pruningalgo, :≥)
        OCpruners = typeof(FCpruner)[]
        if false
            model = MathProgBase.ConicModel(solver)
        else
            model = MathProgBase.LinearQuadraticModel(solver)
        end
        nlds = new{S}(W, h, T, K, C, c, S[], nothing, nothing, CutStore{S}[], CutStore{S}[], localFC, localOC, Float64[], BitSet(), Float64[], nothing, AvgCutGenerator(), nx, nθ, nπ, 1:nπ, 0, Int[], Int[], Vector{Int}[], model, false, false, nothing, newcut, pruningalgo, FCpruner, OCpruners)
        addfollower(localFC, (nlds, (:Feasibility, 0)))
        addfollower(localOC, (nlds, (:Optimality, 1)))
        nlds
    end
end

function NLDS{S}(W::AbstractMatrix, h::AbstractVector, T::AbstractMatrix, K, C, c::AbstractVector, solver, pruningalgo::AbstractCutPruningAlgo, newcut::Symbol=:AddImmediately) where S
    NLDS{S}(AbstractMatrix{S}(W), AbstractVector{S}(h), AbstractMatrix{S}(T), K, C, AbstractVector{S}(c), solver, pruningalgo, newcut)
end

function add_childT!(nlds, childT)
    if childT === nothing
        @assert nlds.childT === nothing
    else
        # If there isn't any child yet, node.childT is null
        if nlds.childT === nothing
            nlds.childT = [childT]
        else
            push!(nlds.childT, childT)
        end
    end
end

function add_scenario_transition!(nlds::NLDS{S}, childFC, childOC, proba, childT) where {S}
    push!(nlds.proba, proba)
    oldnθ = nlds.nθ
    nlds.nθ = nθ(nlds.cutgen, nlds.proba)
    Δθ = nlds.nθ - oldnθ
    push!(nlds.θlb, 0)
    push!(nlds.θfree, length(nlds.θlb))
    if !iszero(Δθ)
        @assert Δθ == 1
        push!(nlds.nρ, 0)
        push!(nlds.ρs, Int[])
        # θ ≧ max β - ⟨a, x⟩
        # so we give :Max and set lazy_minus to true because it is "- ⟨a, x⟩" and not "+ ⟨a, x⟩"
        push!(nlds.OCpruners, CutPruner{nlds.nx, S}(nlds.pruningalgo, :Max, true))
    end
    @assert length(childFC.b) == 0
    addfollower(childFC, (nlds, (:Feasibility, length(nlds.childFC)+1)))
    push!(nlds.childFC, childFC)
    @assert length(childOC.b) == 0
    addfollower(childOC, (nlds, (:Optimality, length(nlds.childOC)+1)))
    push!(nlds.childOC, childOC)
    add_childT!(nlds, childT)
end

function setprobability!(nlds::NLDS, i, proba)
    nlds.proba[i] = proba
end

function getobjectivebound(nlds::NLDS)
    if !isempty(nlds.θfree)
        return -Inf
    end
    for (cone, idxs) in nlds.C
        for i in idxs
            if (cone == :NonNeg && nlds.c[i]  < 0) ||
               (cone == :NonPos && nlds.c[i]  > 0) ||
               (cone == :Free   && nlds.c[i] != 0)
               return -Inf
           end
        end
    end
    return Eθlb(nlds)
end
function setθbound!(nlds::NLDS, i, θlb)
    if isfinite(θlb)
        delete!(nlds.θfree, i)
        nlds.θlb[i] = θlb
    end
end
# Returns the lower bounds and cones for θ
function θC(nlds::NLDS)
    if nlds.nθ == length(nlds.θlb)
        # Multi cut
        θlb = nlds.θlb
        if length(nlds.θfree) < nlds.nθ
            bounded = collect(setdiff(BitSet(1:nlds.nθ), nlds.θfree))
            θC = [(:NonNeg, collect(nlds.nx .+ bounded))]
            if !isempty(nlds.θfree)
                free = collect(nlds.θfree)
                push!(nlds.θC, (:Free, collect(nlds.nx .+ free)))
            end
        else
            θC = [(:Free, collect(nlds.nx .+ (1:nlds.nθ)))]
        end
    elseif nlds.nθ == 1
        # Averaged cut
        if isempty(nlds.θfree)
            θlb = [dot(nlds.proba, nlds.θlb)]
            θC = [(:NonNeg, [nlds.nx + 1])]
        else
            θlb = [0.]
            θC = [(:Free, [nlds.nx + 1])]
        end
    else
        θlb = Float64[]
        θC = Tuple{Symbol, Vector{Int}}[]
    end
    θlb, θC
end

function updatemaxncuts!(nlds::NLDS, maxncuts, i)
    if i == 0
        nlds.FCpruner.maxncuts = maxncuts
    else
        nlds.OCpruners[i].maxncuts = maxncuts
    end
end

# .=== doesn't work :(
function veceqeqeq(v::Vector, x)
    map(el->el === x, v)
end

function padfcut(nlds::NLDS, D::AbstractMatrix{S}) where S
    if nlds.nθ > 0
        [D spzeros(size(D, 1), nlds.nθ)]
    else
        D
    end
end

function getfeasibilitycuts(nlds::NLDS)
    function f(i)
        D = nlds.childFC[i].A
        if nlds.childT !== nothing
            D = D * nlds.childT[i]
        end
        D
    end
    cuts_D = Compat.reduce(vcat, map(i -> f(i), 1:length(nlds.childFC)), init=nlds.localFC.A)
    cuts_d = Compat.reduce(vcat, map(x -> x.b, nlds.childFC), init=nlds.localFC.b)
    mycut = Compat.reduce(vcat, map(x -> veceqeqeq(x.authors, nlds), nlds.childFC), init=veceqeqeq(nlds.localFC.authors, nlds))
    (cuts_D, cuts_d, mycut)
end

# see https://github.com/JuliaLang/julia/issues/16661
function padmultiocutaux(nlds::NLDS, E::AbstractMatrix{S}, i, onesvec) where S
    nrows = size(E, 1)
    [E spzeros(S, nrows, i-1) onesvec spzeros(S, nrows, nlds.nθ-i)]
end
function padmultiocut(nlds::NLDS, E::AbstractMatrix{S}, i) where S
    nrows = size(E, 1)
    padmultiocutaux(nlds, E, i, ones(S, nrows, 1))
end
function padmultiocut(nlds::NLDS, E::AbstractSparseMatrix{S}, i) where S
    nrows = size(E, 1)
    padmultiocutaux(nlds, E, i, sparse(ones(S, nrows, 1)))
end

function padavgocut(nlds::NLDS, E::AbstractMatrix{S}) where S
    nrows = size(E, 1)
    [E ones(S, nrows)]
end
function padavgocut(nlds::NLDS, E::AbstractSparseMatrix{S}) where S
    nrows = size(E, 1)
    [E sparse(ones(S, nrows))]
end

function getoptimalitycuts(nlds::NLDS{S}) where S
    function f(i)
        E = nlds.childOC[i].A
        if nlds.childT !== nothing
            E = E * nlds.childT[i]
        end
        padmultiocut(nlds, E, i)
    end
    if nlds.nθ == 1
        cuts_E = padavgocut(nlds, nlds.localOC.A)
    else
        cuts_E = spzeros(S, 0, nlds.nx + nlds.nθ)
    end
    if nlds.nθ == length(nlds.childOC)
        cuts_E = Compat.reduce(vcat, map(f, 1:length(nlds.childOC)), init=cuts_E)
    end
    cuts_e = Compat.reduce(vcat, map(x -> x.b, nlds.childOC), init=nlds.localOC.b)
    mycut = Compat.reduce(vcat, map(x -> veceqeqeq(x.authors, nlds), nlds.childOC), init=veceqeqeq(nlds.localOC.authors, nlds))
    (cuts_E, cuts_e, mycut)
end

function notifynewcuts(nlds::NLDS{S}, A::AbstractMatrix{S}, b::AbstractVector{S}, attrs, authors::Vector{NLDS{S}}) where S
    @assert attrs[1] in [:Feasibility, :Optimality]
    isfc = attrs[1] == :Feasibility
    nnewcuts = size(A, 1)
    i = attrs[2]
    if i > 0 && nlds.childT !== nothing
        A = A * nlds.childT[i]
    end
    # No .=== :(
    mine = [authors[i] === nlds for i in 1:length(authors)]
    if isfc
        pruner = nlds.FCpruner
    else
        pruner = nlds.OCpruners[i]
    end
    ncur = ncuts(pruner)
    addstatus = addcuts!(pruner, A, b, mine)
    npushed = sum(addstatus .> ncur)
    @assert ncur + npushed == ncuts(pruner)
    cur = nlds.nσ + sum(nlds.nρ)
    if isfc
        nlds.nσ += npushed
        is = nlds.σs
    else
        nlds.nρ[i] += npushed
        is = nlds.ρs[i]
    end
    append!(is, zeros(Int, npushed))
    for j in 1:nnewcuts
        if addstatus[j] > ncur
            cur += 1
            is[addstatus[j]] = cur
        end
        if addstatus[j] > 0 && nlds.loaded
            if 1 <= addstatus[j] <= ncur || nlds.newcut == :InvalidateSolver
                # Need to update cut
                nlds.loaded = false
                nlds.solved = false
            elseif nlds.newcut == :AddImmediately
                idx = collect(1:nlds.nx)
                a = A[j,:]
                if !isfc
                    if i == 0
                        push!(idx, nlds.nx+1)
                    else
                        push!(idx, nlds.nx+i)
                    end
                    a = _veccat(a, one(S), true)
                end
                _addconstr!(nlds.model, idx, a, b[j], :NonPos)
                nlds.solved = false
            else
                error("Invalid newcut option $(nlds.newcut)")
            end
        end
    end
    @assert nlds.newcut == :InvalidateSolver || cur == nlds.nσ + sum(nlds.nρ)
    checkconsistency(nlds)
end

function checkconsistency(nlds)
    @assert length(nlds.πs) == nlds.nπ
    @assert length(nlds.σs) == nlds.nσ
    for i in 1:nlds.nθ
        @assert length(nlds.ρs[i]) == nlds.nρ[i]
    end
    ρs = Compat.reduce(append!, nlds.ρs, init=Int[])
    @assert sort([nlds.πs; nlds.nπ .+ nlds.σs; nlds.nπ .+ ρs]) == collect(1:(nlds.nπ + nlds.nσ + sum(nlds.nρ)))
end

function getrhs(nlds::NLDS{S}) where S
    checkconsistency(nlds)
    bs = [nlds.h - nlds.T * nlds.x_a]
    Ks = [nlds.K]
    Kcut = []
    cur = 0
    b = S[]
    if !isempty(nlds.FCpruner)
        append!(b, nlds.FCpruner.b)
        nlds.nσ = ncuts(nlds.FCpruner)
        nlds.σs = cur .+ (1:nlds.nσ)
        cur += ncuts(nlds.FCpruner)
        push!(Kcut, (:NonPos, nlds.σs))
    end
    for i in 1:nlds.nθ
        if !isempty(nlds.OCpruners[i])
            append!(b, nlds.OCpruners[i].b)
            nlds.nρ[i] = ncuts(nlds.OCpruners[i])
            nlds.ρs[i] = cur .+ (1:nlds.nρ[i])
            cur += ncuts(nlds.OCpruners[i])
            push!(Kcut, (:NonPos, nlds.ρs[i]))
        end
    end
    if !isempty(Kcut)
        push!(bs, b)
        push!(Ks, Kcut)
    end
    checkconsistency(nlds)
    bs, Ks
end

function setparentx(nlds::NLDS, x_a::AbstractVector, xuray_a, objvalxuray_a)
    nlds.x_a = x_a
    unbounded_a = xuray_a !== nothing
    if nlds.xuray_a !== nothing || unbounded_a
        # FIXME do better when delvars!, ... are available in MPB
        nlds.loaded = false
        nlds.solved = false
    end
    if unbounded_a
        nlds.xuray_a = xuray_a
        nlds.objvalxuray_a = objvalxuray_a
    elseif nlds.xuray_a !== nothing
        nlds.xuray_a = nothing
        nlds.objvalxuray_a = nothing
    end
    if nlds.loaded
        if unbounded_a
            nlds.loaded = false
        else
            bs, Ks = getrhs(nlds)
            _setconstrB!(nlds.model, bs, Ks)
        end
        nlds.solved = false
    end
end

function computecuts!(nlds::NLDS)
    nlds.nπ = length(nlds.h)
    nlds.πs = collect(1:nlds.nπ)
    cuts_D, cuts_d, mycut_d = getfeasibilitycuts(nlds)
    cuts_D = padfcut(nlds, cuts_D)
    cuts_E, cuts_e, mycut_e = getoptimalitycuts(nlds)
    # I will add the cuts in notifynewcuts! now
    noneedstored!(nlds.localFC, nlds)
    for store in nlds.childFC
        noneedstored!(store, nlds)
    end
    noneedstored!(nlds.localOC, nlds)
    for store in nlds.childOC
        noneedstored!(store, nlds)
    end
end

function getcutsDE(nlds::NLDS{S}) where S
    checkconsistency(nlds)
    nc = nlds.nσ + sum(nlds.nρ)
    A = spzeros(S, nc, nlds.nx + nlds.nθ)
    if !isempty(nlds.FCpruner)
        A[nlds.σs, 1:nlds.nx] = nlds.FCpruner.A
    end
    for i in 1:nlds.nθ
        if !isempty(nlds.OCpruners[i])
            A[nlds.ρs[i], 1:nlds.nx] = nlds.OCpruners[i].A
            A[nlds.ρs[i], nlds.nx .+ i] .= 1
        end
    end
    A
end

function load!(nlds::NLDS{S}) where S
    if !nlds.loaded
        bigA = nlds.W
        if nlds.nθ > 0
            bigA = [bigA spzeros(size(bigA, 1), nlds.nθ)]
        end
        if nlds.xuray_a !== nothing
            bigA = [bigA nlds.T * nlds.xuray_a]
        end

        # Needs to be done before getcutsDE
        bs, Ks = getrhs(nlds)

        computecuts!(nlds) # FIXME what is the use of this ?
        A = getcutsDE(nlds)
        if nlds.xuray_a !== nothing
            A = [A spzeros(size(A, 1), 1)]
        end
        bigA = [bigA; A]

        bigC = nlds.C
        bigc = nlds.c
        if nlds.nθ > 0
            bigC = [bigC; θC(nlds)[2]]
            if nlds.nθ == 1
                bigc = [bigc; 1]
            else
                bigc = [bigc; nlds.proba]
            end
        end
        if nlds.xuray_a !== nothing
            bigC = [bigC; (:NonNeg, [nlds.nx+nlds.nθ+1])]
            bigc = [bigc; nlds.objvalxuray_a]
        end

        _load!(nlds.model, bigc, bigA, bs, Ks, bigC)
        nlds.loaded = true
    end
end

function checkstatus(status::Symbol)
    if status == :Error
        error("The solver reported an error")
    elseif status == :UserLimit
        error("The solver reached iteration limit or timed out")
    end
end

function Eθlb(nlds::NLDS)
    if nlds.nθ == 0
        0.
    else
        if nlds.nθ == length(nlds.proba)
            dot(nlds.proba, nlds.θlb)
        else
            θlb, _ = θC(nlds)
            @assert length(θlb) == 1
            θlb[1]
        end
    end
end

function solve!(nlds::NLDS{S}) where S
    load!(nlds)
    if !nlds.solved
        MathProgBase.optimize!(nlds.model)
        status = MathProgBase.status(nlds.model)
        checkstatus(status)
        sol = Solution(status, _getobjval(nlds.model) + Eθlb(nlds))

        if status == :Unbounded
            uray = _getunboundedray(nlds.model)
            sol.xuray = uray[1:nlds.nx]
            sol.objvalxuray = dot(nlds.c, sol.xuray)
            sol.θuray = uray[nlds.nx .+ (1:nlds.nθ)]

            # See https://github.com/JuliaOpt/Gurobi.jl/issues/80
            c = MathProgBase.getobj(nlds.model)

            MathProgBase.setobj!(nlds.model, zero(c))
            MathProgBase.optimize!(nlds.model)
            newstatus = MathProgBase.status(nlds.model)
            @assert newstatus != :Unbounded # Now the objective is 0
            checkstatus(newstatus)
            if newstatus == :Infeasible
                # We discard unbounded ray, infeasibility is more important
                status = newstatus
                sol.xuray = nothing
                sol.objvalxuray = nothing
                sol.θuray = nothing
            else
                @assert newstatus == :Optimal
            end
            if newstatus == :Infeasible
                sol.status = :Infeasible
            end
        else
            newstatus = status
        end

        if newstatus == :Infeasible || status != :Unbounded
            # if infeasible dual + λ iray is dual feasible for any λ >= 0
            # here I just take iray to build the feasibility cut
            dual = _getdual(nlds.model)
            if isempty(dual)
                error("Dual vector returned by the solver is empty while the status is $status")
            end
            @assert length(dual) == nlds.nπ + nlds.nσ + sum(nlds.nρ)

            π = dual[nlds.πs]
            σρ = @view dual[nlds.nπ+1:end]

            σ = σρ[nlds.σs]
            addusage!(nlds.FCpruner, σ)
            sol.σd = dot(σ, nlds.FCpruner.b)

            sol.ρe = zero(S)
            for i in 1:nlds.nθ
                ρ = σρ[nlds.ρs[i]]
                addusage!(nlds.OCpruners[i], ρ)
                sol.ρe += dot(ρ, nlds.OCpruners[i].b)
            end

            sol.πT = vec(π' * nlds.T)
            sol.πh = dot(π, nlds.h)
        end

        # Needs to be done after since if status is unbounded I do a resolve
        if newstatus != :Infeasible
            primal = _getsolution(nlds.model)
            sol.x = primal[1:nlds.nx]
            addposition!(nlds.FCpruner, sol.x)
            for i in 1:nlds.nθ
                addposition!(nlds.OCpruners[i], sol.x)
            end
            sol.objvalx = dot(nlds.c, sol.x)
            sol.θ = θC(nlds)[1] .+ primal[nlds.nx .+ (1:nlds.nθ)]
        end

        if status == :Unbounded
            # It needs to be done *after* getsolution for some solver (e.g. CPLEX)
            MathProgBase.setobj!(nlds.model, c)
        end

        nlds.prevsol = sol
    end
end

function getsolution(nlds::NLDS)
    solve!(nlds)
    nlds.prevsol
end


text:
using Test
using QuantumLattices.Prerequisites.CompositeStructures

import QuantumLattices.Prerequisites.Traits: contentnames, getcontent

struct CT{S, N, T} <: CompositeNTuple{N, T}
    info::S
    contents::NTuple{N, T}
end
contentnames(::Type{<:CT}) = (:info, :contents)

@testset "CompositeTuple" begin
    @test contentnames(CompositeTuple) == (:contents,)

    t = CT("Info", (1, 2, 3, 4))
    @test contentnames(typeof(t)) == (:info, :contents)
    @test length(t) == 4
    @test length(typeof(t)) == 4
    @test eltype(t) == Int
    @test eltype(typeof(t)) == Int
    @test t == deepcopy(t)
    @test isequal(t, deepcopy(t))
    @test hash(t) == hash(("Info", (1, 2, 3, 4)))
    @test t[1] == 1
    @test t[end] == 4
    @test t[1:3] == CT("Info", (1, 2, 3))
    @test collect(t) == [1, 2, 3, 4]
    @test collect(t|>Iterators.reverse) == [4, 3, 2, 1]
    @test keys(t) == keys((1, 2, 3, 4))
    @test values(t) == values((1, 2, 3, 4))
    @test pairs(t)|>collect == pairs((1, 2, 3, 4))|>collect
    @test reverse(t) == CT("Info", (4, 3, 2, 1))
    @test convert(Tuple, t) == (1, 2, 3, 4)
end

struct CV{S, T} <: CompositeVector{T}
    info::S
    contents::Vector{T}
end
contentnames(::Type{<:CV}) = (:info, :contents)

@testset "CompositeVector" begin
    @test contentnames(CompositeVector) == (:contents,)

    v = CV("Info", [1, 3, 2, 4])
    @test contentnames(typeof(v)) == (:info, :contents)
    @test size(v) == (4,)
    @test size(v, 1) == 4
    @test length(v) == 4
    @test v == deepcopy(v)
    @test isequal(v, deepcopy(v))
    @test v[end] == 4
    @test v[2] == 3
    @test v[CartesianIndex(2)] == 3
    @test v[[1, 3, 4]] == CV("Info", [1, 2, 4])
    @test v[1:3] == CV("Info", [1, 3, 2])
    @test (v[1] = 10; v[2:3] = 11:12; v == CV("Info", [10, 11, 12, 4]))
    @test push!(v, 1) == CV("Info", [10, 11, 12, 4, 1])
    @test push!(v, 2, 3) == CV("Info", [10, 11, 12, 4, 1, 2, 3])
    @test pushfirst!(v, 20, 21) == CV("Info", [20, 21, 10, 11, 12, 4, 1, 2, 3])
    @test insert!(v, 2, 30) == CV("Info", [20, 30, 21, 10, 11, 12, 4, 1, 2, 3])
    @test append!(v, [0, -1]) == CV("Info", [20, 30, 21, 10, 11, 12, 4, 1, 2, 3, 0, -1])
    @test prepend!(v, [8, 9]) == CV("Info", [8, 9, 20, 30, 21, 10, 11, 12, 4, 1, 2, 3, 0, -1])
    @test (splice!(v, 2) == 9) && (v == CV("Info", [8, 20, 30, 21, 10, 11, 12, 4, 1, 2, 3, 0, -1]))
    @test (splice!(v, 1, 9) == 8) && (v == CV("Info", [9, 20, 30, 21, 10, 11, 12, 4, 1, 2, 3, 0, -1]))
    @test (splice!(v, 1:3) == CV("Info", [9, 20, 30])) && (v == CV("Info", [21, 10, 11, 12, 4, 1, 2, 3, 0, -1]))
    @test (splice!(v, 1:3, [8, 7, 6]) == CV("Info", [21, 10, 11])) && (v == CV("Info", [8, 7, 6, 12, 4, 1, 2, 3, 0, -1]))
    @test deleteat!(v, 4) == CV("Info", [8, 7, 6, 4, 1, 2, 3, 0, -1])
    @test deleteat!(v, [1, 2]) == CV("Info", [6, 4, 1, 2, 3, 0, -1])
    @test deleteat!(v, 5:6) == CV("Info", [6, 4, 1, 2, -1])
    @test (pop!(v) == -1) && (v == CV("Info", [6, 4, 1, 2]))
    @test (popfirst!(v) == 6) && (v == CV("Info", [4, 1, 2]))
    @test (empty!(v) == CV("Info", Int[])) && (v == CV("Info", Int[]))

    v = CV("Info", [1, 3, 2, 2, 4])
    @test empty(v) == CV("Info", Int[])
    @test collect(v) == [1, 3, 2, 2, 4]
    @test keys(v) == keys([1, 3, 2, 2, 4])
    @test values(v) == values([1, 3, 2, 2, 4])
    @test pairs(v) == pairs([1, 3, 2, 2, 4])
    @test convert(Vector, v) == [1, 3, 2, 2, 4]
    @test reverse(v) == CV("Info", [4, 2, 2, 3, 1])
    @test (sort(v) == CV("Info", [1, 2, 2, 3, 4])) && (v == CV("Info", [1, 3, 2, 2, 4]))
    @test (sort!(v) == CV("Info", [1, 2, 2, 3, 4])) && (v == CV("Info", [1, 2, 2, 3, 4]))
    @test (filter(<=(3), v) == CV("Info", [1, 2, 2, 3])) && (v == CV("Info", [1, 2, 2, 3, 4]))
    @test (filter!(<=(3), v) == CV("Info", [1, 2, 2, 3])) && (v == CV("Info", [1, 2, 2, 3]))
    @test findfirst(isequal(2), v) == 2
    @test findlast(isequal(2), v) == 3
    @test findall(isequal(2), v) == [2, 3]
end

struct CD{S, P, I} <: CompositeDict{P, I}
    info::S
    newcontents::Dict{P, I}
end
contentnames(::Type{<:CD}) = (:info, :contents)
getcontent(m::CD, ::Val{:contents}) = getfield(m, :newcontents)

@testset "CompositeDict" begin
    @test contentnames(CompositeDict) == (:contents,)

    d = CD("Info", Dict("a"=>1, "b"=>2))
    @test contentnames(typeof(d)) == (:info, :contents)
    @test d == deepcopy(d)
    @test isequal(d, deepcopy(d))
    @test isempty(d) == false
    @test length(d) == 2
    @test (haskey(d, "a") == true) && (haskey(d, "d") == false)
    @test (Pair("a", 1) ∈ d) && (Pair("d", 4) ∉ d)
    @test get(d, "a", 2) == 1
    @test get(d, "d", 4) == 4
    @test get(()->4, d, "d") == 4
    @test (get!(d, "d", 4) == 4) && (d == CD("Info", Dict("a"=>1, "b"=>2, "d"=>4)))
    @test (get!(()->5, d, "d") == 4) && (d == CD("Info", Dict("a"=>1, "b"=>2, "d"=>4)))
    @test getkey(d, "e", "e") == "e"
    @test d["d"] == 4
    @test (push!(d, Pair("e", 4)) == CD("Info", Dict("a"=>1, "b"=>2, "d"=>4, "e"=>4))) && (d == CD("Info", Dict("a"=>1, "b"=>2, "d"=>4, "e"=>4)))
    @test (d["d"] = 4; d["e"] = 5; d == CD("Info", Dict("a"=>1, "b"=>2, "d"=>4, "e"=>5)))
    @test (pop!(d) == Pair("e", 5)) && (d == CD("Info", Dict("a"=>1, "b"=>2, "d"=>4)))
    @test (pop!(d, "a") == 1) && (d == CD("Info", Dict("b"=>2, "d"=>4)))
    @test (pop!(d, "a", 1) == 1) && (d == CD("Info", Dict("b"=>2, "d"=>4)))
    @test (delete!(d, "b") == CD("Info", Dict("d"=>4))) && (d == CD("Info", Dict("d"=>4)))
    @test (empty!(d) == CD("Info", Dict{String, Int}())) && (d == CD("Info", Dict{String, Int}()))

    d = CD("Info", Dict("a"=>1, "b"=>2, "c"=>3, "d"=>4))
    @test merge(CD("Info", Dict("a"=>1, "b"=>2)), CD("Info", Dict("c"=>3, "d"=>4))) == d
    @test merge(+, CD("Info", Dict("a"=>1, "b"=>2, "c"=>1)), CD("Info", Dict("c"=>2, "d"=>4))) == d
    @test empty(d) == CD("Info", Dict{String, Int}())
    @test collect(d) == collect(Dict("a"=>1, "b"=>2, "c"=>3, "d"=>4))
    @test collect(keys(d)) == collect(keys(Dict("a"=>1, "b"=>2, "c"=>3, "d"=>4)))
    @test collect(values(d)) == collect(values(Dict("a"=>1, "b"=>2, "c"=>3, "d"=>4)))
    @test collect(pairs(d)) == collect(pairs(Dict("a"=>1, "b"=>2, "c"=>3, "d"=>4)))
    @test convert(Dict, d) == Dict("a"=>1, "b"=>2, "c"=>3, "d"=>4)
    @test (filter(p->p.second<=3, d) == CD("Info", Dict("a"=>1, "b"=>2, "c"=>3))) && (d == CD("Info", Dict("a"=>1, "b"=>2, "c"=>3, "d"=>4)))
    @test (filter!(p->p.second<=3, d) == CD("Info", Dict("a"=>1, "b"=>2, "c"=>3))) && (d == CD("Info", Dict("a"=>1, "b"=>2, "c"=>3)))
end

@testset "NamedContainer" begin
    @test NamedContainer{(:a, :b)}((1, 'h')) == (a=1, b='h')
    @test NamedContainer{()}(()) == NamedTuple()
end


text:
include("./computer.jl")
include("./utils.jl")
input_data = read_numbers("input11.txt")
using Test
using OffsetArrays

function paint_move!(painted, current_coordinates, color, new_direction)
  painted[current_coordinates.x, current_coordinates.y] = color
  new_coordinates = (x=current_coordinates.x + new_direction.x_dir, y=current_coordinates.y + new_direction.y_dir)
  return new_coordinates
end
function angle2direction(angle)
  x_dir = -Int(round(cos(angle)))
  y_dir = Int(round(sin(angle)))
  return (x_dir=x_dir, y_dir=y_dir)
end

function count_painted(arr)
  count_black = count(i -> (i == 0), arr)
  count_white = count(i -> (i == 1), arr)
  return count_black + count_white
end

function robot_controller(view_channel, control_channel, initial_color=0)
  dims = (100,100)
  offset_limits = (-50:49, -50:49)
  println(offset_limits)

  painted = OffsetArray(ones(Int8, dims...) * 255, offset_limits)
  current_coordinates = (x=0,y=0)
  current_angle = 0
  current_direction = angle2direction(current_angle)
  color_map = Dict(255 => 0, 0 => 0, 1 => 1)
  function turn_right!()
    current_angle = mod2pi(current_angle + π/2)
    current_direction = angle2direction(current_angle)
  end
  function turn_left!()
    current_angle = mod2pi(current_angle - π/2)
    current_direction = angle2direction(current_angle)
  end

  put!(view_channel, initial_color)
  while true
    @debug(current_coordinates)
    color = take!(control_channel)
    @debug("color: $color")
    if color == 255
      break
    end
    direction = take!(control_channel)
    @debug("direction: $direction")
    if direction == 0
      turn_left!()
    elseif direction == 1
      turn_right!()
    end
    current_coordinates = paint_move!(painted, current_coordinates, color, current_direction)
    @debug("current_coordinates: $current_coordinates")
    view_color = color_map[painted[current_coordinates.x, current_coordinates.y]]
    @debug("view_color: $view_color")
    put!(view_channel, view_color)
  end
  println("Painted pixels: $(count_painted(painted))")
  return painted
end
function robot_program(input_data, view_channel, control_channel)
    run_program_async(input_data, view_channel, control_channel)
    put!(control_channel, 255) #break
end

function robot_state_machine(input_data, initial_color=0)

  view_channel = Channel{Int64}(Inf)
  control_channel = Channel{Int64}(Inf)
  output = undef

  @sync begin
    @async robot_program(input_data, view_channel, control_channel)
    @async output = robot_controller(view_channel, control_channel, initial_color)
  end

  return output

end

function test_robot_program()

  input_data = [3,0,104,1,104,1, 99]
  robot_state_machine(input_data)

end
#run_program_async(input_data, view_channel, control_channel)
test_map = ones(Int8, 3,3) * 255
@test paint_move!(test_map, (x=2, y=2), 1, (x_dir=1, y_dir=0)) == (x=3, y=2)
@test test_map[2,2] == 1
println(test_map)
@test count_painted(test_map) == 1

#test_robot_program()
#display(robot_state_machine(input_data))
offset_array = robot_state_machine(input_data, 1)
normal_array = offset_array.parent
drawArray(normal_array)

max_stars_repo_name: buotex/aoc

max_issues_repo_name: buotex/aoc

max_forks_repo_name: buotex/aoc