text:
"""
    _highlight(x)

Apply style to `x` based on its type.
"""
function _highlight end

_highlight(x::Union{Symbol,AbstractString}) = "[$(theme.symbol)]$x[/$(theme.symbol)]"
_highlight(x::Number) = "[$(theme.number)]$x[/$(theme.number)]"
_highlight(x::DataType) = "[$(theme.type) dim]::$(x)[/$(theme.type) dim]"
_highlight(x::UnitRange) = _highlight(string(x))
function _highlight(x::Union{AbstractArray,AbstractVector,AbstractMatrix})
    return "[$(theme.number)]$x[/$(theme.number)]"
end
_highlight(x::Function) = "[$(theme.function)]$x[/$(theme.function)]"

function _highlight(x)
    return _highlight(string(x))
end

"""
    _highlight_with_type(x)

Apply style to x and and mark its type.
"""
_highlight_with_type(x) = "$(_highlight(x))$(_highlight(typeof(x)))"

"""
    _highlight_numbers(x::AbstractString)

Add style to each number in a string
"""
function _highlight_numbers(x::AbstractString)
    for match in collect(eachmatch(r"[0-9]", x))
        x = replace(x, match.match => "[blue]$(match.match)[/blue]")
    end
    return x
end

"""
    style_error(io::IO, er)

Create a style error panel.

Creates a `Panel` with an error message and optional hints.
"""
function style_error(io::IO, er)
    if haskey(ErrorsExplanations, typeof(er))
        info_msg = ErrorsExplanations[typeof(er)]
    else
        info_msg = nothing
    end

    WIDTH = _width()
    main_message, message = error_message(io, er)
    # @info "Got styled error" main_message info_msg isnothing(info_msg) typeof(message)

    # create panel and text
    if !isnothing(info_msg) &&  Measure(message).w > 0
        panel = Panel(
            main_message / message,
            hLine(WIDTH - 4; style = "red dim"),
            RenderableText(
                "[bold yellow italic underline]hint:[/bold yellow italic underline] [bright_red]$(typeof(er))[/bright_red] " *
                info_msg;
                width = WIDTH - 4,
            );
            title = "ERROR: [bold indian_red]$(typeof(er))[/bold indian_red]",
            title_style = "red",
            style = "dim red",
            width = WIDTH,
            title_justify = :left,
        )
    else
        panel = Panel(
            main_message;
            title = "ERROR: [bold indian_red]$(typeof(er))[/bold indian_red]",
            title_style = "red",
            style = "dim red",
            width = WIDTH,
            title_justify = :left,
        )
    end
    text = RenderableText(
        "\n[bold indian_red]$(typeof(er)):[/bold indian_red] $main_message"
    )

    return panel, text
end

"""
    backtrace_subpanel(line::String, WIDTH::Int, title::String)

Create a subpanel for stacktrace, showing source code.
"""
function backtrace_subpanel(line::String, WIDTH::Int, title::String)
    # @info "getting subpanel"
    # get path to file and line number
    file = split(split_lines(line)[2], " [bold dim]")[1][13:end]
    file, lineno = split(file, ":")
    lineno = parse(Int, lineno)
    # @info "got file and lines" file lineno

    # read and highlight text
    code = ""
    try
        code = load_code_and_highlight(file, lineno; δ = 2)
        if length(code) > 0
            code = TextBox(code; width = WIDTH - 26)
            code = Spacer(8, code.measure.h) * code
        end
    catch SystemError  # file not found
        # @warn "Failed to get code"
        code = ""
    end

    # @info "line&code" line code
    panel = Panel(
        "\n",
        chomp(line),
        code;
        title = title,
        width = :fit,
        style = "dim blue",
        title_style = "bold bright_yellow",
    )
    # @info "Created backtrace subpanel"
    return panel
end

"""
    style_backtrace(io::IO, t::Vector)

Create a Panel with styled error backtrace information.
"""
function style_backtrace(io::IO, t::Vector)
    # @info "styling backtrace"
    WIDTH = _width()

    # create text
    stack_lines::Vector{String} = []
    for (n, frame) in enumerate(t)
        if typeof(frame) ∉ (Ptr, InterpreterIP)
            func_line = "[light_yellow3]($n)[/light_yellow3] [sky_blue3]$(frame.func)[/sky_blue3]"
            file_line = "       [dim]$(frame.file):$(frame.line) [bold dim](line: $(frame.line))[/bold dim][/dim]"
            push!(stack_lines, func_line * "\n" * file_line)
        end
    end

    # trim excedingly long stack traces
    if length(stack_lines) > 10
        stack_lines = vcat(
            stack_lines[1:5],
            "\n[dim bright_blue]        ... skipped $(length(stack_lines)-11) levels in stack trace ...\n",
            stack_lines[(end - 5):end],
        )
    end

    # create layout
    # @info "creating stack panels" length(stack_lines)
    if length(stack_lines) > 0
        if length(stack_lines) > 1
            error_line = backtrace_subpanel(stack_lines[1], WIDTH - 8, "error in")
        else
            error_line = ""
        end

        # @info "error line ready" error_line length(stack_lines)
        if length(stack_lines) > 3
            above = TextBox(join(stack_lines[2:(end - 1)], "\n"); width = WIDTH - 8)
        elseif length(stack_lines) > 2
            above = TextBox(stack_lines[2]; width = WIDTH - 8)
        else
            above = ""
        end
        # @info "above ready"

        offending = backtrace_subpanel(stack_lines[end], WIDTH - 8, "caused by")
        # @info "offending ready" offending
        stack = error_line / above / offending
        # @info "stacked" stack
    else
        stack = "[dim]No stack trace[/dim]"
    end

    # create output
    # @info "creating panel"
    try
        panel = Panel(
            stack;
            title = "StackTrace",
            style = "yellow dim",
            title_style = "yellow",
            title_justify = :left,
            width = WIDTH,
        )
        # @info "panel ready"
        return panel
    catch err
        @warn "failed to crate panel" err
        println(stack)
        return stack
    end
end

"""
    style_stacktrace_simple(stack::Vector)

Simply style a stacktrace. Just adding a miniumum of color.
"""
function style_stacktrace_simple(stack::Vector)
    lines::Vector{String} = []
    for (n, frame) in enumerate(stack)
        if typeof(frame) ∉ (Ptr.Ptr{Nothing}, InterpreterIP)
            push!(
                lines,
                apply_style(
                    "[bright_blue dim]($n)[/bright_blue dim] [yellow]$(frame.func)[/yellow] - [dim]$(line.file):$(line.line)[/dim]",
                ),
            )
        end
    end

    return join(lines, "\n")
end


text:
#EasyPlotQwt: Render EasyPlot-plots with Qwt guiqwt's interface
#-------------------------------------------------------------------------------
module EasyPlotQwt

import EasyPlot #Import only - avoid collisions
using MDDatasets
using Colors

import EasyPlot: render

#Type used to dispatch on a symbol & minimize namespace pollution:
#-------------------------------------------------------------------------------
struct DS{Symbol}; end; #Dispatchable symbol
DS(v::Symbol) = DS{v}()

include("pybase.jl")
include("base.jl")
include("display.jl")

#==Initialization
===============================================================================#
function __init__()
	EasyPlot.registerdefaults(:Qwt,
		maindisplay = PlotDisplay(guimode=true),
		renderdisplay = EasyPlot.NullDisplay() #No support at the moment.
	)
	return
end

end

#Last line

ext: jl


text:
@testset "integral, constant" begin
    N = (10,20)
    v₀ = rand()
    val = fill(v₀, N)

    p₀ = [rand(), rand()]
    r = 5.0
    xₙ, xₚ = p₀[1]-r, p₀[1]+r
    yₙ, yₚ = p₀[2]-r, p₀[2]+r

    xbound = range(xₙ, xₚ, length=N[1]+1)
    ybound = range(yₙ, yₚ, length=N[2]+1)
    vxlbound = (xbound, ybound)

    vc = VoxelwiseConstant(val, vxlbound)

    @test integral(vc) ≈ v₀ * (xₚ-xₙ) * (yₚ-yₙ)
    @test (s = integrate(vc, p->norm(p .- p₀)≤r); abs(s.val - v₀ * π * r^2) ≤ s.err)
end  # @testset "integral, constant"

@testset "integral, non-constant" begin
    N = (10,20)
    val = rand(N...)

    p₀ = [rand(), rand()]
    r = 5.0
    xₙ, xₚ = p₀[1]-r, p₀[1]+r
    yₙ, yₚ = p₀[2]-r, p₀[2]+r

    xbound = range(xₙ, xₚ, length=N[1]+1)
    ybound = range(yₙ, yₚ, length=N[2]+1)
    vxlbound = (xbound, ybound)

    ∆x = (xₚ-xₙ) / N[1]
    ∆y = (yₚ-yₙ) / N[2]

    vc = VoxelwiseConstant(val, vxlbound)

    @test integral(vc) ≈ sum(val) * ∆x * ∆y
end  # @testset "integral, non-constant"


text:
ith = 1
target = 100010
n = 3
while ith < target
    flag = true
    for i = 2:floor(Int64, sqrt(n))
        if n % i == 0
            flag = false
        end
    end

    if flag
        global ith+=1
    end

    if ith == target
        println(n)
        break
    end
    global n += 1
end


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

using Markdown
using InteractiveUtils

# ╔═╡ d941d6c2-55bf-11eb-0002-35c7474e4050
using NiLang, Test

# ╔═╡ 2061b434-0ad1-46eb-a0c7-1a5f432bfa62
begin
	twocol(left, right; llabel="forward", rlabel="backward") = HTML("

	$llabel
	$rlabel

$(html(left))

$(html(right))

")
	example(str) = HTML("""
$str
""")
	title1(str) = HTML("""

'_'))>$str

""")
	title2(str) = HTML("""
'_'))>$str
""")
	titleref(str) = HTML("""'_'))">$str""")
	using PlutoUI: TableOfContents
	using Pkg
	pkgversion(m::Module) = Pkg.TOML.parsefile(NiLang.project_relative_path("Project.toml"))["version"]
end;

# ╔═╡ 8c2c4fa6-172f-4dde-a279-5d0aecfdbe46
module M
using NiLang

# define two functions
function new_forward(x)
	if x > 0
		return x * 2
	elseif x < 0
		return x / 2
	end
end

function new_backward(x)
	if x > 0
		return x / 2
	elseif x < 0
		return x * 2
	end
end

# declare them as reversible to each other
@dual new_forward new_backward

# The following is need only when your function is differentiable
using NiLang.AD: GVar
function new_backward(x::GVar)
	if x.x > 0
		GVar(new_backward(x.x), x.g * 2)
	elseif x.x < 0
		GVar(new_backward(x.x), x.g / 2)
	end
end

function new_forward(x::GVar)
	if x.x > 0
		GVar(new_forward(x.x), x.g / 2)
	elseif x.x < 0
		GVar(new_forward(x.x), x.g * 2)
	end
end
end

# ╔═╡ 0e1ba158-a6bc-401c-9ba7-ed78020ad068
using Base.Threads

# ╔═╡ 3199a048-7b39-40f8-8183-6a54cccd91b6
using BenchmarkTools

# ╔═╡ a4e76427-f051-4b29-915a-fdfce3a299bb
html"""

  <svg class="octicon octicon-mark-github v-align-middle" height="32" viewBox="0 0 16 16" version="1.1" width="32" aria-hidden="true"><path fill-rule="evenodd" d="M8 0C3.58 0 0 3.58 0 8c0 3.54 2.29 6.53 5.47 7.59.4.07.55-.17.55-.38 0-.19-.01-.82-.01-1.49-2.01.37-2.53-.49-2.69-.94-.09-.23-.48-.94-.82-1.13-.28-.15-.68-.52-.01-.53.63-.01 1.08.58 1.23.82.72 1.21 1.87.87 2.33.66.07-.52.28-.87.51-1.07-1.78-.2-3.64-.89-3.64-3.95 0-.87.31-1.59.82-2.15-.08-.2-.36-1.02.08-2.12 0 0 .67-.21 2.2.82.64-.18 1.32-.27 2-.27.68 0 1.36.09 2 .27 1.53-1.04 2.2-.82 2.2-.82.44 1.1.16 1.92.08 2.12.51.56.82 1.27.82 2.15 0 3.07-1.87 3.75-3.65 3.95.29.25.54.73.54 1.48 0 1.07-.01 1.93-.01 2.2 0 .21.15.46.55.38A8.013 8.013 0 0016 8c0-4.42-3.58-8-8-8z"></path></svg>

Download this notebook

"""

# ╔═╡ c2c7b4d4-f8c9-4ebf-8da2-0103f03136e7
md"# NiLang's (v$(pkgversion(NiLang))) Documentation

NiLang is a embeded domain specific language (eDSL) in Julia, so one need to install [Julia](https://julialang.org/) first. Before reading this documentation, you need to know basic Julia grammar, and how to install and use packages.
Also, it might be good to read the [README](https://github.com/GiggleLiu/NiLang.jl) first. In this tutorial, we focus on

* NiLang grammar and design patterns.
* Automatic differentiation based on reversible programming and time-space tradeoff
"

# ╔═╡ 605872cf-f3fd-462e-a2b1-7d1c5ae45efd
title1("Getting started")

# ╔═╡ fb3dee44-5fa9-4773-8b7f-a83c44358545
md"
After installing NiLang in a Julia REPL by typing `]add NiLang`, one can use NiLang and use the macro `@i` to define a reversible function ``f_1: (x, y) → (x+5y, y)``."

# ╔═╡ 70088425-6779-4a2d-ba6d-b0a34c8e93a6
@i function f1(x, y; constant)
	x += y * constant
end

# ╔═╡ f0e94247-f615-472b-8218-3fa287b38aa1
md"The function arguments contains keyword and non-keyword arguments, where non-keyword arguments can only be used as constants. There is no `return` statement because this function returns input non-keyword variables automatically. `return` is both not needed and not allowed."

# ╔═╡ af738f89-3214-429c-9c7d-18a6ea0d9401
f1(1.0, 2.0; constant=5.0)

# ╔═╡ 4d98d529-1e05-49be-9209-f0d9fcc206f7
md"The inverse call (or backward execution) is"

# ╔═╡ 48d7ebc1-5def-4a57-9ec1-3fc370a4543f
 (~f1)(11.0, 2.0; constant = 5.0)

# ╔═╡ e839547b-f47e-4a4e-b4d9-2c22921d80e4
md"Now we see the input arguments are restored"

# ╔═╡ 60575978-081a-4bca-a3ed-2b51cd6abc92
md"One can also differentiating the function. This is a two-in two-out function (the keyword arguments are not included), we need to specify which variable is the loss when asking for the gradients."

# ╔═╡ f98305cb-4ba2-404a-a5c3-65510e059504
NiLang.AD.gradient(f1, (1.0, 2.0); iloss=1, constant=5.0)

# ╔═╡ e8cd6667-597f-458b-8465-1822e09a7891
md"Here, we specify the first variable as the one that stores the loss. We get
```math
\begin{cases}
\frac{\partial x+5y}{\partial x}=1\\
\frac{\partial x+5y}{\partial y}=5\\
\end{cases}
```
"

# ╔═╡ 20145d75-004a-4c2f-b7ff-c400ca846d42
let
	content = md"""The above macro generates two functions, one is `f` and another is `~f` (or `Inv(f)`). The `x += y * constant` is translated to function `(x, y, constant) = PlusEq(*)(x, y, constant)`, where the function `PlusEq(*)` is bijective.

```julia
julia> using MacroTools

julia> MacroTools.prettify(@macroexpand @i function f(x, y; constant)
           x += y * constant
       end)
quote
    $(Expr(:meta, :doc))
    function $(Expr(:where, :(f(x, y; constant))))
        hare = wrap_tuple(((PlusEq)(*))(x, y, constant))
        x = hare[1]
        y = hare[2]
        constant = hare[3]
        (x, y)
    end
    if (NiLangCore)._typeof(f) != _typeof(~f)
        function $(Expr(:where, :((newt::_typeof(~f))(x, y; constant))))
            boar = wrap_tuple(((MinusEq)(*))(x, y, constant))
            x = boar[1]
            y = boar[2]
            constant = boar[3]
            (x, y)
        end
    end
end
```
"""
	HTML("If you know macros$(html(content))")
end

# ╔═╡ b6dcd18c-606f-4340-b2ec-163e8bad03f5
title1("Variables")

# ╔═╡ 2581aa33-1dc5-40b1-aa9f-6a11cc750c93
md"In NiLang, a variable is mutable. After running an instruction that defines the mapping ``\mathbb{R}^n\rightarrow\mathbb{R}^n`` (notice inputs and outputs have the same shapes), the outputs are directly assigned back. In other words, every instruction changes the variable inplace."

# ╔═╡ a1a29f34-f8a9-4e9f-9afe-7d0096771440
title2("Allocate and deallocate a variable")

# ╔═╡ 90bd6ad4-3dd8-4e7c-b445-aed1d248a2ec
md"""
One can allocate a new variable `x` like `x ← constant` and deallocate a variable with known value with `x → known value`. They are reversible to each other with the following relation.

$(
twocol(md"
```julia
x ← constant
```", md"
```julia
x → constant
```
")
)

For example
"""

# ╔═╡ c0259e48-1973-486c-a828-1fcd3e4331c6
@i function f2a()
	tmp ← 1
	# some code that uses `tmp` for computing and restores it to `1`
	tmp → 1
end

# ╔═╡ 7bcc1bc6-5b2c-4aae-b71d-52482ae130e9
md"In reversible computing, it is forbidden to nullify (zero clear or assignment) the content in a variable."

# ╔═╡ fd5d65e7-3d83-45dd-87d7-ff1e1e769430
#`x` is in the input argument list, the content in `x` should not be cleared and re-allocated to another variable.
@test_throws LoadError macroexpand(NiLang, :(@i function f2a(x)
	x ← 1
	x → 1
end))

# ╔═╡ 4601df35-679f-465d-9191-c18748b2fd83
title2("Avoid shared read-write")

# ╔═╡ 4fc72b9d-19a2-40f1-a4a8-5e97d3d5e529
md"Shared read-write is not allowed"

# ╔═╡ a0fde16f-8454-4f5c-a29c-a9e415c0c311
# `y -= y` effectively clears the content in `y`, this is why shared read-write is so dangerous.
@test_throws LoadError macroexpand(NiLang, :(@i function f2b(y)
	y -= y
end))

# ╔═╡ 633ff8f3-8d93-4f73-bec2-c42070e6ece9
md"NiLang is more restrictive, it forbids shared read too. This is in purpose, shared read will become shared write to the gradient field in the back-propagation of gradients - the main goal of NiLang."

# ╔═╡ 57f2d890-b5a0-47b7-9e3d-af4d03b10605
# `shared read is also forbidden`
@test_throws LoadError macroexpand(NiLang, :(@i function f2b(x, y)
	x -= y * y  # should be written as `x -= y^2`
end))

# ╔═╡ 34063cd0-171e-46ce-80dd-52a341fa50a1
md"The correct way of avoiding shared read is renaming one of the variable."

# ╔═╡ 5d5d01db-8ff9-434c-8771-1fec6393e1fb
@i function f2c(x, y)
	tmp ← zero(y)
	tmp += y
	x -= y * tmp
	tmp -= y
	tmp → zero(y)
end

# ╔═╡ 10d85a50-f2f9-403e-8f6c-baef61cf702a
f2c(0.0, 3.0)

# ╔═╡ b52648bf-a28a-48af-8912-31729d943ce0
md"Shared read-write issue is more tricky when one uses NiLang to write kernel function in a parallel program (multi-threading, MPI and CUDA). See $(titleref(\"Multi-threading and CUDA\")) for details."

# ╔═╡ f3503e6d-ded7-4797-9c2f-cb5a7d429949
md"#### Manipulating allocations on stack and dictionary"

# ╔═╡ 7ba86798-daba-49ef-a259-da4db353fdd8
md"""
The forward pass push a variable to the stack and zero-empty the variable, the backward pass. The backward pass checks the variable is empty and pop an variable from the stack and assign it to the variable.

$(
twocol(md"
```julia
PUSH!(stack, variable)
```", md"
```julia
POP!(stack, variable)
```
")
)

Or similarly, there is a `COPYPUSH!` instruction that copies the content in the variable to the stack without zero clearing the variable. In the backward pass, it asserts the target variable is the same as the popped variable.

$(
twocol(md"
```julia
COPYPUSH!(stack, variable)
```", md"
```julia
COPYPOP!(stack, variable)
```
")
)
"""

# ╔═╡ 6e88ede9-879e-42e9-9a7e-2aa3d4837e5d
@i function f2e(stack1::AbstractVector, stack2::AbstractVector)
	var1 ← 3.14
	var2 ← 2.73

	# copy a new variable to a stack
	COPYPUSH!(stack1, var1)
	# add a new variable to a stack and empty the original variable
	PUSH!(stack2, var2)

	# deallocate the original variables
	var2 → 0.0
	var1 → 3.14
end

# ╔═╡ 94cc42f7-6beb-453e-90cb-952d6f497a8e
f2e([], [])

# ╔═╡ 93936612-1447-4114-b864-aba43adef4bd
md"""
The forward pass of dictionary allocation adds an entry to the dictionary (cause key error if the key already exists), the backward pass. The backward pass checks the variable in the dictionary is consistent with the asserted variable, and delete the key.

$(
twocol(md"
```julia
dict[key] ← variable
```", md"
```julia
dict[key] → variable
```
")
)
"""

# ╔═╡ b2add70c-c5d6-4e0f-a153-43e21a197181
@i function f2f(dict::Dict)
	var1 ← 3.14

	# copy a new variable to a dict
	dict["data"] ← var1

	# deallocate the original variables
	var1 → dict["data"]
end

# ╔═╡ 5c03d5a5-99f0-4efd-9a32-ce6d7c2b266c
f2f(Dict())

# ╔═╡ f45db10f-a836-40f3-9d8d-054ea6540e87
title2("Protect special variables")

# ╔═╡ 1903563e-ccc2-44d9-9dbe-e5dede275b3c
md"""
By default, NiLang assigns the value back to the variable after a call. Sometime it can cause issues for special variables like function, type parameters and constants (including the results generated from a function call). One can use `@const` (assert a variable is a constant) or `@skip!` (skip assigning back) to avoid such complications.
"""

# ╔═╡ 7ec577d4-1e19-484e-ac1f-5da3f47d1ec4
md"""
###### A function variable
"""

# ╔═╡ 583f2585-15a3-47c6-a70e-e2f002754028
md"When using a function (e.g. `exp` as shown bellow) as an variable, one should be careful about the scope issue."

# ╔═╡ 60e6ff80-3593-4ae4-a273-914847f692db
@i function f2g(y, f, x)
	y += f(x)
end

# ╔═╡ 9e5cfd68-b58d-4d83-aae2-447e5f805c97
@i function f2h(y, x)
	f2g(y, exp, x)
end

# ╔═╡ dc85a942-cf52-4405-ad03-32a768e1b6e7
@test_throws UndefVarError f2h(0.0, 3.0)

# ╔═╡ 5cdd346b-10a5-485c-ba78-4c0b3cb0e02f
md"We see an error, but why calling `f2g` causes an error? If one check the generated code with `macroexpand`, one will see the `exp` is assigned in the local scope. The compiler takes it as a local variable and compaints that `exp` is not defined."

# ╔═╡ af9287b7-6131-46f6-beb8-6885e55e1975
macroexpand(NiLang, :(@i function f2h(y, x)
	f2g(y, exp, x)
end)) |> NiLangCore.rmlines

# ╔═╡ e20eeabf-1c80-431e-8cfc-4d1b79c52b5a
@i function f2i(y, x)
	f2g(y, (@const exp), x)
end; @test f2i(0.0, 3.0)[1] == exp(3)

# ╔═╡ e45d53b9-6f02-4293-ab4f-85bd751993ad
md"""
###### Type parameters and functions
"""

# ╔═╡ 90d30eea-53de-48a0-9700-ff35681fdf38
md"Type parameter can not be assigned back too."

# ╔═╡ 390f58a5-6f5f-4d3a-bb16-ba04e43a07e7
@test_throws ErrorException Core.eval(NiLang, :(@i function f2j(t::Type{T}, x) where T
	x += one(T)
end))

# ╔═╡ 2b57443e-a516-434b-be86-80616a98e2f5
@i function f2k(t::Type{T}, x) where T
	x += one(@const T)
end; @test f2k(Float64, 0.0)[2] == 1.0

# ╔═╡ fc2e27f9-b7ba-44cc-a953-6745548ad733
md"A function call that returning a constant should also be decorated with the `@const` or `@skip` macro."

# ╔═╡ fc744931-360b-4478-9f77-c50f048de243
@test_throws LoadError macroexpand(NiLang, :(@i function f2l(y, x::Matrix) where T
	y += size(x, 1) * size(x, 2)
end))

# ╔═╡ 9a152b36-f377-44da-9700-ca9e05e365ff
@i function f2m(y, x::Matrix) where T
	y += (@const size(x, 1)) * (@const size(x, 2))
end; @test f2m(0, randn(3,4))[1] == 12

# ╔═╡ 0863bd06-cc70-4dde-b3b2-0a466805a356
md"""$(title1("Integers, floating-point numbers, fixed-point numbers and logarithmic numbers"))

A fixed-point zero with 43 fraction bits can be declared as `x ← Fixed43(0.0)` or `x ← zero(Fixed43)`, a logarithmic one can be declared as `x ← ULogarithmic(1.0)`, `x ← one(ULogarithmic{Float64})` or `x ← ULogarithmic(Fixed43(1.0))`. A complex number zero can be defined as `x ← Complex(0.0, 0.0)`.

| Number Type | +=  | *= | ⊻= | Source |
| ---- | ---- | ---- | ---- | ---- |
| boolean | - |  - | ✓ | JuliaLang |
| integer | ✓ |  × | ✓ | JuliaLang |
| floating-point number | ✓ (rounding error) |  × | - | JuliaLang |
| fixed-point number | ✓ |  × | - | [FixedPointNumbers.jl](https://github.com/JuliaMath/FixedPointNumbers.jl)
| logarithmic number | × |  ✓ | - |  [LogarithmicNumbers.jl](https://github.com/cjdoris/LogarithmicNumbers.jl)

* `✓`: an operation has its reverse when operating on a number type.
* `×`: an operation does not have a reverse when operating on a number type.
* `-`: an operation does not apply on a number type.

The `+=` operation is not regoriously reversible on floating point numbers, but we ignore the rounding errors in NiLang and use the reversibility check to detect the potential too-large rounding errors. Whether logarithmic number has rounding errors depends on its content type. If it uses floating point numbers as storage, then yes, otherwise if it uses fixed point number as the content type, then no.

One can use the `y ⊙= convert(x)` statement to convert `x` to the target type `typeof(y)` and accumulate it to `y`. Here `⊙=` can be one of `+=`, `*=` and `⊻=` that has its reverse on type `typeof(y)`.
"""

# ╔═╡ a0bae195-04e1-4642-9e14-fe4691e0906b
md"Fixed point numbers and Floating point numbers are reversible under the `+=` operation"

# ╔═╡ 20d6e8a0-2cf5-48ad-9549-60506b42b970
@i function f3c(y1::T, x::T) where T<:Union{Fixed43, AbstractFloat}
	y1 += x
end; @test f3c(Fixed43(0.5), Fixed43(0.6))[1] === Fixed43(1.1) && f3c(0.5, 0.6)[1] === 1.1

# ╔═╡ 7dee5748-ed73-4e13-aa80-7a50efbc8449
@i function f3d(y1::T, x::T) where T<:Union{Fixed43, AbstractFloat}
	y1 *= x
end; @test_throws MethodError f3d(Fixed43(1.0), Fixed43(2.0))

# ╔═╡ 4c719e9b-641e-404e-9ab7-59e89135f3ba
md"Logarithmic numbers are reversible under the `*=` operation"

# ╔═╡ 77947e00-42c3-4c9e-b62a-b4b29489db43
@i function f3a(y1::ULogarithmic{T}, x::ULogarithmic{T}) where T
	y1 += x
end; @test_throws MethodError f3a(ULogarithmic(1.0), ULogarithmic(3.0))

# ╔═╡ 2614127d-34fb-4c3d-b678-42693f3c9341
@i function f3b(y1::ULogarithmic{T}, x::ULogarithmic{T}) where T
	y1 *= x
end; @test f3b(ULogarithmic(1.0), ULogarithmic(3.0))[1] === ULogarithmic(3.0)

# ╔═╡ 8f169235-3bd1-4cc4-a083-79736d306ad5
example("computing x ^ 3 with logarithmic numbers")

# ╔═╡ dfc9d305-5bce-4555-bfa3-d8d61fe4ca09
@i function power3(y::ULogarithmic{T}, x::T) where T
	for i=1:3
		y *= convert(x)
	end
end; @test power3(ULogarithmic(1.0), 3.0)[1] ≈ 27

# ╔═╡ c4cd9f88-9cd6-4364-b016-78f90aba6a66
md"""## Mathematical Operations and Elementary Functions
A mathematical operation can be defined in the form `⊙=`, where `⊙` can be `+`, `-`, `*`, `/` or `⊻`. See section $(titleref("Integers, floating-point numbers, fixed-point numbers and logarithmic numbers")) for their supported number types.

When we say an elementary function is supported, we mean its diffrule is defined. A list of supported elementary differentiable functions are

| instruction | output   |
| ----------- | ---------- |
| ``{\rm FLIP}(y)`` | ``\sim y`` |
| ``{\rm NEG}(y)`` | ``-y`` |
| ``{\rm INC}(y)`` | ``y+1`` |
| ``{\rm DEC}(y)`` | ``y-1`` |
| ``{\rm INV}(y)`` | ``y ^ {-1}`` |
| ``{\rm HADAMARD}(x, y)`` | ``\frac{1}{\sqrt{2}}(x+y), \frac{1}{\sqrt{2}}(x-y)``
| ``{\rm SWAP}(a, b)`` | ``b, a`` |
| ``{\rm ROT}(a, b, \theta)`` | ``a \cos\theta - b\sin\theta, b \cos\theta + a\sin\theta, \theta`` |
| ``{\rm IROT}(a, b, \theta)`` | ``a \cos\theta + b\sin\theta, b \cos\theta - a\sin\theta, \theta`` |
| ``y \mathrel{\{+,-\}}= f_{+-}(args...)`` | ``y\{+, -\}f_{+-}(args...), args...`` |
| ``y \mathrel{\{*, /\}}= f_{*/}(args...)`` | ``y\{*, /\}f_{*/}(args...), args...`` |

Functions ``f_{+-} ∈ \rm \{identity, +, -, *, /, ^\wedge, abs, abs2,\rm sqrt, exp, log, sin, sinh, asin, cos, cosh,`` ``acos, tan, tanh, atan, sincos, convert\}`` and ``f_{*/}∈\rm \{identity, +, -, *, /, ^\wedge, convert\}.``
Functions `FLIP`, `NEG`, `INV`, `HADAMARD`, `SWAP` and `y ⊻= f_{⊻}(args...)` are self-reversible (or reflexive). {`ROT`, `IROT`} and {`INC`, `DEC`}, {`y += f_{+-}(args...)`, `y -= f_{+-}(args...)`} and {`y *= f_{*/}(args...)`, `y /= f_{*/}(args...)`} are pair-wise reversible.

For Jacobians and Hessians defined on these instructions, please check this [blog post](https://giggleliu.github.io/2020/01/18/jacobians.html).
"""
# They are translated to `y += f(args...)` is translated to `PlusEq(f)(y, args...)`, `y -= f(args...)` is translated to `MinusEq(f)(y, args...)`, `y *= f(args...)` is translated to `MulEq(f)(y, args...)`, `y /= f(args...)` is translated to `DivEq(f)(y, args...)` and `y ⊻= f(args...)` is translated to `XorEq(f)(y, args...)`.

# ╔═╡ f6049c78-7468-47ce-a4a5-84fab34d115a
title2("How to create a new elementary reversible function")

# ╔═╡ b0f73825-bbb1-448c-b491-bf634fdd398a
md"To define a pair of elementary functions that **reverse to each other**,
1. declare two functions `f` and `g` that each of them defines a mapping ``\mathbb{R}^n \rightarrow \mathbb{R}^n``
2. use `@dual f g` to tell NiLang they are reversible to each other.
3. if you want to make `f` and `g` differentiable, you can specify backward rules on these two function by defining two mappings on ``\mathbb{G}^n\rightarrow \mathbb{G}^n``, where ``\mathbb{G}`` is a 2-tuple of ``\mathbb{R}`` (or `NiLang.AD.GVar`) in NiLang. It is similar to `ForwardDiff.Dual` (check [ForwardDiff](https://github.com/JuliaDiff/ForwardDiff.jl)) but defined for the backward pass.

To define a **self-reversible** elementary function
1. declare a functions `f` that defines a mapping ``\mathbb{R}^n \rightarrow \mathbb{R}^n``
2. use `@selfdual f`
3. define the backward rule on `f` to make it differentiable.
"

# ╔═╡ 648cdcd6-f4f5-461f-a525-4b350cae9eb0
example("defining a new elementary function")

# ╔═╡ d6b1abd6-749d-4591-99e8-64aaa9199ab5
md"""
One can use the invertibility checker to check if the function is really reversible (under a certain tolerance `NiLangCore.GLOBAL_ATOL[]` = $(NiLangCore.GLOBAL_ATOL[])).
"""

# ╔═╡ f502b8c1-9b80-4e67-80e8-a64ddb88fb0f
@test NiLang.check_inv(M.new_forward, (3.0,))

# ╔═╡ 0bce342e-9a8e-4005-8b88-82da2d2c7163
md"""
To check of the gradients are properly defined, one can use `NiLang.AD.check_grad`
"""

# ╔═╡ fd9cf757-2698-4886-9f0a-c6c23ff0d331
@test NiLang.AD.check_grad(M.new_forward, (3.0,); iloss=1)

# ╔═╡ dda6652a-d063-4511-8041-e869bb88ca26
@test NiLang.AD.check_grad(M.new_backward, (3.0,); iloss=1)

# ╔═╡ f6bfa015-c101-45e8-995c-2bb6a3b7dc7d
title1("Types and arrays")

# ╔═╡ 8651d7ec-6bcd-4dbe-a062-c4bde32e5e91
md"
A Julia type can be accessed in NiLang if its default constructor is not overloaded, because NiLang requires the default constructor to \"modify\" a field of a immutable type. For example, one can modify the real and imaginary part of a complex number directly in NiLang."

# ╔═╡ 9f5f9de3-9558-4c18-9d98-b77d19b570ec
example("Complex valued log")

# ╔═╡ 6dfcfa19-f78f-4dac-89f7-d3c5dbe17987
@i function complex_log(y!::Complex{T}, x::Complex{T}) where T
    n ← zero(T)
    n += abs(x)

    y!.re += log(n)
    y!.im += angle(x)

    n -= abs(x)
    n → zero(T)
end; @test complex_log(0.0im, 3.0+2.0im)[1] ≈ log(3.0+2.0im)

# ╔═╡ 99d6fe7b-d704-48f3-b115-2b3159a78068
md"One can also modify the `Array` directly"

# ╔═╡ c21d81c3-981f-4472-ad61-d1661bfe5c4e
example("Implementing \"axpy\" function")

# ╔═╡ 1950ff70-54eb-4ece-a26d-a23fd0e90f5a
@i function arrayaxpy!(y!::Vector{T}, a::T, x::Vector{T}) where T
    for i=1:length(x)
		y![i] += a * x[i]
	end
end; @test arrayaxpy!(zeros(10), 2.0, collect(1.0:10.0))[1] ≈ collect(2.0:2.0:20.0)

# ╔═╡ 21458f81-9007-46f8-92e0-7a17c60beb36
md"To modify an element of a `Tuple`, we need to use a different style to avoid confusion with array"

# ╔═╡ 7813f4ce-6e98-45f3-94a8-7f5981129f2b
@i function tupleaxpy!(y!::NTuple{N,T}, a::T, x::NTuple{N,T}) where {N, T}
    for i=1:length(x)
		(y! |> tget(i)) += a * (x |> tget(i))
	end
end; @test tupleaxpy!((0,0,0), 2, (1,2,3))[1] == (2, 4, 6)

# ╔═╡ 59ec7cb7-6011-456d-9f57-a55bb8ea51a0
md"Here `data |> function` is called a *dataview* of an object, it defines a modifiable view for a field of a data."

# ╔═╡ c7f2f786-d9ff-43b4-baef-5a48d611aa1e
title1("Functions")

# ╔═╡ d6519029-231d-4c63-b47d-684462dab287
md"""
Macros on functions are partly supported, NiLang will render the body of the function first and pass it to the next macro. One can define an inline function like
```julia
@i @inline function f(args...)
	...
end
```

However, the generated functions are not yet supported.
"""

# ╔═╡ b3147f49-c0f8-40cc-a794-282f6950b392
md"
Functions can not be chained.
"

# ╔═╡ 3e1a9006-3930-44f6-9eed-584d78937a57
@test_throws LoadError macroexpand(NiLang, :(@i function f6a(x)
	INV(DEC(x))
end))

# ╔═╡ b5386441-617a-4cc9-a3e8-5d4067e9554e
md"For single input single out function, one can use `-`, `'` and the pipeline `var |> reversible function or data field getter` to make the function shorter."

# ╔═╡ bc19ff2a-d3ba-4a4b-959a-fe381015f2ef
@i function f6b(x)
	INV(-x |> SubConst(4.0))
end

# ╔═╡ 47c95524-cef5-4654-b9ed-324472707ef0
f6b(3.0)

# ╔═╡ 03b31baa-9925-4fcb-a95c-7e48bd7708ca
md"Why the result is not `1/(-3-4)`? This is because the `INV` operation operates on the *dataview* `-x |> SubConst(4.0)`,  rather than on `x` directly. Which mean when it computes the result `-1/7`, it tried to assign this value back to `-x |> SubConst(4.0)`. So it adds 4 to the result and negate it to make the change propagate to the original variable."

# ╔═╡ aacf63a2-9708-40db-8928-049621a7bbc4
md"## Control flows"

# ╔═╡ ad0097e7-c8ad-457a-82a9-18b998a9e9fb
md"""
#### If statement

The condition expression in `if` statement contains two parts, one is precondition and another is postcondition.

$(
twocol(md"
```julia
if (precondition1[, postcondition1])
	...
end
```
", md"
```julia
if (postcondition1[, precondition1])
	~(...)
end
```
")
)

where `...` are statements and `~(...)` are the backward execution of them.
"""

# ╔═╡ 94cd1345-3132-4882-86fe-d2429f610d1d
md"""If no postcondition is provided, it means the precondition is the same as the postcondition. It is translated to `if (cond, cond) ... end`. `elseif` is also supported, to avoid potential complications, we do not discuss it here."""

# ╔═╡ 4a558bd3-6e42-4c61-bd23-888b7f33ae25
@i function f7a(x)
	if x > 1
		x -= 1
	end
end; @test_throws InvertibilityError f7a(1.2)

# ╔═╡ 004f727a-e0c8-49cb-8858-dfdf4d3ac57a
@i function f7b(x, branch_keeper)
	branch_keeper ⊻= x > 1
	if (x > 1, branch_keeper)
		x += 1
	end
end; @test f7b(1.2, false)[1] == 2.2

# ╔═╡ 4c03cde9-b643-40ff-b275-f1795f88949e
title2("For statement")

# ╔═╡ fae7c74e-d25e-4c1e-ac97-199e6dae3365
md"""
The reversible `for` statement is similar to its irreversible counterpart.

$(
twocol(md"
```julia
for iter = start:step:stop
	...
end
```
", md"
```julia
for iter = stop:-step:start
	~(...)
end
```
")
)
"""

# ╔═╡ 2f2b24ea-66d0-4b3b-a460-53b6b3f28ef0
md"The iterator length should not be changed during the iteration."

# ╔═╡ e9dbb64a-27b9-443a-b917-69d55c290235
@i function f7c(x, y)
	for i=1:length(x)
		POP!(x, y[i])
	end
end; @test_throws InvertibilityError f7c([1,2,3], [0,0,0])

# ╔═╡ 0d56ce96-81a5-4102-acbc-7d88f80adcb3
md"There is an `InvertibilityError` because the length of the `x` has been changed. The inverse execution will give incorrect result."

# ╔═╡ 95c41bd1-e50a-42e8-93c3-3b754a458c13
title2("While statement")

# ╔═╡ b8629aeb-6c9a-44ed-87a1-9ab22d9485ed
md"""
The reversible `while` statement starts with the `@from` macro.

$(
twocol(md"
```julia
@from condition1 while condition2
	...
end
```
", md"
```julia
@from !(condition2) while !(condition1)
	...
end
```
")
)
"""

# ╔═╡ 3b211406-041f-4b41-acae-3958e4a37224
md"where the `condition1` in the forward pass is a condition that holds before entering the loop body, but broken at the first iteration, while `condition2` is just a normal while condition. In the backward pass, `!(condition1)` becomes the criteria to break the loop."

# ╔═╡ 62522772-cb59-4d13-acdd-d5067b223910
@i function f7d(x, i)
	@from i==0 while i<10
		i += 1
		x += 1
	end
end

# ╔═╡ 2ba68a0f-6e36-4ea2-a91d-6af43741bad1
f7d(1, 0)

# ╔═╡ 75cebaf1-38de-475f-892e-346fd2b46f6f
@test (~f7d)(11, 10) == (1, 0)

# ╔═╡ 72dcf2fe-eb48-4dee-8121-efafc87637e3
title2("Compute-copy-uncompute statement")

# ╔═╡ 84321198-93d4-4d22-8c0f-a5a10b884e1f
md"The *compute-copy-uncompute* statement is a widely used design pattern in reversible programming.
We compute the forward pass for the result, then we copy the result to the output variable, and run the backward pass to erase intermediate results.

For example, to compute `y = x * exp(k)`, we might write the following code"

# ╔═╡ 32244789-afbf-4215-97cd-15483f438eee
@i function f7e(y, x, k)
	expk ← zero(k)
	expk += exp(k)
	y += x * expk
	# uncompute the ancilla and deallocate it
	expk -= exp(k)
	expk → zero(k)
end; @test f7e(0.0, 2.0, 3.0)[1] ≈ 2.0 * exp(3.0)

# ╔═╡ 4b7f0baf-0316-4da7-9ded-50c064ddbaa3
md"It is equivalent to the following statement that generates the backward pass automatically for you."

# ╔═╡ 0e02952c-7589-4606-b006-16a9f3e52ae1
@i function f7f(y, x, k)
	# record the forward pass
	@routine begin
		expk ← zero(k)
		expk += exp(k)
	end
	y += x * expk
	# reverse execute the recorded the program
	~@routine
end; @test f7f(0.0, 2.0, 3.0)[1] ≈ 2.0 * exp(3.0)

# ╔═╡ f0904d3f-1bf1-459c-9959-b53c0f774e3f
example("Computing Fibonacci numbers")

# ╔═╡ 19bb2af5-2a67-453d-82b0-7d3059b1fa47
md"The sequence of Fibonacci numbers are: 1, 1, 2, 3, 5, 8"

# ╔═╡ 5b5858bf-63ac-4e31-a516-055a9cd18ffe
@i function rfib(out!, n::T) where T
    n1 ← zero(T)
    n2 ← zero(T)
    @routine begin
        n1 += n - 1
        n2 += n - 2
    end
    if (value(n) <= 2, ~)
        out! += 1
    else
        rfib(out!, n1)
        rfib(out!, n2)
    end
    ~@routine
end

# ╔═╡ 95060588-f24b-4eeb-9b0b-ed7159962a3c
@test rfib(0, 6)[1] == 8

# ╔═╡ 91f8cfc6-e261-4945-8506-eed8caa607c2
title1("Multi-threading and CUDA")

# ╔═╡ c82b3b5c-c4e2-4bf6-b4ec-0d05ba9a669b
@i function f8a(y::Vector, x::Vector)
	# check the size of `x` and `y`. `@assert` is not a valid statement in NiLang, so one should decorate it with `@safe` to tell the compiler, doing this is safe, do not check this statement.
	@safe @assert length(x) == length(y)
	@threads for i=1:length(y)
		y[i] += exp(x[i])
	end
end; @test f8a(zeros(3), [1.0, 2.0, 3.0])[1] ≈ [exp(1.0), exp(2.0), exp(3.0)]

# ╔═╡ 8c93a773-edc0-4ec2-88ef-1b58b7deddc5
title2("Shared read-write in parallel computing and autodiff")

# ╔═╡ 16d08950-0575-4a4b-afc8-11ddca3198c7
md"Let's take a look at a correct prallel code that compute `exp(x)` and broadcast it to the output."

# ╔═╡ 7c594d19-59fc-433a-bffa-c63bad46869e
@i function f8c(y::Vector, x::Real, z::Vector)
	@safe @assert length(z) == length(y)
	@threads for i=1:length(y)
		y[i] += x * z[i]
	end
end; @test f8c(zeros(3), 2.0, [1.0, 2.0, 3.0])[1] ≈ [2.0, 4.0, 6.0]

# ╔═╡ 24421e48-a5fa-464e-b4ea-b55dabc7c6f7

# ╔═╡ a7d47e83-7f44-49d0-a43d-e01316fc6eba
title1("Performance Tips")

# ╔═╡ 45985244-adbf-4d6d-9732-a963cca62212
title2("Remove reversibility check")

# ╔═╡ 83d7e75f-7273-4c6a-bec1-a2180ebc3fb9
md"This can be done by putting an `@invcheckoff` before a code block."

# ╔═╡ 7fb05c65-f47c-430a-b588-c2f9bade40a9
example("computing the exp function by Taylor expansion")

# ╔═╡ 14c0caa1-51ea-448c-a7dc-d06e34dd0895
md"Note: this is not a clever implementation. There is an approach of defining it without allocation."

# ╔═╡ 457d07bb-e999-413e-8f29-58714670296f
@i function exp_with_reversibility_check(y::T, x::T) where T
	@routine begin
		N ← 1_000
		anc ← zeros(T, N)
		anc[1] += 1
		anc_y ← T(2.0)
		for i=2:N
			@routine begin
				temp ← zero(T)
				temp += x * anc[i-1]
			end
			anc[i] += temp / i
			anc_y += anc[i]
			~@routine
		end
	end
	y += anc_y
	~@routine
end

# ╔═╡ ac53eac0-1a59-4407-8bf6-3d8b966a9bff
@benchmark exp_with_reversibility_check(0.0, 1.0) seconds=0.3

# ╔═╡ 85c8ac7b-54f5-47dc-bd50-e78ffd6cf1cf
@i function exp_without_reversibility_check(y::T, x::T) where T
	@routine @invcheckoff begin
		N ← 1_000
		anc ← zeros(T, N)
		anc[1] += 1
		anc_y ← T(2.0)
		for i=2:N
			@routine begin
				temp ← zero(T)
				temp += x * anc[i-1]
			end
			anc[i] += temp / i
			anc_y += anc[i]
			~@routine
		end
	end
	y += anc_y
	~@routine
end

# ╔═╡ 95c55847-0591-4f7f-b9a1-aa974ccfef69
@benchmark exp_without_reversibility_check(0.0, 1.0) seconds=0.3

# ╔═╡ f80353d6-0dfe-4b0a-a1af-655d344473bf
title1("Resources and libraries")

# ╔═╡ 7ce31932-0447-4445-99aa-7ebced7d0bad
TableOfContents()

# ╔═╡ Cell order:
# ╟─2061b434-0ad1-46eb-a0c7-1a5f432bfa62
# ╟─a4e76427-f051-4b29-915a-fdfce3a299bb
# ╟─c2c7b4d4-f8c9-4ebf-8da2-0103f03136e7
# ╟─605872cf-f3fd-462e-a2b1-7d1c5ae45efd
# ╟─fb3dee44-5fa9-4773-8b7f-a83c44358545
# ╠═d941d6c2-55bf-11eb-0002-35c7474e4050
# ╠═70088425-6779-4a2d-ba6d-b0a34c8e93a6
# ╟─f0e94247-f615-472b-8218-3fa287b38aa1
# ╠═af738f89-3214-429c-9c7d-18a6ea0d9401
# ╟─4d98d529-1e05-49be-9209-f0d9fcc206f7
# ╠═48d7ebc1-5def-4a57-9ec1-3fc370a4543f
# ╟─e839547b-f47e-4a4e-b4d9-2c22921d80e4
# ╟─60575978-081a-4bca-a3ed-2b51cd6abc92
# ╠═f98305cb-4ba2-404a-a5c3-65510e059504
# ╟─e8cd6667-597f-458b-8465-1822e09a7891
# ╟─20145d75-004a-4c2f-b7ff-c400ca846d42
# ╟─b6dcd18c-606f-4340-b2ec-163e8bad03f5
# ╟─2581aa33-1dc5-40b1-aa9f-6a11cc750c93
# ╟─a1a29f34-f8a9-4e9f-9afe-7d0096771440
# ╟─90bd6ad4-3dd8-4e7c-b445-aed1d248a2ec
# ╠═c0259e48-1973-486c-a828-1fcd3e4331c6
# ╟─7bcc1bc6-5b2c-4aae-b71d-52482ae130e9
# ╠═fd5d65e7-3d83-45dd-87d7-ff1e1e769430
# ╟─4601df35-679f-465d-9191-c18748b2fd83
# ╟─4fc72b9d-19a2-40f1-a4a8-5e97d3d5e529
# ╠═a0fde16f-8454-4f5c-a29c-a9e415c0c311
# ╟─633ff8f3-8d93-4f73-bec2-c42070e6ece9
# ╠═57f2d890-b5a0-47b7-9e3d-af4d03b10605
# ╟─34063cd0-171e-46ce-80dd-52a341fa50a1
# ╠═5d5d01db-8ff9-434c-8771-1fec6393e1fb
# ╠═10d85a50-f2f9-403e-8f6c-baef61cf702a
# ╟─b52648bf-a28a-48af-8912-31729d943ce0
# ╟─f3503e6d-ded7-4797-9c2f-cb5a7d429949
# ╟─7ba86798-daba-49ef-a259-da4db353fdd8
# ╠═6e88ede9-879e-42e9-9a7e-2aa3d4837e5d
# ╠═94cc42f7-6beb-453e-90cb-952d6f497a8e
# ╟─93936612-1447-4114-b864-aba43adef4bd
# ╠═b2add70c-c5d6-4e0f-a153-43e21a197181
# ╠═5c03d5a5-99f0-4efd-9a32-ce6d7c2b266c
# ╟─f45db10f-a836-40f3-9d8d-054ea6540e87
# ╟─1903563e-ccc2-44d9-9dbe-e5dede275b3c
# ╟─7ec577d4-1e19-484e-ac1f-5da3f47d1ec4
# ╟─583f2585-15a3-47c6-a70e-e2f002754028
# ╠═60e6ff80-3593-4ae4-a273-914847f692db
# ╠═9e5cfd68-b58d-4d83-aae2-447e5f805c97
# ╠═dc85a942-cf52-4405-ad03-32a768e1b6e7
# ╟─5cdd346b-10a5-485c-ba78-4c0b3cb0e02f
# ╠═af9287b7-6131-46f6-beb8-6885e55e1975
# ╠═e20eeabf-1c80-431e-8cfc-4d1b79c52b5a
# ╟─e45d53b9-6f02-4293-ab4f-85bd751993ad
# ╟─90d30eea-53de-48a0-9700-ff35681fdf38
# ╠═390f58a5-6f5f-4d3a-bb16-ba04e43a07e7
# ╠═2b57443e-a516-434b-be86-80616a98e2f5
# ╟─fc2e27f9-b7ba-44cc-a953-6745548ad733
# ╠═fc744931-360b-4478-9f77-c50f048de243
# ╠═9a152b36-f377-44da-9700-ca9e05e365ff
# ╟─0863bd06-cc70-4dde-b3b2-0a466805a356
# ╟─a0bae195-04e1-4642-9e14-fe4691e0906b
# ╠═20d6e8a0-2cf5-48ad-9549-60506b42b970
# ╠═7dee5748-ed73-4e13-aa80-7a50efbc8449
# ╟─4c719e9b-641e-404e-9ab7-59e89135f3ba
# ╠═77947e00-42c3-4c9e-b62a-b4b29489db43
# ╠═2614127d-34fb-4c3d-b678-42693f3c9341
# ╟─8f169235-3bd1-4cc4-a083-79736d306ad5
# ╠═dfc9d305-5bce-4555-bfa3-d8d61fe4ca09
# ╟─c4cd9f88-9cd6-4364-b016-78f90aba6a66
# ╟─f6049c78-7468-47ce-a4a5-84fab34d115a
# ╟─b0f73825-bbb1-448c-b491-bf634fdd398a
# ╟─648cdcd6-f4f5-461f-a525-4b350cae9eb0
# ╠═8c2c4fa6-172f-4dde-a279-5d0aecfdbe46
# ╟─d6b1abd6-749d-4591-99e8-64aaa9199ab5
# ╠═f502b8c1-9b80-4e67-80e8-a64ddb88fb0f
# ╟─0bce342e-9a8e-4005-8b88-82da2d2c7163
# ╠═fd9cf757-2698-4886-9f0a-c6c23ff0d331
# ╠═dda6652a-d063-4511-8041-e869bb88ca26
# ╟─f6bfa015-c101-45e8-995c-2bb6a3b7dc7d
# ╟─8651d7ec-6bcd-4dbe-a062-c4bde32e5e91
# ╟─9f5f9de3-9558-4c18-9d98-b77d19b570ec
# ╠═6dfcfa19-f78f-4dac-89f7-d3c5dbe17987
# ╟─99d6fe7b-d704-48f3-b115-2b3159a78068
# ╟─c21d81c3-981f-4472-ad61-d1661bfe5c4e
# ╠═1950ff70-54eb-4ece-a26d-a23fd0e90f5a
# ╟─21458f81-9007-46f8-92e0-7a17c60beb36
# ╠═7813f4ce-6e98-45f3-94a8-7f5981129f2b
# ╟─59ec7cb7-6011-456d-9f57-a55bb8ea51a0
# ╟─c7f2f786-d9ff-43b4-baef-5a48d611aa1e
# ╟─d6519029-231d-4c63-b47d-684462dab287
# ╟─b3147f49-c0f8-40cc-a794-282f6950b392
# ╠═3e1a9006-3930-44f6-9eed-584d78937a57
# ╟─b5386441-617a-4cc9-a3e8-5d4067e9554e
# ╠═bc19ff2a-d3ba-4a4b-959a-fe381015f2ef
# ╠═47c95524-cef5-4654-b9ed-324472707ef0
# ╟─03b31baa-9925-4fcb-a95c-7e48bd7708ca
# ╟─aacf63a2-9708-40db-8928-049621a7bbc4
# ╟─ad0097e7-c8ad-457a-82a9-18b998a9e9fb
# ╟─94cd1345-3132-4882-86fe-d2429f610d1d
# ╠═4a558bd3-6e42-4c61-bd23-888b7f33ae25
# ╠═004f727a-e0c8-49cb-8858-dfdf4d3ac57a
# ╟─4c03cde9-b643-40ff-b275-f1795f88949e
# ╟─fae7c74e-d25e-4c1e-ac97-199e6dae3365
# ╟─2f2b24ea-66d0-4b3b-a460-53b6b3f28ef0
# ╠═e9dbb64a-27b9-443a-b917-69d55c290235
# ╟─0d56ce96-81a5-4102-acbc-7d88f80adcb3
# ╟─95c41bd1-e50a-42e8-93c3-3b754a458c13
# ╟─b8629aeb-6c9a-44ed-87a1-9ab22d9485ed
# ╟─3b211406-041f-4b41-acae-3958e4a37224
# ╠═62522772-cb59-4d13-acdd-d5067b223910
# ╠═2ba68a0f-6e36-4ea2-a91d-6af43741bad1
# ╠═75cebaf1-38de-475f-892e-346fd2b46f6f
# ╟─72dcf2fe-eb48-4dee-8121-efafc87637e3
# ╟─84321198-93d4-4d22-8c0f-a5a10b884e1f
# ╠═32244789-afbf-4215-97cd-15483f438eee
# ╟─4b7f0baf-0316-4da7-9ded-50c064ddbaa3
# ╠═0e02952c-7589-4606-b006-16a9f3e52ae1
# ╟─f0904d3f-1bf1-459c-9959-b53c0f774e3f
# ╟─19bb2af5-2a67-453d-82b0-7d3059b1fa47
# ╠═5b5858bf-63ac-4e31-a516-055a9cd18ffe
# ╠═95060588-f24b-4eeb-9b0b-ed7159962a3c
# ╟─91f8cfc6-e261-4945-8506-eed8caa607c2
# ╠═0e1ba158-a6bc-401c-9ba7-ed78020ad068
# ╠═c82b3b5c-c4e2-4bf6-b4ec-0d05ba9a669b
# ╟─8c93a773-edc0-4ec2-88ef-1b58b7deddc5
# ╠═16d08950-0575-4a4b-afc8-11ddca3198c7
# ╠═7c594d19-59fc-433a-bffa-c63bad46869e
# ╠═24421e48-a5fa-464e-b4ea-b55dabc7c6f7
# ╟─a7d47e83-7f44-49d0-a43d-e01316fc6eba
# ╟─45985244-adbf-4d6d-9732-a963cca62212
# ╟─83d7e75f-7273-4c6a-bec1-a2180ebc3fb9
# ╟─7fb05c65-f47c-430a-b588-c2f9bade40a9
# ╟─14c0caa1-51ea-448c-a7dc-d06e34dd0895
# ╠═457d07bb-e999-413e-8f29-58714670296f
# ╠═3199a048-7b39-40f8-8183-6a54cccd91b6
# ╠═ac53eac0-1a59-4407-8bf6-3d8b966a9bff
# ╠═85c8ac7b-54f5-47dc-bd50-e78ffd6cf1cf
# ╠═95c55847-0591-4f7f-b9a1-aa974ccfef69
# ╟─f80353d6-0dfe-4b0a-a1af-655d344473bf
# ╟─7ce31932-0447-4445-99aa-7ebced7d0bad


text:
"""
Given a remote repo URL and a git tree spec, get a `Project` object
for the project file in that tree and a hash string for the tree.
"""
# function get_project(remote_url::AbstractString, tree_spec::AbstractString)
#     # TODO?: use raw file downloads for GitHub/GitLab
#     mktempdir(mkpath("packages")) do tmp
#         # bare clone the package repo
#         @debug("bare clone the package repo")
#         repo = LibGit2.clone(remote_url, joinpath(tmp, "repo"), isbare=true)
#         tree = try
#             LibGit2.GitObject(repo, tree_spec)
#         catch err
#             err isa LibGit2.GitError && err.code == LibGit2.Error.ENOTFOUND || rethrow(err)
#             error("$remote_url: git object $(repr(tree_spec)) could not be found")
#         end
#         tree isa LibGit2.GitTree || (tree = LibGit2.peel(LibGit2.GitTree, tree))
#
#         # check out the requested tree
#         @debug("check out the requested tree")
#         tree_path = abspath(tmp, "tree")
#         GC.@preserve tree_path begin
#             opts = LibGit2.CheckoutOptions(
#                 checkout_strategy = LibGit2.Consts.CHECKOUT_FORCE,
#                 target_directory = Base.unsafe_convert(Cstring, tree_path)
#             )
#             LibGit2.checkout_tree(repo, tree, options=opts)
#         end
#
#         # look for a project file in the tree
#         @debug("look for a project file in the tree")
#         project_file = Pkg.Types.projectfile_path(tree_path)
#         project_file !== nothing && isfile(project_file) ||
#             error("$remote_url: git tree $(repr(tree_spec)) has no project file")
#
#         # parse the project file
#         @debug("parse the project file")
#         project = Pkg.Types.read_project(project_file)
#         project.name === nothing &&
#             error("$remote_url $(repr(tree_spec)): package has no name")
#         project.uuid === nothing &&
#             error("$remote_url $(repr(tree_spec)): package has no UUID")
#         project.version === nothing &&
#             error("$remote_url $(repr(tree_spec)): package has no version")
#
#         return project, string(LibGit2.GitHash(tree))
#     end
# end

const julia_uuid = "1222c4b2-2114-5bfd-aeef-88e4692bbb3e"

struct RegBranch
    name::String
    version::VersionNumber
    branch::String

    metadata::Dict{String,Any} # "error", "warning", kind etc.

    function RegBranch(pkg::Pkg.Types.Project, branch::AbstractString)
        new(pkg.name, pkg.version, branch, Dict{String,Any}())
    end
end

get_backtrace(ex) = sprint(Base.showerror, ex, catch_backtrace())

function write_registry(registry_path::AbstractString, reg::RegistryData)
    open(registry_path, "w") do io
        TOML.print(io, reg)
    end
end

# error in regbr.metadata["errors"]
# warning in regbr.metadata["warning"]
# version labels for the PR in in regbr.metadata["labels"]
function check_version!(regbr::RegBranch, existing::Vector{VersionNumber})
    ver = regbr.version
    if ver == v"0"
        regbr.metadata["error"] = "Package version must be greater than 0.0.0"
        return regbr
    end

    @assert issorted(existing)
    if isempty(existing)
        push!(get!(regbr.metadata, "labels", String[]), "new package")
        if !(ver in [v"0.0.1", v"0.1", v"1"])
            regbr.metadata["warning"] =
                """This looks like a new registration that registers version $ver.
                Ideally, you should register an initial release with 0.0.1, 0.1.0 or 1.0.0 version numbers"""
        end
        return regbr
    else
        idx = searchsortedlast(existing, ver)
        if idx <= 0
            regbr.metadata["error"] = "Version $ver less than least existing version $(existing[1])"
            return regbr
        end

        prv = existing[idx]
        if ver == prv
            regbr.metadata["error"] = "Version $ver already exists"
            return regbr
        end
        nxt = if ver.major != prv.major
            push!(get!(regbr.metadata, "labels", String[]), "major release")
            VersionNumber(prv.major+1, 0, 0)
        elseif ver.minor != prv.minor
            push!(get!(regbr.metadata, "labels", String[]), "minor release")
            VersionNumber(prv.major, prv.minor+1, 0)
        else
            push!(get!(regbr.metadata, "labels", String[]), "patch release")
            VersionNumber(prv.major, prv.minor, prv.patch+1)
        end
        if ver > nxt
            regbr.metadata["warning"] = "Version $ver skips over $nxt"
            return regbr
        end
    end

    return regbr
end

findpackageerror(name::AbstractString, uuid::Base.UUID, regdata::Array{RegistryData}) =
    findpackageerror(name, string(uuid), regdata)

function findpackageerror(name::AbstractString, u::AbstractString, regdata::Array{RegistryData})
    for _registry_data in regdata
        if haskey(_registry_data.packages, u)
            name_in_reg = _registry_data.packages[u]["name"]
            if name_in_reg != name
                return "Error in (Julia)Project.toml: UUID $u refers to package '$name_in_reg' in registry but Project.toml has '$name'"
            end
            return nothing
        end
    end

    if haskey(BUILTIN_PKGS, name)
        if BUILTIN_PKGS[name] != u
            return "Error in (Julia)Project.toml: UUID $u for package $name should be $(BUILTIN_PKGS[k])"
        end
    else
        return "Error in (Julia)Project.toml: Package '$name' with UUID: $u not found in registry or stdlib"
    end

    nothing
end

import Pkg.Types: VersionRange, VersionBound, VersionSpec

function versionrange(lo::VersionBound, hi::VersionBound)
    lo.t == hi.t && (lo = hi)
    return VersionRange(lo, hi)
end

# Code copied from Pkg found in dev julia (Current version is 1.1)
# TODO: Remove after moving to julia 1.2
"""
    compress_versions(pool::Vector{VersionNumber}, subset::Vector{VersionNumber})
Given `pool` as the pool of available versions (of some package) and `subset` as some
subset of the pool of available versions, this function computes a `VersionSpec` which
includes all versions in `subset` and none of the versions in its complement.
"""
function compress_versions(pool::Vector{VersionNumber}, subset::Vector{VersionNumber})
    # Explicitly drop prerelease/build numbers, as those can confuse this.
    # TODO: Rewrite all this to use VersionNumbers instead of VersionBounds
    drop_build_prerelease(v::VersionNumber) = VersionNumber(v.major, v.minor, v.patch)
    pool = drop_build_prerelease.(pool)
    subset = sort!(drop_build_prerelease.(subset))

    complement = sort!(setdiff(pool, subset))
    ranges = VersionRange[]
    @label again
    isempty(subset) && return VersionSpec(ranges)
    a = first(subset)
    for b in reverse(subset)
        a.major == b.major || continue
        for m = 1:4
            lo = VersionBound((a.major, a.minor, a.patch, a.build)[1:m]...)
            for n = 1:4
                hi = VersionBound((b.major, b.minor, b.patch, b.build)[1:n]...)
                r = versionrange(lo, hi)
                if !any(v in r for v in complement)
                    filter!(!in(r), subset)
                    push!(ranges, r)
                    @goto again
                end
            end
        end
    end
end

function compress_versions(pool::Vector{VersionNumber}, subset)
    compress_versions(pool, filter(in(subset), pool))
end

import Pkg.Compress.load_versions

function compress(path::AbstractString, uncompressed::Dict,
    versions::Vector{VersionNumber} = load_versions(path))
    inverted = Dict()
    for (ver, data) in uncompressed, (key, val) in data
        val isa TOML.TYPE || (val = string(val))
        push!(get!(inverted, key => val, VersionNumber[]), ver)
    end
    compressed = Dict()
    for ((k, v), vers) in inverted
        for r in compress_versions(versions, sort!(vers)).ranges
            get!(compressed, string(r), Dict{String,Any}())[k] = v
        end
    end
    return compressed
end

function save(path::AbstractString, uncompressed::Dict,
    versions::Vector{VersionNumber} = load_versions(path))
    compressed = compress(path, uncompressed)
    open(path, write=true) do io
        TOML.print(io, compressed, sorted=true)
    end
end

# ---- End of code copied from Pkg

function find_package_in_registry(pkg::Pkg.Types.Project,
                                  package_repo::AbstractString,
                                  registry_file::AbstractString,
                                  registry_path::AbstractString,
                                  registry_data::RegistryData,
                                  regbr::RegBranch)
    uuid = string(pkg.uuid)
    if haskey(registry_data.packages, uuid)
        package_data = registry_data.packages[uuid]
        if package_data["name"] != pkg.name
            err = "Changing package names not supported yet"
            @debug(err)
            regbr.metadata["error"] = err
            return nothing, regbr
        end
        package_path = joinpath(registry_path, package_data["path"])
        repo = TOML.parsefile(joinpath(package_path, "Package.toml"))["repo"]
        if repo != package_repo
            err = "Changing package repo URL not allowed, please submit a pull request with the URL change to the target registry and retry."
            @debug(err)
            regbr.metadata["error"] = err
            return nothing, regbr
        end
        regbr.metadata["kind"] = "New version"
    else
        @debug("Package with UUID: $uuid not found in registry, checking if UUID was changed")
        for (k, v) in registry_data.packages
            if v["name"] == pkg.name
                err = "Changing UUIDs is not allowed"
                @debug(err)
                regbr.metadata["error"] = err
                return nothing, regbr
            end
        end

        @debug("Creating directory for new package $(pkg.name)")
        package_path = joinpath(registry_path, package_relpath(registry_data, pkg))
        mkpath(package_path)

        @debug("Adding package UUID to registry")
        push!(registry_data, pkg)
        write_registry(registry_file, registry_data)
        regbr.metadata["kind"] = "New package"
    end

    return package_path, regbr
end

function update_package_file(pkg::Pkg.Types.Project,
                             package_repo::AbstractString,
                             package_path::AbstractString)
    package_info = Dict("name" => pkg.name,
                        "uuid" => string(pkg.uuid),
                        "repo" => package_repo)
    package_file = joinpath(package_path, "Package.toml")
    open(package_file, "w") do io
        TOML.print(io, package_info; sorted=true,
            by = x -> x == "name" ? 1 : x == "uuid" ? 2 : 3)
    end
    nothing
end

function update_versions_file(pkg::Pkg.Types.Project,
                              package_path::AbstractString,
                              regbr::RegBranch,
                              tree_hash::AbstractString)
    versions_file = joinpath(package_path, "Versions.toml")
    versions_data = isfile(versions_file) ? TOML.parsefile(versions_file) : Dict()
    versions = sort!([VersionNumber(v) for v in keys(versions_data)])

    check_version!(regbr, versions)
    if get(regbr.metadata, "error", nothing) !== nothing
        return regbr
    end

    version_info = Dict{String,Any}("git-tree-sha1" => string(tree_hash))
    versions_data[string(pkg.version)] = version_info

    open(versions_file, "w") do io
        TOML.print(io, versions_data; sorted=true, by=x->VersionNumber(x))
    end
    nothing
end

function update_deps_file(pkg::Pkg.Types.Project,
                          package_path::AbstractString,
                          regbr::RegBranch,
                          regdata::Vector{RegistryData})
    if pkg.name in keys(pkg.deps)
        err = "Package $(pkg.name) mentions itself in `[deps]`"
        @debug(err)
        regbr.metadata["error"] = err
        return regbr
    end

    deps_file = joinpath(package_path, "Deps.toml")
    if isfile(deps_file)
        deps_data = Pkg.Compress.load(deps_file)
    else
        deps_data = Dict()
    end

    @debug("Verifying package name and uuid in deps")
    for (k, v) in pkg.deps
        err = findpackageerror(k, v, regdata)
        if err !== nothing
            @debug(err)
            regbr.metadata["error"] = err
            return regbr
        end
    end

    deps_data[pkg.version] = pkg.deps
    save(deps_file, deps_data)
    nothing
end

function update_compat_file(pkg::Pkg.Types.Project,
                            package_path::AbstractString,
                            regbr::RegBranch,
                            regdata::Vector{RegistryData},
                            regpaths::Vector{String})
    err = nothing
    for (p, v) in pkg.compat
        try
            ver = Pkg.Types.semver_spec(v)
            if p == "julia" && any(map(x->!isempty(intersect(Pkg.Types.VersionRange("0-0.6"),x)), ver.ranges))
                err = "Julia version < 0.7 not allowed in `[compat]`"
                @debug(err)
                break
            end
        catch ex
            if isdefined(ex, :msg)
                err = "Error in `[compat]`: $(ex.msg)"
                @debug(err)
                break
            else
                rethrow(ex)
            end
        end
    end

    if err !== nothing
        regbr.metadata["error"] = err
        return regbr
    end

    @debug("update package data: compat file")
    compat_file = joinpath(package_path, "Compat.toml")
    if isfile(compat_file)
        compat_data = Pkg.Compress.load(compat_file)
    else
        compat_data = Dict()
    end

    d = Dict()
    err = nothing
    for (n,v) in pkg.compat
        spec = Pkg.Types.semver_spec(v)
        if n == "julia"
            uuidofdep = julia_uuid
        else
            indeps = haskey(pkg.deps, n)
            inextras = haskey(pkg.extras, n)

            if indeps
                uuidofdep = string(pkg.deps[n])
            elseif inextras
                uuidofdep = string(pkg.extras[n])
            else
                err = "Package $n mentioned in `[compat]` but not found in `[deps]` or `[extras]`"
                @debug(err)
                break
            end

            err = findpackageerror(n, uuidofdep, regdata)
            if err !== nothing
                @debug(err)
                break
            end

            if inextras && !indeps
                @debug("$n is a test-only dependency; omitting from Compat.toml")
                continue
            end
        end

        versionsfileofdep = nothing
        for i=1:length(regdata)
            if haskey(regdata[i].packages, uuidofdep)
                pathofdep = regdata[i].packages[uuidofdep]["path"]
                versionsfileofdep = joinpath(regpaths[i], pathofdep, "Versions.toml")
                break
            end
        end
        # the call to map(versionrange, ) can be removed
        # once Pkg is updated to a version including
        # https://github.com/JuliaLang/Pkg.jl/pull/1181
        ranges = map(r->versionrange(r.lower, r.upper), spec.ranges)
        ranges = VersionSpec(ranges).ranges # this combines joinable ranges
        d[n] = length(ranges) == 1 ? string(ranges[1]) : map(string, ranges)
    end

    if err !== nothing
        regbr.metadata["error"] = err
        return regbr
    end

    compat_data[pkg.version] = d

    save(compat_file, compat_data)
    nothing
end

function get_registrator_tree_sha()
    regtreesha = nothing
    reg_pkgs = Pkg.Display.status(Pkg.Types.Context(),
                   [Pkg.PackageSpec("Registrator", Base.UUID("4418983a-e44d-11e8-3aec-9789530b3b3e"))])
    if !isempty(reg_pkgs)
        regtreesha = reg_pkgs[1].new.hash
    end
    if regtreesha === nothing
        regpath = abspath(joinpath(@__DIR__, "..", ".."))
        if isdir(joinpath(regpath, ".git"))
            regtreesha = LibGit2.head(regpath)
        else
            regtreesha = "unknown"
        end
    end

    return regtreesha
end

"""
    register(package_repo, pkg, tree_hash; registry, registry_deps, push, gitconfig)

Register the package at `package_repo` / `tree_hash` in `registry`.
Returns a `RegEdit.RegBranch` which contains information about the registration and/or any
errors or warnings that occurred.

# Arguments

* `package_repo::AbstractString`: the git repository URL for the package to be registered
* `pkg::Pkg.Types.Project`: the parsed (Julia)Project.toml file for the package to be registered
* `tree_hash::AbstractString`: the tree hash (not commit hash) of the package revision to be registered

# Keyword Arguments

* `registry::AbstractString="$DEFAULT_REGISTRY_URL"`: the git repository URL for the registry
* `registry_deps::Vector{String}=[]`: the git repository URLs for any registries containing
    packages depended on by `pkg`
* `push::Bool=false`: whether to push a registration branch to `registry` for consideration
* `gitconfig::Dict=Dict()`: dictionary of configuration options for the `git` command
"""
function register(
    package_repo::AbstractString, pkg::Pkg.Types.Project, tree_hash::AbstractString;
    registry::AbstractString = DEFAULT_REGISTRY_URL,
    registry_deps::Vector{<:AbstractString} = AbstractString[],
    push::Bool = false,
    force_reset::Bool = true,
    branch::String = registration_branch(pkg),
    cache::RegistryCache=REGISTRY_CACHE,
    gitconfig::Dict = Dict()
)
    # get info from package registry
    @debug("get info from package registry")
    package_repo = GitTools.normalize_url(package_repo)

    # return object
    regbr = RegBranch(pkg, branch)

    # get up-to-date clone of registry
    @debug("get up-to-date clone of registry")
    registry = GitTools.normalize_url(registry)
    registry_repo = get_registry(registry; gitconfig=gitconfig, force_reset=force_reset, cache=cache)
    registry_path = LibGit2.path(registry_repo)

    isempty(registry_deps) || @debug("get up-to-date clones of registry dependencies")
    registry_deps_paths = map(registry_deps) do registry
        LibGit2.path(get_registry(GitTools.normalize_url(registry); gitconfig=gitconfig, force_reset=force_reset, cache=cache))
    end

    clean_registry = true
    err = nothing
    try
        # branch registry repo
        @debug("branch registry repo")
        git = gitcmd(registry_path, gitconfig)
        run(pipeline(`$git checkout -f master`; stdout=devnull))
        if branch != "master"
            run(pipeline(`$git branch -f $branch`; stdout=devnull))
            run(pipeline(`$git checkout -f $branch`; stdout=devnull))
        end

        # find package in registry
        @debug("find package in registry")
        registry_file = joinpath(registry_path, "Registry.toml")
        registry_data = parse_registry(registry_file)
        package_path, regbr = find_package_in_registry(pkg, package_repo,
                                                       registry_file, registry_path,
                                                       registry_data, regbr)
        package_path === nothing && return regbr

        # update package data: package file
        @debug("update package data: package file")
        update_package_file(pkg, package_repo, package_path)

        # update package data: versions file
        @debug("update package data: versions file")
        r = update_versions_file(pkg, package_path, regbr, tree_hash)
        r === nothing || return r

        # update package data: deps file
        @debug("update package data: deps file")
        registry_deps_data = map(registry_deps_paths) do registry_path
            parse_registry(joinpath(registry_path, "Registry.toml"))
        end
        regdata = [registry_data; registry_deps_data]
        r = update_deps_file(pkg, package_path, regbr, regdata)
        r === nothing || return r

        # update package data: compat file
        @debug("check compat section")
        regpaths = [registry_path; registry_deps_paths]
        r = update_compat_file(pkg, package_path, regbr, regdata, regpaths)
        r === nothing || return r

        regtreesha = get_registrator_tree_sha()

        # commit changes
        @debug("commit changes")
        message = """
        $(regbr.metadata["kind"]): $(pkg.name) v$(pkg.version)

        UUID: $(pkg.uuid)
        Repo: $(package_repo)
        Tree: $(string(tree_hash))

        Registrator tree SHA: $(regtreesha)
        """
        run(pipeline(`$git add -- $package_path`; stdout=devnull))
        run(pipeline(`$git add -- $registry_file`; stdout=devnull))
        run(pipeline(`$git commit -m $message`; stdout=devnull))

        # push -f branch to remote
        if push
            @debug("push -f branch to remote")
            run(pipeline(`$git push -f -u origin $branch`; stdout=devnull))
        else
            @debug("skipping git push")
        end

        clean_registry = false
    catch ex
        println(get_backtrace(ex))
        regbr.metadata["error"] = "Unexpected error in registration"
    finally
        if clean_registry
            @debug("cleaning up possibly inconsistent registry", registry_path=showsafe(registry_path), err=showsafe(err))
            rm(registry_path; recursive=true, force=true)
        end
    end
    return regbr
end

struct RegisterParams
    package_repo::AbstractString
    pkg::Pkg.Types.Project
    tree_sha::AbstractString
    registry::AbstractString
    registry_deps::Vector{<:AbstractString}
    push::Bool
    gitconfig::Dict

    function RegisterParams(package_repo::AbstractString,
                            pkg::Pkg.Types.Project,
                            tree_sha::AbstractString;
                            registry::AbstractString=DEFAULT_REGISTRY_URL,
                            registry_deps::Vector{<:AbstractString}=[],
                            push::Bool=false,
                            gitconfig::Dict=Dict())
        new(package_repo, pkg, tree_sha,
            registry, registry_deps,
            push, gitconfig)
    end
end

register(regp::RegisterParams) = register(regp.package_repo, regp.pkg, regp.tree_sha;
                                          registry=regp.registry, registry_deps=regp.registry_deps,
                                          push=regp.push, gitconfig=regp.gitconfig)


text:
struct ObservableDiagram
    intervals::Observable{Vector{Point2f0}}
    ids::Observable{Vector{Int}}
    t_min::Observable{Float64}
    t_max::Observable{Float64}
    infinity::Observable{Float64}
end

function ObservableDiagram(t_min, t_max)
    return ObservableDiagram(Observable(Point2f0[]), Observable(Int[]), t_min, t_max, t_max)
end

function ObservableDiagram(diagram::Vector{<:PersistenceInterval})
    return ObservableDiagram([diagram])
end

function ObservableDiagram(diagrams)
    limits = PersistenceDiagrams.limits(diagrams)
    ids = reduce(vcat, fill(i, length(d)) for (i, d) in enumerate(diagrams))
    intervals = collect(Iterators.flatten(
        convert_arguments(Scatter, diag)[1] for diag in diagrams
    ))
    return ObservableDiagram(intervals, ids, limits...)
end

function Base.show(io::IO, od::ObservableDiagram)
    print(io, "ObservableDiagram with $(length(od)) intervals")
    return nothing
end

Base.length(od::ObservableDiagram) = length(od.intervals[])

function Base.empty!(od::ObservableDiagram)
    od.intervals[] = Point2f0[]
    od.ids[] = Int[]
    return od
end

function Base.append!(od::ObservableDiagram, intervals)
    od.intervals[] = append!(od.intervals[], convert_arguments(Scatter, intervals)[1])
    if !isempty(od.ids[])
        od.ids[] = append!(od.ids[], fill(od.ids[][end] + 1, length(intervals)))
    else
        od.ids[] = append!(od.ids[], fill(1, length(intervals)))
    end
    return od
end

function Base.push!(od::ObservableDiagram, interval, id=od.ids[][end] + 1)
    od.intervals[] = push!(od.intervals[], Point2f0(interval[1], interval[2]))
    od.ids[] = push!(od.ids[], id)
    return od
end

function edit_last!(od::ObservableDiagram, interval, id=od.ids[][end])
    od.intervals[][end] = Point2f0(interval[1], interval[2])
    od.intervals[] = od.intervals[]
    od.ids[][end] = id
    od.ids[] = od.ids[]
    return od
end

function Base.pop!(od::ObservableDiagram)
    pop!(od.intervals[])
    od.intervals[] = od.intervals[]
    pop!(od.ids[])
    return od.ids[] = od.ids[]
end

function AbstractPlotting.default_theme(scene::SceneLike, ::Type{<:Plot(ObservableDiagram)})
    return Theme(;
        palette=DEFAULT_PALETTE, color=[], infinity=nothing, persistence=false, gapwidth=0.1
    )
end

function transform_intervals(intervals, infinity, persistence)
    map(intervals) do interval
        if persistence
            interval -= Point2f0(0, interval[1])
        end
        if isinf(interval[2])
            interval = Point2f0(interval[1], infinity)
        end
        return interval
    end
end

function AbstractPlotting.plot!(p::Plot(ObservableDiagram))
    diagram = to_value(p[1])
    color = lift(p[:palette], p[:color], diagram.ids) do palette, color, ids
        scheme = PlotUtils.get_colorscheme(palette)
        if isempty(ids)
            RGB{Float64}[]
        else
            [scheme[mod1(get(color, id, id), length(scheme))] for id in ids]
        end
    end
    infinity = @lift isnothing($(p[:infinity])) ? $(diagram.infinity) : $(p[:infinity])
    diagrambackground!(
        p,
        diagram.t_min,
        diagram.t_max,
        infinity;
        gapwidth=p[:gapwidth],
        persistence=p[:persistence],
    )
    points = lift(transform_intervals, diagram.intervals, infinity, p[:persistence])
    return scatter!(p, points; color)
end


text:
# ---
# title: 765. Couples Holding Hands
# id: problem765
# author: Tian Jun
# date: 2020-10-31
# difficulty: Hard
# categories: Greedy, Union Find, Graph
# link: <https://leetcode.com/problems/couples-holding-hands/description/>
# hidden: true
# ---
#
# N couples sit in 2N seats arranged in a row and want to hold hands. We want to
# know the minimum number of swaps so that every couple is sitting side by side.
# A _swap_ consists of choosing **any** two people, then they stand up and
# switch seats.
#
# The people and seats are represented by an integer from `0` to `2N-1`, the
# couples are numbered in order, the first couple being `(0, 1)`, the second
# couple being `(2, 3)`, and so on with the last couple being `(2N-2, 2N-1)`.
#
# The couples' initial seating is given by `row[i]` being the value of the
# person who is initially sitting in the i-th seat.
#
# **Example 1:**
#
#
#
#     Input: row = [0, 2, 1, 3]
#     Output: 1
#     Explanation: We only need to swap the second (row[1]) and third (row[2]) person.
#
#
# **Example 2:**
#
#
#
#     Input: row = [3, 2, 0, 1]
#     Output: 0
#     Explanation: All couples are already seated side by side.
#
#
# **Note:**
#
#   1. `len(row)` is even and in the range of `[4, 60]`.
#   2. `row` is guaranteed to be a permutation of `0...len(row)-1`.
#
#
## @lc code=start
using LeetCode

## add your code here:
## @lc code=end


text:
using Distributions
using StatsBase

"""
Gillespie algorithm for ABC SMC simulations

Input:
	- params : parameters values
	- init_values : inital values for the species
	- stoch_matrix : stochiometry matrix for the model
	- hazards : hazard functions of the model
	- times: time points for which simulated data is generated

Output:
	simulated data for all requested time points (in an array where each item is an array with data for all species)
"""
function Gillespie_abc(params::Array{Float64,1},init_values::Array{Float64,1},stoch_matrix::Array{Array{Int64,1},1},hazards::Array{Function,1},times::FloatRange{Float64})

	times = collect(times)
	time = times[1]
	time_counter = 1
	flag = true
	concentrations = Array{Array{Float64,1},1}(length(times))

	current_values = init_values

	time_counter_runs = 0

	#while simulation time hasn't exceeded the last time point
	while flag

		#update the next values and time point
		time_counter_runs +=1

		(next_values,next_time,total_hazard) = Gillespie_abc_one_step(params,current_values,stoch_matrix,hazards,time)
		concentrations[time_counter] = next_values
		time = next_time
		current_values = next_values

		#println("time ", time)
		#println("values ", current_values)

		if sum(current_values) > (1e4*length(current_values))
			#need to end simulation
			time = times[end]
		end

		if (total_hazard <= 0.0)
			#need to end simulation
			time = times[end]
		end

		#don't allow negative concentrations
		if (sum(current_values.<0)>0)
			#one of the species has a negative concentration
			time = times[end]
		end

		#as long as the next time point is larger than the next wanted time point, need to save the same concentration for all the time points in between
		while (time >=times[time_counter])
			time_counter +=1
			time_counter_runs=0
			concentrations[time_counter] = next_values
			if time_counter >= length(times)
				flag = false
				break
			end
		end

	end

	return concentrations
end

"""
Gillespie algorithm for ABC SMC simulations - one step only

Input:
	- params : parameters values
	- values : current values for the species
	- stoch_matrix : stochiometry matrix for the model
	- hazards : hazard functions of the model
	- time : current time point

Output:
	- values : the next values for all species
	- t_next : next time point
"""
function Gillespie_abc_one_step(params::Array{Float64,1},values::Array{Float64,1},stoch_matrix::Array{Array{Int64,1},1},hazards::Array{Function,1},time::Float64)

	#get the hazard function for each reaction
	h_vals = zeros(length(hazards))
	for i in 1:length(hazards)
		#run the hazard function with its corresponding rate
		h_vals[i] = hazards[i](params,values)
	end

	#get the combined hazard function
	h0 = sum(h_vals)

	if h0<=0.0
		return (values,time,h0)
	end

	#get the next time step
	#t_inc = pdf(Exponential(),h0)
	t_inc = (-1.0/h0) * log(rand())
	t_next = time + t_inc

	#get the probability for each reaction
	h_probs = h_vals ./ h0

	#randomly choose a hazard function according to their probabilities
	chosen_hazard = StatsBase.sample(hazards,WeightVec(h_probs))
	chosen_hazard_index = findfirst(hazards,chosen_hazard)

	#update values according to the chosen hazard
	#values is a one dimentional array of the current concentrations
	#we add to it the values of the array in the row that has the same index as the chosen hazard

	values = values .+ stoch_matrix[chosen_hazard_index]

	return (values,t_next,h0)

end


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

using Markdown
using InteractiveUtils

# This Pluto notebook uses @bind for interactivity. When running this notebook outside of Pluto, the following 'mock version' of @bind gives bound variables a default value (instead of an error).
macro bind(def, element)
    quote
        local el = $(esc(element))
        global $(esc(def)) = Core.applicable(Base.get, el) ? Base.get(el) : missing
        el
    end
end

# ╔═╡ d2ff2918-657b-452d-84d7-c26aee730c8c
using PlutoUI, Plots, LaTeXStrings, Latexify

# ╔═╡ df54c070-9a40-11eb-3ff7-c782f6937b66
md"""
### Polinômios de Bernstein

Um exemplo usando os polinômios de bernstein

"""

# ╔═╡ f2819f6f-80fa-4bfc-a911-cc45160079cc
begin
	B0(t) =(1-t)^2
	B1(t) = 2(1-t)t
	B2(t) = t^2
end

# ╔═╡ 0bc5ad31-900f-4116-b1a2-e9f46b4e3615
A=[0,0];

# ╔═╡ 99460aa9-d7e5-4a3a-8f88-6c6298b8277d
@bind x1 Slider(0:5; default = 0, show_value=true )

# ╔═╡ 1d3b3ce0-b629-4307-a400-5efed3b529aa
@bind y1 Slider(0:5; default = 1, show_value=true )

# ╔═╡ fb369ea9-90ea-43b2-9d13-b52f72f98557
B=[x1,y1]

# ╔═╡ 13b439a5-8cf6-415a-af86-0fefaad292d9
md" Agora o ponto de controle"

# ╔═╡ 7ddb0fa7-19c7-47f8-af6e-eeac257f2051
md"""
Xc= $(@bind xc Slider(0:0.2:5; default = 1, show_value=true ))

Yc= $(@bind yc Slider(0:0.2:5; default = 1, show_value=true ))
"""

# ╔═╡ 66a4c751-7fc6-4c1e-8384-d5621e8d99b2
C =[xc,yc]

# ╔═╡ dc41d708-5404-4059-9f18-c5caebeec10a
Γ1(t)= B0(t)*A[1] + B1(t)*C[1] +B2(t)*B[1]

# ╔═╡ cc2b7839-06a3-47a7-98a3-b1f7ed2f6911
Γ2(t) = B0(t)*A[2] + B1(t)*C[2] +B2(t)*B[2]

# ╔═╡ 274c84f7-36a8-4269-b269-8d93d04da8cb
t=range(0,1,length=101)

# ╔═╡ 9e00cca8-1e8b-4ce6-9894-d35377581f4e
x=[Γ1(s) for s in t]

# ╔═╡ 98af52fe-46c3-462c-90f7-d1d0ce875fad
y=[Γ2(s) for s in t]

# ╔═╡ 613e494f-5849-4904-bf5d-0802e55f93f4
begin
	plot(x,y, aspect_ratio=:equal, label=false)
	scatter!((xc,yc), label=false)
end

# ╔═╡ Cell order:
# ╟─df54c070-9a40-11eb-3ff7-c782f6937b66
# ╟─d2ff2918-657b-452d-84d7-c26aee730c8c
# ╟─f2819f6f-80fa-4bfc-a911-cc45160079cc
# ╠═0bc5ad31-900f-4116-b1a2-e9f46b4e3615
# ╟─99460aa9-d7e5-4a3a-8f88-6c6298b8277d
# ╟─1d3b3ce0-b629-4307-a400-5efed3b529aa
# ╠═fb369ea9-90ea-43b2-9d13-b52f72f98557
# ╟─13b439a5-8cf6-415a-af86-0fefaad292d9
# ╠═7ddb0fa7-19c7-47f8-af6e-eeac257f2051
# ╠═66a4c751-7fc6-4c1e-8384-d5621e8d99b2
# ╠═dc41d708-5404-4059-9f18-c5caebeec10a
# ╠═cc2b7839-06a3-47a7-98a3-b1f7ed2f6911
# ╠═274c84f7-36a8-4269-b269-8d93d04da8cb
# ╠═9e00cca8-1e8b-4ce6-9894-d35377581f4e
# ╠═98af52fe-46c3-462c-90f7-d1d0ce875fad
# ╠═613e494f-5849-4904-bf5d-0802e55f93f4


text:
@testset "GammaLikelihood" begin
    for args in ((), (1.0,), (exp,), (ExpLink(),), (1.0, exp), (1.0, ExpLink()))
        lik = GammaLikelihood(args...)
        @test lik isa GammaLikelihood{Float64,ExpLink}
        @test lik.α == 1
    end

    lik = GammaLikelihood(1.0)
    test_interface(lik, SqExponentialKernel(), rand(10); functor_args=(:α, :invlink))
end


text:
mutable struct Sprite
    x::Int
    y::Int

    width::Int
    height::Int

    offsetX::Int
    offsetY::Int
    realWidth::Int
    realHeight::Int

    surface::Cairo.CairoSurface
    filename::String

    function Sprite(x, y, width, height, offsetX, offsetY, realWidth, realHeight, surface, filename)
        return finalizer(s -> deleteSurface(s.surface), new(x, y, width, height, offsetX, offsetY, realWidth, realHeight, surface, filename))
    end
end

mutable struct SpriteHolder
    sprite::Sprite

    readtime::Float64
    path::String
end

function prepareSurfaceForQuad(atlasSurface::Cairo.CairoSurface, x, y, width, height, offsetX, offsetY, realWidth, realHeight)
    surface = CairoARGBSurface(width, height)
    ctx = getSurfaceContext(surface)

    drawImage(ctx, atlasSurface, 0, 0, x, y, width, height)

    return surface
end

function loadSprites(metaFn::String, spritesFn::String)
    splitAtlas = get(config, "split_atlas_into_smaller_surfaces", true)

    surface = open(Cairo.read_from_png, spritesFn)
    buffer = IOBuffer(open(read, metaFn))

    res = Dict{String, SpriteHolder}()

    # Mostly useless information
    read(buffer, Int32)
    Maple.readString(buffer)
    read(buffer, Int32)

    count = read(buffer, Int16)
    for i in 1:count
        dataFile = Maple.readString(buffer)
        sprites = read(buffer, Int16)

        for j in 1:sprites
            path = replace(Maple.readString(buffer), "\\" => "/")

            x = read(buffer, Int16)
            y = read(buffer, Int16)

            width = read(buffer, Int16)
            height = read(buffer, Int16)

            offsetX = read(buffer, Int16)
            offsetY = read(buffer, Int16)
            realWidth = read(buffer, Int16)
            realHeight = read(buffer, Int16)

            quadSurface = splitAtlas ? prepareSurfaceForQuad(surface, x, y, width, height, offsetX, offsetY, realWidth, realHeight) : surface

            sprite = Sprite(
                splitAtlas ? 0 : x,
                splitAtlas ? 0 : y,

                width,
                height,

                offsetX,
                offsetY,

                realWidth,
                realHeight,

                quadSurface,
                spritesFn
            )

            res[path] = SpriteHolder(
                sprite,

                time(),
                spritesFn
            )
        end
    end

    if splitAtlas
        deleteSurface(surface)
    end

    return res
end

function findTextureAnimations(texture::String, sprites::Dict{String, SpriteHolder}; maxPadding=7)::Array{String, 1}
    textures = String[]

    for i in 1:maxPadding
        testFrame = texture * "0"^i

        if haskey(sprites, testFrame) && sprites[testFrame].sprite.width != 0 && sprites[testFrame].sprite.height != 0
            for j in 0:10^i - 1
                frame = texture * lpad(j, i, "0")

                if haskey(sprites, frame) && sprites[frame].sprite.width != 0 && sprites[frame].sprite.height != 0
                    push!(textures, frame)

                else
                    break
                end
            end
        end
    end

    return textures
end

max_stars_repo_name: CelestialCartographers/Ahorn

max_issues_repo_name: CelestialCartographers/Ahorn

max_forks_repo_name: CelestialCartographers/Ahorn


text:
using Pkg
Pkg.activate(".")
Pkg.instantiate()

using Plots
using StatsPlots
using RDatasets
using LinearAlgebra
using Query
using Statistics

iris = dataset("datasets", "iris")

@df iris scatter(
    :PetalLength,
    :PetalWidth;
    label="",
    xlabel = "Petal length",
    ylabel = "Petal width"
)

max_stars_repo_name: rzehumat/JuliaCourse

max_issues_repo_name: rzehumat/JuliaCourse

max_forks_repo_name: rzehumat/JuliaCourse


text:
# This file is a part of JuliaFEM.
# License is MIT: see https://github.com/JuliaFEM/FEMBasis.jl/blob/master/LICENSE

function get_reference_element_coordinates end
function eval_basis! end
function eval_dbasis! end

function calculate_interpolation_polynomials(p, V)
    basis = []
    first(p.args) == :+ || error("Use only summation between terms of polynomial")
    args = p.args[2:end]
    n = size(V, 1)
    b = zeros(n)
    for i in 1:n
        fill!(b, 0.0)
        b[i] = 1.0
        # TODO: Think about numerical stability with
        # inverting Vandermonde matrix?
        solution = V \ b
        N = Expr(:call, :+)
        for (ai, bi) in zip(solution, args)
            isapprox(ai, 0.0) && continue
            push!(N.args, Calculus.simplify( :( $ai * $bi ) ))
        end
        push!(basis, N)
    end
    return basis
end

function calculate_interpolation_polynomial_derivatives(basis, D)
    vars = [:u, :v, :w]
    dbasis = Matrix(undef, D, length(basis))
    for (i, N) in enumerate(basis)
        partial_derivatives = []
        for j in 1:D
            dbasis[j, i] = Calculus.simplify(Calculus.differentiate(N, vars[j]))
        end
    end
    return dbasis
end

function create_basis(name, description, X::Vector{<:Vecish{D}}, p::Expr) where D
    @debug "create basis given antsatz polynomial" name description X p
    V = vandermonde_matrix(p, X)
    basis = calculate_interpolation_polynomials(p, V)
    return create_basis(name, description, X, basis)
end

function create_basis(name, description, X::Vector{<:Vecish{D}}, basis::Vector) where D
    @assert length(X) == length(basis)
    @debug "create basis given basis functions" name description X basis
    dbasis = calculate_interpolation_polynomial_derivatives(basis, D)
    return create_basis(name, description, Vec.(X), basis, dbasis)
end

function create_basis(name, description, X::Vector{<:Vecish{D, T}}, basis, dbasis) where {D, T}
    N = length(X)
    @debug "create basis given basis functions and derivatives" name description X basis dbasis

    Q = Expr(:block)
    for i=1:N
        push!(Q.args, :(N[$i] = $(basis[i])))
    end

    V = Expr(:block)
    for i=1:N
        push!(V.args, :(dN[$i] = Vec(float.(tuple($(dbasis[:, i]...))))))
    end

    if D == 1
        unpack = :((u,) = xi)
    elseif D == 2
        unpack = :((u, v) = xi)
    else
        unpack = :((u, v, w) = xi)
    end

    code = quote
        struct $name <: FEMBasis.AbstractBasis{$D}
        end

        Base.@pure function Base.size(::Type{$name})
            return ($D, $N)
        end

        function Base.size(::Type{$name}, j::Int)
            j == 1 && return $D
            j == 2 && return $N
        end

        Base.@pure function Base.length(::Type{$name})
            return $N
        end

        function FEMBasis.get_reference_element_coordinates(::Type{$name})
            return $X
        end

        @inline function FEMBasis.eval_basis!(::Type{$name}, N::Vector{<:Number}, xi::Vec)
            @assert length(N) == $N
            $unpack
            @inbounds $Q
            return N
        end

        @inline function FEMBasis.eval_dbasis!(::Type{$name}, dN::Vector{<:Vec{$D}}, xi::Vec)
            @assert length(dN) == $N
            $unpack
            @inbounds $V
            return dN
        end
    end
    return code
end

create_basis_and_eval(args...) = eval(create_basis(args...))


text:
using Zygote
using Optim
using DiffImages
using NLsolve

"""
    fit_overpotential(model, k; kwargs...)

Given values for current/rate constant and specified model parameters, find the overpotentials that must have resulted in it.
"""
function fit_overpotential(model::KineticModel, k; kT=.026, kwargs...)
    function compare_k!(storage, V)
        storage .= compute_k(V, model; kT=kT, kwargs...) .- k
    end

    # this isn't quite right yet...
    # function grad!(storage, V)
    #     gs = gradient(V -> compute_k(V, model; kT=kT, kwargs...) .- k, V)
    #     for i in 1:length(V)
    #         storage[i] = gs[i]
    #     end
    # end
    # Vs = nlsolve(compare_k!, grad!, repeat([0.1], length(k)))
    Vs = nlsolve(compare_k!, repeat([0.1], length(k)))
    if !converged(Vs)
        @warn "Overpotential fit not fully converged...you may have fed in an unreachable reaction rate!"
    end
    Vs.zero
end

fitting_params(t::Type{<:KineticModel}) = fieldnames(t)
fitting_params(::Type{MarcusHushChidsey}) = (:A, :λ)
fitting_params(::Type{MarcusHushChidseyDOS}) = (:A, :λ)
const default_param_bounds = Dict(:A => (0.1, 50000), :λ => (0.01, 0.5), :α => (0.01, 0.99))

"""
    fit_model(exp_data, model_type; kwargs...)

# Arguments
* `exp_data::Matrix`: two columns, first with voltage values, second with current
* `model_type::Type{<:KineticModel}`

# Keyword Arguments
Requirements differ by model type...
* `ButlerVolmer`, `AsymptoticMarcusHushChidsey`, `Marcus`: none
* `MarcusHushChidsey`: average_dos::Float64 OR dos::DOSData OR dos_file::String
* `MarcusHushChidseyDOS`: dos::DOSData OR dos_file
Some are always options...
* `param_bounds::Dict{Symbol,Any}`: ranges of guesses for relevant model parameters. (must include all necessary keys, but defaults to some sensible ranges if not provided, see `default_param_bounds`...note that you should provide this for faster fitting if you know bounds)
* E_min and E_max for integral models...defaults to +/- 100kT or in case of MarcusHushChidseyDOS, to energy bounds on DOS data
"""
function fit_model(
    exp_data::Matrix,
    model_type::Type{<:KineticModel};
    param_bounds::Dict = default_param_bounds,
    kT::Real = 0.026,
    kwargs...
)
    V_vals = exp_data[:, 1]
    eval_model(model) = [compute_k(V, model; kT = kT) for V in V_vals]
    _fit_model(exp_data, model_type, param_bounds, eval_model)
end

# TODO: check if this works with Cq
function fit_model(
    exp_data::Matrix,
    model_type::Type{<:IntegralModel};
    param_bounds::Dict = default_param_bounds,
    kT::Real = 0.026,
    E_min = -100 * kT,
    E_max = 100 * kT,
    kwargs...,
)
    V_vals = exp_data[:, 1]
    eval_model(model) = compute_k(V_vals,
                                  model;
                                  kT = kT,
                                  E_min = E_min,
                                  E_max = E_max,
                                  kwargs...)
    _fit_model(exp_data, model_type, param_bounds, eval_model; kwargs...)
end

function _fit_model(
    exp_data,
    model_type::Type{<:KineticModel},
    param_bounds,
    model_evaluator;
    kwargs...,
)
    I_vals = exp_data[:, 2]

    model_builder = _get_model_builder(model_type, param_bounds; kwargs...)
    sq_error(I_pred) = sum((I_pred .- I_vals) .^ 2)

    # WTF, why do I need these here again...just when I think I understand scope
    fitting_params(t::Type{<:KineticModel}) = fieldnames(t)
    fitting_params(::Type{MarcusHushChidsey}) = (:A, :λ)
    fitting_params(::Type{MarcusHushChidseyDOS}) = (:A, :λ)

    # find best-fitting params
    # Zygote is tripped up by QuadGK, so that one has to be done with black-box optimization, but the non-integral models work with autodiff
    opt_func = params -> sq_error(model_evaluator(model_builder(params)))
    local best_params
    function grad!(s, x)
        gs = gradient(params -> opt_func(params), x)[1]
        for i in 1:length(x)
            s[i] = gs[i]
        end
    end
    lower = Float64.([param_bounds[p][1] for p in fitting_params(model_type)])
    upper = Float64.([param_bounds[p][2] for p in fitting_params(model_type)])
    init_guess = 0.5 .* (lower .+ upper)

    # set optimisers based on the model type
    # IntegralModels performed better with GradientDescent
    optimizer, opts = if model_type <: IntegralModel
      opts = Optim.Options(show_trace = false,
                           iterations = 3,
                           outer_iterations = 4,
                           x_tol = 1.,
                           f_tol = 1e3)
      Fminbox(GradientDescent()), opts
    else
      Fminbox(LBFGS()), Optim.Options()
    end
    opt = optimize(opt_func, grad!, lower, upper, init_guess, optimizer, opts)
    best_params = Optim.minimizer(opt)

    # construct and return model
    model_builder(best_params)
end

"""
    is_dosmodel(type{<:KineticModel})

Returns true if the model needs DOS information.
"""
is_dosmodel(t::Type{<:KineticModel}) = false
is_dosmodel(t::Type{<:IntegralModel}) = t == MarcusHushChidsey || t == MarcusHushChidseyDOS

function _get_model_builder(model_type, param_bounds; kwargs...)
    # check that all necessary params are provided
    @assert Set(keys(param_bounds)) >= Set(fitting_params(model_type))

    # check for kwargs, build DOS object if needed, etc.
    arg_names = keys(kwargs)
    local model_builder
    if is_dosmodel(model_type)
        local dos_arg
        if :dos in arg_names
            dos_arg = :dos
        elseif :dos_file in arg_names
            dos_arg = :dos_file
        elseif :average_dos in arg_names
            @assert model_type == MarcusHushChidsey "You haven't provided a DOSData object or a file from which to build it"
            dos_arg = :average_dos
        else
            throw(ArgumentError("I need DOS information for this model type!"))
        end
        model_builder = param_tuple -> model_type(param_tuple..., kwargs[dos_arg])
    else
        model_builder = param_tuple -> model_type(param_tuple...)
    end
    return model_builder
end


text:
function abstract_call_gf_by_type(interp::CustomInterpreter, @nospecialize(f), argtypes::Vector{Any}, @nospecialize(atype), sv::InferenceState,
                                  max_methods::Int = InferenceParams(interp).MAX_METHODS)
    ret = @invoke abstract_call_gf_by_type(interp::AbstractInterpreter, f, argtypes::Vector{Any}, atype, sv::InferenceState,
                                           max_methods::Int)

    return ret
end

function abstract_call_method_with_const_args(interp::CustomInterpreter, @nospecialize(rettype), @nospecialize(f), argtypes::Vector{Any}, match::MethodMatch, sv::InferenceState, edgecycle::Bool)
    ret = @invoke abstract_call_method_with_const_args(interp::AbstractInterpreter, rettype, f, argtypes::Vector{Any}, match::MethodMatch, sv::InferenceState, edgecycle::Bool)

    return ret
end

function abstract_call_method(interp::CustomInterpreter, method::Method, @nospecialize(sig), sparams::SimpleVector, hardlimit::Bool, sv::InferenceState)
    ret = @invoke abstract_call_method(interp::AbstractInterpreter, method::Method, sig, sparams::SimpleVector, hardlimit::Bool, sv::InferenceState)

    return ret
end

function precise_container_type(interp::CustomInterpreter, @nospecialize(itft), @nospecialize(typ), sv::InferenceState)
    ret = @invoke precise_container_type(interp::AbstractInterpreter, itft, typ, sv::InferenceState)

    return ret
end

function abstract_iteration(interp::CustomInterpreter, @nospecialize(itft), @nospecialize(itertype), sv::InferenceState)
    ret = @invoke abstract_iteration(interp::AbstractInterpreter, itft, itertype, sv::InferenceState)

    return ret
end

function abstract_apply(interp::CustomInterpreter, @nospecialize(itft), @nospecialize(aft),    aargtypes::Vector{Any}, sv::InferenceState,
                        max_methods::Int = InferenceParams(interp).MAX_METHODS)
    ret = @invoke abstract_apply(interp::AbstractInterpreter, itft, aft, aargtypes::Vector{Any}, sv::InferenceState,
                                 max_methods::Int)

    return ret
end

function abstract_call_builtin(interp::CustomInterpreter, f::Builtin, fargs::Union{Nothing,Vector{Any}},
                               argtypes::Vector{Any}, sv::InferenceState, max_methods::Int)
    ret = @invoke abstract_call_builtin(interp::AbstractInterpreter, f::Builtin, fargs::Union{Nothing,Vector{Any}},
                                        argtypes::Vector{Any}, sv::InferenceState, max_methods::Int)

    return ret
end

function abstract_call_known(interp::CustomInterpreter, @nospecialize(f),
                             fargs::Union{Nothing,Vector{Any}}, argtypes::Vector{Any},
                             sv::InferenceState,
                             max_methods::Int = InferenceParams(interp).MAX_METHODS)
    ret = @invoke abstract_call_known(interp::AbstractInterpreter, f,
                                      fargs::Union{Nothing,Vector{Any}}, argtypes::Vector{Any},
                                      sv::InferenceState,
                                      max_methods::Int)
end

function abstract_call(interp::CustomInterpreter, fargs::Union{Nothing,Vector{Any}}, argtypes::Vector{Any},
                       sv::InferenceState, max_methods::Int = InferenceParams(interp).MAX_METHODS)
    ret = @invoke abstract_call(interp::AbstractInterpreter, fargs::Union{Nothing,Vector{Any}}, argtypes::Vector{Any},
                                sv::InferenceState, max_methods::Int)

    return ret
end

function abstract_eval_cfunction(interp::CustomInterpreter, e::Expr, vtypes::VarTable, sv::InferenceState)
    ret = @invoke abstract_eval_cfunction(interp::AbstractInterpreter, e::Expr, vtypes::VarTable, sv::InferenceState)

    return ret
end

function abstract_eval_value_expr(interp::CustomInterpreter, e::Expr, vtypes::VarTable, sv::InferenceState)
    ret = @invoke abstract_eval_value_expr(interp::AbstractInterpreter, e::Expr, vtypes::VarTable, sv::InferenceState)

    return ret
end

function abstract_eval_special_value(interp::CustomInterpreter, @nospecialize(e), vtypes::VarTable, sv::InferenceState)
    ret = @invoke abstract_eval_special_value(interp::AbstractInterpreter, e, vtypes::VarTable, sv::InferenceState)

    return ret
end

function abstract_eval_value(interp::CustomInterpreter, @nospecialize(e), vtypes::VarTable, sv::InferenceState)
    ret = @invoke abstract_eval_value(interp::AbstractInterpreter, e, vtypes::VarTable, sv::InferenceState)

    return ret
end

function abstract_eval_statement(interp::CustomInterpreter, @nospecialize(e), vtypes::VarTable, sv::InferenceState)
    ret = @invoke abstract_eval_statement(interp::AbstractInterpreter, e, vtypes::VarTable, sv::InferenceState)

    return ret
end

function typeinf_local(interp::CustomInterpreter, frame::InferenceState)
    ret = @invoke typeinf_local(interp::AbstractInterpreter, frame::InferenceState)

    return ret
end

function typeinf_nocycle(interp::CustomInterpreter, frame::InferenceState)
    ret = @invoke typeinf_nocycle(interp::AbstractInterpreter, frame::InferenceState)

    return ret
end

max_stars_repo_name: aviatesk/AbstractInterpreterPlayground

max_issues_repo_name: aviatesk/AbstractInterpreterPlayground

max_forks_repo_name: aviatesk/AbstractInterpreterPlayground


text:
using MetaPkg
using Base.Test

MetaPkg.dry_run()

# just check that they run without error
spec = MetaPkg.load_meta("MetaGL")
spec = MetaPkg.load_meta("MetaGL", Pkg.dir("MetaPkg","requires"))
meta_add("MetaGL")
meta_add("MetaGL", dir = Pkg.dir("MetaPkg","requires"))
meta_rm("MetaGL")
meta_free("MetaGL")
meta_checkout("MetaGL")
meta_checkout("MetaGL", "dev")


text:
# Diag is treated in a special way.

# f(X)=diagm(X)
import Base.diagm
diagm(A::ADnode)=ADdiag(A)
export diagm

#TODO define diagonalm to construct explicitly a node

#ftranspose(A::ADnode)=ADnode(Ftranspose,A)
#export ftranpose

#Ftranspose(X)=(X',nothing)

#function Ftranspose_inplace(value::Array,auxvalue,X::Array)
#    copy!(value,X')
#end

#function Dtranspose(derivativeIDX,f_c,faux_c,grad_c,grad_n,X::Array)
#    axpy!(1.0,grad_c',grad_n)
#end

#if PROC=="GPU"
#    function Ftranspose_inplace(value::CudaArray,auxvalue,X::CudaArray)
#        tmp=CudaArray(Float64,size(X))
#        CUBLAS.geam!('T','N',1.0,X,0.0,tmp,value) # value=X'
#        free(tmp)
#    end
#
#    function Dtranspose(derivativeIDX,f_c,faux_c,grad_c,grad_n,X::CudaArray)
#        tmp=CudaArray(Float64,size(grad_n))
#        CUBLAS.geam!('T','N',1.0,grad_c,0.0,tmp,tmp) # tmp=grad_c'
#        axpy!(1.0,tmp,grad_n); free(tmp)
#    end
#
#end

#Derivative[Ftranspose]=Dtranspose
#Inplace[Ftranspose]=Ftranspose_inplace
#
#export Ftranspose

import Base.*

function *(A::ADdiag,B::ADnode)
    node[A.index]=nothing
    return ADnode(FdiagAmultX,[node[A.parent] B])
end

function *(A::ADdiag,B::Real)
    node[A.index]=nothing
    return diagm(node[A.parent]*ADconst(B))
end

function *(A::ADnode,B::ADdiag)
    node[B.index]=nothing
    return ADnode(FAmultdiagX,[A node[B.parent]])
end

function *(A::Real,B::ADdiag)
    node[B.index]=nothing
    return diagm(ADconst(A)*node[B.parent])
end

function *(A::ADdiag,B::ADdiag)
    node[A.index]=nothing
    node[B.index]=nothing
    return diagm(node[A.parent].*node[B.parent])
end

export *

#import Base.+#

#function +(A::ADtrans,B::ADnode)
#    node[A.index]=nothing
#    return ADnode(FXtransposePY,[node[A.parent] B])
#end

#function +(A::ADnode,B::ADtrans)
#    node[B.index]=nothing
#    return ADnode(FXPYtranspose,[A node[B.parent]])
#end

#function +(A::ADtrans,B::ADtrans)
#    node[A.index]=nothing
#    node[B.index]=nothing
#    return ADnode(FXtransposePYtranspose,[node[A.parent] node[B.parent]])
#end

#export +

#import Base.-

#function -(A::ADtrans,B::ADnode)
#    node[A.index]=nothing
#    return ADnode(FXtransposeMY,[node[A.parent] B])
#end

#function -(A::ADnode,B::ADtrans)
#    node[B.index]=nothing
#    return ADnode(FXMYtranspose,[A node[B.parent]])
#end

#function -(A::ADtrans,B::ADtrans)
#    node[A.index]=nothing
#    node[B.index]=nothing
#    return ADnode(FXtransposeMYtranspose,[node[A.parent] node[B.parent]])
#end

#trans=ftranspose
#export trans

#export -

max_stars_repo_name: davidbarber/AutoDiff

max_issues_repo_name: davidbarber/AutoDiff

max_forks_repo_name: davidbarber/AutoDiff