text:
using JuMP

"""
n = define(;
    numStates::Int = 0,
    numControls::Int = 0,
    X0 = fill(NaN,numStates),
    XF = fill(NaN,numStates),
    XL = fill(NaN,numStates),
    XU = fill(NaN,numStates),
    CL = fill(NaN,numControls),
    CU = fill(NaN,numControls),
    XS = fill(1.0,numStates),
    CS = fill(1.0,numControls)
)::NLOpt
--------------------------------------------------------------------------------------\n
Author: Huckleberry Febbo, Graduate Student, University of Michigan
Date Create: 1/1/2017, Last Modified: 4/13/2018 \n
Citations: \n
----------\n
Initially Influenced by: S. Hughes.  steven.p.hughes@nasa.gov
Source: DecisionVector.m [located here](https://sourceforge.net/p/gmat/git/ci/264a12acad195e6a2467cfdc68abdcee801f73fc/tree/prototype/OptimalControl/LowThrust/@DecisionVector/)
-------------------------------------------------------------------------------------\n
"""
function define(;
    numStates::Int = 0,
    numControls::Int = 0,
    X0 = fill(NaN,numStates),
    XF = fill(NaN,numStates),
    XL = fill(NaN,numStates),
    XU = fill(NaN,numStates),
    CL = fill(NaN,numControls),
    CU = fill(NaN,numControls),
    XS = fill(1.0,numStates),
    CS = fill(1.0,numControls)
)::NLOpt

    n = NLOpt()

    # validate input
    if numControls <= 0
        error("Controls ($numControls) must be > 0")
    end
    if numStates <= 0
        error("States ($numStates) must be > 0")
    end
    if length(X0) != numStates
        error("Length of X0 ($(length(X0))) must match number of states ($numStates)")
    end
    if length(XF) != numStates
        error("Length of XF ($(length(XF))) must match number of states ($numStates)")
    end
    if length(XL) != numStates
        error("Length of XL ($(length(XL))) must match number of states ($numStates)")
    end
    if length(XU) != numStates
        error("Length of XU ($(length(XU))) must match number of states ($numStates)")
    end
    if length(XS) != numStates
        error("Length of XS ($(length(XS))) must match number of states ($numStates)")
    end
    if length(CL) != numControls
        error("Length of CL ($(length(CL))) must match number of controls ($numControls)")
    end
    if length(CU) != numControls
        error("Length of CU ($(length(CU))) must match number of controls ($numControls)")
    end
    if length(CS) != numControls
        error("Length of CS ($(length(CS))) must match number of controls ($numControls)")
    end

    n.ocp.state     = initState(numStates)
    n.ocp.control   = initControl(numControls)
    n.ocp.X0        = X0
    n.ocp.X0_tol    = fill(NaN, size(X0))
    n.ocp.XF        = XF
    n.ocp.XF_tol    = fill(NaN, size(XF))
    n.ocp.XL        = XL
    n.ocp.XU        = XU
    n.ocp.CL        = CL
    n.ocp.CU        = CU
    n.ocp.XS        = XS
    n.ocp.CS        = CS
    n.f.ocp.defined = true

    return n

end

"""
defineSolver!(n::NLOpt, kw::Dict{Symbol, Symbol})
--------------------------------------------------------------------------------------\n
Author: Huckleberry Febbo, Graduate Student, University of Michigan
Date Create: 2/9/2017, Last Modified: 12/02/2019 \n
Last Modified by: Alexander Buck
-------------------------------------------------------------------------------------\n
"""
function defineSolver!(n, kw)

    if typeof(kw)!=Dict{Symbol,Symbol}
      kw = Dict(kw)
    end

    # Default solver name is :Ipopt
    if haskey(kw, :name)
        n.s.ocp.solver.name = kw[:name]
    else
        @warn "No solver specified: using default solver (IPOPT)"
        n.s.ocp.solver.name = :Ipopt
    end

    # Import solver into Julia environment
    try_import(n.s.ocp.solver.name)

    # Check if user provided MPC defaults
    mpc_defaults = get(kw, :mpc_defaults, false)

    # Setup MPC Solver Settings (if they exist)
    if mpc_defaults
        if n.s.ocp.solver.name == :Ipopt
            n.s.ocp.solver.settings = NLOptControl._Ipopt_MPC
        else
            error("Solver $(n.s.ocp.solver.name) not defined.")
        end
    else
        # NOTE this should already have been done by default, but could get messed up is user is playing with options
        if n.s.ocp.solver.name == :Ipopt
            n.s.ocp.solver.settings = NLOptControl._Ipopt_defaults;
        else
            error("Solver $(n.s.ocp.solver.name) not defined.")
        end
    end

    # Modify additional defaults individually
    for (key, value) in kw
        if haskey(n.s.ocp.solver.settings, key)
            n.s.ocp.solver.settings[key] = value
        elseif key != :name && key != :mpc_defaults # ignore the name and default settings option # TODO could remove them from the Dict
            error("Unknown key: $kw for $(n.s.ocp.solver.name) used in defineSolver!()")
        end
    end

    # Set solver
    if n.s.ocp.solver.name == :Ipopt
        setsolver(
            n.ocp.mdl,
            Ipopt.IpoptSolver(;
                max_cpu_time               = 10000.0, #n.s.ocp.solver.settings[:max_cpu_time],
                print_level                = n.s.ocp.solver.settings[:print_level],
                warm_start_init_point      = n.s.ocp.solver.settings[:warm_start_init_point],
                max_iter                   = 100000, #.s.ocp.solver.settings[:max_iter],
                tol                        = n.s.ocp.solver.settings[:tol],
                dual_inf_tol               = n.s.ocp.solver.settings[:dual_inf_tol],
                constr_viol_tol            = n.s.ocp.solver.settings[:constr_viol_tol],
                compl_inf_tol              = n.s.ocp.solver.settings[:compl_inf_tol],
                acceptable_tol             = n.s.ocp.solver.settings[:acceptable_tol],
                acceptable_constr_viol_tol = n.s.ocp.solver.settings[:acceptable_constr_viol_tol],
                acceptable_dual_inf_tol    = n.s.ocp.solver.settings[:acceptable_dual_inf_tol],
                acceptable_compl_inf_tol   = n.s.ocp.solver.settings[:acceptable_compl_inf_tol],
                acceptable_obj_change_tol  = n.s.ocp.solver.settings[:acceptable_obj_change_tol],
                diverging_iterates_tol     = n.s.ocp.solver.settings[:diverging_iterates_tol]
            )
        )
    else
        error("Solver $(n.s.ocp.solver.name) not defined")
    end

    return nothing

end

"""
OCPdef!(n::NLOpt)
--------------------------------------------------------------------------------------\n
Author: Huckleberry Febbo, Graduate Student, University of Michigan
Date Create: 1/14/2017, Last Modified: 12/06/2019 \n
--------------------------------------------------------------------------------------\n
"""
function OCPdef!(n::NLOpt{T}) where { T <: Number }

    # State variables
    @variable(n.ocp.mdl, x[1:n.ocp.state.pts, 1:n.ocp.state.num])
    #n.r.ocp.x = [ @NLexpression(n.ocp.mdl, [j=1:n.ocp.state.pts], n.ocp.XS[st] * x[j,st] ) for st in 1:n.ocp.state.num ]
    n.r.ocp.x = Matrix{Any}(undef, n.ocp.state.pts,n.ocp.state.num)
    for st in 1:n.ocp.state.num
        n.r.ocp.x[:,st] = @NLexpression(n.ocp.mdl, [j=1:n.ocp.state.pts], n.ocp.XS[st]*x[j,st])
    end

    for st in 1:n.ocp.state.num

        # Lower state constraints
        if !isnan(n.ocp.XL[st])
            if !n.ocp.mXL[st]
                for j in 1:n.ocp.state.pts
                    setlowerbound(x[j,st], n.ocp.XL[st]/n.ocp.XS[st])
                end
            else
                for j in 1:n.ocp.state.pts
                    setlowerbound(x[j,st], n.ocp.XL_var[st][j]/n.ocp.XS[st])
                end
            end
        end

        # Upper state constraint
        if !isnan(n.ocp.XU[st])
            if !n.ocp.mXU[st]
                for j in 1:n.ocp.state.pts
                    setupperbound(x[j,st], n.ocp.XU[st]/n.ocp.XS[st])
                end
            else
                for j in 1:n.ocp.state.pts
                    setlowerbound(x[j,st], n.ocp.XU_var[st,j]/n.ocp.XS[st])
                end
            end
        end
    end

    # control variables
    @variable(n.ocp.mdl,u[1:n.ocp.control.pts,1:n.ocp.control.num])
    #n.r.ocp.u = [ @NLexpression(n.ocp.mdl, [j=1:n.ocp.control.pts], n.ocp.CS[ctr]*u[j,ctr]) for ctr in 1:n.ocp.control.num ]
    n.r.ocp.u = Matrix{Any}(undef,n.ocp.control.pts,n.ocp.control.num)
    for ctr in 1:n.ocp.control.num
        n.r.ocp.u[:,ctr] = @NLexpression(n.ocp.mdl, [j=1:n.ocp.control.pts], n.ocp.CS[ctr]*u[j,ctr])
    end

    for ctr in 1:n.ocp.control.num
        if !isnan.(n.ocp.CL[ctr])
            for j in 1:n.ocp.control.pts
                setlowerbound(u[j,ctr], n.ocp.CL[ctr]/n.ocp.CS[ctr])
            end
        end
        if !isnan.(n.ocp.CU[ctr])
            for j in 1:n.ocp.control.pts
                setupperbound(u[j,ctr], n.ocp.CU[ctr]/n.ocp.CS[ctr])
            end
        end
    end

    # boundary constraints
    if n.s.ocp.x0slackVariables   # create handles for constraining the entire initial state
        n.r.ocp.x0Con = Array{Any}(undef, n.ocp.state.num,2) # this is so they can be easily reference when doing MPC
    else
        n.r.ocp.x0Con = []
    end

    if n.s.ocp.xFslackVariables # create handles for constraining the entire final state
        n.r.ocp.xfCon = Array{Any}(undef, n.ocp.state.num,2) # this is so they can be easily reference when doing MPC
    else
        n.r.ocp.xfCon = []
    end

    if n.s.ocp.x0slackVariables # currently adding slack variables for all initial states even if there are no constraints on them
        @variable(n.ocp.mdl, x0s[1:n.ocp.state.num])
        n.ocp.x0s = x0s
        n.r.ocp.x0sCon = Array{Any}(undef, n.ocp.state.num)
    end
    if n.s.ocp.xFslackVariables # currently adding slack variables for all final states even if there are no constraints on them
        @variable(n.ocp.mdl, xFs[1:n.ocp.state.num])
        n.ocp.xFs = xFs
        n.r.ocp.xfsCon = Array{Any}(undef, n.ocp.state.num)
    end

    for st in 1:n.ocp.state.num
        if !isnan(n.ocp.X0[st]) # could have a bool for this
            if n.s.ocp.x0slackVariables
                n.r.ocp.x0Con[st,1] = @constraint(n.ocp.mdl, x[1,st]*n.ocp.XS[st] -n.ocp.x0s[st] <=  n.ocp.X0[st])
                n.r.ocp.x0Con[st,2] = @constraint(n.ocp.mdl,-x[1,st]*n.ocp.XS[st] -n.ocp.x0s[st] <= -n.ocp.X0[st])
                if !isnan(n.ocp.X0_tol[st])
                    n.r.ocp.x0sCon[st] = @constraint(n.ocp.mdl, n.ocp.x0s[st] <= n.ocp.X0_tol[st])
                end
            else
                n.r.ocp.x0Con = [ n.r.ocp.x0Con ; @constraint(n.ocp.mdl, x[1,st]*n.ocp.XS[st] == n.ocp.X0[st]) ]
            end
        end
        if !isnan(n.ocp.XF[st])
            if n.s.ocp.xFslackVariables
                n.r.ocp.xfCon[st,1] = @constraint(n.ocp.mdl, x[end,st]*n.ocp.XS[st] -n.ocp.xFs[st] <=  n.ocp.XF[st] )
                n.r.ocp.xfCon[st,2] = @constraint(n.ocp.mdl,-x[end,st]*n.ocp.XS[st] -n.ocp.xFs[st] <= -n.ocp.XF[st] )
                if !isnan(n.ocp.XF_tol[st])
                    n.r.ocp.xfsCon[st] = @constraint(n.ocp.mdl, n.ocp.xFs[st] <= n.ocp.XF_tol[st])
                end
            else
                n.r.ocp.xfCon = [n.r.ocp.xfCon ; @constraint(n.ocp.mdl, x[end,st]*n.ocp.XS[st] == n.ocp.XF[st]) ]
            end
        end
    end

    @NLparameter(n.ocp.mdl,t0_param == 0.0);   # for now we just start at zero
    n.ocp.t0 = t0_param  # NOTE consider making a constraint that t0 < tf

    if n.s.ocp.integrationMethod == :ps
        n.r.ocp.dynCon = [Array{Any}(undef, n.ocp.Nck[int],n.ocp.state.num) for int in 1:n.ocp.Ni]
        dynamics_expr = [Array{Any}(undef, n.ocp.Nck[int],n.ocp.state.num) for int in 1:n.ocp.Ni]

        if n.s.ocp.finalTimeDV
            @variable(n.ocp.mdl, n.s.ocp.tfMin <= tf <=  n.s.ocp.tfMax)
            n.ocp.tf = tf
            create_tV!(n)          # make a time vector
        end

        for int in 1:n.ocp.Ni
            x_int,u_int = intervals(n,int,n.r.ocp.x,n.r.ocp.u)

            # dynamics
            L = size(x_int)[1]-1;
            dx = Matrix{Any}(undef, L,n.ocp.state.num)
            for st in 1:n.ocp.state.num
                dx[:,st] = DiffEq(n,x_int,u_int,L,st)
            end

            for st in 1:n.ocp.state.num # TODO consider multiplying X*D to reduce computations (i.e. remove this for loop for the states)
                if n.s.ocp.integrationScheme == :lgrExplicit
                    dynamics_expr[int][:,st] = @NLexpression(n.ocp.mdl, [j in 1:n.ocp.Nck[int]], sum(n.ocp.DMatrix[int][j,i]*x_int[i,st] for i in 1:n.ocp.Nck[int]+1) - ((n.ocp.tf)/2)*dx[j,st]  )
                elseif n.s.ocp.integrationScheme==:lgrImplicit
                    dynamics_expr[int][:,st] = @NLexpression(n.ocp.mdl, [j in 1:n.ocp.Nck[int]], x_int[j+1,st] - x_int[1,st] - ((n.ocp.tf)/2)*sum(n.ocp.IMatrix[int][j,i]*dx[i,st] for i in 1:n.ocp.Nck[int]) )
                end
                for j in 1:n.ocp.Nck[int]
                    n.r.ocp.dynCon[int][j,st] = @NLconstraint(n.ocp.mdl, 0. == dynamics_expr[int][j,st])
                end
            end
            # additional constraints
            for num in 1:length(n.ocp.NLcon)
                ch = addCon(n,x_int,u_int,L,num)
                newConstraint!(n,ch,Symbol(string("ch",num))) # TODO could let the user name these
            end

        end

    elseif n.s.ocp.integrationMethod == :tm
        n.r.ocp.dynCon = Array{Any}(undef,n.ocp.N,n.ocp.state.num)

        if n.s.ocp.finalTimeDV
            @variable(n.ocp.mdl, n.s.ocp.tfMin <= tf <= n.s.ocp.tfMax)
            n.ocp.tf = tf
            create_tV!(n) # make a time vector
        end
        n.ocp.dt = n.ocp.tf / n.ocp.N * ones(n.ocp.N,)
        L = size(n.r.ocp.x)[1]
        #dx = [ DiffEq(n, n.r.ocp.x, n.r.ocp.u, L, st) for st in 1:n.ocp.state.num ]
        dx = Array{Any}(undef,L,n.ocp.state.num)
        for st in 1:n.ocp.state.num
          dx[:,st] = DiffEq(n,n.r.ocp.x,n.r.ocp.u,L,st)
        end
        if n.s.ocp.integrationScheme==:bkwEuler
          for st in 1:n.ocp.state.num
            n.r.ocp.dynCon[:,st] = @NLconstraint(n.ocp.mdl, [j in 1:n.ocp.N], n.r.ocp.x[j+1,st] - n.r.ocp.x[j,st] ==  dx[j+1,st]*n.ocp.tf/(n.ocp.N) );
          end
        elseif n.s.ocp.integrationScheme==:trapezoidal
          for st in 1:n.ocp.state.num
            n.r.ocp.dynCon[:,st] = @NLconstraint(n.ocp.mdl, [j in 1:n.ocp.N], n.r.ocp.x[j+1,st] - n.r.ocp.x[j,st] == 0.5*(dx[j,st] + dx[j+1,st])*n.ocp.tf/(n.ocp.N) )
          end
        end

        # additional constraints
        for num in 1:length(n.ocp.NLcon)
            ch = addCon(n,n.r.ocp.x,n.r.ocp.u,L,num)
            newConstraint!(n,ch,Symbol(string("ch",num))) # TODO could let the user name these
        end
    end

    # If evaluate costates, make sure using Pseudospectral Method
    if n.s.ocp.evalCostates

        if n.s.ocp.integrationMethod != :ps
            error("Costates are only implmented for :ps methods")
        else
            n.s.ocp.evalConstraints = true
        end

    end

    # @info "OCPdef!: size of n.r.ocp.x0Con  : $(size(  n.r.ocp.x0Con ))"
    # @info "OCPdef!: type of n.r.ocp.x0Con  : $(typeof(n.r.ocp.x0Con ))"
    # @info "OCPdef!: size of n.r.ocp.xfCon  : $(size(  n.r.ocp.xfCon ))"
    # @info "OCPdef!: type of n.r.ocp.xfCon  : $(typeof(n.r.ocp.xfCon ))"
    # @info "OCPdef!: size of n.r.ocp.dynCon : $(size(  n.r.ocp.dynCon))"
    # @info "OCPdef!: type of n.r.ocp.dynCon : $(typeof(n.r.ocp.dynCon))"

    # save constraint data
    newConstraint!(n, n.r.ocp.x0Con , :x0_con )
    newConstraint!(n, n.r.ocp.xfCon , :xf_con )
    newConstraint!(n, n.r.ocp.dynCon, :dyn_con)

    # save the current working directory for navigation purposes
    n.r.mainDir = pwd()

    return nothing

end

"""
configure!(n::NLOpt; kwargs...)
--------------------------------------------------------------------------------------\n
Author: Huckleberry Febbo, Graduate Student, University of Michigan
Date Create: 1/1/2017, Last Modified: 12/06/2019 \n
-------------------------------------------------------------------------------------\n
"""
function configure!(n::NLOpt{T}; kwargs... ) where { T <: Number }

    # Gather keyword arguments
    kw = Dict(kwargs)

    # Determine whether final time is a design variable (:finalTimeDV)
    # defaults to not a design variable
    if !haskey(kw,:finalTimeDV);n.s.ocp.finalTimeDV=false;
    else; n.s.ocp.finalTimeDV=get(kw,:finalTimeDV,0);
    end

    # Determine final time (:tf)
    if !haskey(kw,:tf) && !n.s.ocp.finalTimeDV
      error("\n If the final is not a design variable pass it as: (:tf=>Float64(some #)) \n
          If the final time is a design variable, indicate that as: (:finalTimeDV=>true)\n")
    elseif haskey(kw,:tf) && !n.s.ocp.finalTimeDV
      n.ocp.tf = get(kw,:tf,0)
    end

    # Initial State Slack Variables (:x0slackVariables)
    # TODO: figure out a better default value than "true/false"
    # TODO: check to see if user is using tolerances and passed false to configure
    n.s.ocp.x0slackVariables = get(kw, :x0slackVariables, !any(isnan.(n.ocp.X0_tol)))

    # Final State Slack Variables (:xFslackVariables)
    # TODO: figure out a better default value than "true/false"
    n.s.ocp.xFslackVariables = get(kw, :xFslackVariables, !any(isnan.(n.ocp.XF_tol)))

    # Integration Scheme (:integrationScheme)
    # Default = :lgrExplicit
    n.s.ocp.integrationScheme = get(kw, :integrationScheme, :lgrExplicit)

    # Integration Method (:integrationMethod)
    if in( n.s.ocp.integrationScheme,  NLOptControl.pseudoSpectralMethods )

        # Use Pseudospectral Integration Method
        n.s.ocp.integrationMethod = :ps

        if haskey(kw, :N)
            error(":N is not an appropriate keyword argument for :ps method :$(n.s.ocp.integrationScheme), use :Nck")
        else

            n.ocp.Nck = get(kw, :Nck, [10 , 10 , 10 , 10] )

            if any(n.ocp.Nck .< 0)
                error(":Nck must be > 0")
            end

            # Determine number of points
            n.ocp.Ni            = length(n.ocp.Nck)
            n.ocp.state.pts     = sum(n.ocp.Nck) + 1
            n.ocp.control.pts   = sum(n.ocp.Nck)
            n.ocp.Nck_full      = [ 0 ; cumsum(n.ocp.Nck .+ 1) ]
            n.ocp.Nck_cum       = [ 0 ; cumsum(n.ocp.Nck) ]

            # Initialize node data
            taus_and_weights    = [ gaussradau(n.ocp.Nck[int]) for int in 1:n.ocp.Ni ]
            n.ocp.tau           = [taus_and_weights[int][1] for int in 1:n.ocp.Ni]
            n.ocp.w             = [taus_and_weights[int][2] for int in 1:n.ocp.Ni]

            # Create Intervals
            createIntervals!(n)

            #? Create DMatrix
            DMatrix!(n)

        end

    elseif in( n.s.ocp.integrationScheme, NLOptControl.trapezoidalMethods )

        # Use Trapezoidal Method
        n.s.ocp.integrationMethod = :tm

        if haskey(kw, :Nck)
            error(":Nck is not an appropriate keyword argument for :tm method :$(n.s.ocp.integrationScheme), use :N")
        else

            # Default number of points is 100
            n.ocp.N = get(kw, :N, 100)

            # Set number of integration points for state and control
            n.ocp.state.pts = n.ocp.N + 1
            n.ocp.control.pts = n.ocp.N + 1

        end

    else
        error("The :integrationScheme that you specified ($(n.s.ocp.integrationScheme)) is not currently implemeted.")
    end
    n.ocp.mXL = falses(n.ocp.state.num)
    n.ocp.mXU = falses(n.ocp.state.num)
    #n.ocp.XL_var = Vector{Vector{JuMP.JuMPTypes}}([ [ @variable(JuMP.Model(), tmp) for i = 1:n.ocp.state.num ] for j = 1:n.ocp.state.pts] )
    #n.ocp.XU_var = Vector{Vector{JuMP.JuMPTypes}}([ [ @variable(JuMP.Model(), tmp) for i = 1:n.ocp.state.num ] for j = 1:n.ocp.state.pts] )
    n.ocp.XL_var = Matrix{Float64}(undef,n.ocp.state.num,n.ocp.state.pts)
    n.ocp.XU_var = Matrix{Float64}(undef,n.ocp.state.num,n.ocp.state.pts)

    # Define Solver Settings
    defineSolver!(n, Dict(get(kw, :solverSettings, (:name => :Ipopt) )))

    # Definie Optimal Control Problem
    OCPdef!(n)

    return nothing

end


text:
lines = readlines("input.txt")
re = r"(\w+)"

function main()
    #=
    Solve using frequency analysis (i.e.: for digits 0 to 9, we see segment 'a' 8 times)
    For 1: we need c and f, i.e. 8 + 9 = 17. for 7, a + c + f = 8 + 8 + 9 = 25
    This works, because the distribution for the number of times each segment lights
    up for digits 0 to 9 is unique
    =#
    segmentmap = Dict(
        0 => collect("abcefg"),
        1 => collect("cf"),
        2 => collect("acdeg"),
        3 => collect("acdfg"),
        4 => collect("bcdf"),
        5 => collect("abdfg"),
        6 => collect("abdefg"),
        7 => collect("acf"),
        8 => collect("abcdefg"),
        9 => collect("abcdfg")
    )
    freq = X -> Dict([(x, count(==(x), X)) for x = unique(X)])
    X = freq(reduce(vcat, values(segmentmap)))
    numbers = Dict([(i,sum([X[c] for c = set])) for (i,set) = pairs(segmentmap)])
    inverse = Dict(value => key for (key, value) = numbers)
    sumdigits = N -> sum(n * 10^e for (n,e) = zip(N, length(N)-1:-1:0))

    Σ = 0
    for line = lines
        input, output = split(line, '|')
        signals = [m.match for m = eachmatch(re, input)]
        X = freq(reduce(vcat, collect.(signals)))

        Y = [m.match for m = eachmatch(re, output)]
        Σ_Y = map(s -> sum([X[c] for c = collect(s)]), Y)
        N = map(n -> inverse[n], Σ_Y)
        Σ += sumdigits(N)
        println(sumdigits(N))
    end
    println("====")
    println(Σ)
end

main()

ext: jl


text:
using Documenter, GtkReactive

makedocs(format   = :html,
         sitename = "GtkReactive",
         pages    = ["index.md", "controls.md", "drawing.md", "zoom_pan.md", "reference.md"]
         )

deploydocs(repo   = "github.com/JuliaGizmos/GtkReactive.jl",
           target = "build",
           julia  = "1.0",
           deps   = nothing,
           make   = nothing
           )


text:
module FooRouter

using HTTP
const r = HTTP.Router()
f = HTTP.Handlers.RequestHandlerFunction((req) -> HTTP.Response(200))
HTTP.@register(r, "/test", f)

end # module


text:
using ExecutableSpecifications: @given, @when, @then, @expect

@given "that the context variable :x has value 1" begin
    context[:x] = 1
end

@when "the context variable :x is compared with 1" begin
    v = context[:x] == 1
    context[:result] = v
end

@then "the comparison is true" begin
    @expect context[:result]
end

@when "this step fails" begin
    @expect 1 == 2
end

@then "this step is skipped" begin
    @assert false "This code is never executed, because the previous step failed."
end

@when "this step throws an exception, the result is \"Unknown exception\"" begin
    throw(ErrorException("This is an unknown exception"))
end

# This step is commented to show what happens when a step definition is not found.
#@given "a step that has no corresponding step definition in steps/spec.jl" begin
#
#end

@when "it is executed" begin
    @assert false "This step should be skipped, because the previous failed."
end

@then "it fails with error \"No matching step\"" begin
    @assert false "This step should be skipped, because a previous step failed."
end


text:
@testset "execution" begin

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

dummy() = return

@testset "@cuda" begin

@test_throws UndefVarError @cuda undefined()
@test_throws MethodError @cuda dummy(1)

@testset "low-level interface" begin
    k = cufunction(dummy)
    k()
    k(; threads=1)

    CUDAnative.version(k)
    CUDAnative.memory(k)
    CUDAnative.registers(k)
    CUDAnative.maxthreads(k)
end

@testset "launch configuration" begin
    @cuda dummy()

    @cuda threads=1 dummy()
    @cuda threads=(1,1) dummy()
    @cuda threads=(1,1,1) dummy()

    @cuda blocks=1 dummy()
    @cuda blocks=(1,1) dummy()
    @cuda blocks=(1,1,1) dummy()

    @cuda config=(kernel)->() dummy()
    @cuda config=(kernel)->(threads=1,) dummy()
    @cuda config=(kernel)->(blocks=1,) dummy()
    @cuda config=(kernel)->(shmem=0,) dummy()
end

@testset "compilation params" begin
    @cuda dummy()

    @test_throws CuError @cuda threads=2 maxthreads=1 dummy()
    @cuda threads=2 dummy()
end

@testset "reflection" begin
    CUDAnative.code_lowered(dummy, Tuple{})
    CUDAnative.code_typed(dummy, Tuple{})
    CUDAnative.code_warntype(devnull, dummy, Tuple{})
    CUDAnative.code_llvm(devnull, dummy, Tuple{})
    CUDAnative.code_ptx(devnull, dummy, Tuple{})
    CUDAnative.code_sass(devnull, dummy, Tuple{})

    @device_code_lowered @cuda dummy()
    @device_code_typed @cuda dummy()
    @device_code_warntype io=devnull @cuda dummy()
    @device_code_llvm io=devnull @cuda dummy()
    @device_code_ptx io=devnull @cuda dummy()
    @device_code_sass io=devnull @cuda dummy()

    mktempdir() do dir
        @device_code dir=dir @cuda dummy()
    end

    @test_throws ErrorException @device_code_lowered nothing

    # make sure kernel name aliases are preserved in the generated code
    @test occursin("ptxcall_dummy", sprint(io->(@device_code_llvm io=io @cuda dummy())))
    @test occursin("ptxcall_dummy", sprint(io->(@device_code_ptx io=io @cuda dummy())))
    @test occursin("ptxcall_dummy", sprint(io->(@device_code_sass io=io @cuda dummy())))

    # make sure invalid kernels can be partially reflected upon
    let
        invalid_kernel() = throw()
        @test_throws CUDAnative.KernelError @cuda invalid_kernel()
        @test_throws CUDAnative.KernelError @grab_output @device_code_warntype @cuda invalid_kernel()
        out, err = @grab_output begin
            try
                @device_code_warntype @cuda invalid_kernel()
            catch
            end
        end
        @test occursin("Body::Union{}", err)
    end

    let
        range_kernel() = (0.0:0.1:100.0; nothing)

        @test_throws CUDAnative.InvalidIRError @cuda range_kernel()
    end

    # set name of kernel
    @test occursin("ptxcall_mykernel", sprint(io->(@device_code_llvm io=io begin
        k = cufunction(dummy, name="mykernel")
        k()
    end)))
end

@testset "shared memory" begin
    @cuda shmem=1 dummy()
end

@testset "streams" begin
    s = CuStream()
    @cuda stream=s dummy()
end

@testset "external kernels" begin
    @eval module KernelModule
        export external_dummy
        external_dummy() = return
    end
    import ...KernelModule
    @cuda KernelModule.external_dummy()
    @eval begin
        using ...KernelModule
        @cuda external_dummy()
    end

    @eval module WrapperModule
        using CUDAnative
        @eval dummy() = return
        wrapper() = @cuda dummy()
    end
    WrapperModule.wrapper()
end

@testset "calling device function" begin
    @noinline child(i) = sink(i)
    function parent()
        child(1)
        return
    end

    @cuda parent()
end

end

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

@testset "argument passing" begin

dims = (16, 16)
len = prod(dims)

@testset "manually allocated" begin
    function kernel(input, output)
        i = (blockIdx().x-1) * blockDim().x + threadIdx().x

        val = unsafe_load(input, i)
        unsafe_store!(output, val, i)

        return
    end

    input = round.(rand(Float32, dims) * 100)
    output = similar(input)

    input_dev = CuArray(input)
    output_dev = CuArray(output)

    @cuda threads=len kernel(pointer(input_dev), pointer(output_dev))
    @test input ≈ Array(output_dev)
end

@testset "scalar through single-value array" begin
    function kernel(a, x)
        i = (blockIdx().x-1) * blockDim().x + threadIdx().x
        max = gridDim().x * blockDim().x
        if i == max
            _val = unsafe_load(a, i)
            unsafe_store!(x, _val)
        end
        return
    end

    arr = round.(rand(Float32, dims) * 100)
    val = [0f0]

    arr_dev = CuArray(arr)
    val_dev = CuArray(val)

    @cuda threads=len kernel(pointer(arr_dev), pointer(val_dev))
    @test arr[dims...] ≈ Array(val_dev)[1]
end

@testset "scalar through single-value array, using device function" begin
    function child(a, i)
        return unsafe_load(a, i)
    end
    @noinline function parent(a, x)
        i = (blockIdx().x-1) * blockDim().x + threadIdx().x
        max = gridDim().x * blockDim().x
        if i == max
            _val = child(a, i)
            unsafe_store!(x, _val)
        end
        return
    end

    arr = round.(rand(Float32, dims) * 100)
    val = [0f0]

    arr_dev = CuArray(arr)
    val_dev = CuArray(val)

    @cuda threads=len parent(pointer(arr_dev), pointer(val_dev))
    @test arr[dims...] ≈ Array(val_dev)[1]
end

@testset "tuples" begin
    # issue #7: tuples not passed by pointer

    function kernel(keeps, out)
        if keeps[1]
            unsafe_store!(out, 1)
        else
            unsafe_store!(out, 2)
        end
        return
    end

    keeps = (true,)
    d_out = CuArray(zeros(Int))

    @cuda kernel(keeps, pointer(d_out))
    @test Array(d_out)[] == 1
end

@testset "ghost function parameters" begin
    # bug: ghost type function parameters are elided by the compiler

    len = 60
    a = rand(Float32, len)
    b = rand(Float32, len)
    c = similar(a)

    d_a = CuArray(a)
    d_b = CuArray(b)
    d_c = CuArray(c)

    @eval struct ExecGhost end

    function kernel(ghost, a, b, c)
        i = (blockIdx().x-1) * blockDim().x + threadIdx().x
        unsafe_store!(c, unsafe_load(a,i)+unsafe_load(b,i), i)
        return
    end
    @cuda threads=len kernel(ExecGhost(), pointer(d_a), pointer(d_b), pointer(d_c))
    @test a+b == Array(d_c)

    # bug: ghost type function parameters confused aggregate type rewriting

    function kernel(ghost, out, aggregate)
        i = (blockIdx().x-1) * blockDim().x + threadIdx().x
        unsafe_store!(out, aggregate[1], i)
        return
    end
    @cuda threads=len kernel(ExecGhost(), pointer(d_c), (42,))

    @test all(val->val==42, Array(d_c))
end

@testset "immutables" begin
    # issue #15: immutables not passed by pointer

    function kernel(ptr, b)
        unsafe_store!(ptr, imag(b))
        return
    end

    arr = CuArray(zeros(Float32))
    x = ComplexF32(2,2)

    @cuda kernel(pointer(arr), x)
    @test Array(arr)[] == imag(x)
end

@testset "automatic recompilation" begin
    arr = CuArray(zeros(Int))

    function kernel(ptr)
        unsafe_store!(ptr, 1)
        return
    end

    @cuda kernel(pointer(arr))
    @test Array(arr)[] == 1

    function kernel(ptr)
        unsafe_store!(ptr, 2)
        return
    end

    @cuda kernel(pointer(arr))
    @test Array(arr)[] == 2
end

@testset "non-isbits arguments" begin
    function kernel1(T, i)
        sink(i)
        return
    end
    @cuda kernel1(Int, 1)

    function kernel2(T, i)
        sink(unsafe_trunc(T,i))
        return
    end
    @cuda kernel2(Int, 1.)
end

@testset "splatting" begin
    function kernel(out, a, b)
        unsafe_store!(out, a+b)
        return
    end

    out = [0]
    out_dev = CuArray(out)
    out_ptr = pointer(out_dev)

    @cuda kernel(out_ptr, 1, 2)
    @test Array(out_dev)[1] == 3

    all_splat = (out_ptr, 3, 4)
    @cuda kernel(all_splat...)
    @test Array(out_dev)[1] == 7

    partial_splat = (5, 6)
    @cuda kernel(out_ptr, partial_splat...)
    @test Array(out_dev)[1] == 11
end

@testset "object invoke" begin
    # this mimics what is generated by closure conversion

    @eval struct KernelObject{T} <: Function
        val::T
    end

    function (self::KernelObject)(a)
        unsafe_store!(a, self.val)
        return
    end

    function outer(a, val)
       inner = KernelObject(val)
       @cuda inner(a)
    end

    a = [1.]
    a_dev = CuArray(a)

    outer(pointer(a_dev), 2.)

    @test Array(a_dev) ≈ [2.]
end

@testset "closures" begin
    function outer(a_dev, val)
       function inner(a)
            # captures `val`
            unsafe_store!(a, val)
            return
       end
       @cuda inner(pointer(a_dev))
    end

    a = [1.]
    a_dev = CuArray(a)

    outer(a_dev, 2.)

    @test Array(a_dev) ≈ [2.]
end

@testset "conversions" begin
    @eval struct Host   end
    @eval struct Device end

    Adapt.adapt_storage(::CUDAnative.Adaptor, a::Host) = Device()

    Base.convert(::Type{Int}, ::Host)   = 1
    Base.convert(::Type{Int}, ::Device) = 2

    out = [0]

    # convert arguments
    out_dev = CuArray(out)
    let arg = Host()
        @test Array(out_dev) ≈ [0]
        function kernel(arg, out)
            unsafe_store!(out, convert(Int, arg))
            return
        end
        @cuda kernel(arg, pointer(out_dev))
        @test Array(out_dev) ≈ [2]
    end

    # convert tuples
    out_dev = CuArray(out)
    let arg = (Host(),)
        @test Array(out_dev) ≈ [0]
        function kernel(arg, out)
            unsafe_store!(out, convert(Int, arg[1]))
            return
        end
        @cuda kernel(arg, pointer(out_dev))
        @test Array(out_dev) ≈ [2]
    end

    # convert named tuples
    out_dev = CuArray(out)
    let arg = (a=Host(),)
        @test Array(out_dev) ≈ [0]
        function kernel(arg, out)
            unsafe_store!(out, convert(Int, arg.a))
            return
        end
        @cuda kernel(arg, pointer(out_dev))
        @test Array(out_dev) ≈ [2]
    end

    # don't convert structs
    out_dev = CuArray(out)
    @eval struct Nested
        a::Host
    end
    let arg = Nested(Host())
        @test Array(out_dev) ≈ [0]
        function kernel(arg, out)
            unsafe_store!(out, convert(Int, arg.a))
            return
        end
        @cuda kernel(arg, pointer(out_dev))
        @test Array(out_dev) ≈ [1]
    end
end

@testset "argument count" begin
    val = [0]
    val_dev = CuArray(val)
    cuda_ptr = pointer(val_dev)
    ptr = CUDAnative.DevicePtr{Int}(cuda_ptr)
    for i in (1, 10, 20, 35)
        variables = ('a':'z'..., 'A':'Z'...)
        params = [Symbol(variables[j]) for j in 1:i]
        # generate a kernel
        body = quote
            function kernel($(params...))
                unsafe_store!($ptr, $(Expr(:call, :+, params...)))
                return
            end
        end
        eval(body)
        args = [j for j in 1:i]
        call = Expr(:call, :kernel, args...)
        cudacall = :(@cuda $call)
        eval(cudacall)
        @test Array(val_dev)[1] == sum(args)
    end
end

@testset "keyword arguments" begin
    @eval inner_kwargf(foobar;foo=1, bar=2) = nothing

    @cuda (()->inner_kwargf(42;foo=1,bar=2))()

    @cuda (()->inner_kwargf(42))()

    @cuda (()->inner_kwargf(42;foo=1))()

    @cuda (()->inner_kwargf(42;bar=2))()

    @cuda (()->inner_kwargf(42;bar=2,foo=1))()
end

@testset "captured values" begin
    function f(capture::T) where {T}
        function kernel(ptr)
            unsafe_store!(ptr, capture)
            return
        end

        arr = CuArray(zeros(T))
        @cuda kernel(pointer(arr))

        return Array(arr)[1]
    end

    using Test
    @test f(1) == 1
    @test f(2) == 2
end

end

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

@testset "exceptions" begin

@testset "stack traces at different debug levels" begin

script = """
    function kernel(ptr, val)
        unsafe_store!(ptr, Int(val))
        return
    end

    cpu = zeros(Int)
    gpu = CUDAdrv.Mem.alloc(CUDAdrv.Mem.Device, sizeof(cpu))
    gpu_ptr = convert(CUDAdrv.CuPtr{Int}, gpu)
    @cuda kernel(gpu_ptr, 1.2)
    unsafe_copyto!(pointer(cpu), gpu_ptr, 1)
"""

let (code, out, err) = julia_script(script, `-g0`)
    @test code == 1
    @test occursin("ERROR: KernelException: exception thrown during kernel execution on device", err)
    @test isempty(out)
end

let (code, out, err) = julia_script(script, `-g1`)
    @test code == 1
    @test occursin("ERROR: KernelException: exception thrown during kernel execution on device", err)
    @test occursin("ERROR: a exception was thrown during kernel execution", out)
    @test occursin("Run Julia on debug level 2 for device stack traces", out)
end

let (code, out, err) = julia_script(script, `-g2`)
    @test code == 1
    @test occursin("ERROR: KernelException: exception thrown during kernel execution on device", err)
    @test occursin("ERROR: a exception was thrown during kernel execution", out)
    if VERSION < v"1.3.0-DEV.270"
        @test occursin("[1] Type at float.jl", out)
    else
        @test occursin("[1] Int64 at float.jl", out)
    end
    @test occursin("[2] kernel at none:2", out)
end

end

@testset "#329" begin

script = """
    @noinline foo(a, i) = a[1] = i
    bar(a) = (foo(a, 42); nothing)

    ptr = CUDAnative.DevicePtr{Int,AS.Global}(0)
    arr = CuDeviceArray{Int,1,AS.Global}((0,), ptr)

    @cuda bar(arr)
    CUDAdrv.synchronize()
"""

let (code, out, err) = julia_script(script, `-g2`)
    @test code == 1
    @test occursin("ERROR: KernelException: exception thrown during kernel execution on device", err)
    @test occursin("ERROR: a exception was thrown during kernel execution", out)
    @test occursin("foo at none:1", out)
    @test occursin("bar at none:2", out)
end

end

end

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

@testset "shmem divergence bug" begin

@testset "trap" begin
    function trap()
        ccall("llvm.trap", llvmcall, Cvoid, ())
    end

    function kernel(input::Int32, output::CUDAnative.DevicePtr{Int32}, yes::Bool=true)
        i = threadIdx().x

        temp = @cuStaticSharedMem(Cint, 1)
        if i == 1
            yes || trap()
            temp[1] = input
        end
        sync_threads()

        yes || trap()
        unsafe_store!(output, temp[1], i)

        return nothing
    end

    input = rand(Cint(1):Cint(100))
    N = 2

    let output = CuArray(zeros(Cint, N))
        # defaulting to `true` embeds this info in the PTX module,
        # allowing `ptxas` to emit validly-structured code.
        ptr = pointer(output)
        @cuda threads=N kernel(input, ptr)
        @test Array(output) == fill(input, N)
    end

    let output = CuArray(zeros(Cint, N))
        ptr = pointer(output)
        @cuda threads=N kernel(input, ptr, true)
        @test Array(output) == fill(input, N)
    end
end

@testset "unreachable" begin
    function unreachable()
        @cuprintln("go home ptxas you're drunk")
        Base.llvmcall("unreachable", Cvoid, Tuple{})
    end

    function kernel(input::Int32, output::CUDAnative.DevicePtr{Int32}, yes::Bool=true)
        i = threadIdx().x

        temp = @cuStaticSharedMem(Cint, 1)
        if i == 1
            yes || unreachable()
            temp[1] = input
        end
        sync_threads()

        yes || unreachable()
        unsafe_store!(output, temp[1], i)

        return nothing
    end

    input = rand(Cint(1):Cint(100))
    N = 2

    let output = CuArray(zeros(Cint, N))
        # defaulting to `true` embeds this info in the PTX module,
        # allowing `ptxas` to emit validly-structured code.
        ptr = pointer(output)
        @cuda threads=N kernel(input, ptr)
        @test Array(output) == fill(input, N)
    end

    let output = CuArray(zeros(Cint, N))
        ptr = pointer(output)
        @cuda threads=N kernel(input, ptr, true)
        @test Array(output) == fill(input, N)
    end
end

@testset "mapreduce (full)" begin
    function mapreduce_gpu(f::Function, op::Function, A::CuArray{T, N}; dims = :, init...) where {T, N}
        OT = Float32
        v0 =  0.0f0

        threads = 256
        out = CuArray{OT}(undef, (1,))
        @cuda threads=threads reduce_kernel(f, op, v0, A, Val{threads}(), out)
        Array(out)[1]
    end

    function reduce_kernel(f, op, v0::T, A, ::Val{LMEM}, result) where {T, LMEM}
        tmp_local = @cuStaticSharedMem(T, LMEM)
        global_index = threadIdx().x
        acc = v0

        # Loop sequentially over chunks of input vector
        while global_index <= length(A)
            element = f(A[global_index])
            acc = op(acc, element)
            global_index += blockDim().x
        end

        # Perform parallel reduction
        local_index = threadIdx().x - 1
        @inbounds tmp_local[local_index + 1] = acc
        sync_threads()

        offset = blockDim().x ÷ 2
        while offset > 0
            @inbounds if local_index < offset
                other = tmp_local[local_index + offset + 1]
                mine = tmp_local[local_index + 1]
                tmp_local[local_index + 1] = op(mine, other)
            end
            sync_threads()
            offset = offset ÷ 2
        end

        if local_index == 0
            result[blockIdx().x] = @inbounds tmp_local[1]
        end

        return
    end

    A = rand(Float32, 1000)
    dA = CuArray(A)

    @test mapreduce(identity, +, A) ≈ mapreduce_gpu(identity, +, dA)
end

@testset "mapreduce (full, complex)" begin
    function mapreduce_gpu(f::Function, op::Function, A::CuArray{T, N}; dims = :, init...) where {T, N}
        OT = Complex{Float32}
        v0 =  0.0f0+0im

        threads = 256
        out = CuArray{OT}(undef, (1,))
        @cuda threads=threads reduce_kernel(f, op, v0, A, Val{threads}(), out)
        Array(out)[1]
    end

    function reduce_kernel(f, op, v0::T, A, ::Val{LMEM}, result) where {T, LMEM}
        tmp_local = @cuStaticSharedMem(T, LMEM)
        global_index = threadIdx().x
        acc = v0

        # Loop sequentially over chunks of input vector
        while global_index <= length(A)
            element = f(A[global_index])
            acc = op(acc, element)
            global_index += blockDim().x
        end

        # Perform parallel reduction
        local_index = threadIdx().x - 1
        @inbounds tmp_local[local_index + 1] = acc
        sync_threads()

        offset = blockDim().x ÷ 2
        while offset > 0
            @inbounds if local_index < offset
                other = tmp_local[local_index + offset + 1]
                mine = tmp_local[local_index + 1]
                tmp_local[local_index + 1] = op(mine, other)
            end
            sync_threads()
            offset = offset ÷ 2
        end

        if local_index == 0
            result[blockIdx().x] = @inbounds tmp_local[1]
        end

        return
    end

    A = rand(Complex{Float32}, 1000)
    dA = CuArray(A)

    @test mapreduce(identity, +, A) ≈ mapreduce_gpu(identity, +, dA)
end

@testset "mapreduce (reduced)" begin
    function mapreduce_gpu(f::Function, op::Function, A::CuArray{T, N}) where {T, N}
        OT = Int
        v0 = 0

        out = CuArray{OT}(undef, (1,))
        @cuda threads=64 reduce_kernel(f, op, v0, A, out)
        Array(out)[1]
    end

    function reduce_kernel(f, op, v0::T, A, result) where {T}
        tmp_local = @cuStaticSharedMem(T, 64)
        acc = v0

        # Loop sequentially over chunks of input vector
        i = threadIdx().x
        while i <= length(A)
            element = f(A[i])
            acc = op(acc, element)
            i += blockDim().x
        end

        # Perform parallel reduction
        @inbounds tmp_local[threadIdx().x] = acc
        sync_threads()

        offset = blockDim().x ÷ 2
        while offset > 0
            @inbounds if threadIdx().x <= offset
                other = tmp_local[(threadIdx().x - 1) + offset + 1]
                mine = tmp_local[threadIdx().x]
                tmp_local[threadIdx().x] = op(mine, other)
            end
            sync_threads()
            offset = offset ÷ 2
        end

        if threadIdx().x == 1
            result[blockIdx().x] = @inbounds tmp_local[1]
        end

        return
    end

    A = rand(1:10, 100)
    dA = CuArray(A)

    @test mapreduce(identity, +, A) ≈ mapreduce_gpu(identity, +, dA)
end

end

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

@testset "dynamic parallelism" begin

@testset "basic usage" begin
    function hello()
        @cuprint("Hello, ")
        @cuda dynamic=true world()
        return
    end

    @eval function world()
        @cuprint("World!")
        return
    end

    _, out = @grab_output begin
        @cuda hello()
        synchronize()
    end
    @test out == "Hello, World!"
end

@testset "anonymous functions" begin
    function hello()
        @cuprint("Hello, ")
        world = () -> (@cuprint("World!"); nothing)
        @cuda dynamic=true world()
        return
    end

    _, out = @grab_output begin
        @cuda hello()
        synchronize()
    end
    @test out == "Hello, World!"
end

if VERSION >= v"1.1" # behavior of captured variables (box or not) has improved over time
@testset "closures" begin
    function hello()
        x = 1
        @cuprint("Hello, ")
        world = () -> (@cuprint("World $(x)!"); nothing)
        @cuda dynamic=true world()
        return
    end

    _, out = @grab_output begin
        @cuda hello()
        synchronize()
    end
    @test out == "Hello, World 1!"
end
end

@testset "argument passing" begin
    ## padding

    function kernel(a, b, c)
        @cuprint("$a $b $c")
        return
    end

    for args in ((Int16(1), Int32(2), Int64(3)),    # padding
                 (Int32(1), Int32(2), Int32(3)),    # no padding, equal size
                 (Int64(1), Int32(2), Int16(3)),    # no padding, inequal size
                 (Int16(1), Int64(2), Int32(3)))    # mixed
        _, out = @grab_output begin
            @cuda kernel(args...)
            synchronize()
        end
        @test out == "1 2 3"
    end

    ## conversion

    function kernel(a)
        increment(a) = (a[1] += 1; nothing)

        a[1] = 1
        increment(a)
        @cuda dynamic=true increment(a)

        return
    end

    dA = CuArray{Int}(undef, (1,))
    @cuda kernel(dA)
    A = Array(dA)
    @test A == [3]
end

@testset "self-recursion" begin
    @eval function kernel(x::Bool)
        if x
            @cuprint("recurse ")
            @cuda dynamic=true kernel(false)
        else
            @cuprint("stop")
        end
       return
    end

    _, out = @grab_output begin
        @cuda kernel(true)
        synchronize()
    end
    @test out == "recurse stop"
end

@testset "deep recursion" begin
    @eval function kernel_a(x::Bool)
        @cuprint("a ")
        @cuda dynamic=true kernel_b(x)
        return
    end

    @eval function kernel_b(x::Bool)
        @cuprint("b ")
        @cuda dynamic=true kernel_c(x)
        return
    end

    @eval function kernel_c(x::Bool)
        @cuprint("c ")
        if x
            @cuprint("recurse ")
            @cuda dynamic=true kernel_a(false)
        else
            @cuprint("stop")
        end
        return
    end

    _, out = @grab_output begin
        @cuda kernel_a(true)
        synchronize()
    end
    @test out == "a b c recurse a b c stop"
end

@testset "streams" begin
    function hello()
        @cuprint("Hello, ")
        s = CuDeviceStream()
        @cuda dynamic=true stream=s world()
        CUDAnative.unsafe_destroy!(s)
        return
    end

    @eval function world()
        @cuprint("World!")
        return
    end

    _, out = @grab_output begin
        @cuda hello()
        synchronize()
    end
    @test out == "Hello, World!"
end

end

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

if capability(dev) >= v"6.0" && attribute(dev, CUDAdrv.DEVICE_ATTRIBUTE_COOPERATIVE_LAUNCH) == 1

@testset "cooperative groups" begin
    function kernel_vadd(a, b, c)
        i = (blockIdx().x-1) * blockDim().x + threadIdx().x
        grid_handle = this_grid()
        c[i] = a[i] + b[i]
        sync_grid(grid_handle)
        c[i] = c[1]
        return nothing
    end

    a = round.(rand(Float32, (300, 40)) * 100)
    b = round.(rand(Float32, (300, 40)) * 100)
    c = zeros(Float32, (300, 40))
    d_a = CuArray(a)
    d_b = CuArray(b)
    d_c = CuArray(c)  # output array
    @cuda cooperative=true threads=600 blocks=20 kernel_vadd(d_a, d_b, d_c)
    c = Array(d_c)
    @test all(c[1] .== c)
end

end

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

end


text:
function f_gd(x)
  (x[1] - 5.0)^2
end

function g_gd(x, storage)
  storage[1] = 2.0 * (x[1] - 5.0)
end

initial_x = [0.0]

d = DifferentiableFunction(f_gd, g_gd)

results = Optim.gradient_descent(d, initial_x)
@assert isempty(results.trace.states)
@assert results.converged
@assert norm(results.minimum - [5.0]) < 0.01

eta = 0.9

function f_gd(x)
  (1.0 / 2.0) * (x[1]^2 + eta * x[2]^2)
end

function g_gd(x, storage)
  storage[1] = x[1]
  storage[2] = eta * x[2]
end

d = DifferentiableFunction(f_gd, g_gd)

results = Optim.gradient_descent(d, [1.0, 1.0])
@assert isempty(results.trace.states)
@assert results.converged
@assert norm(results.minimum - [0.0, 0.0]) < 0.01


text:
using CUDA
using Adapt
using DataStructures

@enum Dif_equation LiHoChow2005 NaturalNet

struct ContinuousTimeRNN{A,B}
    delta_t::Float32
    number_neurons::Int64
    differential_equation::Dif_equation
    clipping_range_min::Float32
    clipping_range_max::Float32
    set_principle_diagonal_elements_of_W_negative::Bool
    alpha::Float32
    V::A
    W::A
    T::A
    x::B
    input_size::Int64
    output_size::Int64
end

function ContinuousTimeRNN(configuration::OrderedDict, number_inputs::Int, number_outputs::Int, number_individuals::Int)

    if configuration["differential_equation"] == "NaturalNet"
        differential_eq = NaturalNet
    elseif configuration["differential_equation"] == "LiHoChow2005"
        differential_eq = LiHoChow2005
    end

    ContinuousTimeRNN(
        convert(Float32, configuration["delta_t"]),
        configuration["number_neurons"],
        differential_eq,
        convert(Float32, configuration["clipping_range_min"]),
        convert(Float32, configuration["clipping_range_max"]),
        convert(Bool, configuration["set_principle_diagonal_elements_of_W_negative"]),
        convert(Float32, configuration["alpha"]),
        CUDA.fill(0.0f0, (configuration["number_neurons"], number_inputs, number_individuals)),
        CUDA.fill(0.0f0, (configuration["number_neurons"], configuration["number_neurons"], number_individuals)),
        CUDA.fill(0.0f0, (number_outputs, configuration["number_neurons"], number_individuals)),
        CUDA.fill(0.0f0, (configuration["number_neurons"], number_individuals)),
        number_inputs,
        number_outputs,
    )
end

Adapt.@adapt_structure ContinuousTimeRNN

function get_memory_requirements(brains::ContinuousTimeRNN)
    return sizeof(Float32) * brains.number_neurons
end

function initialize(brains::ContinuousTimeRNN, individuals)

    threadID = threadIdx().x
    blockID = blockIdx().x

    v_size = brains.input_size * brains.number_neurons
    w_size = brains.number_neurons * brains.number_neurons

    for i = 1:brains.input_size
        @inbounds brains.V[threadID, i, blockID] = individuals[blockID, threadID+(i-1)*brains.number_neurons]
    end

    for i = 1:brains.number_neurons
        @inbounds brains.W[threadID, i, blockID] = individuals[blockID, v_size+(threadID+(i-1)*brains.number_neurons)]
    end

    sync_threads()

    if brains.set_principle_diagonal_elements_of_W_negative == true
        @inbounds brains.W[threadID, threadID, blockID] = -abs(brains.W[threadID, threadID, blockID])
    end

    for i = 1:brains.output_size
        @inbounds brains.T[i, threadID, blockID] = individuals[blockID, v_size+w_size+(i+(threadID-1)*brains.output_size)]
    end

    sync_threads()

    reset(brains)

    sync_threads()

end

function reset(brains::ContinuousTimeRNN)

    threadID = threadIdx().x
    blockID = blockIdx().x

    if threadID <= brains.number_neurons
        @inbounds brains.x[threadID, blockID] = 0.0
    end

    sync_threads()

end

function step(brains::ContinuousTimeRNN, input, output, offset_shared_memory)

    threadID = threadIdx().x
    blockID = blockIdx().x

    V_value = @cuDynamicSharedMem(Float32, brains.number_neurons, offset_shared_memory)

    if threadID <= brains.number_neurons

        # V_value = V * input (Matrix-Vector-Multiplication)
        V_value[threadID] = 0.0
        for i = 1:brains.input_size
            @inbounds V_value[threadID] += (brains.V[threadID, i, blockID] * input[i])
        end

        sync_threads()

        if brains.differential_equation == NaturalNet

            # W_value = W * tanh(x) (Matrix-Vector-Multiplication)
            W_value = 0.0
            for i = 1:brains.number_neurons
                @inbounds W_value += (brains.W[threadID, i, blockID] * tanh(brains.x[i, blockID]))
            end

            sync_threads()

            # Differential Equation
            dx_dt = -brains.alpha * brains.x[threadID, blockID] + W_value + V_value[threadID]

        elseif brains.differential_equation == LiHoChow2005

            # W_value = W * (tanh(x) + V_value) (Matrix-Vector-Multiplication)
            W_value = 0.0
            for i = 1:brains.number_neurons
                @inbounds W_value += (brains.W[threadID, i, blockID] * (tanh(brains.x[i, blockID] + V_value[i])))
            end

            sync_threads()

            # Differential Equation
            dx_dt = -brains.alpha * brains.x[threadID, blockID] + W_value

        end

        # Euler forward discretization
        @inbounds brains.x[threadID, blockID] += brains.delta_t * dx_dt

        # Clip x to state boundaries
        @inbounds brains.x[threadID, blockID] = clamp(brains.x[threadID, blockID], brains.clipping_range_min, brains.clipping_range_max)

    end

    sync_threads()

    if threadID <= brains.output_size

        # T_value = T * x (Matrix-Vector-Multiplication)
        T_value = 0.0
        for i = 1:brains.number_neurons
            T_value += brains.T[threadID, i, blockID] * brains.x[i, blockID]
        end

        # Calculate outputs
        output[threadID] = tanh(T_value)

    end

    sync_threads()

end

function get_free_parameter_usage(brains)

    usage_dict = Dict()

    usage_dict["V"] = brains.input_size * brains.number_neurons
    usage_dict["W"] = brains.number_neurons * brains.number_neurons
    usage_dict["T"] = brains.number_neurons * brains.output_size

    return usage_dict
end

function get_individual_size(brains)

    usage_dict = get_free_parameter_usage(brains)

    return usage_dict["V"] + usage_dict["W"] + usage_dict["T"]
end


text:
# ```@setup fitnoisemodel
# # the example command line input for this script
# ARGS = ["example.toml", "100", "EE", "--plot"]
# ```
#
# # [Fit Noise Model (fitnoisemodel.jl)](@id fitnoisemodel)
# This script estimates a smooth fit to the noise power spectrum. In this pipeline, we estimate the noise power spectrum
# of each frequency map using the difference between half-mission auto- and cross-spectra,
# ```math
# \begin{split}
# N_{\ell}^{X, \mathrm{hm1} \times \mathrm{hm1}} &=  C_{\ell}^{X, \mathrm{hm1} \times \mathrm{hm1}} -  C_{\ell}^{X, \mathrm{hm1} \times \mathrm{hm2}}, \\
# N_{\ell}^{X, \mathrm{hm2} \times \mathrm{hm2}} &=  C_{\ell}^{X, \mathrm{hm2} \times \mathrm{hm2}} -  C_{\ell}^{X, \mathrm{hm1} \times \mathrm{hm2}}.
# \end{split}
# ```
#
# ```@raw html
# <pre class="shell">
# <code class="language-shell hljs">$ julia fitnoisemodel.jl example.toml 100 EE</code></pre>
# ```

configfile, freq, spec = ARGS
if length(ARGS) > 3
    if ARGS[4] != "--plot"
        maxfuncevals = ARGS[4]
    else
        maxfuncevals = 500_000
    end
end

## setup data
using Plots
using TOML
using BlackBoxOptim
include("util.jl")

config = TOML.parsefile(configfile)
nside = config["general"]["nside"]
run_name = config["general"]["name"]
spectrapath = joinpath(config["scratch"], "rawspectra")
XY = Symbol(spec)
lmax = min(2508, nside2lmax(nside))

# Next, we check to see if we need to render plots for the Documentation.

if "--plot" ∉ ARGS
    Plots.plot(args...; kwargs...) = nothing
    Plots.plot!(args...; kwargs...) = nothing
end

##
# We read in the raw spectra generated from rawspectra.jl.

Cl11 = CSV.read(joinpath(spectrapath,"$(run_name)_P$(freq)hm1xP$(freq)hm1.csv"), DataFrame);
Cl12 = CSV.read(joinpath(spectrapath,"$(run_name)_P$(freq)hm1xP$(freq)hm2.csv"), DataFrame);
Cl22 = CSV.read(joinpath(spectrapath,"$(run_name)_P$(freq)hm2xP$(freq)hm2.csv"), DataFrame);

truncate(vec::AbstractVector, lmax) = SpectralVector(vec[firstindex(vec):(lmax+1)])
truncate(vec::SpectralVector, lmax) = vec[IdentityRange(firstindex(vec):lmax)]

cl11 = truncate(Cl11[!,XY], lmax)
cl12 = truncate(Cl12[!,XY], lmax)
cl22 = truncate(Cl22[!,XY], lmax)

beampath = joinpath(config["scratch"], "beams")
Wl11 = util_planck_beam_Wl(freq, "hm1", freq, "hm1", XY, XY;
    lmax=lmax, beamdir=beampath)
Wl12 = util_planck_beam_Wl(freq, "hm1", freq, "hm2", XY, XY;
    lmax=lmax, beamdir=beampath)
Wl22 = util_planck_beam_Wl(freq, "hm2", freq, "hm2", XY, XY;
    lmax=lmax, beamdir=beampath)

Wl11, Wl12, Wl22 = map(v->truncate(v,lmax), (Wl11, Wl12, Wl22))

cl = cl12 ./ Wl12
nl1 = (cl11 ./ Wl11 .- cl) .* Wl11
nl2 = (cl22 ./ Wl22  .- cl) .* Wl22;

ell = eachindex(cl)
plot(ell, ell.^2 .* cl12, label="unbeamed $(run_name) $(freq) $(spec)",
    xlabel="multipole moment", ylabel="\$\\ell^2 C_{\\ell}^{$spec}\$", xlim=(0,lmax))

#

@. camspec_model(ℓ, α) =
    α[1] * (100. / ℓ)^α[2] + α[3] * (ℓ / 1000.)^α[4] / ( 1 + (ℓ / α[5])^α[6] )^α[7]

#

function fit_bb_model(model, p0, xl, yl, signal; kwargs...)
    lower = map(x -> x-0.9abs(x), p0)
    upper = map(x -> x+0.9abs(x), p0)
    like(α) = (sum((2 .* xl .+ 1) ./ (model(xl, p0).^2 .+ signal.^2) .* (model(xl, α) .- yl).^2))
    println("starting opt ", like(p0))
    res = bboptimize(like; SearchRange=map((i,j)->(i,j), lower, upper), NumDimensions = length(p0),
        MaxFuncEvals=maxfuncevals, TraceInterval=20)
    return best_candidate(res)
end

#
p0_1 = readdlm(joinpath(@__DIR__, "../", "input",
        "initial_planck_noise_coeffs", "$(freq)_hm1_$(spec)_coeff.dat"))[:,1]
p0_2 = readdlm(joinpath(@__DIR__, "../", "input",
        "initial_planck_noise_coeffs", "$(freq)_hm2_$(spec)_coeff.dat"))[:,1]

#
min_ell_ind = 31 # lmin=30
pfit_1 = fit_bb_model(camspec_model, p0_1,
    parent(ell)[min_ell_ind:end], parent(nl1)[min_ell_ind:end], parent(cl)[min_ell_ind:end])
pfit_2 = fit_bb_model(camspec_model, p0_2,
    parent(ell)[min_ell_ind:end], parent(nl2)[min_ell_ind:end], parent(cl)[min_ell_ind:end])

#

mean(x) = sum(x) / length(x)
plot((nl1), alpha=0.5, label="nl hm1")
plot!(2:lmax, [camspec_model(ℓ, p0_1) for ℓ in 2:lmax], ylim=(0.0, 4mean(abs.(nl1))), label="initial model")
plot!(2:lmax, [camspec_model(ℓ, pfit_1) for ℓ in 2:lmax], label="fitted model", linestyle=:dash)

#
##
plot((nl2), alpha=0.5, label="nl hm2")
plot!(2:lmax, [camspec_model(ℓ, p0_2) for ℓ in 2:lmax], ylim=(0.0, 2mean(nl2)), label="initial model")
plot!(2:lmax, [camspec_model(ℓ, pfit_2) for ℓ in 2:lmax], ylim=(0.0, 2mean(nl2)), label="fitted model", linestyle=:dash)

#
coefficientpath = joinpath(config["scratch"], "noise_model_coeffs")
mkpath(coefficientpath)

open(joinpath(coefficientpath, "$(run_name)_$(freq)_$(spec)_hm1.dat"), "w") do io
   writedlm(io, pfit_1)
end
open(joinpath(coefficientpath, "$(run_name)_$(freq)_$(spec)_hm2.dat"), "w") do io
   writedlm(io, pfit_2)
end

max_forks_repo_name: sgiardie/PSpipe


text:
using Plots, Printf, DelimitedFiles, Random, Statistics, ProgressMeter
using MDToolbox

using BSON: @save, @load
@load "afm_generate_test_data.bson" afms models models_template qs myseed sigma_noise nframe

configs = [AfmizeConfig(10.0 * (pi / 180),
    r,
    MDToolbox.Point2D(-250, -200),
    MDToolbox.Point2D(250, 200),
    MDToolbox.Point2D(6.25, 6.25),
    MDToolbox.defaultParameters())
    for r in [10.0, 15.0, 20.0, 25.0, 30.0]]

# single frame
r = getafmposterior(afms[1], models_template, qs, configs)

# all frames with parallel map
r = pmap(x -> getafmposterior(x, models_template, qs, configs), afms[1:20])

@save "afm_analyze_test_data.bson" r configs

max_stars_repo_name: matsunagalab/notebooks

max_issues_repo_name: matsunagalab/notebooks

max_forks_repo_name: matsunagalab/notebooks


text:
using Documenter, BifurcationKit, Setfield
# using DocThemeIndigo
ENV["GKSwstype"] = "100"

makedocs(doctest = false,
	sitename = "Bifurcation Analysis in Julia",
	format = Documenter.HTML(collapselevel = 1,assets = ["assets/indigo.css"]),
	# format = DocumenterLaTeX.LaTeX(),
	authors = "Romain Veltz",
	pages = Any[
		"Home" => "index.md",
		"Overview" => "guidelines.md",
		"Tutorials" => "tutorials/tutorials.md",
		"Functionalities" => [
			"Plotting" => "plotting.md",
			"Periodic Orbits" => [
				"Introduction" => "periodicOrbit.md",
				"Finite Differences" => "periodicOrbitTrapeze.md",
				"Collocation" => "periodicOrbitCollocation.md",
				"Shooting" => "periodicOrbitShooting.md",
				],
			"Symmetries / Waves" => [
				"Introduction" => "intro_wave.md",
				],
			"Continuation methods" => [
					"Introduction" => "IntroContinuation.md",
			"Predictors / correctors" => "Predictors.md",
					"PALC" => "PALC.md",
					"Moore-Penrose Continuation" => "MooreSpence.md",
					"Deflated Continuation" => "DeflatedContinuation.md",
				],
			"Event Handling and Callback" => "EventCallback.md",
			"Bifurcations" => [
				"Bifurcation detection (codim 1)" => "detectionBifurcation.md",
				"Fold / Hopf Continuation (codim 2)" => "codim2Continuation.md",
				],
			"Normal form" =>[
				"Simple branch point" => "simplebp.md",
				"Non-simple branch point" => "nonsimplebp.md",
				"Simple Hopf point" => "simplehopf.md",
			],
			"Branch switching" => "branchswitching.md",
			"Bifurcation diagram" => "BifurcationDiagram.md",
			"Deflated problems" => "deflatedproblem.md",
			"Constrained problem" => "constrainedproblem.md",
			"DiffEq wrapper" => "diffeq.md",
			"Iterator Interface" => "iterator.md",
		],
		"Options" => [
			"Linear Solvers" => "linearsolver.md",
			"Bordered linear solvers" => "borderedlinearsolver.md",
			"Eigen Solvers" => "eigensolver.md",
			"Bordered arrays" => "Borderedarrays.md",
		],
		"Frequently Asked Questions" => "faq.md",
		"Library" => "library.md"
	]
	)

deploydocs(
	repo = "github.com/bifurcationkit/BifurcationKitDocs.jl.git",
	devbranch = "main"
)


text:
using Merly
using JSON
using MySQL
using Tables: rowtable

const host = "127.0.0.1"
const user = "root"
const password = ""

const conn = DBInterface.connect(MySQL.Connection, host, user, password)
DBInterface.execute(conn, "create database if not exists julia_react_todo")

create_table = """
create table if not exists julia_react_todo.todos(
    id int PRIMARY KEY,
    description varchar(255),
    deadline date
)
"""
DBInterface.execute(conn, create_table)

function find_all(conn)
    res = DBInterface.execute(conn, "select * from julia_react_todo.todos")
    return json(rowtable(res))
end

function delete(conn, id)
    DBInterface.execute(conn, "delete from julia_react_todo.todos where id = $id")
end

function insert(conn, todo)
    id = todo["id"]
    description = todo["description"]
    deadline = todo["deadline"]
    DBInterface.execute(conn, "insert into julia_react_todo.todos values($id, \"$description\", \"$deadline\")")
end

useCORS(
AllowOrigins = "*"
, AllowHeaders = "Origin, Content-Type, Accept"
, AllowMethods = "GET,POST,PUT,DELETE"
, MaxAge = "178000")

Get("/todo", (request, HTTP) -> begin
    println(find_all(conn))

    HTTP.Response(201
          , HTTP.mkheaders(["Content-Type" => "application/json"])
          ,body=find_all(conn))
end)

Post("/todo", (request, HTTP) -> begin
    println(request.body)
    todo = JSON.parse(request.body)
    insert(conn, todo)
    HTTP.Response(201, "ok")
end)

Delete("/todo/:id", (request, HTTP) -> begin
    delete(conn, request.params["id"])
    HTTP.Response(201, "ok")
end)

start(host = "127.0.0.1", port = 8086, verbose = true)

max_stars_repo_name: kogad/JuliaReactToDo

max_issues_repo_name: kogad/JuliaReactToDo

max_forks_repo_name: kogad/JuliaReactToDo


text:
module Glimmer

using Blink

using LinearAlgebra
import MeshCat
using CoordinateTransformations
using Rotations
import GeometryBasics
using Colors: RGBA, RGB
import Colors
import ColorSchemes
import FileIO
import MeshIO

import Tables
import JSON, JSONTables
import Dates
using OrderedCollections

using StaticArrays
import UUIDs
import Base64

using Reexport

include("Misc.jl")
include("GridUtils.jl")
include("Material.jl")
include("SceneObject.jl")
include("Scene.jl")

include("JSInterface/JuliaJSBridge.jl")
include("UI/UIVariables.jl")
include("UI/UIControls.jl")
include("UI/WGL.jl")

include("UI/FlexUI.jl")

# reexport all the public stuff from UIControls and UIVariables

@reexport using ..UIControls
@reexport using ..UIVariables
@reexport using ..WGL
@reexport using ..GridUtils

#------------------------------------------------------------------------------
# Exports
#------------------------------------------------------------------------------
export  Scene,
        set_Y_up!,
        set_Z_up!,
        grid!,
        axes!,
        # root, root!,
        parent, parent!,
        material, material!,
        Material,
        transform, lookAt,
        rotation,
        local_tr, local_tr!,
        cameraTransform!, cameraPlanes!,
        EmptySceneObject,
        Box,
        Mesh,
        Arrow,
        Axes,
        LineSegments,
        PointCloud,
        dataPath,
        unitX3,
        unitY3,
        unitZ3,
        zero3,
        cameraTransform!,
        url,
        clear,
        send!,
        root,
        root_transform,
        win,
        vis

export  BasicValidation, on, Variable

export  Container,
        Slider,
        Button,
        MeshCatViewer,
        Label,
        Image,
        PanZoom,
        Field,
        ButtonToggle,
        RadioGroup,
        CheckBox,
        ExpansionPanel,
        Accordion,
        Tabs,
        Tab,
        Divider,
        Card,
        Markdown,

        VContainer,
        HContainer,
        HContainerSpace,
        H1Label,
        H2Label,
        H3Label,
        H4Label

# export the common functions of a UI implementation
export addVariable!, prop, prop!, controls, controls!, renderFunction, renderFunction!

end # module Glimmer


text:
using PeeEmTee
using Test
using StatsBase

const pmt_data = joinpath(@__DIR__, "data", "pmt_data.h5")

@testset "high_voltages()" begin
    @test high_voltages(pmt_data) == ["1000", "1100", "1200"]
end

@testset "WaveSet" begin
    a = [[1, 2, 3] [4, 5, 6] [7, 8, 9]]
    waveset = WaveSet(pmt_data, 1000)
    @test waveset.waveforms == a * waveset.v_gain
    @test waveset.h_int == 1.0
    waveset = WaveSet(pmt_data, 1100)
    @test waveset.waveforms == a * waveset.v_gain
    @test waveset.h_int == 1.0
    waveset = WaveSet(pmt_data, 1200)
    @test waveset.waveforms == a * waveset.v_gain
    @test waveset.h_int == 1.0
end

@testset "calculate_charges" begin
    a = [[0, 1, -45, -53, 0, -1] [-5, 3, -145, -253, 3, -5] [0, 0, -44, -12, -1, 1]]
    @test calculate_charges(a, 1, 2, 3, 4) == [99, 396, 56]
end

@testset "calculate_transit_times" begin
    a = [[0., -1., 1., -53., 0., 1.] [-5., 3., -14., -253., 143., 3.] [100., 102., 85., -9., -1., -15.]]
    @test calculate_transit_times(a, -10.) == [3, 3, 2]
end

@testset "bin_data()" begin
    a = [1, 2, 3, 4, 5, 6, 7, 8, 9, 10]
    @test bin_data(a, 1.:1.:10.) == ([1.5, 2.5, 3.5, 4.5, 5.5, 6.5, 7.5, 8.5, 9.5],
                                     [1, 1, 1, 1, 1, 1, 1, 1, 1])
end

@testset "ChargeDist" begin
    a = [1, 2, 3, 4, 5, 6, 7, 8, 9, 10]
    chargedist = ChargeDist(a, 1.:1.:10.)
    @test chargedist.x == [1.5, 2.5, 3.5, 4.5, 5.5, 6.5, 7.5, 8.5, 9.5]
    @test chargedist.y == [1, 1, 1, 1, 1, 1, 1, 1, 1]
    b = [-.8, -0.3, 0.1, 0.5]
    chargedist = ChargeDist(b, 4)
    @test chargedist.x == [-0.75, -0.25, 0.25, 0.75]
    @test chargedist.y == [1, 1, 1, 1]
end

@testset "norm_gauss()" begin
    @test norm_gauss(0, [0, 1]) == 0.3989422804014327
end

@testset "gauss()" begin
    @test gauss(0, [0, 1, 10]) == 3.9894228040143274
end

@testset "baseline_max" begin
    a = [[0., 0., -1., 0., 0.] [0., 0., -2., 0., -1.] [0., 1., -3., -2., -1.]]
    @test baseline_max(a) == 2
end

@testset "subtract_baseline()" begin
    a = [[-1., 0., -2., 0., 0.] [0., 0., -2., 0., -1.] [0., 1., -3., -2., -1.]]
    a_res1 = [[-0.5, 0.5, -1.5, 0.5, 0.5] [0., 0., -2., 0., -1.] [-0.5, 0.5, -3.5, -2.5, -1.5]]
    a_res2 = [[0, 1., -1., 1., 1.] [0., 0., -2., 0., -1.] [0., 1., -3., -2., -1.]]
    @test subtract_baseline(a) == a_res1
    @test subtract_baseline(a; baseline_range=(1, 1)) == a_res2
end

@testset "pmtresp()" begin
    xs = [0, 0.3]
    params = [0, 0.05, 0.3, 0.2, 0.1, 100.]
    @test pmtresp(xs, params) == [727.8832805768451,
                                  18.4164040244264]
    prffit = PMTRespFit(params...)
    @test pmtresp(xs, prffit) == [727.8832805768451,
                                  18.4164040244264]
end

@testset "pmtresp_uap()" begin
    xs = [0.1, 0.2]
    params = [0, 0.05, 0.3, 0.2, 0.1, 100., 0.1, 0.02, 10.]
    @test pmtresp_uap(xs, params) == [308.25941556084206,
                                      16.40805075212658]
    prffit = PMTRespUapFit(params...)
    @test pmtresp_uap(xs, prffit) == [308.25941556084206,
                                      16.40805075212658]
end

@testset "mean_waveform()" begin
    a = [[0., 0., -1., 0., 0., 0.] [0., 0., 0., -1., 0., 0.] [0., -1., 0., 0., 0., 0.]]
    @test mean_waveform(a) == [0., 0., -1., 0., 0., 0.]
end

@testset "simulate_charges()" begin
    secs = (7, 4.5, 4.5, 4.5, 4.5, 4.5, 4.5, 4.5, 4.5, 4.5)
    nₚₑ = 0.1
    @test isapprox(mean(simulate_charges(nₚₑ, 100000; secs=secs)),
                   prod(secs) * nₚₑ;
                   rtol=0.05)
end


text:
# Put in a separate page so it can be used by SciMLDocs.jl

pages = [
        "Home" => "index.md",
        "Introduction" => "introduction.md",
        "APIs" => "apis.md",
        "References" => "references.md",
    ]


text:
using Images, ImageView, Serialization

# This is pure n-dimensional array form
raw_data = deserialize("imaging_data/neurondata.blob")

# This is normalized gray value form used by Image.jl
# Not mandatory we use it, but plays nicer with tooling in related packages (e.g. ImageView)
img = reinterpret(Gray{N4f12}, raw_data)

imshow(img)

max_stars_repo_name: wsphillips/jlcaimaging

max_issues_repo_name: wsphillips/jlcaimaging

max_forks_repo_name: wsphillips/jlcaimaging


text:
function get_leap(obj::GLONASSScale, arg0::AbsoluteDate)
    return jcall(obj, "getLeap", jdouble, (AbsoluteDate,), arg0)
end

function get_leap(obj::GLONASSScale, arg0::FieldAbsoluteDate)
    return jcall(obj, "getLeap", RealFieldElement, (FieldAbsoluteDate,), arg0)
end

function get_name(obj::GLONASSScale)
    return jcall(obj, "getName", JString, ())
end

function inside_leap(obj::GLONASSScale, arg0::AbsoluteDate)
    return jcall(obj, "insideLeap", jboolean, (AbsoluteDate,), arg0)
end

function inside_leap(obj::GLONASSScale, arg0::FieldAbsoluteDate)
    return jcall(obj, "insideLeap", jboolean, (FieldAbsoluteDate,), arg0)
end

function minute_duration(obj::GLONASSScale, arg0::AbsoluteDate)
    return jcall(obj, "minuteDuration", jint, (AbsoluteDate,), arg0)
end

function minute_duration(obj::GLONASSScale, arg0::FieldAbsoluteDate)
    return jcall(obj, "minuteDuration", jint, (FieldAbsoluteDate,), arg0)
end

function offset_from_tai(obj::GLONASSScale, arg0::AbsoluteDate)
    return jcall(obj, "offsetFromTAI", jdouble, (AbsoluteDate,), arg0)
end

function offset_from_tai(obj::GLONASSScale, arg0::FieldAbsoluteDate)
    return jcall(obj, "offsetFromTAI", RealFieldElement, (FieldAbsoluteDate,), arg0)
end

function offset_to_tai(obj::GLONASSScale, arg0::DateComponents, arg1::TimeComponents)
    return jcall(obj, "offsetToTAI", jdouble, (DateComponents, TimeComponents), arg0, arg1)
end

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


text:
module DBMS # Octo

abstract type AbstractDatabase end

struct SQL <: AbstractDatabase end
struct SQLite <: AbstractDatabase end
struct MySQL <: AbstractDatabase end
struct PostgreSQL <: AbstractDatabase end
struct Hive <: AbstractDatabase end

end # module Octo.DBMS


text:
__precompile__()

module QuantumOptics

export bases, Basis, GenericBasis, CompositeBasis,
        tensor, ⊗,
        states, StateVector, Bra, Ket, basis_bra, basis_ket,
                dagger, normalize, normalize!,
        operators, Operator, DenseOperator, projector,
                expect, dense_identity, ptrace, embed,
        operators_lazy, LazyOperator, LazyTensor, LazySum, LazyProduct,
        operators_sparse, SparseOperator,
        super, DenseSuperOperator, SparseSuperOperator,
                spre, spost, liouvillian,
        fock, FockBasis, number, destroy, create,
                fockstate, coherentstate, qfunc,
        spin, SpinBasis, sigmax, sigmay, sigmaz, sigmap, sigmam, spinup, spindown,
        subspace, SubspaceBasis, projector,
        particle, PositionBasis, MomentumBasis, gaussianstate,
                positionoperator, momentumoperator, laplace_x, laplace_p,
        nparticlebasis, BosonicNParticleBasis, FermionicNParticleBasis, nparticleoperator_1, nparticleoperator_2,
        metrics, tracedistance,
        spectralanalysis, operatorspectrum, operatorspectrum_hermitian,
                eigenstates, eigenstates_hermitian, groundstate,
        timeevolution_simple,
        timeevolution,
        cumulantexpansion,
        steadystate,
        correlations

include("bases.jl")
include("states.jl")
include("operators.jl")
include("sparsematrix.jl")
include("operators_sparse.jl")
include("operators_lazy.jl")
include("superoperators.jl")
include("spin.jl")
include("fock.jl")
include("subspace.jl")
include("particle.jl")
include("nparticles.jl")
include("metrics.jl")
include("ode_dopri.jl")
include("timeevolution_simple.jl")
module timeevolution
    include("master.jl")
    include("schroedinger.jl")
    include("mcwf.jl")
    using .timeevolution_master
    using .timeevolution_schroedinger
    using .timeevolution_mcwf
end
include("cumulantexpansion.jl")
include("steadystate.jl")
include("correlations.jl")
include("spectralanalysis.jl")
include("io.jl")

using .bases
using .states
using .operators
using .operators_sparse
using .operators_lazy
using .superoperators
using .spin
using .fock
using .subspace
using .particle
using .nparticles
using .metrics
using .spectralanalysis

end # module


text:
using JSON
using Test

import StaticArrays
const SA = StaticArrays

using EnsembleKalmanProcesses.ParameterDistributions

using CalibrateEDMF.ModelTypes
using CalibrateEDMF.ReferenceModels
import CalibrateEDMF.ReferenceModels: NameList
using CalibrateEDMF.ReferenceStats
using CalibrateEDMF.DistributionUtils
using CalibrateEDMF.TurbulenceConvectionUtils
import CalibrateEDMF.TurbulenceConvectionUtils: create_parameter_vectors

@testset "TurbulenceConvectionUtils" begin
    @testset "test create_parameter_vectors" begin
        u_names = ["foo", "bar"]
        u = [1.0, 2.0]
        u_names_out, u_out = create_parameter_vectors(u_names, u)
        @test all(u_names_out .∈ [["foo", "bar"]])
        @test only(u_out[u_names_out .== "foo"]) == 1.0
        @test only(u_out[u_names_out .== "bar"]) == 2.0

        u_names = ["bar_{2}", "foo", "bar_{1}"]
        u = [1.0, 2.0, 3.0]
        u_names_out, u_out = create_parameter_vectors(u_names, u)
        @test all(u_names_out .∈ [["foo", "bar"]])
        @test only(u_out[u_names_out .== "foo"]) == 2.0
        @test only(u_out[u_names_out .== "bar"]) == [3.0, 1.0]
    end

    @test get_gcm_les_uuid(1, forcing_model = "model1", month = 1, experiment = "experiment1") ==
          "1_model1_01_experiment1"

    u_names = ["foo", "bar_{1}", "bar_{2}"]
    u = [1.0, 2.0, 3.0]
    u_names_out, u_out = create_parameter_vectors(u_names, u)
    @test length(u_names_out) == 2
    @test isa(u_out, Vector{Any})

    @testset "TC.jl error handling" begin
        # Choose same SCM to speed computation
        data_dir = mktempdir()
        scm_dirs = [joinpath(data_dir, "Output.Rico.000000")]
        # Violate CFL condition for TC.jl simulation to fail
        t_max = 2 * 3600.0
        namelist_args = [
            ("time_stepping", "t_max", t_max),
            ("time_stepping", "dt_max", 200.0),
            ("time_stepping", "dt_min", 200.0),
            ("grid", "dz", 150.0),
            ("grid", "nz", 20),
            ("stats_io", "frequency", 720.0),
        ]

        kwargs_ref_model = Dict(
            :y_names => [["u_mean", "v_mean"]],
            :y_dir => scm_dirs,
            :scm_dir => scm_dirs,
            :case_name => ["Rico"],
            :t_start => [t_max - 3600],
            :t_end => [t_max],
            :Σ_t_start => [t_max - 2.0 * 3600],
            :Σ_t_end => [t_max],
            :namelist_args => repeat([namelist_args], 2),
        )
        ref_models = construct_reference_models(kwargs_ref_model)
        @test_logs (:warn,) match_mode = :any run_reference_SCM.(ref_models, run_single_timestep = false)

        u_names = ["entrainment_factor", "dt_max", "τ_acnv_rai"]
        u = [0.15, 210.0, 2500.0]
        constraints = Dict(
            "entrainment_factor" => [bounded(0.0, 0.5)],
            "dt_max" => [bounded(200.0, 300.0)],
            "τ_acnv_rai" => [no_constraint()],
        )
        prior = construct_priors(constraints)
        ref_stats = ReferenceStatistics(ref_models; y_type = SCM(), Σ_type = SCM())

        res_dir, model_error = run_SCM_handler(ref_models[1], data_dir, u, u_names, namelist_args)
        @test model_error

        @test_logs (:warn,) (:error,) match_mode = :any precondition(
            u,
            prior,
            ref_models,
            ref_stats,
            namelist_args,
            max_counter = 1,
        )
    end

    @testset "Namelist modification" begin

        # Choose same SCM to speed computation
        data_dir = mktempdir()
        scm_dirs = [joinpath(data_dir, "Output.Bomex.000000")]
        case_name = "Bomex"
        t_max = 2 * 3600.0

        kwargs_ref_model = Dict(
            :y_names => [["u_mean", "v_mean"]],
            :y_dir => scm_dirs,
            :scm_dir => scm_dirs,
            :case_name => [case_name],
            :t_start => [t_max - 3600],
            :t_end => [t_max],
            :Σ_t_start => [t_max - 2.0 * 3600],
            :Σ_t_end => [t_max],
        )
        ref_models = construct_reference_models(kwargs_ref_model)
        run_reference_SCM.(ref_models, run_single_timestep = true)

        # ensure namelist generated with `run_reference_SCM` matches default namelist
        init_namelist_path = namelist_directory(scm_dir(ref_models[1]), ref_models[1])
        default_namelist = NameList.default_namelist(case_name, root = scm_dir(ref_models[1]))
        reference_namelist = JSON.parsefile(init_namelist_path)

        namelist_compare_entries = ["microphysics", "time_stepping", "stats_io", "grid", "thermodynamics"]
        for entry in namelist_compare_entries
            @test default_namelist[entry] == reference_namelist[entry]
        end
        # Test nested dictionary of closures with type instabilities from JSON parsing
        for (key, value) in default_namelist["turbulence"]["EDMF_PrognosticTKE"]
            @test Tuple(reference_namelist["turbulence"]["EDMF_PrognosticTKE"][key]) == Tuple(value)
        end

        # ensure namelist in a `run_SCM_handler` call is modified as expected
        u = [0.15, 0.52]
        u_names = ["entrainment_factor", "detrainment_factor"]
        # Test namelist modification for different nesting levels
        namelist_args = [
            ("thermodynamics", "quadrature_type", "gaussian"),
            ("grid", "dz", 150.0),
            ("grid", "nz", 20),
            ("turbulence", "EDMF_PrognosticTKE", "pressure_normalmode_adv_coeff", 0.0),
            ("turbulence", "EDMF_PrognosticTKE", "general_stochastic_ent_params", SA.SVector(0.2, 0.2, 0.01, 0.02)),
        ]

        res_dir, model_error = run_SCM_handler(ref_models[1], data_dir, u, u_names, namelist_args)

        run_scm_namelist_path = namelist_directory(res_dir, ref_models[1])
        run_scm_namelist = JSON.parsefile(run_scm_namelist_path)
        expected_run_scm_namelist = deepcopy(default_namelist)
        expected_run_scm_namelist["turbulence"]["EDMF_PrognosticTKE"]["entrainment_factor"] = 0.15
        expected_run_scm_namelist["turbulence"]["EDMF_PrognosticTKE"]["detrainment_factor"] = 0.52
        expected_run_scm_namelist["turbulence"]["EDMF_PrognosticTKE"]["pressure_normalmode_adv_coeff"] = 0.0
        expected_run_scm_namelist["turbulence"]["EDMF_PrognosticTKE"]["general_stochastic_ent_params"] =
            SA.SVector(0.2, 0.2, 0.01, 0.02)
        expected_run_scm_namelist["thermodynamics"]["quadrature_type"] = "gaussian"
        expected_run_scm_namelist["grid"]["nz"] = 20
        expected_run_scm_namelist["grid"]["dz"] = 150.0
        expected_run_scm_namelist["stats_io"]["calibrate_io"] = true

        for entry in namelist_compare_entries
            @test expected_run_scm_namelist[entry] == run_scm_namelist[entry]
        end
        # Test nested dictionary of closures with type instabilities from JSON parsing
        for (key, value) in expected_run_scm_namelist["turbulence"]["EDMF_PrognosticTKE"]
            @test Tuple(run_scm_namelist["turbulence"]["EDMF_PrognosticTKE"][key]) == Tuple(value)
        end
    end
end


text:
using Plots
using ImplicitEquations
using ForwardDiff

a, b = 1 / √2, 1 / √2
f(x, y) = x^2 + y^2 - 1
f(xy) = f(xy[1], xy[2])
g = xy -> ForwardDiff.gradient(f, xy)

"""
x' = x-a
y' = y-b
"""
function shift(f::Function, a, b)
    t = [
        1 0 a;
        0 1 b;
        0 0 1
    ]
    function tf(xy)
        x, y = xy
        tx, ty, tz = t * [x, y, 1]
        txy = [tx, ty]
        return f(txy)
    end
    tf(x, y) = tf([x, y])
    return tf
end

"""
x'' = x'
y'' = f_x(a,b)*x' + f_y(a,b)*y'
"""
function cvt(f, a, b)
    fx, fy = g([a, b])
    t = [
        1 0 0;
        -fx / fy 1 / fy 0;
        0 0 1
    ]
    t = t / fy
    function tf(xy)
        x, y = xy
        tx, ty, tz = t * [x, y, 1]
        txy = [tx, ty]
        return f(txy)
    end
    tf(x, y) = tf([x, y])
    return tf
end

plot(aspect_ratio = :equal)
plot!(Eq(f, 0))
sf = shift(f, a, b)
plot!(Eq(sf, 0))
plot!(Eq(cvt(sf, a, b), 0))

max_stars_repo_name: terasakisatoshi/juliaExer

max_issues_repo_name: terasakisatoshi/juliaExer

max_forks_repo_name: terasakisatoshi/juliaExer


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

"""

Parameters
----------
distval
    a charasteristic measure to skip element pair, 0..5 => near, 10+ => far
    5 means that distance of slave element midpoint and point to project
    is 5 times larger than length of element
"""
type Contact <: BoundaryProblem
    dimension :: Int
    rotate_normals :: Bool
    finite_sliding :: Bool
    friction :: Bool
    dual_basis :: Bool
    use_forwarddiff :: Bool
    minimum_active_set_size :: Int
    distval :: Float64
    remove_from_set :: Bool # allow removal of non-potential contact pairs
    store_fields :: Vector{AbstractString}
end

function Contact()
    default_fields = ["element area", "contact area", "weighted gap",
        "contact pressure", "active nodes", "inactive nodes", "stick nodes",
        "slip nodes", "complementarity condition", "contact error"]
    return Contact(-1, false, false, false, true, false, 0, 5.0, false, default_fields)
end

function get_unknown_field_name(problem::Problem{Contact})
    return "reaction force"
end

function get_formulation_type(problem::Problem{Contact})
    if problem.properties.use_forwarddiff
        return :forwarddiff
    else
        return :incremental
    end
end

function assemble!(problem::Problem{Contact}, time::Real)
    if problem.properties.dimension == -1
        problem.properties.dimension = dim = size(first(problem.elements), 1)
        info("assuming dimension of mesh tie surface is $dim")
        info("if this is wrong set is manually using problem.properties.dimension")
    end
    dimension = Val{problem.properties.dimension}
    finite_sliding = Val{problem.properties.finite_sliding}
    friction = Val{problem.properties.friction}
    use_forwarddiff = Val{problem.properties.use_forwarddiff}
    assemble!(problem, time, dimension, finite_sliding, friction, use_forwarddiff)
end

typealias ContactElements2D Union{Seg2}


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

name = "CutPursuitBuilder"
version = v"0.0.3"

# Collection of sources required to build CutPursuit
sources = [
    "https://github.com/FugroRoames/cut-pursuit.git" =>
    "f621b1a6850a195e145de2806e2c7fda09815641",
    "https://dl.bintray.com/boostorg/release/1.68.0/source/boost_1_68_0.tar.bz2" =>
    "7f6130bc3cf65f56a618888ce9d5ea704fa10b462be126ad053e80e553d6d8b7",
]

# Bash recipe for building across all platforms
script = raw"""
cd $WORKSPACE/srcdir/boost_1_68_0/
./bootstrap.sh --prefix=$prefix
apk add python-dev
./b2 -j$nproc --layout=system install

cd $WORKSPACE/srcdir/cut-pursuit/src
mkdir build && cd build
cmake -DCMAKE_INSTALL_PREFIX=$prefix -DCMAKE_TOOLCHAIN_FILE=/opt/$target/$target.toolchain -DJlCxx_DIR=../../../destdir/lib/cmake/JlCxx/ ..

# NOTE: This is a hack to remove -std=c++14 which is not supported by gcc-4.8
echo "# CMAKE generated file: DO NOT EDIT!
# Generated by "Unix Makefiles" Generator, CMake Version 3.11

# compile CXX with /opt/x86_64-linux-gnu/bin/x86_64-linux-gnu-g++
CXX_FLAGS =  -Wall -std=c++11 -fopenmp -O3 -fPIC  -std=gnu++1y

CXX_DEFINES = -DJULIA_ENABLE_THREADING -Dcpjl_EXPORTS

CXX_INCLUDES = -I/workspace/srcdir/cut-pursuit/src/../include -isystem /workspace/destdir/include -isystem /workspace/destdir/include/julia" > CMakeFiles/cpjl.dir/flags.make

make -j${nproc}
cp libcpjl.so $prefix/lib
"""

# These are the platforms we will build for by default, unless further
# platforms are passed in on the command line
# platforms = supported_platforms()
platforms = [
    Linux(:x86_64, libc=:glibc)
]

# The products that we will ensure are always built
products(prefix) = [
    LibraryProduct(prefix, "libcpjl", :libcpjl),
    LibraryProduct(prefix, "libboost_math_tr1f", :math_tr1f),
    LibraryProduct(prefix, "libboost_thread", :thread),
    LibraryProduct(prefix, "libboost_unit_test_framework", :unit_test_framework),
    LibraryProduct(prefix, "libboost_type_erasure", :type_erasure),
    LibraryProduct(prefix, "libboost_chrono", :chrono),
    LibraryProduct(prefix, "libboost_math_c99l", :math_c99l),
    LibraryProduct(prefix, "libboost_locale", :locale),
    LibraryProduct(prefix, "libboost_program_options", :program_options),
    LibraryProduct(prefix, "libboost_date_time", :date_time),
    LibraryProduct(prefix, "libboost_graph", :graph),
    LibraryProduct(prefix, "libboost_signals", :signals),
    LibraryProduct(prefix, "libboost_iostreams", :iostreams),
    LibraryProduct(prefix, "libboost_stacktrace_addr2line", :stacktrace_addr2line),
    LibraryProduct(prefix, "libboost_system", :system),
    LibraryProduct(prefix, "libboost_wave", :wave),
    LibraryProduct(prefix, "libboost_wserialization", :wserialization),
    LibraryProduct(prefix, "libboost_math_tr1l", :math_tr1l),
    LibraryProduct(prefix, "libboost_math_tr1", :math_tr1),
    LibraryProduct(prefix, "libboost_filesystem", :filesystem),
    LibraryProduct(prefix, "libboost_random", :random),
    LibraryProduct(prefix, "libboost_coroutine", :coroutine),
    LibraryProduct(prefix, "libboost_serialization", :serialization),
    LibraryProduct(prefix, "libboost_context", :context),
    LibraryProduct(prefix, "libboost_container", :container),
    LibraryProduct(prefix, "libboost_stacktrace_noop", :stacktrace_noop),
    LibraryProduct(prefix, "libboost_contract", :contract),
    LibraryProduct(prefix, "libboost_prg_exec_monitor", :prg_exec_monitor),
    LibraryProduct(prefix, "libboost_regex", :regex),
    LibraryProduct(prefix, "libboost_log_setup", :log_setup),
    LibraryProduct(prefix, "libboost_math_c99", :math_c99),
    LibraryProduct(prefix, "libboost_timer", :timer),
    LibraryProduct(prefix, "libboost_stacktrace_basic", :stacktrace_basic),
    LibraryProduct(prefix, "libboost_math_c99f", :math_c99f),
    LibraryProduct(prefix, "libboost_log", :log),
    LibraryProduct(prefix, "libboost_atomic", :atomic),
    LibraryProduct(prefix, "libcxxwrap_julia", :libcxxwrap)
]

# Dependencies that must be installed before this package can be built
dependencies = [
    "https://github.com/JuliaInterop/libcxxwrap-julia/releases/download/v0.5.1/build_libcxxwrap-julia-1.0.v0.5.1.jl",
    "https://github.com/JuliaPackaging/JuliaBuilder/releases/download/v1.0.0-2/build_Julia.v1.0.0.jl",
    #"https://github.com/twadleigh/BoostBuilder/releases/download/v1.68.0-4/build_Boost.v1.68.0.jl"
]

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

max_stars_repo_name: JoshChristie/CutPursuitBuilder

max_issues_repo_name: JoshChristie/CutPursuitBuilder

max_forks_repo_name: JoshChristie/CutPursuitBuilder