text:
#
# Calculates requests per second from client to server without response, but
# using only naked sockets.
# Results: rps: 8688.474736099924, time: 115.09500002861023
#
type T
  x::Int
  y::Int
  c::UInt32
end
server = listen(2000)
sock   = accept(server)
function run()
  local t::T = T(rand(1:100), rand(1:100), UInt32(1234))
  for i=1:1000000 serialize(sock, t) end
end

max_stars_repo_name: tmptrash/jevo

max_issues_repo_name: tmptrash/jevo

max_forks_repo_name: tmptrash/jevo


text:
require("Ito")
using Base.Test
using Ito.Statistics
using StatsBase

v = [ 3.0, 4.0, 5.0, 2.0, 3.0, 4.0, 5.0, 6.0, 4.0, 7.0 ] #data
w = ones(10) #weights
n=10

@test_approx_eq mean(v, weights(w))  4.3
@test_approx_eq var(v, w)  2.23333333333
@test_approx_eq std(v, w) 1.4944341180973264
#answers verified via excel and quantlib
@test_approx_eq skewness(v, weights(w))*(sqrt(n*(n-1))/(n-2))  0.359543071407
#Limited by excel precision for comparison.
#Note: Wolfram alpha produces unadjusted value, Octave produces strange result
#@test_approx_eq_eps kurtosis(v, w)  -0.151799637209 10e-9

@test_approx_eq regret(v, w, 4) 2.2222222222222223
@test_approx_eq average_shortfall(v, w, 4) 1.3333333333333333

w=[1:10]
@test_approx_eq mean(v, weights(w)) 4.76363636363636

#TODO Lots more tests with realistic data

ext: jl


text:
using OrdinaryDiffEq, Test, LinearAlgebra
prob = ODEProblem((u,p,t)->1.01u,0.5,(0.0,1.0), save_start=false, save_end=false)
function prob_func(prob,i,repeat)
  remake(prob,u0=rand()*prob.u0)
end
ensemble_prob = EnsembleProblem(prob,prob_func=prob_func)
sim = solve(ensemble_prob,Tsit5(),EnsembleThreads(),trajectories=10, save_everystep=false)
@test ndims(sim) == 2
@test length(sim) == 10
@test eltype(sim.u) <: ODESolution
ts = 0.0:0.1:1.0

using SciMLBase.EnsembleAnalysis
sim = solve(ensemble_prob,Tsit5(),EnsembleThreads(),trajectories=10,saveat=0.1)
timeseries_point_meancov(sim,ts)

function prob_sol(_p)
  prob = ODEProblem((u,p,t)->p.*u, _p, (0.0,1.0), _p, save_start=false, save_end=false)
  sim = solve(prob,Tsit5())
end
mapres = SciMLBase.responsible_map(prob_sol, [0.5, diagm([1.0, 1.0])])
@test length(mapres) == 2


text:
@doc raw"""
```
binarize(Yen(), img)
```

Computes the binarization threshold value using Yen's maximum correlation criterion for
bilevel thresholding.

# Output

Returns the binarized image as an `Array{Gray{Bool},2}`.

# Details

This algorithm uses the concept of *entropic correlation* of a gray level histogram to produce a threshold
value.

Let ``f_1, f_2, \ldots, f_I`` be the frequencies in the various bins of the
histogram and ``I`` the number of bins. With ``N = \sum_{i=1}^{I}f_i``, let
``p_i = \frac{f_i}{N}`` (``i = 1, \ldots, I``) denote the probability
distribution of gray levels. From this distribution one derives two additional
distributions. The first defined for discrete values ``1`` to ``s`` and the
other, from ``s+1`` to ``I``. These distributions are

```math
A: \frac{p_1}{P_s}, \frac{p_2}{P_s}, \ldots, \frac{p_s}{P_s}
\quad \text{and} \quad
B: \frac{p_{s+1}}{1-P_s}, \ldots, \frac{p_n}{1-P_s}
\quad \text{where} \quad
P_s = \sum_{i=1}^{s}p_i.
```
The entropic correlations associated with each distribution are

```math
C(A) = -\ln \sum_{i=1}^{s} \left( \frac{p_i}{P_s} \right)^2 \quad \text{and} \quad C(B) = -\ln \sum_{i=s+1}^{I} \left( \frac{p_i}{1 - P_s} \right)^2.
```

Combining these two entropic correlation functions we have

```math
\psi(s) = -\ln \sum_{i=1}^{s} \left( \frac{p_i}{P_s} \right)^2 -\ln \sum_{i=s+1}^{I} \left( \frac{p_i}{1 - P_s} \right)^2.
```
Finding the discrete value ``s`` which maximises the function ``\psi(s)`` produces
the sought-after threshold value (i.e. the bin which determines the threshold).

# Arguments

The function argument is described in more detail below.

##  `img`

An `AbstractArray` representing an image. The image is automatically converted
to `Gray` in order to construct the requisite graylevel histogram.

# Example

Binarize the "cameraman" image in the `TestImages` package.

```julia
using TestImages, ImageBinarization

img = testimage("cameraman")
img_binary = binarize(Yen(), img)
```

# Reference

1. Yen JC, Chang FJ, Chang S (1995), “A New Criterion for Automatic Multilevel Thresholding”, IEEE Trans. on Image Processing 4 (3): 370-378, [doi:10.1109/83.366472](https://doi.org/10.1109/83.366472)
"""
function binarize(algorithm::Yen,  img::AbstractArray{T,2}) where T <: Colorant
  img₀₁ = zeros(Gray{Bool}, axes(img))
  edges, counts = build_histogram(img,  256)
  t = find_threshold(HistogramThresholding.Yen(), counts[1:end], edges)
  for i in CartesianIndices(img)
    img₀₁[i] = img[i] < t ? 0 : 1
  end
  img₀₁
end


text:
@testset "ValuePrecedence2Reification: $(fct_type), dimension $(dim), $(T)" for fct_type in ["vector of variables", "vector affine function"], dim in [2, 3], T in [Int, Float64]
    base_model = if T == Int
        IntReificationEqualToModel{T}()
    elseif T == Float64
        FloatReificationEqualToModel{T}()
    else
        @assert false
    end
    mock = MOIU.MockOptimizer(base_model)
    model = COIB.ValuePrecedence2Reification{T}(mock)

    @test MOI.supports_constraint(model, MOI.VariableIndex, MOI.Integer)
    @test MOI.supports_constraint(
        model,
        MOI.VectorAffineFunction{T},
        CP.Reification{MOI.EqualTo{T}},
    )
    @test MOI.supports_constraint(
        model,
        MOI.ScalarAffineFunction{T},
        MOI.LessThan{T},
    )
    @test MOIB.supports_bridging_constraint(
        model,
        MOI.VectorAffineFunction{T},
        CP.ValuePrecedence{T},
    )

    if T == Int
        x, _ = MOI.add_constrained_variables(model, [MOI.Integer() for _ in 1:dim])
    elseif T == Float64
        x = MOI.add_variables(model, dim)
    end

    fct = if fct_type == "vector of variables"
        MOI.VectorOfVariables(x)
    elseif fct_type == "vector affine function"
        MOIU.vectorize(x)
    else
        @assert false
    end
    c = MOI.add_constraint(model, fct, CP.ValuePrecedence(T(2), T(4), dim))

    for i in 1:dim
        @test MOI.is_valid(model, x[i])
    end
    @test MOI.is_valid(model, c)

    bridge = first(MOIBC.bridges(model))[2]

    @testset "Bridge properties" begin
        @test MOIBC.concrete_bridge_type(typeof(bridge), MOI.VariableIndex, CP.SymmetricAllDifferent) == typeof(bridge)
        @test MOIB.added_constrained_variable_types(typeof(bridge)) == [(MOI.ZeroOne,)]
        @test MOIB.added_constraint_types(typeof(bridge)) == [
            (MOI.VectorAffineFunction{T}, CP.Reification{MOI.EqualTo{T}}),
            (MOI.ScalarAffineFunction{T}, MOI.LessThan{T}),
        ]

        @test MOI.get(bridge, MOI.NumberOfVariables()) == 2 * (dim - 1)
        @test MOI.get(bridge, MOI.NumberOfConstraints{MOI.VariableIndex, MOI.ZeroOne}()) == 2 * (dim - 1)
        @test MOI.get(bridge, MOI.NumberOfConstraints{MOI.VectorAffineFunction{T}, CP.Reification{MOI.EqualTo{T}}}()) == 2 * (dim - 1)
        @test MOI.get(bridge, MOI.NumberOfConstraints{MOI.ScalarAffineFunction{T}, MOI.LessThan{T}}()) == dim - 1

        @test Set(MOI.get(bridge, MOI.ListOfVariableIndices())) == Set([bridge.vars_reif_precv..., bridge.vars_reif_value...])
        @test Set(MOI.get(bridge, MOI.ListOfConstraintIndices{MOI.VariableIndex, MOI.ZeroOne}())) == Set([bridge.vars_reif_precv_bin..., bridge.vars_reif_value_bin...])
        @test Set(MOI.get(bridge, MOI.ListOfConstraintIndices{MOI.VectorAffineFunction{T}, CP.Reification{MOI.EqualTo{T}}}())) == Set([bridge.cons_reif_value..., bridge.cons_reif_precv...])
        @test Set(MOI.get(bridge, MOI.ListOfConstraintIndices{MOI.ScalarAffineFunction{T}, MOI.LessThan{T}}())) == Set(bridge.cons_implication)
    end

    @testset "Set of variables" begin
        @test length(bridge.vars_reif_precv) == dim - 1
        @test length(bridge.vars_reif_value) == dim - 1
        @test length(bridge.vars_reif_precv_bin) == dim - 1
        @test length(bridge.vars_reif_value_bin) == dim - 1

        for i in 1:(dim - 1)
            @test MOI.is_valid(model, bridge.vars_reif_precv[i])
            @test MOI.is_valid(model, bridge.vars_reif_value[i])
            @test MOI.is_valid(model, bridge.vars_reif_precv_bin[i])
            @test MOI.is_valid(model, bridge.vars_reif_value_bin[i])

            @test MOI.get(model, MOI.ConstraintSet(), bridge.vars_reif_precv_bin[i]) == MOI.ZeroOne()
            @test MOI.get(model, MOI.ConstraintFunction(), bridge.vars_reif_precv_bin[i]) == bridge.vars_reif_precv[i]

            @test MOI.get(model, MOI.ConstraintSet(), bridge.vars_reif_value_bin[i]) == MOI.ZeroOne()
            @test MOI.get(model, MOI.ConstraintFunction(), bridge.vars_reif_value_bin[i]) == bridge.vars_reif_value[i]
        end
    end

    @testset "Reify equality to value" begin
        @test length(bridge.cons_reif_value) == dim - 1

        for i in 1:(dim - 1)
            @test MOI.is_valid(model, bridge.cons_reif_value[i])
            @test MOI.get(model, MOI.ConstraintSet(), bridge.cons_reif_value[i]) == CP.Reification(MOI.EqualTo(zero(T)))
            f = MOI.get(model, MOI.ConstraintFunction(), bridge.cons_reif_value[i])

            @test length(f.terms) == 2
            @test f.constants == [zero(T), -T(4)]

            t1 = f.terms[1]
            @test t1.output_index == 1
            @test t1.scalar_term.coefficient === one(T)
            @test t1.scalar_term.variable == bridge.vars_reif_value[i]

            t2 = f.terms[2]
            @test t2.output_index == 2
            @test t2.scalar_term.coefficient === one(T)
            @test t2.scalar_term.variable == x[i + 1]
        end
    end

    @testset "Reify equality to before" begin
        @test length(bridge.cons_reif_precv) == dim - 1

        for i in 1:(dim - 1)
            @test MOI.is_valid(model, bridge.cons_reif_precv[i])
            @test MOI.get(model, MOI.ConstraintSet(), bridge.cons_reif_precv[i]) == CP.Reification(MOI.EqualTo(zero(T)))
            f = MOI.get(model, MOI.ConstraintFunction(), bridge.cons_reif_precv[i])

            @test length(f.terms) == 2
            @test f.constants == [zero(T), -T(2)]

            t1 = f.terms[1]
            @test t1.output_index == 1
            @test t1.scalar_term.coefficient === one(T)
            @test t1.scalar_term.variable == bridge.vars_reif_precv[i]

            t2 = f.terms[2]
            @test t2.output_index == 2
            @test t2.scalar_term.coefficient === one(T)
            @test t2.scalar_term.variable == x[i]
        end
    end

    @testset "Implications" begin
        @test length(bridge.cons_implication) == dim - 1

        for i in 2:dim
            @test MOI.is_valid(model, bridge.cons_implication[i - 1])
            @test MOI.get(model, MOI.ConstraintSet(), bridge.cons_implication[i - 1]) == MOI.LessThan(zero(T))
            f = MOI.get(model, MOI.ConstraintFunction(), bridge.cons_implication[i - 1])
            @test length(f.terms) == i
            @test abs(f.constant) === zero(T)

            t1 = f.terms[1]
            t1.coefficient === one(T)
            t1.variable == bridge.vars_reif_value[i - 1]

            for j in 2:i
                t = f.terms[j]
                @test t.coefficient === -one(T)
                @test t.variable == bridge.vars_reif_precv[j - 1]
            end
        end
    end
end


text:
#!/usr/bin/env julia

# This script sometimes (but not always) crashes Julia v0.6.0, first triggering the warning:
#     WARNING: An error occurred during inference. Type inference is now partially disabled.
#     Base.MethodError(f=typeof(Core.Inference.convert)(), args=(Base.AssertionError, "invalid age range update"), world=0x0000000000000ac2)
# The issue resembles
#     https://github.com/JuliaLang/julia/issues/23200
#     https://github.com/JuliaLang/julia/issues/22355

import Celeste.AccuracyBenchmark
import Celeste.GalsimBenchmark

test_case_names = String[]
for arg in ARGS
    if ismatch(r"^--test_case_names", arg)
        # test case name: function name in test_case_definitions.py, or CL_DESCR
        # field from FITS header
        test_case_names =
            split(replace(arg, "--test_case_names=", ""), ",")
    else
        println("Usage: run_galsim_benchmark.jl " *
                "[--test_case_names=test1,test2,...]")
        error("Invalid argument: ", arg)
    end
end

srand(12345)

truth_catalog, single_predictions = GalsimBenchmark.run_benchmarks(
    test_case_names=test_case_names,
    joint_inference=false,
    verbose=true)
unused, joint_predictions = GalsimBenchmark.run_benchmarks(
    test_case_names=test_case_names,
    joint_inference=true,
    verbose=true)

score_df = AccuracyBenchmark.score_predictions(
    truth_catalog,
    [single_predictions, joint_predictions],
)

println("\n\n=========================\n")
println("overall scores:\n")
display(score_df)
println()


text:
allPerms = Set(([1,2,3], [2,3,1], [3,1,2], [2,1,3], [1,3,2], [3,2,1]))
allSigns = Set(([1,1,1], [-1,1,1], [-1,-1,1], [-1,-1,-1]))

mutable struct Transformation
  offset :: Vector{Int}
  perm   :: Vector{Int}
  sign   :: Vector{Int}
end

# Compose two coordinate transformations such that:
#   1->3 = 2->3 ◌ 1->2
function compose(o12::Transformation, o23::Transformation) :: Transformation
  return Transformation(
    o23.offset .+ transform(o12.offset,o23),
    o12.perm[o23.perm],
    o23.sign .* o12.sign[o23.perm]
  )
end

# Apply transformation to coordinate vector
transform(x::Vector{Int},o::Transformation) :: Vector{Int} = o.offset .+ o.sign .* x[o.perm]

# Get relative transformation between two sets of points
function relativeTransformation(sc1::Vector{Vector{Int}}, sc2::Vector{Vector{Int}}) :: Transformation
  for perm in allPerms
    for sign in allSigns

    end
  end
end

function parseLine(s::String)
  (s == "") && (return nothing)
  (s[1:3] == "---") && (return "new")
  return map(x -> parse(Int,x), split(s,','))
end

function readInput()
  data = map(line -> parseLine(line), eachline("inputs/day19.txt"))

  i = 0
  scanners :: Vector{Vector{Vector{Int}}} = []
  for dd in data
    if dd == "new"
      i += 1
      push!(scanners,[])
    elseif typeof(dd) === Vector{Int}
      push!(scanners[i],dd)
    end
  end
  return scanners
end

scanners = readInput()


text:
# Test the lower bound methods
@testset "Lower Bound" begin
    # initialize model and 3 test variables
    m = InfiniteModel()
    gvref1 = @hold_variable(m, var1)
    gvref2 = @hold_variable(m, var2 >= 0)
    gvref3 = @hold_variable(m, var3 == 0)
    vref1 = dispatch_variable_ref(gvref1)
    vref2 = dispatch_variable_ref(gvref2)
    vref3 = dispatch_variable_ref(gvref3)
    # _temp_constraint_ref
    @testset "_temp_constraint_ref" begin
        cindex = ConstraintIndex(1)
        expected = InfOptConstraintRef(m, cindex)
        @test InfiniteOpt._temp_constraint_ref(m, cindex) == expected
    end
    # is_fixed
    @testset "JuMP.is_fixed" begin
        @test is_fixed(vref3)
        @test !is_fixed(vref2)
    end
    # has_lower_bound
    @testset "JuMP.has_lower_bound" begin
        @test !has_lower_bound(vref1)
        @test has_lower_bound(vref2)
        @test has_lower_bound(gvref2)
    end
    # lower_bound
    @testset "JuMP.lower_bound" begin
        @test lower_bound(vref2) == 0
        @test lower_bound(gvref2) == 0
        @test_throws ErrorException lower_bound(vref1)
    end
    # _lower_bound_index
    @testset "InfiniteOpt._lower_bound_index" begin
        @test InfiniteOpt._lower_bound_index(vref2) == ConstraintIndex(1)
        @test_throws ErrorException InfiniteOpt._lower_bound_index(vref1)
    end
    # _set_lower_bound_index
    @testset "_set_lower_bound_index" begin
        # test function
        @test isa(InfiniteOpt._set_lower_bound_index(vref2, ConstraintIndex(2)),
                  Nothing)
        @test InfiniteOpt._lower_bound_index(vref2) == ConstraintIndex(2)
        # undo changes
        @test isa(InfiniteOpt._set_lower_bound_index(vref2, ConstraintIndex(1)),
                  Nothing)
    end
    # LowerBoundRef
    @testset "JuMP.LowerBoundRef" begin
        # prepare constraint reference
        expected = InfOptConstraintRef(m, ConstraintIndex(1))
        # test for finite variable
        @test LowerBoundRef(vref2) == expected
        @test LowerBoundRef(gvref2) == expected
        # prepare infinite variable
        @infinite_parameter(m, t in [0, 1])
        @infinite_variable(m, ivref(t) >= 0)
        expected = InfOptConstraintRef(m, ConstraintIndex(3))
        @test LowerBoundRef(ivref) == expected
    end
    # set_lower_bound
    @testset "JuMP.set_lower_bound" begin
        # test adding new lower bound
        @test isa(set_lower_bound(vref1, 1), Nothing)
        @test has_lower_bound(vref1)
        @test lower_bound(vref1) == 1
        cref = LowerBoundRef(vref1)
        @test InfiniteOpt._data_object(cref).is_info_constraint
        @test !optimizer_model_ready(m)
        @test InfiniteOpt._core_constraint_object(cref) == ScalarConstraint(gvref1, MOI.GreaterThan(1.))
        # test changing lower bound
        @test isa(set_lower_bound(vref2, 1.5), Nothing)
        @test has_lower_bound(vref2)
        @test lower_bound(gvref2) == 1.5
        cref = LowerBoundRef(vref2)
        @test InfiniteOpt._data_object(cref).is_info_constraint
        @test !optimizer_model_ready(m)
        @test InfiniteOpt._core_constraint_object(cref) == ScalarConstraint(gvref2, MOI.GreaterThan(1.5))
        # test fixed variable error
        @test_throws AssertionError set_lower_bound(vref3, 0)
    end
end

# Test the upper bound methods
@testset "Upper Bound" begin
    # initialize model and 3 test variables
    m = InfiniteModel()
    gvref1 = @hold_variable(m, var1)
    gvref2 = @hold_variable(m, var2 <= 0)
    gvref3 = @hold_variable(m, var3 == 0)
    vref1 = dispatch_variable_ref(gvref1)
    vref2 = dispatch_variable_ref(gvref2)
    vref3 = dispatch_variable_ref(gvref3)
    # is_fixed
    @testset "JuMP.is_fixed" begin
        @test is_fixed(vref3)
        @test !is_fixed(vref2)
    end
    # has_upper_bound
    @testset "JuMP.has_upper_bound" begin
        @test !has_upper_bound(vref1)
        @test has_upper_bound(vref2)
        @test has_upper_bound(gvref2)
    end
    # upper_bound
    @testset "JuMP.upper_bound" begin
        @test upper_bound(vref2) == 0
        @test upper_bound(gvref2) == 0
        @test_throws ErrorException upper_bound(vref1)
    end
    # _upper_bound_index
    @testset "InfiniteOpt._upper_bound_index" begin
        @test InfiniteOpt._upper_bound_index(vref2) == ConstraintIndex(1)
        @test_throws ErrorException InfiniteOpt._upper_bound_index(vref1)
    end
    # _set_upper_bound_index
    @testset "_set_upper_bound_index" begin
        # test function
        @test isa(InfiniteOpt._set_upper_bound_index(vref2, ConstraintIndex(2)),
                  Nothing)
        @test InfiniteOpt._upper_bound_index(vref2) == ConstraintIndex(2)
        # undo changes
        @test isa(InfiniteOpt._set_upper_bound_index(vref2, ConstraintIndex(1)),
                  Nothing)
    end
    # UpperBoundRef
    @testset "JuMP.UpperBoundRef" begin
        # prepare constraint reference
        expected = InfOptConstraintRef(m, ConstraintIndex(1))
        # test for finite variable
        @test UpperBoundRef(vref2) == expected
        @test UpperBoundRef(gvref2) == expected
        # prepare infinite variable
        @infinite_parameter(m, t in [0, 1])
        @infinite_variable(m, ivref(t) <= 0)
        expected = InfOptConstraintRef(m, ConstraintIndex(3))
        @test UpperBoundRef(ivref) == expected
    end
    # set_upper_bound
    @testset "JuMP.set_upper_bound" begin
        # test adding new upper bound
        @test isa(set_upper_bound(vref1, 1), Nothing)
        @test has_upper_bound(vref1)
        @test upper_bound(vref1) == 1
        cref = UpperBoundRef(vref1)
        @test InfiniteOpt._data_object(cref).is_info_constraint
        @test !optimizer_model_ready(m)
        @test InfiniteOpt._core_constraint_object(cref) == ScalarConstraint(gvref1, MOI.LessThan(1.))
        # test changing upper bound
        @test isa(set_upper_bound(vref2, 1.5), Nothing)
        @test has_upper_bound(vref2)
        @test upper_bound(gvref2) == 1.5
        cref = UpperBoundRef(vref2)
        @test InfiniteOpt._data_object(cref).is_info_constraint
        @test !optimizer_model_ready(m)
        @test InfiniteOpt._core_constraint_object(cref) == ScalarConstraint(gvref2, MOI.LessThan(1.5))
        # test fixed variable error
        @test_throws AssertionError set_upper_bound(vref3, 0)
    end
end

# Test the fix methods
@testset "Fix" begin
    # initialize model and 4 test variables
    m = InfiniteModel()
    gvref1 = @hold_variable(m, var1)
    gvref2 = @hold_variable(m, var2 == 0)
    gvref3 = @hold_variable(m, var3)
    gvref4 = @hold_variable(m, var4)
    vref1 = dispatch_variable_ref(gvref1)
    vref2 = dispatch_variable_ref(gvref2)
    vref3 = dispatch_variable_ref(gvref3)
    vref4 = dispatch_variable_ref(gvref4)
    # is_fixed
    @testset "is_fixed" begin
        @test !is_fixed(vref1)
        @test is_fixed(vref2)
        @test is_fixed(gvref2)
    end
    # fix_value
    @testset "JuMP.fix_value" begin
        @test fix_value(vref2) == 0
        @test fix_value(gvref2) == 0
        @test_throws ErrorException fix_value(vref1)
    end
    # _fix_index
    @testset "InfiniteOpt._fix_index" begin
        @test InfiniteOpt._fix_index(vref2) == ConstraintIndex(1)
        @test_throws ErrorException InfiniteOpt._fix_index(vref1)
    end
    # _set_fix_index
    @testset "_set_fix_index" begin
        # test function
        @test isa(InfiniteOpt._set_fix_index(vref2, ConstraintIndex(2)), Nothing)
        @test InfiniteOpt._fix_index(vref2) == ConstraintIndex(2)
        # undo changes
        @test isa(InfiniteOpt._set_fix_index(vref2, ConstraintIndex(1)), Nothing)
    end
    # FixRef
    @testset "JuMP.FixRef" begin
        # prepare constraint reference
        expected = InfOptConstraintRef(m, ConstraintIndex(1))
        # test for finite variable
        @test FixRef(vref2) == expected
        # prepare infinite variable
        @infinite_parameter(m, t in [0, 1])
        @infinite_variable(m, ivref(t) == 0)
        expected = InfOptConstraintRef(m, ConstraintIndex(2))
        # test for infinite variable
        @test FixRef(ivref) == expected
    end
    # fix
    @testset "JuMP.fix" begin
        # test adding new fix
        @test isa(fix(vref1, 1), Nothing)
        @test is_fixed(vref1)
        @test fix_value(vref1) == 1
        cref = FixRef(vref1)
        @test InfiniteOpt._data_object(cref).is_info_constraint
        @test !optimizer_model_ready(m)
        @test InfiniteOpt._core_constraint_object(cref) == ScalarConstraint(gvref1, MOI.EqualTo(1.))
        # test changing fix
        @test isa(fix(gvref2, 1.5), Nothing)
        @test is_fixed(vref2)
        @test fix_value(vref2) == 1.5
        cref = FixRef(vref2)
        @test InfiniteOpt._data_object(cref).is_info_constraint
        @test !optimizer_model_ready(m)
        @test InfiniteOpt._core_constraint_object(cref) == ScalarConstraint(gvref2, MOI.EqualTo(1.5))
        # add lower and upper bounds to vars 3 and 4
        set_lower_bound(vref3, 0.0)
        set_upper_bound(vref4, 0.0)
        # test lower/upper bound  error
        @test_throws ErrorException fix(vref3, 0)
        @test_throws ErrorException fix(vref4, 0)
        # test forcing with lower
        @test isa(fix(vref3, 1, force = true), Nothing)
        @test is_fixed(vref3)
        @test fix_value(vref3) == 1
        cref = FixRef(vref3)
        @test InfiniteOpt._data_object(cref).is_info_constraint
        @test !optimizer_model_ready(m)
        @test InfiniteOpt._core_constraint_object(cref) == ScalarConstraint(gvref3, MOI.EqualTo(1.))
        # test forcing with upper
        @test isa(fix(gvref4, 1.5, force = true), Nothing)
        @test is_fixed(vref4)
        @test fix_value(vref4) == 1.5
        cref = FixRef(vref4)
        @test InfiniteOpt._data_object(cref).is_info_constraint
        @test !optimizer_model_ready(m)
        @test InfiniteOpt._core_constraint_object(cref) == ScalarConstraint(gvref4, MOI.EqualTo(1.5))
    end
end

# Test the start value methods
@testset "Start Value" begin
    # initialize model and 4 test variables
    m = InfiniteModel()
    gvref = @hold_variable(m, var1, start = 0)
    vref = dispatch_variable_ref(gvref)
    gvref2 = @hold_variable(m, var2)
    vref2 = dispatch_variable_ref(gvref2)
    @infinite_parameter(m, t in [0, 1])
    @infinite_variable(m, inf(t), start = 0)
    dinf = dispatch_variable_ref(inf)
    @infinite_variable(m, inf2(t))
    dinf2 = dispatch_variable_ref(inf2)
    # start_value
    @testset "JuMP.start_value" begin
        @test start_value(vref) == 0
        @test start_value(gvref) == 0
        @test start_value(vref2) isa Nothing
        @test start_value(gvref2) isa Nothing
        @test_throws ErrorException start_value(inf)
        @test_throws ErrorException start_value(dinf)
    end
    # set_start_value
    @testset "JuMP.set_start_value" begin
        @test isa(set_start_value(vref, 1.5), Nothing)
        @test start_value(vref) == 1.5
        @test !optimizer_model_ready(m)
        @test isa(set_start_value(gvref, 1), Nothing)
        @test_throws ErrorException set_start_value(inf, 0)
        @test_throws ErrorException set_start_value(dinf, 0)
    end
    # start_value_function
    @testset "start_value_function" begin
        @test start_value_function(inf)([1]) == 0
        @test start_value_function(dinf)([1]) == 0
        @test InfiniteOpt._is_vector_start(dinf)
        @test start_value_function(inf2) isa Nothing
        @test start_value_function(dinf2) isa Nothing
    end
    # set_start_value_function
    @testset "set_start_value_function" begin
        @test set_start_value_function(inf, 1.5) isa Nothing
        @test start_value_function(inf)([0]) == 1.5
        @test InfiniteOpt._is_vector_start(dinf)
        @test !optimizer_model_ready(m)
        func = (a) -> 1
        @test set_start_value_function(inf, func) isa Nothing
        @test start_value_function(inf)(0) == 1
        @test !InfiniteOpt._is_vector_start(dinf)
    end
    # reset_start_value_function
    @testset "reset_start_value_function" begin
        @test reset_start_value_function(inf) isa Nothing
        @test start_value_function(inf) isa Nothing
        @test InfiniteOpt._is_vector_start(dinf)
        @test !optimizer_model_ready(m)
    end
end

# Test the binary methods
@testset "Binary" begin
    # initialize model and 3 test variables
    m = InfiniteModel()
    gvref1 = @hold_variable(m, var1)
    gvref2 = @hold_variable(m, var2, Bin)
    gvref3 = @hold_variable(m, var3, Int)
    vref1 = dispatch_variable_ref(gvref1)
    vref2 = dispatch_variable_ref(gvref2)
    vref3 = dispatch_variable_ref(gvref3)
    # is_binary
    @testset "JuMP.is_binary" begin
        @test !is_binary(vref1)
        @test is_binary(vref2)
        @test is_binary(gvref2)
    end
    # _binary_index
    @testset "InfiniteOpt._binary_index" begin
        @test InfiniteOpt._binary_index(vref2) == ConstraintIndex(1)
        @test_throws ErrorException InfiniteOpt._binary_index(vref1)
    end
    # _set_binary_index
    @testset "_set_binary_index" begin
        # test function
        @test isa(InfiniteOpt._set_binary_index(vref2, ConstraintIndex(2)), Nothing)
        @test InfiniteOpt._binary_index(vref2) == ConstraintIndex(2)
        # undo changes
        @test isa(InfiniteOpt._set_binary_index(vref2, ConstraintIndex(1)), Nothing)
    end
    # BinaryRef
    @testset "BinaryRef" begin
        # prepare constraint reference
        expected = InfOptConstraintRef(m, ConstraintIndex(1))
        # test for finite variable
        @test BinaryRef(vref2) == expected
        # prepare infinite variable
        @infinite_parameter(m, t in [0, 1])
        @infinite_variable(m, ivref(t), Bin)
        expected = InfOptConstraintRef(m, ConstraintIndex(3))
        # test for infinite variable
        @test BinaryRef(ivref) == expected
    end
    # set_binary
    @testset "JuMP.set_binary" begin
        # test adding binary
        @test isa(set_binary(vref1), Nothing)
        @test is_binary(vref1)
        cref = BinaryRef(vref1)
        @test InfiniteOpt._data_object(cref).is_info_constraint
        @test !optimizer_model_ready(m)
        @test InfiniteOpt._core_constraint_object(cref) == ScalarConstraint(gvref1, MOI.ZeroOne())
        # test setting binary again
        @test isa(set_binary(gvref2), Nothing)
        @test is_binary(vref2)
        cref = BinaryRef(vref2)
        @test InfiniteOpt._data_object(cref).is_info_constraint
        @test !optimizer_model_ready(m)
        @test InfiniteOpt._core_constraint_object(cref) == ScalarConstraint(gvref2, MOI.ZeroOne())
        # test integer variable error
        @test_throws ErrorException set_binary(vref3)
    end
end

# Test the integer methods
@testset "Integer" begin
    # initialize model and 3 test variables
    m = InfiniteModel()
    gvref1 = @hold_variable(m, var1)
    gvref2 = @hold_variable(m, var2, Int)
    gvref3 = @hold_variable(m, var3, Bin)
    vref1 = dispatch_variable_ref(gvref1)
    vref2 = dispatch_variable_ref(gvref2)
    vref3 = dispatch_variable_ref(gvref3)
    # is_integer
    @testset "JuMP.is_integer" begin
        @test !is_integer(vref1)
        @test is_integer(vref2)
        @test is_integer(gvref2)
    end
    # _integer_index
    @testset "InfiniteOpt._integer_index" begin
        @test InfiniteOpt._integer_index(vref2) == ConstraintIndex(1)
        @test_throws ErrorException InfiniteOpt._integer_index(vref1)
    end
    # _set_integer_index
    @testset "_set_integer_index" begin
        # test function
        @test isa(InfiniteOpt._set_integer_index(vref2, ConstraintIndex(2)), Nothing)
        @test InfiniteOpt._integer_index(vref2) == ConstraintIndex(2)
        # undo changes
        @test isa(InfiniteOpt._set_integer_index(vref2, ConstraintIndex(1)), Nothing)
    end
    # IntegerRef
    @testset "IntegerRef" begin
        # prepare constraint reference
        expected = InfOptConstraintRef(m, ConstraintIndex(1))
        # test for finite variable
        @test IntegerRef(vref2) == expected
        # prepare infinite variable
        @infinite_parameter(m, t in [0, 1])
        @infinite_variable(m, ivref(t), Int)
        expected = InfOptConstraintRef(m, ConstraintIndex(3))
        # test for infinite variable
        @test IntegerRef(ivref) == expected
    end
    # set_integer
    @testset "JuMP.set_integer" begin
        # test adding integer
        @test isa(set_integer(vref1), Nothing)
        @test is_integer(vref1)
        cref = IntegerRef(vref1)
        @test InfiniteOpt._data_object(cref).is_info_constraint
        @test !optimizer_model_ready(m)
        @test InfiniteOpt._core_constraint_object(cref) == ScalarConstraint(gvref1, MOI.Integer())
        # test setting integer again
        @test isa(set_integer(gvref2), Nothing)
        @test is_integer(vref2)
        cref = IntegerRef(vref2)
        @test InfiniteOpt._data_object(cref).is_info_constraint
        @test !optimizer_model_ready(m)
        @test InfiniteOpt._core_constraint_object(cref) == ScalarConstraint(gvref2, MOI.Integer())
        # test integer variable error
        @test_throws ErrorException set_integer(vref3)
    end
end

# Test relaxation
@testset "JuMP.relax_integrality" begin
    # setup the model
    m = InfiniteModel()
    @infinite_parameter(m, t in [0, 1])
    @infinite_variable(m, y(t), Int)
    @hold_variable(m, 0.25 <= x <= 1.25, Bin)
    @hold_variable(m, z, Bin)
    # test relaxing it
    result_store = []
    @test push!(result_store, relax_integrality(m)) isa Vector
    @test !is_integer(y)
    @test !is_binary(x)
    @test !is_binary(z)
    @test lower_bound(x) == 0.25
    @test upper_bound(x) == 1
    @test lower_bound(z) == 0
    @test upper_bound(z) == 1
    # test unrelaxing it
    @test result_store[1]() isa Nothing
    @test is_integer(y)
    @test is_binary(x)
    @test is_binary(z)
    @test lower_bound(x) == 0.25
    @test upper_bound(x) == 1.25
    @test !has_lower_bound(z)
    @test !has_upper_bound(z)
end


text:
using Documenter, Spectra

makedocs(sitename="Spectra documentation", format = Documenter.HTML(prettyurls = false))

deploydocs(
	repo = "github.com/charlesll/Spectra.jl.git",
	)
#     deps = Deps.pip("mkdocs", "mkdocs-material"),
#     repo = "github.com/charlesll/Spectra.jl.git",
# 	julia  = "1.0",
# 	osname = "linux",
# )


text:
using Test
import sbp
using sbp.Grid:grid_xp, grid_xm, grid_2d_x, grid_2d_y
using sbp.VTK:vtk_write, vtk_read

@testset "Structured Grid 2D data" begin
        nx = 20
        ny = 10
        x1, h = grid_xm(nx)
        y1, h = grid_xp(ny)
        x = grid_2d_x(x1, ny)
        y = grid_2d_y(y1, nx)
        z = x .* y

        vtk_write("test.vtk", x, y, z, nx, ny)
        xr, yr, zr = vtk_read("test.vtk")
        @test size(x, 1) == size(xr, 1)
        @test size(y, 1) == size(yr, 1)
        @test size(z, 1) == size(zr, 1)
        @test isapprox(x, xr, rtol=1e-6)
        @test isapprox(y, yr, rtol=1e-6)
        @test isapprox(z, zr, rtol=1e-6)
        rm("test.vtk")
end

max_stars_repo_name: ooreilly/sbpjl

max_issues_repo_name: ooreilly/sbpjl

max_forks_repo_name: ooreilly/sbpjl


text:
mutable struct Class_kkt_ls <: abstract_ls_info
    step_size_P::Float64
    step_size_D::Float64
    num_steps::Int64
    frac_progress::Float64
    predict_red::Float64
    actual_red::Float64
    cur_merit::Float64

    function Class_kkt_ls(iter::Class_iterate, dir::Class_point, pars::Class_parameters)
        this = new()

        this.step_size_D = 0.0
        this.step_size_P = 0.0
        this.num_steps = 0
        this.predict_red = merit_function_predicted_reduction(iter, dir, 1.0);
        this.do_ls = this.predict_red >= 0.0
        this.frac_progress = NaN
        this.actual_red = NaN
        this.cur_merit = eval_merit_function(iter, pars)

        return this
    end

    function Class_kkt_ls(step_size_P::Float64,step_size_D::Float64,num_steps::Int64,frac_progress::Float64, predict_red::Float64)
        this = new(step_size_P, step_size_D, num_steps, frac_progress, predict_red)
        return this
    end
end

function accept_func_kkt!(accept::abstract_ls_info, iter::Class_iterate, candidate::Class_iterate, dir::Class_point, step_size::Float64, filter::Array{Class_filter,1}, pars::Class_parameters, timer::class_advanced_timer)

    if satisfies_filter!(filter, candidate, accept.step_size_P, pars, timer)
    #if kkt_err_cand / kkt_err < (1.0 - pars.ls.kkt_reduction_factor * accept.step_size_P)
      return :success, accept.step_size_P
    else
      return :not_enough_progress, accept.step_size_P * pars.ls.backtracking_factor
    end
end

function accept_func!(accept::Class_kkt_ls, iter::Class_iterate, candidate::Class_iterate, dir::Class_point, step_size::Float64, filter::Array{Class_filter,1}, pars::Class_parameters, timer::class_advanced_timer)
    return accept_func_kkt!(accept, iter, candidate, dir, step_size, filter, pars, timer)
end


text:
# let Flux discover the trainable parameters

Flux.@functor RxGate
Flux.@functor RyGate
Flux.@functor RzGate
Flux.@functor RotationGate
Flux.@functor PhaseShiftGate
Flux.@functor ControlledGate
Flux.@functor EntanglementXXGate
Flux.@functor EntanglementYYGate
Flux.@functor EntanglementZZGate
Flux.@functor CircuitGate
Flux.@functor Moment
Flux.@functor MeasurementOperator
Flux.@functor Circuit

function collect_gradients(cx::Zygote.Context, q, dq)
    # special cases: circuit gate chain
    # TODO: also support measurement operators
    if typeof(q) <: AbstractVector{<:AbstractCircuitGate} || typeof(q) <: AbstractVector{<:Moment}
        for i in 1:length(q)
            collect_gradients(cx, q[i], dq[i])
        end
    else
        # need to explicitly "accumulate parameters" for Flux Params([...]) to work
        Zygote.accum_param(cx, q, dq)
        for f in fieldnames(typeof(q))
            try
                collect_gradients(cx, getfield(q, f), getfield(dq, f))
            catch UndefRefError
            end
        end
    end
end

# custom adjoint
Zygote.@adjoint apply(ψ::Vector{<:Complex}, moments::Vector{Moment}) = begin
    N = Qaintessent.intlog2(length(ψ))
    length(ψ) == 2^N || error("Vector length must be a power of 2")
    ψ1 = apply(ψ, moments)
    ψ1, function(Δ)
        # factor 1/2 due to convention for Wirtinger derivative with prefactor 1/2
        dmoments, ψbar = Qaintessent.backward(moments, ψ1, 0.5*Δ, N)
        ψbar .*= 2.0
        collect_gradients(__context__, moments, dmoments)
        return (ψbar, dmoments)
    end
end

# custom adjoint
Zygote.@adjoint apply(ψ::AbstractVector, c::Circuit{N}) where {N} = begin
    apply(ψ, c), function(Δ)
        # TODO: don't recompute apply(c, ψ) here
        dc, ψbar = Qaintessent.gradients(c, ψ, Δ)
        collect_gradients(__context__, c, dc)
        return (ψbar, dc)
    end
end

# custom adjoint for applying moments to a density matrix
Zygote.@adjoint apply(ρ::DensityMatrix, moments::Vector{Moment}) = begin
    ρ1 = apply(ρ, moments)
    ρ1, function(Δ)
        if !(Δ isa DensityMatrix)
            Δ = DensityMatrix(Δ[1], ρ.N)
        end
        dmoments, ρbar = Qaintessent.backward_density(moments, ρ1, Δ)
        collect_gradients(__context__, moments, dmoments)
        return (ρbar, dmoments)
    end
end

# custom adjoint for applying a circuit to a density matrix
Zygote.@adjoint apply(ρ::DensityMatrix, c::Circuit{N}) where {N} = begin
    @assert(ρ.N == N)
    apply(ρ, c), function(Δ)
        # TODO: don't recompute apply(c, ρ) here
        dc, ρbar = Qaintessent.gradients(c, ρ, Δ)
        collect_gradients(__context__, c, dc)
        return (ρbar, dc)
    end
end

# custom adjoint for DensityMatrix constructor; setting gradient with respect to integer `N` to zero.
Zygote.@adjoint DensityMatrix(v::AbstractVector{<:Real}, N::Integer) = DensityMatrix(v, N), ρbar -> (ρbar.v, 0)


text:
function merge_pools(pools, names...)
    (names[1] == [] || MLLabelUtils.islabelenc(collect(names), SVDD.learning_pool_enc)) || throw(ArgumentError("$(collect(names)) is not a valid label encoding."))
    return reduce((r, key) -> vcat(r, haskey(pools, key) ? pools[key] : Int64[]), unique(names); init=Int64[])
end

classify(x::Number) = x > 0 ? :outlier : :inlier

function adjust_kernel_matrix(K::Array{T, 2}; tolerance = 1e-15, warn_threshold = 1e-8) where T <: Real
    info(LOGGER, "Adjusting Kernel Matrix.")
    F = eigen(K)
    eltype(F.values) <: Complex && throw(ArgumentError("Matrix K has complex eigenvalues."))
    F.values[F.values .< tolerance] .= 0.0
    K_adjusted::Array{Float64, 2} = F.vectors * diagm(0 => F.values) * inv(F.vectors)
    K_diff = abs.(K_adjusted - K)
    sum_adjustment = sum(K_diff)
    max_adjustment = maximum(K_diff)
    info(LOGGER, "[ADJUST KERNEL] Maximum adjustemt of kernel matrix entry is $max_adjustment. The sum of adjustments is $(sum_adjustment)")
    max_adjustment .> warn_threshold && warn(LOGGER, "[ADJUST KERNEL] Maximum adjustment of kernel matrix entry exceeded the threshold of $(warn_threshold)!.")
    return K_adjusted
end


text:
# functions designed to be used inside the optimization

"""
    get_value(optparam[, x])

Returns the value(s) of the optimization parameter(s) given a scalar/vector of design
variables.  If `x` is not provided, the default parameter values are returned.
"""
get_value

get_value(optparam::OptimizationParameter) = optparam.x0
get_value(optparam::OptimizationParameter{<:Any, <:Any, <:Any, <:Any, false}, x) = optparam.x0
get_value(optparam::OptimizationParameter{<:Any, <:Number, <:Any, <:Any, true}, x) = x[1]/optparam.scaling - design_variable_offset(optparam.lb, optparam.ub)

function get_value(optparam::OptimizationParameter{S, <:AbstractArray{<:Any, N}, <:Any, <:Any, true}, x) where {S, N}

    TFX = eltype(x)

    result = MArray{S, TFX, N}(optparam.x0)

    # find all active design variables for this parameter
    idv = findall(optparam.dv)

    # remove the scaling and offset from the active design variables
    for i = 1:length(idv)
        idx = idv[i]
        scaling = optparam.scaling[idx]
        offset = design_variable_offset(optparam.lb[idx], optparam.ub[idx])
        result[idx] = x[i]/scaling - offset
    end

    return SArray(result)
end

"""
    get_values(optparam[, x])

Returns a named tuple or dictionary of the values of the optimization parameter(s)
given a vector of design variables.  If `x` is not provided, the default parameter
values are returned.
"""
get_values

# same, but for a NamedTuple of OptimizationParameter
function get_values(optparams::NamedTuple{names, types}) where {names, types}
    T = Tuple{eltype.(tuple(types.parameters...))...}
    return NamedTuple{names, T}(optparam.x0 for (name, optparam) in pairs(optparams))
end

function get_values(optparams::NamedTuple{names}, x) where names

    # get design variable indices
    indices, ndv = design_variable_indices(optparams)

    # create tuple of parameter get_values
    get_values = (isempty(indices[name]) ? optparam.x0 : get_value(optparam, x[indices[name][1:length(indices[name])]]) for (name,optparam) in pairs(optparams))

    # return a named tuple of the parameters
    return NamedTuple{names}(get_values)
end

# same, but for a dictionary of OptimizationParameter

function get_values(optparams::AbstractDict)
    return OrderedDict([(name, optparam.x0) for (name, optparam) in pairs(optparams)])
end

function get_values(optparams::AbstractDict, x) where names

    # get design variable indices
    indices, ndv = design_variable_indices(optparams)

    # return dictionary of parameter get_values
    return OrderedDict(
        [isempty(indices[name]) ?
        (name, optparam.x0) :
        (name, get_value(optparam, x[indices[name][1:length(indices[name])]]))
        for (name,optparam) in pairs(optparams)])

end

"""
    design_variable_indices(optparams)

Returns a dictionary which contains the indices of the design variables for
each optimization parameter. Also returns the number of design variables.
"""
# gets the indices of the design variables for each optimization parameter
function design_variable_indices!(indices, optparams)

    ndv = 0
    for (name, optparam) in pairs(optparams)
        idx1 = ndv + 1
        idx2 = ndv + count(optparam.dv)::Int
        indices[name] = idx1:idx2
        ndv = idx2
    end

    return indices, ndv
end

function design_variable_indices(optparams::NamedTuple)
    indices = OrderedDict{Symbol,UnitRange{Int}}()
    return design_variable_indices!(indices, optparams)
end

function design_variable_indices(optparams::AbstractDict)
    indices = OrderedDict{keytype(optparams),UnitRange{Int}}()
    return design_variable_indices!(indices, optparams)
end

"""
    print_design_variables(optparams, x)

Prints design variable get_values for the current iteration
"""
function print_design_variables(optparams, x)

    # get design variable indices
    indices, ndv = design_variable_indices(optparams)

    # print each design variable get_value
    for (name, optparam) in pairs(optparams)
        if !isempty(indices[name])
            if length(indices[name]) == length(optparam.x0)
                println(name, ": ", get_value(optparam, x[indices[name]]))
            else
                dv = findall(optparam.dv)
                println(name, dv,": ", get_value(optparam, x[indices[name]])[dv])
            end
        end
    end

    return nothing
end

"""
    design_variable_offset(lb, ub)

Determines the amount to offset a design variable in order to center it around zero
"""
design_variable_offset(lb, ub) = -(min(ub, 1e20) + max(lb, -1e20))/2


text:
using COSMO, Test, SparseArrays, Random, LinearAlgebra

# This test is precision agnostic and will be run with TestFloat precision
if !@isdefined(UnitTestFloats)
    UnitTestFloats = [Float64] #if not run in full test setup, just do it for one float type
end
@testset "Model" begin

for TestFloat in UnitTestFloats

    rng = Random.MersenneTwister(1872381)

    P_int = 4
    q_int = 2

    P_Float = TestFloat(4.)
    q_Float = TestFloat(2.)

    P_UInt = UInt64(4)
    q_UInt = UInt64(2)

    P_mat = rand(rng, TestFloat, 10, 10)
    q_mat = rand(rng, TestFloat, 10, 1)

    P_sparse = sprand(rng, TestFloat, 10, 10, 0.4)
    q_sparse = rand(rng, TestFloat, 10, 1)

    A_vec = TestFloat[1 2 3 4]
    b_vec = TestFloat(1)

    scalar_c = COSMO.Constraint(TestFloat(1.), TestFloat(2.), COSMO.ZeroSet)
    constr = COSMO.Constraint(diagm(0 => ones(TestFloat, 10)), rand(rng, TestFloat, 10), COSMO.ZeroSet)

    @testset "Model (T = $(TestFloat))" begin
        # default case
        model = COSMO.Model()
        @test typeof(model) == COSMO.Model{Float64}

        # If integers are used, we convert to model parameter type
        model = COSMO.Model{TestFloat}()
        @test assemble!(model, P_int, q_int, [scalar_c]) == nothing

        model = COSMO.Model{TestFloat}()
        @test assemble!(model, P_Float, q_Float, [scalar_c]) == nothing

        model = COSMO.Model{TestFloat}()
        @test assemble!(model, P_Float, q_Float, [scalar_c]) == nothing

        model = COSMO.Model{TestFloat}()
        @test assemble!(model, P_mat, q_mat, [constr]) == nothing

        model = COSMO.Model{TestFloat}()
        @test assemble!(model, P_sparse, q_sparse, [constr]) == nothing

        # Prevent mixing of precisions
        if TestFloat == Float32
            model = COSMO.Model()
            @test_throws ArgumentError assemble!(model, P_sparse, q_sparse, [constr])

            # test that settings and model have consistent types
            model = COSMO.Model{Float64}()
            settings = COSMO.Settings{Float32}()
            @test_throws ArgumentError assemble!(model, rand(Float64, 10, 10), rand(Float64, 10), [constr], settings = settings)
        end

    end #testset
end #UnitTestFloats
end
nothing


text:
# ---
# jupyter:
#   jupytext:
#     cell_metadata_filter: -all
#     formats: jl,ipynb
#     text_representation:
#       extension: .jl
#       format_name: light
#       format_version: '1.5'
#       jupytext_version: 1.3.2
#   kernelspec:
#     display_name: Julia 1.1.0
#     language: julia
#     name: julia-1.1
# ---

# +
using Plots
pyplot();

import Chemostat
const Ch = Chemostat;
# -

function ecToyModel(;Inf_ = 1000)
    # Kinetic params
    k13 = rand()
    k14 = rand()
    k24 = rand()
    k5f = rand()
    k5b = rand()

    S = [# v1     v2     v3    v4/arm  v41    v42    v5f    v5b    v6         e1     e2     e3     e4

           +1      0     -1      0      0      0      0      0      0         0      0      0      0 ;  # M1
            0     +1     -1      0      0      0     -1     +1      0         0      0      0      0 ;  # M2
            0      0     +1     -1      0      0      0      0      0         0      0      0      0 ;  # M3
            0      0      0      0     +1     +1     +1     -1     -1         0      0      0      0 ;  # M4
            0      0      0     +1     -1     -1      0      0      0         0      0      0      0 ;  # PM1

            0      0   -1/k13    0   -1/k14    0      0      0      0        +1      0      0      0 ;  # E1
            0      0      0      0      0   -1/k24    0      0      0         0     +1      0      0 ;  # E2
            0      0      0      0      0      0   -1/k5f -1/k5b    0         0      0     +1      0 ;  # E3
            0      0      0      0      0      0   -2/k5f -2/k5b    0         0      0      0     +1 ;  # E4
        ]

    mets = ["M1", "M2", "M3", "M4", "PM1", "E1", "E2", "E3", "E4"]
    rxns = ["v1", "v2", "v3", "v4/arm", "v41", "v42", "v5f", "v5b", "v6",  "e1" ,  "e2" ,  "e3" ,  "e4" ]
    ub =   [Inf_, Inf_, Inf_,   Inf_  , Inf_ , Inf_ , Inf_ , Inf_ , Inf_, rand(), rand(), rand(), rand()]
    lb = zeros(size(S, 2));
    b = zeros(size(S, 1));
    model = Ch.Utils.MetNet(S, b, lb, ub, rxns, mets);
    return Ch.LP.fva_preprocess(model, verbose = false);
end

# <img src="./ecToyModel.png" width="400" height="700">

model = ecToyModel();

println("Balance equations")
foreach(model.mets) do met
    Ch.Utils.balance_str(model, met) |> println
end

println("Reaction equations")
foreach(model.rxns) do rxn
    eq = Ch.Utils.rxn_str(model, rxn)
    b = Ch.Utils.bounds(model, rxn)
    println(rpad(rxn, 10), " eq: ", eq, "   ", collect(b))
end

# ## FBA

obj_ider = "v6"
fbaout = Ch.LP.fba(model, obj_ider);

let p1 = plot(title = "FBA", xlabel = "flx_id", ylabel = "flux_val"),
    p2 = plot(title = "MODEL", xlabel = "flx_id", ylabel = "UB", yscale = :log10)
    plot([bar!(p1, model.rxns, fbaout.v, label = "", color = :red),
          bar!(p2, model.rxns, model.ub, label = "", yscale = :log10, color = :blue)]...,
    size = [900, 400])
end

# ## EP

epout = Ch.MaxEntEP.maxent_ep(model, verbose = false);

# ## MC-HR

hrout = Ch.MaxEntHR.maxent_hr(model; nsamples = 500_000, drop_samples = true);

let ps = map(model.rxns) do rxn
        p = plot(legend = false, title = rxn, xlim = [0.0, 1.0])
        Ch.Plots.plot_marginal!(p, model, [hrout, epout, fbaout], rxn)
    end
    ncols = 4;
    while mod(length(ps), ncols) != 0
        empty_p = plot([],[],title = "", grid = false, label = "", axis = false)
        push!(ps, empty_p)
    end
    n = length(ps)
    nrows = n / ncols |> Int
    plot(ps...; size = [ncols * 400, nrows * 400],
        layout = grid(nrows, ncols))
end

max_stars_repo_name: josePereiro/Chemostat

max_issues_repo_name: josePereiro/Chemostat

max_forks_repo_name: josePereiro/Chemostat


text:
module MiniProgressBars

export MiniProgressBar, start_progress, end_progress, show_progress, print_progress_bottom

using Printf

Base.@kwdef mutable struct MiniProgressBar
    max::Int = 1.0
    header::String = ""
    color::Symbol = :white
    width::Int = 40
    current::Int = 0.0
    prev::Int = 0.0
    has_shown::Bool = false
    time_shown::Float64 = 0.0
    percentage::Bool = true
    always_reprint::Bool = false
    indent::Int = 4
end

const NONINTERACTIVE_TIME_GRANULARITY = Ref(2.0)
const PROGRESS_BAR_PERCENTAGE_GRANULARITY = Ref(0.1)

function pretend_cursor()
    io = stderr
    bar = MiniProgressBar(; indent=2, header = "Progress", color = Base.info_color(),
    percentage=false, always_reprint=true)
    bar.max = 40
    start_progress(io, bar)
    sleep(0.5)
    for i in 1:40
        bar.current = i
        print_progress_bottom(io)
        println("Downloading ... $i")
        show_progress(io, bar)
        x = randstring(7)
        sleep(0.01)
    end
    end_progress(io, bar)
    print("Hello!")
end

function start_progress(io::IO, _::MiniProgressBar)
    ansi_disablecursor = "\e[?25l"
    print(io, ansi_disablecursor)
end

function show_progress(io::IO, p::MiniProgressBar)
    if p.max == 0
        perc = 0.0
        prev_perc = 0.0
    else
        perc = p.current / p.max * 100
        prev_perc = p.prev / p.max * 100
    end
    # Bail early if we are not updating the progress bar,
    # Saves printing to the terminal
    if !p.always_reprint && p.has_shown && !((perc - prev_perc) > PROGRESS_BAR_PERCENTAGE_GRANULARITY[])
        return
    end
    if !isinteractive()
        t = time()
        if p.has_shown && (t - p.time_shown) < NONINTERACTIVE_TIME_GRANULARITY[]
            return
        end
        p.time_shown = t
    end
    p.prev = p.current
    p.has_shown = true
    n_filled = ceil(Int, p.width * perc / 100)
    n_left = p.width - n_filled
    print(io, " "^p.indent)
    printstyled(io, p.header, color=p.color, bold=true)
    print(io, " [")
    print(io, "="^n_filled, ">")
    print(io, " "^n_left, "]  ", )
    if p.percentage
        @printf io "%2.1f %%" perc
    else
        print(io, p.current, "/",  p.max)
    end
    print(io, "\r")
end

function end_progress(io, p::MiniProgressBar)
    ansi_enablecursor = "\e[?25h"
    ansi_clearline = "\e[2K"
    print(io, ansi_enablecursor)
    print(io, ansi_clearline)
end

# Useful when writing a progress bar in the bottom
# makes the bottom progress bar not flicker
# prog = MiniProgressBar(...)
# prog.end = n
# for progress in 1:n
#     print_progree_bottom(io)
#     println("stuff")
#     prog.current = progress
#     showproress(io, prog)
#  end
#
function print_progress_bottom(io::IO)
    print(io, "\e[S")
    ansi_clearline = "\e[2K"
    ansi_movecol1 = "\e[1G"
    ansi_moveup(n::Int) = string("\e[", n, "A")
    print(io, ansi_moveup(1))
    print(io, ansi_clearline)
    print(io, ansi_movecol1)
end

end


text:
abstract type AbstractDimensionalArray{T,N,D<:Tuple} <: AbstractArray{T,N} end

const AbDimArray = AbstractDimensionalArray

const StandardIndices = Union{AbstractArray,Colon,Integer}

# Interface methods ############################################################

dims(A::AbDimArray) = A.dims
label(A::AbDimArray) = ""

# Array interface methods ######################################################

Base.size(A::AbDimArray) = size(parent(A))
Base.iterate(A::AbDimArray, args...) = iterate(parent(A), args...)
Base.show(io::IO, A::AbDimArray) = begin
    printstyled(io, "\n", label(A), ": "; color=:red)
    show(io, typeof(A))
    show(io, parent(A))
    printstyled(io, "\n\ndims:\n"; color=:magenta)
    show(io, dims(A))
    show(io, refdims(A))
    printstyled(io, "\n\nmetadata:\n"; color=:cyan)
end

Base.@propagate_inbounds Base.getindex(A::AbDimArray, I::Vararg{<:Integer}) =
    getindex(parent(A), I...)
Base.@propagate_inbounds Base.getindex(A::AbDimArray, I::Vararg{<:StandardIndices}) =
    rebuildsliced(A, getindex(parent(A), I...), I)

Base.@propagate_inbounds Base.view(A::AbDimArray, I::Vararg{<:StandardIndices}) =
    rebuildsliced(A, view(parent(A), I...), I)

Base.convert(::Type{Array{T,N}}, A::AbDimArray{T,N}) where {T,N} =
    convert(Array{T,N}, parent(A))

Base.copy(A::AbDimArray) = rebuild(A, copy(parent(A)))
Base.copy!(dst::AbDimArray, src::AbDimArray) = copy!(parent(src), parent(dst))
Base.copy!(dst::AbDimArray, src::AbstractArray) = copy!(parent(src), dst)

Base.BroadcastStyle(::Type{<:AbDimArray}) = Broadcast.ArrayStyle{AbDimArray}()

Base.similar(A::AbDimArray) = rebuild(A, similar(parent(A)))
Base.similar(A::AbDimArray, ::Type{T}) where T = rebuild(A, similar(parent(A), T))
Base.similar(bc::Broadcast.Broadcasted{Broadcast.ArrayStyle{AbDimArray}}, ::Type{ElType}) where ElType = begin
    A = find_dimensional(bc)
    # TODO How do we know what the new dims are?
    rebuildsliced(A, similar(Array{ElType}, axes(bc)), axes(bc))
end

# Need to cover a few type signatures to avoid ambiguity with base
# Don't remove these even though they look redundant
Base.similar(A::AbDimArray, ::Type{T}, I::Vararg{<:Integer}) where T =
    rebuildsliced(A, similar(parent(A), T, I...), I)
Base.similar(A::AbDimArray, ::Type{T}, ::Tuple{Int,Vararg{Int}}) where T =
    rebuild(A, similar(parent(A), T))
Base.similar(A::AbDimArray, ::Type{T}, I::Tuple{Union{Integer,OneTo},Vararg{Union{Integer,OneTo},N}}) where {T,N} =
    rebuildsliced(A, similar(parent(A), T, I...), I)

@inline find_dimensional(bc::Base.Broadcast.Broadcasted) = find_dimensional(bc.args)
@inline find_dimensional(ext::Base.Broadcast.Extruded) = find_dimensional(ext.x)
@inline find_dimensional(args::Tuple{}) = error("dimensional array not found")
@inline find_dimensional(args::Tuple) = find_dimensional(find_dimensional(args[1]), tail(args))
@inline find_dimensional(x) = x
@inline find_dimensional(A::AbDimArray, rest) = A
@inline find_dimensional(::Any, rest) = find_dimensional(rest)

# Concrete implementation ######################################################

"""
A basic DimensionalArray type

Maintains and updates its dimensions through transformations
"""
struct DimensionalArray{T,N,D<:Tuple,R<:Tuple,A<:AbstractArray{T,N}} <: AbstractDimensionalArray{T,N,D}
    data::A
    dims::D
    refdims::R
end
DimensionalArray(A::AbstractArray, dims; refdims=()) =
    DimensionalArray(A, formatdims(A, dims), refdims)

# Getters
refdims(A::DimensionalArray) = A.refdims

# DimensionalArray interface
@inline rebuild(A::DimensionalArray, data, dims, refdims) =
    DimensionalArray(data, dims, refdims)

# Array interface (AbstractDimensionalArray takes care of everything else)
Base.parent(A::DimensionalArray) = A.data

Base.@propagate_inbounds Base.setindex!(A::DimensionalArray, x, I::Vararg{StandardIndices}) =
    setindex!(parent(A), x, I...)


text:
module Benchmark

using Profile
using BSONqs

const BSON = BSONqs

export do_bench, do_profile, profile_load, Foo, Bar, Baz

chars = [x for x in '0':'z']
strings = [String([rand(chars) for _ in 1:20]) for _ in 1:1000]
rstr(n::Int)::String = rand(strings)[1:n]

struct Baz
  going::String
  deeper::String
end

Baz() = Baz(rstr(20), rstr(1))

struct Bar
  level::Int64
  bazes::Vector{Baz}
  salty::AbstractDict{<:AbstractString, <:Unsigned}
end

Bar() = Bar(rand(Int64), [Baz() for _ in 1:50],
            Dict(s => hash(s, UInt64(0xdeadbeef)) for s in (rstr(x) for x in 10:13)))

struct Foo
  agile::String
  software::String
  management::Union{String, Symbol}
  consultant::String
  training::Union{String, Missing}
  projects::Vector{Bar}
end

Foo() = [Foo(rstr(5), rstr(7), rstr(17), rstr(11), rstr(13),
             [Bar() for _ in 1:2000]) for _ in 1:3]

struct Result
  elapsed::Float64
  allocated::Int64
end

const history_file = "./benchmark-history.bson"

macro bench(hist, msg, ex)
  sex = "$ex"
  quote
    GC.gc()
    @info $msg
    local val, t1, bytes, gctime, memallocs = @timed $(esc(ex))
    local mb = ceil(bytes / (1024 * 1024))
    if $sex in keys($(esc(hist)))
      local t0 = $(esc(hist))[$sex].elapsed
      @info $sex elapsed=t1 speedup=t0/t1 allocatedMb=mb gctime
    else
      @info $sex elapsed=t1 allocatedMb=mb gctime
    end
    $(esc(hist))[$sex] = Result(t1, bytes)
    val
  end
end

function do_bench()
  hist = if isfile(history_file)
    BSON.load(history_file)::Dict{String, Result}
  else
    Dict{String, Result}()
  end

  foos = Dict(:Foo => Foo(), :Foo2 => Foo(), :Foo3 => Foo())

  #@bench hist "Roundtrip from cold start (ignore)" BSON.roundtrip(foos)

  io = IOBuffer()

  @bench hist "Bench Save BSON" bson(io, foos)
  seek(io, 0)

  doc = @bench hist "Bench inter load BSON Document" BSON.load_compat(io)
  seek(io, 0)
  GC.gc()
  rfoos = @bench hist "Bench direct load BSON Document" BSON.load(io)

  # Sanity check the results
  rfoos[:Foo][1]::Foo
  rfoos[:Foo][1].projects[1]::Bar
  rfoos[:Foo][1].projects[1].bazes[1]::Baz

  bson(history_file, hist)
end

function do_profile()
  io = IOBuffer()
  foos = Dict(:Foo => Foo(), :Foo2 => Foo(), :Foo3 => Foo())
  minc = parse(Int64, get(ENV, "JULIA_PROFILE_MIN", "0"))
  seek(io, 0)
  GC.gc()
  @info "Profile Save BSON"
  @profile bson(io, foos)
  Profile.print(;noisefloor=2, mincount=minc)
  Profile.clear()
  seek(io, 0)

  # GC.gc()
  # @info "Profile Parse BSON Document"
  # dict = @profile BSON.parse_doc(io)
  # Profile.print(;noisefloor=2, mincount=minc, C=true)
  # Profile.clear()

  # GC.gc()
  # @info "Profile deref"
  # dict = @profile BSON.backrefs!(dict)
  # Profile.print(;noisefloor=2, mincount=minc, C=true)
  # Profile.clear()

  # GC.gc()
  # @info "Profile Raise BSON to Julia types"
  # rfoos = @profile BSON.raise_recursive(dict, IdDict{Any, Any}())
  # Profile.print(;noisefloor=2, mincount=minc, C=true)
  # Profile.clear()

  GC.gc()
  seek(io, 0)
  @info "Profile direct parse"
  dict = @profile BSON.load(io)
  Profile.print(;noisefloor=2, mincount=minc, C=true)
end

function profile_load(obj)
  io = IOBuffer()
  bson(io, Dict(:stuff => obj))

  seek(io, 0)
  BSON.load(io)
  Profile.clear()
  @profile for _ in 1:100
    seek(io, 0)
    BSON.load(io)
  end
end

end


text:
# Use baremodule to shave off a few KB from the serialized `.ji` file
baremodule Arpack_jll
using Base
using Base: UUID
import JLLWrappers

JLLWrappers.@generate_main_file_header("Arpack")
JLLWrappers.@generate_main_file("Arpack", UUID("68821587-b530-5797-8361-c406ea357684"))
end  # module Arpack_jll


text:
using Test
import ACE, JuLIP
import InteractiveUtils

using BenchmarkTools: @belapsed
using LinearAlgebra: eigvals, eigen
using ProgressMeter: @showprogress

import JuLIP.Potentials: F64fun
import JuLIP: Atoms, bulk, rattle!, positions, energy, forces, JVec,
              chemical_symbol, mat, rnn,
              read_dict, write_dict
import JuLIP.MLIPs: combine
import JuLIP.Testing: print_tf

ACE.@aceimports

#---

testfile = @__DIR__() * "/regtest.json"

if isfile(testfile)
   error("regressions file `regtest.json` already exists. Aborting...")
end

# generate the set of parameters
# degrees = Dict(2 => [6, 7, 8],
#                3 => [7, 9, 11],
#                4 => [8, 10, 12])
degrees = Dict(2 => [8, 13, 18],
               4 => [8, 13, 18],
               6 => [10, 16, 18]),
               8 => [10, 13, 16])
species = ["Si"]
r0 = rnn(:Si)
rcuts = [5.5, ]
params = []
for (N, degs) in degrees, maxdeg in degs, rcut in rcuts
   push!( params, Dict( "species" => species, "N" => N,
                        "maxdeg" => maxdeg, "r0" => r0,
                        "rcut" => rcut, "wL" => 1.3 ) )
end
D = Dict("paramsets" => params,
               "at" => write_dict(ACE.Random.rand_config(:Si; repeat=3)))

# generate and save the basis sets if they don't  exist yet
@info("Generate basis and coefficients...")
D["results"] = []
D["labels"] = []
D["versioninfo"] = Dict{String, Any}()
@showprogress for params in D["paramsets"]
   basis = ACE.Testing.test_basis(params)
   coeffs = ACE.Random.randcoeffs(basis; diff=2)
   push!(D["results"], Dict("basis" => write_dict(basis), "coeffs" => coeffs))
end
JuLIP.save_dict(testfile, D)


text:
# This file was generated by the Julia Swagger Code Generator
# Do not modify this file directly. Modify the swagger specification instead.

@doc raw"""AzureDisk represents an Azure Data Disk mount on the host and bind mount to the pod.

    IoK8sApiCoreV1AzureDiskVolumeSource(;
        cachingMode=nothing,
        diskName=nothing,
        diskURI=nothing,
        fsType=nothing,
        kind=nothing,
        readOnly=nothing,
    )

    - cachingMode::String : Host Caching mode: None, Read Only, Read Write.
    - diskName::String : The Name of the data disk in the blob storage
    - diskURI::String : The URI the data disk in the blob storage
    - fsType::String : Filesystem type to mount. Must be a filesystem type supported by the host operating system. Ex. \"ext4\", \"xfs\", \"ntfs\". Implicitly inferred to be \"ext4\" if unspecified.
    - kind::String : Expected values Shared: multiple blob disks per storage account  Dedicated: single blob disk per storage account  Managed: azure managed data disk (only in managed availability set). defaults to shared
    - readOnly::Bool : Defaults to false (read/write). ReadOnly here will force the ReadOnly setting in VolumeMounts.
"""
mutable struct IoK8sApiCoreV1AzureDiskVolumeSource <: SwaggerModel
    cachingMode::Any # spec type: Union{ Nothing, String } # spec name: cachingMode
    diskName::Any # spec type: Union{ Nothing, String } # spec name: diskName
    diskURI::Any # spec type: Union{ Nothing, String } # spec name: diskURI
    fsType::Any # spec type: Union{ Nothing, String } # spec name: fsType
    kind::Any # spec type: Union{ Nothing, String } # spec name: kind
    readOnly::Any # spec type: Union{ Nothing, Bool } # spec name: readOnly

    function IoK8sApiCoreV1AzureDiskVolumeSource(;cachingMode=nothing, diskName=nothing, diskURI=nothing, fsType=nothing, kind=nothing, readOnly=nothing)
        o = new()
        validate_property(IoK8sApiCoreV1AzureDiskVolumeSource, Symbol("cachingMode"), cachingMode)
        setfield!(o, Symbol("cachingMode"), cachingMode)
        validate_property(IoK8sApiCoreV1AzureDiskVolumeSource, Symbol("diskName"), diskName)
        setfield!(o, Symbol("diskName"), diskName)
        validate_property(IoK8sApiCoreV1AzureDiskVolumeSource, Symbol("diskURI"), diskURI)
        setfield!(o, Symbol("diskURI"), diskURI)
        validate_property(IoK8sApiCoreV1AzureDiskVolumeSource, Symbol("fsType"), fsType)
        setfield!(o, Symbol("fsType"), fsType)
        validate_property(IoK8sApiCoreV1AzureDiskVolumeSource, Symbol("kind"), kind)
        setfield!(o, Symbol("kind"), kind)
        validate_property(IoK8sApiCoreV1AzureDiskVolumeSource, Symbol("readOnly"), readOnly)
        setfield!(o, Symbol("readOnly"), readOnly)
        o
    end
end # type IoK8sApiCoreV1AzureDiskVolumeSource

const _property_map_IoK8sApiCoreV1AzureDiskVolumeSource = Dict{Symbol,Symbol}(Symbol("cachingMode")=>Symbol("cachingMode"), Symbol("diskName")=>Symbol("diskName"), Symbol("diskURI")=>Symbol("diskURI"), Symbol("fsType")=>Symbol("fsType"), Symbol("kind")=>Symbol("kind"), Symbol("readOnly")=>Symbol("readOnly"))
const _property_types_IoK8sApiCoreV1AzureDiskVolumeSource = Dict{Symbol,String}(Symbol("cachingMode")=>"String", Symbol("diskName")=>"String", Symbol("diskURI")=>"String", Symbol("fsType")=>"String", Symbol("kind")=>"String", Symbol("readOnly")=>"Bool")
Base.propertynames(::Type{ IoK8sApiCoreV1AzureDiskVolumeSource }) = collect(keys(_property_map_IoK8sApiCoreV1AzureDiskVolumeSource))
Swagger.property_type(::Type{ IoK8sApiCoreV1AzureDiskVolumeSource }, name::Symbol) = Union{Nothing,eval(Base.Meta.parse(_property_types_IoK8sApiCoreV1AzureDiskVolumeSource[name]))}
Swagger.field_name(::Type{ IoK8sApiCoreV1AzureDiskVolumeSource }, property_name::Symbol) =  _property_map_IoK8sApiCoreV1AzureDiskVolumeSource[property_name]

function check_required(o::IoK8sApiCoreV1AzureDiskVolumeSource)
    (getproperty(o, Symbol("diskName")) === nothing) && (return false)
    (getproperty(o, Symbol("diskURI")) === nothing) && (return false)
    true
end

function validate_property(::Type{ IoK8sApiCoreV1AzureDiskVolumeSource }, name::Symbol, val)
end


text:
import Statistics

line = read("input.txt", String)
numbers = parse.(Int, split(line, ','))

# argmin x ∑i=1-N |s_i−x| = ∑i=1-N sign(s_i − x) -> hence, equal number of positive and negative elements = median
position = convert(Int, round(Statistics.median(numbers)))
println(sum(@. abs(numbers - position)))


text:
using .Models

measures = [accuracy, confmat, misclassification_rate]
@test MLJTuning.measure_adjusted_weights([2, 3, 4], measures) == [-2, 0, 4]
@test MLJTuning.measure_adjusted_weights(nothing, measures) == [-1, 0, 0]
@test_throws(DimensionMismatch,
             MLJTuning.measure_adjusted_weights([2, 3], measures))

@testset "losses/scores get minimized/maximimized" begin
    bad_model = KNNClassifier(K=100)
    good_model = KNNClassifier(K=5)

    am = [accuracy, misclassification_rate]
    ma = [misclassification_rate, accuracy]

    # scores when `weights=nothing`
    history = [(model=bad_model, measure=am, measurement=[0, 1]),
               (model=good_model, measure=am, measurement=[1, 0])]
    @test MLJTuning.best(NaiveSelection(), history).model == good_model

    # losses when `weights=nothing`
    history = [(model=bad_model, measure=ma, measurement=[1, 0]),
               (model=good_model, measure=ma, measurement=[0, 1])]
    @test MLJTuning.best(NaiveSelection(), history).model == good_model

    # mixed case favouring the score:
    weights = [2, 1]
    history = [(model=bad_model, measure=am, measurement=[0, 0]),
               (model=good_model, measure=am, measurement=[1, 1])]
    heuristic = NaiveSelection(weights=weights)
    @test MLJTuning.best(heuristic, history).model == good_model

    # mixed case favouring the loss:
    weights = [1, 2]
    history = [(model=bad_model, measure=am, measurement=[1, 1]),
               (model=good_model, measure=am, measurement=[0, 0])]
    heuristic = NaiveSelection(weights=weights)
    @test MLJTuning.best(heuristic, history).model == good_model
end


text:
struct StaticState <: AbstractArray{StaticTeam, 1}
    team1::StaticTeam
    team2::StaticTeam
end

Base.size(s::StaticState) = (2,)
Base.IndexStyle(::Type{<:StaticState}) = IndexLinear()
Base.getindex(s::StaticState, i::UInt8) = i == 0x01 ? s.team1 : s.team2

StaticState(teams::NTuple{6, Union{String, Int}}; league = "great",
    cup = "all") = StaticState(
        StaticTeam((teams[1], teams[2], teams[3]), league = league, cup = cup),
        StaticTeam((teams[4], teams[5], teams[6]), league = league, cup = cup)
    )

struct DynamicState <: AbstractArray{DynamicTeam, 1}
    team1::DynamicTeam
    team2::DynamicTeam
    data::UInt32 # active, pending fast moves, cmp, chance, fast move damage
end

Base.size(d::DynamicState) = (2,)
Base.IndexStyle(::Type{<:DynamicState}) = IndexLinear()
Base.getindex(d::DynamicState, i::UInt8) = i == 0x01 ? d.team1 : d.team2

get_active(state::DynamicState) = UInt8(state.data & UInt32(3)),
    UInt8((state.data >> UInt32(2)) & UInt32(3))
get_fast_moves_pending(state::DynamicState) =
    UInt8((state.data >> UInt32(4)) % UInt32(7)),
    UInt8((state.data ÷ UInt32(112)) % UInt32(7))
get_cmp(state::DynamicState) = UInt8((state.data ÷ UInt32(784)) % UInt32(5))
get_chance(state::DynamicState) = UInt8((state.data ÷ UInt32(3920)) % UInt32(6))
get_fm_damage(state::DynamicState) =
    UInt16((state.data ÷ UInt32(23520)) % UInt32(425)),
    UInt16(state.data ÷ UInt32(9996000))

function DynamicState(s::StaticState)
    fm_dmg_1 = UInt32(calculate_damage(s[0x01][0x01].stats.attack,
            s[0x02][0x01].stats.defense,
            (s[0x02][0x01].primary_type, s[0x02][0x01].secondary_type),
            0x04, 0x04, s[0x01][0x01].fast_move))
    fm_dmg_2 = UInt32(calculate_damage(s[0x02][0x01].stats.attack,
        s[0x01][0x01].stats.defense,
        (s[0x01][0x01].primary_type, s[0x01][0x01].secondary_type),
        0x04, 0x04, s[0x02][0x01].fast_move))
    return DynamicState(DynamicTeam(s[0x01]), DynamicTeam(s[0x02]), UInt32(133) +
        fm_dmg_1 * UInt32(23520) + fm_dmg_2 * UInt32(9996000))
end

DynamicState(team1::DynamicTeam, team2::DynamicTeam, active_1::UInt8,
    active_2::UInt8, fm_pending_1::UInt8, fm_pending_2::UInt8, cmp::UInt8,
    chance::UInt8, fm_dmg_1::UInt16, fm_dmg_2::UInt16) =
    DynamicState(team1, team2, UInt32(active_1) +
        UInt32(active_2) << UInt32(2) + UInt32(fm_pending_1) << UInt32(4) +
        UInt32(112) * UInt32(fm_pending_2) + UInt32(cmp) * UInt32(784) +
        UInt32(chance) * UInt32(3920) + UInt32(fm_dmg_1) * UInt32(23520) +
        UInt32(fm_dmg_2) * UInt32(9996000))


text:
# MIT License
#
# Copyright (c) 2018 Martin Biel
#
# Permission is hereby granted, free of charge, to any person obtaining a copy
# of this software and associated documentation files (the "Software"), to deal
# in the Software without restriction, including without limitation the rights
# to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
# copies of the Software, and to permit persons to whom the Software is
# furnished to do so, subject to the following conditions:
#
# The above copyright notice and this permission notice shall be included in all
# copies or substantial portions of the Software.
#
# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
# IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
# FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
# AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
# LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
# OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
# SOFTWARE.

SubWorker{T,F,I} = RemoteChannel{Channel{Vector{SubProblem{T,F,I}}}}
ScenarioProblemChannel{S} = RemoteChannel{Channel{ScenarioProblems{S}}}
Work = RemoteChannel{Channel{Int}}

function load_subproblems!(subworkers::Vector{SubWorker{T,F,I}},
                           scenarioproblems::DistributedScenarioProblems,
                           decisions::Vector{DecisionChannel},
                           feasibility_strategy::AbstractFeasibilityStrategy,
                           integer_strategy::AbstractIntegerStrategy) where {T <: AbstractFloat,
                                                                             F <: AbstractFeasibilityAlgorithm,
                                                                             I <: AbstractIntegerAlgorithm}
    # Create subproblems on worker processes
    @sync begin
        for w in workers()
            subworkers[w-1] = RemoteChannel(() -> Channel{Vector{SubProblem{T,F,I}}}(1), w)
            prev = map(2:(w-1)) do p
                scenarioproblems.scenario_distribution[p-1]
            end
            start_id = isempty(prev) ? 0 : sum(prev)
            @async remotecall_fetch(initialize_subworker!,
                                    w,
                                    subworkers[w-1],
                                    scenarioproblems[w-1],
                                    feasibility_strategy,
                                    integer_strategy,
                                    start_id)
        end
    end
end

function initialize_subworker!(subworker::SubWorker{T,F,I},
                               scenarioproblems::ScenarioProblemChannel,
                               feasibility_strategy::AbstractFeasibilityStrategy,
                               integer_strategy::AbstractIntegerStrategy,
                               start_id::Integer) where {T <: AbstractFloat,
                                                         F <: AbstractFeasibilityAlgorithm,
                                                         I <: AbstractIntegerAlgorithm}
    sp = fetch(scenarioproblems)
    subproblems = Vector{SubProblem{T,F,I}}(undef, num_subproblems(sp))
    for i in 1:num_subproblems(sp)
        subproblems[i] = SubProblem(
            subproblem(sp, i),
            start_id + i,
            T(probability(scenario(sp, i))),
            feasibility_strategy,
            integer_strategy)
    end
    put!(subworker, subproblems)
    return nothing
end

function mutate_subproblems!(mutator::Function, subworkers::Vector{<:SubWorker})
    @sync begin
        for w in workers()
            @async remotecall_fetch(w, mutator, subworkers[w-1]) do mutator, sw
                for subproblem in fetch(sw)
                    mutator(subproblem)
                end
            end
        end
    end
    return nothing
end

function restore_subproblems!(subworkers::Vector{<:SubWorker})
     @sync begin
        for w in workers()
            @async remotecall_fetch(w, subworkers[w-1]) do sw
                for subproblem in fetch(sw)
                    restore_subproblem!(subproblem)
                end
            end
        end
    end
    return nothing
end

function resolve_subproblems!(subworker::SubWorker{T,F,I},
                              decisions::DecisionChannel,
                              x::AbstractVector,
                              cutqueue::CutQueue{T},
                              aggregator::AbstractAggregator,
                              t::Integer,
                              metadata::MetaDataChannel,
                              worker_metadata::MetaDataChannel) where {T <: AbstractFloat,
                                                                F <: AbstractFeasibilityAlgorithm,
                                                                I <: AbstractIntegerAlgorithm}
    # Fetch all subproblems stored in worker
    subproblems::Vector{SubProblem{T,F,I}} = fetch(subworker)
    if isempty(subproblems)
        # Workers has nothing do to, return.
        return nothing
    end
    # Update subproblems
    update_known_decisions!(fetch(decisions), x)
    # Aggregation policy
    aggregation::AbstractAggregation = aggregator(length(subproblems), T)
    # Solve subproblems
    for subproblem in subproblems
        update_subproblem!(subproblem)
        cut::SparseHyperPlane{T} = subproblem(x, metadata)
        aggregate_cut!(cutqueue, aggregation, worker_metadata, t, cut, x)
    end
    flush!(cutqueue, aggregation, worker_metadata, t, x)
    return nothing
end

function work_on_subproblems!(subworker::SubWorker{T,F,I},
                              decisions::DecisionChannel,
                              work::Work,
                              finalize::Work,
                              cutqueue::CutQueue{T},
                              iterates::RemoteIterates{A},
                              metadata::MetaDataChannel,
                              worker_metadata::MetaDataChannel,
                              aggregator::AbstractAggregator) where {T <: AbstractFloat,
                                                                     A <: AbstractVector,
                                                                     F <: AbstractFeasibilityAlgorithm,
                                                                     I <: AbstractIntegerAlgorithm}
    subproblems::Vector{SubProblem{T,F,I}} = fetch(subworker)
    if isempty(subproblems)
       # Workers has nothing do to, return.
       return nothing
    end
    aggregation::AbstractAggregation = aggregator(length(subproblems), T)
    quit = false
    while true
        t::Int = try
            if isready(finalize)
                quit = true
                take!(finalize)
            else
                wait(work)
                take!(work)
            end
        catch err
            if err isa InvalidStateException
                # Master closed the work/finalize channel. Worker finished
                return nothing
            end
        end
        t == -1 && continue
        x::A = fetch(iterates, t)
        # Update subproblems
        update_known_decisions!(fetch(decisions), x)
        for subproblem in subproblems
            update_subproblem!(subproblem)
            cut::SparseHyperPlane{T} = subproblem(x, metadata)
            !quit && aggregate_cut!(cutqueue, aggregation, worker_metadata, t, cut, x)
        end
        !quit && flush!(cutqueue, aggregation, worker_metadata, t, x)
        if quit
            # Worker finished
            return nothing
        end
    end
end


text:
# a basic create robot type node example

using Caesar, RoME
using CloudGraphs
using IncrementalInference

# # Uncomment out for command line operation
# cloudGraph, addrdict = standardcloudgraphsetup(nparticles=true)
# Nparticles = parse(Int, addrdict["num particles"])
# session = addrdict["session"]

# interactive operation
session = "SESSROX"
Nparticles = 100
include(joinpath(dirname(@__FILE__),"blandauthremote.jl"))

fg = initfg(sessionname=session, cloudgraph=cloudGraph)
println("Working with session $(fg.sessionname)...")

N = 100
initCov = diagm([0.03;0.03;0.001].^2)
odoCov = diagm([3.0;3.0;0.01].^2)

# Some starting position
v1 = addVariable!(fg, :x1, zeros(3,1), diagm([1.0;1.0;0.1]), N=N,labels=["POSE"])
initPosePrior = PriorPose2(MvNormal(zeros(3), initCov) )
f1  = addFactor!(fg,[v1], initPosePrior)

ppc = Pose2Pose2(([50.0;0.0;pi/2]')', odoCov, [1.0])
v2, f2 = addOdoFG!( fg, ppc ,labels=["POSE"] )

v3, f3 = addOdoFG!( fg, ppc ,labels=["POSE"] )

# Some basic graph operations using the database
drawPoses(fg);

# get vertex back from DB
x1neoID = fg.cgIDs[fg.IDs[:x1]]
cv1r = CloudGraphs.get_vertex(fg.cg, x1neoID, false)
# Get neighbors
neighs = CloudGraphs.get_neighbors(fg.cg, cv1r)

# now solve using separate variable
fgs = initfg(sessionname=session, cloudgraph=cloudGraph)

setBackendWorkingSet!(conn, session)
fullLocalGraphCopy!(fgs)
tree = wipeBuildNewTree!(fgs)

inferOverTreeR!(fgs, tree)

@show fgs

# # Robot navigation and inference type stuff
# N=100
# doors = [-100.0;0.0;100.0;300.0]'
# cov = [3.0]
#
# # robot style, add first pose vertex
# v1 = addVariable!(fg,:x1,doors,N=N)
#
# # add a prior for the initial position of the robot
# f0  = addFactor!(fg,[v1], Obsv2(doors, cov', [1.0]))
#
# # add second pose vertex
# tem = 2.0*randn(1,N)+getVal(v1)+50.0
# v2 = addVariable!(fg, :x2, tem, N=N)
#
# # now add the odometry factor between them
# f1 = addFactor!(fg,[v1;v2],Odo([50.0]',[2.0]',[1.0]))

#


text:
# test_config.jl

@testset "Test configuration parameters" begin
    # by default NULLs come first
    con = DBInterface.connect(DuckDB.DB, ":memory:")

    results = DBInterface.execute(con, "SELECT 42 a UNION ALL SELECT NULL ORDER BY a")
    df = DataFrame(results)
    @test names(df) == ["a"]
    @test size(df, 1) == 2
    @test isequal(df.a, [missing, 42])

    DBInterface.close!(con)

    # if we add this configuration flag, nulls should come last
    config = DuckDB.Config()
    DuckDB.set_config(config, "default_null_order", "nulls_last")
    @test_throws DuckDB.QueryException DuckDB.set_config(config, "unrecognized option", "aaa")

    con = DBInterface.connect(DuckDB.DB, ":memory:", config)

    # NULL should come last now
    results = DBInterface.execute(con, "SELECT 42 a UNION ALL SELECT NULL ORDER BY a")
    df = DataFrame(results)
    @test names(df) == ["a"]
    @test size(df, 1) == 2
    @test isequal(df.a, [42, missing])

    DBInterface.close!(config)
    DBInterface.close!(config)
    DBInterface.close!(con)
end

max_stars_repo_name: lokax/duckdb

max_issues_repo_name: lokax/duckdb

max_forks_repo_name: lokax/duckdb


text:
module JuliaFormatter

using CSTParser
import CSTParser.Tokenize.Tokens

export format, format_text, format_file, format_dir

is_str_or_cmd(x::Tokens.Kind) =
    x in (Tokens.CMD, Tokens.TRIPLE_CMD, Tokens.STRING, Tokens.TRIPLE_STRING)

struct Document
    text::AbstractString

    ranges::Vector{UnitRange{Int}}
    line_to_range::Dict{Int,UnitRange{Int}}

    # mapping the offset in the file to the raw literal
    # string and what lines it starts and ends at.
    lit_strings::Dict{Int,Tuple{Int,Int,String}}
    comments::Dict{Int,String}

    # CSTParser does not detect semicolons.
    # It's useful to know where these are for
    # a few node types.
    semicolons::Set{Int}
end
function Document(text::AbstractString)
    ranges = UnitRange{Int}[]
    line_to_range = Dict{Int,UnitRange{Int}}()
    lit_strings = Dict{Int,Tuple{Int,Int,String}}()
    comments = Dict{Int,String}()
    semicolons = Set{Int}()

    for t in CSTParser.Tokenize.tokenize(text)
        if t.kind === Tokens.WHITESPACE
            offset = t.startbyte
            for c in t.val
                if c == '\n'
                    s = length(ranges) > 0 ? last(ranges[end]) + 1 : 1
                    push!(ranges, s:offset+1)
                end
                offset += 1
            end
        elseif t.kind === Tokens.ENDMARKER
            s = length(ranges) > 0 ? last(ranges[end]) + 1 : 1
            push!(ranges, s:t.startbyte)
        elseif is_str_or_cmd(t.kind)
            lit_strings[t.startbyte] = (t.startpos[1], t.endpos[1], t.val)
            if t.startpos[1] != t.endpos[1]
                offset = t.startbyte
                nls = findall(x -> x == '\n', t.val)
                for nl in nls
                    s = length(ranges) > 0 ? last(ranges[end]) + 1 : 1
                    push!(ranges, s:offset+nl)
                end
            end
        elseif t.kind === Tokens.COMMENT
            if t.startpos[1] == t.endpos[1]
                comments[t.startpos[1]] = t.val
            else
                sl = t.startpos[1]
                offset = t.startbyte
                comment_offset = 1
                for (i, c) in enumerate(t.val)
                    if c == '\n'
                        s = length(ranges) > 0 ? last(ranges[end]) + 1 : 1
                        push!(ranges, s:offset+1)
                        comments[sl] = t.val[comment_offset:i-1]
                        sl += 1
                        comment_offset = i + 1
                    end
                    offset += 1
                end
                # last comment
                comments[sl] = t.val[comment_offset:end]
            end
        elseif t.kind === Tokens.SEMICOLON
            push!(semicolons, t.startpos[1])
        end
    end
    for (l, r) in enumerate(ranges)
        line_to_range[l] = r
    end
    Document(text, ranges, line_to_range, lit_strings, comments, semicolons)
end

mutable struct State
    doc::Document
    indent_size::Int
    indent::Int
    offset::Int
    line_offset::Int
    margin::Int
end
State(doc, indent_size, margin) = State(doc, indent_size, 0, 1, 0, margin)

@inline nspaces(s::State) = s.indent
@inline getline(d::Document, line::Int) = d.text[d.line_to_range[line]]
@inline hascomment(d::Document, line::Int) = haskey(d.comments, line)
@inline has_semicolon(d::Document, line::Int) = line in d.semicolons

@inline function cursor_loc(s::State, offset::Int)
    for (l, r) in enumerate(s.doc.ranges)
        if offset in r
            return (l, offset - first(r) + 1, length(r))
        end
    end
    error("Indexing range 1 - $(last(s.doc.ranges[end])), index used = $(offset)")
end
@inline cursor_loc(s::State) = cursor_loc(s, s.offset)

include("pretty.jl")
include("nest.jl")
include("print.jl")

"""
    format_text(
        text::AbstractString;
        indent::Integer = 4,
        margin::Integer = 92,
    ) :: String

Formats a Julia source passed in as a string, returning the formatted
code as another string. The formatting options are:

- `indent` which defaults to 4 spaces
- `margin` which defaults to 92 columns

"""
function format_text(text::AbstractString; indent::Integer = 4, margin::Integer = 92,)
    if isempty(text)
        return text
    end

    x, ps = CSTParser.parse(CSTParser.ParseState(text), true)
    ps.errored && error("Parsing error for input $text")

    # no actual code
    x.args[1].kind === Tokens.NOTHING && length(x) == 1 && return text

    d = Document(text)
    # If "nofmt" occurs in a comment on line 1 do not format
    occursin("nofmt", get(d.comments, 1, "")) && (return text)

    s = State(d, indent, margin)
    t = pretty(x, s)
    hascomment(s.doc, t.endline) && (add_node!(t, InlineComment(t.endline), s))
    nest!(t, s)

    io = IOBuffer()

    # Print comments and whitespace before code.
    if t.startline > 1
        print_tree(io, Notcode(1, t.startline - 1), s)
        print_tree(io, Newline(), s)
    end

    if t.endline < length(s.doc.ranges)
        add_node!(t, Newline(), s)
        add_node!(t, Notcode(t.endline + 1, length(s.doc.ranges)), s)
    end

    print_tree(io, t, s)

    text = String(take!(io))

    _, ps = CSTParser.parse(CSTParser.ParseState(text), true)
    ps.errored && error("Parsing error for formatted $text")
    return text
end

"""
    format_file(
        filename::AbstractString;
        indent::Integer = 4,
        margin::Integer = 92,
        overwrite::Bool = true,
        verbose::Bool = false,
    )

Formats the contents of `filename` assuming it's a Julia source file.

The formatting options are:

- `indent` which defaults to 4 spaces
- `margin` which defaults to 92 columns

If `overwrite` is `true` the file will be reformatted in place, overwriting
the existing file; if it is `false`, the formatted version of `foo.jl` will
be written to `foo_fmt.jl` instead.

If `verbose` is `true` details related to formatting the file will be printed
to `stdout`.
"""
function format_file(
    filename::AbstractString;
    indent::Integer = 4,
    margin::Integer = 92,
    overwrite::Bool = true,
    verbose::Bool = false,
)
    path, ext = splitext(filename)
    if ext != ".jl"
        error("$filename must be a Julia (.jl) source file")
    end
    verbose && println("Formatting $filename with indent = $indent, margin = $margin")
    str = read(filename) |> String
    str = format_text(str, indent = indent, margin = margin)
    overwrite ? write(filename, str) : write(path * "_fmt" * ext, str)
    nothing
end

"""
    format(
        paths; # a path or collection of paths
        indent::Integer = 4,
        margin::Integer = 92,
        overwrite::Bool = true,
        verbose::Bool = false,
    )

Recursively descend into files and directories, formatting and `.jl`
files by calling `format_file` on them.
"""
function format(paths; options...)
    for path in paths
        if isfile(path)
            format_file(path; options...)
        else
            for (root, dirs, files) in walkdir(path)
                for file in files
                    _, ext = splitext(file)
                    ext == ".jl" || continue
                    full_path = joinpath(root, file)
                    ".git" in splitpath(full_path) && continue
                    format_file(full_path; options...)
                end
            end
        end
    end
    nothing
end
format(path::AbstractString; options...) = format((path,); options...)

end


text:
export row, col, reshape
export split, concat
export subtoind, indtosub
export stack, flatten, unstack, vflatten, unflatten
export riffle
export matrix, unmatrix
export lines, unlines, unzip
export findsub
export randsample, flip, flip!, flipdims

#######################################
## row, col, reshape

row(a::Number) = row([a])
row(a) = reshape(a,1,length(a))
row(a::Tuple) = reshape([a...],1,length(a))
row(a...) = row(vcat(a...))

col(a::Number) = col([a])
col(a) = reshape(a,length(a),1)
col(a::Tuple) = reshape([a...],length(a),1)
col(a...) = col(vcat(a...))

import Base.reshape
reshape{T}(a::AbstractArray, siz::Array{T,2}) = reshape(a, siz...)
function reshape(a::AbstractArray)
    if !(ndims(a)==1 || len(a)==1)
        error("FuncionalData.reshape(a): a was expected to be a vector or len(a)==1, actual size was: $(size(a))")
    end
    r2 = sqrt(length(a))
    r3 = cbrt(length(a))
    if round(r2)==r2
        return Base.reshape(a, round(Int,r2), round(Int,r2))
    elseif round(r3)==r3
        return Base.reshape(a, round(Int,r3), round(Int,r3), round(Int,r3))
    else
        error("cannot reshape array of size $(size(a)) to a square or cube")
    end
end

#######################################
##  split, cat

import Base.split
split(a::AbstractArray,x) = split(a, y->y==x)

function split(a::AbstractArray,f::Function)
    start = 1
    r = Any[]
    for i = 1:length(a)
        if f(a[i])
            push!(r, a[start:i-1])
            start = i + 1
        elseif i == length(a)
            push!(r, a[start:i])
        end
    end
    return r
end

concat(a) = concat(a...)
concat(a...) = @p flatten Any[reject(collect(a),x->(!isa(x,Symbol) && isempty(x)))...]

#######################################
## subtoind, indtosub

function subtoind(subs, a)
    strides_ = strides(a)
    r = oneel(subs)
    for i in 1:length(strides_)
        r +=(subs[i]-1) * strides_[i]
    end
    return r
end

indtosub(i::Int, a) = indtosub([i], a)
function indtosub(inds::AbstractArray, a::Union{Array,BitArray})
    @p ind2sub size(a) inds | unstack | map row | col | flatten
end

#######################################
## flatten, stack, unstack, unflatten

function stack{T}(a::Array{T,1})
    r = arraylike(fst(a), len(a))
    for i = 1:len(a)
        setat!(r, i, at(a,i))
    end
    return r
end
stack{T<:Real}(a::DenseArray{T}) = a

typealias StringLike Union{Char, AbstractString}
tostring(a) = string(a)
tostring(a::AbstractString) = a
flatten(a::StringLike) = a
flatten{T<:StringLike}(a::Array{T,1}) = join(a)

function flatten{T<:StringLike}(a::Array{T,2})
    assert(size(a,1)==1)
    join(vec(a))
end

flatten{T<:Real,N}(a::AbstractArray{T,N}) = a
function flatten{T}(a::Array{T,1})
    if isempty(a)
        return similar(a)
    end
    if isa(a[1], StringLike)
        return join(map(a,tostring))
    end
    if VERSION.minor > 3 && !method_exists(ndims,Tuple{typeof(fst(a))})
        return a
    end
    if ndims(fst(a)) == 1
        return vcat(a...)
    end
    r = arraylike(fst(fst(a)), sum([len(x) for x in a]))
    ind = 0
    for i = 1:len(a)
        for j = 1:len(a[i])
            setat!(r, ind+j, at(a[i],j))
        end
        ind += len(a[i])
    end
    return r
end

function flatten{T}(a::Array{T,2})
    if isempty(a)
        return arraylike(a)
    end
    if isa(a[1], StringLike)
        assert(size(a,1) == 1)
        return flatten(vec(a))
    end
    ms = Base.map(x->size(x,1), a)
    ns = Base.map(x->size(x,2), a)
    msum = 0
    nsum = 0
    for mi = 1:size(a,1)
        for ni = 2:size(a,2)
            if ms[mi,ni] != ms[mi,1]
                error("flatten: ms[mi,ni] was $(ms[mi,ni]) but ms[mi,1] was $(ms[mi,1])")
            end
        end
        msum += ms[mi,1]
    end
    for ni = 1:size(a,2)
        for mi = 2:size(a,1)
            if ns[mi,ni] != ns[1,ni]
                error("flatten: ns[mi,ni] was $(ns[mi,ni]) but ns[1,ni] was $(ns[1,ni])")
            end
        end
        nsum += ns[1,ni]
    end

    typ = eltype(a[1])
    for i = 1:length(a)
        typ = promote_type(typ, eltype(a[i]))
    end

    r = Array(typ, msum, nsum)

    ncum = 0
    for n = 1:size(a,2)
        mcum = 0
        for m = 1:size(a,1)
            r[mcum + (1:ms[m,n]), ncum + (1:ns[m,n])] = a[m,n]
            mcum += ms[m,n]
        end
        ncum += ns[1,n]
    end
    r
end

vflatten(a) = @p transpose a | flatten | transpose

unstack(a) = Any[at(a,i) for i in 1:len(a)]

function unflatten(a,template)
    lens = @p mapvec template len
    ends = cumsum(lens)
    starts = @p droplast [1; ends+1]
    @p mapvec2 starts ends (i,j)->part(a,i:j)
end

#######################################
##  riffle

function riffle(a::AbstractString,x::Union{Char,AbstractString})
    join(a,x)
end

function riffle(a::Array{Any,1},x)
    if len(a)==1 return a end

    dims = size(a)
    #dims[end] = dims[end]*2-1
    dims = tuple(dims[1:end-1]..., dims[end]*2-1)
    r = similar(a,dims...)

    for i = 1:len(a)
        r[i*2-1] = a[i]
    end
    for i = 2:2:len(r)-1
        r[i] = x
    end
    return r
end

function riffle(a,x)
    if len(a)==1 return a end

    dims = size(a)
    #dims[end] = dims[end]*2-1
    dims = tuple(dims[1:end-1]..., dims[end]*2-1)
    r = similar(a,dims...)
    # println("dims: $dims   r: $r")
    for i = 1:len(a)
        setat!(r, i*2-1, at(a,i))
    end
    for i = 2:2:len(r)-1
        setat!(r, i, x)
    end
    return r
end

#######################################
##  matrix, unmatrix

matrix{T<:Number}(a::AbstractArray{T}) = @p reshape a div(length(a),len(a)) len(a)
matrix(a) = @p map a col | flatten
function unmatrix(a, example)
    r = @p map a (x->reshape(x,siz(fst(example))))
    if eltype(example)==Any
        r = unstack(r)
    end
    r
end

#######################################
##  lines, unlines

function lines(a)
    r = split(a,r"\r\n|\n")
    if !isempty(r) && last(r)==""
        droplast(r)
    else
        r
    end
end
unlines(a) = join(a,'\n')

#######################################
##  randperm, randsample

import Base.randperm
_randperm = randperm
randperm(a::Number) = _randperm(a)
randperm(a) = part(a, randperm(len(a)))

randsample(a, n = 1) = part(a, rand(1:len(a), round(Int,n)))

flip!(a::Array) = a[:] = flip(a)
flip(a::Pair) = Pair(snd(a),fst(a))
flip(a::Dict) = @p vec a | map flip | Dict
flip(a) = part(a, len(a):-1:1)
function flipdims(a,d1,d2)
    dims = collect(1:ndims(a))
    dims[d1] = d2
    dims[d2] = d1
    permutedims(a, dims)
end

#####################################################
##   unzip

function unzip(a)
    forth(a) = at(a,4)
    fifth(a) = at(a,5)
    if len(fst(a)) == 2
        map(a, fst), map(a, snd)
    elseif len(fst(a)) == 3
        map(a, fst), map(a, snd), map(a, third)
    elseif len(fst(a)) == 4
        map(a, fst), map(a, snd), map(a, third), map(a, forth)
    elseif len(fst(a)) == 5
        map(a, fst), map(a, snd), map(a, third), map(a, forth), map(a, fifth)
    else
        error("FunctionalData.unzip: cannot unzip of len(fst(a))==$(len(fst(a)))")
    end
end

#####################################################
##   find

findsub(a::Union{Array,BitArray}) = indtosub(find(a.!=0), a)


text:
using PackageCompiler
@warn "Manually uncomment MOI test in runtests.jl. This test breaks PackageCompiler."
create_sysimage([:CUDA, :ExaPF]; sysimage_path="exapf.so", precompile_execution_file="test/runtests.jl")


text:
module Test

using DataFrames

df = DataFrame(a=ASCIIString[], b=Int[])

push!(df, ["Hello", 1])
push!(df, ["World", 2])

show(df)

end


text:
using Plots, DelimitedFiles

function data_theoretical(::Type{T}, ns) where T<:Union{Float64,Int64}
    # 0.0581s for 1000 x 1000, Float64
    map(n->n^3/1000^3*0.0862, ns)
end

function data_theoretical(::Type{T}, ns) where T<:Union{Float32}
    # 0.0581s for 1000 x 1000, Float64
    map(n->n^3/1000^3*0.0862/2, ns)
end

function plot_res(::Type{T}; dir=@__DIR__, max_threads=6) where T
    x = 2 .^ (1:12)
    ID1 = if T === Float64
        2
    elseif T === Float32
        1
    elseif T === Int
        3
    else
        error("not defined for type $T")
    end
    ax = plot(T===Int ? [] : x, T===Int ? [] : data_theoretical(T, x), label="theoretical", legendfont=8, legend=:topleft, ylabel="seconds", xlabel="n", xscale=:log10, yscale=:log10)
    data = reshape(readdlm(joinpath(dir, "benchmarks.dat")), 6, 12)
    plot!(ax, x, (data')[:,[2*ID1-1, 2*ID1]] ./ 1e9, label=["Float32-naive" "Float32-serial" "Float64-naive" "Float64-serial" "Int-naive" "Int-serial"][:,[2*ID1-1,2*ID1]])
    for nthreads in 2:max_threads
        data = reshape(readdlm(joinpath(dir, "benchmarks-$nthreads.dat")), 3, 12)
        plot!(ax, x, (data')[:,ID1] .* nthreads ./ 1e9, label=["Float32-$nthreads" "Float64-$nthreads" "Int-$nthreads"][ID1])
    end
    return ax
end