text:
# plotting functions
using Plots, Debugger, Images, Colors
gr()

struct GifPlot
	anim::Animation
	fname::String
end
GifPlot(fname::String) = GifPlot(Animation(), fname)

"""
initialises a plot
"""
function init_plot()
	p = plot(framestyle=:none)
end

"""
Plots the neurons and their static connections
"""
function networkplot(p, nn, connections::SynapseLayers)
	@infiltrate
	ycoords, xcoords = neuronplot!(p, nn)
	for i in 1:length(nn)-1
		y = (ycoords[i], ycoords[i+1])
		x = (xcoords[i], xcoords[i+1])
		plotconnections!(p, y, x, connections[i])
	end

	return p
end
function networkplot(p, nn, connections::SynapseLayers, weights::SynapseLayers; maxneurons=10)

	nn = min.(maxneurons, nn)
	ycoords, xcoords = neuronplot!(p, nn)
	@infiltrate
	for i in 1:length(nn)-1
		y = (ycoords[i], ycoords[i+1])
		x = (xcoords[i], xcoords[i+1])

		viewcons = connections.layers[i].data[1:nn[i],1:nn[i+1]]
		viewweights = weights.layers[i].data[1:nn[i],1:nn[i+1]]

		plotconnections!(p, y, x, viewcons, viewweights)
	end

	return p
end
function networkplot(p, nn, m::MatrixTypes; maxneurons=10)

	nn = min.(maxneurons, nn)
	@infiltrate
	ycoords, xcoords = neuronplot!(p, nn, m, maxneurons)

	for i in 1:length(nn)-1
		y = (ycoords[i], ycoords[i+1])
		x = (xcoords[i], xcoords[i+1])

		viewweights = m.weights.layers[i].data[1:nn[i],1:nn[i+1]] .> 0
		viewspiked = m.spiked.layers[i].data[1:nn[i]] .* viewweights
		plotconnections!(p, y, x, viewweights, viewspiked)
	end

	return p
end

"""
checks the time and shows the plot if appropriate
"""
function shownetwork(p,t, nn, m::MatrixTypes; framerate=1)
	if t%framerate==0
		display(networkplot(p, nn, m))
	end
end
function gifnetwork(gf::GifPlot, p, t, nn, m::MatrixTypes)
	networkplot(p,nn,m)
	frame(gf.anim)
	return gf
end

"""
This returns a plot of circles marking location of neurons.
"""
function neuronplot!(p, nn)

	numlayers = length(nn)
	maxn = maximum(nn)

	xs = fill.([1:numlayers;], nn)
	xcoords = reduce(vcat, xs)

	f(x) = x==1 ? [round(maxn/2)] : [range(1,stop=maxn, length=x);]

	ys = map(x -> f(x),nn)
	ycoords = reduce(vcat, map(x -> f(x),nn))
	markersizes = reduce(vcat, fill.(neuronmarkersize(nn), nn))

	maxmkr = maximum(markersizes)
	maxy = [-maxmkr/5, maxn + maxmkr/5]
	maxx = [1-maxmkr/50, numlayers + maxmkr/50]
	scatter!(p, xcoords, ycoords, markercolor=:white, markersize=markersizes, legend=false, ylims=maxy, xlims=maxx)
	return (ys, xs)
end

#function neuronplot!(p, nn, m::MatrixTypes, maxneurons::Array)
#	for (max_n, num_n) in maxneurons, nn
#		neuronplot!(p, num_n, m, max_n)
#	end
#end

function neuronplot!(p, nn, m::MatrixTypes, maxneurons)

	numlayers = length(nn)
	maxn = maximum(nn)

	xs = fill.([1:numlayers;], nn)
	xcoords = reduce(vcat, xs)

	mcolours = RGB.(vcat(m.input.layers[1][1:maxneurons] .> 0, falses(length(xs[2])), m.spiked.layers[3][1:min(length(m.spiked.layers[3]),maxneurons)]))

	f(x) = x==1 ? [round(maxn/2)] : [range(1,stop=maxn, length=x);]

	ys = map(x -> f(x),nn)
	ycoords = reduce(vcat, map(x -> f(x),nn))
	markersizes = reduce(vcat, fill.(neuronmarkersize(nn), nn))

	maxmkr = maximum(markersizes)
	maxy = [-maxmkr/5, maxn + maxmkr/5]
	maxx = [1-maxmkr/50, numlayers + maxmkr/50]

	scatter!(p, xcoords, ycoords, markercolor=mcolours, markersize=markersizes, legend=false, ylims=maxy, xlims=maxx)
	return(ys, xs)
	# now we need to get just the inputted ones and colour them in white.
end

"""
Plot a line for each synapse, and flash if spiked
"""
function plotconnections!(p, ycoords, xcoords, connections)

	allconnections = vcat.(ycoords[1], ycoords[2]')
	allxs = vcat.(xcoords[1], xcoords[2]')

	linexs = allxs[connections]
	lineys = allconnections[connections]

	for i in 1:length(linexs)
		plot!(p,linexs[i],lineys[i],legend=false, linecolor=RGB(0.5,0.5,0.5), linealpha=0.9)
	end
end
function plotconnections!(p, ycoords, xcoords, weights::BitArray{<:Any}, spikers::BitArray{<:Any})

	allconnections = vcat.(ycoords[1], ycoords[2]')
	allxs = vcat.(xcoords[1], xcoords[2]')

	spikexs = allxs[spikers]
	spikeys= allconnections[spikers]

	linexs = allxs[weights]
	lineys = allconnections[weights]
	linewidths = weights[weights]

	for i in 1:length(linexs)
		iszero(linewidths[i]) && continue
		plot!(p,linexs[i],lineys[i], legend=false, linecolor=RGB(0.5,0.5,0.5), linealpha=0.9, linewidth=linewidths[i])
	end

	plot!(p, spikexs, spikeys, legend=false, linecolor=RGB(0.9,0.1,0.1))
end
function plotconnections!(p, ycoords, xcoords, connections::BitArray{<:Any}, weights::BitArray{<:Any}, spikers::BitArray{<:Any})

	allconnections = vcat.(ycoords[1], ycoords[2]')
	allxs = vcat.(xcoords[1], xcoords[2]')

	spikexs = allxs[spikers]
	spikeys= allconnections[spikers]

	linexs = allxs[connections]
	lineys = allconnections[connections]
	linewidths = weights[connections]

	for i in 1:length(linexs)
		iszero(linewidths[i]) && continue
		plot!(p,linexs[i],lineys[i], legend=false, linecolor=RGB(0.5,0.5,0.5), linealpha=0.9, linewidth=linewidths[i])
	end

	for i in 1:length(spikexs)
		iszero(spikexs[i]) && continue
		plot!(p, spikexs[i], spikeys[i], legend=false, linecolor=RGB(0.9,0.1,0.1))
	end
end
plotconnections!(p, ycoords, xcoords, connections::SynapseLayer, weights::SynapseLayer) = plotconnections!(p, ycoords, xcoords, connections.data .==1, weights.data)
plotconnections!(p, ycoords, xcoords, weights::Array{<:Number,<:Any}, spikers::Array{<:Number,<:Any}) = plotconnections!(p, ycoords, xcoords, weights .> 0.0, spikers)
plotconnections!(p, ycoords, xcoords, connections::Array{<:Number,<:Any}, weights::Array{<:Number,<:Any}, spiked::Array{<:Number, <:Any}) = plotconnections!(p, ycoords, xcoords, connections .==1, weights, spiked .== 1)
plotconnections!(p, ycoords, xcoords, connections::SynapseLayer) = plotconnections!(p, ycoords, xcoords,connections.data .== 1)

"""
for axis with limits 1-n, what is the ideal marker size for n neurons to all show up
"""
neuronmarkersize(numneurons) = 150 ./ (numneurons .+ 2)

"""
given a spikeplot, or some other neuronlayers object, plots a mean value of all neurons in each layer for each layer
"""
function plotmeanlayer(data::Array{<:BrainTypes,N}) where {N}
	p = plot()
	plotmeanlayer!(p,data)
	return p
end
function plotmeanlayer!(p, data::Array{<:BrainTypes,N}) where {N}

	numlayers = length(data[1].layers)
	layercolors = range(HSL(0,0,0), HSL(1,1,1), length=numlayers)
	means = meanlayer(data)
	plot!(p,means)
	return p
end
function meanlayer(data::Array{<:BrainTypes,N}) where {N}

	numsteps = length(data)
	numlayers = length(data[1].layers)
	out = Array{Float64,2}(undef, numsteps,numlayers)

	for i = 1:numsteps
		for l = 1:numlayers
			out[i,l] = mean(data[i].layers[l])
		end
	end
	return out
end

function display_sensory(p, sensory::AbstractInput)
	heatmap!(p, reshape(sensory,length(sensory),1), color=:grays, axis=nothing)
end

function display_sensory(sensory::InputSequence; fname="sensory_sequence.gif")

	p = heatmap(axis=nothing, color=:grays)
	gf = GifPlot(pwd() * "\\output\\figures\\" * fname)
	for s in sensory
		display_sensory(p, s)
		frame(gf.anim)
	end
	gif(gf.anim,gf.fname,fps=2)
end

function bar(ydata::Vector{normdata})

	means = [dp.mean for dp in ydata]
	std = [dp.std for dp in ydata]
	Plots.bar(means, yerror=std)
end

function bar!(ydata::Vector{normdata})

	means = [dp.mean for dp in ydata]
	std = [dp.std for dp in ydata]
	Plots.bar!(means, yerror=std)
end

function scatter(ydata::Array{normdata})

	means = [dp.mean for dp in ydata]
	std = [dp.std for dp in ydata]
	Plots.scatter(means[1], means[2], means[3])
end

function raster(spikes::Array{NeuronLayers{Bool,1},1})

	total_length = length(spikes)
	numlayers = length(spikes[1].layers)
	raster_out = [zeros(Bool, length(layer),total_length) for layer in spikes[1].layers ]

	for t in 1:total_length

		for l in 1:numlayers
			raster_out[l][:,t] .= spikes[t].layers[l]
		end
	end
	return raster_out
end

function raster(spikes::Array{NeuronLayers{Bool,1},1}, layer::Int)

	total_length = length(spikes)
	raster_out = zeros(Bool, length(spikes[1].layers[layer]),total_length)

	for t in 1:total_length

			raster_out[:,t] .= spikes[t].layers[layer]

	end
	return raster_out
end


text:
#=
  All interval problem in Julia ConstraintSolver.jl

  CSPLib problem number 7
  http://www.cs.st-andrews.ac.uk/~ianm/CSPLib/prob/prob007/index.html
  """
  Given the twelve standard pitch-classes (c, c , d, ...), represented by
  numbers 0,1,...,11, find a series in which each pitch-class occurs exactly
  once and in which the musical intervals between neighbouring notes cover
  the full set of intervals from the minor second (1 semitone) to the major
  seventh (11 semitones). That is, for each of the intervals, there is a
  pair of neigbhouring pitch-classes in the series, between which this
  interval appears. The problem of finding such a series can be easily
  formulated as an instance of a more general arithmetic problem on Z_n,
  the set of integer residues modulo n. Given n in N, find a vector
  s = (s_1, ..., s_n), such that (i) s is a permutation of
  Z_n = {0,1,...,n-1}; and (ii) the interval vector
  v = (|s_2-s_1|, |s_3-s_2|, ... |s_n-s_{n-1}|) is a permutation of
  Z_n-{0} = {1,2,...,n-1}. A vector v satisfying these conditions is
  called an all-interval series of size n; the problem of finding such
  a series is the all-interval series problem of size n. We may also be
  interested in finding all possible series of a given size.
  """

  Model created by Hakan Kjellerstrand, hakank@gmail.com
  See also my Julia page: http://www.hakank.org/julia/

=#
using ConstraintSolver, JuMP
using Cbc, GLPK
const CS = ConstraintSolver
include("constraints_utils.jl")

#
# Note: my_abs is defined in constraints_utils.jl
#
function all_interval(n=5,all_solutions=true,print_solutions=true)

    cbc_optimizer = optimizer_with_attributes(Cbc.Optimizer, "logLevel" => 0)
    glpk_optimizer = optimizer_with_attributes(GLPK.Optimizer)

    model = Model(optimizer_with_attributes(CS.Optimizer,   "all_solutions"=> all_solutions,
                                                            # "all_optimal_solutions"=>true,
                                                            "logging"=>[],

                                                            # "traverse_strategy"=>:BFS,
                                                            "traverse_strategy"=>:DFS, # <-
                                                            # "traverse_strategy"=>:DBFS,

                                                            # "branch_split"=>:Smallest,
                                                            # "branch_split"=>:Biggest,
                                                            "branch_split"=>:InHalf, # <-

                                                            # "simplify"=>false,
                                                            # "simplify"=>true, # default

                                                            "time_limit"=>16,

                                                            # "lp_optimizer" => cbc_optimizer, # 45.8s
                                                            "lp_optimizer" => glpk_optimizer, # 21.5s
                                        ))
    n1 = n-1
    @variable(model, 1 <= x[1:n] <= n, Int)
    @variable(model, 1 <= diffs[1:n1] <= n1, Int)

    @constraint(model, x in CS.AllDifferent())
    @constraint(model, diffs in CS.AllDifferent())

    # diffsb = @variable(model, [1:n1], Bin)
    for k in 1:n1
        # These don't work!
        # @constraint(model, diffs[k] == abs(x[k+1] - x[k])) # abs() is not supported!
        # @constraint(model, diffs[k] == (x[k+1] >= x[k]) ? x[k+1] - x[k] : x[k] - x[k+1]) # This don't work!

        # This works!
        # @constraint(model, diffsb[k] := {x[k+1] >= x[k]})
        # @constraint(model, diffsb[k] := {diffs[k] == x[k+1] - x[k]})
        # @constraint(model, !diffsb[k] := {diffs[k] == x[k] - x[k+1]})

        # This works as well and it seems to be about as fast/slow...
        my_abs(model,x[k],x[k+1],diffs[k])
    end

    # symmetry breaking
    @constraint(model, x[1] <= x[n1])
    @constraint(model, diffs[1] <= diffs[2])

    # Solve the problem
    optimize!(model)

    status = JuMP.termination_status(model)
    # println("status:$status")
    if status == MOI.OPTIMAL
        num_sols = MOI.get(model, MOI.ResultCount())
        println("num_sols:$num_sols\n")
        if print_solutions
            for sol in 1:num_sols
                println("solution #$sol")
                xx = convert.(Integer,JuMP.value.(x; result=sol))
                diffsx = convert.(Integer,JuMP.value.(diffs; result=sol))
                println("x:$xx diffs:$diffsx\n")
            end
        end
    else
        println("status:$status")
    end
end

all_solutions = true
print_solutions = true
@time all_interval(9,all_solutions,print_solutions)

# println("n:13")
# @time all_interval(13,false,print_solutions)

# Count the number of solutions
#=
n:3
num_sols:1

  0.001692 seconds (6.08 k allocations: 363.805 KiB)

n:4
num_sols:1

  0.002169 seconds (11.64 k allocations: 652.992 KiB)

n:5
num_sols:3

  0.003961 seconds (26.33 k allocations: 1.381 MiB)

n:6
num_sols:8

  0.013505 seconds (95.76 k allocations: 4.858 MiB)

n:7
num_sols:9

  0.047034 seconds (364.65 k allocations: 17.449 MiB)

n:8
num_sols:15

  0.206969 seconds (1.73 M allocations: 76.011 MiB)

n:9
num_sols:42

  1.171571 seconds (8.01 M allocations: 346.716 MiB, 6.85% gc time)

n:10
num_sols:104

 14.749069 seconds (54.33 M allocations: 2.278 GiB, 8.40% gc time)
=#
for n in 3:10
    println("\nn:$n")
    @time all_interval(n,true,false)
end

ext: jl

max_stars_repo_name: tias/hakank

max_issues_repo_name: tias/hakank

max_forks_repo_name: tias/hakank


text:
using PackageCompiler

"""
    lineedit(editor::Function, filename::String)

Easily edit a file, line by line.  If your `editor` function returns `nothing` the
file is not modified.  If the file does not exist, silently fails.  Usage example:

    lineedit("foo.jl") do lines
        map(lines) do l
            if occursin(r"libblas_name", l)
                return "const libblas_name = \"libfoo.so\"\n"
            end
            return l
        end
    end
"""
function lineedit(editor::Function, filename::String)
    # Silently fail for files that don't exist
    if !isfile(filename)
        return nothing
    end

    lines = open(filename) do io
        readlines(io, keep=true)
    end

    # Run user editor; if something goes wrong, just quit out
    try
        lines = editor(lines)
    catch e
        @error("Error occured while running user line editor:", e)
        return nothing
    end

    # Write it out, if the editor decides something needs to change
    if lines != nothing
        open(filename, "w") do io
            for l in lines
                write(io, l)
            end
        end
    end

    # Return the lines, just for fun
    return lines
end

function replace_libblas(base_dir, name)
    # Replace `libblas` and `liblapack` in build_h.jl
    file = joinpath(base_dir, "build_h.jl")
    lineedit(file) do lines
        @info("Replacing libblas_name in $(file)")
        return map(lines) do l
            if occursin(r"libblas_name", l)
                return "const libblas_name = $(repr(name))\n"
            elseif occursin(r"liblapack_name", l)
                return "const liblapack_name = $(repr(name))\n"
            else
                return l
            end
        end
    end
end

"""
    find_init_func_bounds(lines)

Given a list of `lines`, find the beginning and ending indicies within `lines` that
contains the `__init__()` function, if any.  If none can be found, returns `(nothing, nothing)`.
"""
function find_init_func_bounds(lines)
    start_idx = findfirst(match.(r"^function\s+__init__\(\)\s*$", lines) .!= nothing)
    if start_idx === nothing
        return (nothing, nothing)
    end

    func_len = findfirst(match.(r"^end\s*$", lines[start_idx+1:end]) .!= nothing)
    if func_len === nothing
        return (nothing, nothing)
    end
    return (start_idx, start_idx + func_len)
end

function force_proper_PATH(base_dir, dir_path)
    # We need to be compatible with Julia 1.0+, so we need to search for `__init__()`
    # within both `sysimg.jl` (for Julia 1.0 and 1.1) as well as `Base.jl` (for 1.2+)

    for file in ("sysimg.jl", "Base.jl")
        @info("Checking $(file)")
        lineedit(joinpath(base_dir, file)) do lines
            # Find the start/end of `__init__()` within this file, if possible:
            init_start, init_end = find_init_func_bounds(lines)
            if init_start === nothing
                @info("Could not find init function in $(file)")
                return nothing
            end

            # Scan the function for mentions of our `dir_path`; if it already exists,
            # then call it good, returning `nothing` so the file is not modified.
            if any(match.(r"\s+ENV\[\"PATH\"\] =", lines[init_start+1:init_end-1]) .!= nothing)
                @info("Found ENV already")
                return nothing
            end

            # If we found a function, insert our `ENV["PATH"]` mapping:
            pathsep = Sys.iswindows() ? ';' : ':'
            insert!(
                lines,
                init_start + 1,
                "    ENV[\"PATH\"] = string(ENV[\"PATH\"], $(repr(pathsep)), $(repr(dir_path)))\n",
            )
            @info("Successfully modified $(file)")
            return lines
        end
    end
end

function enable_mkl_startup(libmkl_rt)
    # First, we need to modify a few files in Julia's base directory
    base_dir = joinpath(Sys.BINDIR, Base.DATAROOTDIR, "julia", "base")

    # Replace definitions of `libblas_name`, etc.. to point to MKL
    replace_libblas(base_dir, libmkl_rt)

    # Force-setting `ENV["PATH"]` to include the location of MKL libraries
    # This is required on Windows, where we can't use RPATH
    @info("Checking if we need to update PATH...")
    if Sys.iswindows()
        force_proper_PATH(base_dir, dirname(libmkl_rt))
    end

    # Next, build a new system image
    sysimgpath = PackageCompiler.sysimgbackup_folder("native")
    if ispath(sysimgpath)
        rm(sysimgpath, recursive=true)
    end
    force_native_image!()
end

function enable_openblas_startup(libopenblas = "libopenblas")
    # First, we need to modify a few files in Julia's base directory
    base_dir = joinpath(Sys.BINDIR, Base.DATAROOTDIR, "julia", "base")

    # Replace definitions of `libblas_name`, etc.. to point to MKL
    if Sys.WORD_SIZE == 64
        libopenblas = "$(libopenblas)64_"
    end
    replace_libblas(base_dir, libopenblas)

    # Next, build a new system image
    sysimgpath = PackageCompiler.sysimgbackup_folder("native")
    if ispath(sysimgpath)
        rm(sysimgpath, recursive=true)
    end
    force_native_image!()
end


text:
# VIDEO DEMONSTRATION: https://www.youtube.com/watch?v=iFlwYqjI448

# defining the possible options the user may choose in this game
options = ["rock", "paper", "scissors", "stop", "r", "p", "s"]

#a function allowing the user to define how many rounds they would like in the game
function enterrounds()

    rounds = string(input("Please enter the number of rounds: \n"))

    #tries to convert it to integer
    try
        value = parse(Int64, rounds)

        #when converted, it checks whether the integer is larger than or equal to 1
        (value < 1) ? (

            #restarts the function by recursion if not larger than or equal to 1
            println("Please enter an integer larger than or equal to 1!");
            enterrounds()
        ) : (return value) #else the value is returned
    catch

        #restarts function by recursion if not numeric
        println("Please enter an integer!")
        enterrounds()
    end
end

# function checking whether the input is in the defined array or not
function checkinput(option)

    #asks for input again if the input isn't defined
    while (option in options) == false
        option = lowercase(string(input("Enter Rock (or R), Paper (or P), Scissors (or S), or Stop (to end the game): \n")))
    end
    option
end

# function using ternary operator to check who wins the round
function winner(playermove, aimove)
    (playermove == aimove) ? (return "Tie") : ""
    (playermove == "r" && aimove == "p") ? (return "Computer") : ""
    (playermove == "p" && aimove == "s") ? (return "Computer") : ""
    (playermove == "s" && aimove == "r") ? (return "Computer") : (return "Player")
end

#converts the single-character responses back to full respones to be printed
function myconvert(move)
    (move == "r") ? (return "rock 👊🏻") : ""
    (move == "s") ? (return "scissors ✌🏻") : (return "paper 🖐🏻")
end

#function containing the main gameplay
function main()
    #determining how many rounds the player wants
    rounds = enterrounds()
    playercount = 0
    aicount = 0
    roundscount = 0

    for i in 1:rounds

        roundscount = roundscount + 1
        println()

        #prints the current round number
        println("Round ", string(roundscount))

        #asks for the player's input
        playermove = checkinput(lowercase(string(input("Enter Rock (or R), Paper (or P), Scissors (or S), or Stop (to end the game): \n"))))

        #prints out what the player chose
        println("Player Chose: ", uppercasefirst(myconvert(string(playermove[1]))))

        #check whether the player chose to stop, exits if they did choose to stop
        (playermove == "stop") ? (exit()) : (

            aimove = options[rand(1:3)];

            #prints out what the computer chose
            println("Computer Chose: ", uppercasefirst(myconvert(string(aimove[1]))));

            #determines the winner of the round
            roundwinner = winner(string(playermove[1]), string(aimove[1]));

            #prints out the winner of the round
            println("Round ", string(roundscount), "'s Winner: ",roundwinner);

            #adds 1 to the playercount if the player won
            (roundwinner == "Player") ? (playercount = playercount + 1) : "";

            #adds 1 to the aicount if the computer won
            (roundwinner == "Computer") ? (aicount = aicount + 1) : "";
        )

        #prints out the current statistics so far after each round
        println("Current Statistics:")
        println("   Player Score: ", string(playercount))
        println("   Computer Score: ", string(aicount))
    end
    println()

    #prints out the final result after the whole game is over
    print("RESULT: ")
    (playercount == aicount) ? (println("TIE")) : ""
    (playercount > aicount) ? (println("PLAYER WINS!!!")) : (println("COMPUTER WINS!!!"))

    #asks whether they would like to play again
    playagain()
end

function playagain()

    playagain = lowercase(string(input("Play again? Y/N\n")))

    #ternary operator deciding on whether to continue the game or quit
    (playagain == "y") ? main() : exit()
end

#runs the main function for the first gameplay
main()

max_stars_repo_name: fungktt/JuliaGames

max_issues_repo_name: fungktt/JuliaGames

max_forks_repo_name: fungktt/JuliaGames


text:
@enum_without_prefix CUjit_option CU_
@enum_without_prefix CUjitInputType CU_

function convert_bits(::Type{T}, data::UInt) where T
    if sizeof(data) == sizeof(T)
        return reinterpret(T, data)
    elseif sizeof(data) == 8 && sizeof(T) == 4
        relevant = UInt32(data & 0x00000000ffffffff)
        return reinterpret(T, relevant)
    else
        error("don't know how to extract $(sizeof(T)) bytes out of $(sizeof(data))")
    end
end

function encode(options::Dict{CUjit_option,Any})
    keys = Vector{CUjit_option}()
    vals = Vector{Ptr{Cvoid}}()

    for (opt, val) in options
        push!(keys, opt)
        if opt == JIT_GENERATE_LINE_INFO ||
           opt == JIT_GENERATE_DEBUG_INFO ||
           opt == JIT_LOG_VERBOSE
            push!(vals, convert(Ptr{Cvoid}, convert(Int, val::Bool)))
        elseif opt == JIT_OPTIMIZATION_LEVEL
            push!(vals, convert(Ptr{Cvoid}, convert(Int, val)))
        elseif opt == JIT_INFO_LOG_BUFFER
            buf = val::Vector{UInt8}
            push!(vals, pointer(buf))
            push!(keys, JIT_INFO_LOG_BUFFER_SIZE_BYTES)
            push!(vals, sizeof(buf))
        elseif opt == JIT_ERROR_LOG_BUFFER
            buf = val::Vector{UInt8}
            push!(vals, pointer(buf))
            push!(keys, JIT_ERROR_LOG_BUFFER_SIZE_BYTES)
            push!(vals, sizeof(buf))
        elseif isa(val, Ptr{Cvoid})
            push!(vals, val)
        else
            error("cannot handle option $opt")
        end
    end

    @assert length(keys) == length(vals)
    return keys, vals
end

function decode(keys::Vector{CUjit_option}, vals::Vector{Ptr{Cvoid}})
    @assert length(keys) == length(vals)
    options = Dict{CUjit_option,Any}()

    # decode the raw option value bits to their proper type
    for (opt, val) = zip(keys, vals)
        data = reinterpret(UInt, val)
        if opt == JIT_WALL_TIME
            options[opt] = convert_bits(Cfloat, data)
        elseif opt == JIT_GENERATE_LINE_INFO ||
           opt == JIT_GENERATE_DEBUG_INFO ||
           opt == JIT_LOG_VERBOSE
            options[opt] = Bool(convert_bits(Cint, data))
        elseif opt == JIT_INFO_LOG_BUFFER_SIZE_BYTES ||
           opt == JIT_ERROR_LOG_BUFFER_SIZE_BYTES
            options[opt] = convert_bits(Cuint, data)
        elseif opt == JIT_INFO_LOG_BUFFER ||
           opt == JIT_ERROR_LOG_BUFFER
            options[opt] = reinterpret(Ptr{UInt8}, data)
        else
            options[opt] = data
        end
    end

    # convert some values to easier-to-handle types
    if haskey(options, JIT_INFO_LOG_BUFFER)
        buf = options[JIT_INFO_LOG_BUFFER]
        size = options[JIT_INFO_LOG_BUFFER_SIZE_BYTES]
        delete!(options, JIT_INFO_LOG_BUFFER_SIZE_BYTES)
        options[JIT_INFO_LOG_BUFFER] = unsafe_string(buf, Int(size / sizeof(UInt8)))
    end
    if haskey(options, JIT_ERROR_LOG_BUFFER)
        buf = options[JIT_ERROR_LOG_BUFFER]
        size = options[JIT_ERROR_LOG_BUFFER_SIZE_BYTES]
        delete!(options, JIT_ERROR_LOG_BUFFER_SIZE_BYTES)
        options[JIT_ERROR_LOG_BUFFER] = unsafe_string(buf, Int(size / sizeof(UInt8)))
    end

    return options
end


text:
module DaggerGPU

using Dagger, Requires, Adapt
using Distributed
using KernelAbstractions

macro gpuproc(PROC, T)
    quote
        # Assume that we can run anything
        Dagger.iscompatible_func(proc::$PROC, opts, f) = true
        Dagger.iscompatible_arg(proc::$PROC, opts, x) = true

        # CPUs shouldn't process our array type
        Dagger.iscompatible_arg(proc::Dagger.ThreadProc, opts, x::$T) = false

        # Adapt to/from the appropriate type
        Dagger.move(from_proc::OSProc, to_proc::$PROC, x) = adapt($T, x)
        Dagger.move(from_proc::$PROC, to_proc::OSProc, x) = adapt(Array, x)
    end
end

processor(kind::Symbol) = processor(Val(kind))
processor(::Val) = Dagger.ThreadProc
cancompute(kind::Symbol) = cancompute(Val(kind))
cancompute(::Val) = false

kernel_backend() = kernel_backend(Dagger.Sch.thunk_processor())
kernel_backend(::Dagger.ThreadProc) = CPU()

function __init__()
    @require CUDA="052768ef-5323-5732-b1bb-66c8b64840ba" begin
        include("cu.jl")
    end
    @require ROCArrays="ddf941ca-5d6a-11e9-36cc-a3fed13dd2fc" begin
        include("roc.jl")
    end
end

end


text:
using MuZero
import Flux

#struct TicTacToe <: AbstractEnv
#     history::Vector{<:Integer}
#     rewards::Vector{<:Integer}
#     states::Vector{TicTacToeState}
#     child_visits::Vector
#     root_values::Vector
# end

# @todo figure how to share these constants between game.jl and state.jl
const PLAYER_1 = 1
const PLAYER_2 = 2

struct TicTacToe <: AbstractBoardGame
    history::Vector{<:Integer}
    states::Vector{TicTacToeState}
    rewards::Vector{<:Integer}
    child_visits::Vector
    root_values::Vector
    nodes::Vector{Node}
    backprops::Vector
    function TicTacToe(
        history::Vector{<:Integer},
        states::Vector{TicTacToeState},
        rewards::Vector{<:Integer},
        child_visits=[],
        root_values=[],
        nodes=[],
        backprops=[]
        )
        new(
            history,
            states,
            rewards,
            child_visits,
            root_values,
            nodes,
            backprops,
        )
    end
end

function TicTacToe(state_bits::Integer=0)
    TicTacToe(
        Vector{Int8}(),
        [TicTacToeState(state_bits)],
        Vector{Int8}(),
    )
end

TTT(args...) = TicTacToe(args...)
TicTacToe(state::TicTacToeState) = TicTacToe(key(state))

function TicTacToe(actions::Vector{<:Integer})
    g = TicTacToe()
    for a in actions
        apply!(g, a)
    end
    g
end

MuZero.current_image(game::TicTacToe) = make_image(game, state(game))
MuZero.state(game::TicTacToe) = game.states[end]
MuZero.is_terminal(game::TicTacToe) = game |> state |> is_terminal
MuZero.next_to_play(game::TicTacToe) = just_played(game)
MuZero.legal_actions(game::TicTacToe) = legal_actions(state(game))
MuZero.action_space_size(game::TicTacToe) = 9
MuZero.islegal(game::TicTacToe, action::Integer) = islegal(state(game), action)
MuZero.get_state(game::TicTacToe, state_idx::Integer) = game.states[state_idx]
MuZero.discount(game::AbstractEnv) = Float32(1)

function MuZero.apply!(game::TicTacToe, action::Integer)
    if !MuZero.is_terminal(game) && MuZero.islegal(game, action)
        reward = step!(game, action)
        push!(game.rewards, reward)
        push!(game.history, action)
        return reward
    end
    nothing
end

function MuZero.make_image(game::TicTacToe, image_idx::Integer)
    if 0 < image_idx <= length(game.states)
        return make_image(game, get_state(game, image_idx))
    end
    @assert false "image_idx is $image_idx, but the number of states is $(length(game.states))"
end

function MuZero.make_image_label(game::TicTacToe, image_idx::Integer)
    if 0 < image_idx <= length(game.states)
        actions = game.history[1:(image_idx-1)]
        return join(actions, ",")
        players = player_history(game)[1:(image_idx-1)]
        return join(map(pair->join(pair), zip(players, actions)), ",")
    end
    @assert false "image_idx is $image_idx, but the number of states is $(length(game.states))"
end

function MuZero.reward_matrix(game::TicTacToe)
    rewards_with_initial_zero = [0., game.rewards...]
    m = zeros(Float32, length(rewards_with_initial_zero), 2)
    for (idx,reward) in enumerate(rewards_with_initial_zero)
        m[idx, 1] = iseven(idx) ? reward : -1*reward
        m[idx, 2] = iseven(idx) ? -1*reward : reward
    end
    m
end

function MuZero.policy_matrix(game::TicTacToe)
    # if length(game.child_visits) < length(game.history)
    #     return map(x->uniform_policy(game), 1:(length(game.history)+1))
    # end
    [game.child_visits..., uniform_policy(game)]
end

function MuZero.value_matrix(game::TicTacToe, from::Integer=1; td_steps)
    to = min(length(game.history) + 1, from + td_steps)
    value_matrix(game, from, to)
end

function MuZero.value_matrix(game::TicTacToe, from::Integer=1, to=(length(game.history)+1))
    reward_from = from + 1
    reward_to = to
    rewards = reward_matrix(game)[reward_from:reward_to, :]
    function cumvalue(next_value, row)
        i, reward = row
        reward * discount(game)^i + next_value
    end
    mat = accumulate(
        cumvalue,
        enumerate(eachrow(reverse(rewards; dims=1)));
        init=[0., 0.]) |> reverse
    [mat..., [0.,0.]]
end

MuZero.to_play_at(game::TicTacToe, idx::Integer) = isodd(idx) ? PLAYER_1 : PLAYER_2

function step!(game::TicTacToe, action::Integer)
    @assert !is_terminal(state(game))
    @assert islegal(state(game), action)
    push!(game.states, TicTacToeState(next_state_bits(state(game), action)))
    last_reward(state(game))
end

function Base.copy(game::TicTacToe)
    TicTacToe(
        copy(game.rewards),
        copy(game.history),
        copy.(game.states)
    )
end

to_play(game::TicTacToe) = to_play(state(game))
just_played(arg) = to_play(arg) == 1 ? 2 : 1
calc_cell_states(game::TicTacToe) = calc_cell_states(state(game))
MuZero.next_state(::Type{TicTacToe}, state_key, action) = next_state(TicTacToeState(state_key), action)

Base.show(io::IO, game::TicTacToe) = print_game(game)

function print_game(game::TicTacToe)
    print_board(game)
    player = to_play(game)
    println("\nTo Play: Player $player\n")
    println("History:\n$(join(game.history, ", "))\n")
    println("Rewards:\n$(join(game.rewards, ", "))\n")
    # If search stats created populated, then show.
    if !isempty(game.child_visits)
        println("~~~~~~~~~~~~ Search statistics ~~~~~~~~~~~~\n")
        state_root_values = zip(game.states, game.root_values)
        root_values_str = join(["$(hex_str(state)) => $(root_value)" for (state, root_value) in state_root_values], "\n")
        println("root_values:")
        println("$root_values_str\n")

        state_child_visits = zip(game.states, game.child_visits)
        child_visits_str = join(["$(hex_str(state)) => [$(join(child_visits, ", "))]" for (state, child_visits) in state_child_visits], "\n")

        println("child_visits:")
        println("$child_visits_str\n")

        pos = 1
        println("~~~~~~~~~~~~ Random training sample (from position $pos) ~~~~~~~~~~~~\n")
        sample = TrainingSample(game, pos, 9, 9)
        show(stdout, sample)
    end
    nothing
end

cartesian_coords_to_action(x,y) = (y-1)*3 + x

function print_board(game::TicTacToe)
    function cell_fn(x,y)
        action = cartesian_coords_to_action(x,y)
        cell_player = calc_cell_states(game)[action]
        if cell_player == 1
            return "X"
        elseif cell_player == 2
            return "O"
        end
        nothing
    end
    print_world(Grid(3,3), cell_fn)
end

#= ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ =#

key_to_image(::Type{TicTacToe}, key::UInt32) = TTT(key) |> current_image

function player_history(game::TicTacToe)
    players = ["A","B"]
    map(idx->players[(idx-1)%2+1], 1:length(game.history))
end

function action_player_label()
end

# @todo TicTacToe type is redundant
make_image(::Type{TicTacToe}, state::TicTacToeState) = make_image(state)
make_image(game::TicTacToe, state::TicTacToeState) = make_image(TicTacToe, state)

function make_tictactoe_config()
    make_board_game_config(
        TicTacToe,
        network_type=Incremental,
        weights_file_name="weights-latest.bson",
        max_moves=9,
        root_dirichlet_alpha=0.3,
        action_space_size=9, # Number of all possible actions in the game
        support_size = 10, # Value and reward are scaled (with almost sqrt) and encoded on a vector with a range of -support_size to support_size. Choose it so that support_size <= sqrt(max(abs(discounted reward)))
        observation_shape = (3, 3, 3),
        stacked_observations = 0, # Number of previous observations and previous actions to add to the current observation
        hidden_state_size = 32, # The 1-dimensional size of the encoded game state
        checkpoint_interval = 100, #  Number of training steps before using the model for self-playing
        num_actors = 3,
        weight_decay = 1e-3,
        lr_init=0.002f0,
        save_interval = 10000,
    )
end

mutable struct GameStateIterator
    game_state::TicTacToeState
    by::UInt32
    to::UInt32
end

function all_states(::Type{TicTacToe})
    by = UInt32(1)
    to = UInt32(0b111111111111111111)
    GameStateIterator(TicTacToeState(), by, to)
end

function Base.iterate(i::GameStateIterator, state=UInt32(0))
    i.game_state = TicTacToeState(state)
    if key(i.game_state) > i.to
        return nothing
    end
    i.game_state, state + i.by
end

function MuZero.getstates(t::Type{TicTacToe})
    valid_nonterminal(state) = valid_idx(state) && !is_terminal(state)
    Iterators.filter(valid_nonterminal, all_states(t))
end

function MuZero.getstatesplus(t::Type{TicTacToe})
    valid(state) = valid_idx(state)
    Iterators.filter(valid, all_states(t))
end

function random_game(::Type{TicTacToe})
    g = TTT()
    while !is_terminal(g)
        apply!(g, rand(legal_actions(g)))
    end
    g
end

function check(n::Network, state_key::Vector)
    initial_inference(n, TTT.(state_key)) |> NetworkOutput
end
function check(n::Network, state_key)
    initial_inference(n, TTT(state_key)) |> SingleNetworkOutput
end

check_reward(n, state_key) = check(n, state_key) |> reward
check_value(n, state_key) = check(n, state_key) |> value
check_policy(n::Network, state::AbstractState) = check_policy(n, key(state))
check_policy(n::Network, state_key) = check(n, state_key) |> policy
quick_play() = play_game(config, storage.play_network)

function check_reward(n, state_key, action)
    hid = initial_inference(n, TTT(state_key)) |> SingleNetworkOutput |> hidden_state
    recurrent_inference(n, hid, action) |> SingleNetworkOutput |> reward
end

function check_value(n, state_key, action)
    hid = initial_inference(n, TTT(state_key)) |> SingleNetworkOutput |> hidden_state
    recurrent_inference(n, hid, [action]) |> SingleNetworkOutput |> value
end

function check_policy(n::Network, state_key, action)
    hid = initial_inference(n, TTT(state_key)) |> SingleNetworkOutput |> hidden_state
    recurrent_inference(n, hid, [action]) |> SingleNetworkOutput |> policy
end

function show_policy(state_key)
    show(TTT(state_key))
    check_policy(state_key)
end

function reward_perf()
    pred = check_reward(key.(getstatesplus(TicTacToe)))
    actual = reshape(last_reward.(getstatesplus(TicTacToe)), 1, :)
    Flux.mse(pred, actual)
end

value_fn = value_iteration(TicTacToe)

function value_perf()
    state_keys = key.(getstatesplus(TicTacToe))
    pred = check_value(n, state_keys)
    actual = reshape([value_fn[state_key] for state_key in state_keys], 1, :)
    @assert size(pred) == size(actual)
    Flux.mse(pred, actual)
end

function policy_perf(n::Network; print_failures=false)
    success = Counter(0)
    total = Counter(0)
    pred = check_policy(n, key.(states(TicTacToe)))

    for (idx,state) in enumerate(states(TicTacToe))
        pred_policy = view(pred, :, idx)
        # if rand() < 0.002
        #     a = pred_policy
        #     b =  view(check_policy([state]), :, 1)
        #     c = check_policy(state)
        #     println.([a,b,c])
        #     println(last.(findmax.([a,b,c])))
        #     println()
        # end
        increment!(total)
        legal = legal_actions(state)
        next_vals_inverted = map(action->begin
            state_next = next_state(state, action)
            if in(action, legal)
                return last_reward(state_next) + -1*value_fn[key(state_next)]
            end
            -Inf
        end, 1:9)
        highest_inverted_val = maximum(next_vals_inverted)
        action_choice = (pred_policy |> findmax)[2]

        if in(action_choice, legal)
            if next_vals_inverted[action_choice] == highest_inverted_val
                increment!(success)
            else
                if length(legal) == 1
                    increment!(success)
                else
                    if print_failures
                        println("action $action_choice not highest for $(key(state))")
                    end
                end
            end
        else
            if print_failures
                println("illegal action ($action_choice) for $(key(state))")
            end
        end

    end
    success.n, total.n, success.n/total.n
end

perfect_policy(TicTacToe, state_key) = perfect_policy(TicTacToeState(state_key))

function diagnose_state(state_key)
    state = TicTacToeState(state_key)
    println()
    check_policy(state_key) |> println
    perfect_policy(TicTacToe, state_key) |> println
    print_board(state)
    println()
    println("\nTo Play: Player $(to_play(state))\n")
    for a in sort(collect(legal_actions(state)))
        state_next = next_state(state, a)
        r_pred = check_reward(state_next)
        vn_pred1 = check_value(n, key(state_next))
        vn_pred2 = check_value(n, key(state), a)
        vn_actual = value_fn[key(state_next)]
        println("$a -> $(hex_str(state_next)) (r_pred / vn_pred1 / vn_pred2 / vn_actual: $r_pred / $vn_pred1 / $vn_pred2 / $vn_actual)")
    end
end

function print_board_state(game::TicTacToe)
    println()
    print_board(game)
    println("0x$(key(state(game)))")
    println()
end

print_board(game::TicTacToe) = game |> state |> print_board

reflect(x) = [0 0 1; 0 1 0; 1 0 0] * x

function MuZero.action_image(::Type{TicTacToe}, a::Integer)
    channel = reshape(Flux.onehot(a, 1:9), 3, 3) |> transpose |> reflect
    convert(Array{Float32}, Flux.unsqueeze(channel, 3))
end


text:
using GLM
glm(y ~ 1, DataFrame(y = float64(randbool(10))), Binomial())

n = 2_500_000; srand(1234321)
df2 = DataFrame(x1 = rand(Normal(), n),
                x2 = rand(Exponential(), n),
                ss = pool(rand(DiscreteUniform(50), n)));
beta = unshift!(rand(Normal(),52), 0.5); # "true" parameter values

## Create linear predictor and mean response
df2["eta"] = eta = ModelMatrix(ModelFrame(x1 ~ x1 + x2 + ss, df2)).m * beta;
df2["mu"] = mu = linkinv!(LogitLink(), similar(eta), eta);

df2["y"] = float64(rand(n) .< mu);        # simulate observed responses
head(df2)

gm6 = glm(y ~ x1 + x2 + ss, df2, Binomial())
@time glm(y ~ x1 + x2 + ss, df2, Binomial());


text:
module AlternateFakePackageTests # module name doesn't match file name
using ReTest

@testset "AlternateFakePackageTests" begin
    @test true
end

end # AlternateFakePackageTests

max_issues_repo_name: rfourquet/InlineTest

max_forks_repo_name: rfourquet/InlineTest


text:
module Test01

import Automa
import Automa.RegExp: @re_str
using Test

@testset "Test01" begin
    re = re""
    re.actions[:enter] = [:enter]
    re.actions[:exit] = [:exit]
    machine = Automa.compile(re)
    @test occursin(r"^Automa.Machine\(<.*>\)$", repr(machine))

    for generator in (:table, :inline, :goto), checkbounds in (true, false), clean in (true, false)
        ctx = Automa.CodeGenContext(generator=generator, checkbounds=checkbounds, clean=clean)
        init_code = Automa.generate_init_code(ctx, machine)
        exec_code = Automa.generate_exec_code(ctx, machine, :debug)
        validate = @eval function (data)
            logger = Symbol[]
            $(init_code)
            p_end = p_eof = lastindex(data)
            $(exec_code)
            return cs == 0, logger
        end
        @test validate(b"") == (true, [:enter, :exit])
        @test validate(b"a") == (false, Symbol[])
    end
end

end


text:
import Dagger.Sch: SchedulerOptions, ThunkOptions

@everywhere begin
function inc(x)
    x+1
end
function checkwid(x...)
    @assert myid() == 1
    return 1
end
function checktid(x...)
    @assert Threads.threadid() != 1 || Threads.nthreads() == 1
    return 1
end
end

@testset "Scheduler" begin
    @testset "Scheduler options: single worker" begin
        options = SchedulerOptions(;single=1)
        a = delayed(checkwid)(1)
        b = delayed(checkwid)(2)
        c = delayed(checkwid)(a,b)

        @test collect(Context(), c; options=options) == 1
    end
    @testset "Thunk options: single worker" begin
        options = ThunkOptions(;single=1)
        a = delayed(checkwid; options=options)(1)

        @test collect(Context(), a) == 1
    end
    @static if VERSION >= v"1.3.0-DEV.573"
        if Threads.nthreads() == 1
            @warn "Threading tests running in serial"
        end
        @testset "Scheduler options: threads" begin
            options = SchedulerOptions(;proctypes=[Dagger.ThreadProc])
            a = delayed(checktid)(1)
            b = delayed(checktid)(2)
            c = delayed(checktid)(a,b)

            @test collect(Context(), c; options=options) == 1
        end
        @testset "Thunk options: threads" begin
            options = ThunkOptions(;proctypes=[Dagger.ThreadProc])
            a = delayed(checktid; options=options)(1)

            @test collect(Context(), a) == 1
        end
    end

    @everywhere Dagger.add_callback!(proc->FakeProc())
    @testset "Thunk options: proctypes" begin
        @test Dagger.iscompatible_arg(FakeProc(), nothing, 1) == true
        @test Dagger.iscompatible_arg(FakeProc(), nothing, FakeVal(1)) == true
        @test Dagger.iscompatible_arg(FakeProc(), nothing, 1.0) == false
        @test Dagger.default_enabled(Dagger.ThreadProc(1,1)) == true
        @test Dagger.default_enabled(FakeProc()) == false

        opts = Dagger.Sch.ThunkOptions(;proctypes=[Dagger.ThreadProc])
        as = [delayed(identity; options=opts)(i) for i in 1:5]
        opts = Dagger.Sch.ThunkOptions(;proctypes=[FakeProc])
        b = delayed(fakesum; options=opts)(as...)

        @test collect(Context(), b) == 57
    end
    @everywhere (pop!(Dagger.PROCESSOR_CALLBACKS); empty!(Dagger.OSPROC_CACHE))

end


text:
#!/usr/bin/env julia

function NB_matrix(g)
    edgeidmap, m, aristas = mapa(g)
    B = zeros(Float64, 2*aristas, 2*aristas)
    for (e,u) in edgeidmap
        i, j = src(e), dst(e)
        eles = neighbors(g,j)
        k = j
        for l in eles
            B[edgeidmap[Edge(k,l)],u] = (kron_δ(j,k)*(1-kron_δ(i,l)))
        end
    end
    return B, edgeidmap
end

function mapa(g)
    edgeidmap = Dict{Edge, Int}()
    aristas = ne(g)
    m = 0
    if is_directed(g)
        for e in edges(g)
            m += 1
            edgeidmap[Edge(e.src,e.dst)] = m
            #  edgeidmap[Edge(e.dst,e.src)] = m + aristas
        end
    else
        for e in edges(g)
            m += 1
            edgeidmap[Edge(e.src,e.dst)] = m
            edgeidmap[Edge(e.dst,e.src)] = m + aristas
        end
    end
    edgeidmap, m, aristas
end

function flux_matrix(g)
    edgeidmap, m, aristas = mapa(g)
    B = zeros(Float64, 2*aristas, 2*aristas)
    for (e,u) in edgeidmap
        i, j = src(e), dst(e)
        eles = neighbors(g,j)
        k = j
        for l in eles
            B[edgeidmap[Edge(k,l)],u] = (kron_δ(j,k)*(1-kron_δ(i,l)))*(1/(degree(g,j)))
        end
    end
    return B, edgeidmap
end

function reluctant_matrix(g::SimpleWeightedGraph)
    edgeidmap, m, aristas = mapa(g)
    B = zeros(Float64, 2*aristas, 2*aristas)
    for (e,u) in edgeidmap
        i, j = src(e), dst(e)
        eles = neighbors(g,j)
        k = j
        for l in eles
            B[edgeidmap[Edge(k,l)],u] = (kron_δ(j,k)*(1-kron_δ(l,i)))+(kron_δ(j,k)*kron_δ(l,i)*inv_degree(g,j))
        end
    end
    return B, edgeidmap
end

function normalized_reluctant(g::SimpleWeightedGraph)
    edgeidmap, m, aristas = mapa(g)
    B = zeros(Float64, 2*aristas, 2*aristas)
    for (e,u) in edgeidmap
        i, j = src(e), dst(e)
        eles = neighbors(g,j)
        k = j
        for l in eles
            B[edgeidmap[Edge(k,l)],u] = ((kron_δ(j,k)*(1-kron_δ(i,l)))+(kron_δ(j,k)*kron_δ(l,i)*inv_degree(g,j)))*(1/((degree(g,i)-1)+inv_degree(g,j)))
        end
    end
    return B, edgeidmap
end

function ihara_matrix(g::SimpleWeightedDiGraph)
    v = g.weights
    edgeidmap, m, aristas = mapa(g)
    B = zeros(Float64, aristas, aristas)
    for (e,u) in edgeidmap
        i, j = src(e), dst(e)
        eles = outneighbors(g,j)
        k = j
        for l in eles
            if strengthin(g,l)>0
                B[edgeidmap[Edge(k,l)],u] = (kron_δ(j,k)*(1-kron_δ(i,l)))*(2/(1/v[i,j]+1/v[k,l]))*(1/strengthin(g,l))
            else
                B[edgeidmap[Edge(k,l)],u] = 0
            end
        end
    end
    return B, edgeidmap
end

function ihara_reluctant(g::SimpleWeightedDiGraph)
    v = g.weights
    edgeidmap, m, aristas = mapa(g)
    B = zeros(Float64, aristas, aristas)
    for (e,u) in edgeidmap
        i, j = src(e), dst(e)
        eles = outneighbors(g,j)
        k = j
        for l in eles
            if strengthin(g,l)>0
                B[edgeidmap[Edge(k,l)],u] = (kron_δ(j,k)*(1-kron_δ(i,l)))*((v[i,j]+v[k,l])/2)*(1/strengthin(g,l)) + rel_single(g,k)
            else
                B[edgeidmap[Edge(k,l)],u] = 0
            end
        end
    end
    return B, edgeidmap
end


text:
using ControlToolbox
using Polynomials
using Base.Test

################################################################################
#                              Transfer Function
################################################################################

# Basic Operations
# ================

G = tf([1],[1;1])
H = tf([1],[1;2])

G1 = G-2G

@test G1.num.a == [-1]
@test G1.den.a == [1;1]

G1 = minreal(G + H)

@test G1.num.a == [3;2]
@test G1.den.a == [2;3;1]

G1 = G + H
G2 = parallel(G,H)

@test G1.num == G2.num
@test G1.den == G2.den

G1 = G*H
G2 = series(G,H)

@test G1.num == G2.num
@test G1.den == G2.den

G1 = minreal(G/(1+G*H))
G2 = feedback(G,H)

@test G1.num == G2.num
@test G1.den == G2.den

G = zpk([-1;-1],[+1;+1],4)

@test G.num == 4*poly([-1;-1])
@test G.den == 1*poly([+1,+1])


text:
module Complementarity

using Polyhedra
using StaticArrays
using JuMP, ConditionalJuMP, Cbc
using Test
using Pkg

rot2(θ) = SMatrix{2, 2}(cos(θ), -sin(θ), sin(θ), cos(θ))

struct Obstacle{N, T, H <: HRepresentation{N, T}}
    interior::H
    contact_face::HalfSpace{N, T}
end

struct Environment{N, T, H1 <: HRepresentation{N, T}, H2 <: HRepresentation{N, T}}
    obstacles::Vector{Obstacle{N, T, H1}}
    free_regions::Vector{H2}
end

function contact_basis(face::HalfSpace{2})
    θ = atan(μ)
    R = rot2(θ)
    hcat(R * face.a, R' * face.a)
end

contact_basis(obs::Obstacle) = contact_basis(obs.contact_face)

function ConditionalJuMP.Conditional(op::typeof(in), x::AbstractVector, P::HRepresentation)
    ConditionalJuMP.Conditional(&, [@?(P.A[i, :]' * x <= P.b[i]) for i in 1:length(P)]...)
end

# A simple complementarity-based time-stepping rigid body simulation. All
# notation is taken from Stewart & Trinkle "An Implicit Time-Stepping Scheme for
# Rigid Body Dynamics with Coulomb Friction".

const h = 0.05
const μ = 0.5
const n = [0, 1]
const mass = 1.0
const g = [0, -9.81]

struct ContactResult{T, Tf}
    β::Vector{T}
    λ::T
    c_n::T
    contact_force::Tf
end

JuMP.getvalue(c::ContactResult) = ContactResult(getvalue.((c.β, c.λ, c.c_n, c.contact_force))...)

struct LCPUpdate{T, Tf}
    q::Vector{T}
    v::Vector{T}
    contacts::Vector{ContactResult{T, Tf}}
end

JuMP.getvalue(up::LCPUpdate) =
    LCPUpdate(getvalue.((up.q, up.v))..., getvalue.(up.contacts))

function contact_force(qnext, vnext, obstacle::Obstacle, model::Model)
    n = obstacle.contact_face.a
    D = contact_basis(obstacle)
    k = size(D, 2)

    β = @variable(model,   [1:k], lowerbound=0,   basename="β",     upperbound=100)
    λ = @variable(model,          lowerbound=0,   basename="λ",     upperbound=100)
    c_n = @variable(model,        lowerbound=0,   basename="c_n",   upperbound=100)

    @constraints model begin
        λ .+ D' * vnext .>= 0 # (8)
        μ * c_n .- sum(β) >= 0 # (9)
    end

    @disjunction(model, (n' * qnext - obstacle.contact_face.β == 0), (c_n == 0)) # (10)
    Dtv = D' * vnext
    for j in 1:k
        @disjunction(model, ((λ + Dtv[j]) == 0), β[j] == 0) # (11)
    end
    @disjunction(model, (μ * c_n - sum(β) == 0), (λ == 0)) # (12)

    contact_force = c_n * n .+ D * β
    ContactResult(β, λ, c_n, contact_force)
end

function update(q, v, env::Environment, model::Model)
    qnext = @variable(model, [1:length(q)], lowerbound=-10, basename="qnext", upperbound=10)
    vnext = @variable(model, [1:length(v)], lowerbound=-10, basename="vnext", upperbound=10)

    contacts = [contact_force(qnext, vnext, obs, model) for obs in env.obstacles]
    total_force = mass * g + sum([c.contact_force for c in contacts])

    @constraints model begin
        mass * (vnext - v) .== h * total_force # (5)
        qnext - q .== h .* vnext # (6)
    end

    ConditionalJuMP.disjunction!(
        model,
        [@?(qnext ∈ P) for P in env.free_regions]) # (7)

    LCPUpdate(qnext, vnext, contacts)
end

function simulate(q0, v0, env::Environment, N)
    q, v = q0, v0
    results = LCPUpdate{Float64}[]
    for i in 1:N
        m = Model(solver=CbcSolver())
        up = update(q, v, env, m)
        solve(m)
        push!(results, getvalue(up))
        q = results[end].q
        v = results[end].v
    end
    results
end

function optimize(q0, v0, env::Environment, N::Integer)::Vector{LCPUpdate{Float64}}
    q, v = q0, v0
    m = Model(solver=CbcSolver())
    results = []
    for i in 1:N
        up = update(q, v, env, m)
        push!(results, up)
        q = results[end].q
        v = results[end].v
    end
    solve(m)
    getvalue.(results)
end

function JuMP.setvalue(contact::ContactResult{<:JuMP.AbstractJuMPScalar}, seed::ContactResult{<:Number})
    setvalue(contact.β, seed.β)
    setvalue(contact.λ, seed.λ)
    setvalue(contact.c_n, seed.c_n)
    @assert getvalue(contact.contact_force) ≈ seed.contact_force
end

function JuMP.setvalue(up::LCPUpdate{<:JuMP.AbstractJuMPScalar}, seed::LCPUpdate{<:Number})
    setvalue(up.q, seed.q)
    setvalue(up.v, seed.v)
    setvalue.(up.contacts, seed.contacts)
end

function optimize(q0, v0, env::Environment, seed::Vector{<:LCPUpdate})
    q, v = q0, v0
    m = Model(solver=CbcSolver())
    results = []
    for i in 1:N
        up = update(q, v, env, m)
        setvalue(up, seed[i])
        push!(results, up)
        q = results[end].q
        v = results[end].v
    end
    warmstart!(m, true)
    @assert sum(m.colCat .== :Bin) == 0
    solve(m)
    getvalue.(results)
end

env = Environment(
    [
        Obstacle(
            SimpleHRepresentation{2, Float64}([0 1], [0]),
            HalfSpace{2, Float64}([0, 1], 0)
        ),
        Obstacle(
            SimpleHRepresentation{2, Float64}([-1 0], [-0.2]),
            HalfSpace{2, Float64}([-1, 0], -0.2)
        )
    ],
    [
        SimpleHRepresentation{2, Float64}(
            [0 -1;
             1 0],
            [0, 0.2])
    ]
)

q0 = [-0.5, 0.5]
v0 = [3, -1.5]
N = 20
results1 = simulate(q0, v0, env, N)
results2 = optimize(q0, v0, env, N)
@test all([r1.q ≈ r2.q for (r1, r2) in zip(results1, results2)])

results_seeded = optimize(q0, v0, env, results1)
@test all([r1.q ≈ r2.q for (r1, r2) in zip(results1, results_seeded)])

if haskey(Pkg.installed(), "DrakeVisualizer")
    @eval using DrakeVisualizer;
    @eval using CoordinateTransformations
    DrakeVisualizer.any_open_windows() || DrakeVisualizer.new_window()

    vis = Visualizer()[:block]
    setgeometry!(vis, HyperRectangle(Vec(-0.1, -0.1, 0), Vec(0.2, 0.2, 0.2)))

    q = [r.q for r in results1]
    for qi in q
        settransform!(vis, Translation(qi[1], 0, qi[2]))
        sleep(h)
    end
end

end


text:
module GeometryModule
export Geometry, readGeometryXYZ, angstrom2bohr, angstrom2bohr!

using ..BaseModule
using ..AtomModule

type Geometry
	atoms::Array{Atom,1}
end

function readGeometryXYZ(filename::AbstractString; unit::Symbol=:Bohr)
	geo = Geometry([])
	lines = open(readlines,filename)
	for line in lines[3:end]	# line 1 is n. of atoms, line 2 is comment
		columns = split(line)
		elem = Element(columns[1])
		pos = Position(map(parse,columns[2:4])...)	# float is not recommended for use on tuples
		append!(geo.atoms,[Atom(elem,pos)])
	end
	if (unit == :Bohr) angstrom2bohr!(geo) end	# Should be default as Bohr is an atomic unit
	return geo
end

angstrom2bohr(x::Real) = (return x/0.529177210)
angstrom2bohr(pos::Position) = (return Position(angstrom2bohr(pos.x),angstrom2bohr(pos.y),angstrom2bohr(pos.z)))
function angstrom2bohr!(geo::Geometry)
  for (atom in geo.atoms)
    atom.position = angstrom2bohr(atom.position)
  end
end

end #module

max_stars_repo_name: tkelman/QuantumLab

max_issues_repo_name: tkelman/QuantumLab

max_forks_repo_name: tkelman/QuantumLab


text:
"""
Each entry has the following format:
  - H Accession number
  - D Data description
  - R PMID
  - A Author(s)
  - T Title of the article
  - J Journal reference
  - * Comment or missing
  - C Accession numbers of similar entries with the correlation coefficients
      of 0.8 (-0.8) or more (less).
      Notice: The correlation coefficient is calculated with zeros
      filled for missing values.
  - I Amino acid index data in the following order
      Ala    Arg    Asn    Asp    Cys    Gln    Glu    Gly    His    Ile
      Leu    Lys    Met    Phe    Pro    Ser    Thr    Trp    Tyr    Val
"""
struct Metadata
    id::String
    description::String
    reference::Array{String}
    journal::String
    title::String
    authors::String
    comment::String
end

"""
An amino acid index is a set of 20 numerical values representing various
physico-chemical and biochemical properties of amino acids.

"""
struct Index <: AbstractAAIndex
    data::SVector{20}
    correlation::Dict{String, AbstractFloat}
    metadata::Metadata
end

"""
An amino acid mutation and contact potential matrix is __generally__ 20 x 20
numerical values representing similarity of amino acids.

Note that the row/column identities may not represent standard amino acid
labels depending on the source paper. They should be interpreted by the user.

Each entry has the following format:
  - H Accession number
  - D Data description
  - R PMID
  - A Author
  - T Title of the article
  - J Journal reference
  - * Comment or missing
  - M rows = ARNDCQEGHILKMFPSTWYV, cols = ARNDCQEGHILKMFPSTWYV
      AA
      AR RR
      AN RN NN
      AD RD ND DD
      AC RC NC DC CC
      AQ RQ NQ DQ CQ QQ
      AE RE NE DE CE QE EE
      AG RG NG DG CG QG EG GG
      AH RH NH DH CH QH EH GH HH
      AI RI NI DI CI QI EI GI HI II
      AL RL NL DL CL QL EL GL HL IL LL
      AK RK NK DK CK QK EK GK HK IK LK KK
      AM RM NM DM CM QM EM GM HM IM LM KM MM
      AF RF NF DF CF QF EF GF HF IF LF KF MF FF
      AP RP NP DP CP QP EP GP HP IP LP KP MP FP PP
      AS RS NS DS CS QS ES GS HS IS LS KS MS FS PS SS
      AT RT NT DT CT QT ET GT HT IT LT KT MT FT PT ST TT
      AW RW NW DW CW QW EW GW HW IW LW KW MW FW PW SW TW WW
      AY RY NY DY CY QY EY GY HY IY LY KY MY FY PY SY TY WY YY
      AV RV NV DV CV QV EV GV HV IV LV KV MV FV PV SV TV WV YV VV
"""
struct AMatrix <: AbstractAAIndex
    rowids::String
    columnids::String
    data::Union{SHermitianCompact, SMatrix}
    metadata::Metadata
end


text:
using ..JuliaRL.FeatureCreators
import ..build_algorithm_dict
using Flux

"""
    Agent for learning online in the four rooms environment.
"""
mutable struct NNFourRoomsContAgent{O, P<:Resampling.AbstractPolicy, F} <: JuliaRL.AbstractAgent
    μ::P
    ER::ExperienceReplay
    WER::WeightedExperienceReplay
    opt::O
    normalize::F
    gvf::GVF
    training_gap::Int64
    buffer_size::Int64
    batch_size::Int64
    warm_up::Int64
    algo_dict::Dict{String, Resampling.LearningUpdate}
    sample_dict::Dict{String, String}
    value_dict::Dict{String, Flux.Chain}
    counter::Int64
    state_t::Tuple{Array{Float64, 1}, Bool}
    action_t::Int64
    ones::Array{Array{Float32, 1}, 1}
    to_update::Dict{String, Bool}
    max_is::Float32
    function NNFourRoomsContAgent(μ::P, gvf::GVF, opt::O,
                                  training_gap, buffer_size, batch_size,
                                  warm_up, parsed, env_size;
                                  max_is_ratio=1.0,
                                  rng=Random.GLOBAL_RNG,
                                  normalize::F=identity) where {O, P<:Resampling.AbstractPolicy, F}

        # build tile coder....
        algo_dict, sample_dict, value_type_dict = build_algorithm_dict(parsed; max_is=max_is_ratio)
        value_dict = Dict{String, Flux.Chain}()
        to_update = Dict{String, Bool}()
        init_func=(dims...)->Resampling.glorot_uniform(rng, dims...)
        base_network = Chain(
            Dense(2, 10, relu; initW=init_func),
            Dense(10, 10, relu; initW=init_func),
            Dense(10, 1; initW=init_func)
        )
        for key in keys(algo_dict)
            if value_type_dict[key] == "State"
                value_dict[key] = deepcopy(base_network)
            else
                throw("Not Implemented")
            end
            to_update[key] = true
        end
        er_names = [:s_t, :s_tp1, :a_t, :r, :γ_tp1, :terminal, :mu_t, :pi_t, :ρ]
        state_type = Array{Float64, 1}
        er_types = [state_type, state_type, Int64, Array{Float32,1}, Array{Float32,1}, Bool, Float64, Float64, Array{Float32,1}]
        ER = ExperienceReplay(buffer_size, er_types; column_names=er_names)
        WER = WeightedExperienceReplay(buffer_size, er_types; column_names=er_names)
        new{O, P, F}(μ, ER, WER, opt, normalize, gvf, training_gap, buffer_size, batch_size, warm_up,
                     algo_dict, sample_dict, value_dict, warm_up, ([0.0,0.0], false), 0,
                     [ones(Float32, 1) for i in 1:batch_size], to_update,
                     max_is_ratio)
    end
end

JuliaRL.get_action(agent::NNFourRoomsContAgent, s; rng=Random.GLOBAL_RNG) =
    StatsBase.sample(rng, agent.μ, s[1])

function JuliaRL.start!(agent::NNFourRoomsContAgent, env_s_tp1; rng=Random.GLOBAL_RNG, kwargs...)
    agent.counter = agent.warm_up
    agent.state_t = env_s_tp1
    agent.action_t = JuliaRL.get_action(agent, env_s_tp1; rng=rng)
    return agent.action_t
end

function JuliaRL.step!(agent::NNFourRoomsContAgent, env_s_tp1, r, terminal; rng=Random.GLOBAL_RNG, kwargs...)

    # After taking action at timestep t
    μ_prob = get(agent.μ, agent.state_t[1], agent.action_t, nothing, nothing, nothing)
    c, γ, π_prob = get(agent.gvf, agent.state_t, agent.action_t, env_s_tp1, nothing, nothing)

    experience = (agent.normalize(agent.state_t[1]),
                  agent.normalize(env_s_tp1[1]),
                  agent.action_t, [c], [γ], terminal,
                  μ_prob, π_prob, [π_prob/μ_prob])

    add!(agent.ER, experience)
    add!(agent.WER, experience, π_prob/μ_prob)

    if agent.counter == 0
        agent.counter = agent.training_gap
        train_value_functions(agent; rng=rng)
    end
    agent.counter -= 1

    agent.state_t = env_s_tp1
    agent.action_t = JuliaRL.get_action(agent, agent.state_t; rng=rng)
    return agent.action_t
end

function train_value_functions(agent::NNFourRoomsContAgent; rng=Random.GLOBAL_RNG)

    samp_er = sample(agent.ER, agent.batch_size; rng=rng)
    samp_wer = sample(agent.WER, agent.batch_size; rng=rng)

    action_tp1_er = [StatsBase.sample(rng, agent.μ, s) for s in samp_er[:s_tp1]]
    action_tp1_wer = [StatsBase.sample(rng, agent.μ, s) for s in samp_wer[:s_tp1]]
    # action_tp1 = StatsBase.sample(rng, [1, 2], StatsBase.weights(π), batch_size)

    arg_er = (samp_er[:ρ], samp_er[:s_t], samp_er[:s_tp1],
              samp_er[:r], samp_er[:γ_tp1], samp_er[:terminal],
              samp_er[:a_t], action_tp1_er, agent.gvf.policy)

    # println("ER:", samp_er[:s_t])

    arg_wer = (agent.ones, samp_wer[:s_t], samp_wer[:s_tp1],
               samp_wer[:r], samp_wer[:γ_tp1], samp_wer[:terminal],
               samp_er[:a_t], action_tp1_wer, agent.gvf.policy)
    # println("WER:", samp_wer[:s_t])

    # opt = Descent(α)
    for key in filter((k)->agent.to_update[k], keys(agent.algo_dict))
        try
            if agent.sample_dict[key] == "IR"
                corr_term = 1.0
                if key == "BCIR"
                    corr_term = Resampling.total(agent.WER)/size(agent.WER)[1]
                end
                update!(agent.value_dict[key], agent.opt, agent.algo_dict[key], arg_wer...; corr_term=corr_term)
            elseif agent.sample_dict[key] == "ER"
                if key == "WISBuffer"
                    corr_term = Resampling.total(agent.WER)/size(agent.WER)[1]
                    update!(agent.value_dict[key], agent.opt, agent.algo_dict[key], arg_er...; corr_term=1.0/corr_term)
                elseif key == "NormIS"
                    update!(agent.value_dict[key], agent.opt, agent.algo_dict[key], arg_er[1]./agent.max_is, arg_er[2:end]...;)
                else
                    update!(agent.value_dict[key], agent.opt, agent.algo_dict[key], arg_er...;)
                end
            elseif agent.sample_dict[key] == "Optimal"
                samp_opt = agent.ER.buffer
                arg_opt = (samp_opt[:ρ], samp_opt[:s_t], samp_opt[:s_tp1],
                           samp_opt[:r], samp_opt[:γ_tp1], samp_opt[:terminal])
                update!(agent.value_dict[key], opt, agent.algo_dict[key], arg_opt...;)
            end
        catch ex
            println(ex)
            try
                if ex.msg == "Loss is infinite" || ex.msg == "Loss is NaN"
                    agent.to_update[key] = false
                else
                    throw(ex)
                end
            catch
                throw(ex)
            end
        end
    end
end

function predict(agent::NNFourRoomsContAgent, key::AbstractString, eval_states::Array{Array{Float64, 1}, 1})
    # states = [create_features(agent.tilecoder, eval_state./11.0) for eval_state in eval_states]
    return agent.value_dict[key].(states)
end

function predict!(agent::NNFourRoomsContAgent, eval_states::Array{Array{Float64, 1}, 1}, predict_dict::Dict{String, Array{Float64, 1}})
    states = [agent.normalize(eval_state) for eval_state in eval_states]
    if length(keys(predict_dict)) == 0
        for key in keys(agent.value_dict)
            predict_dict[key] = zeros(length(eval_states))
        end
    end
    for key in keys(agent.value_dict)
        predict_dict[key] .= getindex.(Flux.data.(agent.value_dict[key].(states)), 1)
    end
end


text:
#    Copyright (c) 2021 J.A. Duffek
#    Copyright (c) 2000 D.M. Spink
#
#    This program is free software; you can redistribute it and/or modify
#    it under the terms of the GNU General Public License as published by
#    the Free Software Foundation; either version 2 of the License, or
#    (at your option) any later version.
#
#    This program is distributed in the hope that it will be useful,
#    but WITHOUT ANY WARRANTY; without even the implied warranty of
#    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
#    GNU General Public License for more details.
#
#    You should have received a copy of the GNU General Public License
#    along with this program; if not, see <http://www.gnu.org/licenses/>.

"""
    nrbline()::NURBS1D

or

    nrbline(p1::Vector{F},p2::Vector{F})::NURBS1D where {F<:AbstractFloat}

Construct a straight line.

Constructs NURBS data structure for a straight line. If coordinates are
included the function returns a unit straight line along the x-axis.

# Arguments:
- `p1`: 2D or 3D cartesian coordinate of the start point.
- `p2`: 2D or 3D cartesian coordinate of the end point.

# Output:
- `crv`: NURBS curve for a straight line.

# Examples:
```julia
julia> crv = nrbline()
julia> crv = nrbline(p1,p2)
```
"""
function nrbline(p1::Vector{F},p2::Vector{F})::NURBS1D where {F<:AbstractFloat}
if length(p1) > 3 || length(p2) > 3 || length(p1) != length(p2)
  throw(ArgumentError("Only vectors up to 3D with equal length allowed!"));
end # if
coefs = zeros(F,(4,2));
coefs[4,:] .= 1;
coefs[1:length(p1),1] = p1[:];
coefs[1:length(p2),2] = p2[:];
return nrbmak(coefs, vec([0.0 0.0 1.0 1.0]));
end # nrbline

# overloading of function nrbline with zero arguments
function nrbline()
coefs = zeros(typeof(1.0),(4,2));
coefs[4,:] .= 1;
coefs[1,2] = 1;
return nrbmak(coefs, vec([0.0 0.0 1.0 1.0]));
end # nrbline


text:
"""
** PVT **

- `n_trials`: number of trials
- `trial`: current trial
- `lb`: ISI lower bound
- `ub`: ISI upper bound
- `width`: screen width
- `height`: screen height
- `scheduler`: event scheduler
- `screen`: visual objects on screen
- `canvas`: GTK canvas
- `window`: GTK window
- `visible`: GUI visible
- `speed`: real time speed

Function Signature

````julia
PVT(;n_trials=10, trial=1, lb=2.0, ub=10.0, width=600.0, height=600.0, scheduler=nothing,
    screen=Vector{VisualObject}(), window=nothing, canvas=nothing, visible=false, speed=1.0)
````
"""
mutable struct PVT{T,W,C} <: AbstractTask
    n_trials::Int
    trial::Int
    lb::Float64
    ub::Float64
    width::Float64
    hight::Float64
    scheduler::T
    screen::Vector{VisualObject}
    canvas::C
    window::W
    visible::Bool
    realtime::Bool
    speed::Float64
end

function PVT(;
    n_trials=10,
    trial=1,
    lb=2.0,
    ub=10.0,
    width=600.0,
    height=600.0,
    scheduler=nothing,
    screen=Vector{VisualObject}(),
    window=nothing,
    canvas=nothing,
    visible=false,
    realtime=false,
    speed=1.0,
    )
    visible ? ((canvas,window) = setup_window(width)) : nothing
    visible ? Gtk.showall(window) : nothing
    return PVT(n_trials, trial, lb, ub, width, height, scheduler, screen, canvas, window, visible,
        realtime, speed)
end

function start!(task::PVT, models)
    run_trial!(task, models)
end

function sample_isi(task)
    return rand(Uniform(task.lb, task.ub))
end

function present_stimulus(task, models)
    vo = VisualObject()
    add_to_visicon!(models, vo; stuff=true)
    push!(task.screen, vo)
    w = task.width / 2
    task.visible ? draw_object!(task, "X", w, w) : nothing
    return nothing
end

function run_trial!(task, models)
    isi = sample_isi(task)
    description = "present stimulus"
    register!(task, present_stimulus, after, isi, task, models;
        description)
end

function press_key!(task::PVT, actr, key)
    if key == "sb"
        empty!(task.screen)
        task.visible ? clear!(task) : nothing
        if task.trial < task.n_trials
            task.trial += 1
            run_trial!(task, actr)
        else
            stop!(task.scheduler)
        end
    end
    return nothing
end

function repaint!(task::PVT, actr)
    #clear!(task)
    #draw_attention!(task, actr)
end


text:
#--------------------------------------------------------------------
# Spacetime Discretization methods in Julia
# Soham 03-2018
#--------------------------------------------------------------------

function setcolormap(vec::Array{Float64,1}, map::String, samples::Int)
    nvec = (vec .- minimum(vec))./(maximum(vec) .- minimum(vec))
    clrs = colormap(map, samples+1)
    return clrs[round.(Int, (nvec*samples)+1)]
end

function drawmultipatch(dbase::Dict{Array{Int,1}, Patch}, filename)
    M   = Int(sqrt(length(dbase)))
    AP  = Float64[]

    #-----------------------------------------------
    # Create color pallete
    #-----------------------------------------------
    for mx in 1:M, my in 1:M
        AP = vcat(AP, vec(dbase[[mx, my]].value))
    end
    globalcmap = setcolormap(AP, "Blues", length(AP))

    #-----------------------------------------------
    # Patch sizing
    #-----------------------------------------------
    samplepatch = dbase[[1,1]].value
    (Nx, Ny) = size(samplepatch) .- 1
    (wx, wy) = (chebweights(Nx), chebweights(Ny))
    (lx, ly) = (sort(2.0 - cumsum(wx))*200, sort(2.0 - cumsum(wy))*200)

    #-----------------------------------------------
    # set-up canvas
    #-----------------------------------------------
    canvas   = Drawing(800, 800, "$filename.pdf")
    push!(lx, 400.0)
    push!(ly, 400.0)
    origin(400, 650)
    rotate(-3pi/4)

    setline(0.4)
    for i in 1:Nx+1, j in 1:Ny+1
        rect(lx[i], ly[j], lx[i+1] - lx[i], ly[j+1] - ly[j], :stroke)
    end

    #-----------------------------------------------
    # set ticks
    #-----------------------------------------------
    sethue("black")
    setline(0.4)
    (lx0, ly0) = (lx[1], ly[1])
    (lxe, lye) = (lx[end], ly[end])
    (lxm, lym) = (lx0, ly0)./2 .+ (lxe, lye)./2

    # corners
    line(Point(lx0, ly0), Point(lx0, ly0-10), :stroke)
    line(Point(lx0, ly0), Point(lx0-10, ly0), :stroke)
    line(Point(lxe, lye), Point(lxe, lye+10), :stroke)
    line(Point(lxe, lye), Point(lxe+10, lye), :stroke)
    line(Point(lxe, ly0), Point(lxe, ly0-10), :stroke)
    line(Point(lxe, ly0), Point(lxe+10, ly0), :stroke)
    line(Point(lx0, lye), Point(lx0-10, lye), :stroke)
    line(Point(lx0, lye), Point(lx0, lye+10), :stroke)
    line(Point(lx0, ly0), Point(lx0, ly0-10), :stroke)
    line(Point(lx0, ly0), Point(lx0-10, ly0), :stroke)

    # mid points
    line(Point(lxm, ly0), Point(lxm, ly0-10), :stroke)
    line(Point(lx0, lym), Point(lx0-10, lym), :stroke)
    line(Point(lxm, lye), Point(lxm, lye+10), :stroke)
    line(Point(lxe, lym), Point(lxe+10, lym), :stroke)

    # edges
    line(Point(lx0, ly0), Point(lx0, lye), :stroke)
    line(Point(lx0, ly0), Point(lxe, ly0), :stroke)
    line(Point(lxe, ly0), Point(lxe, lye), :stroke)
    line(Point(lx0, lye), Point(lxe, lye), :stroke)

    #-----------------------------------------------
    # tick labels
    #-----------------------------------------------
    settext("(1, 1)", Point(lx0 - 18, ly0 - 18);
            halign = "center",
            valign = "top")
    settext("(-1, -1)", Point(lxe + 18, lye + 18);
            halign = "center",
            valign = "bottom")
    settext("(-1, 1)", Point(lx0 - 20, lye + 20);
            halign = "center",
            valign = "center")
    settext("(1, -1)", Point(lxe + 20, ly0 - 20);
            halign = "center",
            valign = "center")

    #-----------------------------------------------
    # colorbar
    #-----------------------------------------------

    gsave()
    rotate(3pi/4)
    x   = collect(linspace(-200, 200, length(AP)))
    clr = colormap("Blues", length(AP))

    for i in 1:length(AP)-1
        sethue((clr[i].r, clr[i].g, clr[i].b))
        rect(x[i], 92.5, x[i+1] - x[i], 15, :fill)
    end

    #-----------------------------------------------
    # colobar ticks
    #-----------------------------------------------
    sethue("black")
    setline(0.4)
    line(Point(-200, 107.5), Point(-200, 85), :stroke)
    line(Point(0, 107.5), Point(0, 85), :stroke)
    line(Point(200, 107.5), Point(200, 85), :stroke)

    min = findmin(AP)[1]
    max = findmax(AP)[1]
    med = (min + max)/2

    settext("$min", Point(-200, 80);
                halign = "center",
                valign = "center")
    settext("$med", Point(0, 80);
                halign = "center",
                valign = "center")
    settext("$max", Point(200, 80);
                halign = "center",
                valign = "center")

    grestore()

    #-----------------------------------------------
    # draw patches
    #-----------------------------------------------

    for mx in 1:M, my in 1:M
        gsave()
        Luxor.translate(((my-1)/2M)*(maximum(lx)+400), ((mx-1)/2M)*(maximum(lx)+400))
        scale(1/M)
        patch     = dbase[[mx, my]].value
        localcmap = globalcmap[(sub2ind((M, M), mx, my)-1)*length(patch) + 1:(sub2ind((M, M), mx, my))*length(patch)]
        cmap = reshape(localcmap, size(patch))

        for i in 1:Nx+1, j in 1:Ny+1
           sethue((cmap[i,j].r, cmap[i,j].g, cmap[i,j].b))
           rect(lx[i], ly[j], lx[i+1] - lx[i], ly[j+1] - ly[j], :fill)
        end

        grestore()

    end

    #-----------------------------------------------
    #arrow
    #-----------------------------------------------
    sethue("black")
    setline(0.4)
    origin(700, 100)
    arrow(O, O .+ Point(40, -40))
    settext("v", O .+ Point(40, -40);
                halign = "top",
                valign = "right")
    rotate(-pi/2)
    arrow(O, O .+ Point(40, -40))
    settext("u", O .+ Point(40, -40);
                halign = "bottom",
                valign = "right")

    finish()
    return canvas
end


text:
function std_checkbox()
    map(x->map(RGBA{N0f8}, x),
        (loadasset("checked.png"), loadasset("unchecked.png"))
    )
end

function widget(tick::Signal{T}, window;
        checkbox = std_checkbox(), area=(50, 50),
        kw_args...
    ) where T<:Union{Bool, UInt8}
    icon = map(tick) do t
        checkbox[Bool(t) ? 1 : 2]
    end
    robj = visualize(icon; primitive=SimpleRectangle(0,0,area[2],area[2]), kw_args...)
    robj, togl = toggle_button(tick, robj, window)
    robj, togl
end

function choice_visual(x::Union{Function, DataType, Enum}; kw_args...)
    visualize(string(x); kw_args...)
end
function choice_visual(x; kw_args...)
    visualize(x; kw_args...)
end

function widget(enum::Signal{T}, window;
        kw_args...
    ) where T <: Enum
    all_enums = collect(instances(T))
    i0 = findfirst(x-> x==value(enum), all_enums)
    widget(Signal(all_enums), window; kw_args...)
end

function myscale!(robj, target)
    bb = value(boundingbox(robj))
    w1, w2, w3 = widths(bb)
    w = Vec2f0(w1, w2)
    wt1, wt2, wt3 = widths(target)
    wt = Vec2f0(wt1, wt2)
    t = translationmatrix(-minimum(bb))
    t2 = translationmatrix(minimum(target))
    sm = minimum((1f0 ./ w) .* wt)
    s = scalematrix(Vec3f0(sm, sm, 1))
    m = t2*s*t
    set_arg!(robj, :model, m)
    set_arg!(robj, :boundingbox, Signal(m*bb))
end

function widget(choices::Signal{T}, window;
        text_scale = 4mm, start_idx = 1, area = (150, 30), kw_args...
    ) where T <: AbstractVector
    choices_v = value(choices)
    show_area = SimpleRectangle{Float32}(0, 0, area...)
    vizzes = map(enumerate(choices_v)) do i_c
        i, c = i_c
        vis = choice_visual(c;
            relative_scale = text_scale,
            preferred_camera = :fixed_pixel,
            kw_args...
        )
        vs = vis.children[][:visible]
        vs2 = Signal(i == start_idx)

        vis.children[][:visible] = const_lift(vs, vs2) do a, b
            !a ? false : b # only use parent visibility for hiding, not for showing
        end
        push!(vs2, i == start_idx) # doesn't seem to take the correct value otherwise
        myscale!(vis, AABB{Float32}(Vec3f0(5, 5, 0), Vec3f0(area..., 1) - Vec3f0(5,5,0)))
        vis, vs2
    end
    vis = visualize(show_area, color = RGBA{Float32}(0.95, 0.95, 0.95, 0.4)).children[]
    visual = Context(vis, map(first, vizzes)...)
    selected = foldp(start_idx, toggle(visual, window)) do i0, _
        idx = mod1(i0 + 1, length(vizzes))
        for i = 1:length(vizzes)
            push!(vizzes[i][2], i == idx)
        end
        idx
    end
    visual, map(getindex, choices, selected, typ=eltype(choices_v))
end


text:
function read_nth_lines(stream, num)
  for i = 1:num-1
    readline(stream)
  end
  result = readline(stream)
  print(result != "" ? result : "No such line.")
end

max_stars_repo_name: ethansaxenian/RosettaDecode

max_issues_repo_name: ethansaxenian/RosettaDecode

max_forks_repo_name: ethansaxenian/RosettaDecode


text:
using PeaceVote
using PeaceCypher

using PeaceFounder.Braiders: BraiderConfig
using PeaceFounder.BraidChains: BraidChainConfig

#setnamespace(@__MODULE__)
#uuid = PeaceVote.uuid("Community")

# Cleanup from the previous tests
#dirs = [PeaceVote.keydir(uuid), PeaceVote.datadir(uuid), PeaceVote.communitydir(uuid)]

for dir in [homedir() * "/.peacevote/"]
    isdir(dir) && rm(dir,recursive=true)
end

demespec = PeaceVote.DemeSpec("PeaceDeme",:default,:PeaceCypher,:default,:PeaceCypher,:PeaceFounder)
save(demespec) ### Necessary to connect with Mixer

uuid = demespec.uuid
### I could actually rely on DemeSpec file for creating signers
maintainer = PeaceVote.Signer(uuid,"maintainer")

#mixer = PeaceVote.Signer(uuid,"mixer")
#mixerserver = PeaceFounder.Mixer(1999,deme,mixer)

server = PeaceVote.Signer(uuid,"server")

MIXER_ID = server.id ### Self mixing
SERVER_ID = server.id

certifierconfig = nothing
braiderconfig = PeaceFounder.BraiderConfig(2000,2001,3,SERVER_ID,(uuid,MIXER_ID))
braidchainconfig = PeaceFounder.BraidChainConfig(maintainer.id,[(uuid,maintainer.id),],server.id,2002,2003,2004)
systemconfig = PeaceFounder.SystemConfig(2001,2005,certifierconfig,braiderconfig,braidchainconfig)

PeaceFounder.save(systemconfig,maintainer)

### Starting the server

system = PeaceFounder.System(demespec,server)

sleep(2) # for waiting until server is ready

### Initializing without a ledger. I will change that shortly.
deme = Deme(demespec,nothing)

### Theese are our members

for i in 1:9
    account = "account$i"
    keychain = PeaceVote.KeyChain(deme,account)
    identification = PeaceVote.ID("$i","today",keychain.member.id)
    cert = PeaceVote.Certificate(identification,maintainer)
    @show register(deme,cert)
end

# First braiding
sleep(1)

@sync for i in 1:9
    account = "account$i"
    keychain = PeaceVote.KeyChain(deme,account) ### The issue is perhaps
    @async PeaceVote.braid!(keychain)
end

error("STOP") ### Now I should see members and braids.

pmember = PeaceVote.Member(uuid,"account2")
propose("Found peace for a change?",["yes","no","maybe"],pmember);

# Second braiding

@sync for i in 1:9
    account = "account$i"
    keychain = PeaceVote.KeyChain(uuid,account)
    @async PeaceVote.braid!(keychain)
end

# Now someone sends the proposal

pmember = PeaceVote.Member(uuid,"account1")
propose("Let's vote for a real change",["yes","no"],pmember);

sleep(1)

messages = braidchain()
proposals = PeaceVote.proposals(messages)

# Voting

for i in 1:9
    account = "account$i"
    keychain = PeaceVote.KeyChain(uuid,account)

    # Notice that this is after braiding
    ### We need to also update the registrator
    voter = PeaceVote.Voter(keychain,proposals[1],messages)
    option = PeaceVote.Option(proposals[1],rand(1:3))
    vote(option,voter)

    voter = PeaceVote.Voter(keychain,proposals[2],messages)
    option = PeaceVote.Option(proposals[2],rand(1:2))
    vote(option,voter)
end

max_stars_repo_name: akels/DemocracyCredit

max_issues_repo_name: akels/DemocracyCredit

max_forks_repo_name: akels/DemocracyCredit


text:
@testset "Bold x*" begin # issue 223
    h = raw"**x\***" * J.EOS |> seval
    @test h == "
x*
\n"

    h = raw"_x\__" * J.EOS |> seval
    @test h == "
x_
\n"
end

@testset "Bold code" begin # issue 222
    h = raw"""A **`master`** B.""" |> jd2html
    @test h == "
A master B.
\n"
end

@testset "Tickssss" begin # issue 219
    st = raw"""A `B` C""" * J.EOS
    tokens = J.find_tokens(st, J.MD_TOKENS, J.MD_1C_TOKENS)
    @test tokens[1].name == :CODE_SINGLE
    @test tokens[2].name == :CODE_SINGLE

    st = raw"""A ``B`` C""" * J.EOS
    tokens = J.find_tokens(st, J.MD_TOKENS, J.MD_1C_TOKENS)
    @test tokens[1].name == :CODE_DOUBLE
    @test tokens[2].name == :CODE_DOUBLE

    st = raw"""A ``` B ``` C""" * J.EOS
    tokens = J.find_tokens(st, J.MD_TOKENS, J.MD_1C_TOKENS)
    @test tokens[1].name == :CODE_TRIPLE
    @test tokens[2].name == :CODE_TRIPLE

    st = raw"""A ````` B ````` C""" * J.EOS
    tokens = J.find_tokens(st, J.MD_TOKENS, J.MD_1C_TOKENS)
    @test tokens[1].name == :CODE_PENTA
    @test tokens[2].name == :CODE_PENTA

    st = raw"""A ```b B ``` C""" * J.EOS
    tokens = J.find_tokens(st, J.MD_TOKENS, J.MD_1C_TOKENS)
    @test tokens[1].name == :CODE_LANG
    @test tokens[2].name == :CODE_TRIPLE

    st = raw"""A `````b B ````` C""" * J.EOS
    tokens = J.find_tokens(st, J.MD_TOKENS, J.MD_1C_TOKENS)
    @test tokens[1].name == :CODE_LANG2
    @test tokens[2].name == :CODE_PENTA

    h = raw"""
        A
        `````markdown
        B
        `````
        C
        """ |> jd2html

    @test isapproxstr(h, raw"""

A

B

 C

            """)

    h = raw"""
        A
        `````markdown
        ```julia
        B
        ```
        `````
        C
        """ |> jd2html
    @test isapproxstr(h, raw"""

A

```julia
            B
            ```

 C

            """)

end


text:
function simsys(sys,nt,seed)
	###
	# simulated systems
	###

	# random seed for the simulations
	srand(seed)

	if sys==1
		#
		# system 1

		# ref: "Subspace Identification for Linear Systems"
		# P. Van Overschee and B. De Moor, page 81
		#
		A = [0.6 0.6 0.;-0.6 0.6 0.;0. 0. 0.4]
		K = [0.1706 -0.1507 0.2772]'
		C = [0.7831 0.5351 0.9701]
		R = 6.3663
		Q = K * R * K'

		et = sqrt(R)*randn(1,nt)
		w = K*et

		ns = 3
 		x = zeros(ns,nt)
		x[:,1] = w[:,1]
		for t in 1:nt-1
			x[:,t+1] = A * x[:,t] + w[:,t]
		end
		y = C * x + et

		salida = y, A, K, C, Q, R
		#
	elseif sys==2
		###
		# ref: "Subspace Identification for Linear Systems"
		# P. Van Overschee and B. De Moor, page 115
		###
		A = [0.6 0.6 0.;-0.6 0.6 0.;0. 0. 0.7]
		B = [1.6161 -0.3481 2.6319]'
		K = [-1.1472 -1.5204 -3.1993]'
		C = [-0.4373 -0.5046 -0.0936]
		D = -0.7759
		R = 0.0432
		Q = K * R * K'

		et = sqrt(R)*randn(1,nt)
		w = K*et

		# input = AR(1)
		u = zeros(1,nt)
		a = randn(nt)
		for t in 1:nt-1
			u[t+1] = 0.3*u[t] + a[t]
		end

		# states
		ns = size(A,1)
		x = zeros(ns,nt)
		x[:,1] = w[:,1]
		for t in 1:nt-1
			x[:,t+1:t+1] = A * x[:,t:t] + B * u[t] + w[:,t:t]
		end
		# output
		y = C * x + D * u + et

		salida = y, u, A, B, K, C, D, Q, R
	end

	return salida

end


text:
# using Revise
using ClusterManagers
using Distributed

worker_ids = if "SLURM_JOBID" in keys(ENV)
    wids = ClusterManagers.addprocs_slurm(parse(Int, ENV["SLURM_NTASKS"]))
    wids
else
    addprocs(4)
end

Distributed.@everywhere begin
    using Pkg
    Pkg.activate(".")
    Pkg.instantiate()
end

Distributed.@everywhere begin
    using POMDPs
    using BasicPOMCP
    using POMCPOW
    using VOOTreeSearch
    using ParticleFilters
    using CPUTime
    using DataFrames
    using ARDESPOT
    using Distributions
    using VDPTag2
    using ProgressMeter
    using Dates
    using CSV
    using POMDPSimulators
    using Random
    using Printf
    using LinearAlgebra
    using VDPTag2Experiments
    using LinearAlgebra: diagm, PosDefException, det, eigvals, eigvecs
    using Random: MersenneTwister, seed!
    using Distances: euclidean

    using MCTS
    using StaticArrays
    using POMDPModelTools: GenerativeBeliefMDP
    include("./scripts/VDPTag2/vdp_next_action.jl")
end

@info "Started $(Distributed.nworkers()) workers..."
@info "Precompiling simulation code..."

max_depth = 10
solver_list = Dict("POMCPOW"=>[110.0, 30.0, 30.0, 5.0, 100.0],
                    "VOMCPOW"=>[85.0, 30.0, 30.0, 2.5, 100.0, 0.7])

on_cluster = true
save_data = true
parallelize = true

if on_cluster
    max_time_list = 10.0 .^ (-2:0.25:0)
    max_query = 100_000
    n_iters = 1000
else
    # max_time_list = 10.0 .^ (-2:0.5:0)
    max_time_list = [0.1]
    max_query = 100_000
    n_iters = 50
end

# VDPTag2
pomdp = VDPTagPOMDP(mdp=VDPTagMDP(barriers=CardinalBarriers(0.2, 2.8)))

println(stderr, "Max times: ", max_time_list)
println(stderr, "Max depth: ", max_depth)
println(stderr, "Max query: ", max_query)
println(stderr, "Number of iterations: ", n_iters)
println(stderr, "Solvers: ", solver_list)
println(stderr, "POMDP: ", pomdp)

Distributed.@everywhere begin
    struct SimWrap
        x::Vector{Float64}
        solver_name::String
        i::Int64
        max_time::Float64
        max_depth::Int64
        max_query::Int64
        pomdp::POMDP
    end

    function gen_sim(sim_wrap::SimWrap)
        x = sim_wrap.x
        solver_name = sim_wrap.solver_name
        i = sim_wrap.i
        max_time = sim_wrap.max_time
        max_depth = sim_wrap.max_depth
        max_query = sim_wrap.max_query
        pomdp = sim_wrap.pomdp

        c = x[1]
        k_act = x[2]
        inv_alpha_act = x[3]
        k_obs = x[4]
        inv_alpha_obs = x[5]

        rng = MersenneTwister(7)
        rollout_policy = ToNextMLSolver(rng)

        if solver_name == "POMCPOW"
            solver = POMCPOWSolver(tree_queries=max_query,
                                   criterion=MaxUCB(c),
                                   max_depth=max_depth,
                                   max_time=max_time,
                                   k_action=k_act,
                                   alpha_action=1/inv_alpha_act,
                                   k_observation=k_obs,
                                   alpha_observation=1/inv_alpha_obs,
                                   estimate_value=FORollout(rollout_policy),
                                   next_action=RootToNextMLFirst(rng),
                                   default_action=BasicPOMCP.ReportWhenUsed(TagAction(false, 0.0)),
                                   rng=rng
                                  )

            planner = solve(solver, pomdp)
            filter = SIRParticleFilter(deepcopy(pomdp),
                                       100_000,
                                       rng=MersenneTwister(i+90_000))

            seed!(planner, i+40000)
            sim = Sim(deepcopy(pomdp),
                      planner,
                      filter,
                      rng=MersenneTwister(i+50_000),
                      max_steps=100,
                      metadata=Dict(:i=>i, :solver=>solver_name)
                     )

        elseif solver_name == "VOMCPOW"
            p_voo = x[6]
            sig = diagm([0.1])
            epsilon = 0.1/10.0
            voo_opt = VOOActionGenerator(p_voo, Distributions.Uniform(), MvNormal(sig), rng, epsilon, true, 20)

            solver = POMCPOWSolver(tree_queries=max_query,
                                   criterion=MaxUCB(c),
                                   max_depth=max_depth,
                                   max_time=max_time,
                                   k_action=k_act,
                                   alpha_action=1/inv_alpha_act,
                                   k_observation=k_obs,
                                   alpha_observation=1/inv_alpha_obs,
                                   estimate_value=FORollout(rollout_policy),
                                   next_action = RootToNextMLFirstVOO(voo_opt),
                                   default_action=BasicPOMCP.ReportWhenUsed(TagAction(false, 0.0)),
                                   rng=rng
                                  )

            planner = solve(solver, pomdp)
            filter = SIRParticleFilter(deepcopy(pomdp),
                                       100_000,
                                       rng=MersenneTwister(i+90_000))

            seed!(planner, i+40000)
            sim = Sim(deepcopy(pomdp),
                      planner,
                      filter,
                      rng=MersenneTwister(i+50_000),
                      max_steps=100,
                      metadata=Dict(:i=>i, :solver=>solver_name)
                     )

        end

        result = POMDPs.simulate(sim)
        output = (reward=discounted_reward(result),)
        return merge(sim.metadata, output)
    end
end

datestring = Dates.format(now(), "e_d_u_Y_HH_MM")
results_csv = DataFrame(max_time=Float64[], solver=String[], mean_reward=Float64[], se_reward=Float64[])

try
    for iter_time in max_time_list
        sims = SimWrap[]
        for solver_name in keys(solver_list)
            for j in 1:n_iters
                sim_wrap = SimWrap(solver_list[solver_name], solver_name, j, iter_time, max_depth, max_query, deepcopy(pomdp))
                push!(sims, sim_wrap)
            end
        end
        if parallelize
            results = @showprogress pmap(gen_sim, sims)
        else
            results = @showprogress map(gen_sim, sims)
        end
        println(stderr, "Max time: $iter_time :")
        results = POMDPSimulators.create_dataframe(results)
        println(stderr, results)
        combined = by(results, :solver) do df
            DataFrame(mean_reward=mean(df[:reward]), se_reward=(std(df[:reward])/sqrt(n_iters)))
        end
        for row in eachrow(combined)
            result_row = [iter_time, row.solver, row.mean_reward, row.se_reward]
            println(stderr, @sprintf("- %s reward: %6.3f ± %6.3f", row.solver, row.mean_reward, row.se_reward))
            push!(results_csv, result_row)
        end
        if save_data
            fname_csv = joinpath(@__DIR__, "../data", "VDP_compare_results"*"_"*datestring*".csv")
            CSV.write(fname_csv, results_csv)
            println(stderr,"Saving parameters to ", fname_csv)
        end
    end
finally
    Distributed.rmprocs(worker_ids)
end


text:
using FromFile
@from "Core.jl" import Options, Dataset, Node, copyNode
@from "Utils.jl" import getTime
@from "LossFunctions.jl" import scoreFunc

# Define a member of population by equation, score, and age
mutable struct PopMember{T<:Real}
    tree::Node
    score::T
    birth::Int
end

"""
    PopMember(t::Node, score::T)

Create a population member with a birth date at the current time.

# Arguments

- `t::Node`: The tree for the population member.
- `score::T`: The loss to assign this member.
"""
function PopMember(t::Node, score::T) where {T<:Real}
    PopMember{T}(t, score, getTime())
end

"""
    PopMember(dataset::Dataset{T}, baseline::T,
              t::Node, options::Options)

Create a population member with a birth date at the current time.
Automatically compute the score for this tree.

# Arguments

- `dataset::Dataset{T}`: The dataset to evaluate the tree on.
- `baseline::T`: The baseline loss.
- `t::Node`: The tree for the population member.
- `options::Options`: What options to use.
"""
function PopMember(dataset::Dataset{T},
                   baseline::T, t::Node,
                   options::Options) where {T<:Real}
    PopMember(t, scoreFunc(dataset, baseline, t, options))
end

function copyPopMember(p::PopMember{T}) where {T<:Real}
    tree = copyNode(p.tree)
    score = p.score
    birth = p.birth
    return PopMember{T}(tree, score, birth)
end


text:
# """
# # module ConstructiveSolidGeometry
# """
module ConstructiveSolidGeometry

    # Base Packages
    using LinearAlgebra

    # # Other Packages
    using IntervalSets
    using JSON
    using Parameters
    using PolygonOps
    using Polynomials
    using RecipesBase
    using Rotations
    using StaticArrays
    using Statistics
    using Unitful
    using YAML

    using DataStructures: OrderedDict

    import Base: in, *, +, -, &, size, zero

    abstract type AbstractCoordinateSystem end
    abstract type Cartesian <: AbstractCoordinateSystem end
    abstract type Cylindrical <: AbstractCoordinateSystem end
    const CoordinateSystemType = Union{Type{Cartesian}, Type{Cylindrical}}

    import Base: print, show
    print(io::IO, ::Type{Cartesian}) = print(io, "Cartesian")
    print(io::IO, ::Type{Cylindrical}) = print(io, "Cylindrical")
    show(io::IO, CS::CoordinateSystemType) = print(io, CS)
    show(io::IO,::MIME"text/plain", CS::CoordinateSystemType) = show(io, CS)

    # Tuples with ticks to sample with differently spaced ticks
    const CartesianTicksTuple{T} = NamedTuple{(:x,:y,:z), NTuple{3,Vector{T}}}
    const CylindricalTicksTuple{T} = NamedTuple{(:r,:φ,:z), NTuple{3,Vector{T}}}

    abstract type AbstractGeometry{T <: AbstractFloat} end

    abstract type AbstractPrimitive{T} <: AbstractGeometry{T} end
    abstract type ClosedPrimitive end
    abstract type OpenPrimitive end

    abstract type AbstractVolumePrimitive{T, CO <: Union{ClosedPrimitive, OpenPrimitive}} <: AbstractPrimitive{T} end
    abstract type AbstractSurfacePrimitive{T} <: AbstractPrimitive{T} end
    abstract type AbstractLinePrimitive{T} <: AbstractPrimitive{T} end

    abstract type AbstractConstructiveGeometry{T} <: AbstractGeometry{T} end

    _csg_convert_args(eltype::Type{T}, r::Real) where T = convert(T, r)
    _csg_convert_args(eltype::Type{T}, r::Tuple) where T = broadcast(x -> _csg_convert_args(T, x), r)
    _csg_convert_args(eltype::Type{T}, r::Nothing) where T = nothing

    _csg_get_promoted_eltype(::Type{T}) where {T <: AbstractArray} = eltype(T)
    _csg_get_promoted_eltype(::Type{T}) where {T <: Real} = T
    _csg_get_promoted_eltype(::Type{Nothing}) = Int
    _csg_get_promoted_eltype(::Type{Tuple{T}}) where {T<:Real} = T
    _csg_get_promoted_eltype(::Type{Tuple{T1,T2}}) where {T1<:Real, T2<:Real} = promote_type(T1, T2)
    _csg_get_promoted_eltype(::Type{Tuple{T1,T2}}) where {T1<:Union{Real, Tuple}, T2<:Union{Real, Tuple}} = promote_type(_csg_get_promoted_eltype(T1), _csg_get_promoted_eltype(T2))
    _csg_get_promoted_eltype(::Type{Tuple{Nothing,T2}}) where {T2<:Union{Real, Tuple}} = _csg_get_promoted_eltype(T2)
    _csg_get_promoted_eltype(::Type{Tuple{T1,Nothing}}) where {T1<:Union{Real, Tuple}} = _csg_get_promoted_eltype(T1)

    _handle_phi(φ, rotation) = (φ, rotation)
    _handle_phi(φ::Tuple, rotation) = (abs(φ[2]-φ[1]), rotation*RotZ(φ[1]))

    include("Units.jl")
    include("PointsAndVectors/PointsAndVectors.jl")

    rotation(p::AbstractPrimitive) = p.rotation
    origin(p::AbstractPrimitive) = p.origin
    _transform_into_global_coordinate_system(pt::CartesianPoint, p::AbstractPrimitive) = (rotation(p) * pt) + origin(p)
    _transform_into_global_coordinate_system(pts::AbstractVector{<:CartesianPoint}, p::AbstractPrimitive) =
        broadcast(pt -> _transform_into_global_coordinate_system(pt, p), pts)
    _transform_into_object_coordinate_system(pt::CartesianPoint, p::AbstractPrimitive) = inv(rotation(p)) * (pt - origin(p))
    in(pt::CartesianPoint{T}, p::AbstractPrimitive{T}, csgtol::T = csg_default_tol(T)) where {T} =
        _in(_transform_into_object_coordinate_system(pt, p), p; csgtol = csgtol)
    in(pt::CylindricalPoint{T}, p::AbstractPrimitive{T}, csgtol::T = csg_default_tol(T)) where {T} =
        in(CartesianPoint(pt), p, csgtol)

    """
        extreme_points(es::AbstractPrimitive{T}) where {T}

    Generic fallback of `extreme_points` for any primitive.
    Returns 6 `CartesianPoint`s in both directions of each cartesian axes (x, y and z)
    around the origin of the primitive with distance determined by `extremum(es)`,
    which returns the maximum distance of an primitive to its origin.

    ## Arguments
    * `es::AbstractPrimitive{T}`: Any primitive, e.g. a `Box`.
    """
    function extreme_points(es::AbstractPrimitive{T}) where {T}
        o = origin(es)
        r = extremum(es)
        vX = r * CartesianVector{T}(one(T), zero(T), zero(T))
        vY = r * CartesianVector{T}(zero(T), one(T), zero(T))
        vZ = r * CartesianVector{T}(zero(T), zero(T), one(T))
        SVector{6,CartesianPoint{T}}(
            o - vX, o + vX, o - vY, o + vY, o - vZ, o + vZ,
        )
    end

    get_scale(es::AbstractPrimitive) = extremum(es)

    include("Transformations.jl")
    include("GeometryRounding.jl")
    include("VolumePrimitives/VolumePrimitives.jl")
    include("LinePrimitives/LinePrimitives.jl")
    include("SurfacePrimitives/SurfacePrimitives.jl")
    include("Intervals.jl")
    include("CSG.jl")
    include("IO.jl")

    # Plotting
    include("plotting/plotting.jl")

end


text:
# Use baremodule to shave off a few KB from the serialized `.ji` file
baremodule SPRAL_jll
using Base
using Base: UUID
import JLLWrappers

JLLWrappers.@generate_main_file_header("SPRAL")
JLLWrappers.@generate_main_file("SPRAL", UUID("319450e9-13b8-58e8-aa9f-8fd1420848ab"))
end  # module SPRAL_jll


text:
#
# This file is a part of MolecularGraph.jl
# Licensed under the MIT License http://opensource.org/licenses/MIT
#

@testset "graph.isomorphism.cliquebased" begin

@testset "mcis" begin
    nullg = plaingraph()
    @test mcissize(nullg, nullg) == 0

    g = plaingraph(6, [(1,2), (2,3), (2,4), (2,5), (5,6)])
    h = plaingraph(6, [(1,2), (2,3), (3,4), (3,5), (3,6)])
    @test mcissize(g, h) == 6

    g = pathgraph(4)
    h = cyclegraph(4)
    @test mcissize(g, h) == 3

    g = completegraph(7, mutable=true)
    h = unlinkedges(g, (12,))
    @test mcissize(g, h) == 6
end

@testset "mces" begin
    nullg = plaingraph()
    @test mcessize(nullg, nullg) == 0

    g = plaingraph(6, [(1,2), (2,3), (2,4), (2,5), (5,6)])
    h = plaingraph(6, [(1,2), (2,3), (3,4), (3,5), (3,6)])
    @test mcessize(g, h) == 5

    g = pathgraph(4)
    h = cyclegraph(4)
    @test mcessize(g, h) == 3

    g = completegraph(4, mutable=true)
    h = unlinkedges(g, (6,))
    @test mcessize(g, h) == 5
end

@testset "connected" begin
    g = disjointunion(completegraph(5), completegraph(4))
    h = disjointunion(completegraph(4), completegraph(3))
    @test mcissize(g, h, connected=false) == 7
    @test mcissize(g, h, connected=true) == 4
end

end # graph.isomorphism.cliquebased