text:
function simple_binary_operation(op, a, b, num1, num2)
    a_b = op(a, b)
    interior(a) .= num1
    interior(b) .= num2
    return a_b[2, 2, 2] == op(num1, num2)
end

function three_field_addition(a, b, c, num1, num2)
    a_b_c = a + b + c
    interior(a) .= num1
    interior(b) .= num2
    interior(c) .= num2
    return a_b_c[2, 2, 2] == num1 + num2 + num2
end

function x_derivative(a)
    dx_a = ∂x(a)

    for k in 1:3
        interior(a)[:, 1, k] .= [1, 2, 3]
        interior(a)[:, 2, k] .= [1, 2, 3]
        interior(a)[:, 3, k] .= [1, 2, 3]
    end

    return dx_a[2, 2, 2] == 1
end

function y_derivative(a)
    dy_a = ∂y(a)

    for k in 1:3
        interior(a)[1, :, k] .= [1, 3, 5]
        interior(a)[2, :, k] .= [1, 3, 5]
        interior(a)[3, :, k] .= [1, 3, 5]
    end

    return dy_a[2, 2, 2] == 2
end

function z_derivative(a)
    dz_a = ∂z(a)

    for k in 1:3
        interior(a)[1, k, :] .= [1, 4, 7]
        interior(a)[2, k, :] .= [1, 4, 7]
        interior(a)[3, k, :] .= [1, 4, 7]
    end

    return dz_a[2, 2, 2] == 3
end

function x_derivative_cell(FT, arch)
    grid = RegularCartesianGrid(FT, size=(3, 3, 3), extent=(3, 3, 3))
    a = Field(Cell, Cell, Cell, arch, grid, nothing)
    dx_a = ∂x(a)

    for k in 1:3
        interior(a)[:, 1, k] .= [1, 4, 4]
        interior(a)[:, 2, k] .= [1, 4, 4]
        interior(a)[:, 3, k] .= [1, 4, 4]
    end

    return dx_a[2, 2, 2] == 3
end

function times_x_derivative(a, b, location, i, j, k, answer)
    a∇b = @at location b * ∂x(a)
    return a∇b[i, j, k] == answer
end

function compute_derivative(model, ∂)
    T, S = model.tracers
    S.data.parent .= π

    computation = Computation(∂(S), model.pressures.pHY′)
    compute!(computation)
    result = Array(interior(computation.result))

    return all(result .≈ zero(eltype(model.grid)))
end

function compute_unary(unary, model)
    set!(model; S=π)
    T, S = model.tracers

    computation = Computation(unary(S), model.pressures.pHY′)
    compute!(computation)
    result = Array(interior(computation.result))
    return all(result .≈ unary(eltype(model.grid)(π)))
end

function compute_plus(model)
    set!(model; S=π, T=42)
    T, S = model.tracers

    computation = Computation(S + T, model.pressures.pHY′)
    compute!(computation)
    result = Array(interior(computation.result))

    return all(result .≈ eltype(model.grid)(π+42))
end

function compute_many_plus(model)
    set!(model; u=2, S=π, T=42)
    T, S = model.tracers
    u, v, w = model.velocities

    computation = Computation(@at((Cell, Cell, Cell), u + T + S), model.pressures.pHY′)
    compute!(computation)
    result = Array(interior(computation.result))

    return all(result .≈ eltype(model.grid)(2+π+42))
end

function compute_minus(model)
    set!(model; S=π, T=42)
    T, S = model.tracers

    computation = Computation(S - T, model.pressures.pHY′)
    compute!(computation)
    result = Array(interior(computation.result))

    return all(result .≈ eltype(model.grid)(π-42))
end

function compute_times(model)
    set!(model; S=π, T=42)
    T, S = model.tracers

    computation = Computation(S * T, model.pressures.pHY′)
    compute!(computation)
    result = Array(interior(computation.result))

    return all(result .≈ eltype(model.grid)(π*42))
end

function compute_kinetic_energy(model)
    u, v, w = model.velocities
    set!(u, 1)
    set!(v, 2)
    set!(w, 3)

    kinetic_energy = @at (Cell, Cell, Cell) (u^2 + v^2 + w^2) / 2
    computation = Computation(kinetic_energy, model.pressures.pHY′)
    compute!(computation)
    result = Array(interior(computation.result))

    return all(result .≈ 7)
end

function horizontal_average_of_plus(model)
    S₀(x, y, z) = sin(π*z)
    T₀(x, y, z) = 42*z
    set!(model; S=S₀, T=T₀)
    T, S = model.tracers

    ST = HorizontalAverage(S + T, model)
    computed_profile = ST(model)

    zC = znodes(Cell, model.grid)
    correct_profile = @. sin(π * zC) + 42 * zC

    return all(computed_profile[2:end-1] .≈ correct_profile)
end

function horizontal_average_of_minus(model)
    S₀(x, y, z) = sin(π*z)
    T₀(x, y, z) = 42*z
    set!(model; S=S₀, T=T₀)
    T, S = model.tracers

    ST = HorizontalAverage(S - T, model)
    computed_profile = ST(model)

    zC = znodes(Cell, model.grid)
    correct_profile = @. sin(π * zC) - 42 * zC

    return all(computed_profile[2:end-1] .≈ correct_profile)
end

function horizontal_average_of_times(model)
    S₀(x, y, z) = sin(π*z)
    T₀(x, y, z) = 42*z
    set!(model; S=S₀, T=T₀)
    T, S = model.tracers

    ST = HorizontalAverage(S * T, model)
    computed_profile = ST(model)

    zC = znodes(Cell, model.grid)
    correct_profile = @. sin(π * zC) * 42 * zC

    return all(computed_profile[2:end-1] .≈ correct_profile)
end

function multiplication_and_derivative_ccf(model)
    w₀(x, y, z) = sin(π * z)
    T₀(x, y, z) = 42 * z
    set!(model; w=w₀, T=T₀)

    w = model.velocities.w
    T = model.tracers.T

    wT = HorizontalAverage(w * ∂z(T), model)
    computed_profile = wT(model)

    zF = znodes(Face, model.grid)
    correct_profile = @. 42 * sin(π * zF)

    # Omit both halos and boundary points
    return all(computed_profile[3:end-1] .≈ correct_profile[2:end-1])
end

const C = Cell
const F = Face

function multiplication_and_derivative_ccc(model)
    w₀(x, y, z) = sin(π*z)
    T₀(x, y, z) = 42*z
    set!(model; w=w₀, T=T₀)

    w = model.velocities.w
    T = model.tracers.T

    wT_ccc = @at (C, C, C) w * ∂z(T)
    wT_ccc_avg = HorizontalAverage(wT_ccc, model)
    computed_profile_ccc = wT_ccc_avg(model)

    zF = znodes(Face, model.grid)
    sinusoid = sin.(π * zF)
    interped_sin = [(sinusoid[k] + sinusoid[k+1]) / 2 for k in 1:model.grid.Nz]
    correct_profile = interped_sin .* 42

    # Computed profile includes halos.
    return all(computed_profile_ccc[:][3:end-2] .≈ correct_profile[2:end-1])
end

@testset "Abstract operations" begin
    @info "Testing abstract operations..."

    for FT in float_types
        arch = CPU()
        grid = RegularCartesianGrid(FT, size=(3, 3, 3), extent=(3, 3, 3))
        u, v, w = VelocityFields(arch, grid)
        c = Field(Cell, Cell, Cell, arch, grid, nothing)

        @testset "Unary operations and derivatives [$FT]" begin
            for ψ in (u, v, w, c)
                for op_symbol in Oceananigans.AbstractOperations.unary_operators
                    op = eval(op_symbol)
                    @test typeof(op(ψ)[2, 2, 2]) <: Number
                end

                for d_symbol in Oceananigans.AbstractOperations.derivative_operators
                    d = eval(d_symbol)
                    @test typeof(d(ψ)[2, 2, 2]) <: Number
                end
            end
        end

        @testset "Binary operations [$FT]" begin
            generic_function(x, y, z) = x + y + z
            for (ψ, ϕ) in ((u, v), (u, w), (v, w), (u, c), (generic_function, c), (u, generic_function))
                for op_symbol in Oceananigans.AbstractOperations.binary_operators
                    op = eval(op_symbol)
                    @test typeof(op(ψ, ϕ)[2, 2, 2]) <: Number
                end
            end
        end

        @testset "Multiary operations [$FT]" begin
            generic_function(x, y, z) = x + y + z
            for (ψ, ϕ, σ) in ((u, v, w), (u, v, c), (u, v, generic_function))
                for op_symbol in Oceananigans.AbstractOperations.multiary_operators
                    op = eval(op_symbol)
                    @test typeof(op((Cell, Cell, Cell), ψ, ϕ, σ)[2, 2, 2]) <: Number
                end
            end
        end
    end

    @testset "Fidelity of simple binary operations" begin
        arch = CPU()
        @info "  Testing simple binary operations..."
        for FT in float_types
            num1 = FT(π)
            num2 = FT(42)
            grid = RegularCartesianGrid(FT, size=(3, 3, 3), extent=(3, 3, 3))

            u, v, w = VelocityFields(arch, grid)
            T, S = TracerFields(arch, grid, (:T, :S))

            for op in (+, *, -, /)
                @test simple_binary_operation(op, u, v, num1, num2)
                @test simple_binary_operation(op, u, w, num1, num2)
                @test simple_binary_operation(op, u, T, num1, num2)
                @test simple_binary_operation(op, T, S, num1, num2)
            end
            @test three_field_addition(u, v, w, num1, num2)
        end
    end

    @testset "Derivatives" begin
        arch = CPU()
        @info "  Testing derivatives..."
        for FT in float_types
            grid = RegularCartesianGrid(FT, size=(3, 3, 3), extent=(3, 3, 3),
                                        topology=(Periodic, Periodic, Periodic))

            u, v, w = VelocityFields(arch, grid)
            T, S = TracerFields(arch, grid, (:T, :S))
            for a in (u, v, w, T)
                @test x_derivative(a)
                @test y_derivative(a)
                @test z_derivative(a)
            end
            @test x_derivative_cell(FT, arch)
        end
    end

    @testset "Combined binary operations and derivatives" begin
        @info "  Testing combined binary operations and derivatives..."
        arch = CPU()
        Nx = 3 # Δx=1, xC = 0.5, 1.5, 2.5
        for FT in float_types
            grid = RegularCartesianGrid(FT, size=(Nx, Nx, Nx), extent=(Nx, Nx, Nx))
            a, b = (Field(Cell, Cell, Cell, arch, grid, nothing) for i in 1:2)

            set!(b, 2)
            set!(a, (x, y, z) -> x < 2 ? 3x : 6)

            #                            0   0.5   1   1.5   2   2.5   3
            # x -▶                  ∘ ~~~|--- * ---|--- * ---|--- * ---|~~~ ∘
            #        i Face:    0        1         2        3          4
            #        i Cell:        0         1         2         3         4

            #              a = [    0,       1.5,      4.5,       6,        0    ]
            #              b = [    0,        2,        2,        2,        0    ]
            #          ∂x(a) = [        1.5,       3,       1.5,      -6         ]

            # x -▶                  ∘ ~~~|--- * ---|--- * ---|--- * ---|~~~ ∘
            #        i Face:    0        1         2         3         4
            #        i Cell:        0         1         2         3         4

            # ccc: b * ∂x(a) = [             4.5,      4.5      -4.5,            ]
            # fcc: b * ∂x(a) = [         3,        6,        3,       -6         ]

            @test times_x_derivative(a, b, (C, C, C), 1, 2, 2, 4.5)
            @test times_x_derivative(a, b, (F, C, C), 1, 2, 2, 3)

            @test times_x_derivative(a, b, (C, C, C), 2, 2, 2, 4.5)
            @test times_x_derivative(a, b, (F, C, C), 2, 2, 2, 6)

            @test times_x_derivative(a, b, (C, C, C), 3, 2, 2, -4.5)
            @test times_x_derivative(a, b, (F, C, C), 3, 2, 2, 3)
        end
    end

    for arch in archs
        @testset "Computations [$(typeof(arch))]" begin
            @info "  Testing combined binary operations and derivatives..."
            for FT in float_types
                @info "    Testing computation of abstract operations [$FT, $(typeof(arch))]..."

                model = IncompressibleModel(
                    architecture = arch,
                      float_type = FT,
                            grid = RegularCartesianGrid(FT, size=(16, 16, 16), extent=(1, 1, 1),
                                                        topology=(Periodic, Periodic, Bounded))
                )

                @testset "Derivative computations [$FT, $(typeof(arch))]" begin
                    @info "      Testing compute! derivatives..."
                    @test compute_derivative(model, ∂x)
                    @test compute_derivative(model, ∂y)
                    @test compute_derivative(model, ∂z)
                end

                @testset "Unary computations [$FT, $(typeof(arch))]" begin
                    @info "      Testing compute! unary operations..."
                    for unary in Oceananigans.AbstractOperations.unary_operators
                        @test compute_unary(eval(unary), model)
                    end
                end

                @testset "Binary computations [$FT, $(typeof(arch))]" begin
                    @info "      Testing compute! binary operations..."
                    @test compute_plus(model)
                    @test compute_minus(model)
                    @test compute_times(model)
                end

                @testset "Multiary computations [$FT, $(typeof(arch))]" begin
                    @info "      Testing compute! multiary operations..."
                    @test compute_many_plus(model)

                    @info "      Testing compute! kinetic energy..."
                    @test compute_kinetic_energy(model)
                end

                @testset "Horizontal averages of operations [$FT, $(typeof(arch))]" begin
                    @info "      Testing horizontal averges..."
                    @test horizontal_average_of_plus(model)
                    @test horizontal_average_of_minus(model)
                    @test horizontal_average_of_times(model)

                    @test multiplication_and_derivative_ccf(model)
                    @test multiplication_and_derivative_ccc(model)
                end
            end
        end
    end
end


text:
#!/usr/bin/env julia
using RobotOS
@rosimport obstacle_detector.msg: Obstacles, CircleObstacle
rostypegen()
using obstacle_detector.msg

function setObstacleParams(msg::Obstacles)
  L = length(msg.circles)

  if L > 0
    r=zeros(L);x=r;y=r;vx=r;vy=r;
    r = (); x = (); y = (); vx = (); vy = ();
    for i in 1:L
      r = (r..., msg.circles[i].radius)
      x = (x..., msg.circles[i].center.x)
      y = (y..., msg.circles[i].center.y)
      vx = (vx..., msg.circles[i].velocity.x)
      vy = (vy..., msg.circles[i].velocity.y)
    end

    # update obstacle field parameters
    RobotOS.set_param("obstacle_radius",r)
    RobotOS.set_param("obstacle_x",x)
    RobotOS.set_param("obstacle_y",y)
    RobotOS.set_param("obstacle_vx",vx)
    RobotOS.set_param("obstacle_vy",vy)
  else
    RobotOS.set_param("obstacle_radius",NaN)
    RobotOS.set_param("obstacle_x",NaN)
    RobotOS.set_param("obstacle_y",NaN)
    RobotOS.set_param("obstacle_vx",NaN)
    RobotOS.set_param("obstacle_vy",NaN)
  end

  return nothing
end

init_node("obstacle_params")
sub = Subscriber{Obstacles}("/obstacles", setObstacleParams, queue_size = 10)
spin()

ext: jl


text:
module TestDuplicates
    using Base.Test
    using DataTables

    dt = DataTable(a = [1, 2, 3, 3, 4])
    udt = DataTable(a = [1, 2, 3, 4])
    @test isequal(nonunique(dt), [false, false, false, true, false])
    @test isequal(udt, unique(dt))
    unique!(dt)
    @test isequal(dt, udt)

    pdt = DataTable(a = NullableCategoricalArray(Nullable{String}["a", "a", Nullable(),
                                             Nullable(), "b", Nullable(), "a", Nullable()]),
                    b = NullableCategoricalArray(Nullable{String}["a", "b", Nullable(),
                                                              Nullable(), "b", "a", "a", "a"]))
    updt = DataTable(a = NullableCategoricalArray(Nullable{String}["a", "a", Nullable(), "b", Nullable()]),
                     b = NullableCategoricalArray(Nullable{String}["a", "b", Nullable(), "b", "a"]))
    @test isequal(nonunique(pdt), [false, false, false, true, false, false, true, true])
    @test isequal(nonunique(updt), falses(5) )
    @test isequal(updt, unique(pdt))
    unique!(pdt)
    @test isequal(pdt, updt)

    @testset "missing" begin
        dt = DataTable(A = 1:12, B = repeat('A':'C', inner=4))
        @test DataTables.colmissing(dt) == [0, 0]
    end
end


text:
const shard_type_union = Union{
  GlobalRef, Module, QuoteNode,
  Symbol, AbstractString,
  Integer, AbstractFloat,
  VersionNumber
}

const cur_update_operators = [
  :(+=), :(-=), :(*=), :(/=), :(\=),
  :(÷=), :(%=), :(^=), :(&=), :(|=),
  :($=), :(>>>=), :(>>=), :(<<=)
]

const cur_nested_shards = [
  :->,
  :macrocall, :for,
  :block, :if,
  :call
]

const cur_return_self_shards = [
  :toplevel, :inbounds, :import, :break, :vect,
  :return, :tuple, :where, Symbol("'"),
  :curly, :meta, :ref, :try, :line, :string,
  :(&&), :(||), :(<:), :bitstype, :typealias,
  :(::), :(:), :(.), :(...)
]

const cur_return_nothing_shards = [
  :function, :abstract, :export, :using,
  :const, :type, :stagedfunction,
  :immutable, :module, :macro
]


text:
using DelayDiffEq
using DiffEqDevTools
using LinearAlgebra
using Test

@testset "SIR" begin
    function sir_dde!(du, u, h, p, t)
        S, I, _ = u
        γ, τ = p
        infection = γ * I * S
        Sd, Id, _ = h(p, t-τ)
        recovery = γ * Id * Sd
        @inbounds begin
            du[1] = -infection
            du[2] = infection - recovery
            du[3] = recovery
        end
        nothing
    end
    u0 = [0.99, 0.01, 0.0]
    sir_history(p, t) = [1.0, 0.0, 0.0]
    tspan = (0.0, 40.0)
    p = (γ = 0.5, τ = 4.0)

    prob_dde = DDEProblem(sir_dde!, u0, sir_history, tspan, p; constant_lags = (p.τ,))
    sol_dde = TestSolution(
        solve(prob_dde, MethodOfSteps(Vern9()); reltol=1e-14, abstol=1e-14)
    )

    function sir_ddae!(du, u, h, p, t)
        S, I, R = u
        γ, τ = p
        infection = γ * I * S
        Sd, Id, _ = h(p, t-τ)
        recovery = γ * Id * Sd
        @inbounds begin
            du[1] = -infection
            du[2] = infection - recovery
            du[3] = S + I + R - 1
        end
        nothing
    end

    prob_ddae = DDEProblem(
        DDEFunction{true}(sir_ddae!; mass_matrix=Diagonal([1.0, 1.0, 0.0])),
        u0,
        sir_history,
        tspan,
        p;
        constant_lags = (p.τ,),
    )
    ode_f = DelayDiffEq.ODEFunctionWrapper(prob_ddae.f, prob_ddae.h)
    @test ode_f.mass_matrix == Diagonal([1.0, 1.0, 0.0])

    int = init(prob_ddae, MethodOfSteps(Rosenbrock23()))
    @test int.isdae
    @test int.f.mass_matrix == Diagonal([1.0, 1.0, 0.0])

    for (alg, reltol) in (
        (Rosenbrock23(), nothing),
        (Rodas4(), nothing),
        (QNDF(), 1e-6),
        (Trapezoid(), 1e-6),
    )
        sol_ddae = solve(prob_ddae, MethodOfSteps(alg); reltol=reltol)
        sol = appxtrue(sol_ddae, sol_dde)
        @test sol.errors[:L∞] < 5e-3
    end
end

max_stars_repo_name: JuliaDiffEq/DelayDiffEq

max_issues_repo_name: JuliaDiffEq/DelayDiffEq

max_forks_repo_name: JuliaDiffEq/DelayDiffEq


text:
# special cases for datasets that need more rows for proper column type detection
const Dataset_typedetect_rows = Dict{Tuple{String, String}, Int}(
    ("adehabitatLT", "albatross") => 1000,
    ("adehabitatLT", "bear") => 1200,
    ("adehabitatLT", "buffalo") => 1500,
    ("adehabitatLT", "capreotf") => 600,
    ("boot", "neuro") => 500,
    ("gap", "PD") => 500,
    ("plyr", "baseball") => 1000,
    ("psych", "bfi") => 1000,
    ("survival", "cancer") => 250,
    ("survival", "lung") => 250,
    ("vcd", "Bundesliga") => 15000,
)

function dataset(package_name::AbstractString, dataset_name::AbstractString)
    basename = joinpath(dirname(@__FILE__), "..", "data", package_name)

    rdaname = joinpath(basename, string(dataset_name, ".rda"))
    if isfile(rdaname)
        return load(rdaname)[dataset_name]
    end

    csvname = joinpath(basename, string(dataset_name, ".csv.gz"))
    if isfile(csvname)
        return open(GzipDecompressorStream, csvname, "r") do io
            CSV.read(io, delim=',', quotechar='\"', missingstring="NA",
                     rows_for_type_detect=get(Dataset_typedetect_rows, (package_name, dataset_name), 200))
        end
    end
    error(@sprintf "Unable to locate dataset file %s or %s" rdaname csvname)
end

max_stars_repo_name: SCXFPROTOCOL/GameRewards-System

max_issues_repo_name: SCXFPROTOCOL/GameRewards-System

max_forks_repo_name: SCXFPROTOCOL/GameRewards-System


text:
# ---
# title: 1304. Find N Unique Integers Sum up to Zero
# id: problem1304
# author: Tian Jun
# date: 2020-10-31
# difficulty: Easy
# categories: Array
# link: <https://leetcode.com/problems/find-n-unique-integers-sum-up-to-zero/description/>
# hidden: true
# ---
#
# Given an integer `n`, return **any** array containing `n` **unique**  integers
# such that they add up to 0.
#
#
#
# **Example 1:**
#
#
#
#     Input: n = 5
#     Output: [-7,-1,1,3,4]
#     Explanation: These arrays also are accepted [-5,-1,1,2,3] , [-3,-1,2,-2,4].
#
#
# **Example 2:**
#
#
#
#     Input: n = 3
#     Output: [-1,0,1]
#
#
# **Example 3:**
#
#
#
#     Input: n = 1
#     Output: [0]
#
#
#
#
# **Constraints:**
#
#   * `1 <= n <= 1000`
#
#
## @lc code=start
using LeetCode

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


text:
#=
Dictionary with default values
=#

using SimpleDataStructures

str = "hello hello hello bye"
words = string.(split(str))
sdd = SimpleDefaultDict{String,Int}(0; capacity = 100)

for word in words
    sdd[word] += 1
end

sdd


text:
# Exercises from Chapter 4 of "Think Julia"
using ThinkJulia

🐢 = Turtle()

## There appears to be a problem building Cairo. Here is the last 100 lines of the log:
# Versions and API surfaces shifting like sand is not what I call 'programming'.
# Debugging is important, so I'll get this figured out!
# Maybe there's an opportunity to make a pull request to Cairo or ThinkJulia or Luxor or something.
# Then THAT would be programming :)

# Error: Error building `Cairo`, showing the last 100 of log:
# │ checking for grep that handles long lines and -e... /usr/bin/grep
# │ checking for egrep... /usr/bin/grep -E
# │ checking whether we are using the GNU C Library 2.1 or newer... no
# │ checking whether to enable garbage collector friendliness... no
# │ checking whether to disable memory pools... no
# │ readlink: illegal option -- f
# │ usage: readlink [-n] [file ...]
# │ checking for c++... c++
# │ checking whether we are using the GNU C++ compiler... yes
# │ checking whether c++ accepts -g... yes
# │ checking dependency style of c++... gcc3
# │ checking for gcc... (cached) gcc
# │ checking whether we are using the GNU C compiler... (cached) yes
# │ checking whether gcc accepts -g... (cached) yes
# │ checking for gcc option to accept ISO C89... (cached) none needed
# │ checking dependency style of gcc... (cached) gcc3
# │ checking whether gcc and cc understand -c and -o together... yes
# │ checking for special C compiler options needed for large files... no
# │ checking for _FILE_OFFSET_BITS value needed for large files... no
# │ checking for pkg-config... no
# │ checking for gawk... (cached) awk
# │ checking for perl5... no
# │ checking for perl... perl
# │ checking for indent... indent
# │ checking for perl... /usr/bin/perl
# │ checking for a Python interpreter with version >= 2.5... python
# │ checking for python... /usr/local/anaconda3/bin/python
# │ checking for python version... 3.7
# │ checking for python platform... darwin
# │ checking for python script directory... ${prefix}/lib/python3.7/site-packages
# │ checking for python extension module directory... ${exec_prefix}/lib/python3.7/site-packages
# │ checking for iconv_open... no
# │ checking for libiconv_open in -liconv... no
# │ checking for iconv_open in -liconv... yes
# │ checking whether to cache iconv descriptors... yes
# │ checking for ZLIB... no
# │ checking for ANSI C header files... yes
# │ checking for sys/types.h... yes
# │ checking for sys/stat.h... yes
# │ checking for stdlib.h... yes
# │ checking for string.h... yes
# │ checking for memory.h... yes
# │ checking for strings.h... yes
# │ checking for inttypes.h... yes
# │ checking for stdint.h... yes
# │ checking for unistd.h... yes
# │ checking for inflate in -lz... yes
# │ checking zlib.h usability... yes
# │ checking zlib.h presence... yes
# │ checking for zlib.h... yes
# │ checking for LIBFFI... no
# │ configure: error: in `/Users/xronos/.julia/packages/Cairo/p68X8/deps/builds/gobject':
# │ configure: error: The pkg-config script could not be found or is too old.  Make sure it
# │ is in your PATH or set the PKG_CONFIG environment variable to the full
# │ path to pkg-config.
# │
# │ Alternatively, you may set the environment variables LIBFFI_CFLAGS
# │ and LIBFFI_LIBS to avoid the need to call pkg-config.
# │ See the pkg-config man page for more details.
# │
# │ To get pkg-config, see <http://pkg-config.freedesktop.org/>.
# │ See `config.log' for more details
# │ [ Info: Attempting to create directory /Users/xronos/.julia/packages/Cairo/p68X8/deps/downloads
# │ [ Info: Directory /Users/xronos/.julia/packages/Cairo/p68X8/deps/downloads already exists
# │ [ Info: Downloading file http://ftp.gnome.org/pub/gnome/sources/glib/2.34/glib-2.34.3.tar.xz
# │ [ Info: Done downloading file http://ftp.gnome.org/pub/gnome/sources/glib/2.34/glib-2.34.3.tar.xz
# │ [ Info: Attempting to create directory /Users/xronos/.julia/packages/Cairo/p68X8/deps/src
# │ [ Info: Directory /Users/xronos/.julia/packages/Cairo/p68X8/deps/src already exists
# │ [ Info: Attempting to create directory /Users/xronos/.julia/packages/Cairo/p68X8/deps
# │ [ Info: Directory /Users/xronos/.julia/packages/Cairo/p68X8/deps already exists
# │ [ Info: Path /Users/xronos/.julia/packages/Cairo/p68X8/deps/src/glib-2.34.3 already exists
# │ [ Info: Attempting to create directory /Users/xronos/.julia/packages/Cairo/p68X8/deps/builds/gobject
# │ [ Info: Directory /Users/xronos/.julia/packages/Cairo/p68X8/deps/builds/gobject already exists
# │ [ Info: Changing directory to /Users/xronos/.julia/packages/Cairo/p68X8/deps/builds/gobject
# │ ERROR: LoadError: failed process: Process(setenv(`/Users/xronos/.julia/packages/Cairo/p68X8/deps/src/glib-2.34.3/configure --prefix=/Users/xronos/.julia/packages/Cairo/p68X8/deps/usr`,["XPC_FLAGS=0x0", "PATH=/Users/xronos/.julia/packages/Cairo/p68X8/deps/usr/bin:/Users/xronos/.julia/packages/Homebrew/s09IX/deps/usr/bin:/Users/xronos/.julia/packages/Homebrew/s09IX/deps/usr/sbin:/usr/local/opt/curl/bin:/Users/xronos/opt/bin:/usr/local/cuda/bin:/usr/local/anaconda3/bin:/usr/local/bin:/usr/bin:/bin:/usr/sbin:/sbin:/usr/local/share/dotnet:~/.dotnet/tools:/Library/Frameworks/Mono.framework/Versions/Current/Commands:/Users/xronos/development/adt-bundle-mac-x86_64-20140321/sdk/tools:/Users/xronos/development/adt-bundle-mac-x86_64-20140321/sdk/platform-tools:/Users/xronos/.composer/vendor/bin:/Users/xronos/.julia/packages/Homebrew/s09IX/deps/usr/opt/git/libexec/git-core", "NODE_PATH=/Applications/Atom.app/Contents/Resources/app.asar/exports", "PWD=/", "XPC_SERVICE_NAME=0", "GIT_EXEC_PATH=/Users/xronos/.julia/packages/Homebrew/s09IX/deps/usr/opt/git/libexec/git-core", "SHELL=/bin/bash", "__CF_USER_TEXT_ENCODING=0x1F5:0x0:0x0", "NODE_ENV=production", "TMPDIR=/var/folders/1k/rgrc00jj1j50sbswmmg2f41r0000gn/T/", "HOMEBREW_CACHE=/Users/xronos/Library/Caches/Homebrew.jl/", "CPPFLAGS= -I/Users/xronos/.julia/packages/Cairo/p68X8/deps/usr/include", "LDFLAGS= -L/Users/xronos/.julia/packages/Cairo/p68X8/deps/usr/lib -Wl,-rpath -Wl,/Users/xronos/.julia/packages/Cairo/p68X8/deps/usr/lib", "PS1=\\[\\033[0m\\]\$(git branch &>/dev/null;\\", "SHLVL=1", "CUDA_HOME=/usr/local/cuda", "LOGNAME=xronos", "HOMEBREW_NO_ENV_FILTERING=1", "SSH_AUTH_SOCK=/private/tmp/com.apple.launchd.t1wi7GERbV/Listeners", "ATOM_HOME=/Users/xronos/.atom", "PKG_CONFIG_PATH=/Users/xronos/.julia/packages/Cairo/p68X8/deps/usr/lib/pkgconfig", "JULIA_LOAD_PATH=@:/var/folders/1k/rgrc00jj1j50sbswmmg2f41r0000gn/T/tmpwq880M", "_=/usr/bin/env", "GIT_TEMPLATE_DIR=/Users/xronos/.julia/packages/Homebrew/s09IX/deps/usr/opt/git/share/git-core", "Apple_PubSub_Socket_Render=/private/tmp/com.apple.launchd.5FmVRmq0Rh/Render", "USER=xronos", "JULIA_EDITOR=atom  -a", "HOMEBREW_NO_AUTO_UPDATE=1", "HOME=/Users/xronos", "TERM=xterm-256color", "JULIA_NUM_THREADS=4", "OPENBLAS_MAIN_FREE=1", "HOMEBREW_NO_ANALYTICS=1"]), ProcessExited(1)) [1]
# │ Stacktrace:
# │  [1] error(::String, ::Base.Process, ::String, ::Int64, ::String) at ./error.jl:42
# │  [2] pipeline_error at ./process.jl:785 [inlined]
# │  [3] #run#515(::Bool, ::Function, ::Cmd) at ./process.jl:726
# │  [4] run(::Cmd) at ./process.jl:724
# │  [5] run(::BinDeps.SynchronousStepCollection) at /Users/xronos/.julia/packages/BinDeps/ZEval/src/BinDeps.jl:521
# │  [6] run(::FileRule) at /Users/xronos/.julia/packages/BinDeps/ZEval/src/BinDeps.jl:483
# │  [7] run(::BinDeps.SynchronousStepCollection) at /Users/xronos/.julia/packages/BinDeps/ZEval/src/BinDeps.jl:521 (repeats 2 times)
# │  [8] satisfy!(::BinDeps.LibraryDependency, ::Array{DataType,1}) at ./string:944
# │  [9] satisfy!(::BinDeps.LibraryGroup, ::Array{DataType,1}) at ./string:910
# │  [10] satisfy!(::BinDeps.LibraryGroup) at ./string:874
# │  [11] top-level scope at string:977
# │  [12] include at ./boot.jl:326 [inlined]
# │  [13] include_relative(::Module, ::String) at ./loading.jl:1038
# │  [14] include(::Module, ::String) at ./sysimg.jl:29
# │  [15] include(::String) at ./client.jl:403
# │  [16] top-level scope at none:0
# │ in expression starting at /Users/xronos/.julia/packages/Cairo/p68X8/deps/build.jl:176
# │ ignored: dlclose()
# │ ignored: dlclose()
# │ ignored: dlclose()
# │ ignored: dlclose()
# │ ignored: dlclose()
# │ ignored: dlclose()
# │ ignored: dlclose()
# │ ignored: dlclose()

max_stars_repo_name: ryanorsinger/JuliaPlayground

max_issues_repo_name: ryanorsinger/JuliaPlayground

max_forks_repo_name: ryanorsinger/JuliaPlayground


text:
export sort_pairs, get_corrtype
"""

    sort_pairs(pairs::AbstractArray)

Sort station pairs into their correlation category (auto-correlation, cross-correlation, cross-channel-correlation)

# Arguments
- `pairs::AbstractArray,`    : Array of station pairs, e.g. output of generate_pairs(file["info/stationlist"])

# Output
- `sorted_pairs::Dict{String, Array{String, N}}`    : dictionary containing station pairs grouped by correlation type

"""
function sort_pairs(pairs::AbstractArray)
    # dict containing each correlation category
    # preallocation of dict size assumes all stations have the same number of channels, which may not be true.
    sorted_pairs = Dict{String, Array{String}}("acorr"     => ["", ""],
                                               "xcorr"     => ["", ""],
                                               "xchancorr" => ["", ""])

    # fill dictionary based on detected correlation type
    for stnPair in pairs
        ct = get_corrtype(stnPair)
        sorted_pairs[ct] = hcat(sorted_pairs[ct], stnPair)
    end

    # remove ["",""] used to initialize dict
    sorted_pairs["acorr"]     = sorted_pairs["acorr"][:, 2:end]
    sorted_pairs["xcorr"]     = sorted_pairs["xcorr"][:, 2:end]
    sorted_pairs["xchancorr"] = sorted_pairs["xchancorr"][:, 2:end]

    return sorted_pairs
end

"""

    corrtype(stnPair::Array{String, 1})

Determine correlation type (cross-correlation, auto-correlation, or cross-channel-correlation) using string slicing.
The station string must be: NET.STA.LOC.CHAN
# Arguments
- `stnPair::Array{String, 1},`    : Station pair, e.g. ["BP.SMNB..BP1", "BP.SMNB..BP3"]

# Output
- `corrtype::String`    : correlation type, e.g. "xchancorr"

"""
function get_corrtype(stnPair::Array{String, 1})
    stn1=join(split(stnPair[1],".")[1:3],".")
    stn2=join(split(stnPair[2],".")[1:3],".")
    chan1=split(stnPair[1],".")[4]
    chan2=split(stnPair[2],".")[4]
    # same station, same channel
    if stnPair[1] == stnPair[2]
        ct = "acorr"
    # same station, different channel
    elseif (stn1 == stn2) && (chan1 != chan2) #(stnPair[1][end-3:end] != stnPair[2][end-3:end]) && (stnPair[1][1:end-3] == stnPair[2][1:end-3])
        ct = "xchancorr"
    # different station
    else
        ct = "xcorr"
    end
    return ct
end


text:
@muladd function perform_step!(integrator,
                               cache::Union{ISSEMConstantCache,
                                            ISSEulerHeunConstantCache},
                                            f=integrator.f)
  @unpack t,dt,uprev,u,p = integrator
  @unpack nlsolver = cache
  alg = unwrap_alg(integrator, true)
  theta = alg.theta
  alg.symplectic ? a = dt/2 : a = theta*dt
  OrdinaryDiffEq.markfirststage!(nlsolver)

  # TODO: Stochastic extrapolants?
  u = uprev

  repeat_step = false

  L = integrator.g(uprev,p,t)
  ftmp = integrator.f(uprev,p,t)

  if alg.symplectic
    z = zero(u) # constant extrapolation, justified by ODE IM
  else
    z = dt*ftmp # linear extrapolation
  end
  nlsolver.z = z

  nlsolver.c = a
  if alg.symplectic
    # u = uprev + z then  u = (uprev+u)/2 = (uprev+uprev+z)/2 = uprev + z/2
    #u = uprev + z/2
    tmp = uprev
  else
    tmp = uprev + dt*(1-theta)*ftmp
  end
  nlsolver.tmp = tmp

  z = OrdinaryDiffEq.nlsolve!(nlsolver, integrator, cache, repeat_step)
  OrdinaryDiffEq.nlsolvefail(nlsolver) && return nothing

  if alg.symplectic
    u = tmp + z
  else
    u = tmp + theta*z
  end

  gtmp = L.*integrator.W.dW

  if typeof(cache) <: ISSEulerHeunConstantCache
    utilde = u + gtmp
    gtmp = ((integrator.g(utilde,p,t) + L)/2)*integrator.W.dW
  end

  u += gtmp

  if integrator.opts.adaptive

    if !OrdinaryDiffEq.isnewton(nlsolver)
      is_compos = isa(integrator.alg, StochasticDiffEqCompositeAlgorithm)
      J = OrdinaryDiffEq.calc_J(integrator, nlsolver.cache)
    end
    Ed = dt*(J*ftmp)/2

    if typeof(cache) <: SplitStepEulerConstantCache
      K = @.. uprev + dt * ftmp
      utilde =  @.. K + integrator.sqdt * L
      ggprime = (integrator.g(utilde,p,t).-L)./(integrator.sqdt)
      En = ggprime .* (integrator.W.dW.^2 .- dt)./2
    elseif typeof(cache) <: ISSEulerHeunConstantCache
      utilde = uprev + L*integrator.sqdt
      ggprime = (integrator.g(utilde,p,t).-L)./(integrator.sqdt)
      En = ggprime.*(integrator.W.dW.^2)./2
    end

    resids = calculate_residuals(Ed, En, uprev, u, integrator.opts.abstol,
                                 integrator.opts.reltol, integrator.opts.delta,
                                 integrator.opts.internalnorm,t)
    integrator.EEst = integrator.opts.internalnorm(resids,t)
  end

  integrator.u = u
end

@muladd function perform_step!(integrator, cache::Union{ISSEMCache,
                                                        ISSEulerHeunCache},
                               f=integrator.f)
  @unpack t,dt,uprev,u,p = integrator
  @unpack gtmp,gtmp2,dW_cache,nlsolver = cache
  @unpack z,tmp = nlsolver
  @unpack k,dz = nlsolver.cache # alias to reduce memory
  J = (OrdinaryDiffEq.isnewton(nlsolver) ? nlsolver.cache.J : nothing)
  alg = unwrap_alg(integrator, true)
  alg.symplectic ? a = dt/2 : a = alg.theta*dt
  dW = integrator.W.dW
  mass_matrix = integrator.f.mass_matrix
  theta = alg.theta
  OrdinaryDiffEq.markfirststage!(nlsolver)

  repeat_step = false

  if integrator.success_iter > 0 && !integrator.u_modified && alg.extrapolant == :interpolant
    current_extrapolant!(u,t+dt,integrator)
  elseif alg.extrapolant == :linear
    @.. u = uprev + integrator.fsalfirst*dt
  else # :constant
    copyto!(u,uprev)
  end

  integrator.f(tmp,uprev,p,t)

  if alg.symplectic
    @.. z = zero(eltype(u)) # Justified by ODE solvers, constrant extrapolation when IM
  else
    @.. z = dt*tmp # linear extrapolation
  end

  if alg.symplectic
    #@.. u = uprev + z/2
    @.. tmp = uprev
  else
    #@.. u = uprev + dt*(1-theta)*tmp + theta*z
    @.. tmp = uprev + dt*(1-theta)*tmp
  end
  nlsolver.c = a
  z = OrdinaryDiffEq.nlsolve!(nlsolver, integrator, cache, repeat_step)
  OrdinaryDiffEq.nlsolvefail(nlsolver) && return

  if alg.symplectic
    @.. u = uprev + z
  else
    #@.. u = uprev + dt*(1-theta)*tmp + theta*z
    @.. u = tmp + theta*z
  end

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

  # Handle noise computations

  integrator.g(gtmp,uprev,p,t)

  if is_diagonal_noise(integrator.sol.prob)
    @.. gtmp2 = gtmp*dW
  else
    mul!(gtmp2,gtmp,dW)
  end

  if typeof(cache) <: ISSEulerHeunCache
    gtmp3 = cache.gtmp3
    @.. z = uprev + gtmp2
    integrator.g(gtmp3,z,p,t)
    @.. gtmp = (gtmp3 + gtmp)/2
    if is_diagonal_noise(integrator.sol.prob)
      @.. gtmp2 = gtmp*dW
    else
      mul!(gtmp2,gtmp,dW)
    end
  end

  @.. u += gtmp2

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

  if integrator.opts.adaptive

    if has_Wfact(f)
      # This means the Jacobian was never computed!
      f.jac(J,uprev,p,t)
    end

    mul!(vec(z),J,vec(tmp))
    @.. k = dt*dt*z/2

    # k is Ed
    # dz is En

      if !is_diagonal_noise(integrator.sol.prob)
        g_sized = norm(gtmp,2)
      else
        g_sized = gtmp
      end

      if typeof(cache) <: ISSEMCache
        @.. z = uprev + dt*tmp + integrator.sqdt * g_sized

        if !is_diagonal_noise(integrator.sol.prob)
          integrator.g(gtmp,z,p,t)
          g_sized2 = norm(gtmp,2)
          @.. dW_cache = dW.^2 - dt
          diff_tmp = integrator.opts.internalnorm(dW_cache,t)
          En = (g_sized2-g_sized)/(2integrator.sqdt)*diff_tmp
          @.. dz = En
        else
          integrator.g(gtmp2,z,p,t)
          g_sized2 = gtmp2
          @.. dz = (g_sized2-g_sized)/(2integrator.sqdt)*(dW.^2 - dt)
        end

      elseif typeof(cache) <: ISSEulerHeunCache
        @.. z = uprev + integrator.sqdt * g_sized

        if !is_diagonal_noise(integrator.sol.prob)
          integrator.g(gtmp,z,p,t)
          g_sized2 = norm(gtmp,2)
          @.. dW_cache = dW.^2
          diff_tmp = integrator.opts.internalnorm(dW_cache,t)
          En = (g_sized2-g_sized)/(2integrator.sqdt)*diff_tmp
          @.. dz = En
        else
          integrator.g(gtmp2,z,p,t)
          g_sized2 = gtmp2
          @.. dz = (g_sized2-g_sized)/(2integrator.sqdt)*(dW.^2)
        end

    end

    calculate_residuals!(tmp, k, dz, uprev, u, integrator.opts.abstol,
                         integrator.opts.reltol, integrator.opts.delta,
                         integrator.opts.internalnorm,t)
    integrator.EEst = integrator.opts.internalnorm(tmp,t)

  end
end


text:
struct FFTBasedPoissonSolver{A, G, Λ, S, B, T}
    architecture :: A
            grid :: G
     eigenvalues :: Λ
         storage :: S
          buffer :: B
      transforms :: T
end

function FFTBasedPoissonSolver(arch, grid, planner_flag=FFTW.PATIENT)
    topo = (TX, TY, TZ) =  topology(grid)

    λx = poisson_eigenvalues(grid.Nx, grid.Lx, 1, TX())
    λy = poisson_eigenvalues(grid.Ny, grid.Ly, 2, TY())
    λz = poisson_eigenvalues(grid.Nz, grid.Lz, 3, TZ())

    eigenvalues = (
        λx = arch_array(arch, λx),
        λy = arch_array(arch, λy),
        λz = arch_array(arch, λz)
    )

    storage = arch_array(arch, zeros(complex(eltype(grid)), size(grid)...))

    transforms = plan_transforms(arch, grid, storage, planner_flag)

    # Need buffer for index permutations and transposes.
    buffer_needed = arch isa GPU && Bounded in topo
    buffer = buffer_needed ? similar(storage) : nothing

    return FFTBasedPoissonSolver(arch, grid, eigenvalues, storage, buffer, transforms)
end

"""
    solve_poisson_equation!(solver)

Solves Poisson's equation ∇²ϕ = RHS where `RHS = solver.storage` using periodic or staggered
Neumann boundary conditions as determined by the `solver.grid`. `solver.storage` will be mutated
to contain the solution.
"""
function solve_poisson_equation!(solver::FFTBasedPoissonSolver)
    topo = TX, TY, TZ = topology(solver.grid)
    λx, λy, λz = solver.eigenvalues

    # We can use the same storage for the RHS and the solution ϕ.
    RHS = ϕ = solver.storage

    # Apply forward transforms in order
    [transform!(RHS, solver.buffer) for transform! in solver.transforms.forward]

    # Solve the discrete Poisson equation.
    @. ϕ = -RHS / (λx + λy + λz)

    # Set the volume mean of the solution to be zero.
    # λx[1, 1, 1] + λy[1, 1, 1] + λz[1, 1, 1] = 0 so if ϕ[1, 1, 1] = 0 we get NaNs
    # everywhere after the inverse transform. In eigenspace, ϕ[1, 1, 1] is the
    # "zeroth mode" corresponding to the volume mean of the transform of ϕ, or of ϕ
    # in physical space.
    # Another way of thinking about this: Solutions to Poisson's equation are only
    # unique up to a constant (the global mean of the solution), so we need to pick
    # a constant. ϕ[1, 1, 1] = 0 chooses the constant to be zero so that the solution
    # has zero-mean.
    CUDA.@allowscalar ϕ[1, 1, 1] = 0

    # Apply backward transforms in order
    [transform!(ϕ, solver.buffer) for transform! in solver.transforms.backward]

    return nothing
end


text:
# Testing the polynomial evaluation macro:
let
    p = randn(6)
    x = 0.33
    @test @evalpoly(x, p[1], p[2]) ≈ Psychro.@polyeval(x, p, 2) atol=1e-13
    @test @evalpoly(x, p[1], p[2], p[3]) ≈ Psychro.@polyeval(x, p, 3) atol=1e-13
    @test @evalpoly(x, p[1], p[2], p[3], p[4]) ≈ Psychro.@polyeval(x, p, 4) atol=1e-13
    @test @evalpoly(x, p[1], p[2], p[3], p[4], p[5]) ≈ Psychro.@polyeval(x, p, 5) atol=1e-13
    @test @evalpoly(x, p[1], p[2], p[3], p[4], p[5], p[6]) ≈ Psychro.@polyeval(x, p, 6) atol=1e-13
end


text:
module DataAssimilation

# Write your package code here.
include("AssimilationProblem.jl")
include("EKF.jl")
include("4dvar.jl")

export AssimilationProblem
export timeindex
export update!
export model_forward!
export EKFsolve
export FourDVarSolve
end


text:
mutable struct QRS{T} where T
    P::T
    Q::T
    S::T
    T::T
    form
    noiseLvl::Float32
    validity::Float32
    params::Tuple{Any}
    actions::Tuple{Any}
end
function changeQ(QRS)
end

function changeS(QRS)
end

function newQRS(QRS)
end


text:
module Server
using HttpServer

end


text:
module GridTopologiesTests

using Test
using Gridap.Helpers
using Gridap.Arrays
using Gridap.Fields
using Gridap.ReferenceFEs
using Gridap.Geometry
using Gridap.Geometry: GridTopologyMock

top = GridTopologyMock()
test_grid_topology(top)

@test num_dims(top) == 2
@test num_cell_dims(top) == 2
@test num_point_dims(top) == 2
@test num_faces(top) == 9 + 13 + 5
@test num_faces(top,0) == 9
@test num_faces(top,1) == 13
@test num_faces(top,2) == 5
@test num_cells(top) == 5
@test num_vertices(top) == 9
@test num_edges(top) == 13
@test num_facets(top) == 13
@test get_dimranges(top) == [1:9, 10:22, 23:27]
@test get_dimrange(top,0) == 1:9
@test get_dimrange(top,1) == 10:22
@test get_dimrange(top,2) == 23:27
@test get_offsets(top) == [0, 9, 22]
@test get_offset(top,0) == 0
@test get_offset(top,1) == 9
@test get_offset(top,2) == 22

cell_to_faces = get_cell_faces(top)
@test cell_to_faces == [
  [1,2,4,5,10,11,12,13,23],[2,3,5,14,13,15,24],[3,6,5,16,15,17,25],
  [4,5,7,8,11,18,19,20,26],[5,6,8,9,17,21,20,22,27]]
@test compute_cell_faces(top) == get_cell_faces(top)

@test get_face_vertices(top,0) == get_faces(top,0,0)
@test get_face_vertices(top,1) == get_faces(top,1,0)
@test get_face_vertices(top,2) == get_faces(top,2,0)
@test get_cell_vertices(top) == get_faces(top,2,0)

face_to_vertices = get_face_vertices(top)
@test face_to_vertices == vcat( [ get_face_vertices(top,d) for d in 0:2]... )
@test compute_face_vertices(top) == get_face_vertices(top)

@test is_simplex(top) == false
@test is_n_cube(top) == false
@test is_oriented(top) == false
@test is_regular(top) == true

reffaces, face_to_ftype = compute_reffaces(Polytope{1},top)
@test reffaces == [SEGMENT,]
@test face_to_ftype == [1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1]
@test reffaces == get_reffaces(Polytope{1},top)
@test face_to_ftype == get_face_type(top,1)

@test compute_isboundary_face(top, 0) == Bool[1, 1, 1, 1, 0, 1, 1, 1, 1]
@test compute_isboundary_face(top, 1) == Bool[1, 0, 1, 0, 1, 0, 1, 0, 1, 1, 0, 1, 1]
@test compute_isboundary_face(top, 2) == Bool[1, 1, 1, 1, 1]

@test get_isboundary_face(top, 0) == compute_isboundary_face(top, 0)
@test get_isboundary_face(top, 1) == compute_isboundary_face(top, 1)
@test get_isboundary_face(top, 2) == compute_isboundary_face(top, 2)

@test compute_isboundary_face(top) == Bool[1,1,1,1,0,1,1,1,1,1,0,1,0,1,0,1,0,1,1,0,1,1,1,1,1,1,1]
@test compute_isboundary_face(top) == get_isboundary_face(top)

@test compute_cell_permutations(top,0) == [[1, 1, 1, 1], [1, 1, 1], [1, 1, 1], [1, 1, 1, 1], [1, 1, 1, 1]]
@test compute_cell_permutations(top,1) == [[1, 1, 1, 1], [1, 1, 1], [1, 1, 1], [1, 1, 1, 1], [2, 1, 1, 1]]
@test compute_cell_permutations(top,2) == [[1], [1], [1], [1], [1]]
@test compute_cell_permutations(top) == [
  [1,1,1,1,1,1,1,1,1],[1,1,1,1,1,1,1],[1,1,1,1,1,1,1],[1,1,1,1,1,1,1,1,1],[1,1,1,1,2,1,1,1,1]]

@test get_cell_permutations(top,0) == compute_cell_permutations(top,0)
@test get_cell_permutations(top,1) == compute_cell_permutations(top,1)
@test get_cell_permutations(top,2) == compute_cell_permutations(top,2)
@test get_cell_permutations(top) == compute_cell_permutations(top)

@test get_reffaces(top) == Polytope[VERTEX, SEGMENT, QUAD, TRI]
@test get_face_type(top) == Int8[1,1,1,1,1,1,1,1,1,2,2,2,2,2,2,2,2,2,2,2,2,2,3,4,4,3,3]
@test get_reffaces_offsets(top) == [0, 1, 2]

end # module


text:
@inline function update_batch_edge(mp::MessagePassing, adj, E::Fill{S,2,Axes}, X::CuMatrix, u) where {S,Axes}
    E = fill(E.value, E.axes)
    update_batch_edge(mp, adj, E, X, u)
end

@inline function update_batch_edge(mp::MessagePassing, adj, E::CuMatrix, X::Fill{S,2,Axes}, u) where {S,Axes}
    X = fill(X.value, X.axes)
    update_batch_edge(mp, adj, E, X, u)
end

@inline function update_batch_edge(mp::MessagePassing, adj, E::AbstractMatrix, X::CuMatrix, u)
    E = convert(typeof(X), E)
    update_batch_edge(mp, adj, E, X, u)
end

@inline function update_batch_edge(mp::MessagePassing, adj, E::CuMatrix, X::AbstractMatrix, u)
    X = convert(typeof(E), X)
    update_batch_edge(mp, adj, E, X, u)
end

@inline function update_batch_edge(mp::MessagePassing, adj, E::CuMatrix, X::CuMatrix, u)
    n = size(adj, 1)
    edge_idx = edge_index_table(adj)
    mapreduce(i -> apply_batch_message(mp, i, adj[i], edge_idx, E, X, u), hcat, 1:n)
end

@inline apply_batch_message(mp::MessagePassing, i, js, edge_idx, E::CuMatrix, X::CuMatrix, u) =
    mapreduce(j -> message(mp, _view(X, i), _view(X, j), _view(E, edge_idx[(i,j)])), hcat, js)

@inline function update_batch_vertex(mp::MessagePassing, M::AbstractMatrix, X::CuMatrix, u)
    M = convert(typeof(X), M)
    update_batch_vertex(mp, M, X, u)
end

@inline function update_batch_vertex(mp::MessagePassing, M::CuMatrix, X::AbstractMatrix, u)
    X = convert(typeof(M), X)
    update_batch_vertex(mp, M, X, u)
end

@inline update_batch_vertex(mp::MessagePassing, M::CuMatrix, X::CuMatrix, u) =
    mapreduce(i -> update(mp, _view(M, i), _view(X, i)), hcat, 1:size(X,2))


text:
"""
    AbstractBridge

Represents a bridged constraint or variable in a
[`MathOptInterface.Bridges.AbstractBridgeOptimizer`](@ref). It contains the
indices of the variables and constraints that it has created in the model. These
can be obtained using [`MathOptInterface.NumberOfVariables`](@ref),
[`MathOptInterface.ListOfVariableIndices`](@ref),
[`MathOptInterface.NumberOfConstraints`](@ref) and
[`MathOptInterface.ListOfConstraintIndices`](@ref) using
[`MathOptInterface.get`](@ref) with the bridge in place of the
[`MathOptInterface.ModelLike`](@ref). Attributes of the bridged model such as
[`MathOptInterface.ConstraintDual`](@ref) and
[`MathOptInterface.ConstraintPrimal`](@ref), can be obtained using
[`MathOptInterface.get`](@ref) with the bridge in place of the constraint index.
These calls are used by the
[`MathOptInterface.Bridges.AbstractBridgeOptimizer`](@ref) to communicate with
the bridge so they should be implemented by the bridge.
"""
abstract type AbstractBridge end

"""
    MOI.get(b::AbstractBridge, ::MOI.NumberOfConstraints{F, S}) where {F, S}

The number of constraints of the type `F`-in-`S` created by the bridge `b` in the model.
"""
MOI.get(::AbstractBridge, ::MOI.NumberOfConstraints) = 0

"""
    MOI.get(b::AbstractBridge, ::MOI.ListOfConstraintIndices{F, S}) where {F, S}

A `Vector{ConstraintIndex{F,S}}` with indices of all constraints of
type `F`-in`S` created by the bride `b` in the model (i.e., of length equal to
the value of `NumberOfConstraints{F,S}()`).
"""
function MOI.get(
    ::AbstractBridge,
    ::MOI.ListOfConstraintIndices{F,S},
) where {F,S}
    return MOI.ConstraintIndex{F,S}[]
end

"""
    MOI.supports(model::MOI.ModelLike, attr::MOI.AbstractConstraintAttribute,
                 BT::Type{<:AbstractBridge})

Return a `Bool` indicating whether `BT` supports setting `attr` to `model`.
"""
function MOI.supports(
    ::MOI.ModelLike,
    ::MOI.AbstractConstraintAttribute,
    ::Type{<:AbstractBridge},
)
    return false
end

"""
    function MOI.get(model::MOI.ModelLike, attr::MOI.AbstractConstraintAttribute,
                     bridge::AbstractBridge)

Return the value of the attribute `attr` of the model `model` for the
constraint bridged by `bridge`.
"""
function MOI.get(
    ::MOI.ModelLike,
    attr::MOI.AbstractConstraintAttribute,
    bridge::AbstractBridge,
)
    return throw(
        ArgumentError(
            "Bridge of type `$(typeof(bridge))` does not support accessing the attribute `$attr`.",
        ),
    )
end

function MOI.get(
    model::MOI.ModelLike,
    ::MOI.CanonicalConstraintFunction,
    bridge::AbstractBridge,
)
    return MOI.Utilities.canonical(
        MOI.get(model, MOI.ConstraintFunction(), bridge),
    )
end

"""
    function MOI.set(model::MOI.ModelLike, attr::MOI.AbstractConstraintAttribute,
                     bridge::AbstractBridge, value)

Set the value of the attribute `attr` of the model `model` for the
constraint bridged by `bridge`.
"""
function MOI.set(
    model::MOI.ModelLike,
    attr::MOI.AbstractConstraintAttribute,
    bridge::AbstractBridge,
    value,
)
    if MOI.is_copyable(attr) && !MOI.supports(model, attr, typeof(bridge))
        throw(MOI.UnsupportedAttribute(attr))
    else
        throw(MOI.SetAttributeNotAllowed(attr))
    end
end

"""
    added_constrained_variable_types(BT::Type{<:Variable.AbstractBridge})::Vector{Tuple{DataType}}

Return a list of the types of constrained variables that bridges of concrete
type `BT` add. This is used by the [`LazyBridgeOptimizer`](@ref).
"""
function added_constrained_variable_types end

"""
    added_constraint_types(BT::Type{<:Constraint.AbstractBridge})::Vector{Tuple{DataType, DataType}}

Return a list of the types of constraints that bridges of concrete type `BT`
add. This is used by the [`LazyBridgeOptimizer`](@ref).
"""
function added_constraint_types end

"""
    set_objective_function_type(BT::Type{<:Objective.AbstractBridge})::Type{<:MOI.AbstractScalarFunction}

Return the type of objective function that bridges of concrete type `BT`
set. This is used by the [`LazyBridgeOptimizer`](@ref).
"""
function set_objective_function_type end


text:
module SampleChainsContinuous

using ZigZagBoomerang
using Reexport
@reexport using SampleChains
using NestedTuples
using ElasticArrays
using ConcreteStructs
using TransformVariables
using MappedArrays
using Random
using ZigZagBoomerang

# using SampleChains: chainvec

export ContinuousChain

@concrete struct  ContinuousChain{T} <: AbstractChain{T}
    samples     # :: AbstractVector{T}
    logq        # log-density for distribution the sample was drawn from
    info        # Per-sample metadata, type depends on sampler used
    meta        # Metadata associated with the sample as a whole
    state
    transform
end

# function ContinuousChain(t::TransformVariables.TransformTuple, Q::DynamicHMC.EvaluatedLogDensity, tree_stats, steps)
#     tQq = TransformVariables.transform(t, Q.q)
#     T = typeof(tQq)
#     samples = chainvec(tQq)
#     logq = chainvec(Q.ℓq)
#     info = chainvec(tree_stats)
#     meta = steps
#     transform = t

#     return ContinuousChain{T}(samples, logq, info, meta, Q, transform)
# end

# SampleChains.summarize(ch::ContinuousChain) = summarize(samples(ch))

# function gettransform(chain::ContinuousChain)
#     return getfield(chain, :transform)
# end

# function SampleChains.pushsample!(chain::ContinuousChain, Q::DynamicHMC.EvaluatedLogDensity, tree_stats)
#     push!(samples(chain), transform(gettransform(chain), Q.q))
#     push!(logq(chain), Q.ℓq)
#     push!(info(chain), tree_stats)
# end

# function SampleChains.step!(chain::ContinuousChain)
#     Q, tree_stats = DynamicHMC.mcmc_next_step(getfield(chain, :meta), getfield(chain, :state))
# end

# function SampleChains.initialize!(rng::Random.AbstractRNG, ::Type{ContinuousChain}, ℓ, tr, ad_backend=Val(:ForwardDiff))

#     results = DynamicHMC.mcmc_keep_warmup(
#         rng,
#         ∇P,
#         0;
#         reporter = reporter
#     )

#     steps = DynamicHMC.mcmc_steps(
#         results.sampling_logdensity,
#         results.final_warmup_state,
#     )

#     Q = results.final_warmup_state.Q
#     Q, tree_stats = DynamicHMC.mcmc_next_step(steps, Q)
#     chain = ContinuousChain(tr, Q, tree_stats, steps)
# end

# function SampleChains.initialize!(::Type{ContinuousChain}, ℓ, tr, ad_backend=Val(:ForwardDiff))
#     rng = Random.GLOBAL_RNG
#     return initialize!(rng, ContinuousChain, ℓ, tr, ad_backend)
# end

# function SampleChains.drawsamples!(chain::ContinuousChain, n::Int=1000)
#     @cleanbreak for j in 1:n
#         Q, tree_stats = step!(chain)
#         pushsample!(chain, Q, tree_stats)
#     end
#     return chain
# end

end


text:
function pretty(x, n=8, dx=[1.0, 2.0, 5.0])
    xmin = minimum(x)
    xmax = maximum(x)
    xr = xmax-xmin
    dx0 = xr / n
    logdx = [log10(y) - floor(log10(y))  for y in dx]
    pwr = floor(log10(dx0))
    l0 = log10(dx0) - pwr

    val, idx = findmin(abs.(abs(l0) .- logdx))
    #return val, idx, dx[idx], (l0, pwr)
    dx1 = dx[idx]*10^pwr

    xlim1 = round(xmin/dx1, RoundDown)*dx1
    xlim2 = round(xmax/dx1, RoundUp)*dx1

    return xlim1:dx1:xlim2

end

function solvezero(f, x0, x1, xorder=0.0, EPS=1e-8, MAXITER=100)

    if xorder==0.0
        xorder = 1.0
    end
    y0 = f(x0)
    y1 = f(x1)
    x = 0.0
    for i = 1:MAXITER
        a = (y1-y0) / (x1-x0)
        b = y0 - a*x0
        x = -b/a
        if abs(x-x1) < EPS*abs(xorder)
            return x
        end
        x0 = x1
        y0 = y1
        x1 = x
        y1 = f(x1)
    end

    error("solvezero did not converge")
end


text:
# This file was generated, do not modify it. # hide
using Gnuplot
t = 0:0.001:1
@gp(t, sin.(2π*5*t), "with lines title 'sin' linecolor 'black'")
#save("plt2_ex2.gp") # hide
#save(term="pdf size 5,3", output="plt2_ex2.pdf") # hide
save(term="pngcairo size 600,400", output="plt2_ex2.png") # hide


text:
# Graph
mutable struct Graph
  nN     :: Int
  edges  :: Vector{Vector{Int}}
  big    :: Vector{Bool}
  visits :: Vector{Int}
  tags   :: Dict{String,Int}
  paths  :: Vector{Vector{Int}}
  hasVisitedSmallTwice :: Bool
end

Graph() = Graph(0, [], [], [], Dict{String,Int}(), [], false)

# Adds a new node to the Graph
function addNode!(G::Graph, n::String) :: Int
  haskey(G.tags,n) && return G.tags[n]
  G.nN += 1
  push!(G.edges,[])
  islowercase(n[1]) ? push!(G.big,false) : push!(G.big,true)
  push!(G.visits,0)
  G.tags[n] = G.nN
  return G.tags[n]
end

# Adds a new edge to the Graph
function addEdge!(G::Graph, v::Tuple{String,String})
  idx1 = haskey(G.tags,v[1]) ? G.tags[v[1]] : addNode!(G,v[1])
  idx2 = haskey(G.tags,v[2]) ? G.tags[v[2]] : addNode!(G,v[2])
  push!(G.edges[idx1],idx2)
  push!(G.edges[idx2],idx1)
end

# Recursive function that computes the number of paths that go from 'istart' to 'iend'.
function paths!(G::Graph, istart::Int, iend::Int) :: Int
  (istart == iend) && return 1
  G.visits[istart] += 1
  res = 0
  for e in G.edges[istart]
    if G.visits[e] == 0 || G.big[e]
      res += paths!(G,e,iend)
    end
  end
  G.visits[istart] -= 1
  return res
end

paths!(G::Graph, istart::String, iend::String) = paths!(G,G.tags[istart],G.tags[iend])

# Same as paths, modified for problem 2.
function paths2!(G::Graph, istart::Int, iend::Int) :: Int
  (istart == iend) && return 1
  G.visits[istart] += 1
  res = 0
  for e in G.edges[istart]
    if G.visits[e] == 0 || G.big[e]
      res += paths2!(G,e,iend)
    elseif (G.tags["start"] !== e) && (G.visits[e] == 1) && !(G.hasVisitedSmallTwice)
      G.hasVisitedSmallTwice = true
      res += paths2!(G,e,iend)
      G.hasVisitedSmallTwice = false
    end
  end
  G.visits[istart] -= 1
  return res
end

paths2!(G::Graph, istart::String, iend::String) = paths2!(G,G.tags[istart],G.tags[iend])

# Same as paths, but also saves the actual paths in G.paths
function paths_visual!(G::Graph, istart::Int, iend::Int, p::Vector{Int}) :: Int
  if istart == iend
    push!(G.paths,deepcopy(p))
    return 1
  end

  G.visits[istart] += 1
  res = 0
  for e in G.edges[istart]
    if G.visits[e] == 0 || G.big[e]
      res += paths_visual!(G,e,iend,vcat(p,e))
    end
  end
  G.visits[istart] -= 1
  return res
end

paths_visual!(G::Graph, istart::String, iend::String) = paths_visual!(G,G.tags[istart],G.tags[iend],[G.tags[istart]])

# Prints the paths saved in G.paths
function printPaths(G)
  invMap = Dict(value => key for (key, value) in G.tags)
  for p in G.paths
    println(map(x -> invMap[x], p))
  end
end

parseLine(line::String) = tuple(map(s -> String(s),split(line,'-'))...)

function readInput()
  f = open("./inputs/day12.txt")
  G = Graph()

  edges = map( line -> parseLine(line), eachline(f))
  map(e -> addEdge!(G,e), edges)
  return G
end

println("Problem 1: ")
G = readInput()
res = paths!(G,"start","end")
printPaths(G)
println(res)

println("Problem 2: ")
G = readInput()
res = paths2!(G,"start","end")
println(res)


text:
module Certificate

import MutableArithmetics
const MA = MutableArithmetics
import MultivariatePolynomials
const MP = MultivariatePolynomials
import MultivariateBases
const MB = MultivariateBases
using SemialgebraicSets

using SumOfSquares

include("newton_polytope.jl")

# Certificate              denominator * p = ...                                      domain
# Polya       (1 + x1 + ... + xn)^(2d) * p = ∑ λ_α x.^α                  λ_α ∈ R_+    R_+^n
# Putinar                                p = σ_0 + ∑ σ_i g_i             σ_i SOS      g_i ≥ 0
# Krivine                                p = ∑ λ_{αβ} g.^α (1 .- g).^β   λ_{αβ} ≥ 0   0 ≤ g_i ≤ 1
# Schmüdgen                              p = ∑ σ_α g^α                   σ_α SOS      g_i ≥ 0
# Hilbert                            d * p = σ                           d, σ SOS     R^n

struct PreorderIndex
    value::Int
end
struct IdealIndex
    value::Int
end

const Index = Union{PreorderIndex, IdealIndex}

abstract type Attribute end

# For get
function cone end
function gram_basis end
# FIXME currently, this returns `MB.MonomialBasis` instead of `MB.MonomialBasis{MT, MVT}`
function gram_basis_type end
function reduced_polynomial end
function ideal_certificate end
function preprocessed_domain end

# Only for PreorderIndex
function preorder_indices end
function generator end
function multiplier_basis end
function multiplier_basis_type end

# PreorderCertificate
# IdealCertificate
# FullSpaceCertificate

abstract type AbstractCertificate end

abstract type AbstractPreorderCertificate <: AbstractCertificate end
abstract type AbstractIdealCertificate <: AbstractCertificate end

"""
    struct Putinar{IC <: AbstractIdealCertificate, CT <: SumOfSquares.SOSLikeCone, BT <: MB.AbstractPolynomialBasis} <: AbstractPreorderCertificate
        ideal_certificate::IC
        cone::CT
        basis::Type{BT}
        maxdegree::Int
    end

The `Putinar` certificate ensures the nonnegativity of `p(x)` for all `x` such that
`g_i(x) >= 0` and `h_i(x) = 0` by exhibiting Sum-of-Squares polynomials `σ_i(x)`
such that `p(x) - ∑ σ_i(x) g_i(x)` is guaranteed to be nonnegativity for all `x`
such that `h_i(x) = 0`.
The polynomials `σ_i(x)` are search over `cone` with a basis of type `basis` such that
the degree of `σ_i(x) g_i(x)` does not exceed `maxdegree`.
"""
struct Putinar{IC <: AbstractIdealCertificate, CT <: SumOfSquares.SOSLikeCone, BT <: MB.AbstractPolynomialBasis} <: AbstractPreorderCertificate
    ideal_certificate::IC
    cone::CT
    basis::Type{BT}
    maxdegree::Int
end

cone(certificate::Putinar) = certificate.cone

struct WithVariables{S, V}
    inner::S
    variables::V
end

function MP.variables(v::WithVariables)
    return v.variables
end
function MP.monomials(v::WithVariables)
    return MP.monomials(v.inner)
end

_merge_sorted(a::Vector, ::Tuple{}) = a
function _merge_sorted(a::Vector, b::Vector)
    vars = sort!(vcat(a, b), rev = true)
    unique!(vars)
    return vars
end
_merge_sorted(a::Tuple{}, ::Tuple{}) = a
_merge_sorted(a::Tuple, ::Tuple{}) = a
_merge_sorted(::Tuple{}, b::Tuple) = b
function _merge_sorted(a::Tuple, b::Tuple)
    v = first(a)
    w = first(b)
    if v == w
        return (v, _merge_sorted(Base.tail(a), Base.tail(b))...)
    elseif v > w
        return (v, _merge_sorted(Base.tail(a), b)...)
    else
        return (w, _merge_sorted(a, Base.tail(b))...)
    end
end

_vars(::SemialgebraicSets.FullSpace) = tuple()
_vars(x) = MP.variables(x)

function with_variables(inner, outer)
    return WithVariables(inner, _merge_sorted(_vars(inner), _vars(outer)))
end

function preprocessed_domain(::Putinar, domain::BasicSemialgebraicSet, p)
    return with_variables(domain, p)
end

function preorder_indices(::Putinar, domain::WithVariables)
    return map(PreorderIndex, eachindex(domain.inner.p))
end

function maxdegree_gram_basis(B::Type, variables, maxdegree::Int)
    return MB.maxdegree_basis(B, variables, div(maxdegree, 2))
end
multiplier_maxdegree(maxdegree, q) = maxdegree - MP.maxdegree(q)
function multiplier_basis(certificate::Putinar, index::PreorderIndex, domain::WithVariables)
    q = domain.inner.p[index.value]
    vars = sort!([domain.variables..., MP.variables(q)...], rev = true)
    unique!(vars)
    return maxdegree_gram_basis(certificate.basis, vars, multiplier_maxdegree(certificate.maxdegree, q))
end
function multiplier_basis_type(::Type{Putinar{IC, CT, BT}}) where {IC, CT, BT}
    return BT
end

function generator(::Putinar, index::PreorderIndex, domain::WithVariables)
    return domain.inner.p[index.value]
end

ideal_certificate(certificate::Putinar) = certificate.ideal_certificate
ideal_certificate(::Type{<:Putinar{IC}}) where {IC} = IC

SumOfSquares.matrix_cone_type(::Type{<:Putinar{IC, CT}}) where {IC, CT} = SumOfSquares.matrix_cone_type(CT)

######################
# Ideal certificates #
######################

abstract type SimpleIdealCertificate{CT, BT} <: AbstractIdealCertificate end
reduced_polynomial(::SimpleIdealCertificate, poly, domain) = poly

cone(certificate::SimpleIdealCertificate) = certificate.cone
SumOfSquares.matrix_cone_type(::Type{<:SimpleIdealCertificate{CT}}) where {CT} = SumOfSquares.matrix_cone_type(CT)

# TODO return something else when `PolyJuMP` support other bases.
zero_basis(certificate::SimpleIdealCertificate) = MB.MonomialBasis
zero_basis_type(::Type{<:SimpleIdealCertificate{CT, BT}}) where {CT, BT} = MB.MonomialBasis

"""
    struct MaxDegree{CT <: SumOfSquares.SOSLikeCone, BT <: MB.AbstractPolynomialBasis} <: SimpleIdealCertificate{CT, BT}
        cone::CT
        basis::Type{BT}
        maxdegree::Int
    end

The `MaxDegree` certificate ensures the nonnegativity of `p(x)` for all `x` such that
`h_i(x) = 0` by exhibiting a Sum-of-Squares polynomials `σ(x)`
such that `p(x) - σ(x)` is guaranteed to be zero for all `x`
such that `h_i(x) = 0`.
The polynomial `σ(x)` is search over `cone` with a basis of type `basis` such that
the degree of `σ(x)` does not exceed `maxdegree`.
"""
struct MaxDegree{CT <: SumOfSquares.SOSLikeCone, BT <: MB.AbstractPolynomialBasis} <: SimpleIdealCertificate{CT, BT}
    cone::CT
    basis::Type{BT}
    maxdegree::Int
end
function gram_basis(certificate::MaxDegree, poly) where CT
    return maxdegree_gram_basis(certificate.basis, MP.variables(poly), certificate.maxdegree)
end
function gram_basis_type(::Type{MaxDegree{CT, BT}}) where {CT, BT}
    return BT
end

"""
    struct FixedBasis{CT <: SumOfSquares.SOSLikeCone, BT <: MB.AbstractPolynomialBasis} <: SimpleIdealCertificate{CT, BT}
        cone::CT
        basis::BT
    end

The `FixedBasis` certificate ensures the nonnegativity of `p(x)` for all `x` such that
`h_i(x) = 0` by exhibiting a Sum-of-Squares polynomials `σ(x)`
such that `p(x) - σ(x)` is guaranteed to be zero for all `x`
such that `h_i(x) = 0`.
The polynomial `σ(x)` is search over `cone` with basis `basis`.
"""
struct FixedBasis{CT <: SumOfSquares.SOSLikeCone, BT <: MB.AbstractPolynomialBasis} <: SimpleIdealCertificate{CT, BT}
    cone::CT
    basis::BT
end
function gram_basis(certificate::FixedBasis, poly) where CT
    return certificate.basis
end
function gram_basis_type(::Type{FixedBasis{CT, BT}}) where {CT, BT}
    return BT
end

"""
    struct Newton{CT <: SumOfSquares.SOSLikeCone, BT <: MB.AbstractPolynomialBasis, NPT <: Tuple} <: SimpleIdealCertificate{CT, BT}
        cone::CT
        basis::Type{BT}
        variable_groups::NPT
    end

The `Newton` certificate ensures the nonnegativity of `p(x)` for all `x` such that
`h_i(x) = 0` by exhibiting a Sum-of-Squares polynomials `σ(x)`
such that `p(x) - σ(x)` is guaranteed to be zero for all `x`
such that `h_i(x) = 0`.
The polynomial `σ(x)` is search over `cone` with a basis of type `basis`
chosen using the multipartite Newton polytope with parts `variable_groups`.
If `variable_groups = tuple()` then it falls back to the classical Newton polytope
with all variables in the same part.
"""
struct Newton{CT <: SumOfSquares.SOSLikeCone, BT <: MB.AbstractPolynomialBasis, NPT <: Tuple} <: SimpleIdealCertificate{CT, BT}
    cone::CT
    basis::Type{BT}
    variable_groups::NPT
end
function gram_basis(certificate::Newton{CT, B}, poly) where {CT, B}
    return MB.basis_covering_monomials(B, monomials_half_newton_polytope(MP.monomials(poly), certificate.variable_groups))
end
function gram_basis_type(::Type{<:Newton{CT, BT}}) where {CT, BT}
    return BT
end

"""
    struct Remainder{GCT<:AbstractIdealCertificate} <: AbstractIdealCertificate
        gram_certificate::GCT
    end

The `Remainder` certificate ensures the nonnegativity of `p(x)` for all `x` such that
`h_i(x) = 0` by guaranteeing the remainder of `p(x)` modulo the ideal generated by
`⟨h_i⟩` to be nonnegative for all `x` such that `h_i(x) = 0` using the certificate
`gram_certificate`.
For instance, if `gram_certificate` is [`SumOfSquares.Certificate.Newton`](@ref),
then the certificate `Remainder(gram_certificate)` will take the remainder before
computing the Newton polytope hence might generate a much smaller Newton polytope
hence a smaller basis and smaller semidefinite program.
However, this then corresponds to a lower degree of the hierarchy which might
be insufficient to find a certificate.
"""
struct Remainder{GCT<:AbstractIdealCertificate} <: AbstractIdealCertificate
    gram_certificate::GCT
end

function reduced_polynomial(::Remainder, poly, domain)
    return convert(typeof(poly), rem(poly, ideal(domain)))
end

function gram_basis(certificate::Remainder, poly)
    return gram_basis(certificate.gram_certificate, poly)
end
function gram_basis_type(::Type{Remainder{GCT}}) where GCT
    return gram_basis_type(GCT)
end

cone(certificate::Remainder) = cone(certificate.gram_certificate)
SumOfSquares.matrix_cone_type(::Type{Remainder{GCT}}) where {GCT} = SumOfSquares.matrix_cone_type(GCT)
zero_basis(certificate::Remainder) = zero_basis(certificate.gram_certificate)
zero_basis_type(::Type{Remainder{GCT}}) where {GCT} = zero_basis_type(GCT)

include("Sparsity/Sparsity.jl")
using .Sparsity: SignSymmetry, ChordalCompletion, ClusterCompletion
export Sparsity, SignSymmetry, ChordalCompletion, ClusterCompletion

include("Symmetry/Symmetry.jl")
export Symmetry

end

max_stars_repo_name: blegat/SumOfSquares

max_issues_repo_name: blegat/SumOfSquares

max_forks_repo_name: blegat/SumOfSquares


text:
"""
Data container for collecting vic results.
"""
struct VicRes
    name::String
    regine_main::String
    lat::Array{Float64,1}
    lon::Array{Float64,1}
    gridcel::Array{Int64,1}
    elev::Array{Float64,1}
    data_all::Array{Float64,3}
    data_mean::Array{Float64,2}
    var_names::Array{String,1}
    time::Array{DateTime,1}
end

"""
Data container for collecting vic inputs.
"""
struct VicInput
    name::String
    regine_main::String
    lat::Array{Float64,1}
    lon::Array{Float64,1}
    gridcel::Array{Int64,1}
    elev::Array{Float64,1}
    data_all::Array{Float64,3}
    data_mean::Array{Float64,2}
    var_names::Array{String,1}
    time::Array{DateTime,1}
end

"""
Struct with metadata for watersheds.
"""
struct WatershedData

    name::String
    regine_main::String
    dbk::Float64
    xrange::UnitRange{Int64}
    yrange::UnitRange{Int64}
    gridcel::Array{Int64,2}
    ind_1km::Array{Int64,2}
    ind_5km::Array{Int64,2}
    ind_10km::Array{Int64,2}
    ind_25km::Array{Int64,2}
    ind_50km::Array{Int64,2}
    elev::Array{Float64,2}

end


text:
@testset "ApplyTransform" begin
    func = compile(grad(ApplyTransform(fulltrans, point(2), 2), 2))
    val = ev(func, [0.0, 0.0], squareelt)
    @test val == Matrix(1.0I, 2, 2)
end

@testset "Point" begin
    func = compile(grad(point(2), 2))
    val = ev(func, [0.0, 0.0], squareelt)
    @test val == Matrix(1.0I, 2, 2)
end

@testset "Constant" begin
    Random.seed!(201808081405)
    func = compile(grad(Constant(rand(Float64, 3, 4, 5)), 2))
    val = ev(func, [0.0, 0.0], squareelt)
    @test val == zeros(Float64, 3, 4, 5, 2)
end

@testset "Contract" begin
    Random.seed!(201808081543)

    data = rand(Float64, 3, 2)
    func = grad(Contract(Constant(data), point(2)), 2)
    @test size(func) == (3, 2)
    val = ev(compile(func), [0.0, 0.0], squareelt)
    @test val == data

    data = rand(Float64, 2, 3)
    func = grad(Contract(point(2), Constant(data)), 2)
    @test size(func) == (3, 2)
    val = ev(compile(func), [0.0, 0.0], squareelt)
    @test val == data'

    pt = point(2)
    func = grad(Contract(InsertAxis(pt, 1) .* InsertAxis(pt, 2), pt), 2)
    @test size(func) == (2, 2)
    val = ev(compile(func), [0.3, 0.7], squareelt)
    @test val[:] ≈ [3 * 0.3^2 + 0.7^2, 2 * 0.3 * 0.7, 2 * 0.3 * 0.7, 0.3^2 + 3 * 0.7^2]
end

@testset "Inflate" begin
    index = Constant([3, 5])
    func = grad(Inflate(point(2), (7,), 1 => index), 2)
    @test size(func) == (7, 2)
    val = ev(compile(func), [0.0, 0.0], squareelt)
    @test val == [0.0 0.0; 0.0 0.0; 1.0 0.0; 0.0 0.0; 0.0 1.0; 0.0 0.0; 0.0 0.0]
end

@testset "InsertAxis" begin
    func = compile(grad(InsertAxis(point(2), 1, 2), 2))
    val = ev(func, [0.0, 0.0], squareelt)
    @test size(val) == (1, 2, 1, 2)
    @test val[1,:,1,:] == Matrix(1.0I, 2, 2)
end

@testset "Inv" begin
    pt = point(2)
    func = grad(inv(Add(InsertAxis(pt, 1) .* InsertAxis(pt, 2), Constant(Matrix(1.0I, 2, 2)))), 2)
    @test size(func) == (2, 2, 2)
    val = ev(compile(func), [0.2, 0.9], squareelt)
    det = 1 + 0.2^2 + 0.9^2
    ddet1 = 2 * 0.2
    ddet2 = 2 * 0.9
    @test val[1:8] ≈ [
        -ddet1 * (1 + 0.9^2),
        0.2 * 0.9 * ddet1 - 0.9 * det,
        0.2 * 0.9 * ddet1 - 0.9 * det,
        2 * 0.2 * det - (1 + 0.2^2) * ddet1,
        2 * 0.9 * det - (1 + 0.9^2) * ddet2,
        0.2 * 0.9 * ddet2 - 0.2 * det,
        0.2 * 0.9 * ddet2 - 0.2 * det,
        -ddet2 * (1 + 0.2^2),
    ] / det^2
end

@testset "Monomials" begin
    func = grad(Monomials(point(2), 4), 2)
    @test size(func) == (5, 2, 2)
    val = ev(compile(func), [0.2, 0.5], squareelt)
    @test val[:,1,1] == (0:4) .* 0.2 .^ (-1:3)
    @test val[:,2,1] == fill(0.0, (5,))
    @test val[:,1,2] == fill(0.0, (5,))
    @test val[:,2,2] == (0:4) .* 0.5 .^ (-1:3)
end

@testset "Neg" begin
    func = compile(grad(-point(2), 2))
    val = ev(func, [0.0, 0.0], squareelt)
    @test val == -Matrix(1.0I, 2, 2)
end

@testset "Product" begin
    Random.seed!(201808081416)
    data = rand(Float64, 3, 4, 2)
    func = Constant(data) .* InsertAxis(point(2), 1, 1)
    @test size(func) == (3, 4, 2)
    func = grad(func, 2)
    @test size(func) == (3, 4, 2, 2)

    val = ev(compile(func), [0.0, 0.0], squareelt)
    @test val[:,:,1,1] == data[:,:,1]
    @test val[:,:,2,2] == data[:,:,2]
    @test val[:,:,1,2] == zeros(Float64, 3, 4)
    @test val[:,:,2,1] == zeros(Float64, 3, 4)

    func = compile(grad(InsertAxis(point(2), 1) .* InsertAxis(point(2), 2), 2))
    val = ev(func, [0.0, 0.0], squareelt)
    @test val[:] == fill(0.0, (8,))
    val = ev(func, [0.2, 0.5], squareelt)
    @test val[:] == [0.4, 0.5, 0.5, 0.0, 0.0, 0.2, 0.2, 1.0]
end

@testset "Add" begin
    pt = point(2)
    func = grad(Add(InsertAxis(pt, [1]) .* InsertAxis(pt, [2]), Constant(Matrix(1.0I, 2, 2))), 2)
    @test size(func) == (2, 2, 2)
    func = compile(func)
    val = ev(func, [0.0, 0.0], squareelt)
    @test val[:] == fill(0.0, (8,))
    val = ev(func, [0.2, 0.5], squareelt)
    @test val[:] == [0.4, 0.5, 0.5, 0.0, 0.0, 0.2, 0.2, 1.0]
end

@testset "Zeros" begin
    func = compile(grad(Zeros(Float64, 4, 5), 2))
    val = ev(func, [0.0, 0.0], squareelt)
    @test val == zeros(Float64, 4, 5, 2)
end

max_stars_repo_name: TheBB/Jutils

max_issues_repo_name: TheBB/Jutils

max_forks_repo_name: TheBB/Jutils


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

struct Params
    cache::Vector{InferenceResult}
    world::UInt
    global_cache::Bool

    # optimization
    inlining::Bool
    ipo_constant_propagation::Bool
    aggressive_constant_propagation::Bool
    inline_cost_threshold::Int  # number of CPU cycles beyond which it's not worth inlining
    inline_nonleaf_penalty::Int # penalty for dynamic dispatch
    inline_tupleret_bonus::Int  # extra willingness for non-isbits tuple return types

    # don't consider more than N methods. this trades off between
    # compiler performance and generated code performance.
    # typically, considering many methods means spending lots of time
    # obtaining poor type information.
    # It is important for N to be >= the number of methods in the error()
    # function, so we can still know that error() is always Bottom.
    MAX_METHODS::Int
    # the maximum number of union-tuples to swap / expand
    # before computing the set of matching methods
    MAX_UNION_SPLITTING::Int
    # the maximum number of union-tuples to swap / expand
    # when inferring a call to _apply
    MAX_APPLY_UNION_ENUM::Int

    # parameters limiting large (tuple) types
    MAX_TUPLETYPE_LEN::Int
    TUPLE_COMPLEXITY_LIMIT_DEPTH::Int

    # when attempting to inlining _apply, abort the optimization if the tuple
    # contains more than this many elements
    MAX_TUPLE_SPLAT::Int

    # reasonable defaults
    global function CustomParams(world::UInt,
                    ;
                    inlining::Bool = inlining_enabled(),
                    inline_cost_threshold::Int = DEFAULT_PARAMS.inline_cost_threshold,
                    inline_nonleaf_penalty::Int = DEFAULT_PARAMS.inline_nonleaf_penalty,
                    inline_tupleret_bonus::Int = DEFAULT_PARAMS.inline_tupleret_bonus,
                    max_methods::Int = DEFAULT_PARAMS.MAX_METHODS,
                    tupletype_len::Int = DEFAULT_PARAMS.MAX_TUPLETYPE_LEN,
                    tupletype_depth::Int = DEFAULT_PARAMS.TUPLE_COMPLEXITY_LIMIT_DEPTH,
                    tuple_splat::Int = DEFAULT_PARAMS.MAX_TUPLE_SPLAT,
                    union_splitting::Int = DEFAULT_PARAMS.MAX_UNION_SPLITTING,
                    apply_union_enum::Int = DEFAULT_PARAMS.MAX_APPLY_UNION_ENUM)
        return new(Vector{InferenceResult}(),
                   world, false,
                   inlining, true, false, inline_cost_threshold, inline_nonleaf_penalty,
                   inline_tupleret_bonus, max_methods, union_splitting, apply_union_enum,
                   tupletype_len, tupletype_depth, tuple_splat)
    end
    function Params(world::UInt)
        inlining = inlining_enabled()
        return new(Vector{InferenceResult}(),
                   world, true,
                   #=inlining, ipo_constant_propagation, aggressive_constant_propagation, inline_cost_threshold, inline_nonleaf_penalty,=#
                   inlining, true, false, 100, 1000,
                   #=inline_tupleret_bonus, max_methods, union_splitting, apply_union_enum=#
                   400, 4, 4, 8,
                   #=tupletype_len, tupletype_depth, tuple_splat=#
                   15, 3, 16)
    end
end
const DEFAULT_PARAMS = Params(UInt(0))

max_stars_repo_name: perrutquist/julia

max_issues_repo_name: perrutquist/julia

max_forks_repo_name: perrutquist/julia


text:
# MIT License

# Copyright (c) Microsoft Corporation.

# 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

"""
Either `Rational` or `Euclidean`, used for [`Spline`](@ref)s and [`SplineSurface`](@ref)s.
"""
abstract type CurveType end
abstract type Rational <: CurveType end
abstract type Euclidean <: CurveType end
export CurveType

"""
    Spline{P<:CurveType,S<:Number,N,M}

`M` is the curve order, i.e., the highest power of the parameterizing variable, u.
`P` determines the [`CurveType`](@ref).

All Spline types must implement:
```julia
point(curve,u)
```
and have field `controlpolygon`
"""
abstract type Spline{P<:CurveType,S<:Number,N,M} end

"""
    SplineSurface{P,S,N,M} <: ParametricSurface{S,N}

Curve order, `M`, is the same in the u and v direction and fixed over all spans.
`P` determines the [`CurveType`](@ref).
"""
abstract type SplineSurface{P<:CurveType,S,N,M} <: ParametricSurface{S,N} end
export Spline, SplineSurface

function toeuclidean(point)
    return SVector{length(point) - 1}([point[i] / point[end] for i in 1:(lastindex(point) - 1)])
end

function euclideancontrolpoints(curve::Spline{Rational,S,N,M}) where {S,N,M}
    return toeuclidean.(curve.controlpolygon)
end

function euclideancontrolpoints(curve::Spline{Euclidean,S,N,M}) where {S,N,M}
    return curve.controlpolygon
end

#Bernstein polynomial functions
#B(i,n,u) is broken into Bc,Bu because these two functions will be needed when constructing the matrics to compute moving lines or moving planes
Bc(i, n) = factorial(n) ÷ (factorial(i) * factorial(n - i)) #coefficient of the i,n Bernstein basis polynomial
Bu(i, n, u) = u^i * (1 - u)^(n - i) #polynomial part of Bernstein basis polynomial
B(i, n, u) = Bc(i, n) * Bu(i, n, u)


text:
#InspectDR: IO functionnality for writing plots with Cairo
#-------------------------------------------------------------------------------

#= Comments

 - `show` interface: writing to an IO.
 - `write` interface: writing to a file.
=#

#==Constants
===============================================================================#
const MIMEpng = MIME"image/png"
const MIMEsvg = MIME"image/svg+xml"
const MIMEeps = MIME"image/eps"
const MIMEeps2 = MIME"application/eps" #Apparently this is also a MIME
const MIMEps = MIME"application/postscript" #TODO: support when Cairo.jl supports PSSurface
const MIMEpdf = MIME"application/pdf"
#If an easy way to read Cairo scripts back to a surface is found, we could add:
#const MIMEcairo = MIME"image/cairo"

const MAPEXT2MIME = Dict{String,MIME}(
	".png" => MIMEpng(),
	".svg" => MIMEsvg(),
	".eps" => MIMEeps(),
	".pdf" => MIMEpdf(),
)

#All supported MIMEs:
#EXCLUDE SVG so it can be turnd on/off??
const MIMEall = Union{MIMEpng, MIMEeps, MIMEeps2, MIMEpdf, MIMEsvg}

#=="Constructors"
===============================================================================#
_CairoSurface(io::IO, ::MIMEsvg, w::Float64, h::Float64) =
	Cairo.CairoSVGSurface(io, w, h)
_CairoSurface(io::IO, ::MIMEeps, w::Float64, h::Float64) =
	Cairo.CairoEPSSurface(io, w, h)
_CairoSurface(io::IO, ::MIMEeps2, w::Float64, h::Float64) =
	Cairo.CairoEPSSurface(io, w, h)
_CairoSurface(io::IO, ::MIMEpdf, w::Float64, h::Float64) =
	Cairo.CairoPDFSurface(io, w, h)

#=="withsurf" interface: Make uniform API for all output surfaces
===============================================================================#
function withsurf(fn::Function, stream::IO, mime::MIME, w::Float64, h::Float64)
	surf = _CairoSurface(stream, mime, w, h)
	ctx = CairoContext(surf)
	fn(ctx)
	Cairo.destroy(ctx)
	Cairo.destroy(surf)
end

#There is no PNG stream surface... must write from ARGB surface.
function withsurf(fn::Function, stream::IO, mime::MIMEpng, w::Float64, h::Float64)
	w = round(w); h = round(h)
	surf = Cairo.CairoARGBSurface(w, h)
	ctx = CairoContext(surf)
	fn(ctx)
	Cairo.destroy(ctx)
	Cairo.write_to_png(surf, stream)
	Cairo.destroy(surf)
end

#==MIME interface
===============================================================================#

#w, h: w/h of entire figure.
function _show(stream::IO, mime::MIME, mplot::Multiplot, w::Float64, h::Float64)
	refresh!(mplot.layout) #Copy in values from defaults
	layout = mplot.layout.values
	yoffset = layout.valloc_title
	nrows, ncols = griddims_auto(mplot)
	wplot = w/ncols; hplot = (h-yoffset)/nrows

	withsurf(stream, mime, w, h) do ctx
		_reset(ctx)
		bb = BoundingBox(0,w,0,h)
		Cairo.save(ctx)
			drawrectangle(ctx, bb, layout.frame)
		Cairo.restore(ctx)
		render(ctx, mplot.title, Point2D(w/2, yoffset/2),
			layout.font_title, align=ALIGN_HCENTER|ALIGN_VCENTER
		)

		for (i, plot) in enumerate(mplot.subplots)
			bb = plot.plotbb
			if nothing == bb
				#Use auto-computed grid layout:
				row = div(i-1, ncols) + 1
				col = i - (row-1)*ncols
				xmin = (col-1)*wplot; ymin = yoffset+(row-1)*hplot
				bb = BoundingBox(xmin, xmin+wplot, ymin, ymin+hplot)
			end

			Cairo.save(ctx) #-----
			setclip(ctx, bb) #Avoid accidental overwrites.
			render(ctx, plot, bb)
			Cairo.restore(ctx) #-----
		end
	end
end

#_show() Plot: Leverage write to Multiplot
function _show(stream::IO, mime::MIME, plot::Plot, w::Float64, h::Float64)
	mplot = Multiplot()
	mplot.layout[:valloc_title] = 0
	push!(mplot.subplots, plot)
	_show(stream, mime, mplot, w, h)
end

#_show() Plot2D: Auto-coumpute w/h
function _show(stream::IO, mime::MIME, plot::Plot2D)
	refresh!(plot.layout) #Copy in values from defaults
	bb = plotbounds(plot.layout.values, grid1(plot))
	_show(stream, mime, plot, bb.xmax, bb.ymax)
end

#Default show/MIME behaviour: MethodError
Base.show(io::IO, mime::MIME, mplot::Multiplot) =
	throw(MethodError(show, (io, mime, mplot)))
Base.show(io::IO, mime::MIME, plot::Plot) =
	throw(MethodError(show, (io, mime, plot)))

#Maintain text/plain MIME support (Is this ok?).
Base.show(io::IO, ::MIME"text/plain", plot::Plot) = Base.show(io, plot)
Base.show(io::IO, ::MIME"text/plain", mplot::Multiplot) = Base.show(io, mplot)

#show() Plot/Multiplot: Supported MIMEs:
Base.show(io::IO, mime::MIMEall, mplot::Multiplot) =
	_show(io, mime, mplot, size_auto(mplot)...)
Base.show(io::IO, mime::MIMEall, plot::Plot) =
	_show(io, mime, plot)

#=="showable" interface
===============================================================================#
Base.showable(mime::MIME"text/plain", mplot::Multiplot) = true
Base.showable(mime::MIME, mplot::Multiplot) = false #Default
Base.showable(mime::MIMEall, mplot::Multiplot) = true #Supported
Base.showable(mime::MIMEsvg, mplot::Multiplot) = defaults.rendersvg #depends

Base.showable(mime::MIME"text/plain", p::Plot) = true
Base.showable(mime::MIME, p::Plot) = Base.showable(mime::MIME, Multiplot())

#=="write" interface
===============================================================================#

#_write() Multiplot:
function _write(path::String, mime::MIME, mplot::Multiplot, w::Float64, h::Float64)
	open(path, "w") do io
		_show(io, mime, mplot, w, h)
	end
end

#_write() Plot:
function _write(path::String, mime::MIME, plot::Plot, w::Float64, h::Float64)
	open(path, "w") do io
		_show(io, mime, plot, w, h)
	end
end

#_write() Multiplot:
function _write(path::String, mime::MIME, mplot::Multiplot)
	open(path, "w") do io
		show(io, mime, mplot)
	end
end

#_write() Plot2D:
function _write(path::String, mime::MIME, plot::Plot2D)
	open(path, "w") do io
		_show(io, mime, plot)
	end
end

#==Non-MIME write interface (convenience functions)
===============================================================================#

write_png(path::String, mplot::Multiplot, args...) = _write(path, MIMEpng(), mplot, args...)
write_svg(path::String, mplot::Multiplot, args...) = _write(path, MIMEsvg(), mplot, args...)
write_eps(path::String, mplot::Multiplot, args...) = _write(path, MIMEeps(), mplot, args...)
write_pdf(path::String, mplot::Multiplot, args...) = _write(path, MIMEpdf(), mplot, args...)

write_png(path::String, plot::Plot, args...) = _write(path, MIMEpng(), plot, args...)
write_svg(path::String, plot::Plot, args...) = _write(path, MIMEsvg(), plot, args...)
write_eps(path::String, plot::Plot, args...) = _write(path, MIMEeps(), plot, args...)
write_pdf(path::String, plot::Plot, args...) = _write(path, MIMEpdf(), plot, args...)

#Last line


text:
# ---
# title: 23. Merge k Sorted Lists
# id: problem23
# author: Indigo
# date: 2021-04-14
# difficulty: Hard
# categories: Linked List, Divide and Conquer, Heap
# link: <https://leetcode.com/problems/merge-k-sorted-lists/description/>
# hidden: true
# ---
#
# You are given an array of `k` linked-lists `lists`, each linked-list is sorted
# in ascending order.
#
# _Merge all the linked-lists into one sorted linked-list and return it._
#
#
#
# **Example 1:**
#
#
#
#     Input: lists = [[1,4,5],[1,3,4],[2,6]]
#     Output: [1,1,2,3,4,4,5,6]
#     Explanation: The linked-lists are:
#     [
#       1->4->5,
#       1->3->4,
#       2->6
#     ]
#     merging them into one sorted list:
#     1->1->2->3->4->4->5->6
#
#
# **Example 2:**
#
#
#
#     Input: lists = []
#     Output: []
#
#
# **Example 3:**
#
#
#
#     Input: lists = [[]]
#     Output: []
#
#
#
#
# **Constraints:**
#
#   * `k == lists.length`
#   * `0 <= k <= 10^4`
#   * `0 <= lists[i].length <= 500`
#   * `-10^4 <= lists[i][j] <= 10^4`
#   * `lists[i]` is sorted in **ascending order**.
#   * The sum of `lists[i].length` won't exceed `10^4`.
#
#
## @lc code=start
using LeetCode

Base.isless(l1::ListNode, l2::ListNode) = l1.val < l2.val

function merge_k_lists(lists::Vector{T}) where {T<:Union{ListNode{Int},Nothing}}
    pq = BinaryMinHeap{ListNode{Int}}()
    for l in lists
        isnothing(l) || push!(pq, l)
    end
    res = ListNode{Int}()
    p = res
    while !isempty(pq)
        tp_list = pop!(pq)
        p = p.next = tp_list
        tp_list = tp_list.next
        isnothing(tp_list) || push!(pq, tp_list)
    end
    return res.next
end
## @lc code=end


text:
if basename(pwd()) == "aoc"
    cd("2016/16")
end

bitvec(str::AbstractString) = [c == '1' for c in str]
bitvec(xs::Vector{Int}) = Bool.(xs)

dragonrepeat(str::AbstractString, len::Int) = join(Int.(dragonrepeat(bitvec(str), len)))

function dragonrepeat(data::AbstractVector{Bool}, len::Int)
    while length(data) < len
        append!(data, [0], broadcast(!, reverse(data)))
    end
    data[begin:len]
end

@assert dragonrepeat("1", 3) == "100"
@assert dragonrepeat("0", 3) == "001"
@assert dragonrepeat("11111", 11) == "11111000000"
@assert dragonrepeat("111100001010", 25) == "1111000010100101011110000"

function dragonchecksum(data::AbstractVector{Bool})
    while iseven(length(data))
        data = collect(!(a ⊻ b) for (a, b) in Iterators.partition(data, 2))
    end
    data
end

part1(input::AbstractString, len) = join(Int.(part1(bitvec(input), len)))
part1(input::AbstractVector{Bool}, len) = dragonchecksum(dragonrepeat(input, len))

@assert part1("10000", 20) == "01100"
@time part1("10001001100000001", 272)
@time part1("10001001100000001", 35651584)

max_stars_repo_name: stellartux/advent-of-code

max_issues_repo_name: stellartux/advent-of-code

max_forks_repo_name: stellartux/advent-of-code


text:
for (i, j, k) in zip(('a', 'b', 'c'), ('A', 'B', 'C'), (1, 2, 3))
    println(i, j, k)
end

max_stars_repo_name: ethansaxenian/RosettaDecode

max_issues_repo_name: ethansaxenian/RosettaDecode

max_forks_repo_name: ethansaxenian/RosettaDecode


text:
"pad an integer to get to a 00 format"
"""
:Q: what are, in julia, the pros and cons of declaring this utility
function outside the humanreadable function or nesting it inside it?
"""
function nn(n)
    return lpad(string(n), 2, '0')
end

"convert a number of seconds into a HH:MM:SS string"
function humanreadable(seconds)
    hrs = seconds ÷ 3600;
    mins = (seconds % 3600) ÷ 60;
    secs = seconds % 60;
    return "$(nn(hrs)):$(nn(mins)):$(nn(secs))";
end

max_stars_repo_name: widged/Julia-katas

max_issues_repo_name: widged/Julia-katas

max_forks_repo_name: widged/Julia-katas


text:
# ---
# title: 796. Rotate String
# id: problem796
# author: Tian Jun
# date: 2020-10-31
# difficulty: Easy
# categories:
# link: <https://leetcode.com/problems/rotate-string/description/>
# hidden: true
# ---
#
# We are given two strings, `A` and `B`.
#
# A _shift on`A`_ consists of taking string `A` and moving the leftmost
# character to the rightmost position. For example, if `A = 'abcde'`, then it
# will be `'bcdea'` after one shift on `A`. Return `True` if and only if `A` can
# become `B` after some number of shifts on `A`.
#
#
#
#     **Example 1:**
#     Input: A = 'abcde', B = 'cdeab'
#     Output: true
#
#     **Example 2:**
#     Input: A = 'abcde', B = 'abced'
#     Output: false
#
#
# **Note:**
#
#   * `A` and `B` will have length at most `100`.
#
#
## @lc code=start
using LeetCode

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


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

name = "SuperLU_DIST"
version = v"8.0.0"

# Collection of sources required to complete build
sources = [
    GitSource("https://github.com/xiaoyeli/superlu_dist.git", "4459a89719b982dd47af5bb2494a011bb267195f"),
]

# Bash recipe for building across all platforms
script = raw"""
cd $WORKSPACE/srcdir/superlu_dist*
# allow us to set the name of the shared lib.
sed -i -e 's!OUTPUT_NAME superlu_dist!OUTPUT_NAME "${SUPERLU_OUTPUT_NAME}"!g' SRC/CMakeLists.txt

mkdir build && cd build
if [[ "${target}" == *-mingw* ]]; then
    # This is required to ensure that MSMPI can be found by cmake
    export LDFLAGS="-L${libdir} -lmsmpi"
    PLATFLAGS="-DTPL_ENABLE_PARMETISLIB:BOOL=FALSE -DMPI_C_ADDITIONAL_INCLUDE_DIRS=${includedir}"
else
    PLATFLAGS="-DTPL_PARMETIS_INCLUDE_DIRS=${includedir} -DTPL_PARMETIS_LIBRARIES=${libdir}/libparmetis.${dlext};${libdir}/libmetis.${dlext}"
fi

mkdir ${libdir}/superlu_dist

build_superlu_dist()
{
    if [[ "${1}" == "Int64" ]]; then
        INT=64
    else
        INT=32
    fi
    SUPERLU_PREFIX=${libdir}/superlu_dist/Int${INT}
    mkdir ${SUPERLU_PREFIX}
    cmake -DCMAKE_INSTALL_PREFIX=${SUPERLU_PREFIX} \
    -DCMAKE_TOOLCHAIN_FILE=${CMAKE_TARGET_TOOLCHAIN} \
    -DCMAKE_BUILD_TYPE=Release \
    -DBUILD_SHARED_LIBS=ON \
    -DBUILD_STATIC_LIBS=OFF \
    -DTPL_ENABLE_INTERNAL_BLASLIB=OFF \
    -Denable_tests=OFF \
    -Denable_doc=OFF \
    -Denable_single=ON \
    -Denable_double=ON \
    -Denable_complex16=ON \
    -DTPL_BLAS_LIBRARIES="${libdir}/libopenblas.${dlext}" \
    ${PLATFLAGS} \
    -DCMAKE_C_FLAGS="-std=c99" \
    -DXSDK_INDEX_SIZE=${INT} \
    -DXSDK_ENABLE_Fortran=OFF \
    -DSUPERLU_OUTPUT_NAME="superlu_dist_Int${INT}" \
    -Denable_examples=OFF \
    ..
    make -j${nproc}
    make install
}
build_superlu_dist Int32
build_superlu_dist Int64

"""

# These are the platforms we will build for by default, unless further
# platforms are passed in on the command line
# Excluding Windows due to use of `getline` function, which is non-standard and not provided by MinGW
# per Mose. Will return to it later and attempt to find a solution.
platforms = supported_platforms()
# The products that we will ensure are always built
products = [
    LibraryProduct("libsuperlu_dist_Int32", :libsuperlu_dist_Int32, ["\$libdir/superlu_dist/Int32/lib", "\$libdir/superlu_dist/Int32/bin"]),
    LibraryProduct("libsuperlu_dist_Int64", :libsuperlu_dist_Int64, ["\$libdir/superlu_dist/Int64/lib", "\$libdir/superlu_dist/Int64/bin"])
]

# Dependencies that must be installed before this package can be built
dependencies = [
    Dependency(PackageSpec(name="OpenBLAS32_jll", uuid="656ef2d0-ae68-5445-9ca0-591084a874a2")),
    Dependency(PackageSpec(name="PARMETIS_jll", uuid="b247a4be-ddc1-5759-8008-7e02fe3dbdaa"); platforms=filter(!Sys.iswindows, platforms)),
    Dependency("MPICH_jll"; platforms=filter(!Sys.iswindows, platforms)),
    Dependency("MicrosoftMPI_jll"; platforms=filter(Sys.iswindows, platforms)),
    Dependency("METIS_jll"),
    # For OpenMP we use libomp from `LLVMOpenMP_jll` where we use LLVM as compiler (BSD
    # systems), and libgomp from `CompilerSupportLibraries_jll` everywhere else.
    Dependency(PackageSpec(name="CompilerSupportLibraries_jll", uuid="e66e0078-7015-5450-92f7-15fbd957f2ae"); platforms=filter(!Sys.isbsd, platforms)),
    Dependency(PackageSpec(name="LLVMOpenMP_jll", uuid="1d63c593-3942-5779-bab2-d838dc0a180e"); platforms=filter(Sys.isbsd, platforms)),
]

# Build the tarballs, and possibly a `build.jl` as well.
build_tarballs(ARGS, name, version, sources, script, platforms, products, dependencies; julia_compat="1.6", preferred_gcc_version = v"7")


text:
#=

  Test of of global constraints atmost, atleast, and exactly in Julia ConstraintsSolver.jl

  This Picat model was created by Hakan Kjellerstrand, hakank@gmail.com
  See also my Julia page: http://www.hakank.org/julia/

=#
using ConstraintSolver, JuMP
using Cbc, GLPK, Ipopt
const CS = ConstraintSolver
include("constraints_utils.jl")

# The constraints atmost, atleast, and exactly
# are defined in constraints_utils.jl

function atmost_test(n=4,print_solutions=true,all_solutions=true)

    cbc_optimizer = optimizer_with_attributes(Cbc.Optimizer, "logLevel" => 0)
    glpk_optimizer = optimizer_with_attributes(GLPK.Optimizer)
    ipopt_optimizer = optimizer_with_attributes(Ipopt.Optimizer)

    model = Model(optimizer_with_attributes(CS.Optimizer,   "all_solutions"=> all_solutions,
                                                            # "all_optimal_solutions"=>true,
                                                            "logging"=>[],

                                                            "traverse_strategy"=>:BFS,
                                                            # "traverse_strategy"=>:DFS,
                                                            # "traverse_strategy"=>:DBFS,

                                                            # "branch_split"=>:Smallest,
                                                            # "branch_split"=>:Biggest,
                                                            "branch_split"=>:InHalf,

                                                            # https://wikunia.github.io/ConstraintSolver.jl/stable/options/#branch_strategy-(:Auto)
                                                            # "branch_strategy" => :IMPS, # default
                                                            "branch_strategy" => :ABS, # Activity Based Search
                                                            "activity.decay" => 0.999, # default 0.999
                                                            "activity.max_probes" => 100, # default, 10
                                                            "activity.max_confidence_deviation" => 20, # default 20

                                                            # "simplify"=>false,
                                                            # "simplify"=>true, # default

                                                            "time_limit"=>6,

                                                            # "backtrack" => false, # default true
                                                            # "backtrack_sorting" => false, # default true

                                                            # "lp_optimizer" => cbc_optimizer,
                                                            # "lp_optimizer" => glpk_optimizer,
                                                            # "lp_optimizer" => ipopt_optimizer,
                                        ))

    @variable(model, 1 <= x[1:n] <= n, Int)

    # Here we require that the must be
    #  - max 2 1s
    #  - min 2 2s
    #  - exactly 1 3
    atmost(model,x,1,2)
    atmost(model,x,n,1)
    atleast(model,x,2,2)
    exactly(model,x,3,1)

    # Solve the problem
    optimize!(model)

    status = JuMP.termination_status(model)
    # println("status:$status")
    if status == MOI.OPTIMAL
        num_sols = MOI.get(model, MOI.ResultCount())
        println("num_sols:$num_sols\n")
        if print_solutions
            for sol in 1:num_sols
                # println("solution #$sol")
                xx = convert.(Integer,JuMP.value.(x; result=sol))
                println("x:$xx")

            end
        end
    else
        println("status:$status")
    end

    return status
end

@time atmost_test(5)

max_stars_repo_name: Wikunia/hakank

max_issues_repo_name: Wikunia/hakank

max_forks_repo_name: Wikunia/hakank