text:
"""
	grd2xyz(cmd0::String="", arg1=[], kwargs...)

Reads one 2-D grid and returns xyz-triplets.

Full option list at [`grd2xyz`](http://gmt.soest.hawaii.edu/doc/latest/grd2xyz.html)

Parameters
----------

- $(GMT.opt_J)
- **C** : **row_col** : -- Bool --

    Replace the x- and y-coordinates on output with the corresponding column and row numbers.
    [`-C`](http://gmt.soest.hawaii.edu/doc/latest/grd2xyz.html#c)
- $(GMT.opt_R)
- $(GMT.opt_V)
- **W** : **weight** : -- Str --           Flags = [a|weight]

    Write out x,y,z,w, where w is the supplied weight (or 1 if not supplied) [Default writes x,y,z only].
    [`-W`](http://gmt.soest.hawaii.edu/doc/latest/grd2xyz.html#w)
- **Z** : **flags** : -- Str --

    Write a 1-column table. Output will be organized according to the specified ordering
    convention contained in ``flags``.
    [`-Z`](http://gmt.soest.hawaii.edu/doc/latest/grd2xyz.html#z)
- $(GMT.opt_write)
- $(GMT.opt_append)
- $(GMT.opt_bo)
- $(GMT.opt_d)
- $(GMT.opt_f)
- $(GMT.opt_h)
- $(GMT.opt_o)
- $(GMT.opt_s)
"""
function grd2xyz(cmd0::String="", arg1=[]; kwargs...)

	length(kwargs) == 0 && occursin(" -", cmd0) && return monolitic("grd2xyz", cmd0, arg1)	# Speedy mode

	d = KW(kwargs)

	cmd, = parse_R("", d)
	cmd = parse_V_params(cmd, d)
	cmd, = parse_bo(cmd, d)
	cmd, = parse_d(cmd, d)
	cmd, = parse_f(cmd, d)
	cmd, = parse_h(cmd, d)
	cmd, = parse_o(cmd, d)
	cmd, = parse_s(cmd, d)

	cmd = add_opt(cmd, 'C', d, [:C :row_col])
	cmd = add_opt(cmd, 'W', d, [:W :weight])
	cmd = add_opt(cmd, 'Z', d, [:Z :flags])

	cmd, got_fname, arg1 = find_data(d, cmd0, cmd, 1, arg1)
	if (isa(arg1, Array{<:Number}))		arg1 = mat2grid(arg1)	end
	return common_grd(d, cmd, got_fname, 1, "grd2xyz", arg1)		# Finish build cmd and run it
end

# ---------------------------------------------------------------------------------------------------
grd2xyz(arg1=[], cmd0::String=""; kw...) = grd2xyz(cmd0, arg1; kw...)


text:
module ERFA

using LinearAlgebra: I, cross, dot, norm, normalize

using ERFA_jll

export ERFAException
export DPI, D2PI, DR2D, DD2R, DR2AS, DAS2R, DS2R, TURNAS, DMAS2R, DTY, DAYSEC, DJY, DJC
export DJ00, DJM0, DJM00, DJM77, TTMTAI, DAU, CMPS, AULT, DC, ELG, ELB, TDB0, SRS
export WGS84, GRS80, WGS72

include("exports.jl")
include("erfa_common.jl")

include("a.jl")
include("b.jl")
include("c.jl")
include("d.jl")
include("e.jl")
include("f.jl")
include("g.jl")
include("h.jl")
include("i.jl")
include("j.jl")
include("l.jl")
include("m.jl")
include("n.jl")
include("o.jl")
include("p.jl")
include("r.jl")
include("s.jl")
include("t.jl")
include("u.jl")
include("x.jl")
include("z.jl")

end # module

ext: jl


text:
using DataFrames
using ProgressMeter
using Distributed

"""
    lev(s, t)

Levenshtein distance for any iterable, not just strings. Implemented from the
pseudocode on the Wikipedia [page for Levenshtein distance]
(https://en.wikipedia.org/wiki/Levenshtein_distance).

# Parameters
* **s** The first iterable in the comparison; can be a string, Array, tuple,
    etc., so long as it can be indexed
* **t** The second iterable in the comparison; can be a string, Array, tuple,
    etc., so long as it can be indexed

# Returns
* The calculated Levenshtein distance between `s` and `t`
"""
function lev(s, t)
    m = length(s)
    n = length(t)
    d = Array{Int}(zeros(m+1, n+1))

    for i=2:(m+1)
        @inbounds d[i, 1] = i-1
    end

    for j=2:(n+1)
        @inbounds d[1, j] = j-1
    end

    for j=2:(n+1)
        for i=2:(m+1)
            @inbounds if s[i-1] == t[j-1]
                substitutionCost = 0
            else
                substitutionCost = 1
            end
            @inbounds d[i, j] = min(d[i-1, j] + 1, # Deletion
                            d[i, j-1] + 1, # Insertion
                            d[i-1, j-1] + substitutionCost) # Substitution
        end
    end

    @inbounds return d[m+1, n+1]
end

"""
    pnd(corpus::Array, queries::Array; [progress=true])

Calculate the phonological neighborhood density (pnd) for each item in `queries`
based on the items in `corpus`. This function uses a vantage point tree data
structure to speed up the search for neighbors by pruning the search space. This
function should work regardless of whether the items in `queries` are in
`corpus` or not.

# Parameters
* **corpus** The corpus to be queried for phonological neighbors
* **queries** The items to query phonological neighbors for in `corpus`
* **progress** Whether to display a progress meter or not

# Returns
* A `DataFrame` with the queries in the first column and the phonological
    neighborhood density in the second
"""
function pnd(corpus::Array, queries::Array; progress=true)

    tree = TextVPTree(corpus, lev)
    results = Vector()

    bs = ceil(Int64, length(queries) / nworkers())

    neighbors(query) = length(filter(x -> x != query, radiusSearch(tree, query, 1)))

    if progress
        results = @showprogress pmap(neighbors, queries, batch_size=bs)
    else
        results = pmap(neighbors, queries; batch_size=bs)
    end

    return DataFrame(Query=queries, PND=results)
end


text:
using Pickle
using Pickle: TableBlock, HierarchicalTable

"""
  `load_state(model_name)`

load the state from the given model name as NamedTuple.
"""
function load_state(model_name)
  state_dict = load_state_dict(model_name)
  state = state_dict_to_namedtuple(state_dict)
  return state
end

"""
  `load_state_dict(model_name)`

load the state_dict from the given model name.

See also: [`state_dict_to_namedtuple`](@ref)
"""
function load_state_dict(model_name)
  state_dict = Pickle.Torch.THload(get_registered_weight_path(model_name))
  return state_dict
end

"""
  `state_dict_to_namedtuple(state_dict)`

convert state_dict into NamedTuple.
"""
function state_dict_to_namedtuple(state_dict)
  ht = Pickle.HierarchicalTable()
  foreach(((k, v),)->setindex!(ht, v, k), pairs(state_dict))
  _ht2nt(ht)
end

_ht2nt(x::Some) = something(x)
_ht2nt(x::HierarchicalTable) = _ht2nt(x.head)
function _ht2nt(x::TableBlock)
  if iszero(length(x.entry))
    return ()
  else
    tks = Tuple(keys(x.entry))
    if all(Base.Fix1(all, isdigit), tks)
      inds = Vector(undef, length(tks))
      foreach(tks) do is
        i = parse(Int, is) + 1
        inds[i] = _ht2nt(x.entry[is])
      end
      return inds
    else
      cs = map(_ht2nt, values(x.entry))
      ns = map(Symbol, tks)
      return NamedTuple{ns}(cs)
    end
  end
end


text:
using Plots
"""
    GraphSig_Plot(G::GraphSig; symmetric::Bool = false,
      markersize::Float64 = 2.,
      markercolor::Symbol = :balance,
      markershape::Symbol = :circle,
      markerstrokewidth::Float64 = 1.0,
      markerstrokealpha::Float64 = 1.0,
      markervaluevaries::Bool = true,
      linewidth::Float64 = 1.,
      linecolor::Symbol = :blue,
      linestyle::Symbol = :solid,
      clim::Tuple{Float64,Float64} = (0., 0.),
      notitle::Bool = false, nocolorbar::Bool = false, nolegend::Bool = true,
      stemplot::Bool = false, sortnodes::Symbol = :normal)

Display a plot of the data in a GraphSig object

### Input Argument
*   `G::GraphSig`:        an input GraphSig object
*   `symmetric`           symmetrize the colorbar
*   `markersize`          the size of the nodes
*   `markercolor`         markercolor scheme
*   `markershape`         shape of marker
*   `markerstrokewidth`   width of marker stroke
*   `markerstrokealpha`   capacity of marker stroke
*   `markervaluevaries`   if the marker color depends on the signal value
*   `linewidth`           the width of the lines in gplot
*   `linecolor`           the color of the lines (1D) / graph edges (2D & 3D)
*   `linestyle`:          the style of line
*   `notitle`             display a title
*   `nocolorbar`          display a colorbar
*   `nolegend`            display a legend
*   `stemplot`            use a stem plot
*   `clim`                specify the dynamic display range
*   `sortnodes`           plot the signal values from smallest to largest in magnitude
"""
function GraphSig_Plot(G::GraphSig; symmetric::Bool = false,
      markersize::Any = 2.,
      markercolor::Symbol = :balance,
      markershape::Symbol = :circle,
      markerstrokewidth::Float64 = 1.0,
      markerstrokealpha::Float64 = 1.0,
      markervaluevaries::Bool = true,
      linewidth::Float64 = 1.,
      linecolor::Symbol = :blue,
      linestyle::Symbol = :solid,
      clim::Tuple{Float64,Float64} = (0., 0.),
      notitle::Bool = false, nocolorbar::Bool = false, nolegend::Bool = true,
      stemplot::Bool = false, sortnodes::Symbol = :normal)

    # only run this for 1D signals
    fcols = Base.size(G.f, 2)
    if fcols == 0
        G.f = zeros(G.length, 1)
    elseif fcols != 1
        warn("Cannot handle multiple signals at this time!")
        return
    end

    ## Case 1: the graph does not have spatial coordinates or it is too large
    if G.length > 10^4 || G.dim > 3 || G.dim < 1
        if G.length < 10^6
            spy(G.W, show = true)            # not too clear if we can change the dot size.
            scatter!((1:G.length)', (1:G.length)', marker=(:circle, 4), zcolor = G.f, title = G.name)
        else
            warn("Cannot handle a graph with more than 10^6 nodes at this time!")
        end
        return
    end

    ## Case 2: the graph does have spatial coordinates (G.dim = 1, 2, or 3).

    # extract the plot specifications
    # set(0, 'DefaultFigureVisible', 'on');

    #################################
    ############# 1-D ###############
    #################################
    if G.dim == 1
        if G.length < 65 || stemplot
            # stem(G.xy, G.f, '-o','LineWidth',2,'Color',linecolor,'MarkerFaceColor',linecolor,'MarkerSize',markersize(1));
            # stem(G.xy, G.f, linetype = [:sticks :scatter], linewidth = 2, linecolor = linecolor, markercolor = linecolor, markersize = markersize)
            stem(G.f; x = G.xy, lw = linewidth, c = linecolor, ms = markersize, shape = markershape)
        else
            # plot(G.xy, G.f, '-','Color',linecolor,'LineWidth',linewidth);
            plot(G.xy, G.f, linetype = :line, linecolor = linecolor, linewidth = linewidth);
        end

        xmin = minimum(G.xy)
        xmax = maximum(G.xy)
        dx = xmax - xmin
        if dx < 10 * eps()
            dx = 1
        end
        xmin = xmin - 0.1 * dx
        xmax = xmax + 0.1 * dx

        ymin = minimum(G.f)
        ymax = maximum(G.f)
        dy = ymax - ymin
        if dy < 10 * eps()
            dy = 1
        end
        ymin = ymin - 0.1*dy
        ymax = ymax + 0.1*dy

        # symmetric?
        if symmetric
            ymax = 1.1 * max(ymax, -ymin)
            ymin = -ymax
        end

        xlims!(xmin, xmax)
        ylims!(ymin, ymax)

        # title?
        if ~isempty(G.name) && ~notitle
            title!(G.name)
        end

        # legend?
        if nolegend
            plot!(leg = false)
        end

        #################################
        ######### 2-D  or 3-D ###########
        #################################

    else

        #fig = figure('visible','on');

        # plot the graph
        gplot(G.W, G.xy; style = linestyle, color = linecolor, width = linewidth)
        #hold on

        # plot the nodes
        scatter_gplot!(G.xy;
                       marker = G.f, ms = markersize,
                       shape = markershape, mswidth = markerstrokewidth,
                       msalpha = markerstrokealpha, plotOrder = sortnodes,
                       c = markercolor)

        # colorbar?
        if nocolorbar
            plot!(colorbar = :none)
        end

        # custom dynamic display range?
        if clim != (0., 0.)
            plot!(clim = clim)
        end

        # symmetric?
        if symmetric
            if clim != (0., 0.)
                cmax = maximum(abs.(clim));
            else
                cmax = maximum(abs.(G.f[:]));
            end
            if cmax < 10*eps()
                cmax = 1;
            end
            plot!(cLims = (cmax, cmax))
        end

        # title?
        if !isempty(G.name) && !notitle
            title!(G.name);
        end

        # legend?
        if nolegend
            plot!(leg = false)
        end
    end
end


text:
using OptimalTransport: wasserstein
using Distributions: Categorical
using StringDistances: Levenshtein

function get_T_similarity(snapshot_1, snapshot_2, T)
    distribution_1 = get_genotype_distribution(snapshot_1, T)
    distribution_2 = get_genotype_distribution(snapshot_2, T)

    T_similarity = _wasserstein(distribution_1, distribution_2, Levenshtein())

    return T_similarity
end

function get_genotype_distribution(snapshot, T)
    run_n_timesteps!(snapshot, T)
    return get_genotype_distribution(snapshot)
end

max_stars_repo_name: tochsner/ALifeBenchmark

max_issues_repo_name: tochsner/ALifeBenchmark

max_forks_repo_name: tochsner/ALifeBenchmark


text:
# TODO Make work with multidimensional arrays
immutable EnumerableArray{T,S} <: Enumerable{T}
    source::S
end

function query{TS,N}(source::Array{TS,N})
    return EnumerableArray{TS,Array{TS,N}}(source)
end

function start{T,S}(iter::EnumerableArray{T,S})
    return 1
end

function next{T,S}(iter::EnumerableArray{T,S}, state)
    return iter.source[state], state+1
end

function done{T,S}(iter::EnumerableArray{T,S}, state)
    return state>length(iter.source)
end


text:
#!/usr/bin/env julia

using LombScargle, BenchmarkTools, FFTW

N  = round.(Int, exp10.(range(1, stop=5, length=15)))
for nthreads in (1, 4)
    FFTW.set_num_threads(nthreads)
    info(string(nthreads) * " FFTW thread(s)")
    open(joinpath(@__DIR__, "julia_times-" * string(nthreads) * ".dat"), "w") do file
        for (j, n) in enumerate(N)
            println("Iteration ", j, " (", n, " datapoints)")
            t = collect(range(0.01, stop=10, length=n))
            s = sin.(t) .+ 1.5*cospi.(4t) .+ 3
            plan = LombScargle.plan(t, s, fast = true, flags = FFTW.MEASURE)
            time = minimum(@benchmark(lombscargle($plan)).times)/1e9
            println(" Julia:  ", time, " seconds")
            println(file, n, "\t", time)
        end
    end
end


text:
### A Pluto.jl notebook ###
# v0.11.14

using Markdown
using InteractiveUtils

# ╔═╡ 3da9f1a8-f65a-11ea-2a34-c7eb4db8cd54
using Plots, Distributions, LinearAlgebra,Distributions, ProgressMeter

# ╔═╡ 316f1058-f65a-11ea-33d8-997377fba2b0
md"Preliminary: load packages, setup latex macro"

# ╔═╡ b37d84e0-f659-11ea-0e0e-1d14d3d5aba5
md"# Calculating the ELBO for PFVG
Here we detail the various calculation we need."

# ╔═╡ 5cabfa9e-f65a-11ea-0695-c73c26cf3342
md"## Model
$$\int \mathcal{N}(\mathbf{y}|\mathbf{A}\mathbf{f} - \mathbf{b}, \mathbf{S}) \mathcal{N}(\mathbf{f}|\mathbf{0},\mathbf{K}) \ \mathbf{df} \ \mathbf{dA}\ \mathbf{db}$$
where the matrix $$\mathbf{A}(\mathbf{\alpha})$$ and vector $$\mathbf{b}(\mathbf{\beta})$$ are function of the respective parameters.
"

# ╔═╡ d24a45be-f65b-11ea-2c07-d19988003b56
md"## Lower bound
$$\int q(\mathbf{f}) q(\mathbf{\alpha}) q(\mathbf{\beta}) \log\mathcal{N}(\mathbf{y}|\mathbf{A}\mathbf{f} - \mathbf{b}, \mathbf{S}) \mathcal{N}(\mathbf{f}|\mathbf{0},\mathbf{K}) \ \mathbf{df} \ \mathbf{dA}\ \mathbf{db} + \mathcal{H}[q_f] + \mathcal{H}[q_\alpha] + \mathcal{H}[q_\beta]$$
"

# ╔═╡ 3e54363e-f65c-11ea-0150-7d09f848a86e
md"## Inside the log
$$-\frac{1}{2}(\mathbf{y} - \mathbf{Af}-\mathbf{b})^T\mathbf{S}^{-1} (\mathbf{y} - \mathbf{Af}-\mathbf{b}) - \frac{1}{2}\log|\mathbf{K}| -\frac{1}{2}\mathbf{f}^T\mathbf{K}^{-1}\mathbf{f}$$
"

# ╔═╡ ab41e982-f65c-11ea-26f9-01f9de7053a3
md"## Work out posterior $$q(\mathbf{f})$$ of GP latent variables

Keep only relevant terms:

$$-\frac{1}{2}(-2\mathbf{y}^T\mathbf{S}^{-1}\mathbf{Af} + 2\mathbf{b}^T\mathbf{S}^{-1}\mathbf{Af} + \mathbf{f^T A^T S^{-1} A f})
-\frac{1}{2} \mathbf{f}^T\mathbf{K}^{-1}\mathbf{f}$$

We recognise the quadratic which leads to a Gaussian posterior:

$$-\frac{1}{2}(-2(\mathbf{-y+b})^T\mathbf{S}^{-1}\mathbf{Af}  + \mathbf{f^T A^T S^{-1} A f})
-\frac{1}{2} \mathbf{f}^T\mathbf{K}^{-1}\mathbf{f}$$

$$(\mathbf{y-b})^T\mathbf{S}^{-1}\mathbf{Af}  - \frac{1}{2} \mathbf{f^T (A^T S^{-1} A +K^{-1}) f}$$

The posterior reads:

$$q(\mathbf{f}) = \mathcal{N}(\mathbf{f} |  \mathbf{\Sigma_f} \ \mathbf{A^T S^{-1}(y-b)}, \ \underbrace{\mathbf{A^T S^{-1} A +K^{-1}}}_{\mathbf{\Sigma_f}} )$$

"

# ╔═╡ fe789bd2-f66e-11ea-3400-61f55bbf640b

# ╔═╡ Cell order:
# ╟─316f1058-f65a-11ea-33d8-997377fba2b0
# ╠═3da9f1a8-f65a-11ea-2a34-c7eb4db8cd54
# ╟─b37d84e0-f659-11ea-0e0e-1d14d3d5aba5
# ╟─5cabfa9e-f65a-11ea-0695-c73c26cf3342
# ╟─d24a45be-f65b-11ea-2c07-d19988003b56
# ╟─3e54363e-f65c-11ea-0150-7d09f848a86e
# ╟─ab41e982-f65c-11ea-26f9-01f9de7053a3
# ╠═fe789bd2-f66e-11ea-3400-61f55bbf640b


text:
using Combinatorics

valid_count = 0
for passphrase in eachline(STDIN)
  words = split(passphrase, " ")
  is_valid = true

  for (w1, w2) in combinations(words, 2)
    if w1 == w2
      is_valid = false
      break
    end
  end

  if is_valid
    valid_count += 1
  end
end

println("Valid passphrases for part 1: $valid_count")


text:
@inline Base.Float64(x::DoubleFloat{T}) where {T<:IEEEFloat} = Float64(HI(x))
@inline Base.Float32(x::DoubleFloat{T}) where {T<:IEEEFloat} = Float32(Float64(x))
@inline Base.Float16(x::DoubleFloat{T}) where {T<:IEEEFloat} = Float16(Float64(x))

Base.Int128(x::DoubleFloat{T}) where {T<:IEEEFloat} = Int128(round(BigInt, BigFloat(HI(x)) + BigFloat(LO(x))))
Base.Int64(x::DoubleFloat{T}) where {T<:IEEEFloat} = Int64(round(BigInt, BigFloat(HI(x)) + BigFloat(LO(x))))
Base.Int32(x::DoubleFloat{T}) where {T<:IEEEFloat} = Int32(round(BigInt, BigFloat(HI(x)) + BigFloat(LO(x))))
Base.Int16(x::DoubleFloat{T}) where {T<:IEEEFloat} = Int16(round(BigInt, BigFloat(HI(x)) + BigFloat(LO(x))))
Base.Int8(x::DoubleFloat{T}) where {T<:IEEEFloat} = Int8(round(BigInt, BigFloat(HI(x)) + BigFloat(LO(x))))

@inline Base.BigFloat(x::DoubleFloat{T}) where {T<:AbstractFloat} = BigFloat(HI(x)) + BigFloat(LO(x))
@inline Base.BigInt(x::DoubleFloat{T}) where {T<:AbstractFloat} = BigInt(BigFloat(x))

for (F,D) in ((:Float64, :Double64), (:Float32, :Double32), (:Float16, :Double16))
  @eval begin
    function $D(x::BigFloat)
        hi = $F(x)
        lo = $F(x - hi)
        return $D(hi, lo)
    end
    $D(x::BigInt) = $D(BigFloat(x))
  end
end

# constants for BigFloat precision are > 2*significandbits(DoubleT)

Double64(x::Irrational{S}) where {S} = Double64(BigFloat(x, 250))
Double32(x::Irrational{S}) where {S} = Double32(BigFloat(x, 122))
Double16(x::Irrational{S}) where {S} = Double16(BigFloat(x,  60))

Double64(x::Rational{S}) where {S} = Double64(BigFloat(x, 250))
Double32(x::Rational{S}) where {S} = Double32(BigFloat(x, 122))
Double16(x::Rational{S}) where {S} = Double16(BigFloat(x,  60))


text:
get_pricing_engine{Y}(::Discount, yts::Y) = DiscountingBondEngine(yts)

function apply_termstructure{R <: RateHelper, T <: TermStructure}(rate::R, ts::T)
  newRate = update_termstructure(rate, ts)

  return newRate
end

apply_termstructure{B <: BondHelper, T <: TermStructure}(b::B, ts::T) = update_termstructure(b, ts)
function apply_termstructure{T <: TermStructure}(s::SwapRateHelper, ts::T)
  newSwapHelper = update_termstructure(s, ts)

  return newSwapHelper
end

function apply_termstructure(cds::AbstractCDSHelper, ts::TermStructure)
  # have to clone
  newCDS = clone(cds, ts)
  reset_engine!(newCDS)

  return newCDS
end

# BOOTSTRAPPING
type IterativeBootstrap <: Bootstrap
  firstSolver::BrentSolver
  solver::FiniteDifferenceNewtonSafe

  IterativeBootstrap() = new(BrentSolver(), FiniteDifferenceNewtonSafe())
end

# returns initial bootstrap state to Term Structure
function initialize{T <: TermStructure}(::IterativeBootstrap, ts::T)
  # get the intial data based on trait
  data_initial = initial_value(ts.trait)
  n = length(ts.instruments) + 1

  # get a new pricing engine based on type
  # pe = get_pricing_engine(ts.trait, ts)

  # initialize data
  if data_initial == 1.0
    ts.data = ones(n)
  elseif data_initial == 0.0
    ts.data = zeros(n)
  else
    ts.data = fill(data_initial, n)
  end

  # build times and error vectors (which have the functions for the solver)
  ts.times[1] = time_from_reference(ts, ts.referenceDate)
  ts.dates[1] = ts.referenceDate
  for i = 2:n
    @inbounds ts.times[i] = time_from_reference(ts, maturity_date(ts.instruments[i - 1]))
    @inbounds ts.dates[i] = maturity_date(ts.instruments[i - 1])
    # set yield term Structure
    @inbounds ts.instruments[i - 1] = apply_termstructure(ts.instruments[i - 1], ts)
    # set error function
    @inbounds ts.errors[i] = bootstrap_error(i, ts.instruments[i - 1], ts)
  end

  # initialize interpolation
  QuantLib.Math.initialize!(ts.interp, ts.times, ts.data)
end

function _calculate!{T <: TermStructure}(boot::IterativeBootstrap, ts::T)
  max_iter = max_iterations(ts.trait)
  valid_data = ts.validCurve

  iterations = 0
  # if we get through this loop, we haven't converged
  while iterations < max_iter
    prev_data = copy(ts.data) # need actual copy, not pointer, to check later for convergence
    for i = 2:length(ts.data)

      # bracket root and calculate guess
      min = min_value_after(ts.trait, i, ts, valid_data)
      max = max_value_after(ts.trait, i, ts, valid_data)
      g = guess(ts.trait, i, ts, valid_data)

      # adjust if needed
      if g >= max
        g = max - (max - min) / 5.0
      elseif (g <= min)
        g = min + (max - min) / 5.0
      end

      if !valid_data
        update_idx = i == length(ts.data) ? 1 : i + 1
        QuantLib.Math.update!(ts.interp, update_idx, ts.data[1])
      end

      # put this in a try / catch
      if !valid_data
        # use first solver
        root = solve(boot.firstSolver, ts.errors[i], ts.accuracy, g, min, max)
      else
        root = solve(boot.solver, ts.errors[i], ts.accuracy, g, min, max)
      end
    end

    # let's check for convergence
    change = abs(ts.data[2] - prev_data[2])
    for i=3:length(ts.data)
      change = max(change, abs(ts.data[i] - prev_data[i]))
    end
    if change <= ts.accuracy
      break # bye
    end

    iterations += 1

    valid_data = true
  end
  ts.validCurve = true

  return ts
end

function bootstrap_error{T <: BootstrapHelper, Y <: TermStructure}(i::Int, inst::T, ts::Y)
  function bootstrap_error_inner(g::Float64)
    # update trait
    update_guess!(ts.trait, i, ts, g)
    QuantLib.Math.update!(ts.interp, i, g)
    # qe =
    # if i > 7
    #   println("GUESS: $i : $g")
    #   println("QUOTE ERROR ", qe)
    # end
    return quote_error(inst)
  end

  return bootstrap_error_inner
end

quote_error{B <: BondHelper}(inst::B) = QuantLib.value(inst) - implied_quote(inst) # recalculate
quote_error{R <: RateHelper}(rate::R) = QuantLib.value(rate) - implied_quote(rate)
quote_error(rate::AbstractCDSHelper) = QuantLib.value(rate) - implied_quote(rate)


text:
using NLsolve

λ_p    = 0.15             # net markup
α      = 1-0.60/(1+λ_p)
δ    = 0.025
π_ss = 1.02^(1/4)
γ    = log(1.02)/4
disc = 0.99

h    = 0.5
g_x  = 0.18

s    = 0.995
θ    = 1/4
f_s1 = 3/4
f_h1 = 1/4
A2   = 0.7
v    = 0.7
t_x  = 0.2
τ_x  = 0
#ψ_h1 = θ * f_h1 / 1.1
#ψ_h2 = θ * (1 - f_h1) / 1.1
ψ_h1 = 0
ψ_h2 = 0
H1   = 1
H2   = 1

L1   = ((1-θ) * f_s1 + θ * f_h1) * H1
L2   = ((1-θ) * (1 - f_s1) + θ * (1 - f_h1)) * H2

ρ    = exp(γ) / disc - (1 - δ)

w1   = (α ^ α * (1 - α) ^ (1 - α) / (1 + λ_p) / (ρ ^ α)) ^ (1 / (1 - α))
w2   = A2 * w1

klr1 = w1 / ρ * α / (1 - α)
klr2 = A2 * klr1
k    = klr1 * L1

F1   = λ_p / (1 + λ_p) * L1 * klr1 ^ α
y1   = L1 * klr1 ^ α - F1
y    = y1 / v
y2   = (1 - v) * y
F2   = λ_p * y2 / A2
x    = y

is   = (1 - (1 - δ) * exp(-γ)) * exp(γ) * k / (1 - θ)
g    = g_x * x
c    = y - g - (1 - θ) * is

τ    = τ_x * x
τ_h1 = τ * ψ_h1 / θ / (1 - f_h1)
τ_h2 = τ * ψ_h2 / θ / (1 - f_h1)
τ_s  = τ - θ * (f_h1 * τ_h1 + (1 - f_h1) * τ_h2)

t    = t_x * x
t_h1 = t * ψ_h1 / θ / f_h1
t_h2 = t * ψ_h2 / θ / (1 - f_h1)
t_s  = (t - θ * ( f_h1 * t_h1 + (1 - f_h1) * t_h2)) / (1 - θ)

c_h10 = c / θ / f_h1
c_h20 = c / θ / (1 - f_h1)
c_s0  = c / (1 - θ)
λ_h10 = exp(γ) / (exp(γ) * c_h10 - h * c)
λ_h20 = exp(γ) / (exp(γ) * c_h20 - h * c)
λ_s0  = exp(γ) / (exp(γ) * c_s0 - h * c)
bR0   = (g - t + τ) / (1 - (1.025 / 4) / exp(γ) / π_ss)
R0    = (1.025 / 4)

function distSS(x)
global w1, w2, H1, H2, th1, th2, τ_h1, τ_h2, c, g, t, τ
global s, θ, f_h1, f_s1, γ, π_ss, h, disc

c_h1 = x[1]
c_h2 = x[2]
c_s   = x[3]
λ_h1 = x[4]
λ_h2 = x[5]
λ_s  = x[6]
bR   = x[7]
R    = x[8]

r = zeros(8)

r[1] = w1 * H1 - t_h1 + τ_h1 + R * bR * f_s1 * (1 - s) / θ / f_h1 / exp(γ) / π_ss - c_h1
r[2] = w2 * H2 - t_h2 + τ_h2 + R * bR * (1-f_s1) * (1 - s) / θ / (1 - f_h1) / exp(γ) / π_ss - c_h2
r[3] = c / (1 - θ) - θ * (f_h1 * c_h1 + (1 - f_h1) * c_h2) / (1 - θ) - c_s

r[4] = exp(γ) / (exp(γ) - h) / c_h1 - λ_h1
r[5] = exp(γ) / (exp(γ) - h) / c_h2 - λ_h2
r[6] = exp(γ) / (exp(γ) - h) / c_s  - λ_s

r[7] = (g - t + τ) / (1 - R / exp(γ) / π_ss) - bR

r[8] = exp(γ) * π_ss / disc * λ_s / (s * λ_s + (1 - s) * (f_s1 * λ_h1 + (1 - f_s1) * λ_h2)) - R

return r

end

init = [c_h10 c_h20 c_s0 λ_h10 λ_h20 λ_s0 bR0 R0]
ss = nlsolve(distSS, init)

solution = ss.zero

println(solution)


text:
module Messages

    using Log4jl: Message
    import Log4jl: format, formatted, parameters

    "Message with raw objects"
    type ObjectMessage <: Message
        message::Any
    end
    formatted(msg::ObjectMessage)  = string(msg.message)
    format(msg::ObjectMessage)     = formatted(msg)
    parameters(msg::ObjectMessage) = Any[msg.message]
    show(io::IO, msg::ObjectMessage) = print(io, "ObjectMessage[message=",msg.message,']')

    "Message handles everything as string."
    type SimpleMessage <: Message
        message::AbstractString
        SimpleMessage(msg::AbstractString, params...) = new(msg)
    end
    SimpleMessage(msg::Any) = ObjectMessage(msg)
    formatted(msg::SimpleMessage)  = msg.message
    format(msg::SimpleMessage)     = msg.message
    parameters(msg::SimpleMessage) = nothing
    show(io::IO, msg::SimpleMessage) = print(io, "SimpleMessage[message=",msg.message,']')

    "Message pattern contains placeholders indicated by '{}'"
    type ParameterizedMessage <: Message
        pattern::AbstractString
        params::Vector{Any}
        ParameterizedMessage(ptrn::AbstractString, params...) = new(ptrn, [params...])
    end
    ParameterizedMessage(msg::Any) = ObjectMessage(msg)
    function formatted(msg::ParameterizedMessage)
        length(msg.params) == 0 && return msg.pattern
        offs = map(ss->ss.offset, matchall(r"({})+",msg.pattern))
        @assert length(offs) == length(msg.params) "Pattern does not match parameters"
        sstart = 1
        sformatted = ""
        for (i,send) in enumerate(offs)
            sformatted *= msg.pattern[sstart:send]
            sformatted *= string(msg.params[i])
            sstart = send+3
        end
        sformatted * msg.pattern[sstart:end]
    end
    format(msg::ParameterizedMessage)     = msg.pattern
    parameters(msg::ParameterizedMessage) = msg.params
    show(io::IO, msg::ParameterizedMessage) =
        print(io, "ParameterizedMessage[pattern='",msg.pattern,"', args=",msg.params,']')

    "Message pattern contains 'printf' format string"
    type PrintfFormattedMessage <: Message
        pattern::AbstractString
        params::Vector{Any}
        PrintfFormattedMessage(ptrn::AbstractString, params...) = new(ptrn, [params...])
    end
    PrintfFormattedMessage(msg::Any) = ObjectMessage(msg)
    formatted(msg::PrintfFormattedMessage)  = @eval @sprintf($(msg.pattern), $(msg.params)...)
    format(msg::PrintfFormattedMessage)     = msg.pattern
    parameters(msg::PrintfFormattedMessage) = msg.params
    show(io::IO, msg::PrintfFormattedMessage) =
        print(io, "PrintfFormattedMessage[pattern='",msg.pattern,"', args=",msg.params,']')

    #TODO: MapMessage: XML, JSON, Dict
    #TODO: StructuredDataMessage: RFC 5424

    export ObjectMessage, SimpleMessage, ParameterizedMessage, PrintfFormattedMessage

end


text:
# # [Cantilever with a Uniform Load](@id cantilever)
#
# This example shows how to predict the behavior of a cantilever beam which is partially
# subjected to a uniform distributed load.
#
# ![](../assets/cantilever-drawing.svg)
#
#-
#md # !!! tip
#md #     This example is also available as a Jupyter notebook:
#md #     [`cantilever.ipynb`](@__NBVIEWER_ROOT_URL__/examples/cantilever.ipynb).
#-

using GXBeam, LinearAlgebra

nelem = 12

## create points
a = 0.3
b = 0.7
L = 1.0
n1 = n3 = div(nelem, 3)
n2 = nelem - n1 - n3
x1 = range(0, 0.3, length=n1+1)
x2 = range(0.3, 0.7, length=n2+1)
x3 = range(0.7, 1.0, length=n3+1)
x = vcat(x1, x2[2:end], x3[2:end])
y = zero(x)
z = zero(x)
points = [[x[i],y[i],z[i]] for i = 1:length(x)]

## index of endpoints for each beam element
start = 1:nelem
stop = 2:nelem+1

## create compliance matrix for each beam element
EI = 1e9
compliance = fill(Diagonal([0, 0, 0, 0, 1/EI, 0]), nelem)

## create assembly
assembly = Assembly(points, start, stop, compliance=compliance)

## set prescribed conditions (fixed right endpoint)
prescribed_conditions = Dict(
    nelem+1 => PrescribedConditions(ux=0, uy=0, uz=0, theta_x=0, theta_y=0,
    theta_z=0)
)

## create distributed load
q = 1000
distributed_loads = Dict(
    ielem => DistributedLoads(assembly, ielem; fz = (s) -> q) for ielem in
    n1+1:n1+n2
)

system, converged = static_analysis(assembly;
    prescribed_conditions = prescribed_conditions,
    distributed_loads = distributed_loads,
    linear = true)

state = AssemblyState(system, assembly;
    prescribed_conditions = prescribed_conditions)

#!jl nothing #hide

#
# We can construct the analytical solution for this problem by integrating from the free end of the beam and applying the appropriate boundary conditions.
#

## construct analytical solution
dx = 1e-6
x_a = 0.0:dx:L
q_a = (x) -> a <= x <= b ? -q : 0 # define distributed load
V_a = cumsum(-q_a.(x_a) .* dx) # integrate to get shear
M_a = cumsum(V_a .* dx) # integrate to get moment
slope_a = cumsum(M_a./EI .* dx) # integrate to get slope
slope_a .-= slope_a[end] # apply boundary condition
deflection_a = cumsum(slope_a .* dx) # integrate to get deflection
deflection_a .-= deflection_a[end] # apply boundary condition

## get elastic twist angle
theta_a = -atan.(slope_a)

## switch analytical system frame of reference
M_a = -M_a

#!jl nothing #hide

# Plotting the results reveals that the analytical and computational solutions show
# excellent agreement.

using Plots
#md using Suppressor #hide
pyplot()
#!jl nothing #hide

#-

#md @suppress_err begin #hide
#md local x #hide

## deflection plot
plot(
    xlim = (0.0, 1.0),
    xticks = 0.0:0.2:1.0,
    xlabel = "x (m)",
    ylabel = "Deflection (m)",
    grid = false,
    overwrite_figure=false
    )

x = [assembly.points[ipoint][1] + state.points[ipoint].u[1] for ipoint =
    1:length(assembly.points)]
deflection = [state.points[ipoint].u[3] for ipoint = 1:length(assembly.points)]
plot!(x_a, deflection_a, label="Analytical")
scatter!(x, deflection, label="GXBeam")
#!nb plot!(show=true)
#md savefig("../assets/cantilever-deflection.svg") #hide
#md closeall() #hide
#md end #hide
#md nothing #hide

#md # ![](../assets/cantilever-deflection.svg)

#-

#md @suppress_err begin #hide
#md local x #hide

## elastic twist plot (euler angle)
plot(
    xlim = (0.0, 1.0),
    xticks = 0.0:0.2:1.0,
    xlabel = "x (m)",
    ylabel = "Rotation Angle (rad)",
    grid = false,
    overwrite_figure=false
    )

x = [assembly.points[ipoint][1] + state.points[ipoint].u[1]
    for ipoint = 1:length(assembly.points)]
theta = [4*atan.(state.points[ipoint].theta[2]/4) for ipoint =
    1:length(assembly.points)]
plot!(x_a, theta_a, label="Analytical")
scatter!(x, theta, label="GXBeam")
#!nb plot!(show=true)
#md savefig("../assets/cantilever-twist.svg") #hide
#md closeall() #hide
#md end #hide
#md nothing #hide

#md # ![](../assets/cantilever-twist.svg)

#-

#md @suppress_err begin #hide
#md local x #hide

## bending moment plot
plot(
    xlim = (0.0, 1.0),
    xticks = 0.0:0.2:1.0,
    xlabel = "x (m)",
    ylabel = "Bending Moment (\$Nm\$)",
    grid = false,
    overwrite_figure=false
    )

x = [assembly.elements[ielem].x[1] + state.elements[ielem].u[1] for
    ielem = 1:length(assembly.elements)]
M = [state.elements[ielem].M[2] for ielem = 1:length(assembly.elements)]
plot!(x_a, M_a, label="Analytical")
scatter!(x, M, label="GXBeam")
#!nb plot!(show=true)
#md savefig("../assets/cantilever-moment.svg") #hide
#md closeall() #hide
#md end #hide
#md nothing #hide

#md # ![](../assets/cantilever-moment.svg)

#-

# Note that we could have easily performed a nonlinear analysis for this problem by setting
# `linear=false`.


text:
@testset "18.4sum.jl" begin
    @test four_sum([1, 0, -1, 0, -2, 2], 0) ==
        [(-2, -1, 1, 2), (-2, 0, 0, 2), (-1, 0, 0, 1)]
    @test four_sum(Int[], 0) == NTuple{4,Int}[]
    @test four_sum([-3, -1, 0, 2, 4, 5], 1) == [(-3, -1, 0, 5)]
    @test four_sum([1, 2, 2, 2, 2, 2, 2], 8) == [(2, 2, 2, 2)]
end


text:
# I/O functions

"""
Prints out Node content
"""
function Base.show(io::IO, n::SumNode)
    print(io,"+")
    for i=1:length(n.children)
        print(io, " $(n.children[i]) $(n.weights[i])")
    end
end

function Base.show(io::IO, n::ProductNode)
    print(io,"*")
    for i=1:length(n.children)
        print(io, " $(n.children[i])")
    end
end
function Base.show(io::IO, n::IndicatorFunction)
    print(io, "indicator $(n.scope) $(n.value)")
end
function Base.show(io::IO, n::CategoricalDistribution)
    print(io, "categorical $(n.scope)")
    for v in n.values
        print(io," $v")
    end
end
function Base.show(io::IO, n::GaussianDistribution)
    print(io, "indicator $(n.scope) $(n.mean) $(n.variance)")
end

"""
    summary(io::IO, spn::SumProductNetwork)

Print out information about the network `spn` to stream `io`
"""
function Base.summary(io::IO, spn::SumProductNetwork)
    len, lensum, lenprod, lenleaves = length(spn), length(sumnodes(spn)), length(productnodes(spn)), length(leaves(spn))
    lenvars = length(scope(spn))
    print(io, "Sum-Product Network with ", len, (len==1 ? " node" : " nodes"), " (", lensum, (lensum==1 ? " sum" : " sums"), ", ", lenprod, (lenprod==1 ? " product" : " products"), ", ", lenleaves, (lenleaves==1 ? " leaf" : " leaves"), ") and ", lenvars, (lenvars==1 ? " variable" : " variables"))
    # print(io, "Sum-Product Network with $(length(spn)) nodes ($(length(sumnodes(spn))) sums, $(length(productnodes(spn))) products, $(length(leaves(spn))) leaves), $(nparams(spn)) parameters, and $(length(scope(spn))) variables")
    # #println(io, summary(spn))
    # println(io, "Sum-Product Network with:")
    # println(io, "│\t$(length(spn)) nodes: $(length(sumnodes(spn))) sums, $(length(productnodes(spn))) products, $(length(leaves(spn))) leaves")
    # println(io, "│\t$(nparams(spn)) parameters")
    # println(io, "╰\t$(length(scope(spn))) variables")
    #println(io, "\tdepth = $(length(layers(spn)))")
end
# function Base.summary(spn::SumProductNetwork)
#     io = IOBuffer()
#     print(io, "Sum-Product Network with $(length(spn)) nodes ($(length(sumnodes(spn))) sums, $(length(productnodes(spn))) products, $(length(leaves(spn))) leaves), $(nparams(spn)) parameters, and $(length(scope(spn))) variables.")
#     String(take!(io))
# end

"""
    show(io::IO, spn::SumProductNetwork)

Print the nodes of the network `spn` to stream `io`
"""
function Base.show(io::IO, spn::SumProductNetwork)
    println(io, "SumProductNetwork(IOBuffer(\"\"\"# ", summary(spn))
    for (i, node) in enumerate(spn)
        println(io, i, " ", node)
    end
    print(io,"\"\"\"))")
    # print(io, summary(spn))
end
function Base.show(io::IO, ::MIME"text/plain", spn::SumProductNetwork)
    # recur_io = IOContext(io)
    recur_io = IOContext(io, :SHOWN_SET => spn)
    limit::Bool = get(io, :limit, false)
    if !haskey(io, :compact)
        recur_io = IOContext(recur_io, :compact => true)
    end
    summary(io, spn)
    print(io, ":")
    # print(io, "\n  1: ", spn[1])
    # println(io, "Sum-Product Network with:")
    # println(io, "│\t$(length(spn)) nodes: $(length(sumnodes(spn))) sums, $(length(productnodes(spn))) products, $(length(leaves(spn))) leaves")
    # println(io, "│\t$(nparams(spn)) parameters")
    # println(io, "╰\t$(length(scope(spn))) variables")
    if limit
        sz = displaysize(io)
        rows, cols = sz[1] - 3, sz[2]
        rows < 4   && (print(io, " …"); return)
        cols -= 5 # Subtract the width of prefix "  " and separator " : "
        cols < 12  && (cols = 12) # Minimum widths of 2 for id, 4 for value
        rows -= 1 # Subtract the summary

       # determine max id width to align the output, caching the strings
       ks = Vector{String}(undef, min(rows, length(spn)))
       vs = Vector{String}(undef, min(rows, length(spn)))
       keylen = 0
       vallen = 0
       for (i, n) in enumerate(spn)
           i > rows && break
           ks[i] = sprint(show, i, context=recur_io, sizehint=0)
           vs[i] = sprint(show, n, context=recur_io, sizehint=0)
           keylen = clamp(length(ks[i]), keylen, cols)
           vallen = clamp(length(vs[i]), vallen, cols)
       end
       if keylen > max(div(cols, 2), cols - vallen)
           keylen = max(cld(cols, 3), cols - vallen)
       end
    else
        rows = cols = typemax(Int)
    end
    for (i, node) in enumerate(spn)
        print(io, "\n  ")
        if i == rows < length(spn)
            print(io, rpad("⋮", keylen), " : ⋮")
            # print(io, rpad("⋮", 2))
            break
        end
        if limit
            key = rpad(Base._truncate_at_width_or_chars(ks[i], keylen, "\r\n"), keylen)
        else
            key = sprint(show, i, context=recur_io, sizehint=0)
        end
        print(recur_io, key)
        print(io, " : ")
        if limit
            val = Base._truncate_at_width_or_chars(vs[i], cols - keylen, "\r\n")
            print(io, val)
        else
            show(recur_io, n)
        end
    end
end

"""
    SumProductNetwork(filename::AbstractString; offset=0)::SumProductNetwork
    SumProductNetwork(io::IO=stdin; offset=0)::SumProductNetwork

Reads network from file. Assume 1-based indexing for node ids and values at indicator nodes. Set offset = 1 if these values are 0-based instead.
"""
function SumProductNetwork(filename::String; offset::Integer = 0)
    spn = open(filename) do file
        spn = SumProductNetwork(file, offset=offset)
    end
    spn
end
function SumProductNetwork(io::IO=stdin; offset::Integer = 0)
    # create dictionary of node_id => node (so they can be read in any order)
    nodes = Dict{UInt,Node}()
    # read and create nodes
    for line in eachline(io)
        # remove line break
        line = strip(line)
        # remove comments
        i = findfirst(isequal('#'), line)
        if !isnothing(i)
            line = line[1:i-1]
        end
        if length(line) > 0
            fields = split(line)
            if tryparse(Int, fields[1]) !== nothing
                nodeid = parse(Int,fields[1]) + offset
                nodetype = fields[2][1]
            else
                nodeid = parse(Int,fields[2]) + offset
                nodetype = fields[1][1]
            end
            if nodetype == '+'
                node = SumNode([ parse(Int,ch) + offset for ch in fields[3:2:end] ],
                               [ parse(Float64,w) for w in fields[4:2:end] ])
            elseif nodetype == '*'
                node = ProductNode([ parse(Int,id) + offset for id in fields[3:end] ])
            elseif nodetype == 'c'
                varid = parse(Int,fields[3]) + offset
                node = CategoricalDistribution(varid,[ parse(Float64,value) for value in fields[4:end] ])
            elseif nodetype == 'i' || nodetype == 'l'
                varid = parse(Int,fields[3]) + offset
                value = parse(Float64,fields[4]) + offset
                node = IndicatorFunction(varid,value)
            elseif nodetype == 'g'
                # TODO: read Gaussian leaves
                error("Reading of gaussian nodes is not implemented!")
            end
            nodes[nodeid] = node
        end
    end
    nodelist = Vector{Node}(undef, length(nodes))
    for (id,node) in nodes
        nodelist[id] = node
    end
    spn = SumProductNetwork(nodelist)
    sort!(spn) # ensure nodes are topologically sorted (with ties broken by bfs-order)
    spn
end

"""
    save(spn,filename,[offset=0])

Writes network spn to file. Offset adds constant to node instances (useful for translating to 0 starting indexes).
"""
function save(spn::SumProductNetwork, filename::String, offset = 0)
    open(filename, "w") do io
        for i = 1:length(spn)
            println(io,"$i $(spn[i])")
        end
    end
end

"""
    todot(io, son)

Prints out network structure in graphviz format
"""
function todot(io::IO, spn::SumProductNetwork)
    println(io, "digraph S {")
    for i = 1:length(spn)
	      if isa(spn[i], SumNode)
            if i == 1
		            println(io, "n$i [shape=circle,rank=source,style=filled,color=\"#fed434\",label=\"+\",margin=0.05];")
            else
		            println(io, "n$i [shape=circle,style=filled,color=\"#fed434\",label=\"+\",margin=0.05];")
            end
	      elseif isa(spn[i], ProductNode)
		        println(io, "n$i [shape=circle,style=filled,color=\"#b0db51\",label=\"×\",margin=0.05];")
	      elseif isa(spn[i], LeafNode)
		        println(io, "n$i [shape=circle,rank=sink,style=filled,color=\"#02a1d8\",label=\"X$(spn[i].scope)\",margin=0.05];")
	      end
        if !isa(spn[i],LeafNode)
            print(io, "n$i -> { ")
            for j in children(spn,i)
                print(io, "n$j; ")
            end
            println(io, "};")
        end
    end
    println(io, "}")
end


text:
# do not remove this first line
using PkgPage
using FranklinUtils
using Franklin

#
# Feel free to add whatever custom hfun_* or lx_*
# you might want to use in your site here
#

function env_column(e, _)
    md = Franklin.content(e)
    _, kwargs = lxargs(lxproc(e), "column")
    return _column(md; kwargs...)
end

function _column(md; colclass="col")
    return html("
") * md * html("
")
end

function lx_speakers(_, _)
    include("data/speakers.jl")
    sort!(speakers, by=x->split(x.name)[2]) # sort by surname
    result = ""
    for s in speakers
        result *= """
        \\begin{column}{colclass="col-sm-4"}
        \\figure{path="assets/speakers-img/$(s.imgName)", imgclass="img-organiser"}\\
        **[$(s.name)]($(s.link))**\\
        *$(s.org)*
        \\end{column}
        """
    end
    return result
end


text:
# Exports
export build_automaton

function build_automaton(Lambda::Set{T}) where {T <: Constraint}

    # Find constraint with maximum K value
    lambda_max_k = reduce(Lambda) do x, y
    x.k > y.k ? x : y end

    # Extract necessary parameters
    k = lambda_max_k.k
    n = _get_nbr_cons_hits(Lambda)

    # Decide IntType
    IntType = _word_type(Lambda, k)

    # IF we are not allowed to miss a single deadline, return a one vertex automaton
    if !is_satisfied(Lambda, consword(IntType(k)) << IntType(1))

        # Initialise Automaton
        initvertex              = WordVertex{Int8}(Int8(Hit))
        automaton               = Automaton{Int8}(initvertex)
        initvertex.hit          = initvertex
        automaton[Int8(Hit)]    = initvertex
        return automaton

    # ELSE IF we are allowed to miss infinite number of misses, return another one vertex automaton
    elseif _allmissespermitted(Lambda)
        # Initialise Automaton
        initvertex              = WordVertex{Int8}(Int8(Hit))
        automaton               = Automaton{Int8}(initvertex)
        initvertex.hit          = initvertex
        initvertex.miss         = initvertex
        automaton[Int8(Hit)]    = initvertex
        return automaton

    # ELSE add next vertex to transitions and continue
    else

        # Initialise automaton
        initvertex = WordVertex{IntType}(consword(IntType(n)))
        automaton = Automaton{IntType}(initvertex)

        for lambda in Lambda
            if typeof(lambda) == RowHitConstraint
                # create a curried function in order to be able to call the same function for both rowhit and all other constraints
                return _build_automaton_set!(Lambda,
                                            automaton,
                                            initvertex,
                                            (x1, x2, x3, x4) -> _initiate_vertex_rh_set!(x1, x2, x3, x4, initvertex.w, k))
            end # if
        end # function

        # Since we know something about the approximate size of the graph,
        # improve performance by setting the size to that automatically.
        sizehint!(automaton.data, binomial(k, n))
        return _build_automaton_set!(Lambda, automaton, initvertex, _initiate_vertex_set!)
    end
end # function

###########################
### Auxiliary Functions ###
###########################

#= Initiate vertex by first shifting vertex to acquire children.
 # Then extend automaton with children (if constraint set satisfied) =#
function _initiate_vertex_set!(Lambda::Set{T},
                               automaton::Automaton{S},
                               list::_UninitializedList{S},
                               vertex::WordVertex{S}) where {T <: Constraint, S <: Integer}
    (wm, wh) = _children_set(vertex)
    _extend_automaton_set!(Lambda, automaton, list, vertex, wm, wh)
end
function _initiate_vertex_rh_set!(Lambda::Set{T},
                                  automaton::Automaton{S},
                                  list::_UninitializedList{S},
                                  vertex::WordVertex{S},
                                  wi::S,
                                  k::R) where {T <: Constraint, S <: Integer, R <: Integer}
    (wm, wh) = _children_rh_set(vertex, wi)
    _extend_automaton_rh_set!(Lambda, automaton, list, vertex, wm, wh, k)
end

#= Return the children of vertex, disregarding satisfiability =#
function _children_set(vertex::WordVertex{T}) where {T <: Integer}
    wm = shift(vertex.w, Miss)
    wh = shift(vertex.w, Hit)
    return (wm, wh)
end
function _children_rh_set(vertex::WordVertex{T}, n::T) where {T <: Integer}
    wm = shift_rowhit(vertex.w, Miss, n)
    wh = shift_rowhit(vertex.w, Hit, n)
    return (wm, wh)
end

#= Extend automaton with the children (nm, nh) if they satisfy Lambda =#
function _extend_automaton_set!(Lambda::Set{T},
                                automaton::Automaton{S},
                                list::_UninitializedList{S},
                                vertex::WordVertex{S},
                                wm::S,
                                wh::S) where {T <: Constraint, S <: Integer}

    vertexM = nothing
    vertexH = nothing

    # Check if miss word satisfies constraints and whether it is initialized
    if is_satisfied(Lambda, wm)
        if !haskey(automaton, wm)
            vertexM         = WordVertex{S}(wm)
            automaton[wm]   = vertexM
            _push!(list, vertexM)
        else
            vertexM         = automaton[wm]
        end

        # Add miss path to vertex
        vertex.miss         = vertexM
    end

    # check if hit word is initialised
    if !haskey(automaton, wh)
        vertexH         = WordVertex{S}(wh)
        automaton[wh]   = vertexH
        _push!(list, vertexH)
    else
        vertexH         = automaton[wh]
    end

    # Add hit path to vertex
    vertex.hit = vertexH
end
#= Extend automaton with the children (nm, nh) if they satisfy Lambda =#
function _extend_automaton_rh_set!(Lambda::Set{T},
                                   automaton::Automaton{S},
                                   list::_UninitializedList{S},
                                   vertex::WordVertex{S},
                                   wm::S,
                                   wh::S,
                                   k::R) where {T <: Constraint, S <: Integer, R <: Integer}

    vertexH = nothing
    vertexM = nothing

    # Check if miss word satisfies constraints and whether it is initialized
    if is_satisfied(Lambda, wm)
        wm &= consword(S(k))            # Extracting bits to reduce size of graph
        if !haskey(automaton, wm)
            vertexM         = WordVertex{S}(wm)
            automaton[wm]   = vertexM
            _push!(list, vertexM)
        else
            vertexM         = automaton[wm]
        end

        # Add miss path to vertex
        vertex.miss         = vertexM
    end

    # check if hit word is initialised
    wh &= consword(S(k))            # Extracting bits to reduce size of graph
    if !haskey(automaton, wh)
        vertexH         = WordVertex{S}(wh)
        automaton[wh]   = vertexH
        _push!(list, vertexH)
    else
        vertexH = automaton[wh]
    end

    # Add hit path to vertex
    vertex.hit = vertexH
end

#= Aux function to build automaton =#
function _build_automaton_set!(Lambda::Set{T},
                               automaton::Automaton{S},
                               initvertex::WordVertex{S},
                               f_initiate_vertex!) where {T <: Constraint, S <: Integer}

    # First vertex initialisation
    list = _UninitializedList{S}()
    _push!(list, initvertex)
    automaton[initvertex.w] = initvertex

    # Add the remaining vertices to the while loop
    while !isnothing(list.head) # while we have elements to iterate
        # Extract first element to iterate
        vertex  = _pop!(list)

        # Add new feasible vertices to structures and lists
        f_initiate_vertex!(Lambda, automaton, list, vertex)
    end #while

    return automaton
end # function

@doc """
    build_automaton(Lambda::T) where {T <: Union{Constraint, Set{Constraint}}}

Creates a minimal weakly-hard automaton according to the (set of) weakly-hard
constraint(s) `Lambda`.

NOTE: The function handles both single constraints and sets of constraints.

# Examples
```julia-repl
julia> build_automaton(AnyHitConstraint(1, 3))
Automaton{Int64} with 3 vertices:
{
        WordVertex{Int64}(100 => ---, 001)
        WordVertex{Int64}(010 => 100, 001)
        WordVertex{Int64}(001 => 010, 001)
} with head: WordVertex{Int64}(1 => 10, 1)

julia> build_automaton(Set([AnyHitConstraint(1, 3), RowHitConstraint(2, 6)]))
Automaton{Int64} with 7 vertices:
{
        WordVertex{Int64}(001101 => 011010, 000011)
        WordVertex{Int64}(000110 => 001100, 001101)
        WordVertex{Int64}(011001 => ------, 000011)
        WordVertex{Int64}(011010 => ------, 110101)
        WordVertex{Int64}(001100 => ------, 011001)
        WordVertex{Int64}(000011 => 000110, 000011)
        WordVertex{Int64}(110101 => ------, 000011)
} with head: WordVertex{Int64}(11 => 110, 11)
```
""" build_automaton


text:
module PAMAuth

export authenticate

using Base.Libc
import Libdl

struct PamException <: Exception
    message::String
end

"""
    authenticate(username::AbstractString, password::AbstractString)

Attempts to authenticate a user with the local system via PAM under the given username and password.

# Examples
```julia
using PAM

print("Username: ")
username = readline()
password = read(Base.getpass("Password"), String)

valid = authenticate(username, password)
# Valid is true if authentication was successful
```
"""
function authenticate(username::AbstractString, password::AbstractString)
    libpamjl = Libdl.dlopen(normpath("$(@__DIR__)/../build/libpamjl.so"))
    pamjl_authenticate = Libdl.dlsym(libpamjl, :pamjl_authenticate)
    ccall(pamjl_authenticate, Cint, (Ptr{Cchar}, Ptr{Cchar}), username, password) == 0
end

end # module


text:
# Single branch conditional
if 1 == 1
    print("Of course 1==1 silly")
end

max_stars_repo_name: elucian/julia

max_issues_repo_name: elucian/julia

max_forks_repo_name: elucian/julia


text:
using Gadfly

n = 10
plot(x=rand(n), y=rand(n), color=e.^rand(n),
     Scale.color_log(minvalue=1))


text:
isdefinite(n::Number) = !isnan(n) && !isinf(n)

for n in (1, 1//1, 1.0, 1im, 0)
    d = n / 0
    println("Dividing $n by 0 ", isdefinite(d) ? "results in $d." : "yields an indefinite value ($d).")
end

max_stars_repo_name: ethansaxenian/RosettaDecode

max_issues_repo_name: ethansaxenian/RosettaDecode

max_forks_repo_name: ethansaxenian/RosettaDecode


text:
#=
Peak fitting example using package LUComptonScatter

Author: Miles Cochran-Branson
Date: 2/15/22
=#

using LUComptonScatter, Plots

#function determined from calibration
f(x) = 0.669462 * x + 55.3071

data = read_Txt("examples/exampleData.Txt")
x0, y0 = get_xy_data(data, 100, f)

# Plot raw data
plot(x0, y0, label = "data", xlabel = "Energy [keV]", ylabel = "Counts", legend=:topleft, dpi=400, size=(350, 350))

# extract smoothed curve and peak data
peaks, ySmooth = peak_parameters(x0, y0, 25, 1000.0, 0.5)

#plot this new data
plt = plot(x0, ySmooth, label = "data", xlabel = "Energy [keV]", ylabel = "Counts", legend=:topleft, dpi=400, size=(350, 350))
for i in 1:length(peaks)
    scatter!(x0[peaks[i]], ySmooth[peaks[i]], label = "Peak Data $i")
end
display(plt)

##extract value of peak
result, error = fit_to_gauss(x0, y0, peaks[1])
result[3], error[3]


text:
## path to TRPL folder --> set it before running the code!
path_to_trpl = "D:\\Projects\\TRPL" # some path to the TRPL folder (just an example here)

## activate package environment
cd(@__DIR__)
using Pkg
Pkg.activate("Project.toml")

## loading packages
using CSV, DataFrames

## directory
cd(joinpath(path_to_trpl,"Data\\Real"))

## output results? (true or false)
out = false

## parameters measurement data
name_measurement = "sample1"

if name_measurement == "sample1"
    N_DA = 5e16 # sample1: N_DA = 5e16
elseif name_measurement == "sample2"
    N_DA = 1e17 # sample2: N_DA = 1e17
end

## load measurement data
data = Matrix(CSV.read(name_measurement*".csv",DataFrame,skipto=2))

## laser power P
P = round.(data[1,2:end],sigdigits=3)

## transmission coefficient (calculated in transm_calc.jl)
T = 0.4292

## other experimental parameters
rep_rate = 200e3
lambda   = 520*1e-9
fwhm     = 205*1e-6
d        = 2000*1e-9

## physical constants
h = 6.62607015e-34
c = 299792458

## main calculations
E        = P./rep_rate
N        = T*E/(h*c)*lambda
V        = pi*fwhm^2*d/4
n        = N/V * 1e-6
fwhm_int = 0.754 # from fwhm_int.jl script
eta_0    = round.(fwhm_int*n/N_DA,sigdigits=3)

## print result
println("P = $(P)")
println("η₀ = $(eta_0)")

## write output
if out
    cd(joinpath(path_to_trpl,"Results\\Real"))
    file_name = "eta0_$(name_measurement).txt"
    if isfile(file_name)
        rm(file_name)
    end
    open(file_name,"w") do io
        println(io,"\nP = $(P)")
        println(io,"\nη0 = $(eta_0)")
    end
end

max_stars_repo_name: MaxGrossmann/TRPL

max_issues_repo_name: MaxGrossmann/TRPL

max_forks_repo_name: MaxGrossmann/TRPL


text:
# ************************************************************************
# FAUST Architecture File
# Copyright (C) 2021 GRAME, Centre National de Creation Musicale
# ---------------------------------------------------------------------

# This is sample code. This file is provided as an example of minimal
# FAUST architecture file. Redistribution and use in source and binary
# forms, with or without modification, in part or in full are permitted.
# In particular you can create a derived work of this FAUST architecture
# and distribute that work under terms of your choice.

# This sample code is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
# ************************************************************************

# Architectures files
#include("audio.jl")

# Using PortAudio for audio rendering
using PortAudio

mutable struct portaudio <: audio
    portaudio(buffer_size::Int, sample_rate::Int) = begin
        driver = new()
        driver.buffer_size = buffer_size
        driver.sample_rate = sample_rate
        driver
    end
    dsp::dsp
    buffer_size::Int
    sample_rate::Int
    name::String
end

function init!(driver::portaudio, name::String, dsp::dsp)
    driver.name = name
    driver.dsp = dsp
    init!(dsp, Int32(driver.sample_rate))
end

function run(driver::portaudio)
    PortAudioStream(1, 2) do stream
        outputs = zeros(FAUSTFLOAT, driver.buffer_size, getNumOutputs(driver.dsp))
        while true
            inputs = convert(Matrix{FAUSTFLOAT}, read(stream, driver.buffer_size))
            compute!(driver.dsp, Int32(driver.buffer_size), inputs, outputs)
            write(stream, outputs)
        end
    end
end


text:
# https://github.com/JuliaLang/julia/issues/28679
function rm_stdout(f::Function)
    nullfile = @static Sys.iswindows() ? "nul" : "/dev/null"
    open(nullfile, "w") do io
        redirect_stdout(f, io)
    end
end

rm_stdout() do
    include(joinpath(@__DIR__, "common.jl"))
end

algo_idx = parse(Int, ARGS[1])
problem_idx = parse(Int, ARGS[2])

algo = algos[algo_idx]
prob = problems[problem_idx]

optval, elapsed = algo.f(prob)
println(optval)
println(elapsed)


text:
using DifferentialEquations
using LSODA
using Sobol
using JLD2
using ArgParse

parser = ArgParseSettings(allow_ambiguous_opts=false)
@add_arg_table! parser begin
    "output"
        help = "output JLD2"
        required = true
    "form"
        help = "noise form: add or mul"
        required = true
    "target"
        help = "noise target: rna or both"
        required = true
    "initstate"
        help = "initial state: hi or lo"
        required = true
    "signals"
        help = "signal values"
        required = true
    "--cells"
        help = "number of cells to simulate"
        default = 500
        arg_type = Int
    "--noises"
        help = "noise values"
        default = "0.5:0.5:1.5"
end
args = parse_args(parser)

function model_repressilator(du, u, p, t)
    a1, a2, a3, a, n, dR, dP, y, σR, σP, τ = p

    M1, M2, M3, P1, P2, P3 = u[1], u[2], u[3], u[4], u[5], u[6]

    du[1] = a + a1 / (1 + P3^n) - dR * M1 + u[7]
    du[2] = a + a2 / (1 + P1^n) - dR * M2 + u[8]
    du[3] = a + a3 / (1 + P2^n) - dR * M3 + u[9]

    du[4] = y * M1 - dP * P1 + u[10]
    du[5] = y * M2 - dP * P2 + u[11]
    du[6] = y * M3 - dP * P3 + u[12]

    du[7:9] = - u[7:9]/τ
    du[10:12] = - u[10:12]/τ
end

function stoch_repressilator_add(du,u,p,t)
    a1, a2, a3, a, n, dR, dP, y, σR, σP, τ = p
    du[1:6] .= 0
    du[7:9] .= sqrt((2.0*σR^2.0/τ))
    du[10:12] .= sqrt((2.0*σP^2.0/τ))
end

function stoch_repressilator_mul(du,u,p,t)
    a1, a2, a3, a, n, dR, dP, y, σR, σP, τ = p
    du[1:6] .= 0
    du[7:9] = sqrt(2.0*σR^2.0/τ) .* u[1:3]
    du[10:12] = sqrt(2.0*σP^2.0/τ) .* u[4:6]
end

p_cycle = [
    7, 5, 5, # a1, a2, a3
    0.01, # a
    2, # n
    1, 1, # dR, dP
    1, # y
    1, 1, # σR, σP
    0.1 # τ
];

function filter_sols(sols)
    nvars, ntps, nsols = size(sols)
    good_sols = []
    for i in 1:nsols
        if isassigned(sols[i][1,:], ntps) == true
            push!(good_sols, i)
        end
    end
    println("Obtained $(length(good_sols))/$nsols converged solutions.")
    sols[:, :, good_sols]
end

function sim_pop_ic(model, model_stoch, u0, pars; popsize=1000, selectR::String="lo", saveat=1.0)
    selectR == "lo" || selectR == "hi" || error("invalid selectR")
    prob_ode = ODEProblem(model, u0, (0., 5000.), pars);
    prob_sde = SDEProblem(model, model_stoch, u0, tspan, pars);
    seq = SobolSeq(3)
    init_rnas = [next!(seq) .* 5 .+ [0.01, 0, 0] for _ in 1:100]
    function prob_func_cells_init(prob,i,repeat)
        init_rna = init_rnas[i]
        remake(prob,u0=[init_rna; init_rna; u0[7:end]])
    end
    prob_cells_init = EnsembleProblem(prob_ode, prob_func=prob_func_cells_init)
    @time sim_init = solve(prob_cells_init, lsoda(), saveat=saveat, dt=0.01, reltol=1e-6, abstol=1e-6, trajectories=length(init_rnas));
    sol_init = sim_init[:, 4000:end, :] # The selection of state is sensitive to the start index
    if selectR == "lo"
        sol_init[isnan.(sol_init)] .= Inf
        idsel = argmin(sol_init[1, :, :])
    else
        sol_init[isnan.(sol_init)] .= -Inf
        idsel = argmax(sol_init[1, :, :])
    end
    u0_s = sol_init[:, idsel[1], idsel[2]]
    println("Starting stochastic simulation with init cond ", u0_s)
    function prob_func_cells(prob, i, repeat)
        remake(prob,u0=u0_s)
    end
    prob_cells = EnsembleProblem(prob_sde, prob_func=prob_func_cells)
    @time sim = solve(prob_cells, STrapezoid(), saveat=saveat, maxiter=1e12, dt=0.0001, trajectories=popsize);
    filter_sols(sim)
end

stochs = Dict("add" => stoch_repressilator_add, "mul" => stoch_repressilator_mul)

stoch = stochs[args["form"]]
tspan = (0., 1000.)
u0 = [1.; zeros(11)]
noises = eval(Meta.parse(args["noises"]))
signals = eval(Meta.parse(args["signals"]))

all_results = map(noises) do noise
    results = map(signals) do signal
        p = copy(p_cycle)
        p[1] = signal
        p[9] = noise
        p[10] = args["target"] == "both" ? noise : 0.
        result = sim_pop_ic(model_repressilator, stoch, u0, p; popsize=args["cells"], selectR=args["initstate"])
        signal => result[1:6, :, :]
    end
    noise => Dict(results)
end
save(args["output"], "results", Dict(all_results))


text:
using MacroTools: combinedef

_gradtuple(::Nothing) = nothing
_gradtuple(x::Tuple) = (nothing, x...)
_gradtuple(x) = error("Gradient $x should be a tuple")

function gradm(f, T, args, Ts, body, mut)
  args = esc.(args)
  Ts = esc.(Ts)
  pushfirst!(args, :($(esc(:__context__))::Context), :($f::$T))
  body = quote
    Base.@_inline_meta
    y, back = $(esc(body))
    $(mut ? nothing : :(back2(::Nothing) = nothing))
    # return needed for type inference
    back2(Δ) = return _gradtuple(back(Δ))
    y, back2
  end
  :(Zygote._forward($(args...)) where $(Ts...) = $body; nothing)
end

_gradtuple_kw(::Nothing) = nothing
_gradtuple_kw(x::Tuple) = (nothing, nothing, nothing, x...) # kwfunc, kws, func, args...
_gradtuple_kw(x) = error("Gradient $x should be a tuple")
_untuple_kw(::Nothing) = nothing
_untuple_kw(x::Tuple) = Base.tail(Base.tail(x))

function gradm_kw(f, T, args, Ts, body, mut)
  kws = popfirst!(args)
  Ts = esc.(Ts)
  kT = :(Core.kwftype($T))
  kwargs = [:($(esc(:__context__))::Context), :(::$kT), :kw, f, esc.(args)...]
  kw_wrappers = map(kws.args) do kw
    kw isa Symbol && return :($(esc(kw)) = kw.$kw)
    isexpr(kw, :...) && return :($(esc(kw.args[1])) = kw)
    k, v = kw.args
    :($(esc(k)) = haskey(kw, $(Expr(:quote, k))) ? kw.$k : $(esc(v)))
  end
  body = quote
    Base.@_inline_meta
    $(kw_wrappers...)
    y, back = $(esc(body))
    $(mut ? nothing : :(back2(::Nothing) = nothing))
    # return needed for type inference
    back2(Δ) = return _gradtuple_kw(back(Δ))
    y, back2
  end
  quote
    Zygote._forward($(kwargs...)) where $(Ts...) = $body
    function Zygote._forward(cx::Context, f::$T, $(esc.(args)...)) where $(Ts...)
      y, back = _forward(cx, Core.kwfunc(f), NamedTuple(), f, $(esc.(namify.(args))...)) # TODO unnamed arguments
      return y, Δ -> _untuple_kw(back(Δ))
    end
    nothing
  end
end

function gradm(ex, mut = false)
  @capture(shortdef(ex), (name_(args__) = body_) |
                         (name_(args__) where {Ts__} = body_)) || error("Need a function definition")
  iskw = length(args) > 1 && isexpr(args[1], :parameters)
  name, T = isexpr(name, :(::)) ?
    (length(name.args) == 1 ? (:_, esc(name.args[1])) : esc.(name.args)) :
    (:_, :(typeof($(esc(name)))))
  Ts == nothing && (Ts = [])
  return (iskw ? gradm_kw : gradm)(name, T, args, Ts, body, mut)
end

macro grad(ex)
  gradm(ex)
end

macro grad!(ex)
  gradm(ex, true)
end

macro nograd(ex)
  isexpr(ex, :tuple) || (ex = Expr(:tuple, ex))
  blk = :(;)
  for f in ex.args
    push!(blk.args, :(@inline Zygote._forward(::Context, ::typeof($(esc(f))), args...) = $(esc(f))(args...), Δ -> nothing))
  end
  return blk
end


text: getevennumbers(arr) = filter(iseven,arr)

max_stars_repo_name: UlrichBerntien/Codewars-Katas

max_issues_repo_name: UlrichBerntien/Codewars-Katas

max_forks_repo_name: UlrichBerntien/Codewars-Katas


text:
using Requests
using Base.Test

println("It should respond with index.html")
response = get("http://0.0.0.0:8084")

@test response.status == 200
@test response.data == "<body>hello</body>"

println("It should serve a file")
response = get("http://0.0.0.0:8084/script.js")

@test response.status == 200
@test response.data == "alert;\n"


text:
using Test, YaoFlux
using Yao, Flux
using FiniteDifferences
using LinearAlgebra: norm
using Random

nonzero = xs -> findall(x -> abs2(x) > ϵ, reshape(xs, length(xs)))
ϵ = 1e-7

@testset "Circuit model" begin
    circ = chain(3, put(1 => H), cnot(1, 2), cnot(1, 3), put(2 => Rx(0.)))

    varcirc = Circuit(circ)
    out = varcirc([])
    @test nonzero(out) == [1, 8]

    ps = params(varcirc)
    @test length(ps) == 1
    @test length(ps[1]) == nparameters(circ)
    @test nparameters(varcirc) == nparameters(circ)

    target_ψ = √2/2 * [1, 0, 0, 0, 0, 0, 0, 1]
    loss(x) = norm(varcirc(x) - target_ψ)
    grad = Flux.gradient(ps) do
        train_loss = loss([])
    end
    @test all(grad[varcirc.θ] .< ϵ)
end

@testset "apply_circuit computes gradients correctly" begin
    Random.seed!(3)
    circ = chain(3, put(1 => Rx(0.)), cnot(1, 2), cnot(1, 3), put(2 => Rx(0.)))
    reg = ArrayReg(rand(8, 1) + 1im * rand(8, 1))
    normalize!(reg)
    params = rand(2)
    f = θ -> to_vec(YaoFlux.apply_circuit(reg, θ; circ))[1]
    fdm_jac = jacobian(central_fdm(5, 1), f, params)[1]
    exact_jac = Flux.jacobian(f, params)[1]
    @test norm(fdm_jac - exact_jac) < ϵ
end

max_stars_repo_name: lmondada/YaoFlux

max_issues_repo_name: lmondada/YaoFlux

max_forks_repo_name: lmondada/YaoFlux


text:
"""
    AbstractTuples

Module that provides additional functionally for `Tuple` and `NamedTuple`

### See also
[`AbstractTuple`](@ref), [`AbstractNTuple`](@ref), [`NamedNTuple`](@ref)
"""
module AbstractTuples

export
    AbstractTuple,
    NamedNTuple, AbstractNTuple,
    TypeTuple, NamedTypeTuple, AbstractTypeTuple,
    typetuple_to_type, typetuple,
    replace_in_tuple,
    convert_to_tuple

# ======================================================
# Abstract tuple types
# ======================================================

"""
    AbstractTuple{<: Tuple}

A type `Union` of a `Tuple` and a `NamedTuple`.

Useful for methods that take either a `Tuple` or a `NamedTuple` of some specific
types for parameters.

Note that like `NamedTuple` but unlike `Tuple`, AbstractTuple is not covariant.
I.e. it does not hold that `AbstractTuple{Tuple{Int}} <: AbstractTuple{Tuple{Integer}}`

### See also
[`AbstractNTuple`](@ref), [`Tuple`](@ref), [`NTuple`](@ref)

### Examples

```julia
julia> (1, "xyz") isa AbstractTuple{Tuple{Int, String}}
true

julia> (a=1, b="xyz") isa AbstractTuple{Tuple{Int, String}}
true

julia> (a=1, b="xyz") isa AbstractTuple
true
```
"""
const AbstractTuple{T <: Tuple} = Union{T, NamedTuple{S, T} where S}

"""
    NamedNTuple{N, T}

A `NamedTuple` with `N` arguments of the same type.

Analogue to `NTuple` for `Tuple` but as `NamedTuple` is not covariant so is
also `NamedNTuple` not covariant.

### See also
[`AbstractNTuple`](@ref), [`NTuple`](@ref), [`NamedTuple`](@ref)

### Examples

```julia
julia> (a = 1, b = 2) isa NTuple{2, Int}
true

julia> (a = 1, b = 2) isa NTuple{2}
true

julia> (a = 1, b = "xyz") isa NTuple{2, Int}
false
```
"""
const NamedNTuple{N, T} = NamedTuple{S, NTuple{N, T}} where S

"""
AbstractNTuple{N, T}

A Union of an `NTuple` and a `NamedNTuple`.
### See also
[`NTuple`](@ref), [`NamedNTuple`](@ref)

### Examples

```julia
julia> (1, 2, 3) isa AbstractNTuple{3, Int}
true

julia> (a=1, b=2, c=3) isa AbstractNTuple{3, Int}
true

julia> (1, 2) isa AbstractNTuple{3, Int}
false

julia> (a=1, b=2, c="xyz") isa AbstractNTuple{3, Int}
false
```
"""
const AbstractNTuple{N, T} = Union{NTuple{N, T}, NamedNTuple{N, T}}

# ======================================================
# Type tuples
# ======================================================

const TypeTuple = Tuple{Vararg{Type}}
const NamedTypeTuple = NamedTuple{S, <:Tuple{Vararg{Type}}} where {S}
const AbstractTypeTuple = Union{TypeTuple, NamedTypeTuple}

"""
    typetuple_to_type(t::Tuple)
    typetuple_to_type(t::NamedTuple)

Convert a Tuple or Named tuple containing types to a Tuple or NamedTuple with
these types.

### Examples

```julia
julia> typetuple_to_type((Int64, String))
Tuple{Int64,String}

julia> typetuple_to_type((a = String, b = Float64))
NamedTuple{(:a, :b),Tuple{String,Float64}}
```
"""
function typetuple_to_type end

typetuple_to_type(tup::TypeTuple) = Tuple{ (T for T in tup)... }
typetuple_to_type(tup::NamedTypeTuple) =
    NamedTuple{ Tuple(typeof(tup).names), Tuple{ (T for T in tup)... }}

"""
    typetuple(::Type{<:Tuple})
    typetuple(::Type{<:NamedTuple})

Convert a Tuple or NamedTuple type to a Tuple or NamedTuple
containing the types of the tuple.

### Examples

```julia
julia> typetuple(Tuple{Int64, String})
(Int64, String)

julia> typetuple(typeof(a = "first", b = 2.0))
(a = String, b = Float64)
```
"""
function typetuple end

function typetuple(TT::Type{<:Tuple})

    return Tuple(TT.types)
end

function typetuple(TT::Type{<:NamedTuple})

    return NamedTuple{Tuple(TT.names)}(TT.types)
end

# ======================================================
# replace_in_tuple
# ======================================================

"""
    replace_in_tuple(t::Tuple, key, value)
    replace_in_tuple(t::NamedTuple, key, value)

Return a copy of `t` where the value specified by `key` has been replaced by `value`.
Does not add additional fields if key does not exist.
"""
function replace_in_tuple end

replace_in_tuple(t::Tuple, key, v) = Base.setindex(t, v, key)
replace_in_tuple(nt::NamedTuple, key::Symbol, v) =
    replace_in_tuple(nt, Base.fieldindex(typeof(nt), key), v)

function replace_in_tuple(nt::NamedTuple, key::Integer, v)

    T = typeof(nt)
    return T(replace_in_tuple(Tuple(nt), key, v))
end

# ======================================================
# convert_to_tuple
# ======================================================

# TODO document behaviour
convert_to_tuple(T::Type{<:NamedTuple}, nt::NamedTuple) = convert(T, nt)
convert_to_tuple(T::Type{<:Tuple}, t::Tuple) = convert(T, t)
convert_to_tuple(T::Type{<:NamedTuple}, t::Tuple) = T((t))

end # Module


text:
# Copyright (c) 2022, Eric Sabo
# All rights reserved.
#
# This source code is licensed under the BSD-style license found in the
# LICENSE file in the root directory of this source tree.

include("linearcode.jl")
include("quantumcode.jl")
include("trellis.jl")

#############################
         # General
#############################

# sort print on type here
function show(io::IO, W::WeightEnumerator)
    if length(W.polynomial[1]) == 3
        for (i, term) in enumerate(W.polynomial)
            if term == [1, 0, 0]
                print("1")
            end
            if term[1] != 1
                print(io, term[1])
            end
            if term[2] != 0
                if term[2] != 1
                    print(io, "x^", term[2])
                else
                    print(io, "x")
                end
            end
            if term[3] != 0
                if term[3] != 1
                    print(io, "y^", term[3])
                else
                    print(io, "y")
                end
            end
            if i != length(W.polynomial)
                print(io, " + ")
            # else
            #     print(io, "\n")
            end
        end
    else
        for (i, term) in enumerate(W.polynomial)
            if term == [1, 0, 0, 0]
                print("1")
            end
            if term[1] != 1
                print(io, term[1])
            end
            if term[2] != 0
                if term[2] != 1
                    print(io, "x^", term[2])
                else
                    print(io, "x")
                end
            end
            if term[3] != 0
                if term[3] != 1
                    print(io, "y^", term[3])
                else
                    print(io, "y")
                end
            end
            if term[4] != 0
                if term[4] != 1
                    print(io, "z^", term[4])
                else
                    print(io, "z")
                end
            end
            if i != length(W.polynomial)
                print(io, " + ")
            # else
            #     print(io, "\n")
            end
        end
    end
end

function _weightenumerator(G::fq_nmod_mat)
    E = base_ring(G)
    n = ncols(G)
    weightdist = zeros(Int, n + 1)
    for iter in Base.Iterators.product([0:(Int64(characteristic(E)) - 1) for _ in 1:nrows(G)]...)
        row = E(iter[1]) * G[1, :]
        for r in 2:nrows(G)
            if !iszero(iter[r])
                row += E(iter[r]) * G[r, :]
            end
        end
        wt = weight(row)
        weightdist[wt + 1] += 1
    end

    # weight distribution to weight enumartor
    poly = Vector{Vector{Int64}}()
    for i in 1:n + 1
        if !iszero(weightdist[i])
            push!(poly, [weightdist[i], i - 1, n - i + 1])
        end
    end
    return WeightEnumerator(poly, "Hamming")
end

#############################
        # Classical
#############################

#############################
    # Weight Enumerators
#############################

# speed boost by reducing operation on poly[j][4] at the cost of compile size
function _reducepolyC(poly::Vector{Vector{Int64}})
    reducedpoly = Vector{Vector{Int64}}()
    processed = trues(length(poly))
    for (i, term) in enumerate(poly)
        if processed[i]
            for j in (i + 1):length(poly)
                if processed[j]
                    if term[2] == poly[j][2] && term[3] == poly[j][3]
                        term[1] += poly[j][1]
                        processed[j] = false
                    end
                end
            end
            push!(reducedpoly, term)
            processed[i] = false
        end
    end
    return reducedpoly
end

# for some reason this is broken?
# sorting somehow changes all the elements
function islessPolyC(a::Vector{Int64}, b::Vector{Int64})
    if isless(a[2:3], b[2:3])
        return true
    elseif a[2:3] == b[2:3]
        if isless(a[1], b[1])
            return true
        else
            return false
        end
    else
        return false
    end
end

# set to struct
function weightenumeratorC(T::Trellis, cleanV::Bool=true)
    V = vertices(T)
    E = edges(T)
    for i in 2:length(V)
        for (j, v) in enumerate(V[i])
            outer = Vector{Vector{Int64}}()
            for e in E[i - 1][j]
                inner = deepcopy(V[i - 1][e.outvertex].polynomial)
                for k in e.label
                    if iszero(k)
                        for term in inner
                            term[3] += 1
                        end
                    else
                        for term in inner
                            term[2] += 1
                        end
                    end
                end
                append!(outer, inner)
            end
            v.polynomial = _reducepolyC(outer)
        end
    end
    W = WeightEnumerator(sort!(V[end][1].polynomial, lt=islessPolyC), "Hamming")

    if cleanV
        for i in 2:length(V)
            for v in V[i]
                v.polynomial = [[0, 0, 0]]
            end
        end
    end
    return W
end

function weightenumerator(C::AbstractLinearCode, alg::String="trellis", sect::Bool=false)
    !ismissing(C.weightenum) && return C.weightenum
    alg ∈ ["trellis", "bruteforce"] || error("Algorithm `$alg` is not implemented in weightdistribution.")
    if alg == "trellis"
        T = syndrometrellis(C, sect)
        C.weightenum = weightenumeratorC(T, true)
    elseif alg == "bruteforce"
        C.weightenum = _weightenumerator(generatormatrix(C))
    end
    return C.weightenum
end

# MAGMA returns
# [ <0, 1>, <4, 105>, <6, 280>, <8, 435>, <10, 168>, <12, 35> ]
function weightdistribution(C::AbstractLinearCode, alg::String="trellis", sect::Bool=false)
    ismissing(C.weightenum) && weightenumerator(C, alg, sect)
    temp = zeros(Int64, 1, length(C) + 1)
    for term in C.weightenum.polynomial
        temp[term[2] + 1] = term[1]
    end
    return temp
end

#############################
     # Minimum Distance
#############################

"""
    minimumdistance(C::AbstractLinearCode, alg::String="trellis", sect::Bool=false)

Return the minimum distance of the linear code if known, otherwise computes it
using the algorithm of `alg`. If `alg = "trellis"`, the sectionalization flag
`sect` can be set to true to further compactify the reprsentation.
"""
function minimumdistance(C::AbstractLinearCode, alg::String="trellis", sect::Bool=false)
    !ismissing(C.d) && return C.d
    weightenumerator(C, alg, sect)
    # first element is always the identity y^n, want power of first x term
    C.d = C.weightenum.polynomial[2][2]
    return C.d
end

#############################
         # Quantum
#############################

#############################
    # Weight Enumerators
#############################

function _reducepolyQ(poly::Vector{Vector{Int64}})
    reducedpoly = Vector{Vector{Int64}}()
    processed = trues(length(poly))
    for (i, term) in enumerate(poly)
        if processed[i]
            for j in (i + 1):length(poly)
                if processed[j]
                    if term[2] == poly[j][2] && term[3] == poly[j][3] && term[4] == poly[j][4]
                        term[1] += poly[j][1]
                        processed[j] = false
                    end
                end
            end
            push!(reducedpoly, term)
            processed[i] = false
        end
    end
    return reducedpoly
end

# for some reason this is broken?
# sorting somehow changes all the elements
function islessPolyQ(a::Vector{Int64}, b::Vector{Int64})
    if isless(a[2:4], b[2:4])
        return true
    elseif a[2:4] == b[2:4]
        if isless(a[1], b[1])
            return true
        else
            return false
        end
    else
        return false
    end
end

# right now this only makes sense for qubit codes
# for nonqubit, can add coeff(l, ⋅) instead of 1 but then the concept of Y
# is itself fuzzy
# can just define a qudit version with X, Z tracking only

# does this actually modify T in a way I don't want to keep? should always cleanV?
function Pauliweightenumerator(T::Trellis, cleanV::Bool=true)
    V = vertices(T)
    E = edges(T)
    Int(characteristic(base_ring(E[1][1][1].label))) > 2 &&
        error("The PWE is currently only valid for qubit codes. Run the standard Hamming weight enumerator for qudit codes.")
    for i in 2:length(V)
        for (j, v) in enumerate(V[i])
            outer = Vector{Vector{Int64}}()
            for e in E[i - 1][j]
                inner = deepcopy(V[i - 1][e.outvertex].polynomial)
                for k in e.label
                    if coeff(k, 0) == 1 && iszero(coeff(k, 1))
                        if e.sign == 1
                            for term in inner
                                if term[2] < 0 || term[3] < 0 || term[4] < 0
                                    term[2] -= 1
                                else
                                    term[2] += 1
                                end
                            end
                        else
                            for term in inner
                                term[2] *= -1
                                term[3] *= -1
                                term[4] *= -1
                                if (term[2] < 0 || term[3] < 0 || term[4] < 0) || (term[2] == 0 && term[3] == 0 && term[4] == 0)
                                    term[2] -= 1
                                else
                                    term[2] += 1
                                end
                            end
                        end
                    elseif coeff(k, 0) == 1 && coeff(k, 1) == 1
                        if e.sign == 1
                            for term in inner
                                if term[2] < 0 || term[3] < 0 || term[4] < 0
                                    term[3] -= 1
                                else
                                    term[3] += 1
                                end
                            end
                        else
                            for term in inner
                                term[2] *= -1
                                term[3] *= -1
                                term[4] *= -1
                                if (term[2] < 0 || term[3] < 0 || term[4] < 0) || (term[2] == 0 && term[3] == 0 && term[4] == 0)
                                    term[3] -= 1
                                else
                                    term[3] += 1
                                end
                            end
                        end
                    elseif iszero(coeff(k, 0)) && coeff(k, 1) == 1
                        if e.sign == 1
                            for term in inner
                                if term[2] < 0 || term[3] < 0 || term[4] < 0
                                    term[4] -= 1
                                else
                                    term[4] += 1
                                end
                            end
                        else
                            for term in inner
                                term[2] *= -1
                                term[3] *= -1
                                term[4] *= -1
                                if (term[2] < 0 || term[3] < 0 || term[4] < 0) || (term[2] == 0 && term[3] == 0 && term[4] == 0)
                                    term[4] -= 1
                                else
                                    term[4] += 1
                                end
                            end
                        end
                    else
                        continue
                    end
                end
                append!(outer, inner)
            end
            v.polynomial = _reducepolyQ(outer)
        end
    end
    # sort is somehow broken
    # println(V[end][1].polynomial)
    W = WeightEnumerator(sort!(V[end][1].polynomial, lt=islessPolyQ), "Pauli")
    # W = WeightEnumerator(V[end][1].polynomial, "Pauli")
    # println(W)

    if cleanV
        for i in 2:length(V)
            for v in V[i]
                v.polynomial = [[0, 0, 0, 0]]
            end
        end
    end
    return W
end

# function Pauliweightenumerator(Q::AbstractQuantumCode, Pauli::Char=' ',
#     keeptrellis::Bool=true, cleanV::Bool=true, sect::Bool=false)
function Pauliweightenumerator(Q::AbstractQuantumCode, Pauli::Char=' ',
    sect::Bool=false)

    !ismissing(Q.Pauliweightenum) && return Q.Pauliweightenum
    T = syndrometrellis(Q, "weight", Pauli, sect)
    # if keeptrellis
    #     return Pauliweightenumerator(T, cleanV), T
    # end
    # return PauliweightenumeratorQ(T, cleanV)
    return Pauliweightenumerator(T, false)
end

function PWEtoHWE(PWE::WeightEnumerator)
    poly = deepcopy(PWE)
    for term in poly
        tot = abs(term[2] + term[3] + term[4])
        term[2] = tot
        term[3] = length(S) - tot
        term[4] = 0
    end
    return WeightEnumerator(_reducepolyQ(poly), "Hamming")
end

function PWEtoXWE(PWE::WeightEnumerator)
    poly = deepcopy(PWE)
    for term in poly
        tot = abs(term[2] + term[3])
        term[2] = tot
        term[3] = length(S) - tot
        term[4] = 0
    end
    return WeightEnumerator(_reducepoly(poly), "X")
end

function PWEtoZWE(PWE::WeightEnumerator)
    poly = deepcopy(PWE)
    for term in poly
        tot = abs(term[3] + term[4])
        term[2] = length(S) - tot
        term[3] = 0
        term[4] = tot
    end
    W = WeightEnumerator(_reducepoly(poly), "Z")
    if ismissing(S.Zwtenum)
        S.Zwtenum = W
    end
    return W
end

function HammingweightenumeratorQ(T::Trellis, cleanV::Bool=true)
    V = vertices(T)
    E = edges(T)
    for i in 2:length(V)
        for (j, v) in enumerate(V[i])
            outer = Vector{Vector{Int64}}()
            for e in E[i - 1][j]
                inner = deepcopy(V[i - 1][e.outvertex].polynomial)
                for k in e.label
                    if iszero(k)
                        for term in inner
                            term[3] += 1
                        end
                    else
                        for term in inner
                            term[2] += 1
                        end
                    end
                end
                append!(outer, inner)
            end
            # println(outer)
            v.polynomial = _reducepolyQ(outer)
            # println(v.polynomial)
        end
    end
    W = WeightEnumerator(sort!(V[end][1].polynomial, lt=islessPolyQ), "Hamming")

    if cleanV
        for i in 2:length(V)
            for v in V[i]
                v.polynomial = [[0, 0, 0, 0]]
            end
        end
    end
    return W
end

# function HammingweightenumeratorQ(Q::AbstractQuantumCode, Paui::Char=' ',
#     keeptrellis::Bool=true, cleanV::Bool=true, sect::Bool=false)
function Hammingweightenumerator(Q::AbstractQuantumCode, Paui::Char=' ',
    sect::Bool=false)

    !ismissing(Q.Pauliweightenum) && return PWEtoHWE(Q.Pauliweightenum)
    T = syndrometrellis(Q, "weight", Pauli, sect)
    # if keeptrellis
    #     return HammingweightenumeratorQ(T, cleanV), T
    # end
    # return HammingweightenumeratorQ(T, cleanV)
    return HammingweightenumeratorQ(T, false)
end

# function weightenumerator(Q::AbstractQuantumCode, alg::String="trellis", sect::Bool=false,
#     Pauli::Char=' ', keeptrellis::Bool=true, cleanV::Bool=true, verbose::Bool=false)

# we can organize the S⟂ \ S any way we want as long as we have these elements on the
function weightenumerator(Q::AbstractQuantumCode, alg::String="trellis", sect::Bool=false,
    Pauli::Char=' ', verbose::Bool=true)

    alg ∈ ["trellis"] || error("Algorithm `$alg` is not implemented in weightdistribution.")

    # S⟂ - this is not a stabilizer code
    # if ismissing(Q.dualweightenum)
    #     Q.dualweightenum = weightenumerator(LinearCode(normalizermatrix(Q)))
    # end
    verbose && println("S⟂: ", Q.dualweightenum)

    # S⟂ \ S - this is not a stabilizer code
    # k is always small, so just brute force it
    if ismissing(Q.logsweightenum)
        logspace = LinearCode(logicalspace(Q))
        Q.logsweightenum = _weightenumerator(logspace)
        Q.d = Q.logsweightenum.polynomial[2][2]
    end
    verbose && println("S⟂ / S: ", Q.logsweightenum)

    # S
    if ismissing(Q.Pauliweightenum)
        # this is not valid for qudit codes
        Q.Pauliweightenum = Pauliweightenumerator(Q, Pauli, sect)
    end
    verbose && println("S: ", Q.Pauliweightenum)

    return [Q.dualweightenum, Q.logsweightenum, Q.Pauliweightenum]
end

function weightdistribution(Q::AbstractQuantumCode, alg::String="trellis", sect::Bool=false,
    Pauli::Char=' ', keeptrellis::Bool=true, cleanV::Bool=true, verbose::Bool=true)

    weightenumerator(Q, alg, sect, Pauli, false)
    # S⟂
    temp1 = zeros(Int64, 1, length(Q) + 1)
    for term in Q.dualweightenum
        temp1[term[2] + 1] = term[1]
    end
    verbose && println("S⟂: ", temp1)

    # S⟂ \  S
    temp2 = zeros(Int64, 1, length(Q) + 1)
    for term in Q.logsweightenum
        temp2[term[2] + 1] = term[1]
    end
    verbose && println("S⟂ / S: ", temp2)

    # S
    temp3 = zeros(Int64, 1, length(Q) + 1)
    for term in Q.Pauliweightenum
        temp3[term[2] + 1] = term[1]
    end
    verbose && println("S: ", temp3)

    return [temp1, temp2, temp3]
end

# function weightdistribution(Q::AbstractQuantumCode, alg::String="trellis", sect::Bool=false)
#     alg ∈ ["trellis"] || error("Algorithm `$alg` is not implemented in weightdistribution.")
#     T = syndrometrellis(Q, "weight", ' ', sect)
#     W = PauliweightenumeratorQ(T, true)
#     H = PWEtoHWE(W)
#     temp = zeros(Int64, 1, length(Q))
#     for term in H
#         temp[term[2]] = term[1]
#     end
#
#     if typeof(Q) <: CSSCode
#         X = PWEtoXWE(W)
#         Z = PWEtoZWE(W)
#         dist = Vector{Vector{Int64}}()
#         push!(dist, temp)
#
#         temp = zeros(Int64, 1, length(Q))
#         for term in X
#             temp[term[2] + 1] = term[1]
#         end
#         push!(dist, temp)
#
#         temp = zeros(Int64, 1, length(Q))
#         for term in Z
#             temp[term[4] + 1] = term[1]
#         end
#         push!(dist, temp)
#         return dist
#     end
#     return temp
# end

#############################
     # Minimum Distance
#############################

"""
    minimumdistance(Q::AbstractQuantumCode, alg::String="trellis", sect::Bool=false)

Return the minimum distance of the stabilizer code if known, otherwise computes it.

"""
function minimumdistance(Q::AbstractQuantumCode)
    !ismissing(Q.d) && return Q.d
    if ismissing(Q.logsweightenum)
        logspace = logicalspace(Q)
        Q.logsweightenum = _weightenumerator(logspace)
    end
    Q.d = Q.logsweightenum.polynomial[2][2]
    return Q.d
end

# function minimumdistance(Q::AbstractCSSCode, alg::String="trellis", sect::Bool=false)
#
#
#
#
#
#
#     !ismissing(Q.d) && return Q.d
#     if ismissing(Q.logsweightenum)
#         Q.logsweightenum = weightenumerator(LinearCode(logicalspace(Q)))
#     end
#     # first element is always the identity y^n, want power of first x term
#     Q.d = Q.logsweightenum.polynomial[2][2]
#     return Q.d
# end

max_stars_repo_name: esabo/CodingTheory

max_issues_repo_name: esabo/CodingTheory

max_forks_repo_name: esabo/CodingTheory