text:
using QuasiArrays, ArrayLayouts, LinearAlgebra, Test
import QuasiArrays: apply

@testset "Multiplication" begin
    @testset "Diag * Inclusion" begin
        A = QuasiDiagonal(Inclusion(0:0.1:1))
        @test rowsupport(A,0.1) == colsupport(A,0.1) == 0.1
        b = Inclusion(0:0.1:1)
        Ab = A*b
        @test Ab[0.1] ≈ 0.1^2
        @test_throws DimensionMismatch A*Inclusion(1:2)

        @test A*A isa QuasiDiagonal
        @test (A*A)[0.1,0.1] ≈ A[0.1,0.1]^2
    end
    @testset "Quasi * Quasi" begin
        A = QuasiArray(rand(3,3),(0:0.5:1,0:0.5:1))
        @test *(A) == A
        @test A*A isa QuasiArray
        @test QuasiArray(A*A) == QuasiArray(A.parent*A.parent,A.axes)
    end
    @testset "Quasi * Array" begin
        A = QuasiArray(rand(3,3),(0:0.5:1,Base.OneTo(3)))
        B = rand(3,3)
        @test A*B isa QuasiArray
        @test QuasiArray(A*B) ≈ QuasiArray(A.parent*B,A.axes)
        @test B*A' isa QuasiArray
        @test (A*B)[0.5,1] ≈ (Array(A)*B)[2,1]
    end

    @testset "Triple" begin
        A = QuasiArray(rand(3,3),(0:0.5:1,0:0.5:1))
        @test A isa QuasiArray
        @test Array(A*A) ≈ Array(A)^2
        @test Array(A*A*A) ≈ Array(A)^3
    end

    @testset "absvec" begin
        a = QuasiArray(rand(3),(0:0.5:1,))
        c = QuasiArray(rand(3) .+ im .* randn(3),(0:0.5:1,))
        @test a'a ≈ transpose(a)a ≈ dot(a,a)
        @test c'c ≈ dot(c,c)
        @test transpose(c)c ≈ transpose(c.parent)*c.parent
    end
end


text:
# # Elipse

# ## Model

module elipse   #jl

using ReachabilityAnalysis

A = [3.0 -9.0;
     4.0 -3.0]

function model(X0)
    S = @system(x' = Ax)
    return IVP(S, X0)
end

end #jl

# ## References

ext: jl

max_stars_repo_name: mforets/ReachabilityLibrary

max_issues_repo_name: mforets/ReachabilityLibrary

max_forks_repo_name: mforets/ReachabilityLibrary


text:
"""
    text_output(id)

Render the output from the server as raw text.

This is best used with a formatting tag (such as `p()` or heading `h1()` etc)
"""
function text_output(id)
    (UIElement("""{{ $id }}"""), UIModel(id, "'Loading...'"))
end

"""
    plots_output(id)

Render a plot. Supports plots created by the `Plots` package.

Plots are slow to render, and so should not be directly hooked up to inputs that rapidly
request updates to the server -- like sliders. If you wish to use these input types, add a
button to refresh the plot input and only return the plot when that is `true` (otherwise
return `nothing`).
"""
function plots_output(id::AbstractString)
    UIElement("""<v-container><img v-bind:src = "$id"></v-container>"""), UIModel(id, "null")
end

"""
    datatable_output(id::AbstractString; rows_per_page::Integer

Render a table of data in your UI. The function `id` in your Server module should return a
table that is a subtype of either `AbstractDataFrame` or `Table` (from `Tables.jl`).
"""
function datatable_output(id::AbstractString; rows_per_page::Integer = 10)
    UIElement("""
    <v-data-table
    :headers="$id[0]"
    :items="$id.slice(1)"
    :items-per-page="$rows_per_page"
    no-data-text="No data"
    class="elevation-1"
  ></v-data-table>"""),
    UIModel(id, """[[{text: "Loading...", value: "loading"}], {loading: "Loading..."}]""")
end


text:
module ThinningAndShift

using Distributions
using Random

import LinearAlgebra: dot,diagind
import Statistics: mean,var,cov
import StatsBase: midpoints

abstract type Jittering end

#####
# no jittering
struct NoJitter  <: Jittering end
@inline function jitter!(v::Vector{Float64},::NoJitter)
  return v
end

####
# apply a distribution
struct JitterDistribution{D<:Distribution} <: Jittering
  d::D
end
@inline function jitter!(v::Vector{Float64},j::JitterDistribution{D}) where D<:UnivariateDistribution
  n = length(v)
  if n > 1
    jits = rand(j.d,n)
    v .+= jits
  end
  return v
end
@inline function jitter!(v::Vector{Float64},j::JitterDistribution)
  jits = rand(j.d)
  v .+= jits
  return v
end

# increment by given distribution
struct JitterIncremental{D<:UnivariateDistribution} <: Jittering
  d::D
end

@inline function jitter!(v::Vector{Float64},j::JitterIncremental)
  n = length(v)
  if n>1
    jits = cumsum(rand(j.d,n-1))
    v[2:end] .+= jits
  end
  return v
end

####
## only for sets of two, shifts one spike time w.r.t. the other
struct JitterPaired{D<:UnivariateDistribution} <: Jittering
  d::D
end
@inline function jitter!(v::Vector{Float64},j::JitterPaired)
  ji = rand(j.d)
  v[2] += ji
  return v
end

struct GTAS{R,I,J<:Jittering}
  n::I
  rate_parent::R
  markings::Vector{Vector{I}}
  marking_selector::Categorical{R}
  jitters::Vector{J}
  function GTAS(rate_parent::R,markings::Vector{Vector{I}},
      markings_probs::Vector{R},jitters::Vector{J}) where {R<:Real,I<:Integer,J<:Jittering}
    @assert sum(markings_probs) ≈ 1
    @assert length(markings_probs) == length(markings)
    @assert length(markings_probs) == length(jitters)
    marking_select=Categorical(markings_probs)
    n = maximum(maximum.(markings))
    return new{R,I,J}(n,rate_parent,markings,marking_select,jitters)
  end
end

jitter!(v::Vector{R},k::Integer,g::GTAS) where R = jitter!(v,g.jitters[k])

function make_samples_with_parent(g::GTAS{R,I},t_tot::R) where {R,I}
  ts_ancestor = make_poisson_samples(g.rate_parent,t_tot)
  nt = length(ts_ancestor)
  attributions = rand(g.marking_selector,nt)
  trains = [Vector{R}(undef,0) for _ in 1:g.n]
  # can be optimized by preallocating trains, based on expected rate!
  for (t_k,k) in zip(ts_ancestor,attributions)
    mark_k = g.markings[k]
    n_k = length(mark_k)
    spikes_k = jitter!(fill(t_k,n_k),k,g)
    for (tkk,kk) in zip(spikes_k,mark_k)
      push!(trains[kk],tkk)
    end
  end
  # clean up trains... sorted, and > 0
  sort!.(trains)
  tmin = minimum(first.(trains))
  if tmin < 0
    for tr in trains
      tr .-= tmin
    end
    ts_ancestor .-= tmin
  end
  return trains,ts_ancestor,attributions
end

function make_samples(gtas::GTAS,t_tot::Real)
  return make_samples_with_parent(gtas,t_tot)[1]
end

function make_poisson_samples(rate::R,t_tot::R) where R
  ret = Vector{R}(undef,round(Integer,1.3*rate*t_tot+10)) # preallocate
  t_curr = zero(R)
  k_curr = 1
  while t_curr <= t_tot
    Δt = -log(rand())/rate
    t_curr += Δt
    ret[k_curr] = t_curr
    k_curr += 1
  end
  return keepat!(ret,1:k_curr-2)
end

####
# Cumulants

function get_probability_weight(neu::I,markings::Vector{Vector{I}},
      probs::Vector{<:Real}) where I
  probs_all = map(zip(markings,probs)) do (marks,p)
    fill(p,length(marks))
  end
  probs_all = vcat(probs_all...)
  marks_all = vcat(markings...)
  idxs = findall(==(neu),marks_all)
  return sum(probs_all[idxs])
end

function get_parent_rate(neu::I,rate_target::Real,markings::Vector{Vector{I}},
     probs::Vector{<:Real}) where I<:Integer
  ptot = get_probability_weight(neu,markings,probs)
  return rate_target / ptot
end
function get_parent_rate_target_mean(rate_target::Real,markings::Vector{Vector{I}},
  probs::Vector{<:Real}) where I<:Integer
n = maximum(maximum.(markings))
p_avg = mapreduce(neu->get_probability_weight(neu,markings,probs),+,1:n) / n
return rate_target / p_avg
end

function get_expected_rates(g::GTAS{R,I}) where {R,I}
  mks = g.markings
  probs = g.marking_selector.p
  rat_par = g.rate_parent
  ret = rat_par .* [get_probability_weight(i,mks,probs) for i in 1:g.n]
  return ret
end

# UNDERESTIMATE ! Must account for oder zero ! Chance of events being in the same
# time interval considered. This one converges to numeric only for very small time bins!
function get_expected_Pearson(i::Integer,j::Integer,g::GTAS{R,I}) where {R,I}
  marks = g.markings
  probs = g.marking_selector.p
  uiall = 0.0
  ujall = 0.0
  uijall = 0.0
  for (mark,prob) in zip(marks,probs)
    hasi = i in mark
    hasj = j in mark
    if hasi && hasj
      uijall += prob
    elseif hasi
      uiall += prob
    elseif hasj
      ujall += prob
    end
  end
  return uijall / (uijall+uiall+ujall)
end

function get_expected_Pearson(g::GTAS{R,I}) where {R,I}
  ret = Matrix{R}(undef,g.n,g.n)
  for i in 1:g.n, j in i+1:g.n
    ret[i,j] = get_expected_Pearson(i,j,g)
    ret[j,i] = ret[i,j]
  end
  ret[diagind(ret)] .= one(R)
  return ret
end

#########
# numerical analysis

function bin_spikes(Y::Vector{R},dt::R,Tend::R;
    Tstart::R=0.0) where R
  times = range(Tstart,Tend;step=dt)
  ret = fill(0,length(times)-1)
  for y in Y
    if Tstart < y <= last(times)
      k = searchsortedfirst(times,y)-1
      ret[k] += 1
    end
  end
  return ret
end

function bin_and_rates(Y::Vector{R},dt::R,Tend::R;Tstart::R=0.0) where R
  times = range(Tstart,Tend;step=dt)
  return midpoints(times),bin_spikes(Y,dt,Tend;Tstart=Tstart) ./ dt
end

function bin_and_rates(Ys::Vector{Vector{R}},dt::R,Tend::R;Tstart::R=0.0) where R
  times = range(Tstart,Tend;step=dt)
  ret = [bin_spikes(Y,dt,Tend;Tstart=Tstart) ./ dt for Y in Ys]
  return midpoints(times),vcat(transpose.(ret)...)
end

@inline function get_times_strict(dt::R,Tend::R;Tstart::R=0.0) where R<:Real
  return range(Tstart,Tend-dt;step=dt)
end

# the first element (zero lag) is always rate/dτ
function covariance_self_numerical(Y::Vector{R},dτ::R,τmax::R,
     Tend::Union{R,Nothing}=nothing) where R
  τtimes,ret = covariance_density_numerical([Y,],dτ,τmax;verbose=false,Tend=Tend)
  return  τtimes, ret[:,1,1]
end

function covariance_density_numerical(Ys::Vector{Vector{R}},dτ::Real,τmax::R;
   Tend::Union{R,Nothing}=nothing,verbose::Bool=false) where R
  Tend = something(Tend, maximum(last,Ys)- dτ)
  ndt = round(Integer,τmax/dτ)
  n = length(Ys)
  ret = Array{Float64}(undef,ndt,n,n)
  if verbose
      @info "The full dynamical iteration has $(round(Integer,Tend/dτ)) bins ! (too many?)"
  end
  for i in 1:n
    binnedi = bin_spikes(Ys[i],dτ,Tend)
    fmi = length(Ys[i]) / Tend # mean frequency
    ndt_tot = length(binnedi)
    _ret_alloc = Vector{R}(undef,ndt)
    for j in 1:n
      if verbose
        @info "now computing cov for pair $i,$j"
      end
      binnedj =  i==j ? binnedi : bin_spikes(Ys[j],dτ,Tend)
      fmj = length(Ys[j]) / Tend # mean frequency
      binnedj_sh = similar(binnedj)
      @inbounds @simd for k in 0:ndt-1
        circshift!(binnedj_sh,binnedj,k)
        _ret_alloc[k+1] = dot(binnedi,binnedj_sh)
      end
      @. _ret_alloc = (_ret_alloc / (ndt_tot*dτ^2)) - fmi*fmj
      ret[:,i,j] = _ret_alloc
    end
  end
  return get_times_strict(dτ,τmax), ret
end

function covariance_density_ij(Ys::Vector{Vector{R}},i::Integer,j::Integer,dτ::R,τmax::R;
   Tend::Union{R,Nothing}=nothing) where R
  return covariance_density_ij(Ys[i],Ys[j],dτ,τmax;Tend=Tend)
end

function covariance_density_ij(X::Vector{R},Y::Vector{R},dτ::Real,τmax::R;
    Tend::Union{R,Nothing}=nothing) where R
  times = get_times_strict(dτ,τmax)
  ndt = length(times)
  times_ret = vcat(-reverse(times[2:end]),times)
  Tend = something(Tend, max(X[end],Y[end])- dτ)
  ret = Vector{Float64}(undef,2*ndt-1)
  binnedx = bin_spikes(X,dτ,Tend)
  binnedy = bin_spikes(Y,dτ,Tend)
  fx = length(X) / Tend # mean frequency
  fy = length(Y) / Tend # mean frequency
  ndt_tot = length(binnedx)
  binned_sh = similar(binnedx)
  # 0 and forward
  @simd for k in 0:ndt-1
    circshift!(binned_sh,binnedy,-k)
    ret[ndt-1+k+1] = dot(binnedx,binned_sh)
  end
  # backward
  @simd for k in 1:ndt-1
    circshift!(binned_sh,binnedy,k)
    ret[ndt-k] = dot(binnedx,binned_sh)
  end
  @. ret = (ret / (ndt_tot*dτ^2)) - fx*fy
  return times_ret, ret
end

## include the part with anti-spikes and negative correlations
include("A-GTAS.jl")

end # of Module


text:
using Test, Random, CategoricalArrays, Dates, PooledArrays
const ≅ = isequal
@testset "general usage" begin

    ds = Dataset(x1 = [1,2,3,4], x2 = [1,4,9,16])
    dst =Dataset(_variables_ =["x1","x2"], _c1 = [1,1], _c2 = [2,4],
                    _c3 = [3,9], _c4 = [4,16])
    @test transpose(ds, [:x1, :x2]) == dst

    ds = Dataset(rand(10,5), :auto)
    ds2 = transpose(ds, r"x")

    @test Matrix(ds) == permutedims(Matrix(ds2[!,r"_c"]))

    ds = Dataset(rand(10,5), :auto)
    # allowmissing!(ds)
    ds2 = transpose(ds, r"x")

    @test transpose(ds, r"x") == transpose(ds, :)
    @test transpose(ds, r"x") == transpose(ds, [1,2,3,4,5])
    @test Matrix(ds) == permutedims(Matrix(ds2[!,2:end]))

    ds = Dataset(rand(10,5), :auto)
    # allowmissing!(ds)
    ds[1, 4] = missing
    ds[4, 5] = missing
    ds2 = transpose(ds, r"x")
    @test Matrix(ds) ≅ permutedims(Matrix(ds2[!,Not(:_variables_)]))

    ds = InMemoryDatasets.hcat!(Dataset(rand(Bool, 15,2),[:b1,:b2]), Dataset(rand(15,5),:auto))
    ds2 = transpose(ds, :)
    @test Matrix(ds) == permutedims(Matrix(ds2[!,Not(:_variables_)]))

    ds = Dataset(foo = ["one", "one", "one", "two", "two","two"],
                    bar = ['A', 'B', 'C', 'A', 'B', 'C'],
                    baz = [1, 2, 3, 4, 5, 6],
                    zoo = ['x', 'y', 'z', 'q', 'w', 't'])
    ds2 = transpose(groupby(ds, :foo, stable = true), :baz, id = :bar)
    @test sort(ds2, :) == sort(transpose(gatherby(ds, :foo), :baz, id = :bar), :)
    ds3 = Dataset(foo = ["one", "two"], _variables_ = ["baz", "baz"], A = [1, 4], B = [2, 5], C = [3, 6])

    @test transpose(groupby(ds, :foo, stable = true), :baz, id = :bar) == transpose(groupby(ds, :foo, stable = true), [:baz], id = :bar)
    @test transpose(groupby(ds, :foo, stable = true), :baz, id = :bar) == transpose(groupby(ds, :foo, stable = true), [:baz], id = :bar)
    @test transpose(groupby(ds, :foo, stable = true), :baz, id = :bar) == transpose(groupby(ds, :foo, stable = true), :baz, id = :bar)
    @test ds2 == ds3

    ds = Dataset(id = [1, 2, 3, 1], x1 = rand(4), x2 = rand(4))
    @test_throws AssertionError transpose(ds, r"x", id = :id)

    ds = Dataset(rand(1:100, 1000, 10), :auto)
    insertcols!(ds, 1, :g => repeat(1:100, inner = 10))
    insertcols!(ds, 2, :id => repeat(1:10, 100))
    # duplicate and id within the last group
    ds[1000, :id] = 1
    @test_throws AssertionError transpose(groupby(ds, :g, stable = true), r"x", id = :id, threads = false)

    ds = Dataset(rand(1000, 100), :auto)
    mds = Matrix(ds)

    @test sum.(eachcol(transpose(ds, r"x")[!, r"_c"])) ≈ sum(mds, dims = 2)

    ds = Dataset([[1, 2], [1.1, 2.0],[1.1, 2.1],[1.1, 2.0]]
                    ,[:person, Symbol("11/2020"), Symbol("12/2020"), Symbol("1/2021")])

    groupby!(ds, :person, stable = true)
    dst = transpose(ds, Not(:person),
                        renamerowid = x -> Date(x, dateformat"m/y"),
                        variable_name = "Date",
                         renamecolid = x -> "measurement")
   @test sort(dst, :) == sort(transpose(gatherby(ds, :person), Not(:person),
                     renamerowid = x -> Date(x, dateformat"m/y"),
                     variable_name = "Date",
                      renamecolid = x -> "measurement"),:)
    dstm = Dataset(person = [1,1,1,2,2,2],
                Date = Date.(repeat(["2020-11-01","2020-12-01","2021-01-01"], 2)),
                measurement = [1.1, 1.1, 1.1, 2.0, 2.1, 2.0])
    @test dst == dstm
    dst = transpose(groupby(ds, :person, stable = true), Not(:person),
                        renamerowid = x -> Date(x, dateformat"m/y"),
                        variable_name = "Date",
                         renamecolid = x -> "measurement")
    dstm = Dataset(person = [1,1,1,2,2,2],
                Date = Date.(repeat(["2020-11-01","2020-12-01","2021-01-01"], 2)),
                measurement = [1.1, 1.1, 1.1, 2.0, 2.1, 2.0])
    @test dst == dstm

    ds = Dataset(rand(100,5),:auto)
    insertcols!(ds, 1, :id => repeat(1:20, 5))
    insertcols!(ds, 1, :g => repeat(1:5, inner = 20))

    dst = transpose(groupby(ds, [:g, :id], stable = true), r"x" , variable_name = "variable", renamecolid = x -> "value")

    @test sort(dst, :) == sort(transpose(gatherby(ds, [:g, :id]), r"x" , variable_name = "variable", renamecolid = x -> "value"), :)

    ds = Dataset(group = repeat(1:3, inner = 2),
                                 b = repeat(1:2, inner = 3),
                                 c = repeat(1:1, inner = 6),
                                 d = repeat(1:6, inner = 1),
                                 e = string.('a':'f'))
    # allowmissing!(ds)
    ds[2, :b] = missing
    ds[4, :d] = missing
    ds2 = transpose(groupby(ds, :group, stable = true), 2:4, id = :e, default = 0)

    ds3 = Dataset(group = repeat(1:3, inner = 3),
                    _variables_ = string.(repeat('b':'d', 3)),
                    a = [1,1,1,0,0,0,0,0,0],
                    b = [missing,1,2,0,0,0,0,0,0],
                    c = [0,0,0, 1,1,3, 0,0,0],
                    d = [0,0,0, 2,1,missing, 0,0,0],
                    e = [0,0,0,0,0,0,2,1,5],
                    f = [0,0,0,0,0,0,2,1,6] )
    @test ds2 ≅ ds3

    ds = Dataset(id = ["r3", "r1", "r2" , "r4"], x1 = [1,2,3,4], x2 = [1,4,9,16])
    ds2 = transpose(ds, [:x1,:x2], id = :id)
    ds3 = Dataset([["x1", "x2"],
                        [1, 1],
                        [2, 4],
                        [3, 9],
                        [4, 16]],
                        [:_variables_, :r3, :r1, :r2, :r4])
    @test ds3 == ds2

    pop = Dataset(country = ["c1","c1","c2","c2","c3","c3"],
                            sex = ["male", "female", "male", "female", "male", "female"],
                            pop_2000 = [100, 120, 150, 155, 170, 190],
                            pop_2010 = [110, 120, 155, 160, 178, 200],
                            pop_2020 = [115, 130, 161, 165, 180, 203])

    popt = transpose(groupby(pop, :country, stable = true), r"pop_",
                            id = :sex, variable_name = "year",
                            renamerowid = x -> match(r"[0-9]+",x).match, renamecolid = x -> x * "_pop")
    poptm = Dataset([["c1", "c1", "c1", "c2", "c2", "c2", "c3", "c3", "c3"],
            SubString{String}["2000", "2010", "2020", "2000", "2010", "2020", "2000", "2010", "2020"],
            Union{Missing, Int64}[100, 110, 115, 150, 155, 161, 170, 178, 180],
            Union{Missing, Int64}[120, 120, 130, 155, 160, 165, 190, 200, 203]],
            [:country,:year,:male_pop,:female_pop])

    @test popt == poptm
    popt = transpose(groupby(pop, r"cou", stable = true), r"pop", id = r"sex",  variable_name = "year",
                            renamerowid = x -> match(r"[0-9]+",x).match, renamecolid = x -> x * "_pop")
    @test popt == poptm
    pop.country = PooledArray(pop.country)
    popt = transpose(groupby(pop, r"cou", stable = true), r"pop", id = r"sex",  variable_name = "year",
                            renamerowid = x -> match(r"[0-9]+",x).match, renamecolid = x -> x * "_pop")
    @test popt.country == PooledArray(["c1", "c1", "c1", "c2", "c2", "c2", "c3", "c3", "c3"])
    ds =  Dataset(region = repeat(["North","North","South","South"],2),
                 fuel_type = repeat(["gas","coal"],4),
                 load = [.1,.2,.5,.1,6.,4.3,.1,6.],
                 time = [1,1,1,1,2,2,2,2],
                 )

    ds2 = transpose(groupby(ds, :time, stable = true), :load, id = 1:2)
    ds3 = Dataset([ Union{Missing, Int64}[1, 2],
                 ["load", "load"],
                 Union{Missing, Float64}[0.1, 6.0],
                 Union{Missing, Float64}[0.2, 4.3],
                 Union{Missing, Float64}[0.5, 0.1],
                 Union{Missing, Float64}[0.1, 6.0]], ["time", "_variables_", "(\"North\", \"gas\")", "(\"North\", \"coal\")", "(\"South\", \"gas\")", "(\"South\", \"coal\")"])

     @test ds2 == ds3
     ds = Dataset(A_2018=1:4, A_2019=5:8, B_2017=9:12,
                             B_2018=9:12, B_2019 = [missing,13,14,15],
                              ID = [1,2,3,4])
      f(x) =  match(r"[0-9]+",x).match
      dsA = transpose(groupby(ds, :ID, stable = true), r"A", renamerowid = f, variable_name = "Year", renamecolid = x->"A");
      dsB = transpose(groupby(ds, :ID, stable = true), r"B", renamerowid = f, variable_name = "Year", renamecolid = x->"B");
      ds2 = outerjoin(dsA, dsB, on = [:ID, :Year])
      ds3 = Dataset([[1, 1, 2, 2, 3, 3, 4, 4, 1, 2, 3, 4],
                 SubString{String}["2018", "2019", "2018", "2019", "2018", "2019", "2018", "2019", "2017", "2017", "2017", "2017"],
                 Union{Missing, Int64}[1, 5, 2, 6, 3, 7, 4, 8, missing, missing, missing, missing],
                 Union{Missing, Int64}[9, missing, 10, 13, 11, 14, 12, 15, 9, 10, 11, 12]], [:ID,:Year,:A,:B])
        @test ds2 ≅ ds3
        ds = Dataset(paddockId= [0, 0, 1, 1, 2, 2],
                                color= ["red", "blue", "red", "blue", "red", "blue"],
                                count= [3, 4, 3, 4, 3, 4],
                                weight= [0.2, 0.3, 0.2, 0.3, 0.2, 0.2])
        ds2 = transpose(groupby(transpose(groupby(ds, [:paddockId,:color], stable = true), [:count,:weight]), :paddockId, stable = true),
                             :_c1, id = 2:3)
        ds3 = Dataset([Union{Missing, Int64}[0, 1, 2],
             Union{Missing, String}["_c1", "_c1", "_c1"],
             Union{Missing, Float64}[4.0, 4.0, 4.0],
             Union{Missing, Float64}[0.3, 0.3, 0.2],
             Union{Missing, Float64}[3.0, 3.0, 3.0],
             Union{Missing, Float64}[0.2, 0.2, 0.2]], ["paddockId", "_variables_", "(\"blue\", \"count\")", "(\"blue\", \"weight\")", "(\"red\", \"count\")", "(\"red\", \"weight\")"])

         @test ds2 == ds3

        ds = Dataset(x1 = [9,2,8,6,8], x2 = [8,1,6,2,3], x3 = [6,5,3,10,8])
        ds2 = transpose(ds, r"x", renamerowid = x -> match(r"[0-9]+",x).match,renamecolid = x -> "_column_" * string(x))
        ds3 = Dataset([SubString{String}["1", "2", "3"],
                         [9, 8, 6],
                         [2, 1, 5],
                         [8, 6, 3],
                         [6, 2, 10],
                         [8, 3, 8]],[:_variables_,:_column_1,:_column_2,:_column_3,:_column_4,:_column_5])
         @test ds2 == ds3
         ds = Dataset(a=["x", "y"], b=[1, "two"], c=[3, 4], d=[true, false])
         ds2 = transpose(ds, [:b, :c, :d], id = :a, variable_name = "new_col")
         ds3 = Dataset([["b", "c", "d"],
                         Any[1, 3, true],
                         Any["two", 4, false]],[:new_col,:x,:y])
         @test ds2 == ds3

         ds = Dataset(g1 = ["g2", "g1", "g1", "g2", "g1", "g3"], c1 = ["c1", "c1", "c1", "c1", "c1", "c1"],
                id1 = ["id2_g2", "id1_g1", "id2_g1", "id1_g2", "id3_g1", "id4_g3"], id2 = [2, 3, 1, 5, 4, 6],
                val1 = [1, 2, 3, 4, 5, 6], val2 = [1.1, 2.1, 0.5, 1.3, 1.0, 2])
        ds2 = transpose(ds, [:val1, :val2])
        ds3 = transpose(ds, [:val1, :val2], id = :id2)
        ds4 = transpose(ds, [:val1, :val2], id = [:id1, :id2])
        ds5 = transpose(groupby(ds, 1, stable = true), [:val1, :val2])
        ds6 = transpose(groupby(ds, 1, stable = true), [:val1, :val2], id = :id1)
        ds7 = transpose(groupby(ds, 1, stable = true), [:val1, :val2], id = [:id1, :id2])
        ds8 = transpose(groupby(ds, 1:2, stable = true), [:val1, :val2], id = :id2)
        ds9 = transpose(groupby(ds, 1:2, stable = true), [:val1], id = [:id2, :id1])

        tds2 = Dataset([Union{Missing, String}["val1", "val2"],
                 Union{Missing, Float64}[1.0, 1.1],
                 Union{Missing, Float64}[2.0, 2.1],
                 Union{Missing, Float64}[3.0, 0.5],
                 Union{Missing, Float64}[4.0, 1.3],
                 Union{Missing, Float64}[5.0, 1.0],
                 Union{Missing, Float64}[6.0, 2.0]], ["_variables_", "_c1", "_c2", "_c3", "_c4", "_c5", "_c6"])
        tds3 = Dataset([Union{Missing, String}["val1", "val2"],
                 Union{Missing, Float64}[1.0, 1.1],
                 Union{Missing, Float64}[2.0, 2.1],
                 Union{Missing, Float64}[3.0, 0.5],
                 Union{Missing, Float64}[4.0, 1.3],
                 Union{Missing, Float64}[5.0, 1.0],
                 Union{Missing, Float64}[6.0, 2.0]], ["_variables_", "2", "3", "1", "5", "4", "6"])
        tds4 = Dataset([Union{Missing, String}["val1", "val2"],
                 Union{Missing, Float64}[1.0, 1.1],
                 Union{Missing, Float64}[2.0, 2.1],
                 Union{Missing, Float64}[3.0, 0.5],
                 Union{Missing, Float64}[4.0, 1.3],
                 Union{Missing, Float64}[5.0, 1.0],
                 Union{Missing, Float64}[6.0, 2.0]], ["_variables_", "(\"id2_g2\", 2)", "(\"id1_g1\", 3)", "(\"id2_g1\", 1)", "(\"id1_g2\", 5)", "(\"id3_g1\", 4)", "(\"id4_g3\", 6)"])
        tds5 = Dataset([Union{Missing, String}["g1", "g1", "g2", "g2", "g3", "g3"],
                 Union{Missing, String}["val1", "val2", "val1", "val2", "val1", "val2"],
                 Union{Missing, Float64}[2.0, 2.1, 1.0, 1.1, 6.0, 2.0],
                 Union{Missing, Float64}[3.0, 0.5, 4.0, 1.3, missing, missing],
                 Union{Missing, Float64}[5.0, 1.0, missing, missing, missing, missing]], ["g1", "_variables_", "_c1", "_c2", "_c3"])
        tds6 = Dataset([Union{Missing, String}["g1", "g1", "g2", "g2", "g3", "g3"],
                 Union{Missing, String}["val1", "val2", "val1", "val2", "val1", "val2"],
                 Union{Missing, Float64}[2.0, 2.1, missing, missing, missing, missing],
                 Union{Missing, Float64}[3.0, 0.5, missing, missing, missing, missing],
                 Union{Missing, Float64}[5.0, 1.0, missing, missing, missing, missing],
                 Union{Missing, Float64}[missing, missing, 1.0, 1.1, missing, missing],
                 Union{Missing, Float64}[missing, missing, 4.0, 1.3, missing, missing],
                 Union{Missing, Float64}[missing, missing, missing, missing, 6.0, 2.0]], ["g1", "_variables_", "id1_g1", "id2_g1", "id3_g1", "id2_g2", "id1_g2", "id4_g3"])
        tds7 = Dataset([Union{Missing, String}["g1", "g1", "g2", "g2", "g3", "g3"],
                 Union{Missing, String}["val1", "val2", "val1", "val2", "val1", "val2"],
                 Union{Missing, Float64}[2.0, 2.1, missing, missing, missing, missing],
                 Union{Missing, Float64}[3.0, 0.5, missing, missing, missing, missing],
                 Union{Missing, Float64}[5.0, 1.0, missing, missing, missing, missing],
                 Union{Missing, Float64}[missing, missing, 1.0, 1.1, missing, missing],
                 Union{Missing, Float64}[missing, missing, 4.0, 1.3, missing, missing],
                 Union{Missing, Float64}[missing, missing, missing, missing, 6.0, 2.0]], ["g1", "_variables_", "(\"id1_g1\", 3)", "(\"id2_g1\", 1)", "(\"id3_g1\", 4)", "(\"id2_g2\", 2)", "(\"id1_g2\", 5)", "(\"id4_g3\", 6)"])
        tds8 = Dataset([Union{Missing, String}["g1", "g1", "g2", "g2", "g3", "g3"],
                 Union{Missing, String}["c1", "c1", "c1", "c1", "c1", "c1"],
                 Union{Missing, String}["val1", "val2", "val1", "val2", "val1", "val2"],
                 Union{Missing, Float64}[2.0, 2.1, missing, missing, missing, missing],
                 Union{Missing, Float64}[3.0, 0.5, missing, missing, missing, missing],
                 Union{Missing, Float64}[5.0, 1.0, missing, missing, missing, missing],
                 Union{Missing, Float64}[missing, missing, 1.0, 1.1, missing, missing],
                 Union{Missing, Float64}[missing, missing, 4.0, 1.3, missing, missing],
                 Union{Missing, Float64}[missing, missing, missing, missing, 6.0, 2.0]], ["g1", "c1", "_variables_", "3", "1", "4", "2", "5", "6"])
        tds9 = Dataset([Union{Missing, String}["g1", "g2", "g3"],
                 Union{Missing, String}["c1", "c1", "c1"],
                 Union{Missing, String}["val1", "val1", "val1"],
                 Union{Missing, Int64}[2, missing, missing],
                 Union{Missing, Int64}[3, missing, missing],
                 Union{Missing, Int64}[5, missing, missing],
                 Union{Missing, Int64}[missing, 1, missing],
                 Union{Missing, Int64}[missing, 4, missing],
                 Union{Missing, Int64}[missing, missing, 6]], ["g1", "c1", "_variables_", "(3, \"id1_g1\")", "(1, \"id2_g1\")", "(4, \"id3_g1\")", "(2, \"id2_g2\")", "(5, \"id1_g2\")", "(6, \"id4_g3\")"])
        @test ds2 == tds2
        @test ds3 == tds3
        @test ds4 == tds4
        @test ds5 == tds5
        @test ds6 == tds6
        @test ds7 == tds7
        @test ds8 == tds8
        @test ds9 == tds9
        @test_throws AssertionError transpose(ds, :val1, id = :g1)
        ds = Dataset(g = [1,1,2,2], x1=1:4, y=[454,65,65,65])
        setformat!(ds, 2=>isodd)
        t1 = transpose(groupby(ds, 1), :y, id = :x1)
        t1_t = Dataset("g" => [1,2], "_variables_" => ["y","y"], "true" => [454, 65], "false"=>[65,65])
        @test t1 == t1_t
        t2 = transpose(groupby(ds, 1), :y, id = :x1, mapformats = false)
        t2_t = Dataset("g" => [1,2], "_variables_" => ["y","y"], "1"=>[454, missing], "2"=>[65, missing], "3"=>[missing, 65], "4"=>[missing, 65])
        @test t2 == t2_t
        @test_throws AssertionError transpose(ds, :y, id = :x1, mapformats = true)
        t3 = transpose(ds, :y, id = :x1, mapformats = false)
        t3_t = Dataset("_variables_" => ["y"], "1"=>[454], "2"=>[65], "3"=>[65], "4"=>[65])
        @test t3 == t3_t
        setformat!(ds, 1=>isodd)
        t4 = transpose(ds, :y, id = 1:2)
        t4_t = Dataset([Union{Missing, String}["y"],
                 Union{Missing, Int64}[454],
                 Union{Missing, Int64}[65],
                 Union{Missing, Int64}[65],
                 Union{Missing, Int64}[65]], ["_variables_", "(true, true)", "(true, false)", "(false, true)", "(false, false)"])
        @test t4 == t4_t
        t5 = transpose(ds, :y, id = 1:2, mapformats = false)
        t5_t = Dataset([Union{Missing, String}["y"],
                 Union{Missing, Int64}[454],
                 Union{Missing, Int64}[65],
                 Union{Missing, Int64}[65],
                 Union{Missing, Int64}[65]], ["_variables_", "(1, 1)", "(1, 2)", "(2, 3)", "(2, 4)"])
        @test t5 == t5_t
        ds = Dataset(g1 = [1, 2, 3, 4], x = [1,2,3,4], g2 = [false, false, true, true])
        setformat!(ds, 1=>isodd)
        groupby!(ds, [:g2, :g1])
        dst = transpose(ds, :x, variable_name = nothing)
        dst_t = compare(ds, dst, cols =[:g1=>:g1, :g2=>:g2, :x=>:_c1], mapformats = true)
        @test all(byrow(dst_t, all, :))

        ds = Dataset(x = [1,2,1,2,3], y = [missing, missing, missing, missing, missing])
        @test transpose(groupby(ds,1), :y) == Dataset(x = [1,2,3], _variables_=["y","y","y"], _c1=[missing, missing, missing], _c2=[missing, missing, missing])
        @test transpose(groupby(ds,1), :y, default = 0) == Dataset(x = [1,2,3], _variables_=["y","y","y"], _c1=[missing, missing, missing], _c2=[missing, missing, 0])
        @test transpose(gatherby(ds,1), :y, default = 0) == Dataset(x = [1,2,3], _variables_=["y","y","y"], _c1=[missing, missing, missing], _c2=[missing, missing, 0])

end

@testset "Outputs - Checking types" begin
    ds = Dataset(Fish = CategoricalArray{Union{String, Missing}}(["Bob", "Bob", "Batman", "Batman"]),
                   Key = CategoricalArray{Union{String, Missing}}(["Mass", "Color", "Mass", "Color"]),
                   Value = Union{String, Missing}["12 g", "Red", "18 g", "Grey"])
    levels!(ds[!, 1].val, ["XXX", "Bob", "Batman"])
    levels!(ds[!, 2].val, ["YYY", "Color", "Mass"])
#     Not sure if it is relevant, however, we are doing it here
    ds2 = transpose(groupby!(ds, :Fish, stable = true), [:Value], id = :Key)
    @test levels(ds[!, 1].val) == ["XXX", "Bob", "Batman"] # make sure we did not mess ds[!, 1] levels
    @test levels(ds[!, 2].val) == ["YYY", "Color", "Mass"] # make sure we did not mess ds[!, 2] levels
#     duplicates
    @test_throws AssertionError transpose(ungroup!(ds), [:Value], id = :Key)

    ds = Dataset(Fish = CategoricalArray{Union{String, Missing}}(["Bob", "Bob", "Batman", "Batman"]),
                   Key = CategoricalArray{Union{String, Missing}}(["Mass", "Color", "Mass", "Color"]),
                   Value = Union{String, Missing}["12 g", "Red", "18 g", "Grey"])
    levels!(ds[!, 1].val, ["XXX", "Bob", "Batman"])
    levels!(ds[!, 2].val, ["YYY", "Color", "Mass"])
    ds2 = transpose(groupby(ds, :Fish, stable = true), [:Value], id = :Key, renamecolid=x->string("_", uppercase(string(x)), "_"))
    ds4 = Dataset(Fish = Union{String, Missing}["Bob", "Batman"],
                _variables_ = String["Value", "Value"],
                _MASS_ = Union{String, Missing}["12 g", "18 g"],
                _COLOR_ = Union{String, Missing}["Red", "Grey"])
    @test ds2 == ds4
    # without categorical array
    ds = Dataset(Fish = ["Bob", "Bob", "Batman", "Batman"],
                   Key = ["Mass", "Color", "Mass", "Color"],
                   Value = ["12 g", "Red", "18 g", "Grey"])
    ds2 = transpose(groupby(ds, :Fish, stable = true), [:Value], id = :Key)
    ds4 = Dataset(Fish = ["Batman", "Bob"],
                    _variables_ = ["Value", "Value"],
                    Mass = ["18 g", "12 g"],
                    Color = ["Grey", "Red"]
                    )
    @test ds2 ≅ ds4
    @test eltype(ds2[!, :Fish]) <: Union{String, Missing}
    #Make sure transpose works with missing values at the start of the value column
    # allowmissing!(ds, :Value)
    ds[1, :Value] = missing
    ds2 = transpose(groupby(ds, :Fish, stable = true), [:Value], id =  :Key)
    #This changes the expected result
    # allowmissing!(ds4, :Mass)
    ds4[2, :Mass] = missing
    @test ds2 ≅ ds4

    # test missing value in grouping variable
    mds = Dataset(RowID = 4:-1:1, id=[missing, 1, 2, 3], a=1:4, b=1:4)
    @test select(transpose(groupby(transpose(groupby(mds, [:RowID, :id], stable = true), [:a,:b]), [:RowID, :id]), [:_c1], id = :_variables_), :RowID, :id, :a, :b) ≅ sort!(mds,1)
    @test select(transpose(groupby(transpose(groupby(mds, [:RowID, :id], stable = true), Not(1,2)),[:RowID, :id]), [:_c1], id = :_variables_), :RowID, :id, :a, :b) ≅ sort!(mds,1)

    # test more than one grouping column
    wide = Dataset(id = 1:12,
                     a  = repeat([1:3;], inner = [4]),
                     b  = repeat([1:4;], inner = [3]),
                     c  = randn(12),
                     d  = randn(12))
    w2 = wide[:, [1, 2, 4, 5]]
    InMemoryDatasets.rename!(w2, [:id, :a, :_C_, :_D_])
    long = transpose(groupby(wide, [:id, :a, :b], stable = true), [:c, :d])
    wide3 = transpose(groupby(long, [:id, :a, :b], stable = true), [:_c1], id = :_variables_)
    @test select(wide3, Not(:_variables_)) == wide
    ds = Dataset([repeat(1:2, inner=4), repeat('a':'d', outer=2), collect(1:8)],
                       [:id, :variable, :value])

    uds = transpose(groupby(ds, :id, stable = true), [:value], id = :variable)
    @test select(uds, Not(:_variables_)) == Dataset([Union{Int, Missing}[1, 2], Union{Int, Missing}[1, 5],
                                Union{Int, Missing}[2, 6], Union{Int, Missing}[3, 7],
                                Union{Int, Missing}[4, 8]], [:id, :a, :b, :c, :d])

    @test isa(uds[!, 1].val, Vector{Union{Missing,Int}})
    @test isa(uds[!,:a].val, Vector{Union{Int,Missing}})
    @test isa(uds[!,:b].val, Vector{Union{Int,Missing}})
    @test isa(uds[!,:c].val, Vector{Union{Int,Missing}})
    @test isa(uds[!,:d].val, Vector{Union{Int,Missing}})
    ds = Dataset(Any[categorical(repeat(2:-1:1, inner=4)),
                           categorical(repeat('a':'d', outer=2)), categorical(1:8)],
                       [:id, :variable, :value])

    uds = transpose(groupby!(ds, :id, stable = true), [:value], id = :variable)
    @test isa(uds[!, 1].val, CategoricalArray{Union{Missing, Int64},1,UInt32})
    @test isa(uds[!, :a].val, CategoricalArray{Union{Missing, Int64},1, UInt32})
    @test isa(uds[!, :b].val, CategoricalArray{Union{Missing, Int64},1, UInt32})
    @test isa(uds[!, :c].val, CategoricalArray{Union{Missing, Int64},1, UInt32})
    @test isa(uds[!, :d].val, CategoricalArray{Union{Missing, Int64},1, UInt32})

    ds1 = Dataset(a=["x", "y"], b=rand(2), c=[1, 2], d=rand(Bool, 2))

    @test_throws MethodError transpose(ds1)
    @test_throws ArgumentError transpose(ds1, :bar)

    ds1_pd = transpose(ds1, 2:4, id = 1)
    @test size(ds1_pd, 1) == ncol(ds1) - 1
    @test size(ds1_pd, 2) == nrow(ds1) + 1
    @test names(ds1_pd) == ["_variables_", "x", "y"]
    @test names(transpose(ds1, 2:4, id = 1, variable_name = "foo")) == ["foo", "x", "y"]

end

@testset "passing Tuple as cols" begin
        ds = Dataset([Union{Missing, Int64}[1, 1, 1, 2, 2, 2],
                 Union{Missing, String}["foo", "bar", "monty", "foo", "bar", "monty"],
                 Union{Missing, String}["a", "b", "c", "d", "e", "f"],
                 Union{Missing, Int64}[1, 2, 3, 4, 5, 6]], [:g, :key, :foo, :bar])
        dst = transpose(groupby(ds, :g), (:foo, :bar), id = :key, variable_name = "_variables_")
        dst_t = Dataset([Union{Missing, Int64}[1, 2],
                 Union{Missing, String}["foo", "foo"],
                 Union{Missing, String}["bar", "bar"],
                 Union{Missing, String}["a", "d"],
                 Union{Missing, String}["b", "e"],
                 Union{Missing, String}["c", "f"],
                 Union{Missing, Int64}[1, 4],
                 Union{Missing, Int64}[2, 5],
                 Union{Missing, Int64}[3, 6]], ["g", "_variables_", "_variables__1", "foo", "bar", "monty", "foo_1", "bar_1", "monty_1"])
        @test dst == dst_t
        ds = Dataset(paddockId= [0, 0, 1, 1, 2, 2],
                                color= ["red", "blue", "red", "blue", "red", "blue"],
                                count= [3, 4, 3, 4, 3, 4],
                                weight= [0.2, 0.3, 0.2, 0.3, 0.2, 0.2])
        dst = transpose(groupby(ds, 1), (:count, :weight), id = :color, variable_name = :_variables_)
        dst_t = Dataset([Union{Missing, Int64}[0, 1, 2],
                 Union{Missing, String}["count", "count", "count"],
                 Union{Missing, String}["weight", "weight", "weight"],
                 Union{Missing, Int64}[3, 3, 3],
                 Union{Missing, Int64}[4, 4, 4],
                 Union{Missing, Float64}[0.2, 0.2, 0.2],
                 Union{Missing, Float64}[0.3, 0.3, 0.2]], ["paddockId", "_variables_", "_variables__1", "red", "blue", "red_1", "blue_1"])
        @test dst == dst_t
        ds = Dataset([Union{Missing, Int64}[1, 2, 3, 4, 5, 6, 7, 8, 9, 10],
                 Union{Missing, Int64}[1, 1, 2, 5, 5, 5, 5, 2, 3, 3],
                 Union{Missing, Int64}[5, 5, 5, 5, 5, 5, 5, 5, 5, 5],
                 Union{Missing, Int64}[3, 3, 3, 3, 2, 2, 2, 2, 1, 1],
                 Union{Missing, Int64}[1, 1, 5, 5, missing, missing, missing, missing, missing, missing],
                 Union{Missing, Int64}[2, 2, 5, 5, missing, missing, missing, missing, missing, missing],
                 Union{Missing, Int64}[3, 3, 3, 3, missing, missing, missing, missing, missing, missing]], ["id", "a1_l1", "a2_l1", "a3_l1", "a1_l2", "a2_l2", "a3_l2"])
        ds_sum = combine(ds, 2:ncol(ds) => x->count(isequal(5),x)/IMD.n(x))
        dst = select!(transpose(ds_sum, (r"l1", r"l2")), r"_c1")
        dst_t = Dataset(_c1 = [.4, 1.0, 0], _c1_1 = [.5, .5, 0])
        @test dst == dst_t
end

@testset "flatten  column" begin
    ds_vec = Dataset(a = [1, 2], b = [[1, 2], [3, 4]])
    ds_tup = Dataset(a = [1, 2], b = [(1, 2), (3, 4)])
    ref = Dataset(a = [1, 1, 2, 2], b = [1, 2, 3, 4])
    @test flatten(ds_vec, :b) == flatten(ds_tup, :b) == ref
    @test flatten(ds_vec, "b") == flatten(ds_tup, "b") == ref
    @test flatten!(copy(ds_vec), :b) == flatten(ds_tup, :b) == ref
    @test flatten!(copy(ds_vec), "b") == flatten(ds_tup, "b") == ref
    ds_mixed_types = Dataset(a = [1, 2], b = [[1, 2], ["x", "y"]])
    ref_mixed_types = Dataset(a = [1, 1, 2, 2], b = [1, 2, "x", "y"])
    @test flatten(ds_mixed_types, :b) == ref_mixed_types
    @test flatten!(copy(ds_mixed_types), :b) == ref_mixed_types

    ds_three = Dataset(a = [1, 2, 3], b = [[1, 2], [10, 20], [100, 200, 300]])
    ref_three = Dataset(a = [1, 1, 2, 2, 3, 3, 3], b = [1, 2, 10, 20, 100, 200, 300])
    @test flatten(ds_three, :b) == ref_three
    @test flatten(ds_three, "b") == ref_three
    @test flatten!(copy(ds_three), :b) == ref_three
    @test flatten!(copy(ds_three), "b") == ref_three
    ds_gen = Dataset(a = [1, 2], b = [(i for i in 1:5), (i for i in 6:10)])
    ref_gen = Dataset(a = [fill(1, 5); fill(2, 5)], b = collect(1:10))
    @test flatten(ds_gen, :b) == ref_gen
    @test flatten(ds_gen, "b") == ref_gen
    @test flatten!(copy(ds_gen), :b) == ref_gen
    @test flatten!(copy(ds_gen), "b") == ref_gen
    ds_miss = Dataset(a = [1, 2], b = [Union{Missing, Int}[1, 2], Union{Missing, Int}[3, 4]])
    ref = Dataset(a = [1, 1, 2, 2], b = [1, 2, 3, 4])
    @test flatten(ds_miss, :b) == ref
    @test flatten(ds_miss, "b") == ref
    @test flatten!(copy(ds_miss), :b) == ref
    @test flatten!(copy(ds_miss), "b") == ref
    v1 = [[1, 2], [3, 4]]
    v2 = [[5, 6], [7, 8]]
    v = [v1, v2]
    ds_vec_vec = Dataset(a = [1, 2], b = v)
    ref_vec_vec = Dataset(a = [1, 1, 2, 2], b = [v1 ; v2])
    @test flatten(ds_vec_vec, :b) == ref_vec_vec
    @test flatten(ds_vec_vec, "b") == ref_vec_vec
    @test flatten!(copy(ds_vec_vec), :b) == ref_vec_vec
    @test flatten!(copy(ds_vec_vec), "b") == ref_vec_vec
    ds_cat = Dataset(a = [1, 2], b = [CategoricalArray([1, 2]), CategoricalArray([1, 2])])
    ds_flat_cat = flatten(ds_cat, :b)
    ref_cat = Dataset(a = [1, 1, 2, 2], b = [1, 2, 1, 2])
    @test ds_flat_cat == ref_cat
    @test ds_flat_cat.b.val isa CategoricalArray
    flatten!(ds_cat, :b)
    ref_cat = Dataset(a = [1, 1, 2, 2], b = [1, 2, 1, 2])
    @test ds_cat == ref_cat
    @test ds_cat.b.val isa CategoricalArray

    ds = Dataset(a = [1, 2], b = [[1, 2], [3, 4]], c = [[5, 6], [7, 8]])
    @test flatten(ds, []) == ds
    @test flatten!(copy(ds), []) == ds

    ref = Dataset(a = [1, 1, 2, 2], b = [1, 2, 3, 4], c = [5, 6, 7, 8])
    @test flatten(ds, [:b, :c]) == ref
    @test flatten(ds, [:c, :b]) == ref
    @test flatten(ds, ["b", "c"]) == ref
    @test flatten(ds, ["c", "b"]) == ref
    @test flatten(ds, 2:3) == ref
    @test flatten(ds, r"[bc]") == ref
    @test flatten(ds, Not(:a)) == ref
    @test flatten(ds, Between(:b, :c)) == ref

    @test flatten!(copy(ds), [:b, :c]) == ref
    @test flatten!(copy(ds), [:c, :b]) == ref
    @test flatten!(copy(ds), ["b", "c"]) == ref
    @test flatten!(copy(ds), ["c", "b"]) == ref
    @test flatten!(copy(ds), 2:3) == ref
    @test flatten!(copy(ds), r"[bc]") == ref
    @test flatten!(copy(ds), Not(:a)) == ref
    @test flatten!(copy(ds), Between(:b, :c)) == ref

    ds_allcols = Dataset(b = [[1, 2], [3, 4]], c = [[5, 6], [7, 8]])
    ref_allcols = Dataset(b = [1, 2, 3, 4], c = [5, 6, 7, 8])
    @test flatten(ds_allcols, :) == ref_allcols
    @test flatten!(copy(ds_allcols), :) == ref_allcols
    ds_bad = Dataset(a = [1, 2], b = [[1, 2], [3, 4]], c = [[5, 6], [7]])
    @test_throws ArgumentError flatten(ds_bad, [:b, :c])
    @test_throws ArgumentError flatten!(copy(ds_bad), [:b, :c])
    ds_vec = Dataset(a = [1, missing], b = [[1, missing], [3, 4]])
    ds_tup = Dataset(a = [1, missing], b = [(1, missing), (3, 4)])
    ref = Dataset(a = [1, 1, missing, missing], b = [1, missing, 3, 4])
    @test flatten(ds_vec, :b) == flatten(ds_tup, :b) == ref
    @test flatten(ds_vec, "b") == flatten(ds_tup, "b") == ref
    @test flatten!(copy(ds_vec), :b) == flatten(ds_tup, :b) == ref
    @test flatten!(copy(ds_vec), "b") == flatten(ds_tup, "b") == ref

    ds_cat = Dataset(a = [1, 2], b = [PooledArray([1, 2]), PooledArray([1, 2])])
    repeat!(ds_cat, 1000)
    ds_flat_cat = flatten(ds_cat, :b)
    ref_cat = Dataset(a = repeat([1, 1, 2, 2],1000), b = repeat([1, 2, 1, 2],1000))
    @test ds_flat_cat == ref_cat
    flatten!(ds_cat, :b)
    @test ds_cat == ref_cat

    ds_cat = Dataset(a = [1, 2], b = [CategoricalArray([1, 2]), CategoricalArray([1, 2])])
    repeat!(ds_cat, 1000)
    ds_flat_cat = flatten(ds_cat, :b)
    ref_cat = Dataset(a = repeat([1, 1, 2, 2],1000), b = repeat([1, 2, 1, 2],1000))
    @test ds_flat_cat == ref_cat
    flatten!(ds_cat, :b)
    @test ds_cat == ref_cat

    ds = Dataset(x=1:4, y=["ab,bc","d","ef,gh",missing])
    fmt(x) = split(x, ",")
    fmt(::Missing) = missing
    setformat!(ds, 2 => fmt)
    @test flatten(ds, :y) == Dataset([[1, 1, 1, 1, 1, 2, 3, 3, 3, 3, 3, 4], Union{Missing, Char}['a', 'b', ',', 'b', 'c', 'd', 'e', 'f', ',', 'g', 'h', missing]], [:x, :y])
    @test flatten(ds, :y, mapformats = true) == Dataset([Union{Missing, Int64}[1, 1, 2, 3, 3, 4], Union{Missing, SubString{String}}["ab", "bc", "d", "ef", "gh", missing]], [:x, :y])
    @test flatten(view(ds, :, [2,1]), :y, mapformats = true) == Dataset(reverse([Union{Missing, Int64}[1, 1, 2, 3, 3, 4], Union{Missing, SubString{String}}["ab", "bc", "d", "ef", "gh", missing]]), reverse([:x, :y]))
    flatten!(ds, :y, mapformats = true)
    @test ds == Dataset([Union{Missing, Int64}[1, 1, 2, 3, 3, 4], Union{Missing, SubString{String}}["ab", "bc", "d", "ef", "gh", missing]], [:x, :y])

    for i in 1:10
            ds = Dataset(x=rand(10000),y=[rand(1:100, rand(1:5)) for _ in 1:10000],z=rand(["12,34", "2312,343","32423,,343", missing], 10000))
            fmt__(x) = split(x, ",")
            fmt__(::Missing) = missing
            setformat!(ds, :z => fmt__)
            @test flatten(ds, :y) == flatten!(copy(ds), :y)
            @test flatten(ds, :z, mapformats = true) == flatten!(copy(ds), :z, mapformats = true)
            @test flatten(view(ds, [1,2,5,10], [2,3,1]), :y) == flatten!(ds[[1,2,5,10], [2,3,1]], :y)
            @test flatten(view(ds, [1,2,5,10], [2,3,1]), :z, mapformats = true) == flatten!(ds[[1,2,5,10], [2,3,1]], :z, mapformats = true)
            @test flatten(ds, :y) == flatten!(copy(ds), :y, threads = false)
            @test flatten(ds, :z, mapformats = true, threads = false) == flatten!(copy(ds), :z, mapformats = true)
            @test flatten(view(ds, [1,2,5,10], [2,3,1]), :y, threads = true) == flatten!(ds[[1,2,5,10], [2,3,1]], :y, threads = false)
            @test flatten(view(ds, [1,2,5,10], [2,3,1]), :z, mapformats = true, threads = false) == flatten!(ds[[1,2,5,10], [2,3,1]], :z, mapformats = true, threads = true)
    end
    ds = Dataset(x=1:3, y=[[1,2],[],[[],[]]])
    @test flatten(flatten(ds,2), 2) == Dataset(x=[1,1], y=[1,2])
    flatten!(ds, 2)
    flatten!(ds, :y)
    @test ds == Dataset(x=[1,1], y=[1,2])
end

@testset "transpose - views" begin
        pop = Dataset(country = ["c1","c1","c2","c2","c3","c3"],
                                sex = [1, 2, 1, 2, 1, 2],
                                pop_2000 = [100, 120, 150, 155, 170, 190],
                                pop_2010 = [110, 120, 155, 160, 178, 200],
                                pop_2020 = [115, 130, 161, 165, 180, 203])
        fmt(x) = x == 1 ? "male" : "female"
        setformat!(pop, :sex=>fmt)
        modify!(pop, [:country, :sex] => byrow(x->join(x))=>:c_s)
        pop_v1 = view(pop, :, 5:-1:1)
        pop_v2 = view(pop, [1,2,3,4], [1,2,3,6])

        @test transpose(pop_v2, :pop_2000) == Dataset(_variables_ = "pop_2000", _c1 = 100, _c2=120, _c3=150, _c4=155)
        @test transpose(pop_v2, :pop_2000, id = :c_s) == Dataset(_variables_ = "pop_2000", c11 = 100, c12=120, c21=150, c22=155)
        @test transpose(pop_v2, :pop_2000, id = [:country, :sex], renamecolid = x->join(x), mapformats = false) == Dataset(_variables_ = "pop_2000", c11 = 100, c12=120, c21=150, c22=155)
        @test transpose(copy(pop_v2), :pop_2000, id = :c_s) == Dataset(_variables_ = "pop_2000", c11 = 100, c12=120, c21=150, c22=155)

        poptm = Dataset([["c1", "c1", "c1", "c2", "c2", "c2", "c3", "c3", "c3"],
                ["2000", "2010", "2020", "2000", "2010", "2020", "2000", "2010", "2020"],
                Union{Missing, Int64}[100, 110, 115, 150, 155, 161, 170, 178, 180],
                Union{Missing, Int64}[120, 120, 130, 155, 160, 165, 190, 200, 203]],
                [:country,:year,:male,:female])

        @test transpose(gatherby(pop_v1, :country), [3,2,1], id = :sex, variable_name = "year", renamerowid = x->replace(x, "pop_"=>"")) == poptm
        @test transpose(gatherby(pop_v1, :country), [3,2,1], id = :sex, variable_name = "year", renamerowid = x->replace(x, "pop_"=>""), mapformats = false) == rename(poptm, ["male"=>"1", "female"=>"2"])
        @test transpose(groupby(pop_v1, :country), (1,2), id = :sex, renamecolid = (x,y)->x * "|" * y[1]) == Dataset([Union{Missing, String}["c1", "c2", "c3"], Union{Missing, Int64}[115, 161, 180], Union{Missing, Int64}[130, 165, 203], Union{Missing, Int64}[110, 155, 178], Union{Missing, Int64}[120, 160, 200]], ["country", "male|pop_2020", "female|pop_2020", "male|pop_2010", "female|pop_2010"])
end


text:
const AbstractTensor3{T} = AbstractArray{T,3}
const AbstractTensor5{T} = AbstractArray{T,5}

abstract type AbstractElementShape end
struct Line <: AbstractElementShape end
struct Quad <: AbstractElementShape end
struct Tri <: AbstractElementShape end
struct Hex <: AbstractElementShape end
struct Wed <: AbstractElementShape end
struct Pyr <: AbstractElementShape end
struct Tet <: AbstractElementShape end


text:
function compute(n::Int, digits::AbstractString)::Int
    digits = [parse(Int, digit) for digit ∈ replace(digits, "\n" => "")]
    max_product = 0
    for i ∈ 1:length(digits) - n + 1
        if 0 ∈ digits[i:i+n-1] continue end
        cur_product = prod(digits[i:i+n-1])
        if cur_product > max_product max_product = cur_product end
    end
    return max_product
end

max_stars_repo_name: Dynortice/Project-Euler

max_issues_repo_name: Dynortice/Project-Euler

max_forks_repo_name: Dynortice/Project-Euler


text:
using AdventOfCodeSolutions
using Test

function input(puzzle::Puzzle{2020, 6, n}) where n
    io = openInput(puzzle)
    return read(io, String)
end

function parseInput(input)
    groups = split(input, "\n\n")
    return map(group -> split(group, '\n', keepempty = false), groups)
end

function solve(::Puzzle{2020, 6, 1}, input)
    sum(
        answers -> string(answers...) |> unique |> length,
        input
    )
end

function solve(::Puzzle{2020, 6, 2}, input)
    return sum(
        answers -> intersect(answers...) |> length,
        input
    )
end

@testset "Day 6" begin
    @test(solve(Puzzle(6, 1), parseInput(input(Puzzle(6, 1)))) == 7128)
    @test(solve(Puzzle(6, 2), parseInput(input(Puzzle(6, 2)))) == 3640)
end


text:
__precompile__()

"""
Module for identifying low-rank subsets and subgraphs
"""
module LowRank
using Combinatorics, LinearAlgebra, StatsBase, Parameters
import Base.show
import Base.findmin, Base.findmax

include("generate_datasets.jl")
include("utils.jl")
include("trajectory.jl")
include("low_rank_subset.jl")

export find_low_rank_subset_checkall, find_low_rank_subset_sample_rep, FastLowRank

end # module

# TODO:
# allow for maximization of ∑ σᵢ² for i ≦ k instead of ∑ σᵢ² for i ≧ k+1
# allow for normalization in FastLowRank
# write a test that introduces Gaussian noise to each entry
# write spatial methods
# methods for objectives with normalization?
# write a method to find many low-rank subsets, as opposed to just one


text:
"""
Log info msg
"""
function log_info(msg)
    #st = string(stacktrace()[2])
    #file_name = lstrip(st[last(findlast("at ",st)):end])
    #file_name = split(file_name,":")
    # file_name = lstrip(st[last(findlast("at ",st)):last(findlast(":",st))-1])
    #method_name = lstrip(st[1:last(findfirst("(",st))-1])
    # logger = get_logger(file_name * "|" * method_name)
    # @info(file_name * "]" * msg)
    #@info msg _module="" _file= replace(first(file_name),".jl"=>"") _line=parse(Int64,last(file_name))
    @info msg
end

"""
logs an info msg for tick and agent_id passed
"""
function log_info(t,aid,msg)
    # if t > 50
        @info "[$(t)]($(aid)) $msg"
    # end
end

"""
logs only for a given agent
"""
function log_info(t,aid,only_id,msg)
        if aid in only_id
            @info "[$(t)]($(aid)) $msg"
        end
end

"""
    logs an info msg for tick passed
"""
function log_info(t,msg)
    # if t > 50
        @info "[$(t)] $msg"
    # end
end

function log_debug(t,aid,msg)
    @debug "[$(t)]($(aid)) $msg"
end

function log_debug(t,msg)
    @debug "[$(t)] $msg"
end

function log_debug(msg)
    @debug msg
end


text:
# imports
using Pkg; Pkg.activate(".")
using GRASS
using Statistics

# plotting imports
using LaTeXStrings
import PyPlot; plt = PyPlot; mpl = plt.matplotlib; plt.ioff()
using PyCall; animation = pyimport("matplotlib.animation");
mpl.style.use(GRASS.moddir * "figures/fig.mplstyle")

# define some functions
include(GRASS.moddir * "figures/fig_functions.jl")

# get command line args and output directories
run, plot = parse_args(ARGS)
grassdir, plotdir, datadir = check_plot_dirs()

# figure 2 -- cleaned + extrapolated input data
function main()
    # set number of curves to plot
    ncurves = 25

    # get input data
    bisinfo = GRASS.SolarData(relative=true, extrapolate=true)
    key = (:c, :mu10)
    bis = bisinfo.bis[key]
    wav = bisinfo.wav[key]
    dep = bisinfo.dep[key]
    wid = bisinfo.wid[key]

    # create colormap and colorbar
    nstp = round(Int, size(bis,2) / ncurves)
    iter = 1:nstp:size(bis,2)
    iter = 1:2:30
    cmap = plt.get_cmap("Blues", length(iter))
    cols = [cmap(1*i/length(iter)) for i in 1:length(iter)]
    sm = plt.cm.ScalarMappable(cmap=cmap, norm=mpl.colors.Normalize(vmin=0, vmax=15 * iter[end]/60.0))

    # plot the bisectors
    fig = plt.figure()
    ax1 = fig.add_subplot()
    for (i, t) in enumerate(iter)
        ax1.plot(wav[:,t][2:end].*1000, bis[:,t][2:end], c=cols[i], alpha=0.75)
    end

    # shade upper region
    xlims = ax1.get_xlim()
    ax1.fill_between(range(xlims..., step=0.1), 0.8, 1.0, hatch="/",
                     fc="black", ec="white", alpha=0.15, zorder=0)

    # set axis limits + labels
    ax1.set_xlim(xlims...)
    ax1.set_ylim(0.1, 1.0)
    ax1.set_xlabel(L"{\rm Relative\ Wavelength\ (m\AA)}")
    ax1.set_ylabel(L"{\rm Normalized\ Intensity}")

    # set color bar
    cb = fig.colorbar(sm)
    cb.set_label(L"{\rm Time\ from\ first\ observation\ (min)}")

    # save the plot
    fig.savefig(plotdir * "fig2a.pdf")
    plt.clf(); plt.close()
    println(">>> Figure written to: " * plotdir * "fig2a.pdf")

    # plot the widths
    fig = plt.figure()
    ax1 = fig.add_subplot()
    for (i, t) in enumerate(iter)
        ax1.plot(dep[:,t], wid[:,t], c=cols[i], alpha=0.75)
    end

    # set axis limits + labels
    ax1.set_xlabel(L"{\rm Normalized\ Intensity}")
    ax1.set_ylabel(L"{\rm Width\ across\ line\ (\AA)}")

    # set color bar
    cb = fig.colorbar(sm)
    cb.set_label(L"{\rm Time\ from\ first\ observation\ (min)}")

    # save the plot
    fig.savefig(plotdir * "fig2b.pdf")
    plt.clf(); plt.close()
    println(">>> Figure written to: " * plotdir * "fig2b.pdf")
    return nothing
end

if (run | plot)
    main()
end

max_stars_repo_name: palumbom/GRASS

max_issues_repo_name: palumbom/GRASS

max_forks_repo_name: palumbom/GRASS


text:
using Plots

"""
    default_plotlyjs!()

Set plotting default to that used in Physical Biology of the Cell, 2nd edition.
"""
function default_plotlyjs!()
    plotlyjs(
        background_color="#E3DCD0",
        background_color_outside="white",
        foreground_color_grid="#ffffff",
        gridlinewidth= 0.5,
        guidefontfamily="Lucida Sans Unicode",
        guidefontsize=8,
        tickfontfamily="Lucida Sans Unicode",
        tickfontsize=7,
        titlefontfamily="Lucida Sans Unicode",
        titlefontsize=8,
        dpi=300,
        linewidth=1.25,
        legendtitlefontsize=8,
        legendfontsize=8,
        foreground_color_legend="#E3DCD0",
        color_palette=:seaborn_colorblind
    )
end


text:
using Test
using LogicCircuits

@testset "CNF file parser tests" begin

   circuit = zoo_cnf("iscas89/s386.scan.cnf")

   @test circuit isa PlainLogicCircuit
   @test circuit isa Plain⋀Node
   @test circuit.children[1] isa Plain⋁Node
   @test circuit.children[end] isa Plain⋁Node
   @test num_variables(circuit) == 172
   @test num_children(circuit) == 506

   circuit = zoo_cnf("8.cnf")

   function test8(circuit)
      @test circuit isa PlainLogicCircuit
      @test circuit isa Plain⋀Node
      @test circuit.children[1] isa Plain⋁Node
      @test circuit.children[end] isa Plain⋁Node
      @test num_variables(circuit) == 227
      @test num_nodes(circuit) == 1168
   end

   test8(circuit)

   mktempdir() do tmp

      temp_path = "$tmp/test8.cnf"
      write(temp_path, circuit)
      circuit2 = read(temp_path, LogicCircuit)
      test8(circuit2)

      cnf_orig = zoo_cnf("easy/C17_mince.cnf")
      @test_nowarn write("$tmp/temp.cnf", cnf_orig)
      cnf_new = read("$tmp/temp.cnf", LogicCircuit)

      @test num_variables(cnf_orig) == num_variables(cnf_new)
      @test all(x isa Plain⋁Node || x isa PlainLiteralNode for x in cnf_new.children)
      @test all(all(x isa PlainLiteralNode for x in or_node.children) for or_node in or_nodes(cnf_new))
      @test tree_formula_string(cnf_orig) == tree_formula_string(cnf_new)

      @test_nowarn write("$tmp/temp.cnf.gz", cnf_orig)
      cnf_new = read("$tmp/temp.cnf.gz", LogicCircuit)
      @test tree_formula_string(cnf_orig) == tree_formula_string(cnf_new)

  end

end

@testset "DNF file parser tests" begin

   circuit = zoo_dnf("8.dnf")

   function test8(circuit)
      @test circuit isa PlainLogicCircuit
      @test circuit isa Plain⋁Node
      @test circuit.children[1] isa Plain⋀Node
      @test circuit.children[end] isa Plain⋀Node
      @test num_variables(circuit) == 227
      @test num_nodes(circuit) == 1168
   end

   test8(circuit)

   mktempdir() do tmp

      temp_path = "$tmp/test8.dnf"
      write(temp_path, circuit)
      circuit2 = read(temp_path, LogicCircuit)
      test8(circuit2)

      write("$temp_path.gz", circuit)
      circuit2 = read("$temp_path.gz", LogicCircuit)
      test8(circuit2)

  end

end

max_issues_repo_name: Juice-jl/LogicCircuits

max_forks_repo_name: Juice-jl/LogicCircuits


text:
# Copyright 2020 Google LLC
#
# Licensed under the Apache License, Version 2.0 (the "License");
# you may not use this file except in compliance with the License.
# You may obtain a copy of the License at
#
#     http://www.apache.org/licenses/LICENSE-2.0
#
# Unless required by applicable law or agreed to in writing, software
# distributed under the License is distributed on an "AS IS" BASIS,
# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
# See the License for the specific language governing permissions and
# limitations under the License.

""" Returns a small test LP.
The LP:
  min 5 x_1 + 2 x_2 + x_3 +   x_4 - 14 s.t.
  2 x_1 +   x_2 + x_3 + 2 x_4  = 12
    x_1 +         x_3         >=  7
                  x_3 -   x_4 >=  1
  0 <= x_1 <= 2
  0 <= x_2 <= 4
  0 <= x_3 <= 6
  0 <= x_4 <= 3

Optimum solutions:
  Primal: x_1 = 1, x_2 = 0, x_3 = 6, x_4 = 2. Value: 5 + 0 + 6 + 2 - 14 = -1.
  Dual: [0.5, 4.0, 0.0]  Value: 6 + 28 - 3.5*6 - 14 = -1
"""
function example_lp()
  return FirstOrderLp.linear_programming_problem(
    [0.0, 0.0, 0.0, 0.0],  # variable_lower_bound
    [2.0, 4.0, 6.0, 3.0],  # variable_upper_bound
    [5.0, 2.0, 1.0, 1.0],  # objective_vector
    -14.0,                 # objective_constant
    [
      2.0 1.0 1.0 2.0
      1.0 0.0 1.0 0.0
      0.0 0.0 1.0 -1.0
    ],                     # constraint_matrix
    [12.0, 7.0, 1.0],      # right_hand_side
    1,                     # num_equalities
  )
end

""" Returns a one-variable LP with a no variable bounds.
The LP:
  min -x_1 s.t.
  -x_1 >= -2 (as a linear constraint)

Optimum solution:
  Primal: x_1 = 2, Value: -2.
  Dual: [1.0]  Value: -2.
"""
function example_lp_without_bounds()
  return FirstOrderLp.linear_programming_problem(
    [-Inf],                # variable_lower_bound
    [Inf],                 # variable_upper_bound
    [-1.0],                # objective_vector
    0.0,                   # objective_constant
    reshape([-1.0], 1, 1), # constraint_matrix
    [-2.0],                # right_hand_side
    0,                     # num_equalities
  )
end

""" Returns a small test QP.
The QP:
  min 2 x_1^2 + 0.5 x_2^2 - x_1 - x_2 s.t.
  x_1 + x_2  <= 1

  0 <= x_1 <= 1
  0 <= x_2 <= 1

Optimum solutions:
  Primal: x_1 = 0.2, x_2 = 0.8.
  Dual: [0.2]  Value: 0.04 + 0.32 - 0.2 - 0.8 = -0.6
"""
function example_qp()
  return FirstOrderLp.QuadraticProgrammingProblem(
    [0.0, 0.0],     # variable_lower_bound
    [1.0, 1.0],     # variable_upper_bound
    [
      4.0 0.0
      0.0 1.0
    ],              # objective_matrix
    [-1.0, -1.0],   # objective_vector
    -0.0,           # objective_constant
    [-1.0 -1.0],    # constraint_matrix
    [-1.0],         # right_hand_side
    0,              # num_equalities
  )
end

""" Returns a small test QP.
The QP:
  min 2 x_1^2 + 0.5 x_2^2 - x_1 - x_2 s.t.
  x_1 + x_2  <= 1

  0 <= x_1 <= 1
  0 <= x_2 <= 1

Optimum solutions:
  Primal: x_1 = 0.25, x_2 = 0.0.
  Dual: [0.0]  Value: 0.125 + 0 - 0.25 - 0 = -0.125
"""
function example_qp2()
  return FirstOrderLp.QuadraticProgrammingProblem(
    [0.0, 0.0],     # variable_lower_bound
    [1.0, 1.0],     # variable_upper_bound
    [
      4.0 0.0
      0.0 1.0
    ],              # objective_matrix
    [-1.0, 1.0],    # objective_vector
    -0.0,           # objective_constant
    [-1.0 -1.0],    # constraint_matrix
    [-1.0],         # right_hand_side
    0,              # num_equalities
  )
end

"""Returns a correlation clustering LP.
This is the LP for minimizing disagreements for correlation clustering for the
4-vertex graph
   1 - 3 - 4
   | /
   2
In integer solutions x_ij is 1 if i and j are in the same cluster and 0
otherwise. The 6 variables are in the order
 x_12, x_13, x_14, x_23, x_24, x_34.
For any distinct i,j,k there's a triangle inequality
  (1-x_ik) <= (1-x_ij) + (1-x_jk)
i.e.
  -x_ij - x_jk + x_ik >= -1.
For brevity we only include 3 out of the 12 possible triangle inequalities: two
needed in the optimal solution and 1 other."""
function example_cc_lp()
  return FirstOrderLp.linear_programming_problem(
    [0.0, 0.0, 0.0, 0.0, 0.0, 0.0],     # variable_lower_bound
    [1.0, 1.0, 1.0, 1.0, 1.0, 1.0],     # variable_upper_bound
    [-1.0, -1.0, 1.0, -1.0, 1.0, -1.0], # objective_vector
    4.0,                                # objective_constant
    [
      0.0 -1.0 1.0 0.0 0.0 -1.0
      0.0 0.0 0.0 -1.0 1.0 -1.0
      -1.0 -1.0 0.0 1.0 0.0 0.0
    ],                                  # constraint_matrix
    [-1.0, -1.0, -1.0],                 # right_hand_side
    0,                                  # num_equalities
  )
end

"""Returns another 4-vertex correlation clustering LP.

The variables are x_12, x_13, x_14, x_23, x_24, and x_34.
This time the graph is a star centered at vertex 1.
Only the three triangle inequalities that are needed are included."""
function example_cc_star_lp()
  return FirstOrderLp.linear_programming_problem(
    [0.0, 0.0, 0.0, 0.0, 0.0, 0.0],    # variable_lower_bound
    [1.0, 1.0, 1.0, 1.0, 1.0, 1.0],    # variable_upper_bound
    [-1.0, -1.0, -1.0, 1.0, 1.0, 1.0], # objective_vector
    3.0,                               # objective_constant
    [
      -1.0 -1.0 0.0 1.0 0.0 0.0
      -1.0 0.0 -1.0 0.0 1.0 0.0
      0.0 -1.0 -1.0 0.0 0.0 1.0
    ],      # constraint_matrix
    [-1.0, -1.0, -1.0],                # right_hand_side
    0,                                 # num_equalities
  )
end

""" Returns a small test LP with rows that are linearly dependent.
The LP:
  min x_1 + 2 x_2 + 3 x_3 +  4 x_4 s.t.
    x_1 +     x_2 +   x_3 +   x_4  == 2
    x_1 +     x_2 +   x_3 +   x_4  == 2
    x_1 +                     x_4  == 1
  0 <= x_1
  0 <= x_2
  0 <= x_3
  0 <= x_4

Optimum solution:
  Primal: x_1 = 1, x_2 = 1, x_3 = 0, x_4 = 0.
  Dual solution is not unique.
  Value: 1 + 2 = 3
"""
function example_lp_dependent_rows()
  return FirstOrderLp.linear_programming_problem(
    [0.0, 0.0, 0.0, 0.0],  # variable_lower_bound
    [Inf, Inf, Inf, Inf],  # variable_upper_bound
    [1.0, 2.0, 3.0, 4.0],  # objective_vector
    0.0,                   # objective_constant
    [
      1.0 1.0 1.0 1.0
      1.0 1.0 1.0 1.0
      1.0 0.0 0.0 1.0
    ],                     # constraint_matrix
    [2.0, 2.0, 1.0],       # right_hand_side
    3,                     # num_equalities
  )
end

""" Returns a small test LP that is primal infeasible.
The LP:
  min -x_1 + 0.5 x_2 s.t.
    -x_1 -    x_2   == 1
  0 <= x_1
  0 <= x_2
"""
function example_lp_easy_primal_infeasible()
  return FirstOrderLp.linear_programming_problem(
    [0.0, 0.0],           # variable_lower_bound
    [Inf, Inf],           # variable_upper_bound
    [1.0, 0.5],           # objective_vector
    0.0,                  # objective_constant
    [-1.0 -1.0],          # constraint_matrix
    [1.0],                # right_hand_side
    1,                    # num_equalities
  )
end

""" Returns a small test LP that is primal infeasible. Compared with
example_lp_easy_primal_infeasible() this is a hard problem. Note that
the nonnegativity constraints have no impact on feasibility.
The LP:
  min x_1 + 2 x_2 + 3 x_3 +  4 x_4 s.t.
    x_1 +     x_2                   == 1
              x_2 +   x_3           == 1
                      x_3 +    x_4  == 1
   x_1  +     x_2 +   x_3 +    x_4  == 2 + tol
  0 <= x_1
  0 <= x_2
  0 <= x_3
  0 <= x_4
"""
function example_lp_hard_primal_infeasible(tol::Float64)
  @assert tol > 0.0
  return FirstOrderLp.linear_programming_problem(
    [0.0, 0.0, 0.0, 0.0],      # variable_lower_bound
    [Inf, Inf, Inf, Inf],      # variable_upper_bound
    [1.0, 2.0, 3.0, 4.0],      # objective_vector
    0.0,                       # objective_constant
    [
      1.0 1.0 0.0 0.0
      0.0 1.0 1.0 0.0
      0.0 0.0 1.0 1.0
      1.0 1.0 1.0 1.0
    ],                         # constraint_matrix
    [1.0, 1.0, 1.0, 2 + tol],  # right_hand_side
    4,                         # num_equalities
  )
end

""" Returns a small test LP that is dual infeasible.
The LP:
  min -x_1 + 0.5 x_2 s.t.
    x_1 -     x_2   == 1
  0 <= x_1
  0 <= x_2
"""
function example_lp_dual_infeasible()
  return FirstOrderLp.linear_programming_problem(
    [0.0, 0.0],            # variable_lower_bound
    [Inf, Inf],            # variable_upper_bound
    [-1.0, 0.4],           # objective_vector
    0.0,                   # objective_constant
    [1.0 -2.0],            # constraint_matrix
    [1.0],                 # right_hand_side
    1,                     # num_equalities
  )
end


text:
module Guillotine

include("./../../constraints/geometric/constants.jl")
include("./1_equalities.jl")
include("./2_range.jl")
include("./4_and.jl")
include("./5_find_interval.jl")
using JuMP
using .Constants
using .Equalities
using .Range
using .GuillotineIfAnd
using .FindInterval

export cons_ds_guillotine

function cons_ds_guillotine(m)

    m = first_dim_equality(m)
    m = second_dim_range(m)
    m = all_dimensions_and(m)
    m = find_interval(m)

    @expression(m, guillotine_objective,
        sum(m[:guillotine_equalities_objective]) +
        #sum(m[:guillotine_range_objective]) +
        sum(m[:guillotine_and_objective]) +
        0.0001 * sum(m[:guillotine_interval_objective])
    )

    return m
end

end


text:
#=
juliacaller:
- Julia version: 1.4
- Author: Mehmet Hakan Satman
- Date: 2020-04-03 =#

using Sockets
using Pkg

a = try
    v = Pkg.installed()["JSON"]
	catch
    Pkg.add("JSON")
end

using JSON

function writeln(client, str)
	write(client, str)
	write(client, "\r\n")
end

"""
Handles the client connection.

Command list:
	execute: Evaluates the given statement or expression, possible multiple items seperated with semicolon.

	get: Returns the value of a given variable in JSON format.

	exit: Closes the client connection.

	shutdown: Closes the server.

Examples:

	execute a = 99

	get a

	(Result is {"a":99})

	execute using Distributions

	execute norm = Normal(0, 1)

	execute numbers = rand(norm, 5)

	get numbers

	(Result is {"numbers":[-1.1734648123994915,-0.5577848419508594,0.5628431386136749,0.42211832659070825,-0.4402667291521316]})
"""
function handle_client(server, client)
	while true
		__line__ = readline(client)
		println(__line__)
			if startswith(__line__, "execute ")
				__command__ = __line__[9:end]
				try
					@debug Evaling __command__
					eval(Meta.parse(__command__))
				catch mth_err
					@error mth_err
					close(client)
					close(server)
				end
				# writeln(client, "eval okay")
			elseif startswith(__line__, "get ")
				__varname__ = __line__[5:end]
				__D__ = Dict(__varname__ => eval(Meta.parse(__varname__)))
				writeln(client, json(__D__))
			elseif startswith(__line__, "exit")
				break
			elseif startswith(__line__, "install ")
				__pkg__ = __line__[9:end]
				Pkg.add(__Pkg__)
			elseif startswith(__line__, "shutdown")
				close(client)
				close(server)
				break
			end
	end
	close(client)
end

"""
Creates a TCP server socket and listens on a given port.
# Arguments
- `PORT::Integer`: The port number for the server socket, default is 8000.
"""
function serve(PORT=8000)
	server = listen(PORT)
	println("Listening JuliaCaller on port $PORT")
	while true
		try
			client = accept(server)
			if (@isdefined client)
				handle_client(server, client)
    			else
				close(client)
			end
		catch err
			close(client)
			close(server)
			writeln("Closed connection and server, exiting")
			break
		end
	end
end

max_stars_repo_name: subes/juliacaller

max_issues_repo_name: subes/juliacaller

max_forks_repo_name: subes/juliacaller


text:
using Documenter, MetaGraphs

# index is equal to the README for the time being
cp(normpath(@__FILE__, "../../README.md"), normpath(@__FILE__, "../src/index.md"); force=true)

# same for license
cp(normpath(@__FILE__, "../../LICENSE.md"), normpath(@__FILE__, "../src/license.md"); force=true)

makedocs(
    modules = [MetaGraphs],
    format = :html,
    sitename = "MetaGraphs",
    pages    = Any[
        "Overview"             => "index.md",
        "MetaGraphs Functions" => "metagraphs.md",
        "License Information"  => "license.md"
    ]
)

deploydocs(
    deps=nothing,
    make=nothing,
    repo="github.com/JuliaGraphs/MetaGraphs.jl.git",
    target="build",
    julia="0.6",
    osname = "linux"
)

rm(normpath(@__FILE__, "../src/index.md"))
rm(normpath(@__FILE__, "../src/license.md"))


text:
export
    DriverModel,

    StaticDriver,

    get_name,
    action_type,
    set_desired_speed!,
    observe!,
    reset_hidden_state!,
    prime_with_history!

abstract DriverModel{DriveAction}

get_name(::DriverModel) = "???"
action_type{A}(::DriverModel{A}) = A
set_desired_speed!(model::DriverModel, v_des::Float64) = model # do nothing by default
reset_hidden_state!(model::DriverModel) = model # do nothing by default
observe!{S,D,I,R}(model::DriverModel, scene::EntityFrame{S,D,I}, roadway::R, egoid::Int) = model  # do nothing by default
Base.rand(model::DriverModel) = error("rand not implemented for model $model")
Distributions.pdf{A}(model::DriverModel{A}, a::A) = error("pdf not implemented for model $model")
Distributions.logpdf{A}(model::DriverModel{A}, a::A) = error("logpdf not implemented for model $model")

function prime_with_history!{S,D,I,R}(
    model::DriverModel,
    trajdata::ListRecord{Entity{S,D,I}},
    roadway::R,
    frame_start::Int,
    frame_end::Int,
    egoid::I,
    scene::EntityFrame{S,D,I} = allocate_frame(trajdata),
    )

    reset_hidden_state!(model)

    for frame in frame_start : frame_end
        get!(scene, trajdata, frame)
        observe!(model, scene, roadway, egoid)
    end

    return model
end
function prime_with_history!{S,D,I,R}(model::DriverModel, rec::EntityQueueRecord{S,D,I}, roadway::R, egoid::I;
    pastframe_start::Int=1-nframes(rec),
    pastframe_end::Int=0,
    )

    reset_hidden_state!(model)

    for pastframe in pastframe_start : pastframe_end
        scene = rec[pastframe]
        observe!(model, scene, roadway, egoid)
    end

    model
end

####

type StaticDriver{A,P<:ContinuousMultivariateDistribution} <: DriverModel{A}
    distribution::P
end

get_name(::StaticDriver) = "StaticDriver"
function Base.rand{A,P}(model::StaticDriver{A,P})
    a = rand(model.distribution)
    return convert(A, a)
end
Distributions.pdf{A}(model::StaticDriver{A}, a::A) = pdf(model.distribution, convert(Vector{Float64}, a))
Distributions.logpdf{A}(model::StaticDriver{A}, a::A) = logpdf(model.distribution, convert(Vector{Float64}, a))


text:
module SetPyPlot

# setting pyplot parameters to default values
using PyPlot

export setrcparam, plotyy

"""
setrcparam(;tick_font_size = 8,
legend_font_size = 8,
axis_label_font_size = 9,
fig_width = 9,
fig_height = 7)

sets the font size parameters to my favorite default values
"""
function setrcparam(;tick_font_size = 8,
    legend_font_size = 8,
    axis_label_font_size = 9,
    fig_width = 9,
    fig_height = 7)

    rc("xtick", labelsize=tick_font_size)
    rc("ytick", labelsize=tick_font_size)
    rc("legend", fontsize=legend_font_size)
    rc("axes", labelsize=axis_label_font_size)
    rc("figure", dpi = 300)
    rc("figure", figsize = (fig_width, fig_height)./2.54)
    rc("savefig", bbox = "tight")
end

"""
function plotyy(x1, y1, x2, y2; x_label = "", y1_label = "", y2_label = "", y1_style = "-bo", y2_style = "--vr")
"""
function plotyy(x1, y1, x2, y2; fig_size = (3.5, 2.5), x_label = "", y1_label = "", y2_label = "", y1_style = "-bo", y2_style = "--vr")
    fig, ax1 = subplots(figsize = fig_size)
    ax1.plot(x1, y1, y1_style)
    ax1.set_ylabel(y1_label, color="b")
    ax1.set_xlabel(x_label, color="b")
    ax1.tick_params("y", colors = "b")

    ax2 = ax1.twinx()
    ax2.plot(x2, y2, y2_style)
    ax2.set_ylabel(y2_label, color="r")
    ax2.tick_params("y", colors = "r")

    tight_layout()

    return fig, ax1, ax2
end
# more to come later...

end # module

max_stars_repo_name: simulkade/geothermal

max_issues_repo_name: simulkade/geothermal

max_forks_repo_name: simulkade/geothermal


text:
using Random
using Distributions
using Statistics
using SparseArrays
using MinimalRLCore

mutable struct TabularTMaze <: MinimalRLCore.AbstractEnvironment
    world::Array{Array{String}}
    start_state::Array{Int64}
    current_state::Array{Int64}
    goal_states::Array{String}
    feature_size::Int64
    step::Int
    exploring_starts::Bool
    cumulant_schedule::CumulantSchedule
    extrinsic_reward::Float64

    function TabularTMaze(exploring_starts, cumulant_schedule; extrinsic_reward = -0.01)
        feature_size = 21
        num_actions = 4
        step = 1

        world = [["G1", "0", "0", "0", "0", "0", "G3"],
                 ["1", "0", "0", "0", "0", "0", "1"],
                 ["1", "1", "1", "1", "1", "1", "1"],
                 ["1", "0", "0", "1", "0", "0", "1"],
                 ["G2", "0", "0", "1", "0", "0", "G4"],
                 ["0", "0", "0", "1", "0", "0", "0"],
                 ["0", "0", "0", "1", "0", "0", "0"],
                 ["0", "0", "0", "1", "0", "0", "0"],
                 ["0", "0", "0", "1", "0", "0", "0"]]
        goal_states = ["G1", "G2", "G3", "G4"]
        start_state = [9,4]
        new(world, start_state, [9,4], goal_states, feature_size,step, exploring_starts, cumulant_schedule, extrinsic_reward)
    end
end

Base.size(e::TabularTMaze) = e.feature_size

function generate_obs(state::Array{Int64})
    obs = zeros(Int64, 9, 7)
    obs[state[1],state[2]] = 1
    toReturn = findfirst(!iszero, obs[valid_state_mask()])
    return toReturn
end

function valid_state_mask()
    """
    return a mask of valid states that is 9x7
    """
    world = [["G1", "0", "0", "0", "0", "0", "G3"],
             ["1", "0", "0", "0", "0", "0", "1"],
             ["1", "1", "1", "1", "1", "1", "1"],
             ["1", "0", "0", "1", "0", "0", "1"],
             ["G2", "0", "0", "1", "0", "0", "G4"],
             ["0", "0", "0", "1", "0", "0", "0"],
             ["0", "0", "0", "1", "0", "0", "0"],
             ["0", "0", "0", "1", "0", "0", "0"],
             ["0", "0", "0", "1", "0", "0", "0"]]

    world = permutedims(hcat(world...))
    valid_states = findall(x-> x!="0", world)
    return valid_states
end

function valid_state_action_mask()
"""
returns a mask for the state-action mask (63,4)
"""
    num_actions = 4
    world = [["G1", "0", "0", "0", "0", "0", "G3"],
                      ["1", "0", "0", "0", "0", "0", "1"],
                      ["1", "1", "1", "1", "1", "1", "1"],
                      ["1", "0", "0", "1", "0", "0", "1"],
                      ["G2", "0", "0", "1", "0", "0", "G4"],
                      ["0", "0", "0", "1", "0", "0", "0"],
                      ["0", "0", "0", "1", "0", "0", "0"],
                      ["0", "0", "0", "1", "0", "0", "0"],
                      ["0", "0", "0", "1", "0", "0", "0"]]

    flattened_world = reshape(permutedims(hcat(world...)),(63))
    valid_states = findall(x-> x!="0", flattened_world)
    return valid_states
end

MinimalRLCore.get_reward(env::TabularTMaze) = env.extrinsic_reward
is_terminal(env::TabularTMaze, pos) = env.world[pos[1]][pos[2]] in env.goal_states
MinimalRLCore.is_terminal(env::TabularTMaze) = is_terminal(env, env.current_state)
MinimalRLCore.get_actions(env::TabularTMaze) = 1:4

function MinimalRLCore.get_state(env::TabularTMaze)
    obs = generate_obs(env.current_state)
    cumulants = get_cumulants(env, env.cumulant_schedule, env.current_state)

    return vcat(obs,cumulants)
end

function get_cumulants(env::TabularTMaze, cs::CumulantSchedule, pos)
    num_cumulants = 4
    cumulants = zeros(num_cumulants)
    if env.world[pos[1]][pos[2]] == "G1"
        cumulants[1] = get_cumulant(cs, "G1")
    elseif env.world[pos[1]][pos[2]] == "G2"
        cumulants[2] = get_cumulant(cs, "G2")
    elseif env.world[pos[1]][pos[2]] == "G3"
        cumulants[3] = get_cumulant(cs, "G3")
    elseif env.world[pos[1]][pos[2]] == "G4"
        cumulants[4] = get_cumulant(cs, "G4")
    end
    return cumulants
end

function MinimalRLCore.reset!(environment::TabularTMaze, rng::AbstractRNG=Random.GLOBAL_RNG)
    if environment.exploring_starts == false
        environment.current_state = environment.start_state
    elseif environment.exploring_starts == true
        possible_start_states = findall(x -> x == "1", hcat(environment.world...))
        start_state = possible_start_states[rand(rng, 1:length(possible_start_states))]
        # The start states are flipped since hcat "transposes" a list of list when converting it to a matrix
        environment.current_state = [start_state[2],start_state[1]]
    end
end

function MinimalRLCore.reset!(environment::TabularTMaze,start_state::Int)
    s = valid_state_mask()[start_state]
    environment.current_state = [s[1],s[2]]
end

function MinimalRLCore.reset!(environment::TabularTMaze, start_state::CartesianIndex)
    # throw("Implement env_start with a start_state")
    # environment.current_state = [start_state[2], start_state[1]]
    environment.current_state = [start_state[1], start_state[2]]
end

function MinimalRLCore.environment_step!(environment::TabularTMaze, action, rng::AbstractRNG=Random.GLOBAL_RNG)
    actions = [(-1, 0), (0, 1), (1, 0), (0, -1)] # up, right, down, left
    reward = 0.0
    terminal = false

    potential_state = [environment.current_state[1] + actions[action][1],
                       environment.current_state[2] + actions[action][2]]

    if potential_state[1] in [0, 10] || potential_state[2] in [0, 8] # No nothing, agent doesn't move
    elseif environment.world[potential_state[1]][potential_state[2]] == "1"
        environment.current_state = potential_state
    elseif is_terminal(environment, potential_state)
        environment.current_state = potential_state
        terminal = true
    end
    # update!(environment, environment.cumulant_schedule, environment.current_state)
    update!(environment.cumulant_schedule, environment.current_state)
end

using RecipesBase, Colors
@recipe function f(env::TabularTMaze)
    ticks := nothing
    foreground_color_border := nothing
    grid := false
    legend := false
    aspect_ratio := 1
    xaxis := false
    yaxis := false
    yflip := false

    SIZE=20
    BG = Colors.RGB(1.0, 1.0, 1.0)
    BORDER = Colors.RGB(0.0, 0.0, 0.0)
    WALL = Colors.RGB(0.3, 0.3, 0.3)
    AC = Colors.RGB(0.69921875, 0.10546875, 0.10546875)
    GOAL = Colors.RGB(0.796875, 0.984375, 0.76953125)

    cell = fill(BG, SIZE, SIZE)
    cell[1, :] .= BORDER
    cell[end, :] .= BORDER
    cell[:, 1] .= BORDER
    cell[:, end] .= BORDER

    s_y = length(env.world)
    s_x = length(env.world[1])

    screen = fill(BG, (s_y + 2)*SIZE, (s_x + 2)*SIZE)

    screen[:, 1:SIZE] .= WALL
    screen[1:SIZE, :] .= WALL
    screen[end-(SIZE-1):end, :] .= WALL
    screen[:, end-(SIZE-1):end] .= WALL

    for j ∈ 1:s_x
        for i ∈ 1:s_y
            sqr_i = ((i)*SIZE + 1):((i+1)*SIZE)
            sqr_j = ((j)*SIZE + 1):((j+1)*SIZE)
            if env.current_state[1] == i && env.current_state[2] == j
                v = @view screen[sqr_i, sqr_j]
                v .= cell
                v[Int(SIZE/2)-4:Int(SIZE/2)+5, Int(SIZE/2)-4:Int(SIZE/2)+5] .= AC
            elseif env.world[i][j] == "0"
                screen[sqr_i, sqr_j] .= WALL
            elseif env.world[i][j] == "1"
                screen[sqr_i, sqr_j] .= cell
            elseif env.world[i][j][1] == 'G'
                v = @view screen[sqr_i, sqr_j]
                v .= GOAL
                v[1, :] .= BORDER
                v[:, 1] .= BORDER
                v[end, :] .= BORDER
                v[:, end] .= BORDER
            end
        end
    end
    screen[end:-1:1,:]
end

max_stars_repo_name: MatthewMcLeod/curiosity

max_issues_repo_name: MatthewMcLeod/curiosity

max_forks_repo_name: MatthewMcLeod/curiosity


text:
"""
	isholiday(calendar, dt)

Checks if `dt` is a holiday based on a given `calendar` of holidays.

`calendar` can be an instance of `HolidayCalendar`,  a `Symbol` or an `AbstractString`.

Returns boolean values.
"""
function isholiday(hc::HolidayCalendar, dt::Date)
    error("isholiday for $(hc) not implemented.")
end

isholiday(::NullHolidayCalendar, dt::Date) = false
isholiday(calendar, dt::Date) = isholiday(convert(HolidayCalendar, calendar), dt)

include("calendars/brazil.jl")
include("calendars/uk.jl")
include("calendars/us.jl")


text:
pos(b::Buffer)  = b.cursor.pos
poss(b::Buffer) = "$(pos(b))"
y(b::Buffer)    = b.cursor.pos[1]
ys(b::Buffer)   = "$(y(b))"
x(b::Buffer)    = b.cursor.pos[2]
xs(b::Buffer)   = "$(x(b))"
xu(b::Buffer)   = chr2ind(line(b), x(b))
posu(b::Buffer) = [y(b), xu(b)]

function after(b::Buffer)
    if mode(b)==normal_mode
        after_normal(b)
    elseif mode(b)==insert_mode
        after_insert(b)
    elseif mode(b)==delete_mode
        after_delete(b)
    elseif mode(b)==yank_mode
        after_yank(b)
    elseif mode(b)==find_mode
        after_find(b)
    elseif mode(b)==search_mode
        after_search(b)
    elseif mode(b)==go_mode
        after_go(b)
    elseif mode(b)==replace_mode
        #after_replace(b)
    end
end

top(b::Buffer) = parse(Int64, b.state[:top])
bottom(b::Buffer) = parse(Int64, b.state[:bottom])
left(b::Buffer) = parse(Int64, b.state[:left])
right(b::Buffer) = parse(Int64, b.state[:right])
function height(b::Buffer)
    length(b.text)
end
function line(b::Buffer, line::Integer)
    b.text[line]
end
line(b::Buffer) = line(b, y(b))
function setline(b::Buffer, line::Integer, s::String)
    b.text[line] = s
end
setline(b::Buffer, s::String) = setline(b, y(b), s)

function width(b::Buffer)
    length(line(b))
end

function escape(b::Buffer)
    setmode(b, normal_mode)
    clamp!(b)
end

function clamp!(b::Buffer, edgecase=false)
    if height(b) == 0
        settext(b, [""])
    end
    b.cursor.pos[1] = Base.clamp(y(b), 1, height(b))
    if width(b) > 0 || edgecase
        b.cursor.pos[2] = Base.clamp(x(b), 1, width(b)+edgecase)
    else
        b.cursor.pos[2] = 0
    end
end

function clamp_range(r::UnitRange, low::Number, high::Number)
    Base.clamp(r.start, low, high):Base.clamp(r.stop, low, high)
end
function clamp_range(r::StepRange, low::Number, high::Number)
    Base.clamp(r.start, low, high):step(r):Base.clamp(r.stop, low, high)
end

function findsymbol(s::String, c::Char, n::Integer)
    p = findn([s...] .== c)
    if n>length(p)
        return 0
    else
        return p[n]
    end
end
function findsymbol(b::Buffer, c::Char, pos, n::Integer)
    pos[2] + findsymbol(line(b,pos[1])[pos[2]+1:end], c, n)
end
findsymbol(b::Buffer, c::Char, n::Integer) = findsymbol(b, c, posu(b), n)
findsymbol(b::Buffer, c::Char) = findsymbol(b, c, parse_n(b))

isint(c::Integer) = 47 < c < 58
function parse_n(args::Array)
    n = length(args)
    il = args[indexin([false], isint.(Int.(args)))[1]+1:end]
    if !isempty(il)
        parse(Int, join(il))
    else
        return 1
    end
end
function parse_n(b::Buffer)
    parse_n(b.args)
end

# Should be cleaned up
function set_boundaries(b::Buffer, ind::Integer, first::Symbol, last::Symbol)
    #remove parsing from here
    f, l = parse.([b.state[first], b.state[last]])
    if f > ind
        d = ind - f
    elseif l < ind
        d = ind - l
    else
        d = 0
    end
    b.state[first], b.state[last] = string.([f, l]+d)
end
function set_boundaries(b::Buffer)
    set_boundaries(b, y(b), :top, :bottom)
    df = right(b)-left(b)
    b.state[:left], b.state[:right] = string.([max(1, x(b)-df), max(df+1, x(b))])
end
function resize(b::Buffer, height, width)
    b.state[:top], b.state[:bottom] = string.([top(b), top(b) + height])
    b.state[:left], b.state[:right] = string.([left(b), left(b) + width])
end

function resize(b::Buffer)
    try
        term = Base.REPL.Terminals.TTYTerminal(get(ENV, "TERM", "dumb"), STDIN, STDOUT, STDERR)
        c = Base.Terminals.width(term)
        r = Base.Terminals.height(term)
        resize(b, r-3, c-1)
    catch z
        b.state[:log] = "$z"
    end
end

function clear_arg(b::Buffer)
    resize!(b.args,0)
end

function deleteat(b::Buffer, pos, n)
    pr = b.text[pos[1]]
    b.text[pos[1]] = string(unirange(pr, 1:pos[2]-1), unirange(pr, pos[2]+n:length(pr)))
end

function paste_lines(b::Buffer, ln,
                     s::Array{T,1} where T <: Union{String, SubString{String}})
    settext(b, [b.text[1:ln-1];
                s;
                b.text[ln:end]])
end
paste_lines(b::Buffer, ln, s::String) = paste_lines(b, ln, split(s, '\n'))
function paste_lines(b::Buffer, s::Union{Array{String, 1}, String})
    paste_lines(b, y(b), s)
end
function pastea_lines(b::Buffer, s::String)
    paste(lines, b, y(b)+1, s)
end

function paste_single(b::Buffer, pos, s::String)
    setline(b,
            pos[1],
            string((ln -> unirange(ln, 1:pos[2]-1))(line(b, pos[1])),
                   s,
                   (ln -> unirange(ln, pos[2]:length(ln)))(line(b, pos[1])),
                  )
           )
end
function paste(b::Buffer, pos, s::String)
    if s[1] == '\e' #Line copies add an escape symbol before the text
        b.state[:log] = "$pos"
        paste_lines(b, pos[1], s[2:end])
    else
        s_array = split(s, '\n')
        if length(s_array) == 1
            paste_single(b::Buffer, pos, s)
        else
            ln1 = string((t->unirange(t, 1:pos[2]))(b.text[pos[1]]), s_array[1])
            ln2 = string(s_array[end], (t->unirange(t, pos[2]+1:length(t)))(b.text[pos[1]]))
            settext(b, [b.text[1:pos[1]-1];
                        ln1;
                        s_array[2:end-1];
                        ln2;
                        b.text[pos[1]+1:end]
                       ])
        end
    end
end
function paste(b::Buffer, pos, c::Char)
    if c in keys(b.registers)
        paste(b, pos, b.registers[c])
    else
        self.state[:log] = "Register $c empty"
    end
end
paste(b::Buffer, s) = paste(b, pos(b), s)
function pastea(b::Buffer, pos, s::String)
    if s[1] =='\e'
        paste(b, pos+[1, 0], s)
    else
        paste(b, pos+[0, 1], s)
    end
end
function pastea(b::Buffer, pos, c::Char)
    if c in keys(b.registers)
        pastea(b, pos, b.registers[c])
    else
        self.state[:log] = "Register $c empty"
    end
end
pastea(b::Buffer, s) = pastea(b::Buffer, pos(b), s)

function join_lines(b::Buffer, range::Range, delimiter=" ")
    range = clamp_range(range, 1, height(b))
    jl = join(b.text[range], delimiter)
    deleteat!(b.text, range[2:end])
    b.text[range[1]] = jl
end
function join_lines(b::Buffer, interval::Tuple, delimiter=" ")
    join_lines(b, interval[1]:interval[2], delimiter)
end

function splitlines(b::Buffer, pos)
    nl1 = unirange(b.text[pos[1]], 1:pos[2]-1)
    nl2 = unirange(b.text[pos[1]], max(1, pos[2]):length(b.text[pos[1]]))
    b.text[pos[1]] = nl1
    insert!(b.text, pos[1]+1, nl2)
end

function delete_lines(b::Buffer, range::Range)
    range = clamp_range(range, 1, height(b))
    yank_lines(b, range)
    deleteat!(b.text, range)
end

function order_pos(pos1, pos2)
    if (pos1[1] < pos2[1]) || ((pos1[1] == pos2[1]) && (pos1[2] < pos2[2]))
        return (pos1, pos2)
    else
        return (pos2, pos1)
    end
end

#Caution: Doesn't work exactly like vim (line delete etc)
function delete_between(b::Buffer, pos1, pos2, reg='"')
    #no unicode support for this right now
    sp, ep = order_pos(pos1, pos2)
    yank_between(b, sp, ep, reg)
    nl = string(unirange(b.text[sp[1]], 1:sp[2]-1),
                unirange(b.text[ep[1]], ep[2]:length(b.text[ep[1]])))
    deleteat!(b.text, sp[1]:ep[1]-1)
    b.text[sp[1]] = nl
    b.cursor.pos .= sp
end

function yank_between(b::Buffer, pos1, pos2, reg='"')
    sp, ep = order_pos(pos1, pos2)
    if sp[1] == ep[1]
        b.registers[reg] = string(b.text[sp[1]][sp[2]:ep[2]-1])
    else
        b.registers[reg] = join([unirange(b.text[sp[1]], sp[2]:length(b.text[sp[1]]));
                                b.text[sp[1]+1:ep[1]-1];
                                b.text[ep[1]][1:ep[2]-1]], '\n')
    end
end

function undo(b::Buffer)
    if b.state[:undo]=="undo"
        setredo(b, b.text)
        settext(b, b.undo)
    else
        b.state[:log] = "Multiple undos not supported"
    end
    after(b)
end
function redo(b::Buffer)
    if b.state[:undo]=="redo"
        settext(b, b.redo)
    else
        b.state[:log] = "Already at latest state"
    end
    after(b)
end

function replay(b::Buffer, actions::String, n=1)
    for i in 1:n
        for c in actions
            handle_raw(b, c)
        end
    end
end
replay(b::Buffer, s::Array{Char, 1}, n=1) = replay(b, join(s), n)

source(b::Buffer) = evalcmd(b, join(b.text, '\n'))

function reconfigure(b::Buffer)
    pkdir = Pkg.Dir.path("Joy")
    files = ["include(\"$pkdir/src/base.jl\")",
     "include(\"$pkdir/src/string_utils.jl\")",
     "include(\"$pkdir/src/utils.jl\")",
     "include(\"$pkdir/src/syntax.jl\")",
     "include(\"$pkdir/src/words.jl\")",
     "include(\"$pkdir/src/modes/movements.jl\")",
     "include(\"$pkdir/src/modes/normal.jl\")",
     "include(\"$pkdir/src/modes/insert.jl\")",
     "include(\"$pkdir/src/modes/calls.jl\")",
     "include(\"$pkdir/src/modes/command.jl\")",
     "include(\"$pkdir/src/modes/search.jl\")",
     "include(\"$pkdir/src/config.jl\")",
     "include(\"$pkdir/src/editor.jl\")",
    ]
    (f -> evalcmd(b, f)).(files)
end

function quit(b::Buffer)
    b.state[:running] = "false"
end


text:
struct PairwisePotential
    V
    params
    max_radius
end

@doc raw"""
Pairwise terms: Pairwise potential between nuclei, e.g., Van der Waals potentials, such as
Lennard—Jones terms.
The potential is dependent on the distance between to atomic positions and the pairwise
atomic types:
For a distance `d` between to atoms `A` and `B`, the potential is `V(d, params[(A, B)])`.
The parameters `max_radius` is of `100` by default, and gives the maximum (reduced) distance
between nuclei for which we consider interactions.
"""
function PairwisePotential(V, params; max_radius=100)
    params = Dict(minmax(key[1], key[2]) => value for (key, value) in params)
    PairwisePotential(V, params, max_radius)
end
function (P::PairwisePotential)(basis::PlaneWaveBasis{T}) where {T}
    E = energy_pairwise(basis.model, P.V, P.params; P.max_radius)
    TermPairwisePotential(P.V, P.params, T(P.max_radius), E)
end

struct TermPairwisePotential{TV, Tparams, T} <:Term
    V::TV
    params::Tparams
    max_radius::T
    energy::T
end

function ene_ops(term::TermPairwisePotential, basis::PlaneWaveBasis, ψ, occ; kwargs...)
    (E=term.energy, ops=[NoopOperator(basis, kpt) for kpt in basis.kpoints])
end

@timing "forces: Pairwise" function compute_forces(term::TermPairwisePotential,
                                                   basis::PlaneWaveBasis{T}, ψ, occ;
                                                   kwargs...) where {T}
    forces = zero(basis.model.positions)
    energy_pairwise(basis.model, term.V, term.params; term.max_radius, forces)
    forces
end

"""
Compute the pairwise interaction energy per unit cell between atomic sites. If `forces` is
not nothing, minus the derivatives of the energy with respect to `positions` is computed.
The potential is expected to decrease quickly at infinity.
"""
function energy_pairwise(model::Model{T}, V, params; kwargs...) where {T}
    isempty(model.atoms) && return zero(T)
    symbols = Symbol.(atomic_symbol.(model.atoms))
    energy_pairwise(model.lattice, symbols, model.positions, V, params; kwargs...)
end

# This could be factorised with Ewald, but the use of `symbols` would slow down the
# computationally intensive Ewald sums. So we leave it as it for now.
function energy_pairwise(lattice, symbols, positions, V, params;
                         max_radius=100, forces=nothing)
    T = eltype(lattice)
    @assert length(symbols) == length(positions)

    if forces !== nothing
        @assert size(forces) == size(positions)
        forces_pairwise = copy(forces)
    end

    # The potential V(dist) decays very quickly with dist = ||A (rj - rk - R)||,
    # so we cut off at some point. We use the bound  ||A (rj - rk - R)|| ≤ max_radius
    # where A is the real-space lattice, rj and rk are atomic positions.
    poslims = [maximum(rj[i] - rk[i] for rj in positions for rk in positions) for i in 1:3]
    Rlims = estimate_integer_lattice_bounds(lattice, max_radius, poslims)

    # Check if some coordinates are not used.
    is_dim_trivial = [norm(lattice[:,i]) == 0 for i=1:3]
    max_shell(n, trivial) = trivial ? 0 : n
    Rlims = max_shell.(Rlims, is_dim_trivial)

    #
    # Energy loop
    #
    sum_pairwise::T = zero(T)
    # Loop over real-space
    for R1 in -Rlims[1]:Rlims[1], R2 in -Rlims[2]:Rlims[2], R3 in -Rlims[3]:Rlims[3]
        R = Vec3(R1, R2, R3)
        for i = 1:length(positions), j = 1:length(positions)
            # Avoid self-interaction
            R == zero(R) && i == j && continue
            ai, aj = minmax(symbols[i], symbols[j])
            param_ij = params[(ai, aj)]
            Δr = lattice * (positions[i] - positions[j] - R)
            dist = norm(Δr)
            energy_contribution = V(dist, param_ij)
            sum_pairwise += energy_contribution
            if forces !== nothing
                # We use ForwardDiff for the forces
                dE_ddist = ForwardDiff.derivative(d -> V(d, param_ij), dist)
                dE_dti = lattice' * dE_ddist / dist * Δr
                forces_pairwise[i] -= dE_dti
                forces_pairwise[j] += dE_dti
            end
        end # i,j
    end # R
    energy = sum_pairwise / 2  # Divide by 2 (because of double counting)
    if forces !== nothing
        forces .= forces_pairwise ./ 2
    end
    energy
end


text:
export State, StateError
export @stateidx

import Base: <, >, ==, !=

# ==========================================================================
# Constants for the State

const STATES = [:error, :uninitialized, :initialized, :training, :trystop, :trained, :converged]
const STATES_LOOKUP = Dict{Symbol,Int}()
for i in 1:length(STATES)
  STATES_LOOKUP[STATES[i]] = i
end

macro stateidx(state)
  idx = STATES_LOOKUP[eval(state)]
  :($idx)
end

# ==========================================================================
# State type for convenience and comparison

immutable State
  _idx::Int
  _symb::Symbol

  State(idx::Int) = new(idx, STATES[idx])
  State(symb::Symbol) = new(STATES_LOOKUP[symb], symb)
end

Base.show(io::IO, state::State) = print(io, state._idx, ": ", state._symb)
Base.convert(::Type{State}, idx::Int) = State(idx)
Base.convert(::Type{State}, symb::Symbol) = State(symb)

for op = (:<, :>, :(==), :(!=), :(<=), :(>=))
  @eval begin
    ($op)(s1::State, s2::State) = ($op)(s1._idx, s2._idx)
    ($op)(s::State,  i::Int)    = ($op)(s._idx, i)
    ($op)(i::Int,    s::State)  = ($op)(i, s._idx)
    ($op)(s::State,  k::Symbol) = ($op)(s._idx, STATES_LOOKUP[k])
    ($op)(k::Symbol, s::State)  = ($op)(STATES_LOOKUP[k], s._idx)
  end
end

# ==========================================================================
# State exception used to signal that a function was called in an illegal state

type StateError <: Exception
  msg::String
end


text:
println("------------------------------------")
println("Generic SUNIrrep tests")
println("------------------------------------")
@timedtestset "Basic tests for SUNIrrep{$N}:" for N = 2:5
    I1 = SUNIrrep(tuple(sort(rand(1:9, N); rev=true)..., 1))
    I2 = SUNIrrep(tuple(sort(rand(1:9, N); rev=true)..., 1))
    @constinferred dim(I1)
    d = 0
    for (I, NI) in @constinferred directproduct(I1, I2)
        d += NI*dim(I)
    end
    @test d == dim(I1)*dim(I2)
    if dim(I2) < dim(I1)
        I1, I2 = I2, I1
    end
    b = @constinferred basis(I1)
    (v, s) = @constinferred Nothing iterate(b)
    @constinferred Nothing iterate(b, s)
    v = @constinferred collect(b)
    for i = 1:length(v)-1
        @test isless(v[i], v[i+1])
    end
    @test v[end] == highest_weight(I1)
    s = sprint(show, v[1])
    for k = 1:N
        @test parse(Int, s[17 + 2*k]) == weight(I1)[k]
    end
end