text:
@testset "Bad Variable" begin
    data = LPWriter.newdatastore()
    @test_throws Exception LPWriter.getvariableindex!(data, "x y")
end

@testset "Bounds" begin
    @testset "Infinity" begin
        @test LPWriter.parsefloat("-inf") == -Inf
        @test LPWriter.parsefloat("-iNf") == -Inf
        @test LPWriter.parsefloat("-iNfinity") == -Inf
        @test LPWriter.parsefloat("+inf") == Inf
        @test LPWriter.parsefloat("+iNf") == Inf
        @test LPWriter.parsefloat("+iNfinity") == Inf
    end
    data = LPWriter.newdatastore()
    for (line, lb, ub) in [
            ("x free", -Inf, Inf),
            ("x <= 1.1", -Inf, 1.1),
            ("x < 1.2", -Inf, 1.2),
            ("x >= 1.3", 1.3, Inf),
            ("x > 1.4", 1.4, Inf),
            ("x = 1.5", 1.5, 1.5),
            ("x == 1.6", 1.6, 1.6),
            ("0 < x < 1.7", 0.0, 1.7),
            ("0 <= x <= 1.7", 0.0, 1.7),
            ("0 < x <= 1.7", 0.0, 1.7),
            ("0 <= x < 1.7", 0.0, 1.7),
            ("1.2 > x > 0.1", 0.1, 1.2),
            ("1.2 >= x >= 0.1", 0.1, 1.2),
            ("1.2 > x >= 0.1", 0.1, 1.2),
            ("1.2 >= x > 0.1", 0.1, 1.2),
            ("5.5 <= x <= +inf", 5.5, Inf)
        ]
        LPWriter.parsesection!(Val{:bounds}, data, line)
        @test data[:collb][1] == lb
        @test data[:colub][1] == ub
    end

    for line in [
            "x not free",
            "x notfree",
            "1.1 > x",
            "> x 1.1",
            "x free < 1",
            "x <> 1.4",
            "x",
            "0 < 1 < x < 1.6",
            "0 < x > 1.6"
        ]
        @test_throws Exception LPWriter.parsesection!(Val{:bounds}, data, line)
    end
end

@testset "Variable type" begin
    for (val, cat) in [(Val{:integer}, :Int), (Val{:binary}, :Bin)]
        @testset "$cat" begin
            data = LPWriter.newdatastore()
            LPWriter.addnewvariable!(data, "A")
            LPWriter.addnewvariable!(data, "B")
            LPWriter.addnewvariable!(data, "C")

            @test all(data[:colcat] .== :Cont)

            LPWriter.parsesection!(val, data, "A")
            @test data[:colcat][1] == cat

            LPWriter.parsesection!(val, data, "B C")
            @test all(data[:colcat]  .== cat)

            LPWriter.parsesection!(val, data, "xy")
            @test all(data[:colcat]  .== cat)
            @test length(data[:colcat]) == 4
        end
    end
end

@testset "Objective" begin
    data = LPWriter.newdatastore()
    LPWriter.parsesection!(Val{:obj}, data, "obj: -1 x + 1 y")
    @test data[:c] == [1.0, -1.0]

    data = LPWriter.newdatastore()
    LPWriter.parsesection!(Val{:obj}, data, "-1 x + 1 y")
    @test data[:c] == [1.0, -1.0]

    data = LPWriter.newdatastore()
    LPWriter.parsesection!(Val{:obj}, data, "+ 1 x + 1 y")
    @test data[:c] == [1.0, 1.0]

    data = LPWriter.newdatastore()
    LPWriter.parsesection!(Val{:obj}, data, "+ 1 x - 2.3 y")
    @test data[:c] == [-2.3, 1.0]

    data = LPWriter.newdatastore()
    LPWriter.parsesection!(Val{:obj}, data, "")
    @test data[:c] == Float64[]

    data = LPWriter.newdatastore()
    @test_throws Exception LPWriter.parsesection!(Val{:obj}, data, "-1 x + x y")

    data = LPWriter.newdatastore()
    @test_throws Exception LPWriter.parsesection!(Val{:obj}, data, "-1 x * 1 y")
end
@testset "Constraints" begin
    data = LPWriter.newdatastore()
    LPWriter.parsesection!(Val{:constraints}, data, "C1: 1.2 x - 2.4 y <= 1")
    @test data[:A].i == [1, 1]
    @test data[:A].j == [1, 2]
    @test data[:A].v == [-2.4, 1.2]
    @test data[:rowlb][1] == -Inf
    @test data[:rowub][1] == 1.0

    data = LPWriter.newdatastore()
    @test_throws Exception LPWriter.parsesection!(Val{:constraints}, data, "C1: a x - 2.4 y <= 1")

    data = LPWriter.newdatastore()
    @test_throws Exception LPWriter.parsesection!(Val{:constraints}, data, "C1: 1 x * 2.4 y <= 1")

    data = LPWriter.newdatastore()
    LPWriter.parsesection!(Val{:constraints}, data, "C1: 1 x + 2.4 y")
    @test_throws Exception LPWriter.parsesection!(Val{:constraints}, data, "C1: 1 x + 2.4 y <= 1")

    data = LPWriter.newdatastore()
    LPWriter.parsesection!(Val{:constraints}, data, "C1: ")
    LPWriter.parsesection!(Val{:constraints}, data, "1 x + 2.4 y")
    LPWriter.parsesection!(Val{:constraints}, data, "== 1")
    @test data[:A].i == [1, 1]
    @test data[:A].j == [1, 2]
    @test data[:A].v == [2.4, 1.0]
    @test data[:rowlb][1] == 1.0
    @test data[:rowub][1] == 1.0

    data = LPWriter.newdatastore()
    LPWriter.parsesection!(Val{:constraints}, data, "C1: -1.2 x - 2.4 y = 1")
    @test data[:A].i == [1, 1]
    @test data[:A].j == [1, 2]
    @test data[:A].v == [-2.4, -1.2]
    @test data[:rowlb][1] == 1.0
    @test data[:rowub][1] == 1.0

end

@testset "SOS" begin
    data = LPWriter.newdatastore()
    LPWriter.parsesos!(data, "c5: S1:: x:1 y:2")
    @test length(data[:sos]) == 1
    @test data[:sos][1] == (1, [1,2], [1.0,2.0])

    @test_throws Exception LPWriter.parsesos!(data, "c5: S1::")
    @test_throws Exception LPWriter.parsesos!(data, "c5: S3:: x:1 y:2")
    @test_throws Exception LPWriter.parsesos!(data, "c5: S2:: x:y:1 y:2")
end

@testset "Model 1" begin
    A, collb, colub, c, rowlb, rowub, sense, colcat, sos, Q, modelname, colnames, rownames = LPWriter.read("model1.lp")

    variable_permutation = [5, 4, 1, 2, 3, 7, 6, 8]

    @test A == permute(sparse([1, 2, 3, 4, 4, 4],[1, 2, 3, 5,6,7], [1, 1, 1, 1, 1, 1], 4, 8), 1:4, variable_permutation)
    @test collb == [-Inf, -Inf, -Inf, 5.5, 0, 0, 0, 0][variable_permutation]
    @test colub == [3, 3, 3, Inf, 1, 1, 1, 1][variable_permutation]
    @test c == [0,0,0,-1,1,0,0,0][variable_permutation]
    @test rowlb == [0, 2, -Inf, -Inf]
    @test rowub == [Inf, Inf, 2.5, 1]
    @test sense == :Max
    @test colcat == [:Cont, :Cont, :Cont, :Int, :Cont, :Cont, :Cont, :Bin][variable_permutation]
    @test sos == LPWriter.SOS[
        (1, [3,5,1], [1.0, 2.0, 3.0]),
        (2, [4,2,1], [2.0, 1.0, 2.5])
    ]
    # @test Q == Array{Float64}(0,0)
    # @test modelname == "TestModel"
    @test colnames == ["V$(i)" for i in variable_permutation]
    @test rownames == ["CON$i" for i in 1:4]
end

@testset "Model 2" begin
    A, collb, colub, c, rowlb, rowub, sense, colcat, sos, Q, modelname, colnames, rownames = LPWriter.read("model2.lp")

    variable_permutation = [5, 4, 1, 2, 3, 7, 6, 8]

    @test A == permute(sparse([1, 2, 3, 4, 4, 4],[1, 2, 3, 5,6,7], [1, 1, 1, 1, 1, 1], 4, 8), 1:4, variable_permutation)
    @test collb == [-Inf, -Inf, -Inf, 5.5, 0, 0, 0, 0][variable_permutation]
    @test colub == [3, 3, 3, Inf, 1, 1, 1, 1][variable_permutation]
    @test c == [0,0,0,-1,1,0,0,0][variable_permutation]
    @test rowlb == [0, 2, -Inf, -Inf]
    @test rowub == [Inf, Inf, 2.5, 1]
    @test sense == :Min
    @test colcat == [:Cont, :Cont, :Cont, :Int, :Cont, :Cont, :Cont, :Bin][variable_permutation]
    @test sos == LPWriter.SOS[]
    # @test Q == Array{Float64}(0,0)
    # @test modelname == "TestModel"
    @test colnames == ["V$(i)" for i in variable_permutation]
    @test rownames == ["CON$i" for i in 1:4]
end

@testset "Tricky" begin
    A, collb, colub, c, rowlb, rowub, sense, colcat, sos, Q, modelname, colnames, rownames = LPWriter.read("model1_tricky.lp")

    variable_permutation = [4, 5, 1, 2, 3, 6, 7, 8]

    @test A == permute(sparse([1, 2, 3, 4, 4, 4],[1, 2, 3, 5,6,7], [1, 1, 1, 1, 1, 1], 4, 8), 1:4, variable_permutation)
    @test collb == [-Inf, -Inf, -3, 5.5, 1, -Inf, 0, 0][variable_permutation]
    @test colub == [3, 3, Inf, Inf, 1, Inf, 1, 1][variable_permutation]
    @test c == [0,0,0,-1,1,0,0,0][variable_permutation]
    @test rowlb == [0, 2, -Inf, -Inf]
    @test rowub == [Inf, Inf, 2.5, 1]
    @test sense == :Max
    @test colcat == [:Cont, :Cont, :Cont, :Int, :Int, :Int, :Cont, :Bin][variable_permutation]
    @test sos == LPWriter.SOS[]
    # @test Q == Array{Float64}(0,0)
    # @test modelname == "TestModel"
    cnames = ["V$(i)" for i in variable_permutation]
    cnames[1] = "Var4"
    @test colnames == cnames
    @test rownames == ["CON1", "R2", "CON3", "CON4"]
end

@testset "Corrrupt" begin
    @test_throws Exception LPWriter.read("corrupt.lp")
end


text:
@doc doc"""
Approximation de la solution du problème ``\min_{x \in \mathbb{R}^{n}} f(x)`` en utilisant l'algorithme de Newton

# Syntaxe
```julia
xk,f_min,flag,nb_iters = Algorithme_de_Newton(f,gradf,hessf,x0,option)
```

# Entrées :
   * **f**       : (Function) la fonction à minimiser
   * **gradf**   : (Function) le gradient de la fonction f
   * **hessf**   : (Function) la Hessienne de la fonction f
   * **x0**      : (Array{Float,1}) première approximation de la solution cherchée
   * **options** : (Array{Float,1})
       * **max_iter**      : le nombre maximal d'iterations
       * **Tol_abs**       : la tolérence absolue
       * **Tol_rel**       : la tolérence relative

# Sorties:
   * **xmin**    : (Array{Float,1}) une approximation de la solution du problème  : ``\min_{x \in \mathbb{R}^{n}} f(x)``
   * **f_min**   : (Float) ``f(x_{min})``
   * **flag**     : (Integer) indique le critère sur lequel le programme à arrêter
      * **0**    : Convergence
      * **1**    : stagnation du xk
      * **2**    : stagnation du f
      * **3**    : nombre maximal d'itération dépassé
   * **nb_iters** : (Integer) le nombre d'itérations faites par le programme

# Exemple d'appel
```@example
using Optinum
f(x)=100*(x[2]-x[1]^2)^2+(1-x[1])^2
gradf(x)=[-400*x[1]*(x[2]-x[1]^2)-2*(1-x[1]) ; 200*(x[2]-x[1]^2)]
hessf(x)=[-400*(x[2]-3*x[1]^2)+2  -400*x[1];-400*x[1]  200]
x0 = [1; 0]
options = []
xmin,f_min,flag,nb_iters = Algorithme_De_Newton(f,gradf,hessf,x0,options)
```
"""
function Algorithme_De_Newton(f::Function,gradf::Function,hessf::Function,x0,options)

    "# Si options == [] on prends les paramètres par défaut"
    if options == []
        max_iter = 100
        Tol_abs = sqrt(eps())
        Tol_rel = 1e-15
    else
        max_iter = options[1]
        Tol_abs = options[2]
        Tol_rel = options[3]
    end

    n = length(x0)
    xmin = x0
    f_min = 0
    flag = 0
    nb_iters = 0
    xk = x0
    xk_1 = x0

    while (norm(gradf(xk)) > max(Tol_rel * norm(gradf(x0)), Tol_abs))
        nb_iters += 1
        xk_1 = xk
        dk = -hessf(xk) \ gradf(xk)
        xk = xk + dk

        if norm(gradf(xk_1)) <= max(Tol_abs, Tol_rel * norm(gradf(x0)))
            flag = 0
            break

        elseif norm(dk) <= max(Tol_abs, Tol_rel * norm(xk))
            flag = 1
            break

        elseif abs(f(xk_1) - f(xk)) <= max(Tol_abs, Tol_rel * abs(f(xk)))
            flag = 2
            break
        elseif nb_iters == max_iter
            flag = 3
            break
        end
    end

    xmin = xk
    f_min = f(xmin)

    return xmin,f_min,flag,nb_iters
end

ext: jl

max_stars_repo_name: MOUDDENEHamza/ENSEEIHT

max_issues_repo_name: MOUDDENEHamza/ENSEEIHT

max_forks_repo_name: MOUDDENEHamza/ENSEEIHT


text:
function process_additional_data!(data)
  mva_base = data["baseMVA"]
  rescale = x -> x/mva_base
  rescale_cost = x -> mva_base * x
  if data["multinetwork"] == false
    if haskey(data, "load")
        for (i, load) in data["load"]
          PowerModels.apply_func(load, "pref", rescale)
          PowerModels.apply_func(load, "qref", rescale)
          PowerModels.apply_func(load, "pmin", rescale)
          PowerModels.apply_func(load, "pmax", rescale)
          PowerModels.apply_func(load, "qmin", rescale)
          PowerModels.apply_func(load, "qmax", rescale)
          PowerModels.apply_func(load, "prated", rescale)
          PowerModels.apply_func(load, "qrated", rescale)
          PowerModels.apply_func(load, "voll", rescale_cost)
        end
      end
      if haskey(data["gen"],"1") && haskey(data["gen"]["1"], "prated")
          for (i, gen) in data["gen"]
            PowerModels.apply_func(gen, "pref", rescale)
            PowerModels.apply_func(gen, "qref", rescale)
            PowerModels.apply_func(gen, "prated", rescale)
            PowerModels.apply_func(gen, "qrated", rescale)
          end
      end
      if haskey(data["branch"], "shiftable")
            for (i, branch) in data["branch"]
                PowerModels.apply_func(branch, "shift_fr", PowerModels.deg2rad)
                PowerModels.apply_func(branch, "shift_to", PowerModels.deg2rad)
                PowerModels.apply_func(branch, "shift_fr_max", PowerModels.deg2rad)
                PowerModels.apply_func(branch, "shift_fr_min", PowerModels.deg2rad)
                PowerModels.apply_func(branch, "shift_to_max", PowerModels.deg2rad)
                PowerModels.apply_func(branch, "shift_to_min", PowerModels.deg2rad)
            end
      end
      if haskey(data, "contingencies")
            for (c, cont) in data["contingencies"]
                data["contingencies"][c]["index"] = parse(Int, c)
            end
      end
  else
    for (n, network) in data["nw"]
        if haskey(data["nw"][n], "load")
          for (i, load) in data["nw"][n]["load"]
            PowerModels.apply_func(load, "pref", rescale)
            PowerModels.apply_func(load, "qref", rescale)
            PowerModels.apply_func(load, "pmin", rescale)
            PowerModels.apply_func(load, "pmax", rescale)
            PowerModels.apply_func(load, "qmin", rescale)
            PowerModels.apply_func(load, "qmax", rescale)
            PowerModels.apply_func(load, "prated", rescale)
            PowerModels.apply_func(load, "qrated", rescale)
            PowerModels.apply_func(load, "voll", rescale_cost)
          end
        end
        if haskey(data["nw"][n]["gen"],"1") && haskey(data["nw"][n]["gen"]["1"], "prated")
            for (i, gen) in data["nw"][n]["gen"]
              PowerModels.apply_func(gen, "pref", rescale)
              PowerModels.apply_func(gen, "qref", rescale)
              PowerModels.apply_func(gen, "prated", rescale)
              PowerModels.apply_func(gen, "qrated", rescale)
            end
        end
        if haskey(data["nw"][n]["branch"], "shiftable")
              for (i, branch) in data["nw"][n]["branch"]
                  PowerModels.apply_func(branch, "shift_fr", PowerModels.deg2rad)
                  PowerModels.apply_func(branch, "shift_to", PowerModels.deg2rad)
                  PowerModels.apply_func(branch, "shift_fr_max", PowerModels.deg2rad)
                  PowerModels.apply_func(branch, "shift_fr_min", PowerModels.deg2rad)
                  PowerModels.apply_func(branch, "shift_to_max", PowerModels.deg2rad)
                  PowerModels.apply_func(branch, "shift_to_min", PowerModels.deg2rad)
              end
        end
        if haskey(data["nw"][n], "contingencies")
              for (c, cont) in data["nw"][n]["contingencies"]
                  data["nw"][n]["contingencies"][c]["index"] = parse(Int, c)
              end
        end
    end
  end
end

max_issues_repo_name: frederikgeth/PowerModelsReliability

max_forks_repo_name: frederikgeth/PowerModelsReliability


text:
"""
    quadratic objective
"""
struct QuadraticObjective <: Objective
    Q
    q
    R
    r
    c
end

function quadratic_tracking_objective(Q, R, x, u)
    T = length(Q)
    q = [-2.0 * Q[t] * x[t] for t = 1:T]
    r = [-2.0 * R[t] * u[t] for t = 1:T-1]
    c = [(t == T  ? x[t]' * Q[t] * x[t]
         : x[t]' * Q[t] * x[t] + u[t]' * R[t] * u[t]) for t = 1:T]

    return QuadraticObjective(Q, q, R, r, c)
end

function objective(Z, obj::QuadraticObjective, idx, T)
    J = 0.0

    for t = 1:T
        x = view(Z, idx.x[t])
        Q = obj.Q[t]
        q = obj.q[t]
        c = obj.c[t]

        J += x' * Q * x + q' * x
        J += c

        t == T && continue

        u = view(Z, idx.u[t])
        R = obj.R[t]
        r = obj.r[t]

        J += u' * R * u + r' * u
    end

    return J
end

function objective_gradient!(∇J, Z, obj::QuadraticObjective, idx, T)

    for t = 1:T
        x = view(Z, idx.x[t])
        Q = obj.Q[t]
        q = obj.q[t]
        c = obj.c[t]

        ∇J[idx.x[t]] += 2.0 * Q * x + q

        t == T && continue

        u = view(Z, idx.u[t])
        R = obj.R[t]
        r = obj.r[t]

        ∇J[idx.u[t]] += 2.0 * R * u + r
    end

    return nothing
end

"""
    quadratic time objective
"""
struct QuadraticTimeObjective{T} <: Objective
    obj::QuadraticObjective
    c::T
end

function quadratic_time_tracking_objective(Q, R, x, u, c)
    return QuadraticTimeObjective(quadratic_tracking_objective(Q, R, x, u),c)
end

function objective(Z, obj::QuadraticTimeObjective, idx, T)
    J = 0.0

    for t = 1:T-1
        x = view(Z, idx.x[t])
        Q = obj.obj.Q[t]
        q = obj.obj.q[t]

        u = view(Z, idx.u[t])
        R = obj.obj.R[t]
        r = obj.obj.r[t]

        h = u[end]

        c = obj.obj.c[t]

        J += (x' * Q * x + q' * x) * h
        J += (u' * R * u + r' * u) * h
        J += (c + obj.c) * h
    end

    x = view(Z, idx.x[T])
    Q = obj.obj.Q[T]
    q = obj.obj.q[T]

    J += (x' * Q * x + q' * x)

    return J
end

function objective_gradient!(∇J, Z, obj::QuadraticTimeObjective, idx, T)

    J = 0.0

    for t = 1:T-1
        x = view(Z, idx.x[t])
        Q = obj.obj.Q[t]
        q = obj.obj.q[t]

        u = view(Z, idx.u[t])
        R = obj.obj.R[t]
        r = obj.obj.r[t]

        h = u[end]

        c = obj.obj.c[t]

        ∇J[idx.x[t]] += (2.0 * Q * x + q) * h
        ∇J[idx.u[t]] += (2.0 * R * u + r) * h
        ∇J[idx.u[t][end]] += x' * Q * x + q' * x + u' * R * u + r' * u + c + obj.c
    end

    x = view(Z, idx.x[T])
    Q = obj.obj.Q[T]
    q = obj.obj.q[T]

    ∇J[idx.x[T]] += (2.0 * Q * x + q)

    return nothing
end

max_issues_repo_name: zacmanchester/DirectMotionPlanning

max_forks_repo_name: zacmanchester/DirectMotionPlanning


text:
using JuMP, EAGO

                     m = Model()

                     EAGO.register_eago_operators!(m)

                     @variable(m, -1 <= x[i=1:3] <= 1)

                     @variable(m, -3.428579430434188 <= q <= 8.991868736935128)

                     add_NL_constraint(m, :(swish(0.758903986598098 + 0.8106379052679875*swish(0.47487042435595583 + 0.4488775132618512*swish(0.1716391672055586 + 0.7301849089692007*swish(0.33202390613731003 + -0.8931954910636541*$(x[1]) + 0.7142724585602056*$(x[2]) + -0.28852147983722487*$(x[3])) + 0.6276077641275091*swish(-0.8935890604396688 + -0.29936605585485676*$(x[1]) + -0.5391102865466948*$(x[2]) + 0.8633503800224367*$(x[3]))) + 0.015398554931854491*swish(0.8366427516313624 + 0.8654420505988725*swish(0.33202390613731003 + -0.8931954910636541*$(x[1]) + 0.7142724585602056*$(x[2]) + -0.28852147983722487*$(x[3])) + -0.07448940361379819*swish(-0.8935890604396688 + -0.29936605585485676*$(x[1]) + -0.5391102865466948*$(x[2]) + 0.8633503800224367*$(x[3])))) + -0.6668215058198093*swish(0.8216722410540171 + 0.5667551555087149*swish(0.1716391672055586 + 0.7301849089692007*swish(0.33202390613731003 + -0.8931954910636541*$(x[1]) + 0.7142724585602056*$(x[2]) + -0.28852147983722487*$(x[3])) + 0.6276077641275091*swish(-0.8935890604396688 + -0.29936605585485676*$(x[1]) + -0.5391102865466948*$(x[2]) + 0.8633503800224367*$(x[3]))) + 0.5912184082121308*swish(0.8366427516313624 + 0.8654420505988725*swish(0.33202390613731003 + -0.8931954910636541*$(x[1]) + 0.7142724585602056*$(x[2]) + -0.28852147983722487*$(x[3])) + -0.07448940361379819*swish(-0.8935890604396688 + -0.29936605585485676*$(x[1]) + -0.5391102865466948*$(x[2]) + 0.8633503800224367*$(x[3]))))) + swish(-0.09143860592972564 + 0.3407854070498746*swish(0.47487042435595583 + 0.4488775132618512*swish(0.1716391672055586 + 0.7301849089692007*swish(0.33202390613731003 + -0.8931954910636541*$(x[1]) + 0.7142724585602056*$(x[2]) + -0.28852147983722487*$(x[3])) + 0.6276077641275091*swish(-0.8935890604396688 + -0.29936605585485676*$(x[1]) + -0.5391102865466948*$(x[2]) + 0.8633503800224367*$(x[3]))) + 0.015398554931854491*swish(0.8366427516313624 + 0.8654420505988725*swish(0.33202390613731003 + -0.8931954910636541*$(x[1]) + 0.7142724585602056*$(x[2]) + -0.28852147983722487*$(x[3])) + -0.07448940361379819*swish(-0.8935890604396688 + -0.29936605585485676*$(x[1]) + -0.5391102865466948*$(x[2]) + 0.8633503800224367*$(x[3])))) + 0.24612879529032217*swish(0.8216722410540171 + 0.5667551555087149*swish(0.1716391672055586 + 0.7301849089692007*swish(0.33202390613731003 + -0.8931954910636541*$(x[1]) + 0.7142724585602056*$(x[2]) + -0.28852147983722487*$(x[3])) + 0.6276077641275091*swish(-0.8935890604396688 + -0.29936605585485676*$(x[1]) + -0.5391102865466948*$(x[2]) + 0.8633503800224367*$(x[3]))) + 0.5912184082121308*swish(0.8366427516313624 + 0.8654420505988725*swish(0.33202390613731003 + -0.8931954910636541*$(x[1]) + 0.7142724585602056*$(x[2]) + -0.28852147983722487*$(x[3])) + -0.07448940361379819*swish(-0.8935890604396688 + -0.29936605585485676*$(x[1]) + -0.5391102865466948*$(x[2]) + 0.8633503800224367*$(x[3]))))) - $q <= 0.0))

                     @objective(m, Min, q)

                     return m

max_stars_repo_name: PSORLab/RSActivationFunctions

max_issues_repo_name: PSORLab/RSActivationFunctions

max_forks_repo_name: PSORLab/RSActivationFunctions


text:
module AoC2021

using BenchmarkTools
using Markdown
using Printf
using Reexport
using Statistics

export benchmark_all, run_all

DAYS_DONE = [1,2,3,4,5,6,7,8,9,10,11,12,13,14]

name(day) = @sprintf("day%02d", day)

for day in DAYS_DONE
    dname = name(day)
    include("$dname.jl")
    @eval using .$(Symbol(uppercasefirst(dname)))
    @eval export $(Symbol(uppercasefirst(dname))), $(Symbol(dname))
end

function benchmark_all(days = DAYS_DONE; showall=false)
    totaltime = 0.0
    for day in days
        println("Benchmarking $day")
        b = @eval @benchmark $(Symbol(name(day)))()
        showall && display(b)
        totaltime += median(b.times)
    end
    println("Total time: ", round(totaltime*1e-6; sigdigits=4), "ms")
    return totaltime*1e-6
end

function run_all(days = DAYS_DONE)
    table = Markdown.Table([["Day", "Time", "Part 1", "Part 2"]], [:c, :c, :c, :c])
    for day in days
        res = @timed @eval $(Symbol(name(day)))()
        row = [ @sprintf("Day %02d", day),
                _format_time(res.time),
                string(res.value[1]),
                string(res.value[2])]
        push!(table.rows, row)
    end
    return Markdown.MD(table)
end

function _format_time(time; sigdigits=4)
    magnitudes = [1e-6, 1e-3, 1.0, 60]
    prefixes = ["ns", "µs", "ms", "s"]
    for (mag, pref) in zip(magnitudes, prefixes)
        if time < mag
            return "$(round(time/mag/1e-3; sigdigits))$pref"
        end
    end
    return "$(div(time,60))min $(round(mod(time, 60); digits=1))s"
end

end # module


text:
#using PlinkPlumber
using Base.Test

# Make sure that plink is installed
@test ismatch(r"PLINK.*",readall(`plink --version`))

function create_ped()
# Example from https://www.cog-genomics.org/plink2/formats
	test_ped = """
  	  1 1 0 0 1  0  G G  2 2  C C
  	  1 2 0 0 1  0  A A  0 0  A C
  	  1 3 1 2 1  2  0 0  1 2  A C
  	  2 1 0 0 1  0  A A  2 2  0 0
  	  2 2 0 0 1  2  A A  2 2  0 0
  	  2 3 1 2 1  2  A A  2 2  A A
	"""

	test_map = """
 	 1 snp1 0 1
 	 1 snp2 0 2
 	 1 snp3 0 3
	"""

	plinkstem=tempname()

	open(string(plinkstem,".ped"), "w") do f
		write(f,test_ped)
	end

	open(string(plinkstem,".map"), "w") do f
		write(f,test_map)
	end

	return plinkstem
end

function create_bed(pedstem)
	plinkstem=tempname()
	run(`plink --file $pedstem --make-bed --out $plinkstem`)
	plinkstem
end

function test_readplink()
	a_ped_file = create_ped()
	a_bed_file = create_bed(a_ped_file)

	pedread = readplink(plink(a_ped_file))
	bedread = readplink(plink(a_bed_file))

	println(pedread)
	println(bedread)
end


text:
## Read matrix by [HS GPA, parental income]
# These are transposed in the raw files. Transpose, but then make into a Matrix.
function read_by_gpa_yp(ds :: DataSettings, fPath :: AbstractString)
	df = read_by_xy(fPath);
	m = Matrix{Float64}(data_matrix(df)');
	return m
end

read_by_gpa_yp(ds :: DataSettings, mName :: Symbol; momentType = nothing) =
	read_by_gpa_yp(ds, raw_file_path(ds, mName; momentType = momentType));

## -------------  Individual moments

function load_entry_gpa_yp(ds :: DataSettings)
	load_fct = mt -> read_by_gpa_yp(ds, :fracEnter_gpM; momentType = mt);
	m, ses, cnts = choice_prob_from_xy(load_fct);
	# m = read_by_gpa_yp(ds, :fracEnter_gpM);
	# cnts = read_by_gpa_yp(ds, :fracEnter_gpM; momentType = :count);
	# ses = ses_from_choice_probs(m, cnts);
	# cnts = round.(Int, cnts);
	@assert all(m .< 1.0)  &&  all(m .> 0.0)
	@assert size(m) == (n_gpa(ds), n_parental(ds))
	@assert all(cnts .> 20)
	return m, ses, cnts
end

## ------  Mass by HS GPA / parental. Sums to 1.

# function dm_mass_gpa_yp()
#     target = :mass_gpM;
#     return DataMoment(target, ModelStatistic(:mass_xyM, :gpaYpS),
#         fn_gpa_yp("mass"), nothing,
#         mass_by_gpa_yp, plot_gpa_yp)
# end

function mass_by_gpa_yp(ds :: DataSettings)
	m = read_by_gpa_yp(ds, :mass_gpM; momentType = MtMean);
	cnts = read_by_gpa_yp(ds, :mass_gpM; momentType = :count);
	ses = ses_from_choice_probs(m, cnts);
	cnts = round.(Int, cnts);
    @assert all(m .< 0.5)  &&  all(m .> 0.0)
	@assert size(m) == (n_gpa(ds), n_parental(ds))
	@check sum(m) ≈ 1.0
	@assert all(cnts .> 100)
	return m, ses, cnts
end

## College graduation rates, by HS GPA, parental income.
# This is conditional on entry b/c that helps the calibration.
# Currently not constructed
# function grad_rate_by_gpa_yp(ds :: DataSettings)
# 	target = :fracGradGpaYp;
# 	fPath = data_file()
# 	m = read_by_gpa_yp(ds, target);
#     @assert all(m .< 1.0)  &&  all(m .> 0.0)
#     @assert size(m) == (n_gpa(ds), n_parental(ds))

# 	# Make NOT conditional on entry
# 	# dev = entry_by_gpa_yp(ds);
# 	# m = m .* dev.dataV;
# 	return Deviation{Double}(name = target, dataV = m, modelV = m,
# 		wtV = 1.0 ./ m,  scalarWt = 6.0 ./ length(m),
# 		shortStr = String(target),
# 		longStr = "Graduation rates by HSgpa/parental background (conditional on entry)",
# 		showPath = "fracGradGpaYp.txt")
# end

# ---------------


text:
abstract type EventStore end

const DEFAULT_MONGO_URI = "mongodb://127.0.0.1:27021,127.0.0.1:27022,127.0.0.1:27023/?replicaSet=rst"
const AUTO_SNAPSHOT = 100

mutable struct MongoEventStore <: EventStore
    client::Union{Mongoc.Client,Nothing}
    uri::String
    dbname::String
    namespace::Union{String,Nothing}
    state::ServiceState
    auto_snapshot::UInt
end

function MongoEventStore(;uri::String=DEFAULT_MONGO_URI, dbname::String="eventstore", namespace::Union{String,Nothing}=nothing, auto_snapshot=AUTO_SNAPSHOT)
    MongoEventStore(nothing, uri, dbname, namespace, stopped, auto_snapshot)
end

function collection_name(store::MongoEventStore, type::String)
    if isnothing(store.namespace)
        type
    else
        `$(store.namespace)_$(type)`
    end
end

events_collection_name(store::MongoEventStore) = collection_name(store, "events")
snapshots_collection_name(store::MongoEventStore) = collection_name(store, "snapshots")
counters_collection_name(store::MongoEventStore) = collection_name(store, "counters")

function run!(store::MongoEventStore)
    if store.state == running
        throw(EventException("The event store is already running"))
    end

    store.state = running

    store.client = Mongoc.Client(store.uri)

    # Initialize required collections
    client = store.client
    database = client[store.dbname]

    events_name = events_collection_name(store)
    snapshots_name = snapshots_collection_name(store)
    counters_name = counters_collection_name(store)

    create_events_indexes_cmd = BSON(
        "createIndexes" => events_name,
        "indexes" => [
            BSON(
                "key" => BSON("id" => 1),
                "name" => "$(events_name)_id",
                "unique" => true,
            ),
            BSON(
                "key" => BSON("aggregate.id" => 1),
                "name" => "$(events_name)_aggregate_id"
            ),
            BSON(
                "key" => BSON("aggregate.id" => 1, "aggregate.revision" => 1),
                "name" => "$(events_name)_aggregate_id_revision",
                "unique" => true,
            ),
            BSON(
                "key" => BSON("metadata.position" => 1),
                "name" => "$(events_name)_position"
            ),
        ]
    )

    create_snapshots_indexes_cmd = BSON(
        "createIndexes" => snapshots_name,
        "indexes" => [
            BSON(
                "key" => BSON("aggregate.id" => 1),
                "name" => "$(snapshots_name)_aggregate_id",
                "unique" => true,
            ),
        ]
    )

    create_counters_indexes_cmd = BSON(
        "createIndexes" => counters_name,
        "indexes" => [
            BSON(
                "key" => BSON("counter_id" => 1),
                "name" => "$(counters_name)_counter_id",
                "unique" => true,
            ),
        ]
    )

    reply = Mongoc.write_command(database, create_events_indexes_cmd)
    @assert reply["ok"] == 1
    reply = Mongoc.write_command(database, create_snapshots_indexes_cmd)
    @assert reply["ok"] == 1
    reply = Mongoc.write_command(database, create_counters_indexes_cmd)
    @assert reply["ok"] == 1

    # Initialize the counter
    @assert !isnothing(create_counter(store, "events"))
    nothing
end

# Drops the collections...intended to be used for testing only
function reset!(store::MongoEventStore)
    client = Mongoc.Client(store.uri)
    database = client[store.dbname]

    events_name = events_collection_name(store)
    snapshots_name = snapshots_collection_name(store)
    counters_name = counters_collection_name(store)

    drop_events_cmd = BSON(
        "drop" => events_name,
        "comment" => "reset event store",
    )

    drop_snapshots_cmd = BSON(
        "drop" => snapshots_name,
        "comment" => "reset event store",
    )

    drop_counters_cmd = BSON(
        "drop" => counters_name,
        "comment" => "reset event store",
    )

    try
        Mongoc.write_command(database, drop_events_cmd)
        Mongoc.write_command(database, drop_snapshots_cmd)
        Mongoc.write_command(database, drop_counters_cmd)
        nothing
    catch _err
        # Ignore
    end

    finalize(client)
    nothing
end

function stop!(store::MongoEventStore)
    if store.state == stopped
        throw(EventException("The event store is not running"))
    end

    finalize(store.client)
    store.client = nothing

    store.state = stopped
end

function create_counter(store::MongoEventStore, counter_id::String)
    client = store.client
    database = client[store.dbname]

    query = BSON(
        "counter_id" => counter_id
    )

    update = BSON(
        "\$setOnInsert" => BSON(
            "seq" => 0,
            "counter_id" => counter_id,
        ),
    )

    command = BSON(
        "findAndModify" => counters_collection_name(store),
        "query" => query,
        "update" => update,
        "new" => true,
        "upsert" => true,
    )

    reply = Mongoc.write_command(
        database,
        command,
    )

    return reply["value"]
end

function next_sequence(store::MongoEventStore, database::Mongoc.DatabaseSession, counter_id::String)
    # This should be done transactionally which is why
    # a database instance is required in addition to
    # the store instance

    query = BSON(
        "counter_id" => counter_id
    )

    update = BSON(
        "\$inc" => BSON("seq" => 1)
    )

    command = BSON(
        "findAndModify" => counters_collection_name(store),
        "query" => query,
        "update" => update,
        "new" => true,
    )

    options = BSON()
    err_ref = Ref{Mongoc.BSONError}()
    ok = Mongoc.mongoc_client_session_append(
        Mongoc.get_session(database).handle,
        options.handle,
        err_ref,
    )

    if !ok
        throw(err_ref[])
    end

    reply = Mongoc.write_command(
        database.database,
        command;
        options,
    )

    doc = reply["value"]
    if isnothing(doc)
        throw(EventException("failed to get next sequence for $counter_id"))
    end

    doc["seq"]
end

function last_sequence(store::MongoEventStore, counter_id::String)
    client = store.client
    database = client[store.dbname]
    collection = database[counters_collection_name(store)]

    query = BSON(
        "counter_id" => counter_id
    )

    doc = Mongoc.find_one(
        collection,
        query,
    )

    if isnothing(doc)
        throw(EventException("failed to get last sequence for $counter_id"))
    end

    doc["seq"]
end

function aggregate_events_iterator(store::MongoEventStore, aggregate_id::UUID; from_revision=1, to_revision=2 ^ 31 - 1)
    client = store.client
    database = client[store.dbname]
    collection = database[events_collection_name(store)]

    query = BSON(
        "aggregate.id" => aggregate_id,
        "aggregate.revision" => BSON(
            "\$gte" => from_revision,
            "\$lte" => to_revision,
        ),
    )

    options = BSON(
        "sort" => BSON("aggregate.revision" => 1),
        "projection" => BSON("_id" => 0),
    )

    cursor = Mongoc.find(collection, query; options)

    EventsIterator(cursor)
end

function unpublished_events_iterator(store::MongoEventStore)
    client = store.client
    database = client[store.dbname]
    collection = database[events_collection_name(store)]

    query = BSON(
        "metadata.published" => false,
    )

    options = BSON(
        "sort" => BSON("metadata.position" => 1),
        "projection" => BSON("_id" => 0),
    )

    cursor = Mongoc.find(collection, query; options)

    EventsIterator(cursor)
end

function save_events(store::MongoEventStore, new_events::Vector{Event})
    committed_events = Vector{Event}()
    client = store.client

    # Run transactionally
    Mongoc.transaction(client) do session
        database = session[store.dbname]
        collection = database[events_collection_name(store)]

        for event in new_events
            if event.aggregate.revision < 1
                throw(ArgumentError("aggregate revision must be greater than 0"))
            end

            position = next_sequence(store, database, "events")
            event.metadata.position = position

            doc = BSON(
                "id" => event.id,
                "name" => event.name,
                "context" => BSON(
                    "name" => event.context.name,
                ),
                "aggregate" => BSON(
                    "name" => event.aggregate.name,
                    "id" => event.aggregate.id,
                    "revision" => event.aggregate.revision,
                ),
                "metadata" => BSON(
                    "datetime" => event.metadata.datetime,
                    "correlation_id" => event.metadata.correlation_id,
                    "causation_id" => event.metadata.causation_id,
                    "published" => event.metadata.published,
                    "position" => event.metadata.position,
                ),
                "custom" => event.custom,
                "data" => event.data,
            )

            try
                push!(collection, doc)
            catch err
                # TODO: Make sure this is an error we expect
                # and map the error appropriately
                io = IOBuffer()
                show(io, err)
                errtext = String(take!(io))
                if err.code == 11000 && !isnothing(match(r"aggregate_id_revision", errtext))
                    throw(EventException("duplicate revision"))
                end

                throw(err)
            end

            # Auto-snapshot
            if event.aggregate.revision % store.auto_snapshot == 0 && !isnothing(event.state)
                save_snapshot(store, database, event)
            end

            push!(committed_events, event)
        end
    end

    committed_events
end

function aggregate_last_event(store::MongoEventStore, aggregate_id::UUID)
    client = store.client
    database = client[store.dbname]
    collection = database[events_collection_name(store)]

    query = BSON(
        "aggregate.id" => aggregate_id,
    )

    options = BSON(
        "sort" => BSON("aggregate.revision" => -1),
        "projection" => BSON("_id" => 0),
    )

    doc = Mongoc.find_one(collection, query; options)

    if isnothing(doc)
        doc
    else
        parse(Event, doc)
    end
end

function mark_published_events(store::MongoEventStore, events::Vector{Event})
    client = store.client
    database = client[store.dbname]
    collection = database[events_collection_name(store)]

    query = BSON(
        "id" => BSON(
            "\$in" => map(events) do event event.id end,
        ),
    )

    update = BSON(
        "\$set" => BSON("metadata.published" => true),
    )

    result = Mongoc.update_many(collection, query, update)
    for event in events
        event.metadata.published = true
    end

    events
end

function get_snapshot(store::MongoEventStore, aggregate_id::UUID)
    client = store.client
    database = client[store.dbname]
    collection = database[snapshots_collection_name(store)]

    query = BSON(
        "aggregate.id" => aggregate_id,
    )

    options = BSON(
        "projection" => BSON("_id" => 0),
    )

    doc = Mongoc.find_one(collection, query; options)

    # TODO: Automatically convert to an aggregate
    doc
end

function save_snapshot(store::MongoEventStore, database::Mongoc.DatabaseSession, event::Event)
    if isnothing(event.state)
        throw(ArgumentError("event must have a state to create a snapshot"))
    end

    collection = database[snapshots_collection_name(store)]

    query = BSON(
        "aggregate.id" => event.aggregate.id,
    )

    setter = BSON(
        "\$set" => BSON(
            "aggregate" => BSON(
                "name" => event.aggregate.name,
                "id" => event.aggregate.id,
                "revision" => event.aggregate.revision,
            ),
            "state" => BSON(event.state),
        )
    )

    options = BSON(
        "upsert" => true,
    )

    Mongoc.update_one(collection, query, setter; options)

    nothing
end

function replay_events(store::MongoEventStore; from=1, to=2 ^ 31 -1)
    client = store.client
    database = client[store.dbname]
    collection = database[events_collection_name(store)]

    query = BSON(
        "metadata.position" => BSON(
            "\$gte" => from,
            "\$lte" => to,
        ),
    )

    options = BSON(
        "sort" => BSON("metadata.position" => 1),
        "projection" => BSON("_id" => 0),
    )

    cursor = Mongoc.find(collection, query; options)

    EventsIterator(cursor)
end

struct EventsIterator
    cursor::Mongoc.Cursor
end

function Base.iterate(iter::EventsIterator, state::Nothing=nothing)
    next = Base.iterate(iter.cursor)
    if isnothing(next)
        nothing
    else
        (parse(Event, next[1]), nothing)
    end
end

function Base.collect(iter::EventsIterator)
    out = Vector{Event}()

    for event in iter
        push!(out, event)
    end

    out
end


text:
# Rational Quadratic Covariance Function

include("rq_iso.jl")
include("rq_ard.jl")

@doc """
# Description
Constructor for the Rational Quadratic kernel

# See also
RQIso, RQArd
""" ->
RQ(ll::Float64, lσ::Float64, lα::Float64) = RQIso(ll, lσ, lα)
RQ(ll::Vector{Float64}, lσ::Float64, lα::Float64) = RQArd(ll, lσ, lα)


text:
################################################################################
#
#    This file contains utility functions
#           1. Functions for parameter containers (redefined NamedTuples)
#               1.1 Adding parameter containers
#               1.2 Converting parameter containers into vectors.
#           2. Positions
#           3. Configurations
#               - A Configuration is a type that holds information typically passed to LAMMPS.
#                   - A Configuration is really just a julia wrapper for information contained in DATA files.
#                   - Configuration contains Positions (depending on the number of atoms).
#               - Configurations are particularly important for SNAP implementation, since
#                   SNAP potential, forces, virial, and learning makes use of LAMMPS (and hence should
#                   use configurations, instead of Positions).
#           4. System (not yet implemented)
#               - A system is a collection of Configurations. (Primarily used for training SNAP potentials)

################################################################################
import Base.+
import Base.-
import Base.*
import Base.copy
import Base.vec
import Base.length
# include("param.jl")
include("nnlist.jl")


text:
# LSL.jl: Julia interface for Lab Streaming Layer
# Copyright (C) 2019 Samuel Powell

# test/sample.jl: test sending and receiving of samples supported datatypes

@static if !Sys.iswindows() && Sys.WORD_SIZE == 64
  testtypes = (Int8, Int16, Int32, Int64, Float32, Float64)
else
  testtypes = (Int8, Int16, Int32, Float32, Float64)
end

@testset "Samples" begin
  @testset "Datatype $T" for T in testtypes

    # Create a stream outlet
    count = 8

    info = StreamInfo(name = "TestName",
                      type = "Test$(string(T))",
                      channel_count = count,
                      channel_format = T,
                      source_id = "SampleTest$(string(T))ID")

    outlet = StreamOutlet(info)

    # Resolve the stream with modest timeout, and create inlet
    streams = resolve_byprop("source_id", "SampleTest$(string(T))ID", timeout = 5.0)
    inlet = StreamInlet(streams[1])
    open_stream(inlet)

    # Make some data, and send it once the inlet has opened
    data_in = rand(T, count)
    push_sample(outlet, data_in)

    # Pull the sample
    timestamp, data_out = pull_sample(inlet, timeout = 10.0)
    @test all(data_out .== data_in)

    # Close the stream
    close_stream(inlet)

  end

  # Force cleanup
  GC.gc()

end


text:
# syntax: proto3
using Compat
using ProtoBuf
import ProtoBuf.meta
import Base: hash, isequal, ==

type __enum_Syntax <: ProtoEnum
    SYNTAX_PROTO2::Int32
    SYNTAX_PROTO3::Int32
    __enum_Syntax() = new(0,1)
end #type __enum_Syntax
const Syntax = __enum_Syntax()

type Option
    name::AbstractString
    value::_Any
    Option(; kwargs...) = (o=new(); fillunset(o); isempty(kwargs) || ProtoBuf._protobuild(o, kwargs); o)
end #type Option
hash(v::Option) = ProtoBuf.protohash(v)
isequal(v1::Option, v2::Option) = ProtoBuf.protoisequal(v1, v2)
==(v1::Option, v2::Option) = ProtoBuf.protoeq(v1, v2)

type EnumValue
    name::AbstractString
    number::Int32
    options::Array{Option,1}
    EnumValue(; kwargs...) = (o=new(); fillunset(o); isempty(kwargs) || ProtoBuf._protobuild(o, kwargs); o)
end #type EnumValue
hash(v::EnumValue) = ProtoBuf.protohash(v)
isequal(v1::EnumValue, v2::EnumValue) = ProtoBuf.protoisequal(v1, v2)
==(v1::EnumValue, v2::EnumValue) = ProtoBuf.protoeq(v1, v2)

type __enum_Field_Kind <: ProtoEnum
    TYPE_UNKNOWN::Int32
    TYPE_DOUBLE::Int32
    TYPE_FLOAT::Int32
    TYPE_INT64::Int32
    TYPE_UINT64::Int32
    TYPE_INT32::Int32
    TYPE_FIXED64::Int32
    TYPE_FIXED32::Int32
    TYPE_BOOL::Int32
    TYPE_STRING::Int32
    TYPE_GROUP::Int32
    TYPE_MESSAGE::Int32
    TYPE_BYTES::Int32
    TYPE_UINT32::Int32
    TYPE_ENUM::Int32
    TYPE_SFIXED32::Int32
    TYPE_SFIXED64::Int32
    TYPE_SINT32::Int32
    TYPE_SINT64::Int32
    __enum_Field_Kind() = new(0,1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18)
end #type __enum_Field_Kind
const Field_Kind = __enum_Field_Kind()

type __enum_Field_Cardinality <: ProtoEnum
    CARDINALITY_UNKNOWN::Int32
    CARDINALITY_OPTIONAL::Int32
    CARDINALITY_REQUIRED::Int32
    CARDINALITY_REPEATED::Int32
    __enum_Field_Cardinality() = new(0,1,2,3)
end #type __enum_Field_Cardinality
const Field_Cardinality = __enum_Field_Cardinality()

type Field
    kind::Int32
    cardinality::Int32
    number::Int32
    name::AbstractString
    type_url::AbstractString
    oneof_index::Int32
    packed::Bool
    options::Array{Option,1}
    json_name::AbstractString
    default_value::AbstractString
    Field(; kwargs...) = (o=new(); fillunset(o); isempty(kwargs) || ProtoBuf._protobuild(o, kwargs); o)
end #type Field
const __fnum_Field = Int[1,2,3,4,6,7,8,9,10,11]
meta(t::Type{Field}) = meta(t, ProtoBuf.DEF_REQ, __fnum_Field, ProtoBuf.DEF_VAL, true, ProtoBuf.DEF_PACK, ProtoBuf.DEF_WTYPES, ProtoBuf.DEF_ONEOFS, ProtoBuf.DEF_ONEOF_NAMES)
hash(v::Field) = ProtoBuf.protohash(v)
isequal(v1::Field, v2::Field) = ProtoBuf.protoisequal(v1, v2)
==(v1::Field, v2::Field) = ProtoBuf.protoeq(v1, v2)

type _Enum
    name::AbstractString
    enumvalue::Array{EnumValue,1}
    options::Array{Option,1}
    source_context::SourceContext
    syntax::Int32
    _Enum(; kwargs...) = (o=new(); fillunset(o); isempty(kwargs) || ProtoBuf._protobuild(o, kwargs); o)
end #type _Enum
hash(v::_Enum) = ProtoBuf.protohash(v)
isequal(v1::_Enum, v2::_Enum) = ProtoBuf.protoisequal(v1, v2)
==(v1::_Enum, v2::_Enum) = ProtoBuf.protoeq(v1, v2)

type _Type
    name::AbstractString
    fields::Array{Field,1}
    oneofs::Array{AbstractString,1}
    options::Array{Option,1}
    source_context::SourceContext
    syntax::Int32
    _Type(; kwargs...) = (o=new(); fillunset(o); isempty(kwargs) || ProtoBuf._protobuild(o, kwargs); o)
end #type _Type
hash(v::_Type) = ProtoBuf.protohash(v)
isequal(v1::_Type, v2::_Type) = ProtoBuf.protoisequal(v1, v2)
==(v1::_Type, v2::_Type) = ProtoBuf.protoeq(v1, v2)

export Syntax, _Type, Field_Kind, Field_Cardinality, Field, _Enum, EnumValue, Option


text:
# Supprime la 1ère ligne de motif pat du txt
function eatfirstlinewithpat(txt, pat)
    lines = split(txt, r"\n")
    for i in 1:length(lines)
        if occursin(pat,lines[i])
            deleteat!(lines, i)
            break
        end
    end
    join(lines, "\n")
end


text:
function simulate(gen::AtDynamic{T,U,K}, args) where {T,U,K}
    (key::K, kernel_args) = args
    subtrace::U = simulate(gen.kernel, kernel_args)
    sub_call = get_call_record(subtrace)
    call = CallRecord{T}(sub_call.score, sub_call.retval, args)
    AtDynamicTrace{T,U,K}(call, subtrace, key)
end

max_stars_repo_name: BoltzmannBrain/Gen

max_issues_repo_name: BoltzmannBrain/Gen

max_forks_repo_name: BoltzmannBrain/Gen


text:
function validate_available_devices(device::Type{D}, devices) where {D <: PSY.Device}
    if isempty(devices)
        @warn(
            "The data doesn't include devices of type $(device), consider changing the device models"
        )
        return false
    end
    return true
end

function validate_services!(
    service::Type{S},
    services::Vector{S},
    incompatible_device_types::Vector{<:DataType},
    sys::PSY.System,
) where {S <: PSY.Service}
    services_ = PSY.get_components(S, sys)
    if isempty(services_)
        @warn(
            "The data doesn't include services of type $(service), consider changing the service models"
        )
        return false
    end

    services_mapping = PSY.get_contributing_device_mapping(sys)
    for s in services_
        contributing_devices_ =
            services_mapping[(type = S, name = PSY.get_name(s))].contributing_devices
        contributing_devices =
            [d for d in contributing_devices_ if typeof(d) ∉ incompatible_device_types]
        if isempty(contributing_devices)
            @warn(
                "The contributing devices for service $(PSY.get_name(s)) is empty, consider removing the service from the system"
            )
        else
            push!(services, s)
        end
    end
    return true
end

function validate_services!(
    service::Type{S},
    services::Vector{S},
    ::Vector{<:DataType},
    sys::PSY.System,
) where {S <: PSY.StaticReserveGroup}
    services_ = PSY.get_components(S, sys)
    if isempty(services_)
        @warn(
            "The data doesn't include services of type $(service), consider changing the service models"
        )
        return false
    end

    for s in services_
        contributing_services = PSY.get_contributing_services(s)
        if isempty(contributing_services)
            @warn(
                "The contributing services for group service $(PSY.get_name(service)) is empty, consider removing the group service from the system"
            )
        else
            push!(services, s)
        end
    end
    return !isempty(services)
end


text:
module DataFramesMeta

using Reexport

using MacroTools

@reexport using DataFrames

@reexport using Chain

# Basics:
export @with,
       @subset, @subset!, @rsubset, @rsubset!,
       @orderby, @rorderby,
       @by, @combine,
       @transform, @select, @transform!, @select!,
       @rtransform, @rselect, @rtransform!, @rselect!,
       @eachrow, @eachrow!,
       @byrow, @passmissing,
       @based_on, @where # deprecated

const DOLLAR = raw"$"

include("parsing.jl")
include("macros.jl")
include("linqmacro.jl")
include("eachrow.jl")

end # module


text:
using Roots
using Polynomials
using Base.Test
using Primes

## can use functions
f = x -> (x-1)*(x-2)^2*(x-3)^3
zs, mults = multroot(f)
@test mults == [1,2,3]

x = poly([0.0])

p = (x-1)*(x-2)^2*(x-3)^3
zs, mults = multroot(p)
@test mults == [1,2,3]

p = (x-1)^2*(x-2)^2*(x-3)^4
zs, mults = multroot(p)
@test mults == [2,2,4]

p = (x-1)^2
zs, mults = multroot(p^14)
@test mults == [28]

## test for roots of polynomial functions
roots(x -> x^5 - x + 1)

## test for real roots of polynomial functions
fzeros(x -> x^5 - 1.5x + 1)

## for polynomials in Z[x], Q[x] can use algorithm to be accurate for higher degree

@test fzeros(x -> x - 1)[1] == 1.0 # linear
f = x -> (x-1)*(x-2)*(x-3)^3*(x^2+1)
rts = fzeros(f)
rts = Float64[r for r in rts]
@test maximum(abs(sort(rts) - [1.0, 2.0, 3.0])) <= 1e-12
x = poly([big(0)])
p = prod([x - i for i in 1:20])
Roots.real_roots(p) ## can find this
f = x -> (x-20)^5 - (x-20) + 1
a = fzeros(f)[1]
@assert abs(f(a)) <= 1e-14

x = poly([0.0])
@test abs(fzeros((x-20)^5 - (x-20) + 1)[1] - (20 + fzeros(x^5 - x + 1)[1])) <= 1/2

fzeros(x -> x^5 - 2x^4 + x^3)

## factor
factor(x -> (x-2)^4*(x-3)^9)
factor(x -> (x-1)^3 * (x-2)^3 * (x^5 - x + 1))
factor(x -> x*(x-1)*(x-2)*(x^2 + x + 1))

factor(x -> (x-1)^2 * (x-.99)^2 * (x-1.01)^2) ## can have issue with nearby roots (or high powers)

factor(x -> (x-1//1)^2 * (x-99//100)^2 * (x-101//100)^2) ## conversion is to Float, not Rational{Int}
delta = 1//10
VERSION >= v"0.5.0-" && Primes.factor(convert(Poly{Rational{Int}}, x -> (x-1//1)^2 * (x-1 - delta)^2 * (x-1 + delta)^2))

## Test conversion of polynomials to Int
f = x -> 3x^3 - 2x
convert(Polynomials.Poly{Int}, f)
f = x -> (3x^3 - 2x)/3
convert(Polynomials.Poly{Int}, f)
f = x ->  x^2 / x # rational functions fails
@test_throws MethodError convert(Polynomials.Poly{Int64}, f)
f = x ->  (x^2 + x)^2
convert(Polynomials.Poly{Int}, f)
f = x ->  (x^2 + x)^(1/2)
@test_throws MethodError convert(Polynomials.Poly{Int64}, f)

## polynomial conversions from functions in practice
fzeros(x -> 265 - 0.65x)
fzeros(x -> -16x^2 + 200x)


text:
module GravitationalAcceleration

# This implements the gravitational acceleration

using Parameters, LaTeXStrings, Unitful
using ..Units
using GeoParams: AbstractMaterialParam
import Base.show

abstract type AbstractGravity <: AbstractMaterialParam end

export  ComputeGravity,        # calculation routines
        ConstantGravity        # constant

# Constant Gravity -------------------------------------------------------
"""
    GravityConstant(g=9.81m/s^2)

Set a constant value for the gravitational acceleration:
```math
    g  = 9.81 m s^{-2}
```
"""
@with_kw_noshow mutable struct ConstantGravity <: AbstractGravity
    equation::LaTeXString   =   L"g = 9.81 m s^{-2}"
    g::GeoUnit              =   9.81m/s^2               # gravitational acceleration
end

# Calculation routine
function ComputeGravity(s::ConstantGravity)
    @unpack g   = s

    return g*1.0   # multiply with 1.0, to return Float64
end

# Print info
function show(io::IO, d::ConstantGravity)
    print(io, "Gravitational acceleration: g=$(d.g.val)")
end
#-------------------------------------------------------------------------

# Help info for the calculation routines
"""
ComputeGravity(s:<AbstractGravity)

Returns the gravitational acceleration

"""
ComputeGravity

end


text:
module IterTools

import Base.Iterators: drop, take

import Base: iterate, eltype, length, size, peek
import Base: tail
import Base: IteratorSize, IteratorEltype
import Base: SizeUnknown, IsInfinite, HasLength, HasShape
import Base: HasEltype, EltypeUnknown

export
    takestrict,
    repeatedly,
    chain,
    product,
    distinct,
    partition,
    groupby,
    imap,
    subsets,
    iterated,
    nth,
    takenth,
    peekiter,
    peek,
    ncycle,
    ivec,
    flagfirst,
    takewhile,
    properties,
    propertyvalues,
    fieldvalues

function has_length(it)
    it_size = IteratorSize(it)

    return isa(it_size, HasLength) || isa(it_size, HasShape)
end

# return the size for methods depending on the longest iterator
longest(::T, ::T) where {T<:IteratorSize} = T()
function longest(::S, ::T) where {T<:IteratorSize, S<:IteratorSize}
    longest(T(), S())
end
longest(::HasShape, ::HasShape) = HasLength()
longest(::HasLength, ::HasShape) = HasLength()
longest(::SizeUnknown, ::HasShape) = SizeUnknown()
longest(::SizeUnknown, ::HasLength) = SizeUnknown()
longest(::IsInfinite, ::HasShape) = IsInfinite()
longest(::IsInfinite, ::HasLength) = IsInfinite()
longest(::IsInfinite, ::SizeUnknown) = IsInfinite()

# return the size for methods depending on the shortest iterator
shortest(::T, ::T) where {T<:IteratorSize} = T()
function shortest(::S, ::T) where {T<:IteratorSize, S<:IteratorSize}
    shortest(T(), S())
end
shortest(::HasShape, ::HasShape) = HasLength()
shortest(::HasLength, ::HasShape) = HasLength()
shortest(::IsInfinite, ::HasShape) = HasLength()
shortest(::IsInfinite, ::HasLength) = HasLength()
shortest(::SizeUnknown, ::HasShape) = SizeUnknown()
shortest(::SizeUnknown, ::HasLength) = SizeUnknown()
shortest(::SizeUnknown, ::IsInfinite) = SizeUnknown()

"""
    IterTools.@ifsomething expr

If `expr` evaluates to `nothing`, equivalent to `return nothing`, otherwise the macro
evaluates to the value of `expr`. Not exported, useful for implementing iterators.

```jldoctest
julia> IterTools.@ifsomething iterate(1:2)
(1, 1)

julia> let elt, state = IterTools.@ifsomething iterate(1:2, 2); println("not reached"); end
```
"""
macro ifsomething(ex)
    quote
        result = $(esc(ex))
        result === nothing && return nothing
        result
    end
end

# Iterate through the first n elements, throwing an exception if
# fewer than n items ar encountered.

struct TakeStrict{I}
    xs::I
    n::Int
end
IteratorSize(::Type{<:TakeStrict}) = HasLength()
IteratorEltype(::Type{TakeStrict{I}}) where {I} = IteratorEltype(I)
eltype(::Type{TakeStrict{I}}) where {I} = eltype(I)

"""
    takestrict(xs, n::Int)

Like `take()`, an iterator that generates at most the first `n` elements of `xs`, but throws
an exception if fewer than `n` items are encountered in `xs`.

```jldoctest
julia> a = :1:2:11
1:2:11

julia> collect(takestrict(a, 3))
3-element Array{Int64,1}:
 1
 3
 5
```
"""
takestrict(xs, n::Int) = TakeStrict(xs, n)

function iterate(it::TakeStrict, state=(it.n,))
    n, xs_state = first(state), tail(state)
    n <= 0 && return nothing
    xs_iter = iterate(it.xs, xs_state...)

    if xs_iter === nothing
        throw(ArgumentError("In takestrict(xs, n), xs had fewer than n items to take."))
    end

    v, xs_state = xs_iter
    return v, (n - 1, xs_state)
end

length(it::TakeStrict) = it.n

# Repeat a function application n (or infinitely many) times.

struct RepeatCall{F<:Base.Callable}
    f::F
    n::Int
end
IteratorSize(::Type{<:RepeatCall}) = HasLength()
IteratorEltype(::Type{<:RepeatCall}) = EltypeUnknown()
length(it::RepeatCall) = it.n

"""
    repeatedly(f)
    repeatedly(f, n)

Call function `f` `n` times, or infinitely if `n` is omitted.

```julia
julia> t() = (sleep(0.1); Dates.millisecond(now()))
t (generic function with 1 method)

julia> collect(repeatedly(t, 5))
5-element Array{Any,1}:
 993
  97
 200
 303
 408
```
"""
repeatedly(f, n) = RepeatCall(f, n)
iterate(it::RepeatCall, state=it.n) = state <= 0 ? nothing : (it.f(), state - 1)

struct RepeatCallForever{F<:Base.Callable}
    f::F
end
IteratorSize(::Type{<:RepeatCallForever}) = IsInfinite()
IteratorEltype(::Type{<:RepeatCallForever}) = EltypeUnknown()

repeatedly(f) = RepeatCallForever(f)
iterate(it::RepeatCallForever, state=nothing) = (it.f(), nothing)

@deprecate chain(xss...) Iterators.flatten(xss)
@deprecate product(xss...) Iterators.product(xss...)

# Filter out reccuring elements.

struct Distinct{I, J}
    xs::I

    # Map elements to the index at which it was first seen, so given an iterator
    # state (index) we can test if an element has previously been observed.
    seen::Dict{J, Int}
end

IteratorSize(::Type{<:Distinct}) = SizeUnknown()

eltype(::Type{Distinct{I, J}}) where {I, J} = J

"""
    distinct(xs)

Iterate through values skipping over those already encountered.

```jldoctest
julia> for i in distinct([1,1,2,1,2,4,1,2,3,4])
           @show i
       end
i = 1
i = 2
i = 4
i = 3
```
"""
distinct(xs::I) where {I} = Distinct{I, eltype(xs)}(xs, Dict{eltype(xs), Int}())

function iterate(it::Distinct, state=(1,))
    idx, xs_state = first(state), tail(state)
    xs_iter = iterate(it.xs, xs_state...)

    while xs_iter !== nothing
        val, xs_state = xs_iter
        get!(it.seen, val, idx) >= idx && return (val, (idx + 1, xs_state))

        xs_iter = iterate(it.xs, xs_state)
        idx += 1
    end

    return nothing
end

# Group output from at iterator into tuples.
# E.g.,
#   partition(count(1), 2) = (1,2), (3,4), (5,6) ...
#   partition(count(1), 2, 1) = (1,2), (2,3), (3,4) ...
#   partition(count(1), 2, 3) = (1,2), (4,5), (7,8) ...

struct Partition{I, N}
    xs::I
    step::Int
end
IteratorSize(::Type{<:Partition}) = SizeUnknown()

eltype(::Type{Partition{I, N}}) where {I, N} = NTuple{N, eltype(I)}

"""
    partition(xs, n, [step])

Group values into `n`-tuples.

```jldoctest
julia> for i in partition(1:9, 3)
           @show i
       end
i = (1, 2, 3)
i = (4, 5, 6)
i = (7, 8, 9)
```

If the `step` parameter is set, each tuple is separated by `step` values.

```jldoctest
julia> for i in partition(1:9, 3, 2)
           @show i
       end
i = (1, 2, 3)
i = (3, 4, 5)
i = (5, 6, 7)
i = (7, 8, 9)

julia> for i in partition(1:9, 3, 3)
           @show i
       end
i = (1, 2, 3)
i = (4, 5, 6)
i = (7, 8, 9)

julia> for i in partition(1:9, 2, 3)
           @show i
       end
i = (1, 2)
i = (4, 5)
i = (7, 8)
```
"""
@inline partition(xs, n::Int) = partition(xs, n, n)

function partition(xs::I, n::Int, step::Int) where I
    if step < 1
        throw(ArgumentError("Partition step must be at least 1."))
    end

    if n < 1
        throw(ArgumentError("Partition size n must be at least 1"))
    end

    Partition{I, n}(xs, step)
end

function iterate(it::Partition{I, N}, state=nothing) where {I, N}
    if state === nothing
        result = Vector{eltype(I)}(undef, N)

        result[1], xs_state = @ifsomething iterate(it.xs)

        for i in 2:N
            result[i], xs_state = @ifsomething iterate(it.xs, xs_state)
        end
    else
        (xs_state, result) = state
        result[end], xs_state = @ifsomething iterate(it.xs, xs_state)
    end

    p = similar(result)
    overlap = max(0, N - it.step)
    p[1:overlap] .= result[it.step .+ (1:overlap)]

    # when step > n, skip over some elements
    for i in 1:max(0, it.step - N)
        xs_iter = iterate(it.xs, xs_state)
        xs_iter === nothing && break
        _, xs_state = xs_iter
    end

    for i in (overlap + 1):(N - 1)
        xs_iter = iterate(it.xs, xs_state)
        xs_iter === nothing && break

        p[i], xs_state = xs_iter
    end

    return (tuple(result...)::eltype(Partition{I, N}), (xs_state, p))
end

# Group output from an iterator based on a key function.
# Consecutive entries from the iterator with the same
# key value will be returned in a single array.
# Inspired by itertools.groupby in python.
# E.g.,
#   x = ["face", "foo", "bar", "book", "baz", "zzz"]
#   groupby(z -> z[1], x) =
#       ["face", "foo"]
#       ["bar", "book", "baz"]
#       ["zzz"]
struct GroupBy{I, F<:Base.Callable}
    keyfunc::F
    xs::I
end
IteratorSize(::Type{<:GroupBy}) = SizeUnknown()

eltype(::Type{<:GroupBy{I}}) where {I} = Vector{eltype(I)}

"""
    groupby(f, xs)

Group consecutive values that share the same result of applying `f`.

```jldoctest
julia> for i in groupby(x -> x[1], ["face", "foo", "bar", "book", "baz", "zzz"])
           @show i
       end
i = ["face", "foo"]
i = ["bar", "book", "baz"]
i = ["zzz"]
```
"""
function groupby(keyfunc::F, xs::I) where {F<:Base.Callable, I}
    GroupBy{I, F}(keyfunc, xs)
end

function iterate(it::GroupBy{I, F}, state=nothing) where {I, F<:Base.Callable}
    if state === nothing
        prev_val, xs_state = @ifsomething iterate(it.xs)
        prev_key = it.keyfunc(prev_val)
        keep_going = true
    else
        keep_going, prev_key, prev_val, xs_state = state
        keep_going || return nothing
    end
    values = Vector{eltype(I)}()
    push!(values, prev_val)

    while true
        xs_iter = iterate(it.xs, xs_state)

        if xs_iter === nothing
            keep_going = false
            break
        end

        val, xs_state = xs_iter
        key = it.keyfunc(val)

        if key == prev_key
            push!(values, val)
        else
            prev_key = key
            prev_val = val
            break
        end
    end

    return (values, (keep_going, prev_key, prev_val, xs_state))
end

"""
    imap(f, xs1, [xs2, ...])

Iterate over values of a function applied to successive values from one or more iterators.
Like `Iterators.zip`, the iterator is done when any of the input iterators have been
exhausted.

```jldoctest
julia> for i in imap(+, [1,2,3], [4,5,6])
            @show i
       end
i = 5
i = 7
i = 9
```
"""
imap(mapfunc, it1, its...) = (mapfunc(xs...) for xs in zip(it1, its...))

# Iterate over all subsets of an indexable collection

struct Subsets{C}
    xs::C
end
IteratorSize(::Type{Subsets{C}}) where {C} = longest(HasLength(), IteratorSize(C))

eltype(::Type{Subsets{C}}) where {C} = Vector{eltype(C)}
length(it::Subsets) = 1 << length(it.xs)

"""
    subsets(xs)
    subsets(xs, k)
    subsets(xs, Val{k}())

Iterate over every subset of the indexable collection `xs`. You can restrict the subsets to a
specific size `k`.

Giving the subset size in the form `Val{k}()` allows the compiler to produce code optimized
for the particular size requested. This leads to performance comparable to hand-written
loops if `k` is small and known at compile time, but may or may not improve performance
otherwise.

```jldoctest
julia> for i in subsets([1, 2, 3])
          @show i
       end
i = Int64[]
i = [1]
i = [2]
i = [1, 2]
i = [3]
i = [1, 3]
i = [2, 3]
i = [1, 2, 3]

julia> for i in subsets(1:4, 2)
          @show i
       end
i = [1, 2]
i = [1, 3]
i = [1, 4]
i = [2, 3]
i = [2, 4]
i = [3, 4]

julia> for i in subsets(1:4, Val{2}())
           @show i
       end
i = (1, 2)
i = (1, 3)
i = (1, 4)
i = (2, 3)
i = (2, 4)
i = (3, 4)
```
"""
function subsets(xs)
    Subsets(xs)
end

# state has one extra bit to indicate that we are at the end
function iterate(it::Subsets, state=fill(false, length(it.xs) + 1))
    state[end] && return nothing

    ss = it.xs[state[1:end-1]]

    state = copy(state)
    state[1] = !state[1]
    for i in 2:length(state)
        if !state[i-1]
            state[i] = !state[i]
        else
            break
        end
    end

    return (ss, state)
end

# Iterate over all subsets of an indexable collection with a given size

struct Binomial{Collection}
    xs::Collection
    n::Int64
    k::Int64
end
Binomial(xs::C, n::Integer, k::Integer) where {C} = Binomial{C}(xs, n, k)

IteratorSize(::Type{<:Binomial}) = HasLength()
IteratorEltype(::Type{Binomial{C}}) where {C} = IteratorEltype(C)

eltype(::Type{Binomial{C}}) where {C} = Vector{eltype(C)}
length(it::Binomial) = binomial(it.n,it.k)

subsets(xs, k) = Binomial(xs, length(xs), k)

mutable struct BinomialIterState
    idx::Vector{Int64}
    done::Bool
end

function iterate(it::Binomial, state=BinomialIterState(collect(Int64, 1:it.k), it.k > it.n))
    state.done && return nothing

    idx = state.idx
    set = it.xs[idx]
    i = it.k
    while i > 0
        if idx[i] < it.n - it.k + i
            idx[i] += 1

            for j in 1:it.k-i
                idx[i+j] = idx[i] + j
            end

            break
        else
            i -= 1
        end
    end

    state.done = i == 0

    return set, state
end

# Iterate over all subsets of an indexable collection with a given *statically* known size

struct StaticSizeBinomial{K, Container}
    xs::Container
end

IteratorSize(::Type{StaticSizeBinomial{K, C}}) where {K, C} = HasLength()
IteratorEltype(::Type{StaticSizeBinomial{K, C}}) where {K, C} = IteratorEltype(C)

eltype(::Type{StaticSizeBinomial{K, C}}) where {K, C} = NTuple{K, eltype(C)}
length(it::StaticSizeBinomial{K}) where {K} = binomial(length(it.xs), K)

subsets(xs::C, ::Val{K}) where {K, C} = StaticSizeBinomial{K, C}(xs)

# Special case for K == 0
iterate(it::StaticSizeBinomial{0}, state=false) = state ? nothing : ((), true)

# Generic case K >= 1
pop(t::NTuple) = reverse(tail(reverse(t))), t[end]

function advance(it::StaticSizeBinomial{K}, idx) where {K}
	xs = it.xs
	lidx, i = pop(idx)
    i += 1
	if i > length(xs) - K + length(idx)
		lidx = advance(it, lidx)
		i = lidx[end] + 1
	end
	return (lidx..., i)
end
advance(it::StaticSizeBinomial, idx::NTuple{1}) = (idx[end]+1,)

function iterate(it::StaticSizeBinomial{K}, idx=ntuple(identity, Val{K}())) where K
    idx[end] > length(it.xs) && return nothing
    return (map(i -> it.xs[i], idx), advance(it, idx))
end

# nth : return the nth element in a collection

"""
    nth(xs, n)

Return the `n`th element of `xs`. This is mostly useful for non-indexable collections.

```jldoctest
julia> mersenne = Set([3, 7, 31, 127])
Set([7, 31, 3, 127])

julia> nth(mersenne, 3)
3
```
"""
function nth(xs, n::Integer)
    n > 0 || throw(BoundsError(xs, n))

    # catch, if possible
    has_length(xs) && (n ≤ length(xs) || throw(BoundsError(xs, n)))

    for (i, val) in enumerate(xs)
        i >= n && return val
    end

    # catch iterators with no length but actual finite size less then n
    throw(BoundsError(xs, n))
end

nth(xs::Union{Tuple, Array}, n::Integer) = xs[n]

function nth(xs::AbstractArray, n::Integer)
    idx = eachindex(xs)[n]
    return @inbounds xs[idx]
end

# takenth(xs,n): take every n'th element from xs

struct TakeNth{I}
    xs::I
    interval::UInt
end
IteratorSize(::Type{TakeNth{I}}) where {I} = longest(HasLength(), IteratorSize(I))
IteratorEltype(::Type{TakeNth{I}}) where {I} = IteratorEltype(I)
eltype(::Type{TakeNth{I}}) where {I} = eltype(I)
length(x::TakeNth) = div(length(x.xs), x.interval)

"""
    takenth(xs, n)

Iterate through every `n`th element of `xs`.

```jldoctest
julia> collect(takenth(5:15,3))
3-element Array{Int64,1}:
  7
 10
 13
```
"""
function takenth(xs, interval::Integer)
    if interval <= 0
        throw(ArgumentError(string("expected interval to be 1 or more, ",
                                   "got $interval")))
    end
    TakeNth(xs, convert(UInt, interval))
end

function iterate(it::TakeNth, state...)
    xs_iter = nothing

    for i = 1:it.interval
        xs_iter = @ifsomething iterate(it.xs, state...)
        state = tail(xs_iter)
    end

    return xs_iter
end

struct Iterated{T, F}
    f::F
    seed::T
end
IteratorSize(::Type{<:Iterated}) = IsInfinite()
IteratorEltype(::Type{<:Iterated}) = EltypeUnknown()

"""
    iterated(f, x)

Iterate over successive applications of `f`, as in `x`, `f(x)`, `f(f(x))`, `f(f(f(x)))`, ...

Use `Base.Iterators.take()` to obtain the required number of elements.

```jldoctest
julia> for i in Iterators.take(iterated(x -> 2x, 1), 5)
           @show i
       end
i = 1
i = 2
i = 4
i = 8
i = 16

julia> for i in Iterators.take(iterated(sqrt, 100), 6)
           @show i
       end
i = 100
i = 10.0
i = 3.1622776601683795
i = 1.7782794100389228
i = 1.333521432163324
i = 1.1547819846894583
```
"""
iterated(f, seed) = Iterated(f, seed)

iterate(it::Iterated) = (it.seed, it.seed)
function iterate(it::Iterated{T, F}, state) where {T, F}
    newval = it.f(state)
    return (newval, newval)
end

# peekiter(iter): possibility to peek the head of an iterator

struct PeekIter{I}
    it::I
end

"""
    peekiter(xs)

Lets you peek at the head element of an iterator without updating the state.

```jldoctest
julia> it = peekiter(["face", "foo", "bar", "book", "baz", "zzz"])
IterTools.PeekIter{Array{String,1}}(["face", "foo", "bar", "book", "baz", "zzz"])

julia> @show peek(it);
peek(it) = Some("face")

julia> @show peek(it);
peek(it) = Some("face")

julia> x, s = iterate(it)
("face", ("foo", 3))

julia> @show x;
x = "face"

julia> @show peek(it, s);
peek(it, s) = Some("foo")
```
"""
peekiter(itr) = PeekIter(itr)

eltype(::Type{PeekIter{I}}) where {I} = eltype(I)
IteratorSize(::Type{PeekIter{I}}) where {I} = IteratorSize(I)
IteratorEltype(::Type{PeekIter{I}}) where {I} = IteratorEltype(I)
length(f::PeekIter) = length(f.it)
size(f::PeekIter) = size(f.it)

function iterate(pit::PeekIter, state=iterate(pit.it))
    val, it_state = @ifsomething state
    return (val, iterate(pit.it, it_state))
end

peek(pit::PeekIter, state=iterate(pit)) = Some{eltype(pit)}(first(@ifsomething state))

#NCycle - cycle through an object N times

struct NCycle{I}
    iter::I
    n::Int
end

"""
    ncycle(iter, n)

Cycle through `iter` `n` times.

```jldoctest
julia> for i in ncycle(1:3, 2)
           @show i
       end
i = 1
i = 2
i = 3
i = 1
i = 2
i = 3
```
"""
ncycle(iter, n::Int) = NCycle(iter, n)

eltype(::Type{NCycle{I}}) where {I} = eltype(I)
length(nc::NCycle) = nc.n*length(nc.iter)
IteratorSize(::Type{NCycle{I}}) where {I} = longest(HasLength(), IteratorSize(I))
IteratorEltype(::Type{NCycle{I}}) where {I} = IteratorEltype(I)

function iterate(nc::NCycle, state=(nc.n,))
    nc.n <= 0 && return nothing  # don't do anything if we aren't iterating

    n, inner_state = first(state), tail(state)
    inner_iter = iterate(nc.iter, inner_state...)

    if inner_iter === nothing
        if n <= 1
            return nothing
        else
            inner_iter = @ifsomething iterate(nc.iter)

            n -= 1
        end
    end

    v, inner_state = inner_iter
    return v, (n, inner_state)
end

# IVec - lazy `vec`-like iterator that drops shape

struct IVec{I}
    iter::I
end

"""
    ivec(iter)

Drops all shape from `iter` while iterating.
Like a non-materializing version of [`vec`](https://docs.julialang.org/en/stable/base/arrays/#Base.vec).

```jldoctest
julia> m = collect(reshape(1:6, 2, 3))
2×3 Array{Int64,2}:
 1  3  5
 2  4  6

julia> collect(ivec(m))
6-element Array{Int64,1}:
 1
 2
 3
 4
 5
 6
```
"""
ivec(iter) = IVec(iter)

eltype(::Type{IVec{I}}) where {I} = eltype(I)
length(iv::IVec) = length(iv.iter)
IteratorSize(::Type{IVec{I}}) where {I} = longest(HasLength(), IteratorSize(I))
IteratorEltype(::Type{IVec{I}}) where {I} = IteratorEltype(I)

iterate(iv::IVec, state...) = iterate(iv.iter, state...)

# FlagFirst: prepend a flag that is `true` iff this is the first element

struct FlagFirst{I}
    iter::I
end

"""
    flagfirst(iter)

An iterator that yields `(isfirst, x)` where `isfirst::Bool` is `true` for the first
element, and `false` after that, while the `x`s are elements from `iter`.

```jldoctest
julia> collect(flagfirst(1:3))
3-element Array{Tuple{Bool,Int64},1}:
 (true, 1)
 (false, 2)
 (false, 3)
```
"""
flagfirst(iter) = FlagFirst(iter)

eltype(::Type{FlagFirst{I}}) where I = Tuple{Bool, eltype(I)}
length(ff::FlagFirst) = length(ff.iter)
size(ff::FlagFirst) = size(ff.iter)
IteratorSize(::Type{FlagFirst{I}}) where {I} = IteratorSize(I)
IteratorEltype(::Type{FlagFirst{I}}) where {I} = IteratorEltype(I)

function iterate(ff::FlagFirst, state = (true, ))
    isfirst, rest = first(state), tail(state)
    elt, nextstate = @ifsomething iterate(ff.iter, rest...)
    (isfirst, elt), (isfirst & false, nextstate)
end

# TakeWhile iterates through values from an iterable as long as a given predicate is true.

struct TakeWhile{I}
    cond::Function
    xs::I
end

"""
    takewhile(cond, xs)

An iterator that yields values from the iterator `xs` as long as the
predicate `cond` is true.

```jldoctest
julia> collect(takewhile(x-> x^2 < 10, 1:100))
3-element Array{Int64,1}:
 1
 2
 3
```
"""
takewhile(cond, xs) = TakeWhile(cond, xs)

function Base.iterate(it::TakeWhile, state=nothing)
    (val, state) = @ifsomething (state === nothing ? iterate(it.xs) : iterate(it.xs, state))
    it.cond(val) || return nothing
    val, state
end

Base.IteratorSize(it::TakeWhile) = Base.SizeUnknown()
eltype(::Type{TakeWhile{I}}) where {I} = eltype(I)
IteratorEltype(::Type{TakeWhile{I}}) where {I} = IteratorEltype(I)

# Properties

struct Properties{T}
    x::T
    n::Int
    names
end

"""
    properties(x)

Iterate through the names and value of the properties of `x`.

```jldoctest
julia> collect(properties(1 + 2im))
2-element Array{Any,1}:
 (:re, 1)
 (:im, 2)
```
"""
function properties(x::T) where T
    names = propertynames(x)
    return Properties{T}(x, length(names), names)
end

function iterate(p::Properties, state=1)
    state > length(p) && return nothing

    name = p.names[state]
    return ((name, getproperty(p.x, name)), state + 1)
end

# PropertyValues

struct PropertyValues{T}
    x::T
    n::Int
    names
end

"""
    propertyvalues(x)

Iterate through the values of the properties of `x`.

```jldoctest
julia> collect(propertyvalues(1 + 2im))
2-element Array{Any,1}:
 1
 2
```
"""

function propertyvalues(x::T) where T
    names = propertynames(x)
    return PropertyValues{T}(x, length(names), names)
end

function iterate(p::PropertyValues, state=1)
    state > length(p) && return nothing

    name = p.names[state]
    return (getproperty(p.x, name), state + 1)
end

length(p::Union{Properties, PropertyValues}) = p.n
IteratorSize(::Type{<:Union{Properties, PropertyValues}}) = HasLength()

# FieldValues

struct FieldValues{T}
    x::T
end

"""
    fieldvalues(x)

Iterate through the values of the fields of `x`.

```jldoctest
julia> collect(fieldvalues(1 + 2im))
2-element Array{Any,1}:
 1
 2
```
"""
fieldvalues(x::T) where {T} = FieldValues{T}(x)
length(fs::FieldValues{T}) where {T} = fieldcount(T)
IteratorSize(::Type{<:FieldValues}) = HasLength()

function iterate(fs::FieldValues, state=1)
    state > length(fs) && return nothing

    return (getfield(fs.x, state), state + 1)
end

end # module IterTools


text:
"""
Read a graph from a text file. Skips any lines that begin with '#' or '%'
characters. If oneindex is true, then nodes are mapped to values 1, ..., n.

Returns (A, v), where
   A is the n x n undirected graph adjacency matrix
   v is a length-n vector such that v[i] is the original node identifier
     (if oneindex is false, then v[i] = i, i = 1, ..., n)
"""
function read_undir_graph_txt(filename::AbstractString, oneindex::Bool=false)
    # index mapping
    index_map = Dict{Int64,Int64}()
    index_map_vec = Int64[]
    function get_mapped_index(x::Int64)
        if !haskey(index_map, x)
            next_index = length(index_map) + 1
            index_map[x] = next_index
            push!(index_map_vec, x)
            return next_index
        end
        return index_map[x]
    end

    # Read data
    I = Int64[]
    J = Int64[]
    open(filename) do f
        for line in eachline(f)
            # Skip lines starting with '#' or '%'
            if line[1] == '#' || line[1] == '%'; continue; end
            edge = split(line)
            u = parse(Int64, edge[1])
            v = parse(Int64, edge[2])
            if oneindex
                push!(I, get_mapped_index(u))
                push!(J, get_mapped_index(v))
            else
                push!(I, u)
                push!(J, v)
            end
        end
    end

    # Form adjacency matrix
    if min(minimum(I), minimum(J)) < 1
        error("Minimum node value is less than 1. Try setting ondeindex parameter to true.")
    end
    n = max(maximum(I), maximum(J))
    A = convert(SparseMatrixCSC{Int64,Int64}, sparse(I, J, ones(length(I)), n, n))
    A = max.(A, A')

    if oneindex; return (A, index_map_vec); end
    return (A, collect(1:n))
end

"""
Load an example file from the data directory.
"""
function load_example_data(filename::AbstractString)
    pathname = joinpath(dirname(dirname(@__FILE__)), "data")
    filename = joinpath(pathname, filename)
    if   isfile(filename); return read_undir_graph_txt(filename, true)[1]
    else error(@sprintf("Could not find file %s", name))
    end
end


text:
using Random
using BenchmarkTools

function row_major(A::Array{Float64})

    m, n = size(A)

    for i in 1:m
        for j in 1:n
            A[i,j] = 1
        end
    end

end

function col_major(A::Array{Float64})

    m, n = size(A)

    for j in 1:n
        for i in 1:m
            A[i,j] = 1
        end
    end

end

function main()

    n = 100

    A = randn(n,n)

    # @benchmark row_major($A)
    @benchmark col_major($A)

end

main()

max_stars_repo_name: nunesmelo/djs-office-hours

max_issues_repo_name: nunesmelo/djs-office-hours

max_forks_repo_name: nunesmelo/djs-office-hours


text:
# --------------------------------------------------------------------------- #
# The "Multiscale_3layers" model from (Pan, Yang, Li, et al.; 2018)

struct Multiscale3{I, BN} <: FCModel{I, BN}
  weights :: Vector{Any}
  bnmoments :: Vector{Any}
end

# --------------------------------------------------------------------------- #
# Network architecture (no batch normalization)

function density{I <: Image}(w, s, x, ::Type{Multiscale3{I, false}})
  y1 = rconv(w[2], rconv(w[1], x))
  y2 = rconv(w[4], pool(rconv(w[3], x)))
  y3 = rconv(w[6], rconv(w[5], pool(y1)))

  y4 = rconv(w[8], rconv(w[7], cat3(pool(y2), pool(y3))))
  y5 = rconv(w[11], rconv(w[10], uconv(w[9], y3, y4)))
  y6 = rconv(w[14], rconv(w[13], uconv(w[12], y1, y5)))

  y = conv4(w[15], y6, padding=pad(w[15]))
  return reshape(y, size(y)[[1,2,4]]...)
end

# --------------------------------------------------------------------------- #
# Initialization (no batch normalization)

function init_weights{I<:Image}(::Type{Multiscale3{I, false}},
                                wr :: Function = wr,
                                wu :: Function = wu;
                                potency :: Integer = 6,
                                seed :: Integer = rand(1:10000))
  Knet.setseed(seed)
  c = imgchannels(I)
  s = 2^potency

  return Any[
    wr(c, s), wr(s, s), wr(c, s, (11,11)), wr(s, s, (9,9)),
    wr(s, s), wr(s, s), wr(2s, s), wr(s, s), wu(s, s),
    wr(2s, s), wr(s, s), wu(s, s), wr(2s, s), wr(s, s),
    wr(s, 1) ]
end

init_state{I<:Image}(::Type{Multiscale3{I, false}}) = Any[]

# --------------------------------------------------------------------------- #
# Network architecture (batch normalization)

function density{I <: Image}(w, s, x, ::Type{Multiscale3{I, true}})
  y1 = rbconv(w[4], s[2], w[3], rbconv(w[2], s[1], w[1], x))
  y2 = rbconv(w[8], s[4], w[7], pool(rbconv(w[6], s[3], w[5], x)))
  y3 = rbconv(w[12], s[6], w[11], rbconv(w[10], s[5], w[9], pool(y1)))

  u4 = cat3(pool(y2), pool(y3))
  y4 = rbconv(w[16], s[8], w[15], rbconv(w[14], s[7], w[13], u4))

  u5 = ubconv(w[18], s[9], w[17], y3, y4)
  y5 = rbconv(w[22], s[11], w[21], rbconv(w[20], s[10], w[19], u5))

  u6 = ubconv(w[24], s[12], w[23], y1, y5)
  y6 = rbconv(w[28], s[14], w[27], rbconv(w[26], s[13], w[25], u6))

  y = conv4(w[29], y6, padding=pad(w[29]))
  return reshape(y, size(y)[[1,2,4]]...)
end

# --------------------------------------------------------------------------- #
# Initialization (batch normalization)

function init_weights{I<:Image}(::Type{Multiscale3{I, true}},
                                wr :: Function = wr,
                                wu :: Function = wu;
                                potency :: Integer = 6,
                                seed :: Integer = rand(1:10000))

  # Need the same convolutional weights as in the case without batch
  # normalization. However, also need to insert the trainable bn-values at
  # the right positions of the weights vector
  s = 2^potency

  weights = init_weights(Multiscale3{I, false}, wr, wu, potency = potency, seed = seed)
  return Any[ odd(i) ? weights[ceil(Int, i/2)] : bn(s) for i in 1:29 ]
end

function init_state{I<:Image}(::Type{Multiscale3{I, true}})
  return Any[ bnmoments() for i in 1:14 ]
end

# --------------------------------------------------------------------------- #
# Convenience-constructor

function Multiscale3(I,
                     wr :: Function = wr,
                     wu :: Function = wu;
                     bn :: Bool = false,
                     potency :: Integer = 6,
                     seed :: Integer = rand(1:10000))

  weights = init_weights(Multiscale3{I, bn}, wr, wu,
                         potency = potency, seed = seed)
  state = init_state(Multiscale3{I, bn})

  return Multiscale3{I, bn}(weights, state)
end

max_stars_repo_name: tscode/DCellC

max_issues_repo_name: tscode/DCellC


text:
# ------------ Hierarchical trace ------------ #

struct HierarchicalTrace <: Trace
    calls::Dict{Address, CallSite}
    choices::Dict{Address, ChoiceSite}
    function HierarchicalTrace()
        new(Dict{Address, CallSite}(),
            Dict{Address, ChoiceSite}())
    end
end

Trace() = HierarchicalTrace()

has_top(tr::HierarchicalTrace, addr::T) where T <: Address = haskey(tr.choices, addr)
function has_top(tr::HierarchicalTrace, addr::T) where T <: Tuple
    isempty(addr) && return false
    length(addr) == 1 && return has_top(tr, addr[1])
    has_sub(tr, addr[1]) && has_top(get_sub(tr, addr[1]), addr[2 : end])
end

has_sub(tr::HierarchicalTrace, addr::T) where T <: Address = haskey(tr.calls, addr)
function has_sub(tr::HierarchicalTrace, addr::T) where T <: Tuple
    isempty(addr) && return false
    length(addr) == 1 && return has_sub(tr, addr[1])
    has_sub(tr, addr[1]) && has_sub(get_sub(tr, addr[1]), addr[2 : end])
end

get_sub(tr::HierarchicalTrace, addr::T) where T <: Address = tr.calls[addr]
function get_sub(tr::HierarchicalTrace, addr::T) where T <: Tuple
    isempty(addr) && error("HierarchicalTrace (get_sub): no call in trace at $addr.")
    length(addr) == 1 && return get_sub(tr, addr[1])
    get_sub(get_sub(tr, addr[1]), addr[2 : end])
end

get_top(tr::HierarchicalTrace, addr::T) where T <: Address = tr.choices[addr]
function get_top(tr::HierarchicalTrace, addr::T) where T <: Tuple
    isempty(addr) && error("HierarchicalTrace (get_top): no choice in trace at $addr.")
    length(addr) == 1 && return get_top(tr, addr[1])
    get_top(get_sub(tr, addr[1]), addr[2 : end])
end

function getindex(tr::HierarchicalTrace, addrs...)
    has_sub(tr, addrs) && return get_sub(tr, addrs)
    has_top(tr, addrs) && return get_top(tr, addrs)
    error("HierarchicalTrace (getindex): no choice or call at $addrs.")
end

function Base.haskey(tr::HierarchicalTrace, addr::T) where T <: Address
    has_top(tr, addr) && return true
    has_sub(tr, addr) && return true
    return false
end
function Base.haskey(tr::HierarchicalTrace, addr::T) where T <: Tuple
    isempty(addr) && return false
    length(addr) == 1 && return haskey(tr, addr[1])
    has_sub(tr, addr[1]) && haskey(get_sub(tr, addr[1]), addr[2 : end])
end

# These methods only work for addresses. You should never be adding a call or choice site, except at the current stack level.
add_call!(tr::HierarchicalTrace, addr::T, cs::K) where {T <: Address, K <: CallSite} = tr.calls[addr] = cs
add_choice!(tr::HierarchicalTrace, addr::T, cs::ChoiceSite) where T <: Address = tr.choices[addr] = cs

# ------------- Utility for pretty printing ------------ $

function collect!(par::T, addrs::Vector{Any}, chd::Dict{Any, Any}, tr::HierarchicalTrace, meta) where T <: Tuple
    for (k, v) in tr.choices
        push!(addrs, (par..., k))
        chd[(par..., k)] = v.val
    end
    for (k, v) in tr.calls
        collect!((par..., k), addrs, chd, v.trace, meta)
    end
end

function collect!(addrs::Vector{Any}, chd::Dict{Any, Any}, tr::HierarchicalTrace, meta)
    for (k, v) in tr.choices
        push!(addrs, (k, ))
        chd[(k, )] = v.val
    end
    for (k, v) in tr.calls
        collect!((k, ), addrs, chd, v.trace, meta)
    end
end

# ------------ Hierarchical call site ------------ #

struct HierarchicalCallSite{J, K} <: CallSite
    trace::HierarchicalTrace
    score::Float64
    fn::Function
    args::J
    ret::K
end

has_top(bbcs::HierarchicalCallSite, addr) = has_top(bbcs.trace, addr)

get_top(bbcs::HierarchicalCallSite, addr) = get_top(bbcs.trace, addr)

has_sub(bbcs::HierarchicalCallSite, addr) = has_sub(bbcs.trace, addr)

get_sub(bbcs::HierarchicalCallSite, addr) = get_sub(bbcs.trace, addr)

get_score(bbcs::HierarchicalCallSite) = bbcs.score

getindex(cs::HierarchicalCallSite, addrs...) = getindex(cs.trace, addrs...)

haskey(cs::HierarchicalCallSite, addr) = haskey(cs.trace, addr)

get_ret(cs::HierarchicalCallSite) = cs.ret


text:
## Copyright (c) 2015, Júlio Hoffimann Mendes <juliohm@stanford.edu>
##
## Permission to use, copy, modify, and/or distribute this software for any
## purpose with or without fee is hereby granted, provided that the above
## copyright notice and this permission notice appear in all copies.
##
## THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL WARRANTIES
## WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED WARRANTIES OF
## MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR
## ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL DAMAGES OR ANY DAMAGES
## WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR PROFITS, WHETHER IN AN
## ACTION OF CONTRACT, NEGLIGENCE OR OTHER TORTIOUS ACTION, ARISING OUT OF
## OR IN CONNECTION WITH THE USE OR PERFORMANCE OF THIS SOFTWARE.

function tau_model(events::AbstractVector, D₁::AbstractArray, Dₙ::AbstractArray)
  nevents = length(events)

  # trivial mass function
  nevents == 1 && return ones(1)

  # sources of information: primary + auxiliary₁ + auxiliary₂ + ...
  nsources = 1 + length(Dₙ)

  # distance matrix
  D = zeros(nevents, nsources)
  D[:,1] = D₁[events]
  for j=2:nsources
    D[:,j] = Dₙ[j-1][events]
  end

  # convert distances to ranks
  idx = mapslices(sortperm, D, 1)
  for j=1:nsources
    r = 0; prevdist = -Inf
    for i in view(idx,:,j)
      D[i,j] > prevdist && (r += 1)
      prevdist = D[i,j]
      D[i,j] = r
    end
  end

  # conditional probabilities
  P = nevents - D + 1
  P = broadcast(/, P, sum(P, 1))

  # prior to data all events are equally probable
  x₀ = (1 - 1/nevents) / (1/nevents)

  # assume no redundancy
  X = (1 - P) ./ P
  x = x₀ * prod(X/x₀, 2)

  p = 1 ./ (1 + x)
end


text:
module PlotlyLight

using Random
using JSON3
using EasyConfig
using DefaultApplication
using Scratch
using Printf
using Serialization

export Plot, Config

#-----------------------------------------------------------------------------# init
const plotlyjs = joinpath(@__DIR__, "..", "deps", "plotly-latest.min.js")

plotdir = Ref("")
plot_number = Ref(0)
n_history = Ref(20)  # number of plots to save as history

struct ScratchPlot
    file::String
end

history(; rev=false) = ScratchPlot.(filter!(endswith(".jls"), sort!(readdir(plotdir[]); rev)))

function clean_history!()
    h = history(rev=true)
    n = n_history[]
    if length(h) > n
        rm.(joinpath.(plotdir[], map(x->x.file, h[n+1:end])))
    end
end
set_history!(n::Int) = (n_history[] = n; clean_history!())

function __init__()
    isfile(plotlyjs) || error("Can't find plotly.js.  Try building PlotlyLight again.")
    plotdir[] = @get_scratch!("PlotlyLightHistory")
    h = history()
    if !isempty(h)
        latest_plot = replace(h[end].file, r"(plot_)|(.jls)" => "")
        plot_number[] = parse(Int, latest_plot)
    end
end

#-----------------------------------------------------------------------------# Plot
"""
    Plot(data, layout, config; src=:cdn, class=String[], style="")

A Plotly.js plot with components `data`, `layout`, and `config`.  Each of the three components are
directly converted to JSON.  See the Plotly Javascript docs here: https://plotly.com/javascript/.

- Specify how the Plotly's javascript is loaded via `src`:
    - `:cdn` → load from `"https://cdn.plot.ly/plotly-latest.min.js"` (requires internet access).
    - `:local` → load from `"deps/plotly-latest.min.js"` (downloaded during `Pkg.build("PlotlyLight")`).
    - `:standalone` → Write the `:local` .js file directly into a script tag in the html output.
    - `:none` → `write(io, "text/html"(), plot)` will not add the script tag.
- Classes in `class` will be added to the `<div class="\$(join(class, ' '))">` tag that holds the plot.
- Similar to `class`, `style` will be included via `<div style="\$style">`.

### Example

    Plot(Config(x=1:10, y=randn(10)))

    data = [
        Config(y = randn(10)),
        Config(y = randn(10))
    ]
    layout = Config()
    layout.title.text = "My Title!"
    Plot(data, layout)
"""
mutable struct Plot
    data::Vector{Config}
    layout::Config
    config::Config
    src::Symbol
    class::Vector{String}
    style::String
    scratchfile::String
    function Plot(data = Config[], layout=Config(), config=Config(displaylogo=false); src=:cdn, style="", class=String[])
        @assert src in [:cdn, :local, :standalone]
        n = @sprintf("%010d", plot_number[] += 1)  # keep files sorted alphabetically
        file = joinpath(plotdir[], "plot_$n.jls")
        p = new(data isa Vector ? data : [data], layout, config, src, class, style, file)
        _save_to_scratchspaces(p)
        return p
    end
end

Plot(sp::ScratchPlot) = deserialize(joinpath(plotdir[], sp.file))

"""
    history(i)

Return the `i`-th latest plot.  E.g. `1`
"""
history(i::Integer) = Plot(history(rev=true)[i])

function delete_history!!()
    n = length(rm.(joinpath.(plotdir, map(x->x.file, history()))))
    @info "Deleted PlotlyLight history of $n plots."
end

function _save_to_scratchspaces(o::Plot)
    serialize(touch(o.scratchfile), o)
    clean_history!()
end

function Base.show(io::IO, o::Plot)
    htmlfile = touch(joinpath(plotdir[], "current_plot.html"))
    open(io -> write_html(io, o), htmlfile, "w")
    DefaultApplication.open(htmlfile)
end

#-----------------------------------------------------------------------------# Show text/html
function Base.show(io::IO, ::MIME"text/html", o::Plot)
    o.src in [:cdn, :standalone, :local] || error("`src` must be :cdn, :standalone, :none, or :local")
    id = randstring(20)
    print(io, """

\n""")

    if o.src === :cdn
        write(io, "  ")
    elseif o.src === :standalone
        write(io, "  ")
    elseif o.src === :none
        # No script added
    else # :local
        write(io, "  ")
    end

    print(io, "\n")
end

#-----------------------------------------------------------------------------# HTML
function write_html(io::IO, p::Plot)
    write(io, "<!DOCTYPE html>\n\n  \n  <title>PlotlyLight Viz</title>\n")
    write(io,"\n  \n  \n")
    show(io, MIME"text/html"(), p)
    write(io, "\n  \n  ")
end

end # module


text:
module Essentials

using Reexport: @reexport

export dtype, kind, update!, reset!

function dtype end
function kind end
function update end
function update! end
function reset! end

include("QuantumOperators.jl")
include("QuantumNumbers.jl")
include("Spatials.jl")
include("DegreesOfFreedom.jl")
include("QuantumSystems.jl")
include("Frameworks.jl")

@reexport using .QuantumOperators
@reexport using .QuantumNumbers
@reexport using .Spatials
@reexport using .DegreesOfFreedom
@reexport using .QuantumSystems
@reexport using .Frameworks

end # module


text:
# Autogenerated wrapper script for Graphene_jll for i686-w64-mingw32
export libgraphene

using Glib_jll
JLLWrappers.@generate_wrapper_header("Graphene")
JLLWrappers.@declare_library_product(libgraphene, "libgraphene-1.0-0.dll")
function __init__()
    JLLWrappers.@generate_init_header(Glib_jll)
    JLLWrappers.@init_library_product(
        libgraphene,
        "bin\\libgraphene-1.0-0.dll",
        RTLD_LAZY | RTLD_DEEPBIND,
    )

    JLLWrappers.@generate_init_footer()
end  # __init__()


text:
"""
Outer construct of `FSimBase.Simulator`.
"""
function Simulator(state0, dyn, p=nothing;
        Problem=:ODE,
        solver=nothing,
        kwargs...)
    if Problem == :ODE && solver == nothing
        solver = Tsit5()
    elseif Problem == :Discrete && solver == nothing
        solver = FunctionMap()
    end
    FSimBase.Simulator(state0, dyn, p;
                       Problem=Problem,
                       solver=solver,
                       kwargs...)
end


text:
function apply_BC!(hydro)

	for kk = 1:2
		if kk==1
			nn = 1
			nn1 = 2
		else
			nn = nz
			nn1 = nz-1
		end
		if boundary[kk]=="reflection"
			if geom=="slab"
				#boundary condition for reflection at wall side
				for i = 1:nspec
					hydro.U1D[nn,3*(i-1)+1] = hydro.U1D[nn1,3*(i-1)+1] - lm * hydro.F1D[nn1,3*(i-1)+1] #continuity
					hydro.U1D[nn,3*(i-1)+2] = 0. #momentum
					hydro.U1D[nn,3*(i-1)+3] = hydro.U1D[nn1,3*(i-1)+3] - lm * hydro.F1D[nn1,3*(i-1)+3] #energy
				end
				#now electrons
				hydro.U1D[nn,neqi+1] =  hydro.U1D[nn1,neqi+1] - lm * hydro.F1D[nn1,neqi+1]
			elseif geom=="spherical" #spherical geometry
				#boundary condition for reflection at wall side
				for i = 1:nspec
					hydro.U1D[nn,3*(i-1)+1] = hydro.U1D[nn1,3*(i-1)+1] - lm * hydro.F1D[nn1,3*(i-1)+1] 	-
								2. * dtm_ / r[nz] * phi * hydro.F1D[nz-1,3*(i-1)+1]
					hydro.U1D[nn,3*(i-1)+2] = 0. #momentum
					hydro.U1D[nn,3*(i-1)+3] = hydro.U1D[nn1,3*(i-1)+3] - lm * hydro.F1D[nn1,3*(i-1)+3]  -
								2. * dtm_ / r[nz] * phi * hydro.F1D[nn1,3*(i-1)+3]
				end
				#now electrons
				hydro.U1D[nz,neqi+1] =  hydro.U1D[nz-1,neqi+1]
			end
		elseif boundary[kk]=="open"
			#boundary condition for reflection at wall side
			for i = 1:nspec
				hydro.U1D[nn,3*(i-1)+1] = hydro.rho[nn,i]
				hydro.U1D[nn,3*(i-1)+2] = hydro.rho[nn,i] * hydro.u[nn,i] #momentum
				hydro.U1D[nn,3*(i-1)+3] = hydro.p[nn,i] / (g-1) + 0.5 * hydro.rho[nn,i] * hydro.u[nn,i]^2
			end
			#now electrons
			hydro.U1D[nn,neqi+1] =  hydro.U1D[nn,neqi+1]
		else
			println("wrong boundary condition at boundary #", kk)
		end

		#update variables
		for i = 1:nspec
			hydro.rho[nn,i] = hydro.U1D[nn,3*(i-1)+1]
			hydro.u[nn,i] = hydro.U1D[nn,3*(i-1)+2] / hydro.rho[nn,i]
			hydro.T[nn,i] = (g-1) * mi[i] * ( hydro.U1D[nn,3*(i-1)+3]/hydro.rho[nn,i] - 0.5*hydro.u[nn,i]^2 )
			hydro.p[nn,i] = hydro.rho[nn,i] * hydro.T[nn,i] / mi[i]
		end
	end

	#BC for predictor and corrector
	hydro.U1D_p[1,:] = hydro.U1D[1,:]
	hydro.U1D_c[1,:] = hydro.U1D[1,:]
	hydro.U1D_p[nz,:] = hydro.U1D[nz,:]
	hydro.U1D_c[nz,:] = hydro.U1D[nz,:]

	return hydro
end

max_stars_repo_name: cbellei/JUPL

max_issues_repo_name: cbellei/JUPL

max_forks_repo_name: cbellei/JUPL


text:
module PacanowskiPhilander

using OceanTurb

import .KPP: ∂B∂z

const nsol = 4
@solution U V T S

Base.@kwdef struct Parameters{T} <: AbstractParameters
    Cν₀ :: T = 1e-4
    Cν₁ :: T = 1e-2
    Cκ₀ :: T = 1e-5
    Cκ₁ :: T = 1e-2
    Cc  :: T = 5.0
    Cn  :: T = 2.0
end

struct Model{TS, G, T} <: AbstractModel{TS, G, T}
    @add_standard_model_fields
    parameters :: Parameters{T}
    constants  :: Constants{T}
end

function Model(; N=10, L=1.0,
          grid = UniformGrid(N, L),
     constants = Constants(),
    parameters = Parameters(),
       stepper = :ForwardEuler,
           bcs = BoundaryConditions((ZeroFluxBoundaryConditions() for i=1:nsol)...)
    )

    solution = Solution((CellField(grid) for i=1:nsol)...)
    K = (U=KU, V=KV, T=KT, S=KS)
    R = (U=RU, V=RV, T=RT, S=RS)
    eqn = Equation(K=K, R=R)
    lhs = OceanTurb.build_lhs(solution)

    timestepper = Timestepper(stepper, eqn, solution, lhs)

    return Model(Clock(), grid, timestepper, solution, bcs, parameters, constants)
end

#
# Equation specification
#

function local_richardson(U, V, T, S, g, α, β, i)
    Bz = ∂B∂z(T, S, g, α, β, i)
    S² = ∂z(U, i)^2 + ∂z(V, i)^2

    if S² == 0 && Bz == 0 # Alistair Adcroft's theorem
        return 0
    else
        return Bz / S²
    end
end

local_richardson(m, i) = local_richardson(m.solution.U, m.solution.V, m.solution.T,
                                          m.solution.S, m.constants.g, m.constants.α,
                                          m.constants.β, i)

#
# Equation specification
#

KU(Ri, ν₀, ν₁, c, n) = ν₀ + ν₁ / (1 + c*Ri)^n
KT(Ri, κ₀, κ₁, c, n) = κ₀ + κ₁ / (1 + c*Ri)^(n+1)

KU(m, i) = KU(local_richardson(m, i), m.parameters.Cν₀, m.parameters.Cν₁,
              m.parameters.Cc, m.parameters.Cn)

KT(m, i) = KT(local_richardson(m, i), m.parameters.Cκ₀, m.parameters.Cκ₁,
              m.parameters.Cc, m.parameters.Cn)

const KV = KU
const KS = KT

RU(m, i) =   m.constants.f * m.solution.V[i]
RV(m, i) = - m.constants.f * m.solution.U[i]
RT(m, i) = 0
RS(m, i) = 0

end # module