text:
@enum ConnectionState Busy Idle Dead

"""
    HTTP.Connection

Represents a persistent client connection to a remote host; only created
when a server response includes the "Connection: keep-alive" header. An open and non-idle connection
will be reused when sending subsequent requests to the same host.
"""
mutable struct Connection{I <: IO}
    id::Int
    socket::I
    state::ConnectionState
end

Connection(tcp::IO) = Connection(0, tcp, Busy)
Connection(id::Int, tcp::IO) = Connection(id, tcp, Busy)
busy!(conn::Connection) = (conn.state == Dead || (conn.state = Busy); return nothing)
idle!(conn::Connection) = (conn.state == Dead || (conn.state = Idle); return nothing)
dead!(conn::Connection) = (conn.state == Dead || (conn.state = Dead; close(conn.socket)); return nothing)

"""
    HTTP.Client([logger::IO]; args...)

A type to facilitate connections to remote hosts, send HTTP requests, and manage state between requests.
Takes an optional `logger` IO argument where client activity is recorded (defaults to `STDOUT`).
Additional keyword arguments can be passed that will get transmitted with each HTTP request:

* `chunksize::Int`: if a request body is larger than `chunksize`, the "chunked-transfer" http mechanism will be used and chunks will be sent no larger than `chunksize`
* `connecttimeout::Float64`: sets a timeout on how long to wait when trying to connect to a remote host; default = 10.0 seconds
* `readtimeout::Float64`: sets a timeout on how long to wait when receiving a response from a remote host; default = 9.0 seconds
* `tlsconfig::TLS.SSLConfig`: a valid `TLS.SSLConfig` which will be used to initialize every https connection
* `maxredirects::Int`: the maximum number of redirects that will automatically be followed for an http request
* `allowredirects::Bool`: whether redirects should be allowed to be followed at all; default = `true`
* `forwardheaders::Bool`: whether user-provided headers should be forwarded on redirects; default = `false`
* `retries::Int`: # of times a request will be tried before throwing an error; default = 3
* `managecookies::Bool`: whether the request client should automatically store and add cookies from/to requests (following appropriate host-specific & expiration rules)
* `statusraise::Bool`: whether an `HTTP.StatusError` should be raised on a non-2XX response status code
* `insecure::Bool`: whether an "https" connection should allow insecure connections (no TLS verification); default = `false`
"""
mutable struct Client
    # connection pools for keep-alive; key is host
    httppool::Dict{String, Vector{Connection{TCPSocket}}}
    httpspool::Dict{String, Vector{Connection{TLS.SSLContext}}}
    # cookies are stored in-memory per host and automatically sent when appropriate
    cookies::Dict{String, Set{Cookie}}
    # buffer::Vector{UInt8} #TODO: create a fixed size buffer for reading bytes off the wire and having http_parser use, this should keep allocations down, need to make sure MbedTLS supports blocking readbytes!
    parser::Parser
    logger::Option{IO}
    # global request settings
    options::RequestOptions
    connectioncount::Int
end

Client(logger::Option{IO}, options::RequestOptions) = Client(Dict{String, Vector{Connection{TCPSocket}}}(),
                                                     Dict{String, Vector{Connection{TLS.SSLContext}}}(),
                                                     Dict{String, Set{Cookie}}(),
                                                     Parser(), logger, options, 1)

const DEFAULT_CHUNK_SIZE = 2^20
const DEFAULT_OPTIONS = :((DEFAULT_CHUNK_SIZE, true, 15.0, 15.0, nothing, 5, true, false, 3, true, true, false))

@eval begin
    Client(logger::Option{IO}; args...) = Client(logger, RequestOptions($(DEFAULT_OPTIONS)...; args...))
    Client(; args...) = Client(nothing, RequestOptions($(DEFAULT_OPTIONS)...; args...))
end

function setclient!(client::Client)
    global const DEFAULT_CLIENT = client
end

Base.haskey(::Type{http}, client, host) = haskey(client.httppool, host)
Base.haskey(::Type{https}, client, host) = haskey(client.httpspool, host)

getconnections(::Type{http}, client, host) = client.httppool[host]
getconnections(::Type{https}, client, host) = client.httpspool[host]

setconnection!(::Type{http}, client, host, conn) = push!(get!(client.httppool, host, Connection[]), conn)
setconnection!(::Type{https}, client, host, conn) = push!(get!(client.httpspool, host, Connection[]), conn)

backtrace() = sprint(Base.show_backtrace, catch_backtrace())

abstract type HTTPError <: Exception end

function Base.show(io::IO, e::HTTPError)
    println(io, "$(typeof(e)):")
    println(io, "Exception: $(e.e)")
    print(io, e.msg)
end

struct ConnectError <: HTTPError
    e::Exception
    msg::String
end

struct SendError <: HTTPError
    e::Exception
    msg::String
end

struct ClosedError <: HTTPError
    e::Exception
    msg::String
end

struct ReadError <: HTTPError
    e::Exception
    msg::String
end

struct RedirectException <: Exception
    maxredirects::Int
end
function Base.show(io::IO, err::RedirectException)
    print(io, "RedirectException: more than $(err.maxredirects) redirects attempted")
end
struct StatusError <: Exception
    status::Int
    response::Response
end
function Base.show(io::IO, err::StatusError)
    print(io, "HTTP.StatusError: received a '$(err.status) - $(Base.get(STATUS_CODES, err.status, "Unknown Code"))' status in response")
end

initTLS!(::Type{http}, hostname, opts, socket) = socket

function initTLS!(::Type{https}, hostname, opts, socket)
    stream = TLS.SSLContext()
    TLS.setup!(stream, get(opts, :tlsconfig, TLS.SSLConfig(!opts.insecure::Bool))::TLS.SSLConfig)
    TLS.associate!(stream, socket)
    TLS.hostname!(stream, hostname)
    TLS.handshake!(stream)
    return stream
end

function stalebytes!(c::TCPSocket)
    !isopen(c) && return
    nb_available(c) > 0 && readavailable(c)
    return
end
stalebytes!(c::TLS.SSLContext) = stalebytes!(c.bio)

function connect(client, sch, hostname, port, opts, verbose)
    logger = client.logger
    if haskey(sch, client, hostname)
        @log "checking if any existing connections to '$hostname' are re-usable"
        conns = getconnections(sch, client, hostname)
        inds = Int[]
        i = 1
        while i <= length(conns)
            c = conns[i]
            # read off any stale bytes left over from a possible error in a previous request
            # this will also trigger any sockets that timed out to be set to closed
            stalebytes!(c.socket)
            if !isopen(c.socket) || c.state == Dead
                @log "found dead connection #$(c.id) to delete"
                dead!(c)
                push!(inds, i)
            elseif c.state == Idle
                @log "found re-usable connection #$(c.id)"
                busy!(c)
                try
                    deleteat!(conns, sort!(unique(inds)))
                end
                return c
            end
            i += 1
        end
        try
            deleteat!(conns, sort!(unique(inds)))
        end
    end
    # if no re-usable connection was found, make a new connection
    try
        # EH: throws DNSError, OutOfMemoryError, or SystemError; retry once, but otherwise, we can't do much
        ip = @retry Base.getaddrinfo(hostname)
        # EH: throws error, ArgumentError for out-of-range port, UVError; retry if UVError
        tcp = @retryif Base.UVError @timeout(opts.connecttimeout::Float64,
                                        Base.connect(ip, Base.parse(Int, port)), error("connect timeout"))
        socket = initTLS!(sch, hostname, opts, tcp)
        conn = Connection(client.connectioncount, socket)
        client.connectioncount += 1
        setconnection!(sch, client, hostname, conn)
        @log "created new connection #$(conn.id) to '$hostname'"
        return conn
    catch e
        throw(ConnectError(e, backtrace()))
    end
end

function addcookies!(client, host, req, verbose)
    logger = client.logger
    # check if cookies should be added to outgoing request based on host
    if haskey(client.cookies, host)
        cookies = client.cookies[host]
        tosend = Set{Cookie}()
        expired = Set{Cookie}()
        for (i, cookie) in enumerate(cookies)
            if Cookies.shouldsend(cookie, scheme(uri(req)) == "https", host, path(uri(req)))
                cookie.expires != DateTime() && cookie.expires < now(Dates.UTC) && (push!(expired, cookie); @log("deleting expired cookie: " * cookie.name); continue)
                push!(tosend, cookie)
            end
        end
        setdiff!(client.cookies[host], expired)
        if length(tosend) > 0
            @log "adding cached cookies for host to request header: " * join(map(x->x.name, tosend), ", ")
            req.headers["Cookie"] = string(Base.get(req.headers, "Cookie", ""), [c for c in tosend])
        end
    end
end

function connectandsend(client, ::Type{sch}, hostname, port, req, opts, verbose) where sch
    logger = client.logger
    conn = connect(client, sch, hostname, port, opts, verbose)
    opts.managecookies::Bool && addcookies!(client, hostname, req, verbose)
    try
        @log "sending request over the wire\n"
        reqstr = string(req, opts)
        verbose && (println(client.logger, "HTTP.Request:\n"); println(client.logger, reqstr))
        # EH: throws ArgumentError if socket is closed, UVError; retry if UVError,
        @retryif Base.UVError write(conn.socket, reqstr)
        !isopen(conn.socket) && throw(CLOSED_ERROR)
    catch e
        @log backtrace()
        typeof(e) <: ArgumentError && throw(ClosedError(e, backtrace()))
        throw(SendError(e, backtrace()))
    end
    return conn
end

function redirect(response, client, req, opts, stream, history, retry, verbose)
    logger = client.logger
    @log "checking for location to redirect"
    key = haskey(response.headers, "Location") ? "Location" :
          haskey(response.headers, "location") ? "location" : ""
    if key != ""
        push!(history, response)
        length(history) > opts.maxredirects::Int && throw(RedirectException(opts.maxredirects::Int))
        newuri = URIs.URL(response.headers[key])
        u = uri(req)
        newuri = !isempty(hostname(newuri)) ? newuri : URIs.URI(scheme=scheme(u), hostname=hostname(u), port=port(u), path=path(newuri), query=query(u))
        if opts.forwardheaders::Bool
            h = headers(req)
            delete!(h, "Host")
            delete!(h, "Cookie")
        else
            h = Headers()
        end
        redirectreq = Request(req.method, newuri, h, req.body)
        @log "redirecting to $(newuri)"
        return request(client, redirectreq, opts, stream, history, retry, verbose)
    end
end

const CLOSED_ERROR = ClosedError(ErrorException(""), "error receiving response; connection was closed prematurely")
function getbytes(socket, tm)
    try
        # EH: returns UInt8[] when socket is closed, error when socket is not readable, AssertionErrors, UVError;
        buffer = @retry @timeout(tm, readavailable(socket), error("read timeout"))
        return buffer, CLOSED_ERROR
    catch e
        return UInt8[], ReadError(e, backtrace())
    end
end

function processresponse!(client, conn, response, host, method, maintask, stream, tm, verbose)
    logger = client.logger
    while true
        buffer, err = getbytes(conn.socket, tm)
        if length(buffer) == 0 && !isopen(conn.socket)
            @log "socket closed before full response received"
            dead!(conn)
            close(response.body)
            # retry the entire request
            return false, err
        end
        @log "received bytes from the wire, processing"
        # EH: throws a couple of "shouldn't get here" errors; probably not much we can do
        errno, headerscomplete, messagecomplete, upgrade = HTTP.parse!(response, client.parser, buffer; host=host, method=method, maintask=maintask)
        if errno != HPE_OK
            dead!(conn)
            throw(ParsingError("error parsing response: $(ParsingErrorCodeMap[errno])\nCurrent response buffer contents: $(String(buffer))"))
        elseif messagecomplete
            http_should_keep_alive(client.parser, response) || (@log("closing connection (no keep-alive)"); dead!(conn))
            # idle! on a Dead will stay Dead
            idle!(conn)
            return true, StatusError(status(response), response)
        elseif stream && headerscomplete
            @log "processing the rest of response asynchronously"
            response.body.task = @async processresponse!(client, conn, response, host, method, maintask, false, tm, false)
            return true, nothing
        end
    end
    # shouldn't ever reach here
    dead!(conn)
    return false
end

function request(client::Client, req::Request, opts::RequestOptions, stream::Bool, history::Vector{Response}, retry::Int, verbose::Bool)
    retry = max(0, retry) # ensure non-negative
    update!(opts, client.options)
    verbose && not(client.logger) && (client.logger = STDOUT)
    logger = client.logger
    @log "using request options:\n\t" * join((s=>getfield(opts, s) for s in fieldnames(typeof(opts))), "\n\t")
    u = uri(req)
    host = hostname(u)
    sch = scheme(u) == "http" ? http : https
    @log "making $(method(req)) request for host: '$host' and resource: '$(resource(u))'"
    # maybe allow retrying for all kinds of errors?
    p = port(u)
    conn = @retryif ClosedError 4 connectandsend(client, sch, host, ifelse(p == "", "80", p), req, opts, verbose)

    response = Response(stream ? DEFAULT_CHUNK_SIZE : FIFOBuffers.DEFAULT_MAX, req)
    reset!(client.parser)
    success, err = processresponse!(client, conn, response, host, HTTP.method(req), current_task(), stream, opts.readtimeout::Float64, verbose)
    if !success
        retry >= opts.retries::Int && throw(err)
        return request(client, req, opts, stream, history, retry + 1, verbose)
    end
    @log "received response"
    opts.managecookies::Bool && !isempty(response.cookies) && (@log("caching received cookies for host: " * join(map(x->x.name, response.cookies), ", ")); union!(get!(client.cookies, host, Set{Cookie}()), response.cookies))
    response.history = history
    if opts.allowredirects::Bool && req.method != HEAD && (300 <= status(response) < 400)
        return redirect(response, client, req, opts, stream, history, retry, verbose)
    end
    if success && ((200 <= status(response) < 300) || !opts.statusraise::Bool)
        return response
    else
        retry >= opts.retries::Int && throw(err)
        return request(client, req, opts, stream, history, retry + 1, verbose)
    end
end

request(req::Request;
            opts::RequestOptions=RequestOptions(),
            stream::Bool=false,
            history::Vector{Response}=Response[],
            retry::Int=0,
            verbose::Bool=false,
            args...) =
    request(DEFAULT_CLIENT, req, RequestOptions(opts; args...), stream, history, retry, verbose)
request(client::Client, req::Request;
            opts::RequestOptions=RequestOptions(),
            stream::Bool=false,
            history::Vector{Response}=Response[],
            retry::Int=0,
            verbose::Bool=false,
            args...) =
    request(client, req, RequestOptions(opts; args...), stream, history, retry, verbose)

# build Request
function request(client::Client, method, uri::URI;
                 headers::Headers=Headers(),
                 body=FIFOBuffers.EMPTYBODY,
                 stream::Bool=false,
                 verbose::Bool=false,
                 args...)
    opts = RequestOptions(; args...)
    not(client.logger) && (client.logger = STDOUT)
    client.logger != STDOUT && (verbose = true)
    req = Request(method, uri, headers, body; options=opts, verbose=verbose, logger=client.logger)
    return request(client, req; opts=opts, stream=stream, verbose=verbose)
end
request(uri::String; verbose::Bool=false, query="", args...) = request(DEFAULT_CLIENT, GET, URIs.URL(uri; query=query); verbose=verbose, args...)
request(uri::URI; verbose::Bool=false, args...) = request(DEFAULT_CLIENT, GET, uri; verbose=verbose, args...)
request(method, uri::String; verbose::Bool=false, query="", args...) = request(DEFAULT_CLIENT, convert(HTTP.Method, method), URIs.URL(uri; query=query); verbose=verbose, args...)
request(method, uri::URI; verbose::Bool=false, args...) = request(DEFAULT_CLIENT, convert(HTTP.Method, method), uri; verbose=verbose, args...)

for f in [:get, :post, :put, :delete, :head,
          :trace, :options, :patch, :connect]
    f_str = uppercase(string(f))
    meth = convert(Method, f_str)
    @eval begin
        @doc """
    $($f)(uri; kwargs...) -> Response
    $($f)(client::HTTP.Client, uri; kwargs...) -> Response

Build and execute an http "$($f_str)" request. Query parameters can be passed via the `query` keyword argument as a `Dict`. Multiple
query parameters with the same key can be passed like `Dict("key1"=>["value1", "value2"], "key2"=>...)`.
Returns a `Response` object that includes the resulting status code (`HTTP.status(r)` and `HTTP.statustext(r)`),
response headers (`HTTP.headers(r)`), cookies (`HTTP.cookies(r)`), response history if redirects were involved
(`HTTP.history(r)`), and response body (`HTTP.body(r)` or `String(r)` or `take!(r)`).

The body or payload for a request can be given through the `body` keyword arugment.
The body can be given as a `String`, `Vector{UInt8}`, `IO`, `HTTP.FIFOBuffer` or `Dict` argument type.
See examples below for how to use an `HTTP.FIFOBuffer` for asynchronous streaming uploads.

If the body is provided as a `Dict`, the request body will be uploaded using the multipart/form-data encoding.
The key-value pairs in the Dict will constitute the name and value of each multipart boundary chunk.
Files and other large data arguments can be provided as values as IO arguments: either an `IOStream` such as returned via `open(file)`,
an `IOBuffer` for in-memory data, or even an `HTTP.FIFOBuffer`. For complete control over the multipart details, an
`HTTP.Multipart` type is provided to support setting the `Content-Type`, `filename`, and `Content-Transfer-Encoding` if desired. See `?HTTP.Multipart` for more details.

Additional keyword arguments supported, include:

* `headers::Dict{String,String}`: headers given as Dict to be sent with the request
* `body`: a request body can be given as a `String`, `Vector{UInt8}`, `IO`, `HTTP.FIFOBuffer` or `Dict`; see example below for how to utilize `HTTP.FIFOBuffer` for "streaming" request bodies; a `Dict` argument will be converted to a multipart form upload
* `stream::Bool=false`: enable response body streaming; depending on the response body size, the request will return before the full body has been received; as the response body is read, additional bytes will be recieved and put in the response body. Readers should read until `eof(response.body) == true`; see below for an example of response streaming
* `chunksize::Int`: if a request body is larger than `chunksize`, the "chunked-transfer" http mechanism will be used and chunks will be sent no larger than `chunksize`
* `connecttimeout::Float64`: sets a timeout on how long to wait when trying to connect to a remote host; default = 10.0 seconds
* `readtimeout::Float64`: sets a timeout on how long to wait when receiving a response from a remote host; default = 9.0 seconds
* `tlsconfig::TLS.SSLConfig`: a valid `TLS.SSLConfig` which will be used to initialize every https connection
* `maxredirects::Int`: the maximum number of redirects that will automatically be followed for an http request
* `allowredirects::Bool`: whether redirects should be allowed to be followed at all; default = `true`
* `forwardheaders::Bool`: whether user-provided headers should be forwarded on redirects; default = `false`
* `retries::Int`: # of times a request will be tried before throwing an error; default = 3
* `managecookies::Bool`: whether the request client should automatically store and add cookies from/to requests (following appropriate host-specific & expiration rules)
* `statusraise::Bool`: whether an `HTTP.StatusError` should be raised on a non-2XX response status code

Simple request example:
```julia
julia> resp = HTTP.get("http://httpbin.org/ip")
HTTP.Response:
\"\"\"
HTTP/1.1 200 OK
Connection: keep-alive
X-Powered-By: Flask
Content-Length: 32
Via: 1.1 vegur
Access-Control-Allow-Credentials: true
X-Processed-Time: 0.000903129577637
Date: Wed, 23 Aug 2017 23:35:59 GMT
Content-Type: application/json
Access-Control-Allow-Origin: *
Server: meinheld/0.6.1
Content-Length: 32

{
  "origin": "50.207.241.62"
}
\"\"\"

julia> String(resp)
"{\n  \"origin\": \"65.130.216.45\"\n}\n"
```

Response streaming example (asynchronous download):
```julia
julia> r = HTTP.get("http://httpbin.org/stream/100"; stream=true)
HTTP.Response:
\"\"\"
HTTP/1.1 200 OK
Connection: keep-alive
X-Powered-By: Flask
Transfer-Encoding: chunked
Via: 1.1 vegur
Access-Control-Allow-Credentials: true
X-Processed-Time: 0.000981092453003
Date: Wed, 23 Aug 2017 23:36:56 GMT
Content-Type: application/json
Access-Control-Allow-Origin: *
Server: meinheld/0.6.1

[HTTP.Response body of 27415 bytes]
Content-Length: 27390

{"id": 0, "origin": "50.207.241.62", "args": {}, "url": "http://httpbin.org/stream/100", "headers": {"Connection": "close", "User-Agent": "HTTP.jl/0.0.0", "Host": "httpbin.org", "Accept": "text/html,application/xhtml+xml,application/xml;q=0.9,*/*;q=0.8,application/json"}}
{"id": 1, "origin": "50.207.241.62", "args": {}, "url": "http://httpbin.org/stream/100", "headers": {"Connection": "close", "User-Agent": "HTTP.jl/0.0.0", "Host": "httpbin.org", "Accept": "text/html,application/xhtml+xml,application/xml;q=0.9,*/*;q=0.8,application/json"}}
{"id": 2, "origin": "50.207.241.62", "args": {}, "url": "http://httpbin.org/stream/100", "headers": {"Connection": "close", "User-Agent": "HTTP.jl/0.0.0", "Host": "httpbin.org", "
⋮
```
\"\"\"

julia> body = HTTP.body(r)
HTTP.FIFOBuffers.FIFOBuffer(27390, 1048576, 27390, 1, 27391, -1, 27390, UInt8[0x7b, 0x22, 0x69, 0x64, 0x22, 0x3a, 0x20, 0x30, 0x2c, 0x20  …  0x6e, 0x2f, 0x6a, 0x73, 0x6f, 0x6e, 0x22, 0x7d, 0x7d, 0x0a], Condition(Any[]), Task (done) @0x0000000112d84250, true)

julia> while true
           println(String(readavailable(body)))
           eof(body) && break
       end
{"id": 0, "origin": "50.207.241.62", "args": {}, "url": "http://httpbin.org/stream/100", "headers": {"Connection": "close", "User-Agent": "HTTP.jl/0.0.0", "Host": "httpbin.org", "Accept": "text/html,application/xhtml+xml,application/xml;q=0.9,*/*;q=0.8,application/json"}}
{"id": 1, "origin": "50.207.241.62", "args": {}, "url": "http://httpbin.org/stream/100", "headers": {"Connection": "close", "User-Agent": "HTTP.jl/0.0.0", "Host": "httpbin.org", "Accept": "text/html,application/xhtml+xml,application/xml;q=0.9,*/*;q=0.8,application/json"}}
{"id": 2, "origin": "50.207.241.62", "args": {}, "url": "http://httpbin.org/stream/100", "headers": {"Connection": "close", "User-Agent": "HTTP.jl/0.0.0", "Host": "httpbin.org", "Accept": "text/html,application/xhtml+xml,application/xml;q=0.9,*/*;q=0.8,application/json"}}
{"id": 3, "origin": "50.207.241.62", "args": {}, "url": "http://httpbin.org/stream/100", "headers": {"Connection": "close", "User-Agent": "HTTP.jl/0.0.0", "Host": "httpbin.org", "Accept": "text/html,application/xhtml+xml,application/xml;q=0.9,*/*;q=0.8,application/json"}}
...
```

Request streaming example (asynchronous upload):
```julia
# create a FIFOBuffer for sending our request body
f = HTTP.FIFOBuffer()
# write initial data
write(f, "hey")
# start an HTTP.post asynchronously
t = @async HTTP.post("http://httpbin.org/post"; body=f)
write(f, " there ") # as we write to f, it triggers another chunk to be sent in our async request
write(f, "sailor")
close(f) # setting eof on f causes the async request to send a final chunk and return the response

resp = wait(t) # get our response by getting the result of our asynchronous task
```
        """ function $(f) end
        ($f)(uri::AbstractString; verbose::Bool=false, query="", args...) = request(DEFAULT_CLIENT, $meth, URIs.URL(uri; query=query, isconnect=$(f_str == "CONNECT")); verbose=verbose, args...)
        ($f)(uri::URI; verbose::Bool=false, args...) = request(DEFAULT_CLIENT, $meth, uri; verbose=verbose, args...)
        ($f)(client::Client, uri::AbstractString; query="", args...) = request(client, $meth, URIs.URL(uri; query=query, isconnect=$(f_str == "CONNECT")); args...)
        ($f)(client::Client, uri::URI; args...) = request(client, $meth, uri; args...)
    end
end

function download(uri::AbstractString, file; threshold::Int=50000000, verbose::Bool=false, query="", args...)
    res = request(GET, uri; verbose=verbose, query=query, stream=true, args...)
    body = HTTP.body(res)
    file = Base.get(HTTP.headers(res), "Content-Encoding", "") == "gzip" ? string(file, ".gz") : file
    threshold_step = threshold
    nbytes = 0
    open(file, "w") do f
        while !eof(body)
            nbytes += write(f, readavailable(body))
            if verbose && nbytes > threshold
                println("[$(now())]: downloaded $nbytes bytes...")
                flush(STDOUT)
                threshold += threshold_step
            end
        end
        length(body) > 0 && write(f, readavailable(body))
    end
    return file
end


text:
function deim_get_Π(Ψ)
    N, k = size(Ψ)
    j = argmax(Ψ[:,1])
    Π = zeros(N, k)
    Π[j,1] = 1.0
    r = zeros(N)

    for i in 2:k
        c  = @views (Π[:,1:(i-1)]' * Ψ[:,1:(i-1)]) \ (Π[:,1:(i-1)]' * Ψ[:,i])
        r .= @views Ψ[:,i] .- Ψ[:,1:(i-1)] * c
        j = argmax(r)

        Π[j,i] = 1.0
    end
    return Π
end

ext: jl


text:
function extract_paths(tour, graph, graph_type)

    X = tour["X"]
    Y = tour["Y"]

    path = DataFrame(
            start_node = Int64[],
            end_node = Int64[],
            nb_passage = Float64[],
            serv = Bool[],
            )

    if nrow(graph["req_arcs"])>0
      for arc in 1:nrow(graph["req_arcs"])
        start_node = graph["req_arcs"][arc, :start_node]
        end_node = graph["req_arcs"][arc, :end_node]
        nb_passage = X["req_arcs"][arc] + Y["req_arcs"][arc]
        serv = X["req_arcs"][arc] >= 1
        if nb_passage > 10^-5
          push!(path, (start_node, end_node, nb_passage, serv))
        end
      end
    end

    if nrow(graph["noreq_arcs"])>0
      for arc in 1:nrow(graph["noreq_arcs"])
        start_node = graph["noreq_arcs"][arc, :start_node]
        end_node = graph["noreq_arcs"][arc, :end_node]
        nb_passage = Y["noreq_arcs"][arc]
        if nb_passage > 10^-5
          push!(path, (start_node, end_node, nb_passage, false))
        end
      end
    end

    if nrow(graph["req_edges"])>0
      for arc in 1:nrow(graph["req_edges"])
        start_node = graph["req_edges"][arc, :start_node]
        end_node = graph["req_edges"][arc, :end_node]
        nb_passage = X["req_edges"][arc] + Y["req_edges"][arc]
        serv = X["req_edges"][arc] >= 1
        if nb_passage > 10^-5
          push!(path, (start_node, end_node, nb_passage, serv))
        end
      end
    end

    if (nrow(graph["req_edges"])>0)&(graph_type=="directed")
      for arc in 1:nrow(graph["req_edges"])
        start_node = graph["req_edges"][arc, :end_node]
        end_node = graph["req_edges"][arc, :start_node]
        nb_passage = X["req_edges_swap"][arc] + Y["req_edges_swap"][arc]
        serv = X["req_edges_swap"][arc] >= 1
        if nb_passage > 10^-5
          push!(path, (start_node, end_node, nb_passage, serv))
        end
      end
    end

    if nrow(graph["noreq_edges"])>0
      for arc in 1:nrow(graph["noreq_edges"])
        start_node = graph["noreq_edges"][arc, :start_node]
        end_node = graph["noreq_edges"][arc, :end_node]
        nb_passage = Y["noreq_edges"][arc]
        if nb_passage > 10^-5
          push!(path, (start_node, end_node, nb_passage, false))
        end
      end
    end

    if (nrow(graph["req_edges"])>0)&(graph_type=="directed")
      for arc in 1:nrow(graph["noreq_edges"])
        start_node = graph["noreq_edges"][arc, :end_node]
        end_node = graph["noreq_edges"][arc, :start_node]
        nb_passage = Y["noreq_edges_swap"][arc]
        if nb_passage > 10^-5
          push!(path, (start_node, end_node, nb_passage, false))
        end
      end
    end

    return(path)
end

function display_results(tours, graph, data, solving_time)

  println("Solving time : ", solving_time)
  println("Total cost : ", sum([tour["cost"] for tour in tours]))
  println("Number of vehicule(s) used : ", length([tour for tour in tours if tour["cost"]!=0]))
  println("--------------")
  for (i,tour) in enumerate(tours)
    if tour["cost"] != 0
      println("Vehicle ", i)
      println("\tCost : ", tour["cost"])
      println("\tCapacity used : ", tour["capacity"])
      if (nrow(graph["req_arcs"])+nrow(graph["noreq_arcs"])==0) & (length(data[1, :depots])==1) & (data[1, :nb_vehicles]==1)
        graph_type = "undirected"
      else
        graph_type = "directed"
      end
      path = extract_paths(tour, graph, graph_type)
      println(path)
    end
  end

end

max_stars_repo_name: TheoGuyard/Groomer-Optimization-Problem

max_issues_repo_name: TheoGuyard/Groomer-Optimization-Problem

max_forks_repo_name: TheoGuyard/Groomer-Optimization-Problem


text:
using Luxor, Test

pi = π

function paramCircle(t,r)
    return Point(r*cos(t),r*sin(t))
end

function drawCircle(radius)
    start = 0
    step = 0.001
    endT = 2π
    prevPoint = paramCircle(start, radius)
    newpath()
    for t in step:step:endT
        currPoint = paramCircle(t,radius)
        line(prevPoint,currPoint)
        prevPoint = currPoint
    end
    strokepath()
end

function involuteCircle(t,r,a)
    xOfT = r*(cos(t)+(t-a)*sin(t))
    yOfT = r*(sin(t)-(t-a)*cos(t))
    return Point(xOfT,yOfT)
end

function drawInvolute(startRadian, radius; clockwise= true)
    start = startRadian
    if clockwise == true
        step = 0.001
        endT = startRadian+2*pi
    else
        step = -0.01
        endT = startRadian-2*pi
    end
    sethue("red")
    a = start
    prevPoint = involuteCircle(start,radius,a)
    newpath()
    for t in start+step:step:endT
        currPoint = involuteCircle(t,radius,a)
        line(prevPoint,currPoint)
        prevPoint = currPoint
    end
    strokepath()
end

function test_circleInvolute(fname)
    pagewidth, pageheight = 1200, 1400
    Drawing(pagewidth, pageheight, fname)
	origin()
	sethue("green")
	radius = 200
	radiusOffset = 50

	drawCircle(radius)
	sethue("blue")
	#drawCircle(radius+radiusOffset)

	rads = 0:pi/10:2*pi

	r = radius+radiusOffset+40

	for myrad in rads
		drawInvolute(myrad,radius)
		drawInvolute(myrad,radius,clockwise=false)
		#newpath()
		#sethue("black")
		#line(Point(0,0),Point(r*cos(myrad),r*sin(myrad)))
		#strokepath()
	end

	@test finish() == true
	println("...finished circle tangent test, saved in $(fname)")
end

fname = "circleInvolute_test.pdf"

test_circleInvolute(fname)
preview()


text:
using DataFrames, CSV;
using Serialization;

function calcEnt(datasets)
    num_samples = size(datasets, 1)
    label_counts = Dict()
    for label in datasets[end]
        if label in keys(label_counts)
            label_counts[label] += 1
        else
            label_counts[label] = 1
        end
    end
    # $ Ent(D) = - ∑_{p ∈ P} p log_2(p) $
    return - sum([p * log2(p) for p in values(label_counts) ./ num_samples ])
end

function calcGain(datasets, class, ent=calcEnt(datasets))
    num_samples = size(datasets, 1)
    groups = groupby(datasets, class)
    entries = [(size(subsets, 1), calcEnt(subsets)) for subsets in groups]
    # $Gain(D,a) = Ent(D) - ∑_{v ∈ V} |D^v| / |D| Ent(D^v)$
    return ent - sum([group_size / num_samples * grout_ent for (group_size, grout_ent) in entries])
end

function calcGainRatio(datasets, class, ent=calcEnt(datasets))
    num_samples = size(datasets, 1)
    groups = groupby(datasets, class)
    entries = [(size(subsets, 1), calcEnt(subsets)) for subsets in groups]
    gain =  ent - sum([group_size / num_samples * grout_ent for (group_size, grout_ent) in entries])
    # $IV = - ∑^V_{v=1} |D^v| / |D| log2(|D^v| / |D|)$
    IV = - sum([group_size / num_samples * log2(group_size / num_samples) for (group_size, grout_ent) in entries])
    # $Gain(D,a) = Gain(D,a) / IV$
    return gain / IV
end

function maxLabel(datasets)
    return by(datasets, :好瓜, df -> DataFrame(N=size(df, 1)))
end

function createTree(datasets, classes=names(datasets))
    groups = groupby(datasets, classes[end])
    if length(groups) == 1          # 取出分类标签（好瓜/坏瓜），如果只有一种，不需要继续分割
        return groups[end][end][end]
    end
    if length(classes) == 1         # 没有别的类别可以用来分割，取概率最大的分类标签
        max_group = groups[1]
        for group in groups[2:end - 1]
            if length(group) > length(max_group)
                max_group = group
            end
        end
        return max_group[end][end]
    end

    best_class = classes[1]         # 取出最好的分类类别，即信息增益最大的类别
    best_gain = calcGain(datasets, best_class)
    for class in classes[1:end - 1]
        curr_gain = calcGain(datasets, class)
        if curr_gain > best_gain
            best_class = class
            best_gain = curr_gain
        end
    end

    new_classes = copy(classes)     # 递归调用，下一个分支需要用到的分类标签
    filter!(x -> x != best_class, new_classes)   # 由于julia传引用，需要拷贝一份
    best_class_label = Set(datasets[best_class])
    tree = Dict(best_class => Dict())

    # 遍历该类别标准下所有子分类，递归求解
    for group in groupby(datasets, best_class)
        tree[best_class][group[best_class][1]] = createTree(group, new_classes)
    end

    return tree
end

datasets = CSV.read("data/2.0.csv", DataFrame);
tree = createTree(datasets)

print(tree)

describe(datasets)
by(datasets, :好瓜, df -> DataFrame(N=size(df, 1)))

stream = open("data/tree.data", "w")
serialize(stream,tree)
close(stream)

stream = open("data/tree.data", "r")
read_tree = deserialize(stream)
close(stream)

max_stars_repo_name: AlongWY/HIT-WorkShop

max_issues_repo_name: AlongWY/HIT-WorkShop

max_forks_repo_name: AlongWY/HIT-WorkShop


text:
"""
    meanfield(em, x0, kwargs)

Solve the mean field equations for the `AbstractEpidemicModel` `em`, from initial condition `x0`.

Keyword argument:
* `tmax=100.0`: Maximum time of the solution.
* `saveat=[]`: Times at which the solution should be computed. See the `DifferentialEquations` documentation on `solve` for more information.
"""
function meanfield(em::AbstractEpidemicModel, state_0; tmax=100.0, saveat=[])
    f! = meanfield_fun(em)
    u0 = init_state_mf(em, state_0)
    prob = ODEProblem(f!, u0, (0.0, tmax))
    sol = solve(prob, Tsit5(), saveat=saveat)
    finalize_meanfield(em, sol)
end

init_state_mf(em::AbstractEpidemicModel, x0) = float(x0)

finalize_meanfield(em::AbstractEpidemicModel, sol) = sol.t, sol.u, sol


text:
function myfunc()
more stuff
if something
print('something')
else
print('not something')
end
for something
print('something')
else
print('not something')
end
more stuff
end


text:
# This file is a part of Julia. License is MIT: https://julialang.org/license

using Test, Libdl, MPFR_jll

@testset "MPFR_jll" begin
    vn = VersionNumber(unsafe_string(ccall((:mpfr_get_version,libmpfr), Cstring, ())))
    @test vn == v"4.1.0"
end

max_stars_repo_name: jonas-schulze/julia

max_issues_repo_name: jonas-schulze/julia

max_forks_repo_name: jonas-schulze/julia


text:
# Bhattacharyya distances. Much like for KLDivergence we assume the vectors to
# be compared are probability distributions, frequencies or counts rather than
# vectors of samples. Pre-calc accordingly if you have samples.

type BhattacharyyaDist <: SemiMetric end

type HellingerDist <: Metric end

# Bhattacharyya coefficient

function bhattacharyya_coeff{T<:Number}(a::AbstractVector{T}, b::AbstractVector{T})
    n = length(a)
    sqab = zero(T)
    # We must normalize since we cannot assume that the vectors are normalized to probability vectors.
    asum = zero(T)
    bsum = zero(T)

    @simd for i = 1:n
        @inbounds ai = a[i]
        @inbounds bi = b[i]
        sqab += sqrt(ai * bi)
        asum += ai
        bsum += bi
    end

    sqab / sqrt(asum * bsum)
end

bhattacharyya_coeff{T <: Number}(a::T, b::T) = throw("Bhattacharyya coefficient cannot be calculated for scalars")

# Faster pair- and column-wise versions TBD...

# Bhattacharyya distance
evaluate{T<:Number}(dist::BhattacharyyaDist, a::AbstractVector{T}, b::AbstractVector{T}) = -log(bhattacharyya_coeff(a, b))
bhattacharyya(a::AbstractVector, b::AbstractVector) = evaluate(BhattacharyyaDist(), a, b)
evaluate{T <: Number}(dist::BhattacharyyaDist, a::T, b::T) = throw("Bhattacharyya distance cannot be calculated for scalars")
bhattacharyya{T <: Number}(a::T, b::T) = evaluate(BhattacharyyaDist(), a, b)

# Hellinger distance
evaluate{T<:Number}(dist::HellingerDist, a::AbstractVector{T}, b::AbstractVector{T}) = sqrt(1 - bhattacharyya_coeff(a, b))
hellinger(a::AbstractVector, b::AbstractVector) = evaluate(HellingerDist(), a, b)
evaluate{T <: Number}(dist::HellingerDist, a::T, b::T) = throw("Hellinger distance cannot be calculated for scalars")
hellinger{T <: Number}(a::T, b::T) = evaluate(HellingerDist(), a, b)


text:
function BoundaryProjection(arg0::Point, arg1::Point, arg2::jdouble)
    return BoundaryProjection((Point, Point, jdouble), arg0, arg1, arg2)
end

function get_offset(obj::BoundaryProjection)
    return jcall(obj, "getOffset", jdouble, ())
end

function get_original(obj::BoundaryProjection)
    return jcall(obj, "getOriginal", Point, ())
end

function get_projected(obj::BoundaryProjection)
    return jcall(obj, "getProjected", Point, ())
end


text:
using aoc15
using Test

@testset "memoryGameTests" begin

    @testset "startingListIsSpoken" begin
    # Write your tests here.
        input = [4,1]
        @test aoc15.spoken(input, 1) == 4
        @test aoc15.spoken(input, 2) == 1
    end

    @testset "newNumberWasAdded" begin
        input = [0,3,6]
        @test aoc15.spoken(input, 4) == 0
        @test aoc15.spoken(input, 8) == 0
        @test aoc15.spoken(input, 10) == 0
    end

    @testset "oldNumberIsFound" begin
        input = [0,3,6]
        @test aoc15.spoken(input, 5) == 3
        @test aoc15.spoken(input, 6) == 3
        @test aoc15.spoken(input, 7) == 1
        @test aoc15.spoken(input, 9) == 4
    end

    @testset "2020Examples" begin
        @test aoc15.spoken([0,3,6], 2020) == 436
        @test aoc15.spoken([1,3,2], 2020) == 1
        @test aoc15.spoken([2,1,3], 2020) == 10
        @test aoc15.spoken([1,2,3], 2020) == 27
        @test aoc15.spoken([2,3,1], 2020) == 78
        @test aoc15.spoken([3,2,1], 2020) == 438
        @test aoc15.spoken([3,1,2], 2020) == 1836
    end

    @testset "30000000Examples" begin
        @test aoc15.spoken([0,3,6], 30000000) == 175594
        @test aoc15.spoken([1,3,2], 30000000) == 2578
        @test aoc15.spoken([2,1,3], 30000000) == 3544142
        @test aoc15.spoken([1,2,3], 30000000) == 261214
        @test aoc15.spoken([2,3,1], 30000000) == 6895259
        @test aoc15.spoken([3,2,1], 30000000) == 18
        @test aoc15.spoken([3,1,2], 30000000) == 362
    end

end


text:
function get_covariance(obj::StatisticalMultivariateSummary)
    return jcall(obj, "getCovariance", RealMatrix, ())
end

function get_dimension(obj::StatisticalMultivariateSummary)
    return jcall(obj, "getDimension", jint, ())
end

function get_geometric_mean(obj::StatisticalMultivariateSummary)
    return jcall(obj, "getGeometricMean", Vector{jdouble}, ())
end

function get_max(obj::StatisticalMultivariateSummary)
    return jcall(obj, "getMax", Vector{jdouble}, ())
end

function get_mean(obj::StatisticalMultivariateSummary)
    return jcall(obj, "getMean", Vector{jdouble}, ())
end

function get_min(obj::StatisticalMultivariateSummary)
    return jcall(obj, "getMin", Vector{jdouble}, ())
end

function get_n(obj::StatisticalMultivariateSummary)
    return jcall(obj, "getN", jlong, ())
end

function get_standard_deviation(obj::StatisticalMultivariateSummary)
    return jcall(obj, "getStandardDeviation", Vector{jdouble}, ())
end

function get_sum(obj::StatisticalMultivariateSummary)
    return jcall(obj, "getSum", Vector{jdouble}, ())
end

function get_sum_log(obj::StatisticalMultivariateSummary)
    return jcall(obj, "getSumLog", Vector{jdouble}, ())
end

function get_sum_sq(obj::StatisticalMultivariateSummary)
    return jcall(obj, "getSumSq", Vector{jdouble}, ())
end


text:
using FrankWolfe
using ProgressMeter
using Arpack
using Plots
using DoubleFloats
using ReverseDiff

using LinearAlgebra

n = Int(1e2);
k = n

f(x) = dot(x, x)
function grad!(storage, x)
    @. storage = 2 * x
end

# pick feasible region
lmo = FrankWolfe.ProbabilitySimplexOracle{Rational{BigInt}}(1); # radius needs to be integer or rational

# compute some initial vertex
x0 = FrankWolfe.compute_extreme_point(lmo, zeros(n));

# benchmarking Oracles
FrankWolfe.benchmark_oracles(f, grad!, () -> rand(n), lmo; k=100)

# the algorithm runs in rational arithmetic even if the gradients and the function itself are not rational
# this is because we replace the descent direction by the directions of the LMO are rational

@time x, v, primal, dual_gap, trajectory = FrankWolfe.frank_wolfe(
    f,
    grad!,
    lmo,
    x0,
    max_iteration=k,
    line_search=FrankWolfe.Agnostic(),
    print_iter=k / 10,
    verbose=true,
    emphasis=FrankWolfe.blas,
);

println("\nOutput type of solution: ", eltype(x))

# you can even run everything in rational arithmetic using the shortstep rule
# NOTE: in this case the gradient computation has to be rational as well

@time x, v, primal, dual_gap, trajectory = FrankWolfe.frank_wolfe(
    f,
    grad!,
    lmo,
    x0,
    max_iteration=k,
    line_search=FrankWolfe.RationalShortstep(),
    L=2,
    print_iter=k / 10,
    verbose=true,
    emphasis=FrankWolfe.blas,
);

println("\nOutput type of solution: ", eltype(x))

println("\nNote: the last step where we exactly close the gap. This is not an error. ")
fract = 1 // n
println(
    "We have *exactly* computed the optimal solution with with the $fract * (1, ..., 1) vector.\n",
)
println("x = $x")

####################################################################################################################
### APPRROXIMATE CARATHEODORY WITH PLANTED SOLUTION
####################################################################################################################

rhs = 1
n = 40
k = 1e5

xpi = rand(big(1):big(100), n)
total = sum(xpi)
xp = xpi .// total

f(x) = norm(x - xp)^2
function grad!(storage, x)
    @. storage = 2 * (x - xp)
end

lmo = FrankWolfe.ProbabilitySimplexOracle{Rational{BigInt}}(rhs)
direction = rand(n)
x0 = FrankWolfe.compute_extreme_point(lmo, direction)

@time x, v, primal, dual_gap, trajectory = FrankWolfe.frank_wolfe(
    f,
    grad!,
    lmo,
    x0,
    max_iteration=k,
    line_search=FrankWolfe.Agnostic(),
    print_iter=k / 10,
    emphasis=FrankWolfe.memory,
    verbose=true,
)
println("\nOutput type of solution: ", eltype(x))

println("Computed solution: x = $x")


text:
module Queller

using TextWrap

export main,
	check_queller_graphs

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

not_unique(v, f=isequal) = filter(e -> (count(f.([e],v)) > 1), v)
strvec2str(v,sep='\n') = isempty(v) ? "" : reduce((s,t) -> s*sep*t, v)
fieldsdefined(obj) = all(i -> isdefined(obj,i), 1:nfields(obj))

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

include("dice_and_strategy.jl")
include("cli.jl")
include("graph.jl")
include("quellerstate.jl")
include("crawler.jl")

const PKG_DIR = abspath(joinpath(@__DIR__, ".."))
const GRAPH_FILES = filter(p-> splitext(p)[2] == ".jl", readdir("$(PKG_DIR)/graphs", join=true))
const GRAPHS = load_graphs(GRAPH_FILES...)

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

mutable struct ProgramState
	phase::Int
	phases::Vector{Function}
	graphs::Dict{String,QuellerGraph}

	available_dice::Vector{Die.Face}

	iop::IOParser

	reset::Bool
	reset_phase::Bool
	exit::Bool
end

function ProgramState()
	phase = 1
	phases = [phase1,phase2,phase3,phase4,phase5]

	available_dice = Vector{Die.Face}()

	iop = IOParser([
		CMD.ResetPhase()
		CMD.Exit()
		CMD.Repeat()
		CMD.Phase(length(phases))
		CMD.Help()
		])

	cmds = Vector{CMD.Command}()
	return ProgramState(phase,phases,GRAPHS,available_dice,iop,false,false,false)
end

greeting_str = """
Queller CLI: War of the Ring Shadow AI

Type 'help' and press enter for more information regarding the operation of the program.
"""

help_str = """
Queller CLI: War of the Ring Shadow AI

-- Usage --

This program will present a number of statements, selections, and actions for each of the 5 game phases. Simply answer/perform these when they are presented to you. All statments refer to the Shadow player's hand, armies etc. unless specifically stated that it refers to the Free People's hand, armies etc.

Read the short manual before starting and have it available. Words between asterisks, e.g., *threat*, refers to the glossary in the manual.

-- Dice --

The program will keep track of the dice available to the Shadow player and instruct you on which to use. However, it is possible, due to both mistakes and legitimate reasons, that the program will get out of sync with reality. There is a menu option in Phase 5 for correcting this.

-- Inputs --

The valid options for a query are given by the input prompt. If the options are separated by "/" only one of the options should be given. If they are separated by "," several options can be given. To select an options, simply type in the selected option and press enter, multiple options should be separated by spaces or nothing at all.

All inputs are case insensitive and can in many cases be shorten to one letter, i.e., "t" for "true", "f" for "false", "u" for "undo" etc.

Examples:

[true/false] > true
[true/false] > f
[1/2/3/4] > 3
[1,2,3,4] > 2 4 1
[C,A,M,H,P] > mappc

-- Commands --

The following commands can (almost) always be used when prompted for input. To use a command type the name of the command and press enter.

help        :: Shows this help message.
exit        :: Exits the program.
undo        :: Undo the latest input and step back. Note, everything can not be undone, use the "reset" or "phase" command in these circumstances.
repeat      :: Repeats the latest query.
reset       :: Reset and restart the current phase.
phase <nbr> :: Jumps to the beginning of a phase, e.g., type "phase 3" to jump to phase 3

----------
"""

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

function print_read_process(state, msg, options, callback, silent_options=Vector{CMD.Command}())
	display_message(state.iop, msg)

	while true
		cmd = read_input(state.iop, options, silent_options)

		if cmd in options || cmd in silent_options
			callback(cmd)
			return true
		end

		cmd isa CMD.Command && handle_general_command(state,cmd)
		cmd isa CMD.Repeat && return print_read_process(state, msg, options, callback, silent_options)
		cmd isa CMD.AbortingCommand && return false
	end
end

function resolve_decision_graph(state, graph)
	gc = GraphCrawler(graph, state.graphs, state.available_dice)

	abort = Ref(false)
	callback(cmd) = proceed!(gc, cmd)
	callback(cmd::CMD.Undo) = !undo!(gc) && (abort[] = !print_read_process(state, "Cannot undo more.", [CMD.Blank()], x->nothing))

	while !at_end(gc)
		msg, options = getinteraction(gc)
		!print_read_process(state, msg, options, callback, [CMD.Undo()]) && return false
		abort[] && return false
	end

	state.available_dice = get_available_dice(gc)
	return true
end

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

function main()::Cint
	try
		state = ProgramState()
		display_message(state.iop, greeting_str)

		while !state.exit
			state.phases[state.phase](state)

			state.reset_phase && (state.reset_phase = false; continue)

			state.phase = 1 + (state.phase % length(state.phases))
		end
	catch e
		println(e)
		return 1
	end
	return 0
end

handle_general_command(state,cmd::CMD.Command) = nothing
handle_general_command(state,cmd::CMD.ResetPhase) = (state.reset_phase = true)
handle_general_command(state,cmd::CMD.Exit) = (state.exit = true)
handle_general_command(state,cmd::CMD.Phase) = (state.reset_phase = true; state.phase = cmd.nbr)
handle_general_command(state,cmd::CMD.Help) = display_message(state.iop, help_str)

phase1(state) = graph_phase(state, "phase_1", "Phase 1")
phase2(state) = graph_phase(state, "phase_2", "Phase 2")
phase3(state) = graph_phase(state, "phase_3", "Phase 3")
phase4(state) = graph_phase(state, "phase_4", "Phase 4")

function graph_phase(state, graph, name)
	display_message(state.iop, name)
	resolve_decision_graph(state, graph)
end

function phase5(state)
	menu = """
	Available Dice: $(strvec2str(sort(Die.char.(state.available_dice)),','))

	Select phase 5 action.

	1. Choose Shadow action
	2. Resolve a card effect
	3. Resolve a battle
	4. Recruit a minion as Shadow's final action with die set aside earlier
	5. Change the available dice (use this if the bot's available dice do not match reality)
	6. End turn and go to phase 1
	"""

	graphs = [
			  "phase_5",
			  "event_cards_resolve_effect",
			  "battle",
			  "muster_minion_selection",
			  "adjust_dice",
			  ]

	choice = Ref(0)
	menu_callback(cmd) = (choice[] = cmd.opt)
	!print_read_process(state, menu, CMD.Option.(1:6), menu_callback) && return

	if 1 <= choice[] <= 5
		!resolve_decision_graph(state, graphs[choice[]]) && return

	elseif choice[] == 6
		# Confirm exit
		msg = "Exit phase 5 and return to phase 1."
		exit = Ref(false)

		exit_callback(cmd::CMD.False) = (exit[] = false)
		exit_callback(cmd::CMD.True) = (exit[] = true)

		!print_read_process(state, msg, [CMD.True(), CMD.False()], exit_callback) && return

		exit[] && return
	end

	phase5(state)
end

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

include("graphviz.jl")

function check_queller_graphs()
	println("\nAll graphs in PostScript can be found in Queller/graph_output.")
	graph_output_file(f) = joinpath(PKG_DIR, "graph_output", splitext(basename(f))[1]*".ps")

	for f in getfield.(values(GRAPHS),:source_file)
		graph2ps(f, graph_output_file(f))
	end

	unjumped = get_graphs_not_jumped_to(values(GRAPHS))
	println("\nGraphs that are not jumped to from another graph:\n$(strvec2str(unjumped))")

	return
end

end

max_stars_repo_name: mvmorin/queller-bot

max_issues_repo_name: mvmorin/queller-bot

max_forks_repo_name: mvmorin/queller-bot


text:
using HMRowUnification
using MLStyle
using Test

@testset "HMRowUnification.jl" begin
    # Write your tests here.

    using HMRowUnification

    tctx = HMT[]

    vids = UInt[]
    st = mk_tcstate(tctx, x -> push!(vids, x))

    tvar1 = st.new_tvar()
    tvar2 = st.new_tvar()
    int_t = Nom(:int)

# int -> int ~ 'tvar1 -> 'tvar2
# =>
# 'tvar1 = int
# 'tvar2 = int
    @test st.unify(Arrow(int_t, int_t), Arrow(tvar1, tvar2))

# true

    @test st.prune.([tvar1, tvar1]) == [int_t, int_t]
# 2-element Array{Nom,1}:
# Nom(:int)
# Nom(:int)

    tv3 = st.new_tvar()
    tv4 = st.new_tvar()

    rec1 = row_of_pairs(
    Dict(:a => int_t, :b => int_t),
    RowPoly(tv3)
)
    rec2 = row_of_pairs(
    Dict(:a => int_t, :c => int_t),
    RowPoly(tv4)
)

    @test Set([st.int_of_tvar(tv4).value]) == ftv(Record(rec2))

    @test st.unify(Record(rec1), Record(rec2))
    println(st.prune(tv3), "  ", st.prune(tv4))
    # @test @match st.prune(tv3), st.prune(tv4) begin
    #     (Record(rec1), Record(rec2)) =>
    #         begin
    #             m1, rowt′ = st.extract_row(rec1)
    #             m2, rowt = st.extract_row(rec2)
    #             @test m1 == m2
    #             sort(collect(m1)) == sort([
    #                 :a => int_t,
    #                 :b => int_t,
    #                 :c => int_t
    #             ])
    #         end
    # end

    t1 = Forall((:a, :b), App(Fresh(:b), Fresh(:a)))
    t2 = Forall((:b, :a), App(Fresh(:b), Fresh(:a)))
    @test t1 ⪯ t2
    @test t2 ⪯ t1

    @test st.unify(t1, t2)

    t3 = Forall((:b, :a), App(Fresh(:a), Fresh(:b)))
    t3′ = Forall((:a, ), App(Fresh(:a), Fresh(:a)))
    @test t3′ ⪯ t1
    @test (t1 ⪯ t3′) == false

    @test st.unify(t1, t3)

    println(t3)
    println(Genvar(0, :a))
    println(rec2)
    @test repr(Record(RowMono)) == "{}"

    tv5 = st.new_tvar()
    arrow_t = Arrow(int_t, Tup((int_t, tv5)))
    tv6 = st.new_tvar()
    @test st.unify(Arrow(tv6, Tup((tv5, tv6))), arrow_t)
    @test st.prune(tv5) == st.prune(tv6) == int_t
    println(st.prune(arrow_t))

    t4 = Forall((:b, :a), App(Fresh(:a), Fresh(:a)))
    @test !(st.unify(t1, t4))

    @test isempty(ftv(t4))
    @test st.unify(st.new_tvar(), Fresh(:a))

    @test int_t ⪯ Forall((:a, ), Fresh(:a))
    @test int_t ⪯ int_t
    @test Arrow(int_t, int_t) ⪯ Forall((:a, ), Arrow(Fresh(:a), Fresh(:a)))
    @test Arrow(
            Forall((:a, ),
                Arrow(Fresh(:a), Fresh(:a))),
            int_t) ⪯
            Arrow(Arrow(int_t, int_t), int_t)
    tv = st.new_tvar()
    @test Tup((int_t, int_t)) ⪯ Forall((:a, ), Tup((Fresh(:a), Fresh(:a))))
    @test st.unify(Fresh(:a), Fresh(:a))
    @test Fresh(:a) ⪯ Fresh(:a)
    println(st.instantiate(t4))
    @test isempty(unbound(t4))
    @test collect(unbound(Tup((int_t,  Fresh(:a))))) == [:a]
    l = st.new_tvar()
    println(st.generalise(Dict(l.var.i => :a), l) == Forall((:a, ), Fresh(:a)))
    x, y, a, b = [st.new_tvar() for i = 1:4]
    @test st.unify(Arrow(Tup((x, y)), x), Arrow(a, b))
    @test st.prune(x) == st.prune(b)
    println(st.genvar_links)
    st.unlink(0)

    a = Fresh(:a)
    v = st.new_tvar()
    level = length(st.genvar_links)
    @test level == 0
    st.type_less(Arrow(v, v), Forall((:a, ), Arrow(a, a)))
    println(st.genvar_links)
    st.unlink(level)
    println(st.genvar_links)

    t1 = Forall((:a, ), Record(
        RowCons(:a, int_t,
        RowCons(:b, int_t,
        RowPoly(a)))))

    t2 = Forall((:a, ), Record(
        RowCons(:a, int_t,
        RowCons(:c, int_t,
        RowPoly(a)))))

    t1_ = st.instantiate(t1)
    t2_ = st.instantiate(t2)
    @test st.unify(t1_, t2_)
    t1_ = st.instantiate(t1)
    t2_ = st.instantiate(t2)
    @test (t1_ ⪯ t2_) == false
    println(st.prune(t2_))
    println(vids)

    a = Fresh(:a)
    s = Forall((:a, ), Arrow(a, a))
    @test Arrow(s, s) ⪯ s
    @test !(s ⪯ Arrow(s, s))
end


text:
"""
Dancing Links
=============

Finds the first solution of an exact cover problem using Donald Knuth's Dancing Links technique.

References
----------

- The "Dancing Links" paper: https://arxiv.org/pdf/cs/0011047.pdf
- A good tutorial: http://garethrees.org/2007/06/10/zendoku-generation/#section-4
- C implementation: http://maskray.me/blog/2009-11-23-sudoku-with-dancing-links-and-algorithm-x
"""
module DancingLinks

export DLinks, init_headers!, insert!, cover!, uncover!, dlx!

mutable struct DLinks
    L::Vector{Int}
    R::Vector{Int}
    U::Vector{Int}
    D::Vector{Int}

    h::Int
    C::Vector{Int}  # points to the column object
    S::Vector{Int}  # the size of column (number of 1s)
end

"Initializes root and column headers"
function init_headers!(d::DLinks, n_cols::Int)
    d.h = n_cols + 1
    for col in 1 : d.h
        d.L[col] = mod1(col - 1, d.h)
        d.R[col] = mod1(col + 1, d.h)
        d.U[col] = col
        d.D[col] = col
    end
end

function insert!(d::DLinks, cur::Int, left::Int, right::Int, col::Int)
    # Row
    d.L[cur] = left
    d.R[left] = cur

    d.R[cur] = right
    d.L[right] = cur

    # Column
    d.U[cur] = d.U[col]
    d.D[d.U[cur]] = cur

    d.D[cur] = col
    d.U[col] = cur

    # Header
    d.C[cur] = col
    d.S[col] += 1
end

"Removes the column."
function cover!(d::DLinks, col::Int)
    d.R[d.L[col]] = d.R[col]
    d.L[d.R[col]] = d.L[col]

    for r in Loop(d.D, col)
        for c in Loop(d.R, r)
            d.D[d.U[c]] = d.D[c]
            d.U[d.D[c]] = d.U[c]

            d.S[d.C[c]] -= 1
        end
    end
end

"Puts the column back."
function uncover!(d::DLinks, col::Int)
    d.R[d.L[col]] = col
    d.L[d.R[col]] = col

    for r in Loop(d.U, col)
        for c in Loop(d.L, r)
            d.D[d.U[c]] = c
            d.U[d.D[c]] = c

            d.S[d.C[c]] += 1
        end
    end
end

"""
The "Algorithm X".

It's a depth-first search with backtracking.
"""
function dlx!(d::DLinks, O::Vector{Int})::Bool
    if d.R[d.h] == d.h
        return true
    end

    # finds column with the minimum size
    col = begin
        min_col = nothing
        s_min = typemax(Int)
        for col in Loop(d.R, d.h)
            if d.S[col] < s_min
                s_min = d.S[col]
                min_col = col
            end
        end
        min_col
    end

    cover!(d, col)
    push!(O, 0)

    for r in Loop(d.D, col)
        O[end] = r

        for c in Loop(d.R, r)
            cover!(d, d.C[c])
        end

        if dlx!(d, O)
            return true  # only the first
        end

        for c in Loop(d.L, r)
            uncover!(d, d.C[c])
        end
    end

    uncover!(d, col)
    false
end

mutable struct Loop
    circle::Vector{Int}
    cur::Int
    stop::Int

    Loop(circle::Vector{Int}, start::Int) = new(circle, start, start)
end

function Base.iterate(iter::Loop, _state=iter)
    next = iter.circle[iter.cur]
    if next == iter.stop
        nothing
    else
        iter.cur = next
        (iter.cur, iter)
    end
end

include("Sudoku.jl")

end  # module

max_stars_repo_name: gyk/TrivialSolutions

max_issues_repo_name: gyk/TrivialSolutions

max_forks_repo_name: gyk/TrivialSolutions


text:
########################################
#              Preamble
########################################
using PyCall
using FFTW
using LinearAlgebra
using JLD2
using YAML
using CosmoMMF

config = YAML.load_file("config.yaml")

########################################
#            Load in Data
########################################

@load config["output_directory"]*config["run_name"]*"_"*config["snapfile_root"]*".jld2" den

########################################
#  Run n' Save NEXUS/NEXUS+ Algorithm(s)
########################################

Rs = (√2) .^ (0:config["num_scales"])
print("\n loading different scales \n")

save_name = config["run_name"]*"_"*config["snapfile_root"]

if ARGS[1] == "NEXUS"
    @time combined_NEXUS = CosmoMMF.maximum_signature(Rs, den; alg=:NEXUS);
    CosmoMMF.save_max_sigs(config["output_directory"], "Sigs_Nexus_"*save_name, combined_NEXUS)
    combined_NEXUS = nothing
elseif ARGS[1] == "NEXUSPLUS"
    @time combined_NEXUSPLUS = CosmoMMF.maximum_signature(Rs, den; alg=:NEXUSPLUS);
    CosmoMMF.save_max_sigs(config["output_directory"], "Sigs_NexusPlus_"*save_name, combined_NEXUSPLUS)
    combined_NEXUSPLUS = nothing
elseif ARGS[1] == "BOTH"
    @time combined_NEXUS = CosmoMMF.maximum_signature(Rs, den; alg=:NEXUS);
    CosmoMMF.save_max_sigs(config["output_directory"], "Sigs_Nexus_"*save_name, combined_NEXUS)
    combined_NEXUS = nothing

    @time combined_NEXUSPLUS = CosmoMMF.maximum_signature(Rs, den; alg=:NEXUSPLUS);
    CosmoMMF.save_max_sigs(config["output_directory"], "Sigs_NexusPlus_"*save_name, combined_NEXUSPLUS)
    combined_NEXUSPLUS = nothing
else
    printstyled("Please give a valid algorithm input argument: \n NEXUS, NEXUSPLUS, or BOTH"; color = :red)

########################################
#            Miscellaneous
########################################

# print(ARGS)

# x = parse(Int32,ARGS[3])
# y = parse(Int32,ARGS[4])

# print("\n")
# print(x+y)
# print("\n")


text:
import GeometryTypes: numtype
@testset "Convex Hulls" begin
    T = Float64
    V = Vec{2, T}
    s = Simplex((Vec(0.0, 0.), Vec(0.0, 1.0), Vec(1.0, 0.0)))
    fs = FlexibleSimplex([Vec(0.0, 0.0), Vec(0.0, 1.0), Vec(1.0, 0.0)])
    fh = FlexibleConvexHull([Vec(0.0, 0.0), Vec(0.0, 1.0), Vec(1.0, 0.0)])
    types = (Simplex, FlexibleSimplex, FlexibleConvexHull)
    objects = (s,fs,fh)

    @testset "conversions" begin
        @test typeof(@inferred FlexibleSimplex(s)) == FlexibleSimplex{V}
        @test typeof(@inferred FlexibleConvexHull(s)) == FlexibleConvexHull{V}
        @test typeof(@inferred FlexibleConvexHull(fs)) == FlexibleConvexHull{V}
        for i1 in objects, i2 in objects
            @test isapprox(FlexibleConvexHull(i1), i2)
        end
    end

    @testset "Utility functions" begin

        v_matrix = [0.0 0.0 1.0; 0.0 1.0 0.0]
        v_mat = Mat{2, 3}(v_matrix)
        for shape in objects
            @test (@inferred eltype(shape)) == V
            @test (@inferred spacedim(shape)) == 2
            @test nvertices(shape) == 3

            @test (@inferred vertexmatrix(shape)) == v_matrix
            @test vertexmat(shape) == v_mat

        end
        @test (@inferred vertexmat(s)) == v_mat
        @test (@inferred nvertices(s)) == 3
    end

    @testset "isapprox" begin
        s2 = rand(Simplex{3, Vec{2, Float64}})

        @test !isapprox(s, s2)
        @test isapprox(s, s2, atol = 100.0)
    end

    @testset "Rects" begin
        c = HyperCube(Vec(1.0, 2.0), 1.0)
        r = HyperRectangle(Vec(1.0, 2.0), Vec(1.0, 1.0))
        fh = FlexibleConvexHull(Vec{2, Float64}[(1.0, 2.0), (1.0, 3.0), (2.0,2.0), (2.0,3.0)])
        objects = (c, r, fh)
        @test (@inferred convert(HyperRectangle, c)) == r
        fh2 = (@inferred convert(FlexibleConvexHull, c))
        @test vertices(fh2) == vertices(fh)
        for x in objects, y in objects
            @test isapprox(x, y)
        end
        for x in objects
            @test nvertices(c) == 4
            @test spacedim(c) == 2
            @test numtype(c) == Float64
        end
    end
end


text:
using ForwardDiff, LinearAlgebra

struct TimingData{T<:AbstractFloat}
    N::Vector{Int64}          # Number of transits for each body
    nt::Vector{Vector{Int64}} # Transit numbers
    tt::Vector{Vector{T}}     # Transit times for each body
    tterr::Vector{Vector{T}}  # Transit time errors
    BJD::T                    # Timing offset
end

function chi_squared(model::Vector{T}, data::Vector{T}, errors::Vector{T}) where T<:Real
    res = (data .- model) ./ errors
    return dot(res,res)
end

function chi_squared(model::T, data::T, error::T) where T<:Real
    res = (data - model) / error
    return res * res
end

function ∇chi_squared(model::Vector{T}, dmodel::Matrix{T}, data::Vector{T}, errors::Vector{T}, dll::Vector{T}) where T<:Real
    res = (data .- model) ./ (errors .* errors)
    # Iterate over each model parameter
    for i in eachindex(dmodel[:,1])
        dll[i] -= @views 2.0 * dot(dmodel[i,:], res)
    end
end

function ∇chi_squared(model::T, dmodel::Vector{T}, data::T, error::T, dll::Vector{T}) where T<:Real
    res = (data - model) / error
    for i in eachindex(dmodel)
        dll[i] -= @views 2.0 .* dmodel[:] * res / error
    end
    return
end

# Hessian of nll wrt model, data, error
chi_squared(p) = chi_squared(p[1], p[2], p[3])
hess_chi_squared = p -> ForwardDiff.hessian(chi_squared, p)

# Bounds for optimizer
function get_upper_open(elements)
    m,P,t0,ecosϖ,esinϖ,I,Ω = unpack_elements(elements)
    N = length(m)
    ms = 1.1 .* m[:]
    Ps = 1.2 .* P[:]
    t0s = t0[:] .+ (1.1 .* Ps[:])
    ecs = ones(N) .* 0.5
    ess = ones(N) .* 0.5
    Is = ones(N) .* (π/2 + π/8)
    Ωs = ones(N) .* 0.1
    return [ms...,Ps...,t0s...,ecs...,ess...,Is...,Ωs...]
end

function get_lower_open(elements)
    m,P,t0,ecosϖ,esinϖ,I,Ω = unpack_elements(elements)
    N = length(m)
    ms = ones(N) .* 1e-6
    Ps = ones(N) .* 0.9
    t0s = t0[:] .- (2.0 .* Ps[:])
    ecs = -ones(N) .* 0.5
    ess = -ones(N) .* 0.5
    Is = ones(N) .* (π/2 - π/8)
    Ωs = ones(N) .* -0.1
    return [ms...,Ps...,t0s...,ecs...,ess...,Is...,Ωs...]
end

unpack_elements(el) = [el[:,i] for i in 1:7]

function logP(lnlike::Function, elements::Matrix{T}, H::Vector{Int64}, data::TimingData, t0::T,
    intr::Integrator; parameter_mask=nothing, grad::Bool=true) where T <: Real

    ## Check things here ##
    if parameter_mask == nothing
        # Assume all planet parameters, and fully-nested orbits
        parameter_mask = zeros(Bool, size(elements))
        parameter_mask[2:end,:] .= true
    end

    # Setup and run simulation
    ic = ElementsIC(t0,H,elements)
    s = State(ic)
    tt = TransitTiming(intr.tmax, ic)
    intr(s,tt,grad=grad)

    # calculate -loglikelihood function
    params = elements[parameter_mask][:]
    if @generated
        nll = :(0.0)
        for i in 1:length(data.N)
            nll = :($nll + lnlike(tt.tt[$i + 1, data.nt[$i]],data.tt[$i],data.tterr[$i],grad))
        end
        return :($nll)
    else
        nll = 0.0
        for i in 1:length(data.N)
            nll += lnlike(tt.tt[i+1,data.nt[i]],data.tt[i],data.tterr[i],grad)
        end
        return nll
    end
end

function logP(lnlike::Function, dlnlike::Function, elements::Matrix{T}, H::Vector{Int64},
    data::TimingData, t0::T, intr::Integrator; parameter_mask=nothing, grad=true) where T <: Real

    ## Check things here ##
    if parameter_mask == nothing
        # Assume all planet parameters, and fully-nested orbits
        parameter_mask = zeros(Bool, size(elements))
        parameter_mask[2:end,:] .= true
    end

    # Check eccentricity
    e = zeros(length(elements[2:end,1]))
    e = sqrt.(elements[:,4].^2 .+ elements[:,5].^2)
    if any(e .> 1.0) || any(e .< 0.0)
        if grad; return 1e30, ones(length(elements[2:end,:][:])) * 1e30; end
        return 1e30
    end

    # Setup and run simulation
    ic = ElementsIC(t0, H, elements)
    s = State(ic)
    tt = TransitTiming(intr.tmax, ic)
    intr(s, tt)

    # calculate -loglikelihood function
    params = elements[parameter_mask][:]
    N_params = sum(parameter_mask)
    ll = 0.0
    for i in 1:length(data.N)
        ll += lnlike(tt.tt[i+1, data.nt[i]], data.tt[i], data.tterr[i])
    end

    # Now gradients
    dll = zeros(T, N_params)
    N = length(data.N); n = sum(parameter_mask[:,1])
    inds = [7,1,2,3,4,5,6]
    for i in 1:N # <- Make this a list of planet indices..
        dmodel = zeros(N_params, length(data.tt[i])) # Find a way to allocate this outside loop...
        for k in 1:data.N[i]
            dmodel[:,k] = tt.dtdelements[i+1, data.nt[i][k], inds, :]'[parameter_mask][:]
        end
        dlnlike(tt.tt[i+1, data.nt[i]], dmodel, data.tt[i], data.tterr[i], dll)
    end

    # Add priors here
    return ll, dll
end

function calc_hessian(hess_func::Function, elements::Matrix{T}, H::Vector{Int64},
    data::TimingData, t0::T, intr::Integrator; parameter_mask=nothing) where T<:Real

    ## Check things here ##
    if parameter_mask == nothing
        # Assume all planet parameters, and fully-nested orbits
        parameter_mask = zeros(Bool, size(elements))
        parameter_mask[2:end,:] .= true
    end

    # Setup and run simulation
    ic = ElementsIC(t0, H, elements)
    s = State(ic)
    tt = TransitTiming(intr.tmax, ic)
    intr(s, tt)

    # Calculate hessian
    inds = [7,1,2,3,4,5,6]
    N_params = sum(parameter_mask); N = length(data.N)
    hessian = zeros(T, N_params, N_params)
    for i in 1:N
        dmodel = zeros(N_params, length(data.tt[i]))
        hmodel = zeros(length(data.tt[i]))
        for k in 1:data.N[i]
            dmodel[:,k] .= tt.dtdelements[i+1, data.nt[i][k], inds, :]'[parameter_mask][:]
        end
        for k in 1:N_params, l in 1:N_params
            for j in 1:length(data.tt[i])
                hmodel[j] = hess_func([tt.tt[i+1, data.nt[i][j]], data.tt[i][j], data.tterr[i][j]])[1]
            end
            hessian[k,l] -= sum(hmodel .* dmodel[k,:] .* dmodel[l,:] ./ (data.tterr[i].^2))
        end
    end
    return 2.0 .* hessian
end


text:
using Test, Distributed, Random, LinearAlgebra, Distances, Pathogen

addprocs(3)

@everywhere using Pathogen
@everywhere include(joinpath(@__DIR__, "risk_functions.jl"))

# Set RNG seed
Random.seed!(4321)

# using Logging
# global_logger(ConsoleLogger(stderr, LogLevel(-10000)))

# Define population
n = 25
risks = DataFrame(x = rand(Uniform(0, 5), n),
                  y = rand(Uniform(0, 5), n),
                  riskfactor1 = rand(Gamma(), n))

pop = Population(risks,
                 [euclidean([risks[i, :x]; risks[i, :y]], [risks[j, :x]; risks[j, :y]]) for i = 1:n, j = 1:n])

@testset "SEIR Model" begin
  # Some commonly used functions/examples provided in helpers/RiskFunctions.jl
  # For SEIR, risk functions and parameters in order of: sparks, susceptibility, infectivity, transmissibility, latency, and removal
  rf = RiskFunctions{SEIR}(_constant,
                           _coefficient,
                           _powerlaw,
                           _one,
                           _constant,
                           _constant)

  rparams = RiskParameters{SEIR}([0.01],
                                 [1.0],
                                 [2.0, 5.0],
                                 Float64[],
                                 [0.1],
                                 [0.05])

  sim = Simulation(pop,
                   rf,
                   rparams)

  simulate!(sim, tmax=100.0)

  state_counts = [Pathogen._count_by_state(sim.events, State_S, 50.0)
                  Pathogen._count_by_state(sim.events, State_E, 50.0)
                  Pathogen._count_by_state(sim.events, State_I, 50.0)
                  Pathogen._count_by_state(sim.events, State_R, 50.0)]

  @test sum(state_counts) == n
  @test length(sim.disease_states) == n
  @test size(sim.transmission_network.internal) == (n, n)

  obs = observe(sim, Uniform(), Uniform())
  @test length(obs.infection) == n

  rpriors = RiskPriors{SEIR}([Uniform(0.0, 0.1)],
                             [Uniform(0.0, 2.0)],
                             [Uniform(0.0, 4.0); Uniform(1.0, 8.0)],
                             UnivariateDistribution[],
                             [Uniform(0.0, 1.0)],
                             [Uniform(0.0, 1.0)])
  ee = EventExtents{SEIR}((2.0, 20.0), 2.0, 2.0)
  mcmc = MCMC(obs, ee, pop, rf, rpriors)
  start!(mcmc, 3, attempts=100)
  @test length(mcmc.markov_chains) == 3
  @test all([length(mcmc.markov_chains[i].risk_parameters[1]) for i=1:3] .== length(rparams))

  # Check bounds of Events initialization
  for j = 1:1, i = 1:n
    if mcmc.markov_chains[j].events[1].exposure[i] > -Inf
      l, u = Pathogen._bounds(i, State_E, ee, obs, mcmc.markov_chains[j].events[1], mcmc.markov_chains[j].transmission_network[1])
      @test l < u
    end
    if mcmc.markov_chains[j].events[1].infection[i] > -Inf
      l, u = Pathogen._bounds(i, State_I, ee, obs, mcmc.markov_chains[j].events[1], mcmc.markov_chains[j].transmission_network[1])
      @test l < u
    end
    if mcmc.markov_chains[j].events[1].removal[i] > -Inf
      l, u = Pathogen._bounds(i, State_R, ee, obs, mcmc.markov_chains[j].events[1], mcmc.markov_chains[j].transmission_network[1])
      @test l < u
    end
  end
  iterate!(mcmc, 100, 1.0, condition_on_network=true)
  @test all([mcmc.markov_chains[i].iterations for i = 1:3] .== 100)
  tnd = TNDistribution(mcmc, burnin=2)
  @test sum(tnd) ≈ sum(obs.infection .!== NaN)
  @test sum(mode(tnd)) ≈ sum(obs.infection .!== NaN)
  @test all(Pathogen._outdegree(tnd) .>= 0.0)
  mcmc = MCMC(obs, ee, pop, rf, rpriors, tnprior=tnd)
  start!(mcmc, attempts=100)
  iterate!(mcmc, 100, 1.0, condition_on_network=false)
  @test length(mcmc.markov_chains) == 1
  @test typeof(summary(mcmc, burnin=50, thin=2)) == DataFrame

  tn = mcmc.markov_chains[1].transmission_network[end]
  for i = 1:n
    pt = Pathogen._pathway_to(i, tn, depth=1)
    if length(pt) >= 2
      @test i in Pathogen._pathway_from(pt[2], tn, depth=1)
    end
  end
end

@testset "SEI Model" begin
  rf = RiskFunctions{SEI}(_constant,
                          _coefficient,
                          _powerlaw,
                          _one,
                          _constant)

  rparams = RiskParameters{SEI}([0.01],
                                [1.0],
                                [2.0, 5.0],
                                Float64[],
                                [0.1])

  sim = Simulation(pop,
                   rf,
                   rparams)

  simulate!(sim, tmax=100.0)

  state_counts = [Pathogen._count_by_state(sim.events, State_S, 50.0)
                  Pathogen._count_by_state(sim.events, State_E, 50.0)
                  Pathogen._count_by_state(sim.events, State_I, 50.0)]

  @test sum(state_counts) == n

  @test length(sim.disease_states) == n
  @test size(sim.transmission_network.internal) == (n, n)

  obs = observe(sim, Uniform())
  @test length(obs.infection) == n

  rpriors = RiskPriors{SEI}([Uniform(0.0, 0.1)],
                            [Uniform(0.0, 2.0)],
                            [Uniform(0.0, 4.0); Uniform(1.0, 8.0)],
                            UnivariateDistribution[],
                            [Uniform(0.0, 1.0)])

  ee = EventExtents{SEI}((2.0, 20.0), 2.0)
  mcmc = MCMC(obs, ee, pop, rf, rpriors)
  start!(mcmc, 3, attempts = 100)
  @test length(mcmc.markov_chains) == 3
  @test all([length(mcmc.markov_chains[i].risk_parameters[1]) for i=1:3] .== length(rparams))

  # Check bounds of Events initialization
  for j = 1:3, i = 1:n
    if mcmc.markov_chains[j].events[1].exposure[i] > -Inf
      l, u = Pathogen._bounds(i, State_E, ee, obs, mcmc.markov_chains[j].events[1], mcmc.markov_chains[j].transmission_network[1])
      @test l < u
    end
    if mcmc.markov_chains[j].events[1].infection[i] > -Inf
      l, u = Pathogen._bounds(i, State_I, ee, obs, mcmc.markov_chains[j].events[1], mcmc.markov_chains[j].transmission_network[1])
      @test l < u
    end
  end
  iterate!(mcmc, 100, 1.0, condition_on_network=true)
  @test all([mcmc.markov_chains[i].iterations for i = 1:3] .== 100)
  tnd = TNDistribution(mcmc, burnin=2)
  @test sum(tnd) ≈ sum(obs.infection .!== NaN)
  @test sum(mode(tnd)) ≈ sum(obs.infection .!== NaN)
  @test all(Pathogen._outdegree(tnd) .>= 0.0)
  mcmc = MCMC(obs, ee, pop, rf, rpriors, tnprior=tnd)
  start!(mcmc, attempts=50)
  iterate!(mcmc, 100, 1.0, condition_on_network=false)
  @test length(mcmc.markov_chains) == 1
  @test typeof(summary(mcmc, burnin=50, thin=2)) == DataFrame

  tn = mcmc.markov_chains[1].transmission_network[end]
  for i = 1:n
    pt = Pathogen._pathway_to(i, tn, depth=1)
    if length(pt) >= 2
      @test i in Pathogen._pathway_from(pt[2], tn, depth=1)
    end
  end
end

@testset "SIR Model" begin
  rf = RiskFunctions{SIR}(_constant,
                          _coefficient,
                          _powerlaw,
                          _one,
                          _constant)

  rparams = RiskParameters{SIR}([0.01],
                                [1.0],
                                [2.0, 5.0],
                                Float64[],
                                [0.05])

  sim = Simulation(pop,
                   rf,
                   rparams)

  simulate!(sim, tmax=100.0)

  state_counts = [Pathogen._count_by_state(sim.events, State_S, 50.0)
                  Pathogen._count_by_state(sim.events, State_I, 50.0)
                  Pathogen._count_by_state(sim.events, State_R, 50.0)]

  @test sum(state_counts) == n

  @test length(sim.disease_states) == n
  @test size(sim.transmission_network.internal) == (n, n)

  obs = observe(sim, Uniform(), Uniform())
  @test length(obs.infection) == n

  rpriors = RiskPriors{SIR}([Uniform(0.0, 0.1)],
                            [Uniform(0.0, 2.0)],
                            [Uniform(0.0, 4.0); Uniform(1.0, 8.0)],
                            UnivariateDistribution[],
                            [Uniform(0.0, 1.0)])

  ee = EventExtents{SIR}(2.0, 2.0)
  mcmc = MCMC(obs, ee, pop, rf, rpriors)
  start!(mcmc, 3, attempts=100)
  @test length(mcmc.markov_chains) == 3
  @test all([length(mcmc.markov_chains[i].risk_parameters[1]) for i=1:3] .== length(rparams))

  # Check bounds of Events initialization
  for j = 1:3, i = 1:n
    if mcmc.markov_chains[j].events[1].infection[i] > -Inf
      l, u = Pathogen._bounds(i, State_I, ee, obs, mcmc.markov_chains[j].events[1], mcmc.markov_chains[j].transmission_network[1])
      @test l < u
    end
    if mcmc.markov_chains[j].events[1].removal[i] > -Inf
      l, u = Pathogen._bounds(i, State_R, ee, obs, mcmc.markov_chains[j].events[1], mcmc.markov_chains[j].transmission_network[1])
      @test l < u
    end
  end
  iterate!(mcmc, 100, 1.0, condition_on_network=true)
  @test all([mcmc.markov_chains[i].iterations for i = 1:3] .== 100)
  tnd = TNDistribution(mcmc, burnin=2)
  @test sum(tnd) ≈ sum(obs.infection .!== NaN)
  @test sum(mode(tnd)) ≈ sum(obs.infection .!== NaN)
  @test all(Pathogen._outdegree(tnd) .>= 0.0)
  mcmc = MCMC(obs, ee, pop, rf, rpriors, tnprior=tnd)
  start!(mcmc, attempts=50)
  iterate!(mcmc, 100, 1.0, condition_on_network=false)
  @test length(mcmc.markov_chains) == 1
  @test typeof(summary(mcmc, burnin=50, thin=2)) == DataFrame

  tn = mcmc.markov_chains[1].transmission_network[end]
  for i = 1:n
    pt = Pathogen._pathway_to(i, tn, depth=1)
    if length(pt) >= 2
      @test i in Pathogen._pathway_from(pt[2], tn, depth=1)
    end
  end
end

@testset "SI Model" begin
  rf = RiskFunctions{SI}(_constant,
                         _coefficient,
                         _powerlaw,
                         _one)

  rparams = RiskParameters{SI}([0.01],
                               [1.0],
                               [2.0, 5.0],
                               Float64[])

  sim = Simulation(pop,
                   rf,
                   rparams)

  simulate!(sim, tmax=100.0)

  state_counts = [Pathogen._count_by_state(sim.events, State_S, 50.0)
                  Pathogen._count_by_state(sim.events, State_I, 50.0)]

  @test sum(state_counts) == n

  @test length(sim.disease_states) == n
  @test size(sim.transmission_network.internal) == (n, n)

  obs = observe(sim, Uniform())
  @test length(obs.infection) == n

  rpriors = RiskPriors{SI}([Uniform(0.0, 0.1)],
                           [Uniform(0.0, 2.0)],
                           [Uniform(0.0, 4.0); Uniform(1.0, 8.0)],
                           UnivariateDistribution[])

  ee = EventExtents{SI}(2.0)
  mcmc = MCMC(obs, ee, pop, rf, rpriors)
  start!(mcmc, 3, attempts=100)
  @test length(mcmc.markov_chains) == 3
  @test all([length(mcmc.markov_chains[i].risk_parameters[1]) for i=1:3] .== length(rparams))
  # Check bounds of Events initialization
  for j = 1:3, i = 1:n
    if mcmc.markov_chains[j].events[1].infection[i] > -Inf
      l, u = Pathogen._bounds(i, State_I, ee, obs, mcmc.markov_chains[j].events[1], mcmc.markov_chains[j].transmission_network[1])
      @test l < u
    end
  end
  iterate!(mcmc, 100, 1.0, condition_on_network=true)
  @test all([mcmc.markov_chains[i].iterations for i = 1:3] .== 100)
  tnd = TNDistribution(mcmc, burnin=2)
  @test sum(tnd) ≈ sum(obs.infection .!== NaN)
  @test sum(mode(tnd)) ≈ sum(obs.infection .!== NaN)
  @test all(Pathogen._outdegree(tnd) .>= 0.0)
  mcmc = MCMC(obs, ee, pop, rf, rpriors, tnprior=tnd)
  start!(mcmc, attempts=50)
  iterate!(mcmc, 100, 1.0, condition_on_network=false)
  @test length(mcmc.markov_chains) == 1
  @test typeof(summary(mcmc, burnin=50, thin=2)) == DataFrame

  tn = mcmc.markov_chains[1].transmission_network[end]
  for i = 1:n
    pt = Pathogen._pathway_to(i, tn, depth=1)
    if length(pt) >= 2
      @test i in Pathogen._pathway_from(pt[2], tn, depth=1)
    end
  end
end


text:
using StanSample, Random

Random.seed!(1233)

ProjDir = @__DIR__
cd(ProjDir)

include(joinpath(@__DIR__, "LBA_functions.jl"))

LBA_stan = "
functions{

     real lba_pdf(real t, real b, real A, real v, real s){
          //PDF of the LBA model

          real b_A_tv_ts;
          real b_tv_ts;
          real term_1;
          real term_2;
          real term_3;
          real term_4;
          real pdf;

          b_A_tv_ts = (b - A - t*v)/(t*s);
          b_tv_ts = (b - t*v)/(t*s);
          term_1 = v*Phi(b_A_tv_ts);
          term_2 = s*exp(normal_lpdf(b_A_tv_ts|0,1));
          term_3 = v*Phi(b_tv_ts);
          term_4 = s*exp(normal_lpdf(b_tv_ts|0,1));
          pdf = (1/A)*(-term_1 + term_2 + term_3 - term_4);

          return pdf;
     }

     real lba_cdf(real t, real b, real A, real v, real s){
          //CDF of the LBA model

          real b_A_tv;
          real b_tv;
          real ts;
          real term_1;
          real term_2;
          real term_3;
          real term_4;
          real cdf;

          b_A_tv = b - A - t*v;
          b_tv = b - t*v;
          ts = t*s;
          term_1 = b_A_tv/A * Phi(b_A_tv/ts);
          term_2 = b_tv/A   * Phi(b_tv/ts);
          term_3 = ts/A     * exp(normal_lpdf(b_A_tv/ts|0,1));
          term_4 = ts/A     * exp(normal_lpdf(b_tv/ts|0,1));
          cdf = 1 + term_1 - term_2 + term_3 - term_4;

          return cdf;

     }

     real lba_lpdf(matrix RT, real k, real A, vector v, real s, real tau){

          real t;
          real b;
          real cdf;
          real pdf;
          vector[rows(RT)] prob;
          real out;
          real prob_neg;

          b = A + k;
          for (i in 1:rows(RT)){
               t = RT[i,1] - tau;
               if(t > 0){
                    cdf = 1;

                    for(j in 1:num_elements(v)){
                         if(RT[i,2] == j){
                              pdf = lba_pdf(t, b, A, v[j], s);
                         }else{
                              cdf = (1-lba_cdf(t, b, A, v[j], s)) * cdf;
                         }
                    }
                    prob_neg = 1;
                    for(j in 1:num_elements(v)){
                         prob_neg = Phi(-v[j]/s) * prob_neg;
                    }
                    prob[i] = pdf*cdf;
                    prob[i] = prob[i]/(1-prob_neg);
                    if(prob[i] < 1e-10){
                         prob[i] = 1e-10;
                    }

               }else{
                    prob[i] = 1e-10;
               }
          }
          out = sum(log(prob));
          return out;
     }
}

data{
     int N;
     int Nc;
     vector[N] rt;
     vector[N] choice;
}

parameters {
     real<lower=0> k;
     real<lower=0> A;
     real<lower=0> tau;
     vector<lower=0>[Nc] v;
}

model {
     real s;
     matrix[N,2] RT;
     s=1;
     RT[:,1] = rt;
     RT[:,2] = choice;
     k ~ normal(.5,1)T[0,];
     A ~ normal(.5,1)T[0,];
     tau ~ normal(.5,.5)T[0,];
     for(n in 1:Nc){
          v[n] ~ normal(2,1)T[0,];
     }
     RT ~ lba(k,A,v,s,tau);
}
";

# First 2 runs are using the standard 2.19.1 version of cmdstan

# This run tests passing a data file name as data in the stan_sample() call

stanmodel = SampleModel("LBA", LBA_stan;
  method = StanSample.Sample(adapt = StanSample.Adapt(delta = 0.95)));

N = 200
v = [1.0, 1.5]
Nc = length(v)
data = simulateLBA(;Nd=N,v=v,A=.8,k=.2,tau=.4)

@time (sample_file, log_file) = stan_sample(stanmodel; data=data, n_chains=4)

if !(sample_file == Nothing)
  # Use StanSamples to read a chain in NamedTupla format
  nt = read_samples(stanmodel.sm; chain = 3)

  # Convert to an MCMCChains.Chains object
  chns = read_samples(stanmodel)

  # Describe the MCMCChains using MCMCChains statistics
  display(describe(chns))

  # Show the same output in DataFrame format
  #sdf = StanSample.read_summary(stanmodel)
  #display(sdf)
end


text:
# Autogenerated wrapper script for oxigraph_server_jll for powerpc64le-linux-gnu
export oxigraph_server

using OpenSSL_jll
JLLWrappers.@generate_wrapper_header("oxigraph_server")
JLLWrappers.@declare_executable_product(oxigraph_server)
function __init__()
    JLLWrappers.@generate_init_header(OpenSSL_jll)
    JLLWrappers.@init_executable_product(
        oxigraph_server,
        "bin/oxigraph_server",
    )

    JLLWrappers.@generate_init_footer()
end  # __init__()


text:
"Places refers to a subset of the axes (aka dimension, column, attributes) of relation"
struct Places{T <: Tuple} end
const Place{I} = Places{Tuple{I}}

# First input is INPUTOFFSETth place
const INPUTOFFSET = 2

places(ts) = Places{Tuple{ts...}}
inputplaces(ids) = places(INPUTOFFSET .+ ids)
# inputplaces(s::Statement) = inputplaces(1:length(1:s.expr.args))

const Θ = Place{1}
const Z = Place{2}
const A = Place{3}
const B = Place{4}
const C = Place{5}
const D = Place{5}

Base.:*(p1::Type{Places{Tuple{I}}}, p2::Type{Places{Tuple{J}}}) where {I, J} = Places{Tuple{I, J}}
Base.:*(p1::Type{Places{Tuple{I, J}}}, p2::Type{Places{Tuple{K}}}) where {I, J, K} = Places{Tuple{I, J, K}}
Base.:*(p1::Type{Places{Tuple{I}}}, p2::Type{Places{Tuple{J, K}}}) where {I, J, K} = Places{Tuple{I, J, K}}

const DIMSTOSYM = [:Θ, :Z, Symbol.(collect('A':'Y'))...]

h(I) = I in 1:length(DIMSTOSYM) ? DIMSTOSYM[I] : I

Base.show(io::IO, ::Type{Places{Tuple{I}}}) where {I} = print(io, "Places{$(h(I))}")
Base.show(io::IO, ::Type{Places{Tuple{I, J}}}) where {I, J} = print(io, "Places{$(h(I)), $(h(J))}")
Base.show(io::IO, ::Type{Places{Tuple{I, J, K}}}) where {I, J, K} = print(io, "Places{$(h(I)), $(h(J)), $(h(K))}")
Base.show(io::IO, ::Type{Places{Tuple{I, J, K, L}}}) where {I, J, K, L} = print(io, "Places{$(h(I)), $(h(J)), $(h(K)), $(h(L))}")
Base.show(io::IO, ::Type{Places{T}}) where {T<:Tuple} = print(io, "Places{$(join([h(i) for i in 1:length(T.types)], ", "))}")

const ZA = Z * A
const ZB = Z * B
const ZC = Z * C
const AB = A * B
const BC = B * C
const ABC = A * B * C
const ZΘ = Z * Θ

# Produces places correpsoinding to argument types
places_from_argtypes(t::Type{Tuple{T1}}) where {T1} = A
places_from_argtypes(t::Type{Tuple{T1, T2}}) where {T1, T2} = A * B
places_from_argtypes(t::Type{Tuple{T1, T2, T3}}) where {T1, T2, T3} = A * B * C
places_from_argtypes(t::Type{<:Tuple}) = inputplaces(1:length(t.types))


text:
function printdict(dict, sequence=false)
    # Find longest key
    longestkey = " "
    maxchar = length(longestkey)
    for (key, value) in dict
        char = length(key)
        if char > maxchar
            maxchar = char
            longestkey = key
        end
    end
    # Print dict in a nice format
    println("\n")
    if sequence
        for i in 1:length(dict)
            key = string(i)
            value = dict[key]
            if key != longestkey
                spacediff = maxchar - length(key)
                space = " "^spacediff
                key = "    "*space*key
            else
                key = "    "*key
            end
            println(key, ": ", value)
        end
    else
        for (key, value) in dict
            if key != longestkey
                spacediff = maxchar - length(key)
                space = " "^spacediff
                key = "    "*space*key
            else
                key = "    "*key
            end
            println(key, ": ", value)
        end
    end
end


text:
module BenchBasic
using BenchmarkTools
using BudAutomata: LevenshteinAutomaton, build_dfa

suite = BenchmarkGroup()

lev = LevenshteinAutomaton("woof", 1)
suite["woof_1"] = @benchmarkable build_dfa($lev)

lev = LevenshteinAutomaton("woof", 2)
suite["woof_2"] = @benchmarkable build_dfa($lev)

lev = LevenshteinAutomaton("banana", 1)
suite["banana_1"] = @benchmarkable build_dfa($lev)

end # module

BenchBasic.suite


text:
"""
    struct Turn{T} <: Segment{T}
        axis::Tensors.Vec{3, Float64}
        α::typeof(1.0°)
        r::T
        p0::T
        α0::typeof(1.0°)
    end
A circular turn is parameterized by the turn angle `α` and turning radius `r`.
It begins initial angle `α0` around the axis `axis`.
The center of the circle is given by:
`cen = p0 + Point(r*cos(α0+sign(α)*π/2), r*sin(α0+sign(α)*π/2))`
The parametric function over `t ∈ [0,1]` describing the turn is given by:
`t -> cen + Point(r*cos(α0-sign(α)*π/2+α*t), r*sin(α0-sign(α)*π/2+α*t))`
"""
struct Turn{T} <: Segment{T}
    α::typeof(1.0°)
    r::T
    p0::T
    axis::Tensors.Vec{3, Float64}
    α0::typeof(1.0°)
end

"""
    Turn{T<:Unitful.Length}(α, r::T, p0::T,
        axis=Tensors.Vec{3}((0., 0., 1.)), α0=0.0°)
Outer constructor for `Turn` segments.
"""
Turn(α, r::T;
     p0::T=zero(T),
     axis=Tensors.Vec{3}((0., 0., 1.)),
     α0=0.0°) where {T <: Unitful.Length} =
    Turn{T}(α, r, p0, axis, α0)

pathlength(s::Turn{T}) where {T} = T(abs(s.r*s.α))
function orientation(s::Turn, z)
    axis = s.axis
    angle = s.α0 + (z-s.p0)/s.r
    return (axis, angle)
end

"""
    add_turn!(w::Waveguide{T}, α, r::Unitful.Length,
        med::Medium=lastmedium(w)) where {T <: Unitful.Length}
Turn a waveguide `w` by angle `α` with turning radius `r` in the
current direction. By default, we take the last continuous medium
in the waveguide
"""
function add_turn!(w::Waveguide{T}, α, r::Unitful.Length,
        med::Medium=lastmedium(w)) where {T <: Unitful.Length}
    Unitful.dimension(T) != Unitful.dimension(typeof(r)) && throw(Unitful.DimensionError(T(1),l))
    @assert r > zero(r) "tried to turn by infinitesimal amount."
    p0 = pathlength(w)
    axis, α0 = orientation(w, p0)
    s = Turn{T}(α, r, p0, axis, α0)
    push!(w.structures, Structure(s, med))
    return nothing
end


text:
@enum(CipherID,
    CIPHER_ID_NONE=0,
    CIPHER_ID_NULL,
    CIPHER_AES,
    CIPHER_DES,
    CIPHER_3DES,
    CIPHER_CAMELLIA,
    CIPHER_BLOWFISH,
    CIPHER_ARC4)

@enum(CipherKind,
    CIPHER_NONE = 0,
    CIPHER_NULL,
    CIPHER_AES_128_ECB,
    CIPHER_AES_192_ECB,
    CIPHER_AES_256_ECB,
    CIPHER_AES_128_CBC,
    CIPHER_AES_192_CBC,
    CIPHER_AES_256_CBC,
    CIPHER_AES_128_CFB128,
    CIPHER_AES_192_CFB128,
    CIPHER_AES_256_CFB128,
    CIPHER_AES_128_CTR,
    CIPHER_AES_192_CTR,
    CIPHER_AES_256_CTR,
    CIPHER_AES_128_GCM,
    CIPHER_AES_192_GCM,
    CIPHER_AES_256_GCM,
    CIPHER_CAMELLIA_128_ECB,
    CIPHER_CAMELLIA_192_ECB,
    CIPHER_CAMELLIA_256_ECB,
    CIPHER_CAMELLIA_128_CBC,
    CIPHER_CAMELLIA_192_CBC,
    CIPHER_CAMELLIA_256_CBC,
    CIPHER_CAMELLIA_128_CFB128,
    CIPHER_CAMELLIA_192_CFB128,
    CIPHER_CAMELLIA_256_CFB128,
    CIPHER_CAMELLIA_128_CTR,
    CIPHER_CAMELLIA_192_CTR,
    CIPHER_CAMELLIA_256_CTR,
    CIPHER_CAMELLIA_128_GCM,
    CIPHER_CAMELLIA_192_GCM,
    CIPHER_CAMELLIA_256_GCM,
    CIPHER_DES_ECB,
    CIPHER_DES_CBC,
    CIPHER_DES_EDE_ECB,
    CIPHER_DES_EDE_CBC,
    CIPHER_DES_EDE3_ECB,
    CIPHER_DES_EDE3_CBC,
    CIPHER_BLOWFISH_ECB,
    CIPHER_BLOWFISH_CBC,
    CIPHER_BLOWFISH_CFB64,
    CIPHER_BLOWFISH_CTR,
    CIPHER_ARC4_128,
    CIPHER_AES_128_CCM,
    CIPHER_AES_192_CCM,
    CIPHER_AES_256_CCM,
    CIPHER_CAMELLIA_128_CCM,
    CIPHER_CAMELLIA_192_CCM,
    CIPHER_CAMELLIA_256_CCM)

@enum(CipherMode,
    CIPHER_MODE_NONE = 0,
    CIPHER_MODE_ECB,
    CIPHER_MODE_CBC,
    CIPHER_MODE_CFB,
    CIPHER_MODE_OFB,
    CIPHER_MODE_CTR,
    CIPHER_MODE_GCM,
    CIPHER_MODE_STREAM,
    CIPHER_MODE_CCM)

@enum(Padding,
    PADDING_PKCS7 = 0,
    PADDING_ONE_AND_ZEROS,
    PADDING_ZEROS_AND_LEN,
    PADDING_ZEROS,
    PADDING_NONE)

@enum(Operation,
    OPERATION_NONE = -1,
    DECRYPT = 0,
    ENCRYPT)

type CipherInfo
    data::Ptr{Void}
end

type Cipher
    data::Ptr{Void}

    function Cipher()
        ctx = new()
        ctx.data = Libc.malloc(200) # 88
        ccall((:mbedtls_cipher_init, MBED_CRYPTO), Void,
            (Ptr{Void},), ctx.data)

        finalizer(ctx, ctx->begin
            ccall((:mbedtls_cipher_free, MBED_CRYPTO), Void,
                (Ptr{Void},), ctx.data)
            Libc.free(ctx.data)
        end)

        ctx
    end
end

function CipherInfo(name::AbstractString)
    ptr = ccall((:mbedtls_cipher_info_from_string, MBED_CRYPTO), Ptr{Void},
        (Cstring,), String(name))
    ptr == C_NULL && error("No cipher for $name found")
    CipherInfo(ptr)
end

function CipherInfo(kind::CipherKind)
    ptr = ccall((:mbedtls_cipher_info_from_type, MBED_CRYPTO), Ptr{Void},
        (Cint,), Int(kind))
    ptr == C_NULL && error("No cipher for $kind found")
    CipherInfo(ptr)
end

"""
`CipherInfo(id::CipherID, key_bitlength, mode::CipherMode) -> CipherInfo`

Construct a custom cipher info object.

- `id`: A CipherID, such as CIPHER_AES of CIPHER_BLOWFIHS
- `key_bitlength`: The bit length of the secret key. Available options depend
on the specific cipher id.
- `mode`: Either `Encrypt` or `Decrypt` to indicate which cipher operation will
be performed with this cipher info.
"""
function CipherInfo(id::CipherID, key_bitlen, mode::CipherMode)
    ptr = ccall((:mbedtls_cipher_info_from_values, MBED_CRYPTO), Ptr{Void},
        (Cint, Cint, Cint), Int(id), key_bitlen, Int(mode))
    ptr == C_NULL && error("No cipher for ($id, $(key_bitlen), $mode) found")
    CipherInfo(ptr)
end

"""
`CipherInfo(id::CipherID) -> CipherInfo`

Build a CipherInfo for the given cipher id using the strongest available
key size and the CBC block mode.

*Warning*: In CBC block mode, it is imperative that you use a unique IV (initial value)
for each encryption operation to maintain security.
"""
function CipherInfo(id::CipherID)
    canonical_map = Dict(
        CIPHER_AES => CIPHER_AES_256_CBC,
        CIPHER_DES => CIPHER_DES_CBC,
        CIPHER_CAMELLIA => CIPHER_CAMELLIA_256_CBC,
        CIPHER_BLOWFISH => CIPHER_BLOWFISH_CBC,
        CIPHER_ARC4 => CIPHER_ARC4_128
    )
    if haskey(canonical_map, id)
        canonical_map[id]
    else
        error("No default cipher found for $id. Use the three-argument constructor
              of CipherInfo to explicitly create a cipher info object.")
    end
end

function Cipher(info::CipherInfo)
    cipher = Cipher()
    @err_check ccall((:mbedtls_cipher_setup, MBED_CRYPTO), Cint,
        (Ptr{Void}, Ptr{Void}), cipher.data, info.data)
    cipher
end

"""
`Cipher(info::Union{CipherID, CipherKind}) -> Cipher`

Construct a cipher object and set it to use the specified cipher algorithm.

The algorithm can either be specific (ie, `CIPHER_AES_256_CBC`), or general
(ie, `CIPHER_AES`). In the latter case, a default choice of specific cipher will
be used. See `?CipherInfo` for more details.
"""
Cipher(info::Union{CipherID, CipherKind}) = Cipher(CipherInfo(info))

function Base.show(io::IO, cipher::Cipher)
    print(io, "Cipher($(cipher.data)")
end

function get_key_bitlen(cipher::Cipher)
    ret = ccall((:mbedtls_cipher_get_key_bitlen, MBED_CRYPTO), Cint,
        (Ptr{Void},), cipher.data)
    Int(ret)
end

tobytes(x::Vector{UInt8}) = x
tobytes(x) = String(x).data

function set_key!(cipher::Cipher, key, op::Operation)
    key_b = tobytes(key)
    keysize = 8sizeof(key_b)  # Convert key size from bytes to bits
    @err_check ccall((:mbedtls_cipher_setkey, MBED_CRYPTO), Cint,
        (Ptr{Void}, Ptr{Void}, Cint, Cint),
        cipher.data, pointer(key_b), keysize, Int(op))
    key
end

function set_padding_mode!(cipher::Cipher, padding::Padding)
    @err_check ccall((:mbedtls_cipher_set_padding_mode, MBED_CRYPTO), Cint,
        (Ptr{Void}, Cint), cipher.data, Int(padding))
end

function set_iv!(cipher::Cipher, iv)
    iv_b = tobytes(iv)
    @err_check ccall((:mbedtls_cipher_set_iv, MBED_CRYPTO), Cint,
        (Ptr{Void}, Ptr{Void}, Csize_t),
        cipher.data, pointer(iv_b), sizeof(iv_b))
end

"""
`update!(cipher::Cipher, buf_in, buf_out::Vector{UInt8}) -> Int`

Run the given cipher on `buf_in` (a String or `Vector{UInt8}`) and
store the result of the cipher in `buf_out` (a `Vector{UInt8}`).

It is your responsibility to ensure that `buf_out` is at least as large as necessary
to hold the result. It should be at least as big as the size of `buf_in` plus the
block size associated with `cipher`.
"""
function update!(cipher::Cipher, buf_in, buf_out)
    buf_in_b = tobytes(buf_in)
    out_ref = Ref{Csize_t}(sizeof(buf_out))
    @err_check ccall((:mbedtls_cipher_update, MBED_CRYPTO), Cint,
        (Ptr{Void}, Ptr{Void}, Csize_t, Ptr{Void}, Ptr{Csize_t}),
        cipher.data, buf_in, sizeof(buf_in), buf_out, out_ref)
    Int(out_ref[])
end

function finish!(cipher::Cipher, buf_out)
    out_ref = Ref{Csize_t}(sizeof(buf_out))
    @err_check ccall((:mbedtls_cipher_finish, MBED_CRYPTO), Cint,
        (Ptr{Void}, Ptr{Void}, Ptr{Csize_t}),
        cipher.data, buf_out, out_ref)
    Int(out_ref[])
end

function process_iv(iv, cipher)
    if isempty(iv)
        process_iv(Nullable())
    else
        iv_b = tobytes(iv)
        iv_b, sizeof(iv_b)
    end
end

function process_iv(iv::Nullable, cipher)
    if isnull(iv)  # Return a default IV
        # todo: Don't hard-code a block size (this assumes 128-bit, as for AES)
        # todo: Think about what appropriate default (if any) should be used here
        zeros(Int8, 16), 16
    else
        process_iv(get(iv))
    end
end

function crypt!(cipher::Cipher, iv, buf_in, buf_out)
    olen_ref = Ref{Csize_t}(sizeof(buf_out))
    iv_b, iv_size = process_iv(iv, cipher)
    buf_in_b = tobytes(buf_in)
    @err_check ccall((:mbedtls_cipher_crypt, MBED_CRYPTO), Cint,
        (Ptr{Void}, Ptr{Void}, Csize_t, Ptr{Void}, Csize_t, Ptr{Void}, Ptr{Csize_t}),
        cipher.data, iv_b, iv_size, buf_in_b, sizeof(buf_in_b),
        buf_out, olen_ref)
    Int(olen_ref[])
end

crypt!(cipher::Cipher, buf_in, buf_out) = crypt!(cipher, C_NULL, buf_in, buf_out)

function crypt(cipher_info, op::Operation, key, iv, msg)
    const max_block_size = 256  # todo: obtain this dynamically from mbedtls
    buf = tobytes(msg)
    cipher = Cipher(cipher_info)
    set_key!(cipher, key, op)
    buf_out = Array(UInt8, sizeof(buf) + max_block_size)
    olen = crypt!(cipher, iv, buf, buf_out)
    resize!(buf_out, olen)
    buf_out
end

"""
`encrypt(cipher, key, msg, [iv]) -> Vector{UInt8}`

Encrypt a message using the given cipher. The cipher can be specified as

- a generic cipher (like CIPHER_AES)
- a specific cipher (like CIPHER_AES_256_CBC)
- a Cipher object

`key` is the symmetric key used for cryptography, given as either a String or a
`Vector{UInt8}`. It must be the right length for the chosen cipher; for example,
CIPHER_AES_256_CBC requires a 32-byte (256-bit) key.

`msg` is the message to be encoded. It should either be convertible to a String
or be a `Vector{UInt8}`.

`iv` is the initialization vector, whose size must match the block size of the
cipher (eg, 16 bytes for AES). By default, it will be set to all zeros, which
is not secure. For security reasons, it should be set to a different value for each
encryption operation.
"""
encrypt(cipher, key, msg, iv=Nullable()) = crypt(cipher, ENCRYPT, key, iv, msg)

"""
`decrypt(cipher, key, msg, [iv]) -> Vector{UInt8}`

Decrypt a message using the given cipher. The cipher can be specified as

- a generic cipher (like CIPHER_AES)
- a specific cipher (like CIPHER_AES_256_CBC)
- a Cipher object

`key` is the symmetric key used for cryptography, given as either a String or a
`Vector{UInt8}`. It must be the right length for the chosen cipher; for example,
CIPHER_AES_256_CBC requires a 32-byte (256-bit) key.

`msg` is the message to be encoded. It should either be convertible to a String
or be a `Vector{UInt8}`.

`iv` is the initialization vector, whose size must match the block size of the
cipher (eg, 16 bytes for AES) and correspond to the iv used by the encryptor.
By default, it will be set to all zeros.
"""
decrypt(cipher, key, msg, iv=Nullable()) = crypt(cipher, DECRYPT, key, iv, msg)