text:
# A modification of Tang's table-driven logarithm for different base.

const t_log2_64 = Array((Float64,Float64),129)
N=39 # (can be up to N=42, which appears to be what Apple's libm uses).
sN = 2.0^N
isN = 1.0/sN
s7 = 2.0^7
is7 = 1.0/s7

for j=0:128
    l_big = Base.log2(big(1.0+j*is7))
    l_hi = isN*float64(round(sN*l_big))
    l_lo = float64(l_big-l_hi)
    t_log2_64[j+1] = (l_hi,l_lo)
    # c_invF[j+1] = 1.0+is7*j
end

const inv_ln2 = 0x1.71547652b82fep+0
const inv_ln2_hi = 0x1.7154766p+0 # 28 bits
const inv_ln2_lo = -0x1.a8fa03d1105ep-29

# const t_log_32 = Array((Float32,Float32),129)
# N=16
# sN = 2f0^N
# isN = 1f0/sN
# for j=0:128
#     l_big = Base.log(big(1.0+j*is7))
#     l_hi = isN*float32(round(sN*l_big))
#     l_lo = float32(l_big-l_hi)
#     t_log_32[j+1] = (l_hi,l_lo)
# end

# Procedure 1
@inline function log2_proc1(y::Float64,mf::Float64,F::Float64,f::Float64,jp::Int)
    ## Steps 1 and 2
    @inbounds hi,lo = libm.t_log2_64[jp]
    l_hi = mf + hi
    l_lo = lo

    ## Step 3' (alternative)
    u = (2.0f)/(y+F)
    us = inv_ln2*u
    v = u*u
    q = us*v*@horner(v,
                    0x1.5_5555_5555_0286p-4,
                    0x1.9_99a0_bc71_2416p-7)

    ## Step 4
    l_hi + (us + (q + l_lo))
end

@inline function log2_proc2(f::Float64)
    ## Step 1
    g = 1.0/(2.0+f)
    u = 2.0*f*g
    us = 1.4426950408889634*u
    v = u*u

    ## Step 2
    q = us*v*@horner(v,
                    0x1.5_5555_5555_54e6p-4,
                    0x1.9_9999_99ba_c6d4p-7,
                    0x1.2_4923_07f1_519fp-9,
                    0x1.c_8034_c85d_fff0p-12)

    ## Step 3
    u1 = float64(float32(u)) # round to 24 bits
    f1 = float64(float32(f))
    f2 = f-f1
    u2 = ((2.0*(f-u1)-u1*f1)-u1*f2)*g

    ## Step 4
    u1*inv_ln2_hi + (u2*inv_ln2 + (u1*inv_ln2_lo + q))
end

function log2_tang(x::Float64)
    if x > 0.0
        x == Inf && return x

        # Step 2
        if 0x1.e_0fab_fbc7_02a3p-1 < x < 0x1.1_082b_577d_34eep0
            f = x-1.0
            return log2_proc2(f)
        end

        # Step 3
        xu = reinterpret(Uint64,x)
        m = int(xu >> 52) & 0x07ff
        if m == 0 # x is subnormal
            x *= 0x1p54 # normalise significand
            xu = reinterpret(Uint64,x)
            m = int(xu >> 52) & 0x07ff - 54
        end
        m -= 1023
        y = reinterpret(Float64,(xu & 0x000f_ffff_ffff_ffff) | 0x3ff0_0000_0000_0000)

        mf = Float64(m)
        F = (y + 0x1p45) - 0x1p45 # 0x1p-7*round(0x1p7*y)
        f = y-F
        jp = itrunc(0x1p7*F)-127

        log2_proc1(y,mf,F,f,jp)
    elseif x == 0.0
        -Inf
    elseif isnan(x)
        NaN
    else
        throw(DomainError())
    end
end

@vectorize_1arg Real log2_tang


text:
#  Copyright 2016, Iain Dunning, Joey Huchette, Miles Lubin, and contributors
#  This Source Code Form is subject to the terms of the Mozilla Public
#  License, v. 2.0. If a copy of the MPL was not distributed with this
#  file, You can obtain one at http://mozilla.org/MPL/2.0/.
using JuMP

# clnlbeam
# Based on AMPL model
# Copyright (C) 2001 Princeton University
# All Rights Reserved
#
# Permission to use, copy, modify, and distribute this software and
# its documentation for any purpose and without fee is hereby
# granted, provided that the above copyright notice appear in all
# copies and that the copyright notice and this
# permission notice appear in all supporting documentation.

#   Source:
#   H. Maurer and H.D. Mittelman,
#   "The non-linear beam via optimal control with bound state variables",
#   Optimal Control Applications and Methods 12, pp. 19-31, 1991.
let
    N     = 1000
    ni    = N
    h     = 1/ni
    alpha = 350

    m = Model()

    @variable(m, -1 <= t[1:(ni+1)] <= 1)
    @variable(m, -0.05 <= x[1:(ni+1)] <= 0.05)
    @variable(m, u[1:(ni+1)])

    @NLobjective(m, Min, sum{ 0.5*h*(u[i+1]^2 + u[i]^2) + 0.5*alpha*h*(cos(t[i+1]) + cos(t[i])), i = 1:ni})

    # cons1
    for i in 1:ni
        @NLconstraint(m, x[i+1] - x[i] - (0.5h)*(sin(t[i+1])+sin(t[i])) == 0)
    end
    # cons2
    for i in 1:ni
        @constraint(m, t[i+1] - t[i] - (0.5h)*u[i+1] - (0.5h)*u[i] == 0)
    end

    solve(m)

end

ext: jl


text:
#
#  The following functions allow us to do math on vectors and NaNs
#
#   dependencies: StatsBase
#
#   nanmean(x::Vector{Float64})
#   nansum(x::Vector{Float64})
#   nanmat(r, c)
#   nan_elementwise_multiply(x::Vector{Float64}, y::Vector{Float64})
#   nan_elementwise_sum(x::Vector{Float64}, y::Vector{Float64})
#   mean(y::Vector)
#   std(x::Vector; force_nan_zero=false)
#   checksign(p1, p2)

using StatsBase

function nanmean(x::Vector{Float64})
    xx = filter(!isnan, x)
    sum(xx)/length(xx)
end
function nansum(x::Vector{Float64})
    xx = filter(!isnan, x)
    sum(xx)
end;
function nan_elementwise_multiply(x::Vector{Float64}, y::Vector{Float64})
    xx = findall(x->isnan(x), x)
    x[xx] .= 1;
    yy = findall(y->isnan(y), y)
    y[yy] .= 1;
    x.*y
end
function nan_elementwise_sum(x::Vector{Float64}, y::Vector{Float64})
    xx = findall(x->isnan(x), x)
    x[xx] .= 0;
    yy = findall(y->isnan(y), y)
    y[yy] .= 0;
    x.+y
end
function nanmat(r,c)
    NaN.*ones(r,c)
end
# function mean(y)
#     sum(y)/length(y)
# end
function nanstd(x; force_nan_zero=false)
    y = StatsBase.std(x)
    if isnan(y) && force_nan_zero
        y=0
    else
        return y
    end
end

function checksign(d1,d2)
   sign(d1)==sign(d2)
end


text:
module TestQueue
using Test
using StructDatabaseMapping
const Queue = StructDatabaseMapping.Queue
function test_queue()
    a = Queue(3)
    push!(a, 5)
    push!(a, 5)
    push!(a, 1)
    @test_throws StructDatabaseMapping.Full StructDatabaseMapping.push!_nowait(a, 12)
end
end

max_stars_repo_name: lucianolorenti/StructDatabaseMapping

max_issues_repo_name: lucianolorenti/StructDatabaseMapping

max_forks_repo_name: lucianolorenti/StructDatabaseMapping


text:
module TestReplica
using Omega

function testreplica()
  x = normal(0.0, 1.0)
  y = normal(10.0, 3.0)
  whi = bernoulli(0.5) + 1
  z = randarray([x, y])
  a = z[whi]
  o = normal(a, 1.0)
  rand(a, o ==ₛ 23.0)
end

# testreplica()

function testreplica2()
  x = normal(0.0, 1.0)
  rand(x, x >ₛ 0.0, 1000; alg = Replica, swapevery = 10)
end

function testreplica3()
  x = normal(0.0, 1.0)
  y = normal(0.0, 1.0)
  samples = rand(x, x ==ₛ y, 1000; alg = Replica, swapevery = 10)
  (samples = samples, x = x, y = y)
  s = randtuple((x, y)).(samples)
  xs = [sa[2] for sa in s]
  ys = [sa[2] for sa in s]
  scatter(xs, ys)
end

function testreplica4()
  x = normal(0.0, 1.0)
  y = normal(0.0, 1.0)
  sam = withkernel(kseα(9.7448e9)) do
    rand(x, x ==ₛ y, 1000; alg = SSMH)
  end
  s = randtuple((x, y)).(sam)
  xs = [sa[2] for sa in s]
  ys = [sa[2] for sa in s]
  scatter(xs, ys)
end

function testreplica5()
  x = uniform(0.0, 1.0)
  c = ((x >ₛ 0.05) & (x <ₛ 0.15)) | ((x >ₛ 0.75) & (x <ₛ 0.85))
  # samples = withkernel(kseα(9.7448e9)) do
  #   samples = rand(x, c, 1000; alg = SSMH)
  # end
  samples = rand(x, c, 10000; alg = Replica, swapevery = 10)
  # (samples = samples, x = x, y = y)
  @show x.(samples)
  histogram(x.(samples), xlim = (0, 1), nbins=50)
end

function testreplica7(n; t = 1000, kwargs...)
  x = uniform(-1.0, 1.0)
  y = uniform(-1.0, 1.0)
  xy = randtuple((x, y))
  c = abs(x) ==ₛ abs(y)
  s = withkernel(kseα(t)) do
    rand(xy, c, n; kwargs...)
  end
  # s = rand(randtuple((x, y)), c, 10000; alg = Replica, swapevery = 10, kwargs...)
  xs = [sa[1] for sa in s]
  ys = [sa[2] for sa in s]
  @show length(xs)
  scatter(xs, ys, xlim=(-1, 1), ylim = (-1, 1))
end

end


text:
# Resampling

# Taken from Turing.jl
function resample_systematic(w::AbstractVector{<:Real}, num_particles::Integer)

    Q, N = cumsum(w), num_particles

    T = collect(range(0, stop = maximum(Q)-1/N, length = N)) .+ rand()/N
    push!(T, 1)

    indx, i, j = Array{Int}(undef, N), 1, 1
    while i <= N
        if T[i] < Q[j]
            indx[i] = j
            i += 1
        else
            j += 1
        end
    end
    return indx
end

function resample_stratified(w::AbstractVector{<:Real}, num_particles::Integer)

    Q, N = cumsum(w), num_particles

    T = Array{Float64}(undef, N + 1)
    for i=1:N,
        T[i] = rand() / N + (i - 1) / N
    end
    T[N+1] = 1

    indx, i, j = Array{Int}(undef, N), 1, 1
    while i <= N
        if T[i] < Q[j]
            indx[i] = j
            i += 1
        else
            j += 1
        end
    end
    return indx
end


text:
#!/usr/bin/env julia

using ArgParse

function compute_pi(a::Float64, b::Float64, n::Int64)::Float64
    result = zeros(Threads.nthreads())
    delta = (b - a)/n
    Threads.@threads for i in 0:(n - 1)
        result[Threads.threadid()] += sqrt(1.0 - (a + i*delta)^2)
    end
    return 2.0*delta*sum(result)
end

function main()
    arg_parser = ArgParseSettings()
    arg_parser.description = "Compute pi using a quadrature method on multiple threads"
    @add_arg_table arg_parser begin
        "n"
            help = "number of point for the quadrature method"
            arg_type = Int64
            required = true
    end
    options = parse_args(ARGS, arg_parser)
    nr_threads = Threads.nthreads()
    println("running with $nr_threads threads")
    nr_points = options["n"]
    println(compute_pi(-1.0, 1.0, nr_points))
end

main()


text:
# Conjugate

"""
  Conjugate(f::ProximableFunction)

Returns the conjugate function of `f`, that is `f*(x) = sup{y'x - f(y)}`.
"""

immutable Conjugate{T <: ProximableFunction} <: ProximableFunction
  f::T
end

fun_dom(f::Conjugate) = fun_dom(f.f)

# only prox! is provided here, call method would require being able to compute
# an element of the subdifferential of the conjugate

function prox!{R <: Real}(g::Conjugate, x::AbstractArray{R}, y::AbstractArray{R}, gamma::Real=1.0)
  # need to make a copy, doing this in place is probably not possible in general:
  # if object_id(x) == object_id(y), then prox! causes the loss of x (which is needed afterwards)
  x_copy = copy(x)
  # Moreau identity
  v = prox!(g.f, x/gamma, y, 1.0/gamma)
  v = vecdot(x_copy,y) - gamma*vecdot(y,y) - v
  for k in eachindex(y)
    y[k] = x_copy[k] - gamma*y[k]
  end
  return v
end

# complex case, need to cast inner products to real

function prox!{R <: Real}(g::Conjugate, x::AbstractArray{Complex{R}}, y::AbstractArray{Complex{R}}, gamma::Real=1.0)
  x_copy = copy(x)
  v = prox!(g.f, x/gamma, y, 1.0/gamma)
  v = real(vecdot(x_copy,y)) - gamma*real(vecdot(y,y)) - v
  for k in eachindex(y)
    y[k] = x_copy[k] - gamma*y[k]
  end
  return v
end

# special case for indicator of convex cones (the conjugate is also an indicator function)

prox!{R <: Real, C <: IndicatorConvexCone}(g::Conjugate{C}, x::AbstractArray{R}, y::AbstractArray{R}, gamma::Real=1.0) =
  prox_conjugate_convex_cone!(g.f, x, y, gamma)

prox!{R <: Real, C <: IndicatorConvexCone}(g::Conjugate{C}, x::AbstractArray{Complex{R}}, y::AbstractArray{Complex{R}}, gamma::Real=1.0) =
  prox_conjugate_convex_cone!(g.f, x, y, gamma)

function prox_conjugate_convex_cone!{T <: RealOrComplex, C <: IndicatorConvexCone}(f::C, x::AbstractArray{T}, y::AbstractArray{T}, gamma::Real=1.0)
  x_copy = copy(x)
  prox!(f, x/gamma, y, 1.0/gamma)
  for k in eachindex(y)
    y[k] = x_copy[k] - gamma*y[k]
  end
  return 0.0
end

# naive implementation

function prox_naive{T <: RealOrComplex}(g::Conjugate, x::AbstractArray{T}, gamma::Real=1.0)
  y, v = prox_naive(g.f, x/gamma, 1.0/gamma)
  return x - gamma*y, real(vecdot(x,y)) - gamma*real(vecdot(y,y)) - v
end

is_prox_accurate(f::Conjugate) = is_prox_accurate(f.f)


text:
using ConstraintSolver, JuMP
CS = ConstraintSolver

function read_puzzle(pname)
    dir = pkgdir(ConstraintSolver)
    lines = readlines(pname)
    width, height = parse.(Int, split(lines[1]))
    lines = lines[2:end]
    npieces = length(lines)
    puzzle = zeros(Int, (npieces, 5))
    for i in 1:npieces
        puzzle[i, 1] = i
        parts = split(lines[i], " ")
        puzzle[i, 2:end] = parse.(Int, parts)
    end
    return puzzle, width, height
end

function get_rotations(puzzle)
    npieces = size(puzzle)[1]
    rotations = zeros(Int, (npieces * 4, 5))
    rotation_indices = [[2, 3, 4, 5], [3, 4, 5, 2], [4, 5, 2, 3], [5, 2, 3, 4]]
    for i in 1:npieces
        j = 1
        for rotation in rotation_indices
            rotations[(i - 1) * 4 + j, 1] = i
            rotations[(i - 1) * 4 + j, 2:end] = puzzle[i, rotation]
            j += 1
        end
    end

    return rotations
end

function main(pname; time_limit = 1800)
    puzzle, width, height = read_puzzle(pname)
    rotations = get_rotations(puzzle)
    npieces = size(puzzle)[1]
    ncolors = maximum(puzzle[:, 2:end])

    m = Model(optimizer_with_attributes(CS.Optimizer, "time_limit" => time_limit))

    @variable(m, 1 <= p[1:height, 1:width] <= npieces, Int)
    @variable(m, 0 <= pu[1:height, 1:width] <= ncolors, Int)
    @variable(m, 0 <= pr[1:height, 1:width] <= ncolors, Int)
    @variable(m, 0 <= pd[1:height, 1:width] <= ncolors, Int)
    @variable(m, 0 <= pl[1:height, 1:width] <= ncolors, Int)

    @constraint(m, p[:] in CS.AllDifferent())
    for i in 1:height, j in 1:width
        @constraint(
            m,
            [p[i, j], pu[i, j], pr[i, j], pd[i, j], pl[i, j]] in CS.TableSet(rotations)
        )
    end

    # borders
    # up and down
    for j in 1:width
        @constraint(m, pu[1, j] == 0)
        @constraint(m, pd[height, j] == 0)

        if j != width
            @constraint(m, pr[1, j] == pl[1, j + 1])
            @constraint(m, pr[height, j] == pl[height, j + 1])
        end
    end

    # right and left
    for i in 1:height
        @constraint(m, pl[i, 1] == 0)
        @constraint(m, pr[i, width] == 0)

        if i != height
            @constraint(m, pd[i, 1] == pu[i + 1, 1])
            @constraint(m, pd[i, width] == pu[i + 1, width])
        end
    end

    for i in 1:(height - 1), j in 1:(width - 1)
        @constraint(m, pd[i, j] == pu[i + 1, j])
        @constraint(m, pr[i, j] == pl[i, j + 1])
    end

    if width == height
        start_piece = findfirst(i -> count(c -> c == 0, puzzle[i, :]) == 2, 1:npieces)
        @constraint(m, p[1, 1] == start_piece)
    end

    optimize!(m)

    status = JuMP.termination_status(m)
    if status == MOI.OPTIMAL
        print("$status, $(JuMP.objective_value(m)), $(JuMP.solve_time(m))")
    else
        print("$status, NaN, $(time_limit)")
    end
end


text:
export select_device

import MPI
@static if ENABLE_CUDA
    using CUDA
end

"""
    select_device()

Select the device (GPU) corresponding to the node-local MPI rank and return its ID. This function only needs to be called when using nodes with more than one device.
"""
function select_device()
    check_initialized();
    @static if ENABLE_CUDA
        @assert CUDA.functional(true)
        comm_l = MPI.Comm_split_type(comm(), MPI.MPI_COMM_TYPE_SHARED, me())
    	if (MPI.Comm_size(comm_l) > length(CUDA.devices())) error("More processes have been launched per node than there are GPUs available."); end
    	me_l = MPI.Comm_rank(comm_l)
        CUDA.device!(me_l)
        return me_l
    else
        error("Cannot select a GPU because ImplicitGlobalGrid was not precompiled for GPU usage (as CUDA was not functional).")
    end
end


text:
export rand, Ginibre
import Base.rand

#Samples a matrix from the Ginibre ensemble
#This ensemble lives in GL(N, F), the set of all invertible N x N matrices
#over the field F
#For beta=1,2,4, F=R, C, H respectively
struct Ginibre <: ContinuousMatrixDistribution
   beta::Float64
   N::Integer
end

function rand(W::Ginibre)
    beta, n = W.beta, W.N
    if beta==1
        randn(n,n)
    elseif beta==2
        randn(n,n)+im*randn(n,n)
    elseif beta==4
        Q0=randn(n,n)
        Q1=randn(n,n)
        Q2=randn(n,n)
        Q3=randn(n,n)
        [Q0+im*Q1 Q2+im*Q3;-Q2+im*Q3 Q0-im*Q1]
    else
        error(string("beta = ", beta, " not implemented"))
    end
end

function jpdf(Z::AbstractMatrix{z}) where {z<:Complex}
    pi^(size(Z,1)^2)*exp(-trace(Z'*Z))
end

#Dyson Circular orthogonal, unitary and symplectic ensembles
struct CircularOrthogonal
    N :: Int64
end

function rand(M::CircularOrthogonal)
    U = rand(Ginibre(2, M.N))
    U * U'
end

struct CircularUnitary
    N :: Int64
end

rand(M::CircularUnitary) = rand(Ginibre(2, M.N))

struct CircularSymplectic
    N :: Int64
end

rand(M::CircularSymplectic) = error("Not implemented")


text:
function CostGradient(X::Array,y::Vector;λ=0, Norm=false)

   #Size of the X Parameters
    m=size(X,1)
    n=size(X,2)
    #if Normalize

    if Norm==true
       μ=mean(X,dims=1)
       σ=std(X,dims=1)
       X=(X.-μ)./σ
    else
        μ=zeros(1,n)
        σ=ones(1,n)
    end

    #Add vector of ones to X parameters
    X=[ones(m) X]

    CostFunction=function (θ)
                    se=(X*θ).-y
                    J=(0.5./m).*((se'*se)+λ*(θ[2:end]'*θ[2:end]))
                    return J
                end
    θr=zeros(n+1)
    Gradient=function (θ)
             θr.=θ
             θr[1]=0
             g=(1/m).*(X'*((X*θ).-y)) .+ ((λ./m).*θr)
        return g
            end

    return CostFunction, Gradient, μ, σ
end


text:
using Plots
using Random
Random.seed!(0)
A = randn(40,30)
heatmap(A)


text:
module Dash
using DashBase
import HTTP, JSON2, CodecZlib, MD5
using Sockets
const ROOT_PATH = realpath(joinpath(@__DIR__, ".."))
const RESOURCE_PATH = realpath(joinpath(ROOT_PATH, "resources"))
include("exceptions.jl")
include("Components.jl")
include("Front.jl")
include("HttpHelpers/HttpHelpers.jl")

using .HttpHelpers

export dash, Component, Front, callback!,
enable_dev_tools!, ClientsideFunction,
run_server, PreventUpdate, no_update, @var_str,
Input, Output, State, make_handler, callback_context,
ALL, MATCH, ALLSMALLER

include("Contexts/Contexts.jl")
include("env.jl")
include("utils.jl")
include("app.jl")
include("resources/application.jl")
include("handlers.jl")
include("server.jl")

@doc """
    module Dash

Julia backend for [Plotly Dash](https://github.com/plotly/dash)

# Examples
```julia

using Dash
using DashHtmlComponents
using DashCoreComponents
app = dash(external_stylesheets=["https://codepen.io/chriddyp/pen/bWLwgP.css"]) do
    html_div() do
        dcc_input(id="graphTitle", value="Let's Dance!", type = "text"),
        html_div(id="outputID"),
        dcc_graph(id="graph",
            figure = (
                data = [(x = [1,2,3], y = [3,2,8], type="bar")],
                layout = Dict(:title => "Graph")
            )
        )

    end
end
callback!(app, Output("outputID", "children"), Input("graphTitle","value"), State("graphTitle","type")) do value, type
    "You've entered: '\$(value)' into a '\$(type)' input control"
end
callback!(app, Output("graph", "figure"), Input("graphTitle", "value")) do value
    (
        data = [
            (x = [1,2,3], y = abs.(randn(3)), type="bar"),
            (x = [1,2,3], y = abs.(randn(3)), type="scatter", mode = "lines+markers", line = (width = 4,))
        ],
        layout = (title = value,)
    )
end
run_server(app, HTTP.Sockets.localhost, 8050)
```

""" Dash

function __init__()
    DashBase.main_registry().dash_dependency = (
        dev = ResourcePkg(
            "dash_renderer",
            RESOURCE_PATH, version = "1.5.0",
            [
                Resource(
                relative_package_path = "react@16.13.0.js",
                external_url = "https://unpkg.com/react@16.13.0/umd/react.development.js"
                ),
                Resource(
                    relative_package_path = "react-dom@16.13.0.js",
                    external_url = "https://unpkg.com/react-dom@16.13.0/umd/react-dom.development.js"
                ),
                Resource(
                    relative_package_path = "polyfill@7.8.7.min.js",
                    external_url = "https://unpkg.com/@babel/polyfill@7.8.7/dist/polyfill.min.js"
                ),
                Resource(
                    relative_package_path = "prop-types@15.7.2.js",
                    external_url = "https://unpkg.com/prop-types@15.7.2/prop-types.js",
                ),
            ]
        ),
        prod = ResourcePkg(
            "dash_renderer",
            RESOURCE_PATH, version = "1.2.2",
            [
                Resource(
                relative_package_path = "react@16.13.0.min.js",
                external_url = "https://unpkg.com/react@16.13.0/umd/react.production.min.js"
                ),
                Resource(
                    relative_package_path = "react-dom@16.13.0.min.js",
                    external_url = "https://unpkg.com/react-dom@16.13.0/umd/react-dom.production.min.js"
                ),
                Resource(
                    relative_package_path = "polyfill@7.8.7.min.js",
                    external_url = "https://unpkg.com/@babel/polyfill@7.8.7/dist/polyfill.min.js"
                ),
                Resource(
                    relative_package_path = "prop-types@15.7.2.min.js",
                    external_url = "https://unpkg.com/prop-types@15.7.2/prop-types.min.js"
                ),
            ]
        )
    )

    DashBase.main_registry().dash_renderer = ResourcePkg(
        "dash_renderer",
        RESOURCE_PATH, version = "1.5.0",
        [
            Resource(
                relative_package_path = "dash_renderer.min.js",
                dev_package_path = "dash_renderer.dev.js",
                external_url = "https://unpkg.com/dash-renderer@1.5.0/dash_renderer/dash_renderer.min.js"
            ),
            Resource(
                relative_package_path = "dash_renderer.min.js.map",
                dev_package_path = "dash_renderer.dev.js.map",
                dynamic = true,
            ),
        ]
    )

end

end # module


text:
using skalman
using Test

@testset "skalman.jl" begin
    # Write your tests here.
    @test legendrepoly(0.0,0) == [1.0]

end


text:
"""
    AbstractPermutation

Abstract type representing a compile-time permutation.

Subtypes are [`Permutation`](@ref) and [`NoPermutation`](@ref).
"""
abstract type AbstractPermutation end

"""
    Permutation{p} <: AbstractPermutation

Describes a compile-time dimension permutation.

The type parameter `p` should be a valid permutation such as `(3, 1, 2)`.

---

    Permutation(perm::Vararg{Int})
    Permutation(perm::NTuple{N,Int})

Constructs a `Permutation`.

# Example

Both are equivalent:

```julia
p1 = Permutation(2, 3, 1)
p2 = Permutation((2, 3, 1))
```
"""
struct Permutation{p,N} <: AbstractPermutation
    @inline Permutation(perm::Vararg{Int}) = new{perm, length(perm)}()
end

@inline Permutation(perm::Tuple) = Permutation(perm...)

@inline Base.getindex(::Permutation{p}, ::Val{i}) where {p,i} = p[i]
@inline Base.getindex(p::Permutation, i::Integer) = p[Val(i)]

Base.show(io::IO, ::Permutation{p}) where {p} = print(io, "Permutation", p)

"""
    NoPermutation <: AbstractPermutation

Represents an identity permutation.

It is functionally equivalent to `Permutation(1, 2, 3, ...)`, and is included
for convenience.
"""
struct NoPermutation <: AbstractPermutation end

@inline Base.getindex(::NoPermutation, ::Val{i}) where {i} = i
@inline Base.getindex(::NoPermutation, i::Integer) = i

Base.show(io::IO, ::NoPermutation) = print(io, "NoPermutation()")


text:
module DesignStructureMatrix

using Luxor
using LinearAlgebra

export plotDSM
export toReachableMatrix
export OrderReachable
export Sequencing
export Clustering

include("plotDSM.jl")
include("sequencing.jl")
include("clustering.jl")

end


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

include("model_Resource.jl")
include("model_RouteListResult.jl")
include("model_RoutePropertiesFormat.jl")
include("model_RouteTableListResult.jl")
include("model_RouteTablePropertiesFormat.jl")
include("model_SubResource.jl")
include("model_Route.jl")
include("model_RouteTable.jl")
include("model_Subnet.jl")


text:
using Oceananigans.BoundaryConditions: fill_halo_regions!

import Oceananigans.TimeSteppers: update_state!

"""
    update_state!(model::ShallowWaterModel)

Fill halo regions for `model.solution` and `model.tracers`.
"""
function update_state!(model::ShallowWaterModel)

    # Fill halos for velocities and tracers
    fill_halo_regions!(merge(model.solution, model.tracers),
                       model.architecture,
                       model.clock,
                       fields(model))

    return nothing
end


text:
module PlanarPoisson

import FinElt: VariationalProblem, EigenProblem, Mesh
import FinElt: TRIANGLE, LINE
using LinearAlgebra: rmul!

export add_bilin_form!, add_lin_functnl!
export barycentric, barycentric!, centroids_of
export grad_dot_grad!, func_times_func!, source_times_func!
export bdry_func_times_func!, bdry_source_times_func!
export shape_params

function add_bilin_form!(vp::VariationalProblem, name::String,
                         elm_mat!::Function, coef=1.0)
    if !(name in keys(vp.mesh.elmtype))
        error("$name: unknown physical name")
    end
    if !(elm_mat! in keys(BILIN2GEOMTYPE))
        error("$name: unknown element matrix routine")
    end
    if BILIN2GEOMTYPE[elm_mat!] != vp.mesh.elmtype[name]
        error("Element type does not match bilinear form")
    end
    push!(vp.bilin_form, (name, elm_mat!,coef))
    return
end

function add_lin_functnl!(vp::VariationalProblem, name::String,
                          elm_vec!::Function, f::Function)
    if !(name in keys(vp.mesh.elmtype))
        error("$name: unknown physical name")
    end
    if !(elm_vec! in keys(LIN2GEOMTYPE))
        error("$name: unknown element matrix routine")
    end
    if LIN2GEOMTYPE[elm_vec!] != vp.mesh.elmtype[name]
        error("Element type does not match bilinear form")
    end
    push!(vp.lin_functnl, (name, elm_vec!, f))
    return
end

function add_bilin_form!(ep::EigenProblem, name::String,
                         elm_mat!::Function, side::Symbol,
                         coef=1.0)
    if !(name in keys(ep.mesh.elmtype))
        error("$name: unknown physical name")
    end
    if !(elm_mat! in keys(BILIN2GEOMTYPE))
        error("$name: unknown element matrix routine")
    end
    if BILIN2GEOMTYPE[elm_mat!] != ep.mesh.elmtype[name]
        error("Element type does not match bilinear form")
    end
    if side == :LHS
        push!(ep.LHS_bilin_form, (name, elm_mat!,coef))
    elseif side == :RHS
        push!(ep.RHS_bilin_form, (name, elm_mat!,coef))
    else
        error("$side: illegal argument")
    end
    return
end

"""
b, area = barycentric(z)

For a planar triangle with vertices z_1, z_2, z_3,
compute vectors b1, b2, b3 such that if

x = lambda_1 z_1 + lambda_2 z_2 + lambda_3 z_3

then the barycentric coordinates are given by

lambda_p = 1/3 + b_p dot ( x - centroid )

Also compute the area of the triangle.
"""
function barycentric(z::Matrix)
    b = zeros(2, 3)
    b[1,1] = z[1,1] - z[1,3]
    b[2,1] = z[2,1] - z[2,3]
    b[1,2] = z[1,2] - z[1,3]
    b[2,2] = z[2,2] - z[2,3]
    # overwrite B[1:2,1:2] with its inverse transpose
    area = b[1,1] * b[2,2] - b[1,2] * b[2,1]
    b[1,1], b[2,2] =  b[2,2]/area,  b[1,1]/area
    b[1,2], b[2,1] = -b[2,1]/area, -b[1,2]/area

    b[1,3] = -b[1,1] - b[1,2]
    b[2,3] = -b[2,1] - b[2,2]
    area = abs(area) / 2
    return b, area
end

function centroid(z::Matrix)
    cntrd = zeros(2)
    cntrd[1] = ( z[1,1] + z[1,2] + z[1,3] ) / 3
    cntrd[2] = ( z[2,1] + z[2,2] + z[2,3] ) / 3
    return cntrd
end

function centroids_of(mesh::Mesh, name::String)
    noelms = count_elms(mesh, name)
    elm = mesh.elms_of[name]
    @assert mesh.elmtype[name] == TRIANGLE
    cntrd = zeros(2, noelms)
    z = zeros(2, 3)
    for k = 1:noelms
        for p = 1:3, i = 1:2
            z[i,p] = mesh.coord[i,elm[p,k]]
        end
#        println(z)
        cntrd[:,k] = centroid(z)
    end
    return cntrd
end

function grad_dot_grad!(A::Matrix,         # output 3x3
                        z::Matrix)         # input
    # A = element stiffness matrix
    b, area = barycentric(z)
    for q=1:3
        for p = q:3
            A[p,q] = area * ( b[1,p] * b[1,q] + b[2,p] * b[2,q] )
        end
    end
    for q = 2:3
        for p = 1:q-1
            A[p,q] = A[q,p]
        end
    end
    return
end

function grad_dot_grad!(A::Matrix,         # output 3x3
                        z::Matrix,         # input
                        coef::Float64      )
    # A = element stiffness matrix
    grad_dot_grad!(A, z)
    rmul!(A, coef)
    return
end

function grad_dot_grad!(A::Matrix,         # output 3x3
                        z::Matrix,         # input
                        coef::Function)
    # A = element stiffness matrix
    cntrd = centroid(z)
    grad_dot_grad!(A, z, coef(cntrd))
    return
end

function func_times_func!(A::Matrix,       # output 3x3
                          z::Matrix)       # input  2x3
    # A = element mass matrix
    b, area = barycentric(z)
    fill!(A, area/12)
    for p=1:3
        A[p,p] *= 2
    end
    return
end

function func_times_func!(A::Matrix,
                          z::Matrix,
                          coef::Float64)
    func_times_func!(A, z)
    rmul!(A, coef)
    return
end

function func_times_func!(A::Matrix,
                          z::Matrix,
                          coef::Function)
    cntrd = centroid(z)
    func_times_func!(A, z, coef(cntrd))
    return
end

function bdry_func_times_func!(A::Matrix,   # output 2x2
                               z::Matrix)   # input  2x2
    # A = element mass matrix for the line segment from z[:,1]
    # to z[:,2]
    h = hypot(z[1,2]-z[1,1], z[2,2]-z[2,1])
    A[1,1] = h / 3
    A[2,1] = h / 6
    A[1,2] = A[2,1]
    A[2,2] = A[1,1]
    return
end

function source_times_func!(v::Vector,      # output 3
                            z::Matrix,      # input  2x3
                            f::Function)
    # v = element load vector, using piecewise-linear interpolation
    # of f
    A = zeros(3, 3)
    func_times_func!(A, z)
    fz = zeros(3)
    for p = 1:3
        fz[p] = f(z[:,p])
    end
    v[:] = A * fz
    return
end

function bdry_source_times_func!(v::Vector,     # output 3
                                z::Matrix,      # input  2x2
                                g::Function)
    # v = element load vector, using piecewise-linear interpolation
    # of g
    A = zeros(2, 2)
    bdry_func_times_func!(A, z)
    gz = zeros(2)
    for p = 1:2
        gz[p] = g(z[:,p])
    end
    v[:] = A * gz
    return
end

"""
Returns the diameter h and shape regularity measure h^2/area
for the triangle with vertices z[:,1], z[:,2], z[:,3]
(assumed to lie in the plane, i.e., only the first two
components of each column are used).
"""
function shape_params(z::Matrix)
    b, area = barycentric(z)
    side = zeros(2,3)
    for p = 1:2
        side[p,1] = z[p,3] - z[p,2]
        side[p,2] = z[p,1] - z[p,3]
        side[p,3] = z[p,2] - z[p,1]
    end
    h = min(norm(side[:,1]), norm(side[:,2]), norm(side[:,3]))
    sr = h^2 / area
    return h, sr
end

function shape_params(mesh::Mesh)
    h  = Dict{String, Vector{Float64}}()
    sr = Dict{String, Vector{Float64}}()
    z = zeros(size(mesh.coord,1), 3)
    for name in keys(mesh.elms_of)
        if mesh.elmtype[name] != TRIANGLE
            continue
        end
        elms = mesh.elms_of[name]
        noelms = size(elms, 2)
        h[name]  = zeros(noelms)
        sr[name] = zeros(noelms)
        for k = 1:noelms
            for i = 1:size(z,1), p = 1:3
                z[i,p] = mesh.coord[i,elms[p,k]]
            end
            h[name][k], sr[name][k] = shape_params(z)
        end
    end
    return h, sr
end

const BILIN2GEOMTYPE = Dict(
    grad_dot_grad!          => TRIANGLE,
    func_times_func!        => TRIANGLE,
    source_times_func!      => TRIANGLE,
    bdry_func_times_func!   => LINE,
    bdry_source_times_func! => LINE)

const LIN2GEOMTYPE = Dict(
    source_times_func!      => TRIANGLE,
    bdry_source_times_func! => LINE)

end # module


text:
"""Structs for defining various terms in a spin Hamiltonian.
"""

abstract type AbstractInteraction end      # Subtype this for user-facing interfaces
abstract type AbstractInteractionCPU end   # Subtype this for actual internal CPU implementations
abstract type AbstractInteractionGPU end   # Subtype this for actual internal GPU implementations

struct QuadraticInteraction <: AbstractInteraction
    J     :: Mat3
    bond  :: Bond
    label :: String
end

struct ExternalField <: AbstractInteraction
    B :: Vec3
end

struct DipoleDipole <: AbstractInteraction
    extent   :: Int
    η        :: Float64
end

"""
    SiteInfo(site::Int, S, g)

Characterizes the degree of freedom located at a given `site` index with
 a spin magnitude `S` and g-tensor `g`. When provided to a `SpinSystem`,
 this information is automatically propagated to all symmetry-equivalent
 sites.
"""
struct SiteInfo
    site  :: Int      # Site index
    S     :: Float64  # Magnitude of the spin
    g     :: Mat3     # Spin g-tensor
end

# Helper constructors
SiteInfo(site::Int, S, g::Number) = SiteInfo(site, S, Mat3(g * I))
SiteInfo(site::Int, S) = SiteInfo(site, S, 2.0)

function Base.show(io::IO, mime::MIME"text/plain", int::QuadraticInteraction)
    b = repr(mime, int.bond)
    J = int.J
    if int.bond.i == int.bond.j && iszero(int.bond.n)  # Catch on-site anisotropies
        @assert J ≈ J'
        # Check if it is easy-axis or easy-plane
        λ, V = eigen(J)
        nonzero_λ = findall(x -> abs(x) > 1e-12, λ)
        if length(nonzero_λ) == 1
            i = nonzero_λ[1]
            dir = V[:, i]
            if count(<(0.0), dir) >= 2
                dir = -dir
            end
            name, D = λ[i] < 0 ? ("easy_axis", -λ[i]) : ("easy_plane", λ[i])
            @printf io "%s(%.4g, [%.4g, %.4g, %.4g], %d)" name D dir[1] dir[2] dir[3] int.bond.i
        else
            @printf io "single_ion_anisotropy([%.4g %.4g %.4g; %.4g %.4g %.4g; %.4g %.4g %.4g], %d)" J[1,1] J[1,2] J[1,3] J[2,1] J[2,2] J[2,3] J[3,1] J[3,2] J[3,3] int.bond.i
        end
    elseif J ≈ -J'                         # Catch purely DM interactions
        x = J[2, 3]
        y = J[3, 1]
        z = J[1, 2]
        @printf io "dm_interaction([%.4g, %.4g, %.4g], %s)" x y z b
    elseif diagm(fill(J[1,1], 3)) ≈ J      # Catch Heisenberg interactions
        @printf io "heisenberg(%.4g, %s)" J[1,1] b
    elseif diagm(diag(J)) ≈ J              # Catch diagonal interactions
        @printf io "exchange(diagm([%.4g, %.4g, %.4g]), %s)" J[1,1] J[2,2] J[3,3] b
    else                                   # Rest -- general exchange interactions
        @printf io "exchange([%.4g %.4g %.4g; %.4g %.4g %.4g; %.4g %.4g %.4g], %s)" J[1,1] J[1,2] J[1,3] J[2,1] J[2,2] J[2,3] J[3,1] J[3,2] J[3,3] b
        # TODO: Figure out how to reenable this depending on context:
        # @printf io "exchange([%.4f %.4f %.4f\n"   J[1,1] J[1,2] J[1,3]
        # @printf io "          %.4f %.4f %.4f\n"   J[2,1] J[2,2] J[2,3]
        # @printf io "          %.4f %.4f %.4f],\n" J[3,1] J[3,2] J[3,3]
        # @printf io "    %s)" b
    end
end

"""
    exchange(J, bond::Bond, label="Exchange")

Creates a quadratic interaction,

```math
    ∑_{⟨ij⟩} 𝐒_i^T J^{(ij)} 𝐒_j
```

where ``⟨ij⟩`` runs over all bonds (not doubly counted) that are symmetry
equivalent to `bond`. The ``3 × 3`` interaction matrix ``J^{(ij)}`` is the
covariant transformation of `J` appropriate for the bond ``⟨ij⟩``.
"""
function exchange(J, bond::Bond, label::String="Exchange")
    QuadraticInteraction(Mat3(J), bond, label)
end

"""
    heisenberg(J, bond::Bond, label::String="Heisen")

Creates a Heisenberg interaction
```math
    J ∑_{⟨ij⟩} 𝐒_i ⋅ 𝐒_j
```
where ``⟨ij⟩`` runs over all bonds symmetry equivalent to `bond`.
"""
heisenberg(J, bond::Bond, label::String="Heisen") = QuadraticInteraction(J*Mat3(I), bond, label)

"""
    dm_interaction(DMvec, bond::Bond, label::String="DMInt")

Creates a DM Interaction
```math
    ∑_{⟨ij⟩} 𝐃^{(ij)} ⋅ (𝐒_i × 𝐒_j)
```
where ``⟨ij⟩`` runs over all bonds symmetry equivalent to `bond`, and
``𝐃^{(ij)}`` is the covariant transformation of the DM pseudo-vector `DMvec`
appropriate for the bond ``⟨ij⟩``.
"""
function dm_interaction(DMvec, bond::Bond, label::String="DMInt")
    J = SA[     0.0   DMvec[3] -DMvec[2]
           -DMvec[3]       0.0  DMvec[1]
            DMvec[2] -DMvec[1]      0.0]
    QuadraticInteraction(J, bond, label)
end

"""
    single_ion_anisotropy(J, site, label="Anisotropy")

Creates a quadratic single-ion anisotropy,
```math
    ∑_i 𝐒_i^T J^{(i)} 𝐒_i
```
where ``i`` runs over all sublattices that are symmetry equivalent to `site`,
and ``J^{(i)}`` is the covariant transformation of the ``3 × 3`` anisotropy
matrix `J` appropriate for ``i``. Without loss of generality, we require that
`J` is symmetric.
"""
function single_ion_anisotropy(J, site::Int, label::String="Anisotropy")
    if !(J ≈ J')
        error("Single-ion anisotropy must be symmetric.")
    end
    QuadraticInteraction(Mat3(J), Bond(site, site, [0,0,0]), label)
end

"""
    easy_axis(D, n, site, label="EasyAxis")

Creates an easy axis anisotropy,
```math
    - D ∑_i (𝐧̂^{(i)}⋅𝐒_i)^2
```
where ``i`` runs over all sublattices that are symmetry equivalent to `site`,
``𝐧̂^{(i)}`` is the covariant transformation of the unit vector `n`, and ``D > 0``
is the interaction strength.
"""
function easy_axis(D, n, site::Int, label::String="EasyAxis")
    if D <= 0
        error("Parameter `D` must be nonnegative.")
    end
    if !(norm(n) ≈ 1)
        error("Parameter `n` must be a unit vector. Consider using `normalize(n)`.")
    end
    QuadraticInteraction(-D*Mat3(n*n'), Bond(site, site, [0,0,0]), label)
end

"""
    easy_plane(D, n, site, label="EasyPlane")

Creates an easy plane anisotropy,
```math
    + D ∑_i (𝐧̂^{(i)}⋅𝐒_i)^2
```
where ``i`` runs over all sublattices that are symmetry equivalent to `site`,
``𝐧̂^{(i)}`` is the covariant transformation of the unit vector `n`, and ``D > 0``
is the interaction strength.
"""
function easy_plane(D, n, site::Int, label::String="EasyAxis")
    if D <= 0
        error("Parameter `D` must be nonnegative.")
    end
    if !(norm(n) ≈ 1)
        error("Parameter `n` must be a unit vector. Consider using `normalize(n)`.")
    end
    QuadraticInteraction(+D*Mat3(n*n'), Bond(site, site, [0,0,0]), label)
end

"""
    external_field(B::Vec3)

Adds an external field ``𝐁`` with Zeeman coupling,

```math
    -∑_i 𝐁 ⋅ 𝐌_i.
```

The magnetic moments are ``𝐌_i = μ_B g 𝐒_i`` where ``g`` is the g-factor or
g-tensor, and the spin magnitude ``|𝐒_i|`` is typically a multiple of 1/2. The
Bohr magneton ``μ_B`` is a physical constant, with numerical value determined by
the unit system.
"""
external_field(B) = ExternalField(Vec3(B))

function Base.show(io::IO, ::MIME"text/plain", int::ExternalField)
    B = int.B
    @printf io "external_field([%.4g, %.4g, %.4g])" B[1] B[2] B[3]
end

"""
    dipole_dipole(; extent::Int=4, η::Float64=0.5)

Includes long-range dipole-dipole interactions,

```math
    -(μ₀/4π) ∑_{⟨ij⟩}  (3 (𝐌_j⋅𝐫̂_{ij})(𝐌_i⋅𝐫̂_{ij}) - 𝐌_i⋅𝐌_j) / |𝐫_{ij}|^3
```

where the sum is over all pairs of spins (singly counted), including periodic
images, regularized using the Ewald summation convention. The magnetic moments
are ``𝐌_i = μ_B g 𝐒_i`` where ``g`` is the g-factor or g-tensor, and the spin
magnitude ``|𝐒_i|`` is typically a multiple of 1/2. The Bohr magneton ``μ_B``
and vacuum permeability ``μ_0`` are physical constants, with numerical values
determined by the unit system.

`extent` controls the number of periodic copies of the unit cell summed over in
the Ewald summation (higher is more accurate, but higher creation-time cost),
while `η` controls the direct/reciprocal-space tradeoff in the Ewald summation.
"""
dipole_dipole(; extent=4, η=0.5) = DipoleDipole(extent, η)

#= Energy and field functions for "simple" interactions that aren't geometry-dependent.
   See Hamiltonian.jl for expectations on `_accum_neggrad!` functions.
=#

struct ExternalFieldCPU
    effBs :: Vector{Vec3}  # |S_b|gᵀB for each basis index b
end

function ExternalFieldCPU(ext_field::ExternalField, sites_info::Vector{SiteInfo}; μB=BOHR_MAGNETON)
    # As E = -∑_i 𝐁^T g 𝐒_i, we can precompute effB = g^T S B, so that
    #  we can compute E = -∑_i effB ⋅ 𝐬_i during simulation.
    # However, S_i may be basis-dependent, so we need to store an effB
    #  per sublattice.
    effBs = [μB * site.g' * site.S * ext_field.B for site in sites_info]
    ExternalFieldCPU(effBs)
end

function energy(spins::Array{Vec3, 4}, field::ExternalFieldCPU)
    E = 0.0
    for b in 1:size(spins, 1)
        effB = field.effBs[b]
        for s in selectdim(spins, 1, b)
            E += effB ⋅ s
        end
    end
    return -E
end

"Accumulates the negative local Hamiltonian gradient coming from the external field"
@inline function _accum_neggrad!(B::Array{Vec3, 4}, field::ExternalFieldCPU)
    for b in 1:size(B, 1)
        effB = field.effBs[b]
        for idx in CartesianIndices(size(B)[2:end])
            B[b, idx] = B[b, idx] + effB
        end
    end
end


text:
include("World.jl")
using Random
function asd()
rng = MersenneTwister(0);
adjacency_list = Vector{Union{Vector{Union{Nothing,World.GraphNode}},Nothing}}(nothing,1)

first_node = World.GraphNode(0,0,World.Terrain(1))
neighboring_nodes = Vector{Union{Nothing,World.GraphNode}}(nothing,1)
for xpos in 1:2
    for ypos in 1:2
        new_node = World.GraphNode(xpos,ypos,World.Terrain(1))
        if (neighboring_nodes[end] == nothing)
            neighboring_nodes[end] = new_node
        else
            push!(neighboring_nodes,new_node)

        end
    end
end
first_list = Vector{Union{Nothing,World.GraphNode}}(nothing,1)
first_list[1] = pop!(neighboring_nodes)
adjacency_list[1] = first_list
for node in neighboring_nodes
    push!(adjacency_list[1],node)
end
number_graph_nodes = length(adjacency_list[1])
am = zeros(Int8,number_graph_nodes,number_graph_nodes)
i::Int = 1
j::Int = 1
for lists in adjacency_list
    if(j > 0) am[j,1] = 1
    for nodes in lists
        if ( i != j)
            am[i,j] = 1
        end
        println(i)
        println(j)
        j += 1
    end
    i += 1
    j = 1
end
print(am)
end

asd()

max_stars_repo_name: VesaIhanainen/node-world

max_issues_repo_name: VesaIhanainen/node-world

max_forks_repo_name: VesaIhanainen/node-world


text:
################################################################################

struct ActiveDie
	die::Die.Face
	use_as::Die.Face
	use_modt::Bool
	use_ring::Bool

	used::Bool
end

ActiveDie(die, use_as, use_modt, use_ring) = ActiveDie(die, use_as, use_modt, use_ring, false)
use_active_die(d) = ActiveDie(d.die, d.use_as, d.use_modt, d.use_ring, true)

article(d::ActiveDie) = Die.article(d.die)

function Base.string(d::ActiveDie)
	d.die == d.use_as && return string(d.die)

	text = "$(string(d.die)) as $(string(d.use_as))"
	d.use_modt && return text*" by using the Messenger of the Dark Tower ability"
	d.use_ring && return text*" by using an Elven Ring"
	return text
end

function get_active_die(die, available_dice, strategy, may_use_modt, may_use_ring)
	die in available_dice && return ActiveDie(die,die,false,false)

	(Die.ArmyMuster in available_dice) && (die == Die.Army || die == Die.Muster) && return ActiveDie(Die.ArmyMuster,die,false,false)
	(Die.Muster in available_dice) && (die == Die.Army && may_use_modt) && return ActiveDie(Die.Muster,die,true,false)

	strategy == Strategy.Military && (discardable_dice = [Die.Character, Die.Event])
	strategy == Strategy.Corruption && (discardable_dice = [Die.Army, Die.Muster, Die.ArmyMuster, Die.Event])

	discardable_dice = intersect(available_dice, discardable_dice)
	may_use_ring && !isempty(discardable_dice) && return ActiveDie(rand(discardable_dice),die,false,true)

	return nothing
end

################################################################################

mutable struct QuellerState
	strategy::Strategy.Choice

	active_die::Union{ActiveDie,Nothing}
	available_dice::Vector{Die.Face}

	ring_available::Bool
	modt_available::Bool
end

function QuellerState(available_dice = Vector{Die.Face}())
	strategy = rand(instances(Strategy.Choice))

	active_die = nothing

	ring_available = false
	modt_available = false

	QuellerState(strategy,active_die,available_dice,ring_available,modt_available)
end

################################################################################

mutable struct CheckStrategy <: StateInteractionNode
	id::String
	strategy::Strategy.Choice
	n_true::Node
	n_false::Node

	CheckStrategy(strategy) = (obj = new(); obj.strategy = StrategyChoice(strategy); obj)
end
getnext!(n::CheckStrategy, state) = (state.strategy == n.strategy ? n.n_true : n.n_false)
setnext!(n::CheckStrategy, name::Symbol, next::Node) = setfield!(n,name,next)

getmsg(n::CheckStrategy) = "The $(string(n.strategy)) strategy is used."
getmsg(n::CheckStrategy, state) = "$(getmsg(n)): $(state.strategy == n.strategy)"

################################################################################

mutable struct SetStrategy <: StateInteractionNode
	id::String
	strategy::Strategy.Choice
	next::Node

	SetStrategy(strategy) = (obj = new(); obj.strategy = StrategyChoice(strategy); obj)
end
setnext!(n::SetStrategy, next::Node) = (n.next = next)
getnext!(n::SetStrategy, state) = (state.strategy = n.strategy; n.next)

getmsg(n::SetStrategy) = "Set Queller strategy to $(string(n.strategy))."

################################################################################

mutable struct CheckActiveDie <: StateInteractionNode
	id::String
	die::Die.Face
	n_true::Node
	n_false::Node

	CheckActiveDie(die) = (obj = new(); obj.die = DieFace(die); obj)
end
setnext!(n::CheckActiveDie, name::Symbol, next::Node) = setfield!(n,name,next)
getnext!(n::CheckActiveDie, state) = (state.active_die.use_as == n.die ? n.n_true : n.n_false)
getmsg(n::CheckActiveDie) = "The active die is $(Die.article(n.die)) $(string(n.die))"
getmsg(n::CheckActiveDie, state) = "$(getmsg(n)): $(state.active_die.use_as == n.die)"

################################################################################

mutable struct SetActiveDie <: StateInteractionNode
	id::String
	next::Node
	no_die::Node

	die::Die.Face
	may_use_ring::Bool

	SetActiveDie(die; may_use_ring=false) = (obj = new(); obj.die = DieFace(die); obj.may_use_ring = may_use_ring; obj)
end
setnext!(n::SetActiveDie, name::Symbol, next::Node) = setfield!(n,name,next)

function get_active_die(n::SetActiveDie,state)
	!isnothing(state.active_die) && state.active_die.used && error("Trying to set an active die when a die already have been used.")
	return get_active_die(n.die, state.available_dice, state.strategy, state.modt_available, state.ring_available && n.may_use_ring)
end

function getnext!(n::SetActiveDie,state)
	state.active_die = get_active_die(n, state)
	return isnothing(state.active_die) ? n.no_die : n.next
end

getmsg(n::SetActiveDie, state) = "The active die is set to: $(string(get_active_die(n,state)))"

function getmsg(n::SetActiveDie)
	prio_list = [string(n.die)]
	if n.die == Die.Army
		push!(prio_list, "$(string(Die.ArmyMuster)) as $(string(Die.Army))")
		push!(prio_list, "$(string(Die.Muster)) and Messenger of the Dark Tower as $(string(Die.Army))")
	elseif n.die == Die.Muster
		push!(prio_list, "$(string(Die.ArmyMuster)) as a $(string(Die.Muster))")
	end

	n.may_use_ring && push!(prio_list, "A random non-*preferred* die and an Elven Ring as $(string(n.die))")

	text = """
	Set the first matching die as the Active Die:

	"""
	for (i, d) in enumerate(prio_list)
		text *= "$(i). $(d)\n"
	end
	return text
end

################################################################################

mutable struct UseActiveDie <: StateInteractionNode
	id::String
	next::Node

	UseActiveDie() = new()
end
setnext!(n::UseActiveDie, next::Node) = (n.next = next)
function getnext!(n::UseActiveDie, state)
	isnothing(state.active_die) && error("No active die have been set before it is used.")

	if !state.active_die.used
		state.active_die = use_active_die(state.active_die)
		i = findfirst(==(state.active_die.die), state.available_dice)
		state.available_dice = deleteat!(state.available_dice, i)
	end
	return n.next
end

getmsg(n::UseActiveDie) = "Use the Active Die to:"
function getmsg(n::UseActiveDie, state)
	die = state.active_die
	return "Use $(article(die)) $(string(die)) to:"
end

################################################################################

mutable struct SetRingAvailable <: StateInteractionNode
	id::String
	next::Node
	ring_available::Bool

	SetRingAvailable(ring_available) = (obj = new(); obj.ring_available = ring_available; obj)
end
setnext!(n::SetRingAvailable, next::Node) = (n.next = next)
getnext!(n::SetRingAvailable, state) = (state.ring_available = n.ring_available; n.next)
getmsg(n::SetRingAvailable) = n.ring_available ? "Set an Elven ring as available." : "Set an Elven ring as not available."

################################################################################

mutable struct SetMoDTAvailable <: StateInteractionNode
	id::String
	next::Node
	modt_available::Bool

	SetMoDTAvailable(modt_available) = (obj = new(); obj.modt_available = modt_available; obj)
end
setnext!(n::SetMoDTAvailable, next::Node) = (n.next = next)
getnext!(n::SetMoDTAvailable, state) = (state.modt_available = n.modt_available; n.next)
getmsg(n::SetMoDTAvailable) = n.modt_available ? "Set 'Messenger of the Dark Tower' as available." : "Set 'Messenger of the Dark Tower' as not available."

################################################################################

mutable struct GetAvailableDice <: InteractiveNode
	id::String
	next::Node
	action::String

	GetAvailableDice(action = """
					 Input the available action dice (not counting dice set aside for later use).
					 """) = (obj = new(); obj.action = strip(action); obj)
end
setnext!(n::GetAvailableDice, next::Node) = (n.next = next)
getopt(n::GetAvailableDice) = [CMD.Dice()]
getnext!(n::GetAvailableDice, state, opt::CMD.Dice) = (state.available_dice = opt.dice; n.next)

function getmsg(n::GetAvailableDice)
	diefaces = instances(Die.Face)
	legends = Die.char.(diefaces).*": ".*string.(diefaces)

	return """
	$(n.action)

	$(strvec2str(legends))
	"""
end

max_stars_repo_name: mvmorin/queller-bot

max_issues_repo_name: mvmorin/queller-bot

max_forks_repo_name: mvmorin/queller-bot


text:
# Author: Ali Siahkoohi, alisk@gatech.edu
# Date: September 2020
# Copyright: Georgia Institute of Technology, 2020

export logpdf, gradlogpdf

function logpdf(μ::Float32, σ::Float32, X)

    f = -.5f0*((X .- μ)/σ).^2
    f = f .+ -5f-1 * log(2f0π) .- log(σ)

    return f
end

function gradlogpdf(μ::Float32, σ::Float32, X)

    g = -(X .- μ)/σ^2

    return g
end


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

@testset "graph.triangle" begin

@testset "triangles" begin
    graph1 = pathgraph(5)
    @test isempty(triangles(graph1))
    graph2 = plaingraph(5, [(1, 2), (2, 3), (3, 1)])
    @test issetequal(collect(triangles(graph2))[1], 1:3)
    graph3 = plaingraph(8, [
        (1, 2), (2, 3), (1, 3), (3, 4), (4, 5),
        (5, 6), (4, 6), (5, 7), (7, 8), (8, 6)
    ])
    @test length(triangles(graph3)) == 2
    graph4 = plaingraph(10, [
        (1, 2), (2, 3), (3, 4), (3, 5), (5, 6),
        (6, 7), (7, 5), (5, 8), (8, 9), (9, 5), (5, 10)
    ])
    @test length(triangles(graph4)) == 2
end

end # graph.triangle


text:
using Flux
using Base: @get!
using MacroTools: @forward

const ϵ = 1e-8

# TODO: should use weak refs

"""
  Descent(η)

Classic gradient descent optimiser with learning rate `η`.
For each parameter `p` and its gradient `δp`, this runs `p -= η*δp`

## Parameters
  - Learning Rate (η): The amount by which the gradients are discounted before updating the weights. Defaults to `0.1`.

## Example
```julia-repl
opt = Descent() # uses default η (0.1)

opt = Descent(0.3) # use provided η

ps = params(model)

gs = gradient(ps) do
  loss(x, y)
end

Flux.Optimise.update!(opt, ps, gs)
```
"""
mutable struct Descent
  eta::Float64
end

Descent() = Descent(0.1)

function apply!(o::Descent, x, Δ)
  Δ .*= o.eta
end

"""
    Momentum(η, ρ)

Gradient descent with learning rate `η` and momentum `ρ`.

## Parameters
  - Learning Rate (`η`): Amount by which gradients are discounted before updating the weights. Defaults to `0.01`.
  - Momentum (`ρ`): Parameter that accelerates descent in the relevant direction and dampens oscillations. Defaults to `0.9`.

## Examples
```julia
opt = Momentum() # uses defaults of η = 0.01 and ρ = 0.9

opt = Momentum(0.01, 0.99)
```
"""
mutable struct Momentum
  eta::Float64
  rho::Float64
  velocity::IdDict
end

Momentum(η = 0.01, ρ = 0.9) = Momentum(η, ρ, IdDict())

function apply!(o::Momentum, x, Δ)
  η, ρ = o.eta, o.rho
  v = get!(o.velocity, x, zero(x))::typeof(x)
  @. v = ρ * v - η * Δ
  @. Δ = -v
end

"""
    Nesterov(η, ρ)

Gradient descent with learning rate  `η` and Nesterov momentum `ρ`.

## Parameters
  - Learning Rate (η): Amount by which the gradients are dicsounted berfore updating the weights. Defaults to `0.001`.
  - Nesterov Momentum (ρ): Paramters controlling the amount of nesterov momentum to be applied. Defaults to `0.9`.

## Examples
```julia
opt = Nesterov() # uses defaults η = 0.001 and ρ = 0.9

opt = Nesterov(0.003, 0.95)
```
"""
mutable struct Nesterov
  eta::Float64
  rho::Float64
  velocity::IdDict
end

Nesterov(η = 0.001, ρ = 0.9) = Nesterov(η, ρ, IdDict())

function apply!(o::Nesterov, x, Δ)
  η, ρ = o.eta, o.rho
  v = get!(o.velocity, x, zero(x))::typeof(x)
  d = @. ρ^2 * v - (1+ρ) * η * Δ
  @. v = ρ*v - η*Δ
  @. Δ = -d
end

"""
    RMSProp(η, ρ)

Implements the RMSProp algortihm. Often a good choice for recurrent networks. Paramters other than learning rate generally don't need tuning.

## Parameters
  - Learning Rate (η): Defaults to `0.001`.
  - Rho (ρ): Defaults to `0.9`.

## Examples
```julia
opt = RMSProp() # uses default η = 0.001 and ρ = 0.9

opt = RMSProp(0.002, 0.95)
```

## References
[RMSProp](https://www.cs.toronto.edu/~tijmen/csc321/slides/lecture_slides_lec6.pdf)
"""
mutable struct RMSProp
  eta::Float64
  rho::Float64
  acc::IdDict
end

RMSProp(η = 0.001, ρ = 0.9) = RMSProp(η, ρ, IdDict())

function apply!(o::RMSProp, x, Δ)
  η, ρ = o.eta, o.rho
  acc = get!(o.acc, x, zero(x))::typeof(x)
  @. acc = ρ * acc + (1 - ρ) * Δ^2
  @. Δ *= η / (√acc + ϵ)
end

"""
    ADAM(η, β::Tuple)

Implements the ADAM optimiser.

## Paramters
  - Learning Rate (`η`): Defaults to `0.001`.
  - Beta (`β::Tuple`): The first element refers to β1 and the second to β2. Defaults to `(0.9, 0.999)`.

## Examples

```julia
opt = ADAM() # uses the default η = 0.001 and β = (0.9, 0.999)

opt = ADAM(0.001, (0.9, 0.8))
```
## References
[ADAM](https://arxiv.org/abs/1412.6980v8) optimiser.
"""
mutable struct ADAM
  eta::Float64
  beta::Tuple{Float64,Float64}
  state::IdDict
end

ADAM(η = 0.001, β = (0.9, 0.999)) = ADAM(η, β, IdDict())

function apply!(o::ADAM, x, Δ)
  η, β = o.eta, o.beta
  mt, vt, βp = get!(o.state, x, (zero(x), zero(x), β))
  @. mt = β[1] * mt + (1 - β[1]) * Δ
  @. vt = β[2] * vt + (1 - β[2]) * Δ^2
  @. Δ =  mt / (1 - βp[1]) / (√(vt / (1 - βp[2])) + ϵ) * η
  o.state[x] = (mt, vt, βp .* β)
  return Δ
end

"""
    RADAM(η, β::Tuple)

Implements the rectified ADAM optimizer.

## Parameters
  - Learning Rate (η): Defaults to `0.001`
  - Beta (β::Tuple): The first element refers to β1 and the second to β2. Defaults to `(0.9, 0.999)`.

## Examples

```julia
opt = RADAM() # uses the default η = 0.001 and β = (0.9, 0.999)

opt = RADAM(0.001, (0.9, 0.8))
```

## References
[RADAM](https://arxiv.org/pdf/1908.03265v1.pdf) optimiser (Rectified ADAM).
"""
mutable struct RADAM
  eta::Float64
  beta::Tuple{Float64,Float64}
  state::IdDict
end

RADAM(η = 0.001, β = (0.9, 0.999)) = RADAM(η, β, IdDict())

function apply!(o::RADAM, x, Δ)
  η, β = o.eta, o.beta
  ρ∞ = 2/(1-β[2])-1
  mt, vt, βp, t = get!(o.state, x, (zero(x), zero(x), β, 1))
  @. mt = β[1] * mt + (1 - β[1]) * Δ
  @. vt = β[2] * vt + (1 - β[2]) * Δ^2
  ρ = ρ∞ - 2t*βp[2]/(1-βp[2])
  if ρ > 4
    r = sqrt((ρ-4)*(ρ-2)*ρ∞/((ρ∞-4)*(ρ∞-2)*ρ))
    @. Δ =  mt / (1 - βp[1]) / (√(vt / (1 - βp[2])) + ϵ) * η * r
  else
    @. Δ =  mt / (1 - βp[1]) * η
  end
  o.state[x] = (mt, vt, βp .* β, t+1)
  return Δ
end

"""
    AdaMax(η, β::Tuple)

Variant of ADAM based on ∞-norm.

## Parameters
  - Learning Rate (η): Defaults to `0.001`
  - Beta (β::Tuple): The first element refers to β1 and the second to β2. Defaults to `(0.9, 0.999)`.

## Examples
```julia
opt = AdaMax() # uses default η and β

opt = AdaMax(0.001, (0.9, 0.995))
```
## References
[AdaMax](https://arxiv.org/abs/1412.6980v9) optimiser.
"""
mutable struct AdaMax
  eta::Float64
  beta::Tuple{Float64,Float64}
  state::IdDict
end

AdaMax(η = 0.001, β = (0.9, 0.999)) = AdaMax(η, β, IdDict())

function apply!(o::AdaMax, x, Δ)
  η, β = o.eta, o.beta
  mt, ut, βp = get!(o.state, x, (zero(x), zero(x), β))
  @. mt = β[1] * mt + (1 - β[1]) * Δ
  @. ut = max(β[2] * ut, abs(Δ))
  @. Δ = (η/(1 - βp[1])) * mt/(ut + ϵ)
  o.state[x] = (mt, ut, βp .* β)
  return Δ
end

"""
    ADAGrad(η)

Implements AdaGrad. It has parameter specific learning rates based on how frequently it is updated.

## Parameters
  - Learning Rate (η): Defaults to `0.1`

## Examples
```julia
opt = ADAGrad() # uses default η = 0.1

opt = ADAGrad(0.001)
```

## References
[ADAGrad](http://www.jmlr.org/papers/volume12/duchi11a/duchi11a.pdf) optimiser.
Parameters don't need tuning.
"""
mutable struct ADAGrad
  eta::Float64
  acc::IdDict
end

ADAGrad(η = 0.1) = ADAGrad(η, IdDict())

function apply!(o::ADAGrad, x, Δ)
  η = o.eta
  acc = get!(o.acc, x, fill!(zero(x), ϵ))::typeof(x)
  @. acc += Δ^2
  @. Δ *= η / (√acc + ϵ)
end

"""
    ADADelta(ρ)

Version of ADAGrad that adapts learning rate based on a window of past gradient updates. Parameters don't need tuning.

## Parameters
  - Rho (ρ): Factor by which gradient is decayed at each time step. Defaults to `0.9`.

## Examples
```julia
opt = ADADelta() # uses default ρ = 0.9
opt = ADADelta(0.89)
```

## References
[ADADelta](https://arxiv.org/abs/1212.5701) optimiser.
"""
mutable struct ADADelta
  rho::Float64
  state::IdDict
end

ADADelta(ρ = 0.9) = ADADelta(ρ, IdDict())

function apply!(o::ADADelta, x, Δ)
  ρ = o.rho
  acc, Δacc = get!(o.state, x, (zero(x), zero(x)))
  @. acc = ρ * acc + (1 - ρ) * Δ^2
  @. Δ *= √Δacc/ (√acc + ϵ)
  @. Δacc = ρ * Δacc + (1 - ρ) * Δ^2
  return Δ
end

"""
    AMSGrad(η, β::Tuple)

Implements AMSGrad version of the ADAM optimiser. Parameters don't need tuning.

## Parameters
  - Learning Rate (η): Defaults to `0.001`.
  - Beta (β::Tuple): The first element refers to β1 and the second to β2. Defaults to `(0.9, 0.999)`.

## Examples
```julia
opt = AMSGrad() # uses default η and β
opt = AMSGrad(0.001, (0.89, 0.995))
```

## References
[AMSGrad](https://openreview.net/forum?id=ryQu7f-RZ) optimiser.
"""
mutable struct AMSGrad
  eta::Float64
  beta::Tuple{Float64, Float64}
  state::IdDict
end

AMSGrad(η = 0.001, β = (0.9, 0.999)) = AMSGrad(η, β, IdDict())

function apply!(o::AMSGrad, x, Δ)
  η, β = o.eta, o.beta
  mt, vt, v̂t = get!(o.state, x, (fill!(zero(x), ϵ), fill!(zero(x), ϵ), fill!(zero(x), ϵ)))
  @. mt = β[1] * mt + (1 - β[1]) * Δ
  @. vt = β[2] * vt + (1 - β[2]) * Δ ^ 2
  @. v̂t = max(v̂t, vt)
  @. Δ = η * mt / (√v̂t + ϵ)
end

"""
    NADAM(η, β::Tuple)

Nesterov variant of ADAM. Parameters don't need tuning.

## Parameters
  - Learning Rate (η): Defaults to `0.001`.
  - Beta (β::Tuple): The first element refers to β1 and the second to β2. Defaults to `(0.9, 0.999)`.

## Examples
```julia
opt = NADAM() # uses default η and β
opt = NADAM(0.002, (0.89, 0.995))
```

## References
[NADAM](http://cs229.stanford.edu/proj2015/054_report.pdf) optimiser.
"""
mutable struct NADAM
  eta::Float64
  beta::Tuple{Float64, Float64}
  state::IdDict
end

NADAM(η = 0.001, β = (0.9, 0.999)) = NADAM(η, β, IdDict())

function apply!(o::NADAM, x, Δ)
  η, β = o.eta, o.beta
  mt, vt, (β1p, β2p) = get!(o.state, x, (zero(x), zero(x), o.beta))
  @. mt = β[1] * mt + (1 - β[1]) * Δ
  @. vt = β[2] * vt + (1 - β[2]) * Δ^2
  @. Δ = (β[1] * mt / (1 - β[1] * β1p) + (1 - β[1]) * Δ / (1 - β1p)) / (√(vt * β[2] / (1 - β2p)) + ϵ) * η
  o.state[x] = (mt, vt, (β1p * β[1], β2p * β[2]))
  return Δ
end

"""
    ADAMW(η, β::Tuple, decay)

Variant of ADAM defined by fixing weight decay regularization.

## Parameters
  - Learning Rate (η): Defaults to `0.001`.
  - Beta (β::Tuple): The first element refers to β1 and the second to β2. Defaults to (0.9, 0.999).
  - decay: Decay applied to weights during optimisation. Defaults to 0.

## Examples
```julia
opt = ADAMW() # uses default η, β and decay
opt = ADAMW(0.001, (0.89, 0.995), 0.1)
```

## References
[ADAMW](https://arxiv.org/abs/1711.05101)
"""
ADAMW(η = 0.001, β = (0.9, 0.999), decay = 0) =
  Optimiser(ADAM(η, β), WeightDecay(decay))

# Compose optimizers

"""
    Optimiser(a, b, c...)

Combine several optimisers into one; each optimiser produces a modified gradient
that will be fed into the next, and this is finally applied to the parameter as
usual.
"""
mutable struct Optimiser
  os::Vector{Any}
end

Optimiser(o...) = Optimiser(Any[o...])

@forward Optimiser.os Base.getindex, Base.first, Base.last, Base.lastindex, Base.push!, Base.setindex!
@forward Optimiser.os Base.iterate

Base.getindex(c::Optimiser, i::AbstractArray) = Optimiser(c.os[i]...)

function apply!(o::Optimiser, x, Δ)
  for opt in o.os
    Δ = apply!(opt, x, Δ)
  end
  return Δ
end

"""
  InvDecay(γ)

Applies inverse time decay to an optimiser

## Parameters
  - gamma (γ): Defaults to `0.001`

## Example
```julia
  Optimiser(InvDecay(..), Opt(..))
```
"""
mutable struct InvDecay
  gamma::Float64
  state::IdDict
end

InvDecay(γ = 0.001) = InvDecay(γ, IdDict())

function apply!(o::InvDecay, x, Δ)
  γ = o.gamma
  n = get!(o.state, x, 1)
  Δ .*= 1 / (1 + γ * n)
  o.state[x] = n + 1
  return Δ
end

"""
  ExpDecay(eta, decay, decay_step, clip)

Discount the learning rate `eta` by `decay` every `decay_step` till a minimum of `clip`.

## Parameters
  - Learning Rate (eta): Defaults to `0.001`.
  - decay: Factor by which the learning rate is discounted. Defaults to `0.1`.
  - decay_step: Schedules decay operations by setting number of steps between two decay operations. Defaults to `1000`.
  - clip: Minimum value of learning rate. Defaults to `1e-4`.

## Example
To apply exponential decay to an optimiser:
```julia
  Optimiser(ExpDecay(..), Opt(..))

  opt = Optimiser(ExpDecay(), ADAM())
```
"""
mutable struct ExpDecay
  eta::Float64
  decay::Float64
  step::Int64
  clip::Float64
  current::IdDict
end

ExpDecay(opt = 0.001, decay = 0.1, decay_step = 1000, clip = 1e-4) = ExpDecay(opt, decay, decay_step, clip, IdDict())

function apply!(o::ExpDecay, x, Δ)
  η, s, decay = o.eta, o.step, o.decay
  n = o.current[x] = get(o.current, x, 0) + 1
  if o.current[x]%s == 0 && count(x -> x%s == 0, values(o.current)) == 1
    η = max(η * decay^(s / n), o.clip)
    o.eta = η
  end
  @. Δ *= η
end

"""
  WeightDecay(wd)

Decays the weight by `wd`

## Parameters
  - weight decay (wd): 0
"""
mutable struct WeightDecay
  wd::Real
end

WeightDecay() = WeightDecay(0)

function apply!(o::WeightDecay, x, Δ)
  wd = o.wd
  @. Δ += wd * x
end


text:
# ----------------------- #
# Display-esque functions #
# ----------------------- #
Base.istextmime(::MIME"application/vnd.plotly.v1+json") = true

function html_body(p::Plot)
    """

    """
end

function script_content(p::Plot)
    lowered = JSON.lower(p)
    """
    Plotly.newPlot('$(p.divid)', $(json(lowered[:data])),
                   $(json(lowered[:layout])), {showLink: false});
    """
end

const js_default = Ref(:local)

function stringmime(::MIME"text/html", p::Plot, js::Symbol=js_default[])

    if js == :local
        script_txt = ""
    elseif js == :remote
        script_txt = ""
    elseif js == :embed
        script_txt = ""
    else
        msg = """
        Unknown value for argument js: $js.
        Possible choices are `:local`, `:remote`, `:embed`
            """
        throw(ArgumentError(msg))
    end

    """

         $script_txt

         $(html_body(p))

    """

end

Base.show(io::IO, ::MIME"text/html", p::Plot, js::Symbol=js_default[]) =
    print(io, stringmime(MIME"text/html"(), p, js))

function Base.show(io::IO, ::MIME"text/plain", p::Plot)
    println(io, """
    data: $(json(map(_describe, p.data), 2))
    layout: "$(_describe(p.layout))"
    """)
end

Base.show(io::IO, p::Plot) = show(io, MIME("text/plain"), p)

# ----------------------------------------- #
# SyncPlot -- sync Plot object with display #
# ----------------------------------------- #
struct SyncPlot{TD<:AbstractPlotlyDisplay}
    plot::Plot
    view::TD
end

plot(args...; kwargs...) = SyncPlot(Plot(args...; kwargs...))

## API methods for SyncPlot
for f in [:restyle!, :relayout!, :addtraces!, :deletetraces!, :movetraces!,
          :redraw!, :extendtraces!, :prependtraces!, :purge!, :to_image,
          :download_image]
    @eval function $(f)(sp::SyncPlot, args...; kwargs...)
        $(f)(sp.plot, args...; kwargs...)
        $(f)(sp.view, args...; kwargs...)
    end

    no_!_method = Symbol(string(f)[1:end-1])
    @eval function $(no_!_method)(sp::SyncPlot, args...; kwargs...)
        sp2 = fork(sp)
        $f(sp2.plot, args...; kwargs...)  # only need to update the julia side
        sp2  # return so we display fresh
    end
end

# add extra same extra methods we have on ::Plot for these functions
for f in (:extendtraces!, :prependtraces!)
    @eval begin
        function $(f)(p::AbstractPlotlyDisplay, inds::Vector{Int}=[0],
                      maxpoints=-1; update...)
            ($f)(p, Dict(map(x->(x[1], _tovec(x[2])), update)), inds, maxpoints)
        end

        function $(f)(p::AbstractPlotlyDisplay, ind::Int, maxpoints=-1;
                      update...)
            ($f)(p, [ind], maxpoints; update...)
        end

        function $(f)(p::AbstractPlotlyDisplay, update::Associative, ind::Int,
                      maxpoints=-1)
            ($f)(p, update, [ind], maxpoints)
        end
    end
end

Base.show(io::IO, ::MIME"text/html", sp::SyncPlot, js::Symbol=js_default[]) =
    print(io, stringmime(MIME"text/html"(), sp.plot, js))

function Base.show(io::IO, ::MIME"application/vnd.plotly.v1+json", p::Union{Plot,SyncPlot})
    JSON.print(io, p)
end

# Add some basic Julia API methods on SyncPlot that just forward onto the Plot
Base.size(sp::SyncPlot) = size(sp.plot)
Base.copy(sp::SyncPlot) = fork(sp)  # defined by each SyncPlot{TD}

# ----------------- #
# Display frontends #
# ----------------- #

include("displays/juno.jl")
include("displays/electron.jl")
include("displays/ijulia.jl")

# methods to convert from one frontend to another
let
    all_frontends = [:ElectronPlot, :JupyterPlot]
    for fe_to in all_frontends
        for fe_from in all_frontends
            @eval $(fe_to)(sp::$(fe_from)) = $(fe_to)(sp.plot)
        end
    end
end

# -------- #
# Defaults #
# -------- #


text:
cd(@__DIR__) do
    if !ispath("call_FGb")
        run(`wget http://www-polsys.lip6.fr/\~jcf/FGb/C/@downloads/call_FGb6.maclinux.x64.tar.gz`)
        run(`tar -xzf call_FGb6.maclinux.x64.tar.gz`)
    end

    Cdir = joinpath("call_FGb","nv","maple","C")

    BLOBSDIR = is_linux() ? "x64" : is_apple() ? "macosx" : assert(false)
    LIBDIR = joinpath(Cdir, BLOBSDIR)

    run(`gcc -m64 -c -fPIC dummy-hooks.c`)
    run(`g++ -m64 -shared -o libfgb.so dummy-hooks.o -L$(LIBDIR) -Wl,-zmuldefs -Wl,--whole-archive -lfgb -lfgbexp -lgb -lgbexp -lminpoly -lminpolyvgf -Wl,--no-whole-archive -lgmp -lm -fopenmp`)
end
if nothing == try Pkg.installed("PolynomialRings") end
    Pkg.clone("https://github.com/tkluck/PolynomialRings.jl")
end


text:
# Portions translated from SLICOT-Reference distribution
# Copyright (c) 2002-2020 NICONET e.V.
function run_mb04dz(datfile, io=stdout)
    NIN = 5
    NOUT = 6
    NMAX = 100
    LDA = NMAX
    LDQG = NMAX
    A = Array{ComplexF64,2}(undef, LDA,NMAX)
    QG = Array{ComplexF64,2}(undef, LDQG,NMAX+1)
    SCALE = Array{Float64,1}(undef, NMAX)
    DUMMY = Array{Float64,1}(undef, 1)
    f = open(datfile,"r")
    readline(f)
    vs = split(readline(f))
    N = parse(BlasInt, vs[1])
    JOB =  vs[2][1]
    if ( N<=0 || N>NMAX )
        @error "Illegal N=$N"
    end

    vs = String[]
    _isz,_jsz = (N,N)
    while length(vs) < _isz*_jsz
        append!(vs, replace.(split(readline(f)),'D'=>'E'))
    end
    for i in 1:_isz
       _i0 = (i-1)*_jsz
       A[i,1:_jsz] .= parsex.(ComplexF64, vs[_i0+1:_i0+_jsz])
    end

    vs = String[]
    _isz,_jsz = (N,N+1)
    while length(vs) < _isz*_jsz
        append!(vs, replace.(split(readline(f)),'D'=>'E'))
    end
    for i in 1:_isz
       _i0 = (i-1)*_jsz
       QG[i,1:_jsz] .= parsex.(ComplexF64, vs[_i0+1:_i0+_jsz])
    end
    close(f)
    # interp call 1

    ILO, INFO = SLICOT.mb04dz!(JOB, N, A, QG, SCALE)
    @test INFO == 0
    INFO == 0 || return
    println(io, "ILO = $ILO")

# interp output 1
    println(io, "A:")
    _nc = N
    _nr = N
    show(io, "text/plain", A[1:_nr,1:_nc])
    println(io)

# interp output 2
    println(io, "QG:")
    _nc = N+1
    _nr = N
    show(io, "text/plain", QG[1:_nr,1:_nc])
    println(io)

if ( ILO>1 )
        subd = hypot(norm(tril(view(A,2:N,1:ILO-1))),
                    norm(tril(view(QG,1:N,1:ILO-1))))
        println(io, "subdiagonal norm: $subd")
end # if

end # run_mb04dz()


text:
function pcor(x::VecOrMat, y::VecOrMat, z::Vector)
    n = length(z)
    ρxz = cor(x,z)
    ρzy = cor(z,y)
    ρxy = cor(x,y)
    @inbounds for j in 1:size(y, 2)
        ρzy_j = ρzy[1, j]
        @simd for i in 1:size(x, 2)
            ρxz_i = ρxz[i, 1]
            ρxy[i,j] = (ρxy[i,j] - ρxz_i * ρzy_j)/sqrt((1.0 - ρxz_i * ρxz_i))/sqrt(1.0 - ρzy_j * ρzy_j)
        end
    end
    ρxy
end


text:
using OrdinaryDiffEq, DifferentialEquations, Turing

tspan = (0.0f0, 50.0f0)

"""
Lorenz system with chaotic parameter values
"""
function lorenz_system(du, u, p, t;
                     U_func_time=U_func_time_trivial,
                     U_func_space=U_func_space_trivial)
    ρ, σ, β = p
    x, y, z = u

    Ft = U_func_time(t, u)
    Fs = U_func_space(u)

    du[1] = σ*(y-x)
    du[2] = x*(ρ-z) - y
    du[3] = x*y - β*z
    du .+= Ft .+ Fs
end

# Unforced version for use with Turing.jl
function lorenz_system(u, p, t)
    u = convert.(eltype(p),u)

    ρ, σ, β = p
    x, y, z = u
    du = [σ*(y-x); x*(ρ-z) - y; x*y - β*z]
    return du
end

# Generate data, with a time component
p = [28, 10, 8//3]
u0 = [10.0f0, 0.0f0, 0.0f0]
ts = range(tspan[1], tspan[2], length=5000)

function solve_lorenz_system(U_func_time=U_func_time_trivial,
                              U_func_space=U_func_space_trivial)
    prob = ODEProblem((du, u, p, t)->lorenz_system(du, u, p, t,
                                U_func_time=U_func_time,
                                U_func_space=U_func_space),
                    u0, tspan, p)
    sol = solve(prob, Tsit5(), saveat=ts);
    return sol
end

#####
##### Turing function for Bayesian parameter estimation
#####
function lorenz_system(u, p, t)
    u = convert.(eltype(p),u)

    ρ, σ, β = p
    x, y, z = u
    du = [σ*(y-x); x*(ρ-z) - y; x*y - β*z]
    return du
end

@model lorenz_grad_residual(y, dat) = begin
    # Lorenz parameters
    ρ ~ Normal(10, 10.0)
    σ ~ Normal(10, 1.0)
    β ~ Normal(5, 1.0)
    params = [ρ, σ, β]
    noise ~ Truncated(Normal(5, 5.0), 0, 20)

    t = [0] # Not time dependent
    for i in 1:size(y,2)
        u = dat[:,i]
        du_params = lorenz_system(u, params, t)
        y[:, i] ~ MvNormal(du_params, [noise, noise, noise])
    end
end;

#####
##### True model in SINDy syntax
#####
r, s, b = p
#      x  y z  c xx xy xz yy yz zz
 A = [[-s s 0  0  0  0  0  0  0  0];
      [r -1 0  0  0  0 -1  0  0  0];
      [0 0 -b  0  0  1  0  0  0  0]]
n = size(A, 1)
sindy_library = Dict(
    "cross_terms"=>2,
    "constant"=>nothing
);
core_dyn_true = sindyModel(ts, A,
                            convert_string2function(sindy_library),
                            ["x", "y", "z"])

# Actually export
export ts, solve_lorenz_system, lorenz_system, core_dyn_true,
    lorenz_grad_residual


text:
# Note that this script can accept some limited command-line arguments, run
# `julia build_tarballs.jl --help` to see a usage message.
using BinaryBuilder

name = "CImGui"
version = v"1.77.0"

# Collection of sources required to build CImGui
sources = [
    GitSource("https://github.com/ocornut/imgui.git",
              "9418dcb69355558f70de260483424412c5ca2fce"),

    GitSource("https://github.com/cimgui/cimgui.git",
              "bc3bc4583d0e857c05d0d76e9275c7a0652108f1"),

    DirectorySource("./bundled"),
]

# Bash recipe for building across all platforms
script = raw"""
cd $WORKSPACE/srcdir
rm cimgui/CMakeLists.txt
mv imgui wrapper/helper.c wrapper/helper.h wrapper/CMakeLists.txt cimgui/
mkdir -p cimgui/build && cd cimgui/build
cmake .. -DCMAKE_INSTALL_PREFIX=${prefix} -DCMAKE_TOOLCHAIN_FILE=${CMAKE_TARGET_TOOLCHAIN} -DCMAKE_BUILD_TYPE=Release
make -j${nproc}
make install
install_license ../LICENSE ../imgui/LICENSE.txt
"""

# These are the platforms we will build for by default, unless further
# platforms are passed in on the command line
platforms = supported_platforms()

# The products that we will ensure are always built
products = [
    LibraryProduct("libimgui-cpp", :libimgui),
    LibraryProduct("libcimgui", :libcimgui),
    LibraryProduct("libcimgui_helper", :libcimgui_helper),
    FileProduct("share/compile_commands.json", :compile_commands)
]

# Dependencies that must be installed before this package can be built
dependencies = Dependency[
]

# Build the tarballs, and possibly a `build.jl` as well.
build_tarballs(ARGS, name, version, sources, script, platforms, products, dependencies)


text:
function DECMesh(vertices, edges, triangles)
    tri = make_tri(vertices, edges, triangles)
    vint, dtri = make_dual(tri)
    vint, DECMesh(tri, dtri)
end

function make_tri(vertices, edges, triangles)
    vertices = convert(Vector{SVector{2,Float64}}, vertices)
    edges = convert(Vector{SVector{2,Int}}, edges)
    rectify_orientation!(triangles, vertices)
    triangles = convert(Vector{SVector{3,Int}}, triangles)
    return Triangulation(vertices, edges, triangles)
end

function make_dual(tri)
    vertices = tri.vertices
    edges = tri.edges
    triangles = tri.triangles

    edges_neighbors =
        [[i for i in 1:length(triangles) if edge ⊆ triangles[i]] for edge in edges]
    vertices_neighbors = [findall(x -> (i in x), edges) for i = 1:length(vertices)]

    edges_ext = findall(x -> length(x) == 1, edges_neighbors)
    edges_int = setdiff(1:length(edges), edges_ext)
    vertices_ext = findall(x -> !isempty(intersect(x, edges_ext)), vertices_neighbors)
    vertices_int = setdiff(1:length(vertices), vertices_ext)

    nvertices_ext = length(vertices_ext)

    mid_points = [mid_point(vertices[edges[edge]]...) for edge in edges_ext]
    nmid_points = length(mid_points)

    dtriangles = [circumcenter(vertices[triangle]...) for triangle in triangles]
    ndtriangles = length(dtriangles)

    append!(dtriangles, mid_points)
    append!(dtriangles, vertices[vertices_ext])

    dedges = fill([0, 0], length(edges))
    for i in edges_int
        d = edges_neighbors[i]
        rectify_orientation!(d, dtriangles[d], vertices[edges[i]])
        dedges[i] = d
    end
    j = 1
    for i in edges_ext
        d = [edges_neighbors[i][1], ndtriangles + j]
        rectify_orientation!(d, dtriangles[d], vertices[edges[i]])
        dedges[i] = d
        j += 1
    end

    dvertices = fill(Int[], length(vertices))
    for i in vertices_int
        dvertices[i] = catsort(dedges[vertices_neighbors[i]])
    end
    j = 1
    for i in vertices_ext
        dvertices[i] = catsort(dedges[vertices_neighbors[i]], ndtriangles + nmid_points + j)
        j += 1
    end
    rectify_orientation!(dvertices, dtriangles)

    return vertices_int,
    DualTriangulation(
        (ndtriangles, nmid_points, nvertices_ext),
        dtriangles,
        dedges,
        dvertices,
    )
end

function catsort(dedges, last = 0)
    vsorted = if last == 0
        fill(0, length(dedges))
    else
        fill(0, length(dedges) + 2)
    end

    uniques = reduce(union, dedges)
    graph = Dict(
        [i =>
            dedges |>
            y ->
                filter(x -> (i in x), y) |>
                x -> reduce(union, x) |> x -> setdiff(x, [i]) |> sort for i in uniques],
    )

    if last == 0
        vsorted[begin] = uniques[1]
    else
        vsorted[begin] = findfirst(v -> (length(v) == 1), graph)
        vsorted[end] = last
    end

    i = 2
    len = if last == 0
        length(uniques) - 1
    else
        length(uniques)
    end
    for _ = 1:len
        for n in graph[vsorted[i - 1]]
            if n ∉ vsorted
                vsorted[i] = n
                i += 1
                break
            end
        end
    end

    return vsorted
end


text:
function ElevationFilter(arg0::GroundStation, arg1::jdouble)
    return ElevationFilter((GroundStation, jdouble), arg0, arg1)
end

function equals(obj::Object, arg0::Object)
    return jcall(obj, "equals", jboolean, (Object,), arg0)
end

function filter(obj::ElevationFilter, arg0::ObservedMeasurement, arg1::SpacecraftState)
    return jcall(obj, "filter", void, (ObservedMeasurement, SpacecraftState), arg0, arg1)
end

function get_class(obj::Object)
    return jcall(obj, "getClass", Class, ())
end

function hash_code(obj::Object)
    return jcall(obj, "hashCode", jint, ())
end

function notify(obj::Object)
    return jcall(obj, "notify", void, ())
end

function notify_all(obj::Object)
    return jcall(obj, "notifyAll", void, ())
end

function to_string(obj::Object)
    return jcall(obj, "toString", JString, ())
end

function wait(obj::Object)
    return jcall(obj, "wait", void, ())
end

function wait(obj::Object, arg0::jlong)
    return jcall(obj, "wait", void, (jlong,), arg0)
end

function wait(obj::Object, arg0::jlong, arg1::jint)
    return jcall(obj, "wait", void, (jlong, jint), arg0, arg1)
end


text:
using FileIO

function fileio_load(f::File{format"JPEG"}; kwargs...)
    open(f.filename, "r") do io
        jpeg_decode(io; kwargs...)
    end
end
fileio_load(io::Stream{format"JPEG"}; kwargs...) = jpeg_decode(read(io); kwargs...)

function fileio_save(f::File{format"JPEG"}, img::AbstractArray; kwargs...)
    open(f.filename, "w") do io
        jpeg_encode(io, img; kwargs...)
    end
end
function fileio_save(io::Stream{format"JPEG"}, img::AbstractArray; kwargs...)
    jpeg_encode(io.io, img; kwargs...)
end


text:
using Base64, HTTP, JSON

include("auth.jl")

MEME_API = "https://meme-api.herokuapp.com/"
POST_TWEET_URL = "https://api.twitter.com/1.1/statuses/update.json"
TWITTER_UPLOAD_URL = "https://upload.twitter.com/1.1/media/upload.json"

# This is a deprecated method for media upload
function upload_meme_file(name)
    img_data = ""
    try
        open(name) do file
            img_data = read(file)
        end
    catch
        println("File does not exist.")
        exit(1)
    end
    options = Dict("media" => base64encode(img_data),
                   "media_category" => "tweet_image")
    req = post_oauth(TWITTER_UPLOAD_URL, options)
    print("UPLOAD (deprecated): ", req.status)
    json = JSON.parse(String(req.body))
    return json["media_id"]
end

function upload_by_chunk(name)
    img_data = ""
    img_size = 0
    try
        open(name) do file
            img_data = read(file)
            img_size = stat(file).size
        end
    catch
        println("File does not exist.")
        exit(1)
    end
    init_opt = Dict("command" => "INIT",
                    "total_bytes" => img_size,
                    "media_type" => "image/jpeg",
                    "media_category" => "tweet_image")
    r_init = post_oauth(TWITTER_UPLOAD_URL, init_opt)
    println("INIT: ", r_init.status)
    json = JSON.parse(String(r_init.body))
    return upload_append(img_data, json["media_id"])
end

function upload_append(img_data, media_id)
    append_opt = Dict("command" => "APPEND",
                      "media_id" => media_id,
                      "media" => base64encode(img_data),
                      "segment_index" => 0)
    r_append = post_oauth(TWITTER_UPLOAD_URL, append_opt)
    println("APPEND: ", r_append.status)
    return upload_finalize(media_id)
end

function upload_finalize(media_id)
    final_opt = Dict("command" => "FINALIZE",
                     "media_id" => media_id)
    r_final = post_oauth(TWITTER_UPLOAD_URL, final_opt)
    println("FINALIZE: ", r_final.status)
    return media_id
end

function get_meme()
    try
        channel_list = ["ProgrammerHumor", "dankmemes"]
        selection = rand(1 : 2, 1)[1]
        url = MEME_API * "gimme/" * channel_list[selection] * "/1"
        req = HTTP.request("GET", url)
        println("MEME: ", req.status)
        json = JSON.parse(String(req.body))
        img_url = json["memes"][1]["url"]
        extension = split(img_url, ".")[end]
        img_name = "image.$extension"
        download(img_url, img_name)
        return json, img_name
    catch
        println("Something went wrong. Sorry.")
        exit(1)
    end
end

function post_tweet()
    meme, name = get_meme()
    post_link = meme["memes"][1]["postLink"]
    title = meme["memes"][1]["title"]
    media_id = upload_by_chunk(name)
    options = Dict("status" => "$title\n\nSource: $post_link",
                   "media_ids" => media_id)
    r_tweet = post_oauth(POST_TWEET_URL, options)
    println("TWEET: ", r_tweet.status)
end

post_tweet()

max_stars_repo_name: abhishalya/twitter-bot

max_issues_repo_name: abhishalya/twitter-bot

max_forks_repo_name: abhishalya/twitter-bot


text:
# Views (django-sense)

function index(req, res)
	page = open(readall, "_pages/index")
    tpl = open(readall, "templates/page.html")
    body = Markdown.html(Markdown.parse(page))
    content = render(tpl, { "content" => body })

    basetpl = open(readall, "templates/base.html")
    render(basetpl, { "content" => content })
end


text:
export smooth, forget, propagate_constants, deepcopy, condition, replace_node,
    split, clone, merge, split_candidates, random_split, split_step, struct_learn,
    clone_candidates, standardize_circuit

"""
    smooth(root::Node)::Node

Create an equivalent smooth circuit from the given circuit.
"""
function smooth(root::Node)::Node
    lit_nodes = canonical_literals(root)
    f_con(n) = (n, BitSet())
    f_lit(n) = (n, BitSet(variable(n)))
    f_a(n, call) = begin
        parent_scope = BitSet()
        smooth_children = Vector{Node}(undef, num_children(n))
        map!(smooth_children, children(n)) do child
            (smooth_child, scope) = call(child)
            union!(parent_scope, scope)
            smooth_child
        end
        (conjoin([smooth_children...]; reuse=n), parent_scope)
    end
    f_o(n, call) = begin
        parent_scope = mapreduce(c -> call(c)[2], union, children(n))
        smooth_children = Vector{Node}(undef, num_children(n))
        map!(smooth_children, children(n)) do child
            (smooth_child, scope) = call(child)
            smooth_node(smooth_child, setdiff(parent_scope, scope), lit_nodes)
        end
        return (disjoin([smooth_children...]; reuse=n), parent_scope)
    end
    foldup(root, f_con, f_lit, f_a, f_o, Tuple{Node, BitSet})[1]
end

"""
    smooth_node(node::Node, missing_scope, lit_nodes)

Return a smooth version of the node where
the `missing_scope` variables are added to the scope, using literals from `lit_nodes`
"""
function smooth_node(node::Node, missing_scope, lit_nodes)
    if isempty(missing_scope)
        return node # avoid adding unnecessary nodes
    else
        get_lit(l) = get!(() -> compile(node, l), lit_nodes, l)
        ors = map(collect(missing_scope)) do v
            lit = var2lit(Var(v))
            lit_node = get_lit(lit)
            not_lit_node = get_lit(-lit)
            disjoin([lit_node, not_lit_node]; reuse=node)
        end
        conjoin([node, ors...]; reuse=node)
    end
end

"""
    forget(root::Node, is_forgotten::Function)::Node

Forget variables from the circuit.
Warning: this may or may not destroy the determinism property.
"""
function forget(root::Node, is_forgotten::Function, )::Node
    (_, true_node) = canonical_constants(root)
    if isnothing(true_node)
        true_node = compile(typeof(root), true)
    end
    f_con(n) = n
    f_lit(n) = is_forgotten(variable(n)) ? true_node : n
    f_a(n, cn) = conjoin([cn...]; reuse=n) # convert type of cn
    f_o(n, cn) = disjoin([cn...]; reuse=n)
    foldup_aggregate(root, f_con, f_lit, f_a, f_o, Node)
end

"""
    propagate_constants(root::Node)

Remove all constant leafs from the circuit
"""
function propagate_constants(root::Node; remove_unary=false)
    (false_node, true_node) = canonical_constants(root)
    if isnothing(true_node)
        true_node = compile(typeof(root), true)
    end
    if isnothing(false_node)
        false_node = compile(typeof(root), false)
    end
    f_con(n) = n
    f_lit(n) = n
    f_a(n, cn) = begin
        if any(c -> isfalse(c), cn)
            false_node
        else
            T = promote_type(typeof(n), typeof.(cn)...)
            proped_children = convert(Vector{T}, filter(c -> !istrue(c), cn))
            if isempty(proped_children)
                return true_node
            elseif remove_unary && issingle(proped_children)
                return proped_children[1]
            else
                return conjoin(proped_children; reuse=n)
            end
        end
    end
    f_o(n, cn) = begin
        if any(c -> istrue(c), cn)
            true_node
        else
            T = promote_type(typeof(n), typeof.(cn)...)
            proped_children = convert(Vector{T}, filter(c -> !isfalse(c), cn))
            if isempty(proped_children)
                return false_node
            elseif remove_unary && issingle(proped_children)
                return proped_children[1]
            else
                return disjoin(proped_children; reuse=n)
            end
        end
    end
    foldup_aggregate(root, f_con, f_lit, f_a, f_o, Node)
end

import Base.deepcopy

"""
    deepcopy(n::Node, depth::Int64)

Recursively create a copy circuit rooted at `n` to a certain depth `depth`
"""
function deepcopy(n::Node, depth::Int64, old2new::Dict{Node, Node} = Dict{Node, Node}(); cache=true)
    # TODO: can we remove the Dict and use the data field instead?
    if depth == 0 || isliteralgate(n) || isconstantgate(n)
        return n
    elseif cache && haskey(old2new, n)
        return old2new[n]
    else
        cns = map(children(n)) do c
            deepcopy(c, depth - 1, old2new; cache=cache)
        end
        if is⋀gate(n)
            new_n = conjoin(cns)
        else
            @assert is⋁gate(n)
            new_n = disjoin(cns)
        end
        if cache
            old2new[n] = new_n
        end
        return new_n
    end
end

"""
    condition(root::Node, lit::Lit; callback::Function)::Node

Return the circuit conditioned on given literal constrains
`callback` is called after modifying conjunction node
"""
function condition(root::Node, lit::Lit; callback=noop)::Node
    literals = canonical_literals(root)
    (false_node, ) = canonical_constants(root) # reuse constants when possible
    if isnothing(false_node)
        try
            false_node = compile(typeof(root), false)
        catch e
            false_node = nothing
        end
    end

    if !haskey(literals, -lit)
        root
    elseif haskey(literals, -lit) && !haskey(literals, lit)
        false_node
    else
        f_con(n) = (n, true)
        f_lit(n) = begin
            if variable(n) == lit2var(lit) && literal(n) != lit
                (false_node, false)
            else
                (n, true)
            end
        end
        f_a(n, cv) = begin
            kept = last.(cv)
            new_children = first.(cv)[kept]
            if all(kept)
                (conjoin([new_children...]; reuse=n), true)
            else
                (false_node, false)
            end
        end
        f_o(n, cv) = begin
            kept = last.(cv)
            new_children = first.(cv)[kept]
            if any(kept) && (length(new_children) == 1) && isliteralgate(new_children[1])
                (new_children[1], true)
            elseif any(kept)
                new_n = disjoin([new_children...]; reuse=n)
                callback(new_n, n, kept)
                (new_n, true)
            else
                (false_node, false)
            end
        end
        foldup_aggregate(root, f_con, f_lit, f_a, f_o, Tuple{Union{Nothing,Node}, Bool})[1]
    end
end

import Base.split

"""
    split(root::Node, (or, and)::Tuple{Node, Node}, var::Var; depth=0, sanity_check=true)

Return the circuit after spliting on edge `edge` and variable `var`
"""
function split(root::Node, (or, and)::Tuple{Node, Node}, var::Var; depth=0, sanity_check=true)
    # sanity check
    if sanity_check
        @assert depth >= 0
        @assert is⋁gate(or) && is⋀gate(and) && and in children(or)
        @assert or in or_nodes(root) && and in and_nodes(root)
        literals = canonical_literals(and)
        @assert haskey(literals, var2lit(var)) && haskey(literals, - var2lit(var))
    end

    # split
    new_children1 = map(children(and)) do c
        condition(c, var2lit(var))
    end

    new_children2 = map(children(and)) do c
        condition(c, - var2lit(var))
    end

    new_and1 = deepcopy(conjoin(new_children1), depth; cache=false)
    new_and2 = deepcopy(conjoin(new_children2), depth; cache=false)
    new_or = disjoin([[new_and1, new_and2]; filter(c -> c != and, children(or))])

    replace_node(root, or, new_or), new_or
end

"""
    split_candidates(circuit::Node)::Tuple{Vector{Tuple{Node, Node}}, Dict{Node, BitSet}}

Return the edges and variables which can be splited on
"""
function split_candidates(circuit::Node)::Tuple{Vector{Tuple{Node, Node}}, Dict{Node, BitSet}}
    candidates::Vector{Tuple{Node, Node}} = Vector{Tuple{Node, Node}}(undef, 0)
    variable_scope::Dict{Node, BitSet} = Dict{Node, BitSet}()

    f_con(n) = begin
        variable_scope[n] = BitSet()
        (false, BitSet(), BitSet())
    end
    f_lit(n) = begin
        lit = literal(n)
        var = variable(n)
        variable_scope[n] = BitSet(var)

        lit > 0 ? (false, variable_scope[n], BitSet()) : (false, BitSet(), variable_scope[n])
    end
    f_a(n, cv) = begin
        pos_scope = reduce(union, [item[2] for item in cv])
        neg_scope = reduce(union, [item[3] for item in cv])
        variable_scope[n] = intersect(pos_scope, neg_scope)

        (!isdisjoint(pos_scope, neg_scope), pos_scope, neg_scope)
    end
    f_o(n, cv) = begin
        pos_scope = reduce(union, [item[2] for item in cv])
        neg_scope = reduce(union, [item[3] for item in cv])
        variable_scope[n] = intersect(pos_scope, neg_scope)

        map(zip(first.(cv), children(n))) do (splitable, c)
            if splitable
                push!(candidates, (n, c))
            end
        end

        (any(first.(cv)), pos_scope, neg_scope)
    end

    foldup_aggregate(circuit, f_con, f_lit, f_a, f_o, Tuple{Bool, BitSet, BitSet})

    candidates, variable_scope
end

function clone_candidates(circuit::Node)::Dict{Node, Vector{Node}}
    candidates = Dict{Node, Vector{Node}}()
    parents = Dict{Node, Vector{Node}}()

    f_con(n) = begin
        false
    end
    f_lit(n) = begin
        false
    end
    f_a(n, cv) = begin
        for c in children(n)
            if !is⋁gate(c)
                continue
            end

            if c in keys(parents)
                push!(parents[c], n)
            else
                parents[c] = [n]
            end
        end
        false
    end
    f_o(n, cv) = begin
        if !(n in keys(parents))
            parents[n] = []
        end
        true
    end

    foldup_aggregate(circuit, f_con, f_lit, f_a, f_o, Bool)

    # Find candidates
    candidates = filter(p->(length(last(p)) == 2), parents) # Set of AND gates shared by exactly 2 OR gates
    candidates
end

"""
Randomly picking egde and variable from candidates
"""
function random_split(circuit::Node)
    candidates, variable_scope = split_candidates(circuit)
    or, and = rand(candidates)
    vars = Var.(collect(variable_scope[and]))
    var = rand(vars)
    (or, and), var
end

"""
Split step
"""
function split_step(circuit::Node; loss=random_split, depth=0, sanity_check=true)
    edge, var = loss(circuit)
    split(circuit, edge, var; depth=depth, sanity_check=sanity_check)
end

"""
Structure learning manager
"""
function struct_learn(circuit::Node;
    primitives=[split_step],
    kwargs=Dict(split_step=>(loss=random_split, depth=0)),
    maxiter=typemax(Int), stop::Function=x->false)

    for iter in 1 : maxiter
        primitive_step = rand(primitives)
        kwarg = kwargs[primitive_step]
        c2, _ = primitive_step(circuit; kwarg...)
        if stop(c2)
            return c2
        end
        circuit = c2
    end
    circuit
end

"""
Clone the `or` node and redirect one of its parents to the new copy
"""
function clone(root::Node, and1::Node, and2::Node, or::Node; depth=1)
    # sanity check
    @assert depth >= 1
    @assert is⋀gate(and1) && is⋀gate(and2) && is⋁gate(or)
    @assert or in children(and1) && or in children(and2)
    @assert and1 in and_nodes(root) && and2 in and_nodes(root) && or in or_nodes(root)

    # clone
    new_and2 = nothing
    new_or = deepcopy(or, depth; cache=false)

    if or == children(and1)[1]
        new_and2 = conjoin([[new_or]; filter(c -> c != or, children(and2))])
    else
        new_and2 = conjoin([filter(c -> c != or, children(and2))..., new_or])
    end

    replace_node(root, and2, new_and2)
end

import Base.merge
"""
    merge(root::Node, or1::Node, or2::Node)

Merge two circuits.
"""
function merge(root::Node, or1::Node, or2::Node)
    # sanity check
    @assert is⋁gate(or1) && is⋁gate(or2)
    @assert or1 in or_nodes(root) && or2 in or_nodes(root)

    # Simply replace the occurrencce of the larger OR with the smaller OR
    n_or1 = num_nodes(or1)
    n_or2 = num_nodes(or2)

    if n_or1 >= n_or2
        replace_node(root, or1, or2)
    else
        replace_node(root, or2, or1)
    end
end

"""
Replace node `old` with node `new` in circuit `root`
"""
function replace_node(root::Node, old::Node, new::Node; callback::Function=(x, y, z) -> nothing)::Node
    @assert GateType(old) == GateType(new)
    f_con(n) = old == n ? new : n
    f_lit = f_con
    f_a(n, cns) = old == n ? new : conjoin([cns...]; reuse=n)
    f_o(n, cns) = old == n ? new : begin
            new_n = disjoin([cns...]; reuse=n)
            callback(new_n, n, trues(length(cns)))
            new_n
        end
    foldup_aggregate(root, f_con, f_lit, f_a, f_o, Node)
end

"""
Standardize the circuit:

1. Children of or nodes are and nodes.
2. Children of and nodes are or nodes.
3. Each and node has exactly two children.

Note: for circuits with parameters this function will not keep the
      parameters equavalent.
"""
function standardize_circuit(circuit::LogicCircuit)
    compile_and_node(left_child::LogicCircuit, right_child::LogicCircuit) = begin
        if is⋀gate(left_child)
            left_child = disjoin(left_child)
        end
        if is⋀gate(right_child)
            right_child = disjoin(right_child)
        end
        disjoin(conjoin(left_child, right_child))
    end

    f_con(node) = compile(node, constant(node))
    f_lit(node) = compile(node, literal(node))
    f_a(node, cs) = begin
        node_children::Array{LogicCircuit} = Array{LogicCircuit}(undef, 0)

        @assert length(cs) > 1 "⋀gate $(node) should contain at least two children, but only has one."

        # Flatten the children by ensuring all children of `node` are ⋁gate or leaf gate.
        for child_node in cs
            if is⋀gate(child_node)
                for n in children(child_node)
                    push!(node_children, n)
                end
            else
                push!(node_children, child_node)
            end
        end

        if length(node_children) > 2
            n = compile_and_node(node_children[end - 1], node_children[end])
            for idx = length(node_children) - 2 : -1 : 2
                n = compile_and_node(node_children[idx], n)
            end
            node_children = [first(node_children), n]
        end

        conjoin(node_children...)
    end
    f_o(node, cs) = begin
        node_children::Array{LogicCircuit} = Array{LogicCircuit}(undef, 0)

        for child_node in cs
            if is⋁gate(child_node)
                for n in children(child_node)
                    push!(node_children, n)
                end
            else # is⋀gate(child_node) || isleafgate(child_node)
                push!(node_children, child_node)
            end
        end

        disjoin(node_children...)
    end

    foldup_aggregate(circuit, f_con, f_lit, f_a, f_o, LogicCircuit)
end


text:
    using Revise, Parameters, Logging, Test
    using CUDA
    includet("C:\\GitHub\\GitHub\\NuclearMedEval\\test\\includeAllUseFullForTest.jl")
    using ..CUDAGpuUtils ,..IterationUtils,..ReductionUtils , ..MemoryUtils,..CUDAAtomicUtils
    using ..BitWiseUtils,..ResultListUtils, ..MetadataAnalyzePass,..MetaDataUtils,..WorkQueueUtils,..ProcessMainDataVerB,..HFUtils, ..ScanForDuplicates

numb = UInt32(0)
pos = UInt32(2)
@test !isBit1AtPos(numb,pos)

@setBitTo1(numb,pos)

@test isBit1AtPos(numb,pos)
@test !isBit1AtPos(numb,pos-1)
@test !isBit1AtPos(numb,pos+1)

numb = UInt32(0)
pos = UInt32(2)

@setBitTo(numb,pos,true)
@test isBit1AtPos(numb,pos)

@setBitTo(numb,pos,false)
@test !isBit1AtPos(numb,pos)

numbs=[2]
@setBitTo(numbs[1],pos,true)

##### bitDilatate
numb = UInt32(0)
locArr = UInt32(0)
@setBitTo1UINt32(numb,1)
@setBitTo1UINt32(numb,5)
typeof(numb)

numbB= bitDilatate(numb,locArr)

@test isBit1AtPos(numbB,1)
@test isBit1AtPos(numbB,2)
@test !isBit1AtPos(numbB,3)
@test isBit1AtPos(numbB,4)
@test isBit1AtPos(numbB,5)
@test isBit1AtPos(numbB,6)
@test !isBit1AtPos(numbB,7)

@test !isBit1AtPos(numb,2)

#### bitPassOnes

source = UInt32(0)
target = UInt32(0)
@setBitTo1UINt32(source,1)
@setBitTo1UINt32(source,5)
@setBitTo1UINt32(target,4)
new =llvmPassOnes(source,target)

@test isBit1AtPos(new,1)
@test !isBit1AtPos(new,2)
@test !isBit1AtPos(new,3)
@test isBit1AtPos(new,4)
@test isBit1AtPos(new,5)

nnn = UInt8(1)
isBit1AtPos(nnn,1)

function llvmPassOnesBB(source::UInt32,target::UInt32,array::Vector{UInt32})::UInt32
    Base.llvmcall("""
    %3 = or i32 %0, %1
    %2[1]= %3
    ret i32 %3""", UInt32, Tuple{UInt32,UInt32}, source, target)
end

arr = [UInt32(1), UInt32(2)]
source = UInt32(0)
target = UInt32(0)
@setBitTo1UINt32(source,1)
@setBitTo1UINt32(source,5)
@setBitTo1UINt32(target,4)
new =llvmPassOnesBB(source,target,arr)

# using Revise, Parameters, Logging, Test
# using CUDA
# includet("C:\\GitHub\\GitHub\\NuclearMedEval\\test\\includeAllUseFullForTest.jl")

# using CUDA, ..CUDAGpuUtils
# function simpleShuffle()
#     xx = threadIdxX()
#     xx+=shfl_down_sync(FULL_MASK, xx, 1)
#     CUDA.@cuprint "xx $(xx) \n"

#     return
# end

# @cuda threads=(32,1) blocks=1 simpleShuffle()