text:
# TODO this should contain helper functions for testing
using InteractiveUtils: subtypes

const test_vertex_types = subtypes(Signed) ∪ subtypes(Unsigned)
#
# TODO needs more different types, also tuples and named tuples, unitful units, Union types, Nothing
const test_edgeval_types = [Rational{Int}, Float16, Float64, BigFloat, Int, UInt8, Tuple{Int, Int}, NamedTuple{(:a, :b), Tuple{Int, Int}}]

function make_testgraphs(G::Type{<:AbstractGraph}; kwargs...)
    return Channel(c -> _make_testgraphs(c, G, kwargs...))
end

# Similar to `rand` but for more types.
# We do not want type piracy, so we need this.
function rand_sample(::Type{<:Rational{T}}) where {T}
    # TODO do we want to avoid division by zero?
    return rand(T) // rand(T)
end

rand_sample(T::Type) = rand(T)

rand_sample(::Type{Nothing}) = nothing
rand_sample(::Type{Missing}) = missing

function rand_sample(U::Union)
    i = rand(1:length(propertynames(U)))
    return rand_sample(getproperty(U, propertynames(U)[i]))
end

function rand_sample(T::Type, dims...)
    result = Array{T}(dims)
    for i in eachindex(result)
        result[i] = rand_sample(T)
    end
    return result
end

function rand_sample(T::Type{<:Union{Tuple, NamedTuple}})
    return T(rand_sample(TT) for TT in T.types)
end

function _make_testgraphs(c::Channel, ::Type{SimpleGraph}; kwargs...)
    for V in test_vertex_types
        g = SimpleGraph{V}(0)
        info = """type: $(typeof(g))
                  0-graph"""
        put!(c, (graph=g, info=info))
    end
    for V in test_vertex_types
        g = SimpleGraph{V}(1)
        info = """type: $(typeof(g))
                  1-graph"""
        put!(c, (graph=g, info=info))
    end
    for V in test_vertex_types
        g = SimpleGraph{V}(1)
        add_edge!(g, 1, 1)
        info = """type: $(typeof(g))
                  1-graph with self-loop"""
        put!(c, (graph=g, info=info))
    end
    for V in test_vertex_types
        g = SimpleGraph{V}(2)
        info = """type: $(typeof(g))
                  2-vertex-graph, no edges"""
        put!(c, (graph=g, info=info))
    end
    for V in test_vertex_types
        g = SimpleGraph{V}(2)
        add_edge!(g, 1, 2)
        info = """type: $(typeof(g))
                  2-vertex-graph, 1 edge"""
        put!(c, (graph=g, info=info))
    end
    for V in test_vertex_types
        g = SimpleGraph{V}(2)
        add_edge!(g, 2, 2)
        info = """type: $(typeof(g))
                  2-vertex-graph, no edge, 1 self-loop"""
        put!(c, (graph=g, info=info))
    end
    for V in test_vertex_types
        g = SimpleGraph{V}(2)
        add_edge!(g, 2, 2)
        info = """type: $(typeof(g))
                  2-vertex-graph, no edge, 2 self-loops"""
        put!(c, (graph=g, info=info))
    end
    for V in test_vertex_types
        g = SimpleGraph{V}(2)
        add_edge!(g, 1, 2)
        add_edge!(g, 1, 1)
        info = """type: $(typeof(g))
                  2-vertex-graph, 1 edge, 1 self-loop"""
        put!(c, (graph=g, info=info))
    end
    for V in test_vertex_types
        g = SimpleGraph{V}(2)
        add_edge!(g, 1, 2)
        add_edge!(g, 1, 1)
        add_edge!(g, 2, 2)
        info = """type: $(typeof(g))
                  2-vertex-graph, 1 edge, 2 self-loops"""
        put!(c, (graph=g, info=info))
    end
end

function _make_testgraphs(c::Channel, ::Type{SimpleDiGraph}; kwargs...)
    for V in test_vertex_types
        g = SimpleDiGraph{V}(0)
        info = """type: $(typeof(g))
                  0-graph"""
        put!(c, (graph=g, info=info))
    end
    for V in test_vertex_types
        g = SimpleDiGraph{V}(1)
        info = """type: $(typeof(g))
                  1-graph"""
        put!(c, (graph=g, info=info))
    end
    for V in test_vertex_types
        g = SimpleDiGraph{V}(1)
        add_edge!(g, 1, 1)
        info = """type: $(typeof(g))
                  1-graph with self-loop"""
        put!(c, (graph=g, info=info))
    end
    for V in test_vertex_types
        g = SimpleDiGraph{V}(2)
        info = """type: $(typeof(g))
                  2-vertex-graph, no edges"""
        put!(c, (graph=g, info=info))
    end
    for V in test_vertex_types
        g = SimpleDiGraph{V}(2)
        add_edge!(g, 1, 2)
        info = """type: $(typeof(g))
                  2-vertex-graph, 1 edge"""
        put!(c, (graph=g, info=info))
    end
    for V in test_vertex_types
        g = SimpleDiGraph{V}(2)
        add_edge!(g, 2, 2)
        info = """type: $(typeof(g))
                  2-vertex-graph, 1 self loop"""
        put!(c, (graph=g, info=info))
    end
    for V in test_vertex_types
        g = SimpleDiGraph{V}(2)
        add_edge!(g, 1, 1)
        add_edge!(g, 1, 2)
        info = """type: $(typeof(g))
                  2-vertex-graph, 1 edge, 1 self loop"""
        put!(c, (graph=g, info=info))
    end
    for V in test_vertex_types
        g = SimpleDiGraph{V}(2)
        add_edge!(g, 1, 1)
        add_edge!(g, 2, 2)
        info = """type: $(typeof(g))
                  2-vertex-graph, 2 self loops"""
        put!(c, (graph=g, info=info))
    end
    for V in test_vertex_types
        g = SimpleDiGraph{V}(2)
        add_edge!(g, 1, 2)
        add_edge!(g, 2, 1)
        info = """type: $(typeof(g))
                  2-vertex-graph, 1 edge, 1 reverse edge"""
        put!(c, (graph=g, info=info))
    end
    for V in test_vertex_types
        g = SimpleDiGraph{V}(2)
        add_edge!(g, 1, 2)
        add_edge!(g, 2, 1)
        add_edge!(g, 2, 2)
        info = """type: $(typeof(g))
                  2-vertex-graph, 1 edge, 1 reverse edge, 1 self loop"""
        put!(c, (graph=g, info=info))
    end
    for V in test_vertex_types
        g = SimpleDiGraph{V}(2)
        add_edge!(g, 1, 2)
        add_edge!(g, 1, 1)
        add_edge!(g, 2, 2)
        info = """type: $(typeof(g))
                  2-vertex-graph, 1 edge, 2 self loops"""
        put!(c, (graph=g, info=info))
    end
    for V in test_vertex_types
        g = SimpleDiGraph{V}(2)
        add_edge!(g, 1, 2)
        add_edge!(g, 2, 1)
        add_edge!(g, 1, 1)
        add_edge!(g, 2, 2)
        info = """type: $(typeof(g))
                  2-vertex-graph, 1 edge, 1 reverse edge, 2 self loops"""
        put!(c, (graph=g, info=info))
    end
end

function testgraphA()
    g = SimpleValueGraph(CompleteBipartiteGraph(3, 4))
    add_edge!(g, 1, 4, 10.0)
    return g
end

# TODO only works for undirected graphs
# TODO allow for only non-negative or positive weights
function graph_with_randvals(g::AbstractGraph{T}) where {T}
    U = Float64
    resultg = SimpleValueGraph{T, U}(nv(g))
    for e in edges(g)
        s, d = src(e), dst(e)
        add_edge!(resultg, s, d, convert(U, randn()))
    end
    return resultg
end


text: Dict{String, Any}("a" => Dict{String, Any}("b" => Dict{String, Any}("value" => "false", "type" => "bool"), "a" => Dict{String, Any}("value" => "true", "type" => "bool")))

ext: jl


text:
using IterTools

function generate_cut(layout::RQCLayout, s, cycle_start, cycle_end)
    ncycles = length(gates(layout, 1)) ÷ 2
    vs = sort!(collect(vertices(layout)))
    s1 = filter!(v -> has_vertex(layout, v), s)
    s2 = setdiff(vs, s)
    cycle_end > ncycles && (cycle_end = ncycles)
    cuts = []
    for c in cycle_start:cycle_end
        for v1 in s1
            v2 = partner(layout, v1, c)
            if v2 in s2
                push!(cuts, (v1, 2*(c-1)+1), (v2, 2*(c-1)+1))
            end
        end
    end
    return cuts
end

function gen_cut_info(layout::RQCLayout, vs_cut)
    vs = collect(vertices(layout))
    cuts_info = Dict{Vector{Int}, Tuple{Int, Int}}()
    ncycle = depth_to_cycle(depth(layout))
    for S in subsets(vs_cut)
        L = setdiff(vs, S)
        layout = google_layout_53(53, ncycle)
        cuts = XEB.generate_cut(layout, L, 1, ncycle)
        XEB.simplify!(layout, cuts);

        n_open = length(filter(v -> gates(layout, v)[end] !== XEB.Noise, L))
        cuts_info[S] = (flux(layout, L), n_open)
    end
    return cuts_info
end


text:
# broadcast
Base.Broadcast.broadcasted(f, t::Tensor, args...) = f(t, args...)
Base.Broadcast.broadcasted(::typeof(Base.:*), a::Tensor, b::Tensor) = ThC.mul(a, b)

# operators
Base.:+(r::TorchNumber, t::Tensor) = ThC.add1(t, r)
Base.:+(t::Tensor, r::TorchNumber) = r + t
Base.:+(a::Tensor{T, N}, b::Tensor{T, N}) where {T, N} = ThC.opt_add(a, b)

Base.:-(r::TorchNumber, t::Tensor) = ThC.ones_like(t) * r - t
Base.:-(t::Tensor, r::TorchNumber) = ThC.sub1(t, r)
Base.:-(a::Tensor{T, N}, b::Tensor{T, N}) where {T, N} = ThC.sub(a, b)
Base.:-(a::Tensor) = 0 - a

Base.:*(r::TorchNumber, t::Tensor) = ThC.mul1(t, r)
Base.:*(t::Tensor, r::TorchNumber) = r * t
Base.:*(a::Tensor{T, N}, b::Tensor{T, N}) where {T, N} = ThC.mm(a, b)

Base.:/(n::TorchNumber, t::Tensor) = ThC.ones_like(t) * n / t
Base.:/(t::Tensor, n::TorchNumber) = ThC.div1(t, n)
Base.:/(a::Tensor{T, N}, b::Tensor{T, N}) where {T, N} = ThC.div(a, b)
Base.div(n::TorchNumber, t::Tensor, r::RoundingMode=RoundToZero) = n / t
Base.div(t::Tensor, n::TorchNumber, r::RoundingMode=RoundToZero) = t / n
Base.div(a::Tensor{T, N}, b::Tensor{T, N}, r::RoundingMode=RoundToZero) where {T, N} = div(a, b)

Base.:^(t::Tensor, r::TorchNumber) = ThC.pow(t, r)

Base.:(==)(t1::Tensor, t2::Tensor) = ThArrays.ThC.all(ThArrays.ThC.eq1(t1, t2))[]

function Base.ones(::Type{Tensor{T}}, I::Vararg{Int}; dev::Device=CPU()) where T
    dims = Int64[I...]
    ThC.ones(dims, eltype_id(T), convert(Int, dev))
end

function Base.zeros(::Type{Tensor{T}}, I::Vararg{Int}; dev::Device=CPU()) where T
    dims = Int64[I...]
    ThC.zeros(dims, eltype_id(T), convert(Int, dev))
end

function Base.rand(::Type{Tensor{T}}, I::Vararg{Int}; dev::Device=CPU()) where T
    dims = Int64[I...]
    ThC.rand(dims, eltype_id(T), convert(Int, dev))
end

ThC.eye(::Type{T}, n::Int64; dev::Device=CPU()) where T =
    ThC.eye(n, eltype_id(T), convert(Int, dev))
ThC.eye(::Type{Tensor{T}}, n::Int64; dev::Device=CPU()) where T =
    ThC.eye(T, n, dev=dev)

Base.sum(t::Tensor{T}) where T = ThC.sum(t, eltype_id(T))
Base.view(t::Tensor{T}, I...) where T = error("Not implement yet.")
Base.transpose(t::Tensor{T, 2}) where T = ThC.t(t)
Base.adjoint(t::Tensor) = error("Not implement yet.")
# LinearAlgebra.det(t::Tensor) = error("Not implement yet.")
# LinearAlgebra.logdet(t::Tensor) = error("Not implement yet.")
# LinearAlgebra.logabsdet(t::Tensor) = error("Not implement yet.")
# Base.repeat(t::Tensor; kw...) = error("Not implement yet.")
# Base.reshape(t::Tensor, dims...) = error("Not implement yet.")
# Base.permutedims(t::Tensor, perm) = error("Not implement yet.")
# Base.PermutedDimsArray(t::Tensor, perm) = error("Not implement yet.")
# Base.reverse(t::Tensor; dims) = error("Not implement yet.")
# Base.reverse(t::Tensor) = error("Not implement yet.")
# Base.reverse(t::Tensor, start, stop) = error("Not implement yet.")
# Base.inv(t::Tensor) = error("Not implement yet.")
# Base.:\(a::Tensor, b::Tensor) = error("Not implement yet.")
# Base.prod(t::Tensor, dim) = error("Not implement yet.")
# Base.prod(t::Tensor) = error("Not implement yet.")
# Base.prod(f::Union{Function, Type}, t::Tensor) = error("Not implement yet.")
# Statistics.mean(t::Tensor; dims = :) = error("Not implement yet.")
# Base.maximum(t::Tensor; dims = :) = error("Not implement yet.")
# Base.minimum(t::Tensor; dims = :) = error("Not implement yet.")
# Base.dot(a::Tensor, b::Tensor) = error("Not implement yet.")
# LinearAlgebra.diagm(...)
# NNlib
# softmax, logsoftmax, depthwiseconv, conv, ∇conv_data, maxpool, meanpool


text:
import Conda
import Pkg

Conda.add_channel("r")
Conda.add("r-base")

if Base.Sys.iswindows()
    separator = ";"
else
    separator = ":"
end
original_path = ENV["PATH"]
new_path = string(
    Conda.BINDIR,
    separator,
    original_path,
    )
ENV["PATH"] = new_path

ENV["R_HOME"] = "*"
Pkg.add("RCall")
Pkg.build("RCall")

import RCall


text:
using PyPlot
using JLD
#using Optim

# type Measurements{T}
#     i::Int64                # measurement counter
#     t::Array{Float64}       # time data (when it was received by this recorder)
#     t_msg::Array{Float64}   # time that the message was sent
#     z::Array{T}             # measurement values
# end

function data_analysis(code::AbstractString)

	log_path_record = "$(homedir())/open_loop/output-record-$code.jld"
    d_rec = load(log_path_record)

    cmd_pwm_log = d_rec["cmd_pwm_log"]
    vel_est     = d_rec["vel_est"]
    # t_cmd       = d_rec["t_cmd"]
    # t_no_cmd    = d_rec["t_no_cmd"]

    t0 = max(cmd_pwm_log.t[1],vel_est.t[1])
    t_end = min(cmd_pwm_log.t[end],vel_est.t[end])

    t = t0+0.1:.02:t_end-0.1
    v = zeros(length(t))
    cmd = zeros(length(t))

    for i=1:length(t)
        v[i] = vel_est.z[t[i].>vel_est.t,1][end]
        #v[i] = pos_info.z[t[i].>pos_info.t,15][end]
        cmd[i] = cmd_pwm_log.z[t[i].>cmd_pwm_log.t,1][end]
    end

    v_opt   = v[1:end]
    cmd_opt = cmd[1:end]
    v_opt   = v[cmd.>90]
    cmd_opt = cmd[cmd.>90]

    figure()
    plot(v_opt)

    # half = round(length(v_opt)/2)
    # println(half)

    #v_mean = mean(v_opt[half:end])

    #v_mean = mean(v_opt[length(v_opt)/2:end])

    #a = v_mean/4

    a = v_opt[end]/4

    println("a= ",a)

    figure(2)

    plot(cmd_opt[1],a, "*")

    # println(t_cmd)
    # println(t_no_cmd)
    # println(t0)

    # figure()
    # plot(cmd,v)

    # figure()
    # plot(t,cmd)

    # figure()
    # plot(t,v)

    # figure()
    # plot(v)

    # figure()
    # plot(t)
end

max_stars_repo_name: RoyaFiroozi/barc

max_issues_repo_name: RoyaFiroozi/barc

max_forks_repo_name: RoyaFiroozi/barc


text:
using Test
using ChainedFixes
using ChainedFixes.ChainedCore
using Documenter
using Base: Fix1

empty_named_tuple = NamedTuple{(),Tuple{}}(())

@test @inferred(getfxn(1)) == identity
@test @inferred(getfxn(+)) == +

@test !is_fixed_function(1)

@testset "Fix1" begin
    fix1_fxn = Fix1(<, 1)
    @test is_fixed_function(typeof(fix1_fxn))
    @test @inferred(getargs(fix1_fxn)) == (1,)
    @test @inferred(getfxn(fix1_fxn)) == <
    @test @inferred(getkwargs(fix1_fxn)) == empty_named_tuple
end

@testset "Not" begin
    notfxn = !(+)
    @test isa(notfxn, Not)
    @test !isa(+, Not)
    @test @inferred(is_fixed_function(typeof(notfxn)))
    @test @inferred(getfxn(notfxn)) == !
    @test @inferred(getargs(notfxn)) == (+,)
    @test @inferred(getkwargs(notfxn)) == empty_named_tuple
end

@testset "In" begin
    infxn = in(1)
    @test isa(infxn, In)
    @test !isa(+, In)

    @test @inferred(is_fixed_function(typeof(infxn)))
    @test @inferred(getfxn(infxn)) == in
    @test @inferred(getargs(infxn)) == (1,)
    @test @inferred(getkwargs(infxn)) == empty_named_tuple
end

@testset "NotIn" begin
    notinfxn = !in(1)
    @test isa(notinfxn, NotIn)
    @test !isa(!(==(1)), In)

    @test @inferred(is_fixed_function(typeof(notinfxn)))
    @test @inferred(getfxn(notinfxn)) == !
    @test @inferred(getargs(notinfxn)) == (in(1),)
    @test @inferred(getkwargs(notinfxn)) == empty_named_tuple
end

@testset "Approx" begin
    isapprox_fxn = isapprox(1; atol=2)
    @test isa(isapprox_fxn, Approx)
    @test !isa(!(==(1)), Approx)

    @test @inferred(is_fixed_function(typeof(isapprox_fxn)))
    @test @inferred(getfxn(isapprox_fxn)) == isapprox
    @test @inferred(getargs(isapprox_fxn)) == (1,)
    @test @inferred(getkwargs(isapprox_fxn)) == (atol=2,)
end

@testset "NotApprox" begin
    notisapprox_fxn = !isapprox(1; atol=2)
    @test isa(notisapprox_fxn, NotApprox)
    @test !isa(!(==(1)), NotApprox)

    @test @inferred(is_fixed_function(typeof(notisapprox_fxn)))
    @test @inferred(getfxn(notisapprox_fxn)) == !
    @test @inferred(getargs(notisapprox_fxn)) == (isapprox(1; atol=2),)
end

@testset "Less" begin
    ltfxn = <(1)
    @test isa(ltfxn, Less)
    @test !isa(isapprox(1), Less)

    @test @inferred(is_fixed_function(typeof(ltfxn)))
    @test @inferred(getfxn(ltfxn)) == <
    @test @inferred(getargs(ltfxn)) == (1,)
    @test @inferred(getkwargs(ltfxn)) == empty_named_tuple
end

@testset "Equal" begin
    eqfxn = ==(1)
    @test isa(eqfxn, Equal)
    @test !isa(isapprox(1), Equal)

    @test @inferred(is_fixed_function(typeof(eqfxn)))
    @test @inferred(getfxn(eqfxn)) == ==
    @test @inferred(getargs(eqfxn)) == (1,)
    @test @inferred(getkwargs(eqfxn)) == empty_named_tuple
end

@testset "EndsWith" begin
    endswith_fxn = endswith("i")
    @test endswith("i")("hi")
    @test endswith("i") isa EndsWith

    @test @inferred(is_fixed_function(typeof(endswith_fxn)))
    @test @inferred(getfxn(endswith_fxn)) == endswith
    @test @inferred(getargs(endswith_fxn)) == ("i",)
    @test @inferred(getkwargs(endswith_fxn)) == empty_named_tuple
end

@testset "StartsWith" begin
    startswith_fxn = startswith("h")
    @test startswith("h")("hi")
    @test startswith("h") isa StartsWith

    @test @inferred(is_fixed_function(typeof(startswith_fxn)))
    @test @inferred(getfxn(startswith_fxn)) == startswith
    @test @inferred(getargs(startswith_fxn)) == ("h",)
    @test @inferred(getkwargs(startswith_fxn)) == empty_named_tuple
end

@testset "and" begin
    @test and(true, <(5))(1)
    @test !and(<(5), false)(1)
    @test and(and(<(5), >(1)), >(2))(3)
    @test and(<(5) ⩓ >(1), >(2))(3)  # ⩓ == \\And

    @test @inferred(and(true, true))
    @test !@inferred(and(false, false, ==(true)))
    @test !@inferred(and(false, ==(true), false))
    @test @inferred(and(<=(1), <=(10))) == <=(1)
    @test @inferred(and(<=(10), <=(1))) == <=(1)
    @test @inferred(and(<(1), <(10))) == <(1)
    @test @inferred(and(<(10), <(1))) == <(1)
    @test @inferred(and(>(1), >(10))) == >(10)
    @test @inferred(and(>(10), >(1))) == >(10)
    @test @inferred(and(>=(1), >=(10))) == >=(10)
    @test @inferred(and(>=(10), >=(1))) == >=(10)

    and_fxn = and(true, <(5))
    @test @inferred(is_fixed_function(typeof(and_fxn)))
    @test @inferred(getfxn(and_fxn)) == and
    @test @inferred(getargs(and_fxn)) == (true, <(5))

    @test_throws MethodError and(true, true, true)
end

@testset "or" begin
    @test or(true, <(5))(1)
    @test or(<(5), false)(1)
    @test or(<(5) ⩔ >(1), >(2))(3)  # ⩔ == \\Or
    @test @inferred(or(true, false, ==(true)))
    @test @inferred(or(true, false))
    @test @inferred(or(false, ==(true), true))
    @test @inferred(or(<(1), <(10))) == <(10)
    @test @inferred(or(<(10), <(1))) == <(10)
    @test @inferred(or(<=(1), <=(10))) == <=(10)
    @test @inferred(or(<=(10), <=(1))) == <=(10)
    @test @inferred(or(>(1), >(10))) == >(1)
    @test @inferred(or(>(10), >(1))) == >(1)
    @test @inferred(or(>=(1), >=(10))) == >=(1)
    @test @inferred(or(>=(10), >=(1))) == >=(1)

    or_fxn = or(true, <(5))
    @test getfxn(or_fxn) == or
    @test getargs(or_fxn) == (true, <(5))
    @test_throws MethodError or(true, true, true)
end

@test ChainedFixes.ArgPosition(1) === ChainedFixes.ArgPosition{1}()

fxn1(x::Integer, y::AbstractFloat, z::AbstractString) = Val(1)
fxn1(x::Integer, y::AbstractString, z::AbstractFloat) = Val(2)
fxn1(x::AbstractFloat, y::Integer, z::AbstractString) = Val(3)
fxn1(x::AbstractFloat, y::AbstractString, z::Integer) = Val(4)
fxn1(x::AbstractString, y::Integer, z::AbstractFloat) = Val(5)
fxn1(x::AbstractString, y::AbstractFloat, z::Integer) = Val(6)
fxn2(; x, y, z) = fxn1(x, y, z)
fxn3(args...; kwargs...) = (fxn1(args...), fxn2(; kwargs...))

f = @nfix fxn1(1, 2.0, _)
@test @inferred(f("a")) == Val{1}()
@test @inferred(is_fixed_function(typeof(f)))

f = @nfix fxn1(1, _, 2.0)
@test @inferred(f("a")) == Val{2}()

f = @nfix fxn1(1.0, _, "")
@test @inferred(f(2)) == Val{3}()

f = @nfix fxn2(x=1, y=2.0)
@test @inferred(f(z = "a")) == Val{1}()

f = @nfix fxn2(x=1, z=2.0)
@test @inferred(f(y = "a")) == Val{2}()

f = @nfix fxn3(1, 2.0, _; x = 1.0, z= "")
@test @inferred(f(""; y = 1)) == (Val{1}(), Val{3}())

f = pipe_chain(@nfix(_ * "is "), @nfix(_ * "a "), @nfix(_ * "sentence."))
@test f("This ") == "This is a sentence."

f2 = pipe_chain(f, endswith("sentence."))
@test f2("This ")

f2 = pipe_chain(f, startswith("This"))
@test f2("This ")

f = pipe_chain(and(<=(3), !=(2)), ==(true), in(trues(2)), !in(falses(2)))
@test f(1)

f = pipe_chain(isapprox(0.1), !isapprox(0.2))
@test f(0.1 - 1e-10)

splat_pipe(op, args::Tuple) = op(args...)
splat_pipe(op) = @nfix splat_pipe(op, _...)
f = pipe_chain(extrema, splat_pipe(+))
@test @inferred(f([1 2; 3 4])) == 5

@testset "docs" begin
    doctest(ChainedFixes)
end


text:
struct VectorStrategy <: SearchStrategy
    elements::SearchStrategy
    minsize::UInt
    maxsize::UInt

    VectorStrategy(elements::SearchStrategy, minsize::Nothing, maxsize::Nothing) = new(elements, 0, 9999)
    VectorStrategy(elements::SearchStrategy, minsize::Nothing, maxsize) = new(elements, 0, maxsize)
    VectorStrategy(elements::SearchStrategy, minsize, maxsize::Nothing) = new(elements, minsize, 9999)
    function VectorStrategy(elements::SearchStrategy, minsize, maxsize)
        minsize <= maxsize || error("VectorStrategy has unsatisfiable bounds")
        new(elements, minsize, maxsize)
    end
end

function draw(strategy::VectorStrategy)
    len = rand(strategy.minsize:strategy.maxsize)
    result = Base.eltype(strategy)()
    for _ in 1:len
        push!(result, draw(strategy.elements))
    end
    return result
end

vectors(elements::SearchStrategy; minsize=nothing, maxsize=nothing) = VectorStrategy(elements, minsize, maxsize)
Base.eltype(s::VectorStrategy) = Vector{Base.eltype(s.elements)}


text:
using SLFN
using Test
using Random

const TS = (AlgebraicNetwork, ELM, OSELM, ROSELM, TSELM)
const verbose = true

Random.seed!(42)
x = range(0, stop=1, length=50) + 0.1*randn(50); y = sin.(6x)
for T in TS
    n = T(x, y; s=30)
    verbose && println(typeof(n), " 1")
    @test maximum(abs.(n(x) - y)) < 1e-2
end

# multi-dim version
nx, ny, nz = 10, 8, 9
x = range(0, stop=1, length=nx)
y = range(0, stop=1, length=ny)
z = range(0.1, stop=1, length=nz)

f1(x) = @. sin(2pi/nx * (x-1))
f2(x, y) = @. f1(x) .* cos(4pi/ny * (y-1))
f3(x, y, z) = @. f2(x, y) .* log(z)

# 2d case
X2 =[kron(ones(ny), x) kron(y, ones(nx))]
Y2 = f2(X2[: ,1], X2[:, 2])

for T in TS
    n = T(X2, Y2, s=40)
    verbose && println(typeof(n), " 2")
    @test maximum(abs.(n(X2) - Y2)) < 1e-1
end

X3 =[kron(ones(ny*nz), x) kron(ones(nz), y, ones(nx)) kron(z, ones(nx*ny))]
Y3 = f3(X3[:, 1], X3[:, 2], X3[:, 3])

for T in TS
    n = T(X3, Y3, s=200)
    verbose && println(typeof(n), " 3")
    @test maximum(abs.(n(X3) - Y3)) < 1e-1
end

# test non [0, 1] domain works (i.e. standardizing is working)
x = 150 * rand(5, 5)
y = sin.(range(0, stop=6, length=size(x, 1)))

for T in TS
    n = T(x, y, s=size(x, 1))
    verbose && println(typeof(n), " standardizing")
    @test maximum(abs.(n(x) - y)) < 1e-10
end


text:
module OSMTests

using OpenStreetMap
using Base.Test

tests = [
    "read_data",
    "simcity",
    "crop_map",
    "classes",
    "intersections",
    "routes",
    "plots",
    "examples" ]

for t in tests
    println("testing $t ...")
    @time include("$t.jl")
end

end # Module OSMTests


text:
@testset "State Bound Constraint" begin

    # Test ControlBoundConstraint
    T = Float64
    dt = 0.1
    n = 6
    m = 6
    X = SVector{n,T}([13.0, 1.0, -12.0, 1.0, 2.0, 30.0])
    U = SVector{m,T}([13.0, 1.0, -12.0, 1.0, 2.0, 30.0])
    Z = KnotPoint(X,U,dt)
    x_max = [Inf,   Inf,  -11.0, 15.0, 2.0,  30.0]
    x_min = [-Inf, -10.0, -Inf,  1.0, -2.0, -30.0]

    con = StateBoundConstraint(n,x_max=x_max, x_min=x_min)
    @test TrajectoryOptimization.evaluate(con,Z) == [-1.0, -14.0, 0.0, 0.0, -11.0, 0.0, -4.0, -60.0]
    con.inds == [3,4,5,6,8,10,11,12]
    @test (@ballocated TrajectoryOptimization.evaluate($con,$Z)) == 0

    function easy_jacobian(con,Z)
        X = state(Z)
        function local_evaluate(X)
            Z_ = KnotPoint(X, control(Z), Z.dt)
            return TrajectoryOptimization.evaluate(con,Z_)
        end
        return ForwardDiff.jacobian(local_evaluate, state(Z))
    end
    ∇c_easy = easy_jacobian(con, Z)
    ∇c = zeros(MMatrix{8,n,T,8*n})
    jacobian!(∇c,con,Z)
    @test norm(∇c - ∇c_easy, 1) < 1e-10
    @test (@ballocated jacobian!($∇c,$con,$Z)) == 0

end


text:
using Receive: Break

struct SelectInterrupt <: Exception end

function waiter(ch::Channel, queue::Channel, ping::Condition)
    try
        while true
            wait(ch)
            put!(queue, ch)
            wait(ping)
            yield()
        end
    catch err
        if err isa Union{InvalidStateException, SelectInterrupt}
            return nothing
        end
    end
end

function receive_csp(f::Function, ch1::Channel, ch2::Channel)

    queue = Channel{Channel}()

    ping = Condition()

    t1 = @async waiter(ch1, queue, ping)
    t2 = @async waiter(ch2, queue, ping)

    try

        while true

            ch = take!(queue)
            value = take!(ch)

            res = f(ch, value)

            if res isa Break
                break
            end

            notify(ping)
        end

    finally
        close(queue)

        istaskdone(t1) || Base.schedule(t1, SelectInterrupt(), error=true)
        istaskdone(t2) || Base.schedule(t2, SelectInterrupt(), error=true)
    end
end


text:
# Superselection sectors (quantum numbers):
# for defining graded vector spaces and invariant subspaces of tensor products
#==============================================================================#
"""
    abstract type Sector end

Abstract type for representing the (isomorphism classes of) simple objects in (unitary
and pivotal) (pre-)fusion categories, e.g. the irreducible representations of a finite or
compact group. Subtypes `I<:Sector` as the set of labels of a `GradedSpace`.

Every new `I<:Sector` should implement the following methods:
*   `one(::Type{I})`: unit element of `I`
*   `conj(a::I)`: ``a̅``, conjugate or dual label of ``a``
*   `⊗(a::I, b::I)`: iterable with unique fusion outputs of ``a ⊗ b``
    (i.e. don't repeat in case of multiplicities)
*   `Nsymbol(a::I, b::I, c::I)`: number of times `c` appears in `a ⊗ b`, i.e. the
    multiplicity
*   `FusionStyle(::Type{I})`: `Abelian()`, `SimpleNonAbelian()` or
    `DegenerateNonAbelian()`
*   `BraidingStyle(::Type{I})`: `Bosonic()`, `Fermionic()`, `Anyonic()`, ...
*   `Fsymbol(a::I, b::I, c::I, d::I, e::I, f::I)`: F-symbol: scalar (in case of
    `Abelian`/`SimpleNonAbelian`) or matrix (in case of `DegenerateNonAbelian`)
*   `Rsymbol(a::I, b::I, c::I)`: R-symbol: scalar (in case of
    `Abelian`/`SimpleNonAbelian`) or matrix (in case of `DegenerateNonAbelian`)
and optionally
*   `dim(a::I)`: quantum dimension of sector `a`
*   `frobeniusschur(a::I)`: Frobenius-Schur indicator of `a`
*   `Bsymbol(a::I, b::I, c::I)`: B-symbol: scalar (in case of
    `Abelian`/`SimpleNonAbelian`) or matrix (in case of `DegenerateNonAbelian`)
*   `twist(a::I)` -> twist of sector `a`
and optionally, if `FusionStyle(I) isa DegenerateNonAbelian`
*   `vertex_ind2label(i::Int, a::I, b::I, c::I)` -> a custom label for the `i`th copy of
    `c` appearing in `a ⊗ b`

Furthermore, `iterate` and `Base.IteratorSize` should be made to work for the singleton type
[`SectorValues{I}`](@ref).
"""
abstract type Sector end

# iterator over the values (i.e., elements of representative set of simple objects)
# in the sector
"""
    struct SectorValues{I<:Sector}

Singleton type to represent an iterator over the possible values of type `I`, whose
instance is obtained as `values(I)`. For a new `I::Sector`, the following should be defined
*   `Base.iterate(::SectorValues{I}[, state])`: iterate over the values
*   `Base.IteratorSize(::Type{SectorValues{I}})`: `HasLenght()`, `SizeUnkown()`
    or `IsInfinite()` depending on whether the number of values of type `I` is finite
    (and sufficiently small) or infinite; for a large number of values, `SizeUnknown()` is
    recommend because this will trigger the use of `GenericGradedSpace`.
If `IteratorSize(I) == HasLength()`, also the following must be implemented:
*   `Base.length(::SectorValues{I})`: the number of different values
*   `Base.getindex(::SectorValues{I}, i::Int)`: a mapping between an index `i` and an
    instance of `I`
*   `findindex(::SectorValues{I}, c::I)`: reverse mapping between a value `c::I` and an
    index `i::Integer ∈ 1:length(values(I))`
"""
struct SectorValues{I<:Sector} end
Base.IteratorEltype(::Type{<:SectorValues}) = HasEltype()
Base.eltype(::Type{SectorValues{I}}) where {I<:Sector} = I
Base.values(::Type{I}) where {I<:Sector} = SectorValues{I}()

# Define a sector for ungraded vector spaces
struct Trivial <: Sector
end
Base.show(io::IO, ::Trivial) = print(io, "Trivial()")

Base.IteratorSize(::Type{SectorValues{Trivial}}) = HasLength()
Base.length(::SectorValues{Trivial}) = 1
Base.iterate(::SectorValues{Trivial}, i = false) = return i ? nothing : (Trivial(), true)
Base.getindex(::SectorValues{Trivial}, i::Int) = i == 1 ? Trivial() :
    throw(BoundsError(values(Trivial), i))
findindex(::SectorValues{Trivial}, c::Trivial) = 1

"""
    one(::Sector) -> Sector
    one(::Type{<:Sector}) -> Sector

Return the unit element within this type of sector.
"""
Base.one(a::Sector) = one(typeof(a))
Base.one(::Type{Trivial}) = Trivial()

"""
    dual(a::Sector) -> Sector

Return the conjugate label `conj(a)`.
"""
dual(a::Sector) = conj(a)
Base.conj(::Trivial) = Trivial()

"""
    isreal(::Type{<:Sector}) -> Bool

Return whether the topological data (Fsymbol, Rsymbol) of the sector is real or not (in
which case it is complex).
"""
function Base.isreal(I::Type{<:Sector})
    u = one(I)
    return (eltype(Fsymbol(u, u, u, u, u, u))<:Real) && (eltype(Rsymbol(u, u, u))<:Real)
end
Base.isreal(::Type{Trivial}) = true

Base.isless(::Trivial, ::Trivial) = false

# FusionStyle: the most important aspect of Sector
#---------------------------------------------
"""
    ⊗(a::I, b::I) where {I<:Sector}

Return an iterable of elements of `c::I` that appear in the fusion product `a ⊗ b`.

Note that every element `c` should appear at most once, fusion degeneracies (if
`FusionStyle(I) == DegenerateNonAbelian()`) should be accessed via `Nsymbol(a, b, c)`.
"""
⊗(::Trivial, ::Trivial) = (Trivial(),)

"""
    Nsymbol(a::I, b::I, c::I) where {I<:Sector} -> Integer

Return an `Integer` representing the number of times `c` appears in the fusion product
`a ⊗ b`. Could be a `Bool` if `FusionStyle(I) == Abelian()` or `SimpleNonAbelian()`.
"""
function Nsymbol end
Nsymbol(::Trivial, ::Trivial, ::Trivial) = true

# trait to describe the fusion of superselection sectors
abstract type FusionStyle end
struct Abelian <: FusionStyle
end
abstract type NonAbelian <: FusionStyle end
struct SimpleNonAbelian <: NonAbelian # non-abelian fusion but multiplicity free
end
struct DegenerateNonAbelian <: NonAbelian # non-abelian fusion with multiplicities
end

"""
    FusionStyle(a::Sector) -> ::FusionStyle
    FusionStyle(I::Type{<:Sector}) -> ::FusionStyle

Return the type of fusion behavior of sectors of type I, which can be either
*   `Abelian()`: single fusion output when fusing two sectors;
*   `SimpleNonAbelian()`: multiple outputs, but every output occurs at most one,
    also known as multiplicity free (e.g. irreps of ``SU(2)``);
*   `DegenerateNonAbelian()`: multiple outputs that can occur more than once (e.g. irreps
    of ``SU(3)``).
There is an abstract supertype `NonAbelian` of which both `SimpleNonAbelian` and
`DegenerateNonAbelian` are subtypes.
"""
FusionStyle(a::Sector) = FusionStyle(typeof(a))
FusionStyle(::Type{Trivial}) = Abelian()

# NOTE: the following inline is extremely important for performance, especially
# in the case of Abelian, because ⊗(...) is computed very often
@inline function ⊗(a::I, b::I, c::I, rest::Vararg{I}) where {I<:Sector}
    if FusionStyle(I) isa Abelian
        return a ⊗ first(⊗(b, c, rest...))
    else
        s = Set{I}()
        for d in ⊗(b, c, rest...)
            for e in a ⊗ d
                push!(s, e)
            end
        end
        return s
    end
end

"""
    Fsymbol(a::I, b::I, c::I, d::I, e::I, f::I) where {I<:Sector}

Return the F-symbol ``F^{abc}_d`` that associates the two different fusion orders of sectors
`a`, `b` and `c` into an ouput sector `d`, using either an intermediate sector ``a ⊗ b → e``
or ``b ⊗ c → f``:
```
a-<-μ-<-e-<-ν-<-d                                     a-<-λ-<-d
    ∨       ∨       -> Fsymbol(a,b,c,d,e,f)[μ,ν,κ,λ]      ∨
    b       c                                             f
                                                          v
                                                      b-<-κ
                                                          ∨
                                                          c
```
If `FusionStyle(I)` is `Abelian` or `SimpleNonAbelian`, the F-symbol is a number. Otherwise
it is a rank 4 array of size
`(Nsymbol(a, b, e), Nsymbol(e, c, d), Nsymbol(b, c, f), Nsymbol(a, f, d))`.
"""
function Fsymbol end
Fsymbol(::Trivial, ::Trivial, ::Trivial, ::Trivial, ::Trivial, ::Trivial) = 1

"""
    Rsymbol(a::I, b::I, c::I) where {I<:Sector}

Returns the R-symbol ``R^{ab}_c`` that maps between ``c → a ⊗ b`` and ``c → b ⊗ a`` as in
```
a -<-μ-<- c                                 b -<-ν-<- c
     ∨          -> Rsymbol(a,b,c)[μ,ν]           v
     b                                           a
```
If `FusionStyle(I)` is `Abelian()` or `SimpleNonAbelian()`, the R-symbol is a number.
Otherwise it is a square matrix with row and column size `Nsymbol(a,b,c) == Nsymbol(b,a,c)`.
"""
function Rsymbol end
Rsymbol(::Trivial, ::Trivial, ::Trivial) = 1

# If a I::Sector with `fusion(I) == DegenerateNonAbelian` fusion wants to have custom vertex
# labels, a specialized method for `vertindex2label` should be added
"""
    vertex_ind2label(k::Int, a::I, b::I, c::I) where {I<:Sector}

Convert the index `k` of the fusion vertex (a,b)->c into a label. For
`FusionStyle(I) == Abelian()` or `FusionStyle(I) == NonAbelian()`, where every fusion
output occurs only once and `k == 1`, the default is to suppress vertex labels by setting
them equal to `nothing`. For `FusionStyle(I) == DegenerateNonAbelian()`, the default is to
just use `k`, unless a specialized method is provided.
"""
vertex_ind2label(k::Int, a::I, b::I, c::I) where {I<:Sector}=
    _ind2label(FusionStyle(I), k::Int, a::I, b::I, c::I)
_ind2label(::Abelian, k, a, b, c) = nothing
_ind2label(::SimpleNonAbelian, k, a, b, c) = nothing
_ind2label(::DegenerateNonAbelian, k, a, b, c) = k

"""
    vertex_labeltype(I::Type{<:Sector}) -> Type

Return the type of labels for the fusion vertices of sectors of type `I`.
"""
vertex_labeltype(I::Type{<:Sector}) = typeof(vertex_ind2label(1, one(I), one(I), one(I)))

# combine fusion properties of tensor products of sectors
Base.:&(f::F, ::F) where {F<:FusionStyle} = f
Base.:&(f1::FusionStyle, f2::FusionStyle) = f2 & f1

Base.:&(::SimpleNonAbelian, ::Abelian) = SimpleNonAbelian()
Base.:&(::DegenerateNonAbelian, ::Abelian) = DegenerateNonAbelian()
Base.:&(::DegenerateNonAbelian, ::SimpleNonAbelian) = DegenerateNonAbelian()

# properties that can be determined in terms of the F symbol
# TODO: find mechanism for returning these numbers with custom type T<:AbstractFloat
"""
    dim(a::Sector)

Return the (quantum) dimension of the sector `a`.
"""
function dim(a::Sector)
    if FusionStyle(a) isa Abelian
        1
    elseif FusionStyle(a) isa SimpleNonAbelian
        abs(1/Fsymbol(a, conj(a), a, a, one(a), one(a)))
    else
        abs(1/Fsymbol(a, conj(a), a, a, one(a), one(a))[1])
    end
end

"""
    frobeniusschur(a::Sector)

Return the Frobenius-Schur indicator of a sector `a`.
"""
function frobeniusschur(a::Sector)
    if FusionStyle(a) isa Abelian || FusionStyle(a) isa SimpleNonAbelian
        sign(Fsymbol(a, conj(a), a, a, one(a), one(a)))
    else
        sign(Fsymbol(a, conj(a), a, a, one(a), one(a))[1])
    end
end

"""
    twist(a::Sector)

Return the twist of a sector `a`
"""
twist(a::Sector) = sum(dim(b)/dim(a)*tr(Rsymbol(a,a,b)) for b in a ⊗ a)

"""
    Bsymbol(a::I, b::I, c::I) where {I<:Sector}

Return the value of ``B^{ab}_c`` which appears in transforming a splitting vertex
into a fusion vertex using the transformation
```
a -<-μ-<- c                                                    a -<-ν-<- c
     ∨          -> √(dim(c)/dim(a)) * Bsymbol(a,b,c)[μ,ν]           ∧
     b                                                            dual(b)
```
If `FusionStyle(I)` is `Abelian()` or `SimpleNonAbelian()`, the B-symbol is a number.
Otherwise it is a square matrix with row and column size
`Nsymbol(a, b, c) == Nsymbol(c, dual(b), a)`.
"""
function Bsymbol(a::I, b::I, c::I) where {I<:Sector}
    if FusionStyle(I) isa Abelian || FusionStyle(I) isa SimpleNonAbelian
        sqrt(dim(a)*dim(b)/dim(c))*Fsymbol(a, b, dual(b), a, c, one(a))
    else
        reshape(sqrt(dim(a)*dim(b)/dim(c))*Fsymbol(a, b, dual(b), a, c, one(a)),
            (Nsymbol(a, b, c), Nsymbol(c, dual(b), a)))
    end
end

# Not necessary
function Asymbol(a::I, b::I, c::I) where {I<:Sector}
    if FusionStyle(I) isa Abelian || FusionStyle(I) isa SimpleNonAbelian
        sqrt(dim(a)*dim(b)/dim(c))*
            conj(frobeniusschur(a)*Fsymbol(dual(a), a, b, b, one(a), c))
    else
        reshape(sqrt(dim(a)*dim(b)/dim(c))*
                    conj(frobeniusschur(a)*Fsymbol(dual(a), a, b, b, one(a), c)),
                (Nsymbol(a, b, c), Nsymbol(dual(a), c, b)))
    end
end

# Braiding:
#-------------------------------------------------
# trait to describe type to denote how the elementary spaces in a tensor product space
# interact under permutations or actions of the braid group
abstract type BraidingStyle end # generic braiding
abstract type SymmetricBraiding <: BraidingStyle end # symmetric braiding => actions of permutation group are well defined
struct Bosonic <: SymmetricBraiding end # trivial under permutations
struct Fermionic <: SymmetricBraiding end
struct Anyonic <: BraidingStyle end

Base.:&(b::B, ::B) where {B<:BraidingStyle} = b
Base.:&(B1::BraidingStyle, B2::BraidingStyle) = B2 & B1

Base.:&(::Bosonic, ::Fermionic) = Fermionic()
Base.:&(::Bosonic, ::Anyonic) = Anyonic()
Base.:&(::Fermionic, ::Anyonic) = Anyonic()

"""
    BraidingStyle(::Sector) -> ::BraidingStyle
    BraidingStyle(I::Type{<:Sector}) -> ::BraidingStyle

Return the type of braiding and twist behavior of sectors of type `I`, which can be either
*   `Bosonic()`: symmetric braiding with trivial twist (i.e. identity)
*   `Fermionic()`: symmetric braiding with non-trivial twist (squares to identity)
*   `Anyonic()`: general ``R_(a,b)^c`` phase or matrix (depending on `SimpleNonAbelian` or
    `DegenerateNonAbelian` fusion) and arbitrary twists

Note that `Bosonic` and `Fermionic` are subtypes of `SymmetricBraiding`, which means that
braids are in fact equivalent to crossings (i.e. braiding twice is an identity:
`isone(Rsymbol(b,a,c)*Rsymbol(a,b,c)) == true`) and permutations are uniquely defined.
"""
BraidingStyle(a::Sector) = BraidingStyle(typeof(a))
BraidingStyle(::Type{Trivial}) = Bosonic()

# SectorSet:
#-------------------------------------------------------------------------------
# Custum generator to represent sets of sectors with type inference
struct SectorSet{I<:Sector, F, S}
    f::F
    set::S
end
SectorSet{I}(::Type{F}, set::S) where {I<:Sector, F, S} = SectorSet{I, Type{F}, S}(F, set)
SectorSet{I}(f::F, set::S) where {I<:Sector, F, S} = SectorSet{I, F, S}(f, set)
SectorSet{I}(set) where {I<:Sector} = SectorSet{I}(identity, set)

Base.IteratorEltype(::Type{<:SectorSet}) = HasEltype()
Base.IteratorSize(::Type{SectorSet{I, F, S}}) where {I<:Sector, F, S} = Base.IteratorSize(S)

Base.eltype(::SectorSet{I}) where {I<:Sector} = I
Base.length(s::SectorSet) = length(s.set)
Base.size(s::SectorSet) = size(s.set)

function Base.iterate(s::SectorSet{I}, args...) where {I<:Sector}
    next = iterate(s.set, args...)
    next === nothing && return nothing
    val, state = next
    return convert(I, s.f(val)), state
end

# possible sectors
include("groups.jl")
include("irreps.jl") # irreps of symmetry groups, with bosonic braiding
# include("fermions.jl") # irreps with defined fermionparity and fermionic braiding
include("anyons.jl") # non-group sectors
include("product.jl") # direct product of different sectors


text:
"""
module Iso3D extends PyPlot to easily create isosurface plots from regularly-
spaced 3d arrays of data.
It requires matplotlib, numpy and scikit-image be importable into the default Python.
"""
module Iso3D

using PyCall, PyPlot

export isosurface

@pyimport mpl_toolkits.mplot3d.art3d as art3d
@pyimport skimage.measure as measure

"""
`isosurface(a::Array{T,3}, v)`

Draw a 3D isosurface at value `v` drawn from the 3D array `a`.
"""
function isosurface{T}(a::Array{T,3}, v)
	minimum(a) < v && maximum(a) > v || error("Isosurface value is not contained in data")
	verts, faces = measure.marching_cubes(a, v)
	nf = size(faces, 1)
	triangles = Array{eltype(verts)}(nf, 3, 3)
	x = verts[:,1]'
	y = verts[:,2]'
	z = verts[:,3]'
	@inbounds for i = 1:nf, j = 1:3, k = 1:3
		triangles[i,j,k] = verts[faces[i,j]+1,k]
	end
	collection = art3d.Poly3DCollection(triangles)
	fig = figure(figsize=(10, 10))
	ax = subplot(111, projection="3d")
	ax[:add_collection](collection)
	xlim(minimum(x), maximum(x))
	ylim(minimum(y), maximum(y))
	zlim(minimum(z), maximum(z))
	xlabel("x")
	ylabel("y")
	zlabel("z")
end

end

max_stars_repo_name: anowacki/assorted-julia-modules

max_issues_repo_name: anowacki/assorted-julia-modules

max_forks_repo_name: anowacki/assorted-julia-modules


text:
#TODO: we could first add the VTKPlotter here
#TODO: then matplotlib plotting would be nicer (maybe we can do it via Plots.jl?)


text:
function simpson(x0,xN,r,dx,f,Param)
  n=(xN-x0)/dx+1;
  h=(xN-x0)/n;
  res=0.5*(f(x0,r,Param)+f(xN,r,Param));
  xi=x0;
  for i=1:(n-1)
    xi=xi+h;
    res=res+f(xi,r,Param);
  end
  xi=x0-0.5*h;
  for i=1:n
    xi=xi+h;
    res=res+2.0*f(xi,r,Param);
  end
  res=res*h/3.0;
  return res
end


text:
function handle_drop(files::Vector{String})
    for f in files
        try
            glplot(load(f))
        catch e
            warn(e)
        end
    end
end

toolbar_area(pa, toolbar_width) = IRect(0, 0, toolbar_width, pa.h)

function viewing_area(area_l, area_r)
    IRect(area_l.x+area_l.w, 0, area_r.x-area_l.w, area_r.h)
end

function edit_rectangle(visible, area, tarea)
    w = visible ? 70mm : 1.5mm
    x = area.w - w
    IRect(x, 0, w, area.h)
end

function edit_item_area(la, item_height, left_gap)
    y = la.y - item_height-2
    return IRect(left_gap, y, la.w, item_height)
end

layout_pos_ho(i) = map(icon_size()) do ip
    SimpleRectangle{Float32}(0, i*ip + i*2, ip, ip)
end
layout_pos_ver(i, border) = map(icon_size()) do ip
    SimpleRectangle{Float32}(i*ip + i*border, 0, ip, ip)
end

function save_record(frames)
    path = joinpath(homedir(), "Desktop")
    GLVisualize.create_video(frames, "test.webm", path, 0)
end

const _compute_callbacks = []
register_compute(f) = push!(_compute_callbacks, f)
poll_reactive() = (Base.n_avail(Reactive._messages) > 1) && Reactive.run_till_now()

const plotting_screens = Screen[]

viewing_screen() = (init(); plotting_screens[1])
edit_screen() = (init(); plotting_screens[2])
tool_screen() = (init(); plotting_screens[3])

function glplot_renderloop(window, compute_s, record_s)
    frames = []
    i = 0
    Reactive.stop()
    while isopen(window)
        tic()
        GLWindow.poll_glfw() # GLFW poll
        if Base.n_avail(Reactive._messages) > 0
            GLWindow.poll_reactive() # reactive poll
            record = !value(record_s)
            GLWindow.render_frame(window)
            GLWindow.swapbuffers(window)
        end
        yield()
        diff = (1 / 60) - toq()
        while diff >= 0.001
            tic()
            sleep(0.001) # sleep for the minimal amount of time
            diff -= toq()
        end
        yield()
    end
    GLWindow.destroy!(window)
end

const _icon_size = Signal(10mm)
icon_size() = _icon_size

global set_screenshotpath!, get_screenshotpath

let

    const screenshotpath = Ref(joinpath(homedir(), "Desktop", "glplot.png"))
    """
    Sets the current screen shot/recording path. Default is on users desktop
    """
    function set_screenshotpath!(path)
        screenshotpath[] = path
    end

    """
    Gets the current screen shot/recording path. Default is on users desktop
    """
    function get_screenshotpath()
        screenshotpath[]
    end

end

function init()
    if !isempty(plotting_screens) && isopen(first(plotting_screens))
        return # already initialized
    end
    empty!(plotting_screens)
    w = glscreen("GLPlot")

    preserve(map(handle_drop, w.inputs[:dropped_files]))

    w.inputs[:key_pressed] = const_lift(GLAbstraction.singlepressed,
        w.inputs[:mouse_buttons_pressed],
        GLFW.MOUSE_BUTTON_LEFT
    )
    button_pos = map(w.area) do a
        Point2f0[(0, a.h/2)]
    end
    button_width = 1.5mm
    edit_screen_show_button = visualize(
        (SimpleRectangle{Float32}(0, 0, button_width, button_width*2), button_pos),
        offset = Vec2f0(0, -button_width),
        color = RGBA{Float32}(0.6,0.6,0.6,1)
    )
    tarea = map(toolbar_area, w.area, icon_size())

    show_edit_screen = toggle(edit_screen_show_button, w, false)
    edit_screen_area = map(edit_rectangle,
        show_edit_screen, w.area, tarea
    )

    viewing_screen = Screen(w,
        name = Symbol("Viewing Screen"),
        area = map(viewing_area, tarea, edit_screen_area),
        color = RGBA{Float32}(1,1,1)
    )
    toolbar_screen = Screen(w, area = tarea)
    edit_screen = Screen(
        w, area = edit_screen_area,
        color = RGBA{Float32}(0.9, 0.9, 0.9)
    )
    push!(plotting_screens, viewing_screen)
    push!(plotting_screens, edit_screen)
    push!(plotting_screens, toolbar_screen)

    _view(edit_screen_show_button, edit_screen, camera = :fixed_pixel)
    play_record, record_sig = toggle_button(imload("record.png"), imload("break.png"), w)
    compute_record, compute_sig = toggle_button(imload("play.png"), imload("break.png"), w)
    persp_ortho, persp_ortho_toggle_sig = toggle_button(imload("perspective.png"), imload("ortho.png"), w)

    persp_ortho_sig = map(persp_ortho_toggle_sig) do isp
        isp && return GLAbstraction.PERSPECTIVE
        GLAbstraction.ORTHOGRAPHIC
    end
    cube = cubecamera(viewing_screen, persp_ortho_sig)

    image_names = ["center", "screenshot"]
    tools = Matrix{BGRA{N0f8}}[imload("$name.png") for name in image_names]
    center_b, center_s = button(tools[1], w)
    screenshot_b, screenshot_s = button(tools[2], w)

    tools = [center_b, screenshot_b, play_record, persp_ortho, compute_record]
    tools_robjs = Any[]

    i = 0
    for tool in tools
        robj = layout!(layout_pos_ho(i), visualize(tool))
        _view(robj, toolbar_screen, camera = :fixed_pixel)
        push!(tools_robjs, robj.children[])
        i += 1
    end
    preserve(map(center_s) do pressed
        if pressed
            center!(viewing_screen)
        end
        nothing
    end)
    preserve(map(screenshot_s) do pressed
        if pressed
            screenshot(viewing_screen, path = get_screenshotpath())
        end
        nothing
    end)
    rot = cube.children[][:model]

    cube.children[][:model] = map(rot, icon_size()) do r, ip
        half = ip/2
        translationmatrix(Vec3f0(half,i*ip + i*2 + half,0))*r*scalematrix(Vec3f0(half))
    end
    _view(cube, toolbar_screen, camera=:fixed_pixel)

    GLVisualize.add_screen(viewing_screen)
    @materialize scroll, mouseposition = edit_screen.inputs
    should_scroll = map(mouseposition) do mb
        isinside(value(w.area), mb...)
    end
    scroll = filterwhen(should_scroll, value(scroll), scroll)
    edit_screen.inputs[:menu_scroll] = foldp(0, scroll) do v0, s
        v0 + (ceil(Int, s[2]) * 15)
    end
    global _renderloop_task = @async glplot_renderloop(w, compute_sig, record_sig)
    viewing_screen
end

"""
Blocks the julia process as long as a window is open.
Usefull when running a GLPlot script not within some REPL.
"""
function block()
    global _renderloop_task
    wait(_renderloop_task)
end


text:
import Base: inv

"""
Family of nonlinear transformation parametrized by λ that defaults to Box Cox.
N.B. We assume in derivatives.jl that the nonlinear transformation f is monotonic increasing.
"""

@doc raw"""
    BoxCox(λ::Real)

The Box-Cox power transformation function with parameter λ.

```math
g_\lambda(x) = \begin{cases}
    \frac{x^\lambda - 1}{\lambda} & \lambda \neq 0\\
    \ln(x) & \lambda = 0
\end{cases}
```

External links
* [Power Transform on Wikipedia](https://en.wikipedia.org/wiki/Power_transform)

...
# Arguments
* `x::float64`: evaluation point
* `lambda::float=1.0`: lambda hyperparameter

"""
abstract type AbstractTransform end

abstract type NonlinearTransform <: AbstractTransform end

struct BoxCox <:NonlinearTransform end

struct anotherone<:NonlinearTransform end

#jitter = 1e-6
(::BoxCox)(x, λ::Union{Array{T, 1}, T} where T<:Real) = (@assert (typeof(λ) <:Array{T} where T<:Real && size(λ, 2)==1) || (typeof(λ)<:Real);
                                                 λ[1] == 0 ? log.(x ) : expm1.(log.(x) .* λ[1]) ./ λ[1])
partialx(::BoxCox, x, λ) =  λ[1]==0 ? float(x).^(-1) : float(x).^(λ[1] .-1)
partialxx(::BoxCox, x, λ) = λ[1]==0 ? - float(x).^(- 2) : (λ[1] .- 1) * float(x) .^ (λ[1] .- 2)
partialλ(::BoxCox, x, λ) = λ[1]==0 ? 0 : (λ[1] .* float(x).^λ[1] .* log.(x) .- float(x).^λ[1] .+ 1)/λ^2
partialλλ(::BoxCox, x, λ) = (λ[1]==0 ? 0 : float(x).^λ[1] .* (Base.log.(x)^2)/λ[1] -
                            (2 * x.^λ[1] .* Base.log.(x))/λ[1]^2 + 2*(float(x).^λ[1] - 1)/λ[1]^3)
partialxλ(::BoxCox, x, λ) = λ[1]==0 ? 0 : float(x).^(λ[1] -1) .* Base.log.(x)
inverse(::BoxCox, x, λ) = λ[1]==0 ? Base.exp.(x) : Base.exp.(Base.log.(λ[1].*x.+1)./λ[1])

struct ShiftedBoxCox <:NonlinearTransform end
(::ShiftedBoxCox)(x, λ::Union{Array{T, 1}, T} where T<:Real) = (@assert (typeof(λ) <:Array{T} where T<:Real && size(λ, 2)==1) || (typeof(λ)<:Real);
                                                 λ[1] == 0 ? log.(x .+ 1) : expm1.(log.(x .+ 1) .* λ[1]) ./ λ[1])
partialx(::ShiftedBoxCox, x, λ) =  λ[1]==0 ? float(x .+ 1).^(-1) : float(x .+ 1).^(λ[1] .-1)
partialxx(::ShiftedBoxCox, x, λ) = λ[1]==0 ? - float(x .+ 1).^(- 2) : (λ[1] .- 1) * float(x .+ 1) .^ (λ[1] .- 2)

struct YeoJohnson <:NonlinearTransform end
function (::YeoJohnson)(x, λ::Union{Array{T, 1}, T}) where T<:Real
    @assert (typeof(λ) <:Array{T} where T<:Real && size(λ, 2)==1) || (typeof(λ)<:Real)
    if x >= 0
        g = BoxCox()
        return g(x+1,  λ)
    elseif λ == 2
        return -log.(1 .- x)
    else
        @assert (x<0) && (λ!=2)
        return (expm1.(log.(1-x) .* (2-λ))) /(λ-2)
    end
end
function partialx(::YeoJohnson, x, λ)
    if x >= 0
        g = BoxCox()
        return partialx(g, x+1, λ)
    elseif λ == 2
        return 1 ./ (1 .- float(x))
    else
        @assert (x<0) && (λ!=2)
        (1-float(x)).^(1-λ)
    end
end

struct IdentityTransform <:NonlinearTransform end
(::IdentityTransform)(x, λ::Union{Array{T, 1}, T} where T<:Real) = (@assert (typeof(λ) <:Array{T} where T<:Real && size(λ, 2)==1) || (typeof(λ)<:Real);
                                                            identity.(x))
(::IdentityTransform)(x) = (identity.(x))

partialx(::IdentityTransform, x, λ) = 1
partialxx(::IdentityTransform, x, λ)  = 0
partialλ(::IdentityTransform, x, λ) = 0
partialλλ(::IdentityTransform, x, λ) = 0
inverse(::IdentityTransform, x, λ) = identity.(x)

max_stars_repo_name: dbindel/btg

max_issues_repo_name: dbindel/btg

max_forks_repo_name: dbindel/btg


text:
module TestTwoStage
using Test
using CompressedSensing: gaussian_data, perturb, srr, sp, ompr, sparse_vector

@testset "Two-Stage Algorithms" begin
    n, m, k = 32, 64, 3
    A, x, b = gaussian_data(n, m, k)
    δ = 1e-2
    y = perturb(b, δ/2)

    @testset "Stagewise Regression with Replacement" begin
        # noiseless
        xsrr = srr(A, b, k)
        @test xsrr.nzind == x.nzind
        @test xsrr.nzval ≈ x.nzval

        # noisy
        xsrr = srr(A, y, k)
        @test xsrr.nzind == x.nzind
        @test isapprox(xsrr.nzval, x.nzval, atol = 3δ)

        # special case
        x1 = sparse_vector(m, 1)
        xsrr = srr(A, A*x1, 1)
        @test xsrr.nzind == x1.nzind
        @test xsrr.nzval ≈ x1.nzval

        # l-step version of srr
        l = k
        # noiseless
        xsrr = srr(A, b, k, l = k)
        @test xsrr.nzind == x.nzind
        @test xsrr.nzval ≈ x.nzval

        # noisy
        xsrr = srr(A, y, k, l = k)
        @test xsrr.nzind == x.nzind
        @test isapprox(xsrr.nzval, x.nzval, atol = 3δ)

    end

    @testset "Subspace Pursuit" begin
        xsp = sp(A, b, k)
        # noiseless
        @test xsp.nzind == x.nzind
        @test xsp.nzval ≈ x.nzval

        # noisy
        xssp = sp(A, y, k, δ)
        @test xssp.nzind == x.nzind
        @test isapprox(xssp.nzval, x.nzval, atol = 3δ)
    end

    @testset "OMP with replacement" begin
        # noiseless
        xompr = ompr(A, b, k, 1e-6)
        @test xompr.nzind == x.nzind
        @test xompr.nzval ≈ x.nzval

        # slightly noisy
        xomp = ompr(A, y, k, δ)
        @test xomp.nzind == x.nzind
        @test isapprox(xomp.nzval, x.nzval, atol = 3δ)
    end
end
end # TestTwoStage


text:
"""
    Poisson{N,M}

Composite type for conservative variable coefficient Poisson equations:

    ∮ds β ∂x/∂n = σ

The resulting linear system is

    Ax = [L+D+L']x = b

where A is symmetric, block-tridiagonal and extremely sparse. Implemented on a
structured grid of dimension N, then L has dimension M=N+1 and size(L,M)=N.
Moreover, D[I]=-∑ᵢ(L[I,i]+L'[I,i]). This means matrix storage, multiplication,
ect can be easily implemented and optimized without external libraries.

To help iteratively solve the system above, the Poisson structure holds
helper arrays for inv(D), the error ϵ, and residual r=b-Ax. An iterative
solution method then estimates the error ϵ=̃A⁻¹r and increments x+=ϵ, r-=Aϵ.
"""
abstract type AbstractPoisson{N,M} end
struct Poisson{N,M} <: AbstractPoisson{N,M}
    L :: Array{Float64,M} # Lower diagonal coefficients
    D :: Array{Float64,N} # Diagonal coefficients
    iD :: Array{Float64,N} # 1/Diagonal
    x :: Array{Float64,N} # approximate solution
    ϵ :: Array{Float64,N} # increment/error
    r :: Array{Float64,N} # residual
    n :: Vector{Int16}    # pressure solver iterations
    function Poisson(L::Array{Float64,m}) where m
        M = size(L); N = M[1:end-1]; n = m-1
        @assert M[end] == n
        x,ϵ,r,D,iD = zeros(N),zeros(N),zeros(N),zeros(N),zeros(N)
        set_diag!(D,iD,L)
        new{n,m}(L,D,iD,x,ϵ,r,[])
    end
end
function set_diag!(D::Array{T,n},iD,L) where {T,n}
    @inbounds @simd for I ∈ inside(D)
        D[I] = -sum(@inbounds(L[I,i]+L[I+δ(i,n),i]) for i ∈ 1:n)
        iD[I] = abs2(D[I])<1e-8 ? 0. : inv(D[I])
    end
end
set_diag!(p::Poisson) = set_diag!(p.D,p.iD,p.L)
update!(p::Poisson,L::Array) = (p.L .= L; set_diag!(p))

@fastmath @inline multL(I::CartesianIndex{d},L,x) where {d} =
    sum(@inbounds(x[I-δ(a,I)]*L[I,a]) for a ∈ 1:d)
@fastmath @inline multU(I::CartesianIndex{d},L,x) where {d} =
    sum(@inbounds(x[I+δ(a,I)]*L[I+δ(a,I),a]) for a ∈ 1:d)
@fastmath @inline mult(I,L,D,x) = @inbounds(x[I]*D[I])+multL(I,L,x)+multU(I,L,x)
function mult(p::Poisson{n},x::Array{Float64,n}) where n
    @assert size(p.x)==size(x)
    b = zeros(size(p.x))
    @inside b[I] = mult(I,p.L,p.D,x)
    return b
end

@fastmath residual!(p::Poisson,b::Array{Float64}) =
    @inside p.r[I] = b[I]-mult(I,p.L,p.D,p.x)

@fastmath increment!(p::Poisson) = @inbounds @simd for I ∈ inside(p.x)
    p.x[I] += p.ϵ[I]
    p.r[I] -= mult(I,p.L,p.D,p.ϵ)
end
"""
SOR!(p::Poisson;ω=1.5)

Successive Over Relaxation preconditioner. The routine uses backsubstitution
to compute ϵ=̃A⁻¹r, where ̃A=[D/ω+L], and then increments x,r.
"""
@fastmath function SOR!(p::Poisson; ω=1.5)
    @inbounds for I ∈ inside(p.r) # order matters here
        σ = p.r[I]-multL(I,p.L,p.ϵ)
        p.ϵ[I] = ω*σ*p.iD[I]
    end
    increment!(p)
end
"""
GS!(p::Poisson;it=0)

Gauss-Sidel smoother. When it=0, the function serves as a Jacobi preconditioner.
"""
@fastmath function GS!(p::Poisson;it=0)
    @inside p.ϵ[I] = p.r[I]*p.iD[I]
    for i ∈ 1:it; @inbounds for I ∈ inside(p.r, reverse = i%2==0) # order matters here
        σ = p.r[I]-multL(I,p.L,p.ϵ)-multU(I,p.L,p.ϵ)
        p.ϵ[I] = σ*p.iD[I]
    end;end
    increment!(p)
end

function solver!(x::Array{Float64,n},p::Poisson{n},b::Array{Float64,n};log=false,tol=1e-4,itmx=1e3) where n
    p.x .= x
    residual!(p,b); r₂ = L₂(p.r)
    log && (res = [r₂])
    nᵖ=0
    while r₂>tol && nᵖ<itmx
        SOR!(p,ω=1.8); r₂ = L₂(p.r)
        log && push!(res,r₂)
        nᵖ+=1
    end
    x .= p.x
    push!(p.n,nᵖ)
    log && return res
end


text:
function prepend_prefix!(v::Vector{<:Line},y)
    for x in v
        prepend!(x.prefix.prefix, y)
        filter!(t -> !in(value(t),[ "\n", "—" ]), x.tokens)
    end
    filter!(t -> !isempty(t.tokens), v)
end
trimstring(x) = x
word_string(w) =
    trimstring(
        join(string.(filter(t-> !(t isa TokenPair && t.key==:paren), w[1].tokens))))

function parse_overview(namespace, title, w, t::Nothing)
    images = Line{NamedString,LineContent}[]
    inflections = Pair{String,Token}[]
    language, wordtype, genus = string((filter(t-> t isa Template,w[1].tokens)[end].arguments[1].second)), "",""
    ( # word = Token(:literal, intern(word_string(w))),
      title = Token(namespace, intern(title)),
      language = Token(:language,language),
      wordtype = Token(:wordtype,wordtype),
      inflections = inflections,
      images = images,
      wikitext = [ w[1], w[2]... ] )
end

image_argument_parser =
    Sequence(!!re"Bild[^[:digit:]]*",
             !!re"[[:digit:]]*",
             AtEnd()) do v
                 NamedString[ NamedString("image", v[2]),
                              NamedString("property", v[1]) ]
             end

function parse_overview(namespace, title, w, t::Template)
    images = Line{NamedString,LineContent}[]
    inflections = Pair{String, Token}[]
    overview_parser =
        Sequence(whitespace_maybe, !!re"[-\p{L}]+", whitespace_horizontal,
                 !!Optional(re".*[^ ]"),
                 Optional(re" *Übersicht *"),
                 !!Optional(re".*[^ ]"),
                 whitespace_maybe,
                 Repeat(newline)) do v
                     String[v[2],v[4],v[6]]
                 end

    language, wordtype, genus = parse(overview_parser, t.template)
    args = Pair{String,AbstractToken}[]
    for (i,a) in enumerate(t.arguments)
        if !isempty(a.second) && !isempty(a.second[1].tokens)
            val = a.second[1].tokens[1]
            imp = tryparse(image_argument_parser, a.first)
            if imp !== nothing
                append!(images, prepend_prefix!(a.second, imp))
            elseif a.first==""
                ## push!(args, "$i" => val)
            elseif a.first in ["Genus"]
                push!(args, a.first => val)
            else
                push!(inflections, a.first => Token(:literal,string(val)))
            end
        end
    end
    ( # word = Token(:literal, intern(intern(word_string(w)))),
      title = Token(namespace, intern(title)),
      language = Token(:language,language),
      wordtype = Token(:wordtype,wordtype),
      inflections = inflections,
      images = images,
      wikitext = [ w[1], w[2]..., Line(LineContent[t]) ] )
end

import CombinedParserTools
import CombinedParserTools.Tokens
function parse_overview(namespace, title, w, v::Vector{<:Line})
    parse_overview(namespace, title, w, isempty(v) ? nothing : v[1].tokens[1])
end

is_template(template::String) =
    IteratorParser(Template,template) do x
        x.name==template
    end
is_template_line(x::String,f::Function,T::Type) =
    map(f,is_template_line(x))
is_template_line(template::String) =
    IteratorParser(Line{NamedString,LineContent},template) do x
        !isempty(x.tokens) &&
            x.tokens[1] isa Template &&
            x.tokens[1].template==template
    end
is_template_line(pred::Function, T = Line{NamedString,LineContent}) =
    IteratorParser{T}("template") do x
        (x) isa Line &&
            !isempty(x.tokens) &&
            (x.tokens[1]) isa Template &&
            pred(x.tokens[1])
    end
is_template_line(T::Type, pred::Function) =
    is_template_line(pred, T)

ignore = (v) -> missing
ignore_heading = (v) -> begin
    length(v) == 1 ? [] : v[2:end]
end

wiktionary_de_content=[
    is_template_line("Übersetzungen", ignore, Missing) => (:translations, is_line()),
    is_template_line(t -> match(re"Übersicht", t.template)!==nothing) => (:overview, missing),
    is_template_line("Nebenformen", ignore, Missing) => (:variants, is_line()),
    is_template_line("Entlehnungen", ignore, Missing) => (:variants, is_line()),
    is_template_line("Nicht mehr gültige Schreibweisen") => (:deprecated, is_line()),
    is_template_line("Alternative Schreibweisen") => (:versions, is_line()),
    is_template_line("Anmerkungen", ignore, Missing) => (:annotations, is_line()),
    is_template_line("Anmerkung", ignore, Missing) => (:annotations, is_line()),
    is_template_line("Worttrennung", ignore, Missing) => (:hyphenation, is_line()),
    is_template_line("Aussprache", ignore, Missing) => (:phonetic, is_line()),
    is_template_line("Oberbegriffe", ignore, Missing) => (:superterms, is_line()),
    is_template_line("Unterbegriffe", ignore, Missing) => (:subterms, is_line()),
    is_template_line("Teilbegriffe", ignore, Missing) => (:parts, is_line()),
    is_template_line("Bedeutungen", ignore, Missing) => (:meaning, is_line()),
    is_template_line("Synonyme", ignore, Missing) => (:synonyms, is_line()),
    is_template_line("Weibliche Wortformen", ignore, Missing) => (:synonyms, is_line()),
    is_template_line("Männliche Wortformen", ignore, Missing) => (:synonyms, is_line()), # exists?
    is_template_line("Sinnverwandte Wörter", ignore, Missing) => (:synonyms, is_line()),
    is_template_line("Wortbildungen", ignore, Missing) => (:formations, is_line()),
    is_template_line("Gegenwörter", ignore, Missing) => (:antonyms, is_line()),
    is_template_line("Abkürzungen", ignore, Missing) => (:abbreviations, is_line()),
    is_template_line("Wortfamilie", ignore, Missing) => (:family, is_line()),
    is_template_line("Grammatische Merkmale", ignore, Missing) => (:grammar, is_line()),
    is_template_line("Redewendungen", ignore, Missing) => (:phrases, is_line()),
    is_template_line("Siehe auch", ignore, Missing) => (:etymology, is_line()),
    is_template_line("erweitern", ignore, Missing) => (:etymology, is_line()),
    is_template_line("Herkunft", ignore, Missing) => (:etymology, is_line()),
    is_template_line("Beispiele", ignore, Missing) => (:examples, is_line()),
    is_template_line("Beispiele fehlen") => (:examples, is_line()),
    is_template_line("erweitern") => (:examples, is_line()),
    is_template_line("Charakteristische Wortkombinationen", ignore, Missing) => (:combinations, is_line()),
    is_template_line("Referenzen", ignore, Missing) => (:references, is_line()),
    is_template_line("Quellen", ignore, Missing) => (:sources, is_line())
]

export wiktionary_defs
wiktionary_defs =
    Repeat(Sequence(
        :prefix => Repeat(is_line()),
        :word => Sequence(is_heading(x -> isequal(x.value,"2")),
                          Repeat(is_line())),
        :defs => Repeat(Sequence(is_heading(x -> isequal(x.value,"3")),
                                 greedy(wiktionary_de_content...; alt=[:lines => is_line() ]))
                        )
    ));

number_line = Sequence(
    map(value,IteratorParser(x->variable(x)==:number, Token, "index")),
    Repeat(IteratorParser(!isinformative, Token, "whitespace")),
    Repeat(IteratorParser(LineContent, "rest"))) do v
        v[1] => v[3]
    end

using Epitypes
export wiki_meaning
function wiki_meaning(title,v;namespace = "wikt:de")
    L = Line{NamedString,LineContent}
    getval(val,p) =  p => ( ( haskey(val,p) && val[p]!==missing ) ? val[p] : L[] )
    fields = [ x.second[1] for x in wiktionary_de_content ]
    function inner(wt)
        meaning_data = Dict{String,Dict{Symbol,Vector{L}}}()
        function push_field(num,f,val)
            meandata = get!(meaning_data, num) do
                Dict{Symbol,Vector{Line{NamedString,LineContent}}}()
            end
            data = get!(() -> Line{NamedString,LineContent}[], meandata, f)
            val !== nothing && push!(data, Line(NamedString[], val))
        end
        for (f,v) in wt[2]
            x = parse(
                Repeat(map(is_line()) do l
                       r = tryparse(number_line, l.tokens)
                       r === nothing ? "?" => l.tokens : r
                       end),
                v)
            lastnum = "?"
            for e in x
                num = e.first=="?" ? lastnum : e.first
                let is = tryparse(expand_numbers,InternedStrings.intern(num))
                    if is === nothing
                        push_field(num,f,e.second)
                    else
                        for i in is
                            push_field(i,f,e.second)
                        end
                    end
                end
                lastnum=num
            end
        end
        ## fill meaning_data values with all section defs
        for f in fields
            for num in keys(meaning_data)
                push_field(num,f,nothing)
            end
        end
        common = get(()->Dict{Symbol,Vector{L}}(), meaning_data,"[*]")
        for (m,val) in pairs(meaning_data)
            if !in(m, ["?","[*]"])
                for (f,v) in pairs(common)
                    if haskey(val,f)
                        append!(val[f],v)
                    else
                        val[f]=v
                    end
                end
            end
        end
        val = get(meaning_data,"?",Dict())
        base = parse_overview(namespace, title, v.word, getval(val,:overview).second)
        meanings = [
            NamedStruct{:meaning}(
                ( title=base.title,
                  # , word=base.word,
                  order=m,
                  ( getval(val,p)
                    for p in fields)... ))
            for (m,val) in pairs(meaning_data)
            if !in(m, ["?","[*]"])
        ]

        NamedStruct{:word}(
            merge(
                base,
                (;( getval(val,p)
                    for p in fields
                    if p != :overview)...
                 , meanings = sort(meanings; by=x->x.order)
                 ))
        )
    end
    Any[ inner(wt) for  wt = v.defs ]
end


text:
using Test, PreallocationTools, ForwardDiff, LinearAlgebra

#test for downsizing cache
randmat = rand(5, 3)
sto = similar(randmat)
stod = dualcache(sto)

function claytonsample!(sto, τ, α; randmat=randmat)
    sto = get_tmp(sto, τ)
    sto .= randmat
    τ == 0 && return sto

    n = size(sto, 1)
    for i in 1:n
        v = sto[i, 2]
        u = sto[i, 1]
        sto[i, 1] = (1 - u^(-τ) + u^(-τ)*v^(-(τ/(1 + τ))))^(-1/τ)*α
        sto[i, 2] = (1 - u^(-τ) + u^(-τ)*v^(-(τ/(1 + τ))))^(-1/τ)
    end
    return sto
end

#taking the derivative of claytonsample! with respect to τ only
df1 = ForwardDiff.derivative(τ -> claytonsample!(stod, τ, 0.0), 0.3)
@test size(randmat) == size(df1)

#calculating the jacobian of claytonsample! with respect to τ and α
df2 = ForwardDiff.jacobian(x -> claytonsample!(stod, x[1], x[2]), [0.3; 0.0]) #should give a 15x2 array,
#because ForwardDiff flattens the output of jacobian, see: https://juliadiff.org/ForwardDiff.jl/stable/user/api/#ForwardDiff.jacobian

@test (length(randmat), 2) == size(df2)
@test df1[1:5,2] ≈ df2[6:10,1]

#test for enlarging cache
function rhs!(du, u, p, t)
    A = p
    mul!(du, A, u)
end

function loss(du, u, p, t)
    _du = get_tmp(du, p)
    rhs!(_du, u, p, t)
    l = _du[1]
    return l
end

u = [3.0, 0.0]
A = ones(2,2)

du = similar(u)
_du = dualcache(du)
f = A -> loss(_du, u, A, 0.0)
analyticalsolution = [3.0 0; 0 0]
@test ForwardDiff.gradient(f, A) ≈ analyticalsolution


text:
template_plugins(::PDTN) = nothing

function create_template(t::PDTN{name}, dir, user, interactive) where name
    plugins = template_plugins(t)
    if !isnothing(plugins)
        if isempty(user)
            @warn "git host user name not found, please specify " *
                "using --user to enable GitHub related plugins " *
                "or use -i/--interactive to configure interactively"
            return Template(;
                dir=dir,
                interactive=interactive,
                plugins = [plugins..., !Git, !TagBot, !CompatHelper],
            )
        else
            return Template(;
                dir=dir,
                user=user,
                interactive=interactive,
                plugins = plugins,
            )
        end
    else
        template = search_local_template(t)
        isnothing(template) && error("template $(name) not found")
        return template
    end
end

function search_local_template(::PDTN{name}) where name
    template_dir = PATH.dot_ion("templates")
    ispath(template_dir) || return
    template = string(name)
    for dir in readdir(template_dir)
        if dir == template
            error("read/save local template is not supported yet")
            # return read_template(joinpath(template_dir, dir))
        end
    end
    return
end

# function read_template(path)
#     file = joinpath(path, "README.md")
#     plugins = []
#     if isfile(file)
#         push!(plugins, Readme(;file=file))
#     end

#     file = joinpath(path, "Project.toml")

# end


text:
# Note: Test environment is setup by https://github.com/JuliaCI/PkgEvalHadoopEnv
include("hdfstests.jl")
include("yarntests.jl")

test_hdfs()
test_yarn()


text:
module JungleHelperSnake

using ..Ahorn, Maple
@pardef Snake(x1::Integer, y1::Integer, x2::Integer=x1, y2::Integer=y1+ 16, left::Bool=false, sprite::String="") = Entity("JungleHelper/Snake", x=x1, y=y1, nodes=Tuple{Int, Int}[(x2, y2)], left=left, sprite=sprite)

const placements = Ahorn.PlacementDict(
    "Snake (Jungle Helper)" => Ahorn.EntityPlacement(
        Snake,
        "rectangle",
        Dict{String, Any}(),
        function(entity)
            x, y = Int(entity.data["x"]), Int(entity.data["y"])
            entity.data["x"], entity.data["y"] = x + 64, y
            entity.data["nodes"] = [(x, y)]
        end
    ),
    "Snake (Facing Left) (Jungle Helper)" => Ahorn.EntityPlacement(
        Snake,
        "rectangle",
        Dict{String, Any}(
            "left" => true,
        ),
        function(entity)
            x, y = Int(entity.data["x"]), Int(entity.data["y"])
            entity.data["x"], entity.data["y"] = x + 64, y
            entity.data["nodes"] = [(x, y)]
        end
    )
)

Ahorn.nodeLimits(entity::Snake) = 1, 1

const sprite = "JungleHelper/Snake/IdleAggro/snake_idle00"

function Ahorn.selection(entity::Snake)
    x, y = Ahorn.position(entity)
    nx, ny = Int.(entity.nodes[1])

    return Ahorn.Rectangle[Ahorn.Rectangle(x, y + 8, 64, 8), Ahorn.Rectangle(nx, ny + 8, 64, 8)]
end

function Ahorn.renderSelectedAbs(ctx::Ahorn.Cairo.CairoContext, entity::Snake)
    px, py = Ahorn.position(entity)
    nx, ny = Int.(entity.nodes[1])

    if get(entity.data, "left", false)
        Ahorn.drawSprite(ctx, sprite, nx + 64, ny, sx=-1, jx=0, jy=0)
    else
        Ahorn.drawSprite(ctx, sprite, nx, ny, jx=0, jy=0)
    end

	Ahorn.drawArrow(ctx, px + 32, py + 12, nx + 32, ny + 12, Ahorn.colors.selection_selected_fc, headLength=6)
end

function Ahorn.renderAbs(ctx::Ahorn.Cairo.CairoContext, entity::Snake, room::Maple.Room)
    x, y = Ahorn.position(entity)

    if get(entity.data, "left", false)
        Ahorn.drawSprite(ctx, sprite, x + 64, y, sx=-1, jx=0, jy=0)
    else
        Ahorn.drawSprite(ctx, sprite, x, y, jx=0, jy=0)
    end
end

end

max_stars_repo_name: Jackal-Celeste/JungleHelper

max_issues_repo_name: Jackal-Celeste/JungleHelper

max_forks_repo_name: Jackal-Celeste/JungleHelper


text:
# This implements the DGM operator for Diversity-Guided Mutation as described in
#
#  M. S. Alam, Md. M. Islam, X. Yao, K. Murase, "Diversity Guided Evolutionary
#  Programming: A novel approach for continuous optimization", ASOC, 2013
#

# DGM operator - this is a slight variation from the one in their paper
# were we do not need to set a min distance but rather select the min and non-min
# sets by randomly sampling a subset of individuals in the population.
#
# Input:
#  - population is a matrix of Float64 where each column represents one individual in the population
#  - parent is an Int64 index that points out one column in the population which will act as the parent to be mutated
#  - population_size
# Output:
#  - a child that is a diversity-mutated variant of the parent
#
function diversity_guided_mutation(population::Array{Float64, 2}, parent::Int64;
  sample_size = 10,
  indices = collect(1:size(population, 2)))

  p = population[:,parent]
  population_size = size(population, 2)
  shuffle!(indices)

  # Calculate the distances from sample_size randomly sampled individuals
  # to the parent.
  distances = zeros(sample_size)
  subset = zeros(Int64, sample_size)
  di = i = 1
  while di <= sample_size
    si = indices[i]
    if si == parent
      i += 1
      si = indices[i]
    end
    distances[di] = norm(population[:,si] - p)
    di += 1
    i += 1
  end

  # Sort the distances in the subset so we can have the low 3rd as Nx and the
  # high 3rd as Pop-Nx.
  perm = sortperm(distances)
  third_size = round(sample_size/3)
  subset = subset[perm]
  nx = subset[1:third_size]
  pt_minus_nx = subset[]


text:
using Docile, Lexicon, QuantEcon

const api_directory = "api"
const modules = [QuantEcon, QuantEcon.Models]

# const nb_src_dir = joinpath(dirname(dirname(@__FILE__)), "solutions")
# const nb_output_dir = joinpath(dirname(dirname(@__FILE__)), "solutions")
# const nb_list = filter(x->endswith(x, "ipynb"), readdir(nb_src_dir))

cd(dirname(@__FILE__)) do
    # Generate and save the contents of docstrings as markdown files.
    index  = Index()
    for mod in modules
        Lexicon.update!(index, save(joinpath(api_directory, "$(mod).md"), mod))
    end
    save(joinpath(api_directory, "index.md"), index; md_subheader = :category)

    # Add a reminder not to edit the generated files.
    open(joinpath(api_directory, "README.md"), "w") do f
        print(f, """
        Files in this directory are generated using the `build.jl` script. Make
        all changes to the originating docstrings/files rather than these ones.
        """)
    end

    # for nb in nb_list
    #     nb_path = joinpath(nb_src_dir, nb)
    #     md_path = joinpath(nb_output_dir, replace(nb, ".ipynb", ".md"))
    #     run(`ipython nbconvert --to markdown $nb_path --output=$md_path`)
    # end

    # info("Adding all documentation changes in $(api_directory) to this commit.")
    # success(`git add $(api_directory)`) || exit(1)

end

cd(dirname(dirname(@__FILE__))) do
    yaml = """
    # This is automatically generated by docs/build.jl. Edit that file if you
    # wish to make any permenant changes.
    site_name: QuantEcon.jl
    site_description: QuantEcon.jl Quantitative economics in Julia
    site_author: Spencer Lyon
    repo_name: GitHub
    site_favicon: favicon.ico
    docs_dir: 'docs'
    site_dir: 'site'
    repo_url: https://github.com/QuantEcon/QuantEcon.jl/
    google_analytics: ['UA-64434120-1', 'auto']

    pages:
    - Home: 'index.md'
    - Overview: 'api/index.md'
    - API Docs:
      - QuantEcon: 'api/QuantEcon.md'
      - Models: 'api/QuantEcon.Models.md'
    """

    # TODO: add the solutions if I figure out how to get them to render properly

    open("mkdocs.yml", "w") do f
        write(f, yaml)
    end
end


text:
module PolicyViz

using GridInterpolations, Interact, PGFPlots, Colors, ColorBrewer, HDF5, SparseArrays, Printf

include("./viz_policy_constants.jl")
include("./viz_policy.jl")
include("./nnet_calculations.jl")

end

max_stars_repo_name: cherrywoods/HorizontalCAS

max_issues_repo_name: ahmed-irfan/HorizontalCAS

max_forks_repo_name: ahmed-irfan/HorizontalCAS


text:
### Day 4 - Problem 1 & 2
## Security Through Obscurity
## Author: Thanasis Georgiou

# Read input
f = open("rooms.txt")
input = readlines(f)
close(f)

# Verify room
function getRoomID(room::String)
    # Remove newline
    room = chomp(room)

    # Get checksum
    checksum = room[end - 5:end - 1]
    room = room[1:end - 7]

    # Get ID
    id = room[end - 2:end]
    room = room[1:end - 4]

    # Count each letter
    letters = Dict(map(x -> Pair(x, 0), unique(room)))
    foreach(x -> letters[x] += 1, room)
    delete!(letters, '-') # Drop dash

    # Sort most common letters first by how many times
    # they appear and then alphabetically by key
    comp(a, b) = begin
        if last(a) > last(b)
            return false
        elseif last(a) == last(b)
            return first(a) > first(b)
        else
            return true
        end
    end
    trueChecksumArr = sort(collect(letters), lt = comp, rev = true)
    trueChecksum = ""
    for kv in trueChecksumArr[1:5]
        trueChecksum = "$(trueChecksum)$(kv[1])"
    end

    # Verify checksum
    if trueChecksum == checksum
        # Grab ID
        roomId = parse(Int64, id)

        # Decrypt name
        name = map(room) do c
            if 'a' <= c <= 'z'
                (c - 'a' + roomId) % 26 + 'a'
            elseif c == '-'
                ' '
            end
        end

        if contains(name, "pole")
            println("North Pole Objects: $roomId")
        end

        # Return room ID
        return roomId
    else
        return 0
    end
end

# Get room IDs
result = sum(map(getRoomID, input))


text:
loadbutton = filepicker()
columnbuttons = Observable{Any}(dom"div"())
data = Observable{Any}(Sequence)
#Assigning function of data when load button (filepicker) is changed
map!(f->begin
        global seq = JLD2.load(FileIO.File{FileIO.DataFormat{:JLD2}}(f),"seq")
        print("Loaded! ... $f\n")
    end
    , data, loadbutton)

p = MRIsim.plot_seq(seq)
plt = Observable{Any}(p)

function makebuttons(seq)
    namesseq = ["Sequence","k-space"]
    buttons = button.(string.(namesseq))
    for (btn, name) in zip(buttons, namesseq)
        if name == "Sequence"
            map!(t-> begin
                MRIsim.plot_seq(seq)
                end
                , plt, btn)
        elseif name == "k-space"
            hh, ww = 600, 600
            l = PlotlyJS.Layout(;
                scene=attr(xaxis_title=L"kx [m^-1]",
                        yaxis_title="ky [m^-1]",
                        zaxis_title="kz [m^-1]",),
                modebar=attr(orientation="v"),height=hh,width=ww,hovermode="closest",
                scene_camera_eye=attr(x=0, y=0, z=2),
                scene_camera_up=attr(x=0, y=1., z=0),)
            map!(t->begin
                ACQ = is_DAC_on.(seq)
                #Adding gradients before and after acq, to show the centered k-space
                ACQ = (circshift(ACQ,1)+circshift(ACQ,-1)+ACQ).!=0 #TODO: Change to accumulated moment-0
                seqADC = sum(seq[ACQ])
                kspace = MRIsim.get_designed_kspace(seqADC)
                N = size(kspace,1)-1
                M = size(kspace,2)
                c = ["hsv($((i-1)/N*255),255,25)" for i=1:N-1]
                lines = PlotlyJS.scatter3d(x=kspace[:,1],y=kspace[:,2],z=kspace[:,3],mode="lines",
                        line=attr(color=c),aspectratio=1)
                PlotlyJS.plot(lines,l)
            end
            , plt, btn)
        end
    end
    dom"div"(hbox(buttons))
end

map!(makebuttons, columnbuttons, data)

columnbuttons = makebuttons(seq)
ui = dom"div"(loadbutton, columnbuttons, plt)
content!(w, "div#content", ui)


text:
abstract type AbstractNamedDimDataGraph{VD,ED,V,E} <: AbstractDataGraph{VD,ED,V,E} end

#
# Specializations for AbstractNamedDimDataGraph
#

## function IndexType(graph::AbstractNamedDimDataGraph, v_or_e)
##   (tuple_convert(v_or_e) isa eltype(graph)) && return VertexIndex()
##   return error("$v_or_e doesn't represent a vertex or an edge for graph:\n$graph.")
## end

ElementIndex_to_VertexIndex(x) = x
ElementIndex_to_VertexIndex(::ElementIndex) = VertexIndex()

IndexType(::AbstractNamedDimDataGraph, I...) = ElementIndex_to_VertexIndex(IndexType(I...))

# Ambiguity error with AbstractGraph version
IndexType(::AbstractNamedDimDataGraph{VD,ED,V}, ::V) where {VD,ED,V} = VertexIndex()

# Ambiguity error with AbstractGraph version
IndexType(::AbstractNamedDimDataGraph, ::AbstractEdge) = EdgeIndex()

# Ambiguity error with AbstractGraph version
IndexType(::AbstractNamedDimDataGraph, ::Pair) = EdgeIndex()

function getindex(::SliceIndex, graph::AbstractNamedDimDataGraph, index...)
  underlying_subgraph = getindex(underlying_graph(graph), index...)
  subvertices = vertices(underlying_subgraph)
  subvertex_data = vertex_data(graph)[subvertices]
  subedge_data_indices = filter(
    e -> src(e) ∈ subvertices && dst(e) ∈ subvertices, keys(edge_data(graph))
  )
  subedge_data = getindices(edge_data(graph), subedge_data_indices)
  return typeof(graph)(underlying_subgraph, subvertex_data, subedge_data)
end

# TODO: define VertexNamedDimDataGraph, a graph with only data on the
# vertices, and EdgeNamedDimDataGraph, a graph with only data on the edges.
struct NamedDimDataGraph{VD,ED,V,E,G<:AbstractGraph} <: AbstractNamedDimDataGraph{VD,ED,V,E}
  underlying_graph::G
  vertex_data::MultiDimDictionary{V,VD}
  edge_data::Dictionary{E,ED}
end
underlying_graph(graph::NamedDimDataGraph) = graph.underlying_graph
vertex_data(graph::NamedDimDataGraph) = graph.vertex_data
edge_data(graph::NamedDimDataGraph) = graph.edge_data

function is_directed(::Type{<:NamedDimDataGraph{VD,ED,V,E,G}}) where {VD,ED,V,E,G}
  return is_directed(G)
end

function copy(graph::NamedDimDataGraph)
  # Need to use deepcopy of Dictionaries, see:
  # https://github.com/andyferris/Dictionaries.jl/issues/98
  return NamedDimDataGraph(
    copy(underlying_graph(graph)), deepcopy(vertex_data(graph)), deepcopy(edge_data(graph))
  )
end

function NamedDimDataGraph{VD,ED}(
  underlying_graph::AbstractGraph,
  vertex_data::MultiDimDictionary{V,VD},
  edge_data::Dictionary{E,ED},
) where {VD,ED,V,E}
  G = typeof(underlying_graph)
  return NamedDimDataGraph{VD,ED,V,E,G}(underlying_graph, vertex_data, edge_data)
end

function NamedDimDataGraph{VD,ED}() where {VD,ED}
  return NamedDimDataGraph{VD,ED}(NamedDimGraph())
end

function NamedDimDataGraph{VD}() where {VD}
  return NamedDimDataGraph{VD,Any}(NamedDimGraph())
end

function NamedDimDataGraph()
  return NamedDimDataGraph{Any,Any}(NamedDimGraph())
end

function NamedDimDataGraph{VD,ED}(underlying_graph::NamedDimGraph) where {VD,ED}
  V = eltype(underlying_graph)
  E = edgetype(underlying_graph)
  vertex_data = MultiDimDictionary{V,VD}()
  edge_data = Dictionary{E,ED}()
  return NamedDimDataGraph{VD,ED}(underlying_graph, vertex_data, edge_data)
end

function NamedDimDataGraph{VD}(underlying_graph::NamedDimGraph) where {VD}
  return NamedDimDataGraph{VD,Any}(underlying_graph)
end

function NamedDimDataGraph(underlying_graph::NamedDimGraph)
  return NamedDimDataGraph{Any,Any}(underlying_graph)
end

function NamedDimDataGraph{VD,ED}(underlying_graph::Graph; kwargs...) where {VD,ED}
  return NamedDimDataGraph{VD,ED}(NamedDimGraph(underlying_graph; kwargs...))
end

function NamedDimDataGraph{VD}(underlying_graph::Graph; kwargs...) where {VD}
  return NamedDimDataGraph{VD}(NamedDimGraph(underlying_graph; kwargs...))
end

function NamedDimDataGraph(underlying_graph::Graph; kwargs...)
  return NamedDimDataGraph(NamedDimGraph(underlying_graph; kwargs...))
end

function hvncat(
  dim::Int, graph1::NamedDimDataGraph, graph2::NamedDimDataGraph; new_dim_names=(1, 2)
)
  # Concatenate the underlying graph
  new_underlying_graph = hvncat(
    dim, underlying_graph(graph1), underlying_graph(graph2); new_dim_names
  )

  # The new vertices of the graph.
  # This can introduce new dimensions, shift values, etc.
  new_vertices = vertices(new_underlying_graph)

  # Map the vertices of the vertex data to the new vertices (for the first input graph)
  map_vertices1 = Dictionary(vertices(graph1), new_vertices[1:nv(graph1)])
  new_vertices1 = getindices(map_vertices1, keys(vertex_data(graph1)))
  new_vertex_data1 = MultiDimDictionary(new_vertices1, vertex_data(graph1))

  # Map the vertices of the vertex data to the new vertices (for the second input graph)
  map_vertices2 = Dictionary(vertices(graph2), new_vertices[(nv(graph1) + 1):end])
  new_vertices2 = getindices(map_vertices2, keys(vertex_data(graph2)))
  new_vertex_data2 = MultiDimDictionary(new_vertices2, vertex_data(graph2))

  new_vertex_data = merge(new_vertex_data1, new_vertex_data2)

  # The new edges of the graph.
  # This can introduce new dimensions, shift values, etc.
  new_edges = edges(new_underlying_graph)

  # Map the edges of the edge data to the new vertices (for the first input graph)
  map_edges1_out = Dictionary(edges(graph1), new_edges[1:ne(graph1)])
  map_edges1_in = Dictionary(reverse.(keys(map_edges1_out)), reverse.(map_edges1_out))
  map_edges1 = merge(map_edges1_out, map_edges1_in)
  new_edges1 = getindices(map_edges1, keys(edge_data(graph1)))
  new_edge_data1 = Dictionary(new_edges1, edge_data(graph1))

  # Map the edges of the edge data to the new vertices (for the second input graph)
  map_edges2_out = Dictionary(edges(graph2), new_edges[(ne(graph1) + 1):end])
  map_edges2_in = Dictionary(reverse.(keys(map_edges2_out)), reverse.(map_edges2_out))
  map_edges2 = merge(map_edges2_out, map_edges2_in)
  new_edges2 = getindices(map_edges2, keys(edge_data(graph2)))
  new_edge_data2 = Dictionary(new_edges2, edge_data(graph2))

  # convert is needed because type information is lost by Dictionaries.jl
  new_edge_data = convert(typeof(edge_data(graph1)), merge(new_edge_data1, new_edge_data2))

  return NamedDimDataGraph(new_underlying_graph, new_vertex_data, new_edge_data)
end


text:
"""
Subroutine dmid

This subroutine computes the projection of x
on the n-dimensional interval [xl,xu].

MINPACK-2 Project. March 1999.
Argonne National Laboratory.
Chih-Jen Lin and Jorge J. More'.
"""
function dmid(n,x,xl,xu)
    @inbounds for i=1:n
        x[i] = max(xl[i],min(x[i],xu[i]))
    end
    return
end


text:
using DSGE
using HDF5, Base.Test
using DataFrames

path = dirname(@__FILE__)

# Set up model for testing
custom_settings = Dict{Symbol, Setting}(
    :date_forecast_start  => Setting(:date_forecast_start, quartertodate("2015-Q4")))
m = AnSchorfheide(custom_settings = custom_settings, testing = true)

# Read in the data, mode, and hessian
mode = h5read("$path/../reference/hessian.h5","paramsmode")
data = h5read("$path/../reference/hessian.h5","data")'

# Read in the covariance matrix for Metropolis-Hastings and reference parameter draws
hessian_inv, ref_draws, ref_cov =
    h5open("$path/../reference/metropolis_hastings.h5", "r") do file
        read(file, "hessian_inv"),
        read(file, "mhparams"),
        read(file, "ref_cov")
    end

# Set up and run metropolis-hastings
update!(m, mode)
prop_cov = DegenerateMvNormal(mode, hessian_inv)

metropolis_hastings(prop_cov, m, data, .01, .09, verbose=:none)
compute_parameter_covariance(m)

# Read in the parameter draws and covariance just generated from estimate.
test_draws = h5open(rawpath(m, "estimate", "mhsave.h5"), "r") do file
    read(file, "mhparams")
end
test_cov = h5open(workpath(m, "estimate", "parameter_covariance.h5"), "r") do file
    read(file, "mhcov")
end

# Test that the parameter draws and covariance matrices are equal
@testset "Check equality of parameter draws and cov matrices in MH" begin
    @test @test_matrix_approx_eq ref_draws test_draws
    @test @test_matrix_approx_eq ref_cov test_cov
end

nothing


text:
using FinEtools
using FinEtoolsFlexBeams.FESetCorotBeamModule: FESetL2CorotBeam
using FinEtoolsFlexBeams.CrossSectionModule: CrossSectionCircle, CrossSectionRectangle
using FinEtoolsFlexBeams.MeshFrameMemberModule: frame_member, merge_members
using FinEtoolsFlexBeams.VisUtilModule: plot_nodes, plot_midline, render, plot_space_box, plot_solid, space_aspectratio, default_layout_3d, plot_from_json
using PlotlyJS
using JSON

L = 42
xyz = [0 0 0;
0 L/4 L*1/4;
L/4 L/4 L*2/4;
L/4 0 L*3/4;
0 0 L]
nL = 20

cs = CrossSectionCircle(s -> 5.9910, s -> [0.0, 0.0, 1.0])
fens, fes = frame_member(xyz, nL, cs)
plots = cat(plot_nodes(fens),
    plot_midline(fens, fes; color = "rgb(155, 155, 255)", lwidth = 4),
    dims = 1)
layout = default_layout_3d(width = 500, height = 500)
layout[:title] = ""
config  = PlotConfig(plotlyServerURL="https://chart-studio.plotly.com", showLink=true)

pl = plot(plots, layout; config=config)

display(pl)

savejson(pl, "plot.json")

# The below works when run interactively from the command line:
pl = plot_from_json("plot.json")


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

# LQ Factorizations

struct LQ{T,S<:AbstractMatrix} <: Factorization{T}
    factors::S
    τ::Vector{T}
    LQ{T,S}(factors::AbstractMatrix{T}, τ::Vector{T}) where {T,S<:AbstractMatrix} = new(factors, τ)
end

struct LQPackedQ{T,S<:AbstractMatrix} <: AbstractMatrix{T}
    factors::Matrix{T}
    τ::Vector{T}
    LQPackedQ{T,S}(factors::AbstractMatrix{T}, τ::Vector{T}) where {T,S<:AbstractMatrix} = new(factors, τ)
end

LQ(factors::AbstractMatrix{T}, τ::Vector{T}) where {T} = LQ{T,typeof(factors)}(factors, τ)
LQPackedQ(factors::AbstractMatrix{T}, τ::Vector{T}) where {T} = LQPackedQ{T,typeof(factors)}(factors, τ)

"""
    lqfact!(A) -> LQ

Compute the LQ factorization of `A`, using the input
matrix as a workspace. See also [`lq`](@ref).
"""
lqfact!(A::StridedMatrix{<:BlasFloat}) = LQ(LAPACK.gelqf!(A)...)
"""
    lqfact(A) -> LQ

Compute the LQ factorization of `A`. See also [`lq`](@ref).
"""
lqfact(A::StridedMatrix{<:BlasFloat})  = lqfact!(copy(A))
lqfact(x::Number) = lqfact(fill(x,1,1))

"""
    lq(A; [thin=true]) -> L, Q

Perform an LQ factorization of `A` such that `A = L*Q`. The
default is to compute a thin factorization. The LQ factorization
is the QR factorization of `A.'`. `L` is not extended with
zeros if the full `Q` is requested.
"""
function lq(A::Union{Number, AbstractMatrix}; thin::Bool=true)
    F = lqfact(A)
    F[:L], full(F[:Q], thin=thin)
end

copy(A::LQ) = LQ(copy(A.factors), copy(A.τ))

convert{T}(::Type{LQ{T}},A::LQ) = LQ(convert(AbstractMatrix{T}, A.factors), convert(Vector{T}, A.τ))
convert{T}(::Type{Factorization{T}}, A::LQ{T}) = A
convert{T}(::Type{Factorization{T}}, A::LQ) = convert(LQ{T}, A)
convert(::Type{AbstractMatrix}, A::LQ) = A[:L]*A[:Q]
convert(::Type{AbstractArray}, A::LQ) = convert(AbstractMatrix, A)
convert(::Type{Matrix}, A::LQ) = convert(Array, convert(AbstractArray, A))
convert(::Type{Array}, A::LQ) = convert(Matrix, A)
full(A::LQ) = convert(AbstractArray, A)

ctranspose{T}(A::LQ{T}) = QR{T,typeof(A.factors)}(A.factors', A.τ)

function getindex(A::LQ, d::Symbol)
    m, n = size(A)
    if d == :L
        return tril!(A.factors[1:m, 1:min(m,n)])
    elseif d == :Q
        return LQPackedQ(A.factors,A.τ)
    else
        throw(KeyError(d))
    end
end

function getindex(A::LQPackedQ, i::Integer, j::Integer)
    x = zeros(eltype(A), size(A, 1))
    x[i] = 1
    y = zeros(eltype(A), size(A, 2))
    y[j] = 1
    return dot(x, A*y)
end

getq(A::LQ) = LQPackedQ(A.factors, A.τ)

convert{T}(::Type{LQPackedQ{T}}, Q::LQPackedQ) = LQPackedQ(convert(AbstractMatrix{T}, Q.factors), convert(Vector{T}, Q.τ))
convert{T}(::Type{AbstractMatrix{T}}, Q::LQPackedQ) = convert(LQPackedQ{T}, Q)
convert(::Type{Matrix}, A::LQPackedQ) = LAPACK.orglq!(copy(A.factors),A.τ)
convert(::Type{Array}, A::LQPackedQ) = convert(Matrix, A)
function full{T}(A::LQPackedQ{T}; thin::Bool = true)
    #= We construct the full eye here, even though it seems inefficient, because
    every element in the output matrix is a function of all the elements of
    the input matrix. The eye is modified by the elementary reflectors held
    in A, so this is not just an indexing operation. Note that in general
    explicitly constructing Q, rather than using the ldiv or mult methods,
    may be a wasteful allocation. =#
    if thin
        convert(Array, A)
    else
        A_mul_B!(A, eye(T, size(A.factors,2), size(A.factors,1)))
    end
end

size(A::LQ, dim::Integer) = size(A.factors, dim)
size(A::LQ) = size(A.factors)
function size(A::LQPackedQ, dim::Integer)
    if 0 < dim && dim <= 2
        return size(A.factors, dim)
    elseif 0 < dim && dim > 2
        return 1
    else
        throw(BoundsError())
    end
end

size(A::LQPackedQ) = size(A.factors)

## Multiplication by LQ
A_mul_B!{T<:BlasFloat}(A::LQ{T}, B::StridedVecOrMat{T}) = A[:L]*LAPACK.ormlq!('L','N',A.factors,A.τ,B)
A_mul_B!{T<:BlasFloat}(A::LQ{T}, B::QR{T}) = A[:L]*LAPACK.ormlq!('L','N',A.factors,A.τ,full(B))
A_mul_B!{T<:BlasFloat}(A::QR{T}, B::LQ{T}) = A_mul_B!(zeros(full(A)), full(A), full(B))
function *{TA,TB}(A::LQ{TA},B::StridedVecOrMat{TB})
    TAB = promote_type(TA, TB)
    A_mul_B!(convert(Factorization{TAB},A), copy_oftype(B, TAB))
end
function *{TA,TB}(A::LQ{TA},B::QR{TB})
    TAB = promote_type(TA, TB)
    A_mul_B!(convert(Factorization{TAB},A), convert(Factorization{TAB},B))
end
function *{TA,TB}(A::QR{TA},B::LQ{TB})
    TAB = promote_type(TA, TB)
    A_mul_B!(convert(Factorization{TAB},A), convert(Factorization{TAB},B))
end

## Multiplication by Q
### QB
A_mul_B!{T<:BlasFloat}(A::LQPackedQ{T}, B::StridedVecOrMat{T})   = LAPACK.ormlq!('L','N',A.factors,A.τ,B)
function (*)(A::LQPackedQ,B::StridedVecOrMat)
    TAB = promote_type(eltype(A), eltype(B))
    A_mul_B!(convert(AbstractMatrix{TAB}, A), copy_oftype(B, TAB))
end

### QcB
Ac_mul_B!{T<:BlasReal}(A::LQPackedQ{T}, B::StridedVecOrMat{T})    = LAPACK.ormlq!('L','T',A.factors,A.τ,B)
Ac_mul_B!{T<:BlasComplex}(A::LQPackedQ{T}, B::StridedVecOrMat{T}) = LAPACK.ormlq!('L','C',A.factors,A.τ,B)
function Ac_mul_B(A::LQPackedQ, B::StridedVecOrMat)
    TAB = promote_type(eltype(A), eltype(B))
    if size(B,1) == size(A.factors,2)
        Ac_mul_B!(convert(AbstractMatrix{TAB}, A), copy_oftype(B, TAB))
    elseif size(B,1) == size(A.factors,1)
        Ac_mul_B!(convert(AbstractMatrix{TAB}, A), [B; zeros(TAB, size(A.factors, 2) - size(A.factors, 1), size(B, 2))])
    else
        throw(DimensionMismatch("first dimension of B, $(size(B,1)), must equal one of the dimensions of A, $(size(A))"))
    end
end

### QBc/QcBc
for (f1, f2) in ((:A_mul_Bc, :A_mul_B!),
                 (:Ac_mul_Bc, :Ac_mul_B!))
    @eval begin
        function ($f1)(A::LQPackedQ, B::StridedVecOrMat)
            TAB = promote_type(eltype(A), eltype(B))
            BB = similar(B, TAB, (size(B, 2), size(B, 1)))
            ctranspose!(BB, B)
            return ($f2)(A, BB)
        end
    end
end

### AQ
A_mul_B!{T<:BlasFloat}(A::StridedMatrix{T}, B::LQPackedQ{T}) = LAPACK.ormlq!('R', 'N', B.factors, B.τ, A)
function *{TA,TB}(A::StridedMatrix{TA},B::LQPackedQ{TB})
    TAB = promote_type(TA,TB)
    if size(B.factors,2) == size(A,2)
        A_mul_B!(copy_oftype(A, TAB),convert(AbstractMatrix{TAB},B))
    elseif size(B.factors,1) == size(A,2)
        A_mul_B!( [A zeros(TAB, size(A,1), size(B.factors,2)-size(B.factors,1))], convert(AbstractMatrix{TAB},B))
    else
        throw(DimensionMismatch("second dimension of A, $(size(A,2)), must equal one of the dimensions of B, $(size(B))"))
    end
end

### AQc
A_mul_Bc!{T<:BlasReal}(A::StridedMatrix{T}, B::LQPackedQ{T})    = LAPACK.ormlq!('R','T',B.factors,B.τ,A)
A_mul_Bc!{T<:BlasComplex}(A::StridedMatrix{T}, B::LQPackedQ{T}) = LAPACK.ormlq!('R','C',B.factors,B.τ,A)
function A_mul_Bc{TA<:Number,TB<:Number}( A::StridedVecOrMat{TA}, B::LQPackedQ{TB})
    TAB = promote_type(TA,TB)
    A_mul_Bc!(copy_oftype(A, TAB), convert(AbstractMatrix{TAB},(B)))
end

### AcQ/AcQc
for (f1, f2) in ((:Ac_mul_B, :A_mul_B!),
                 (:Ac_mul_Bc, :A_mul_Bc!))
    @eval begin
        function ($f1)(A::StridedMatrix, B::LQPackedQ)
            TAB = promote_type(eltype(A), eltype(B))
            AA = similar(A, TAB, (size(A, 2), size(A, 1)))
            ctranspose!(AA, A)
            return ($f2)(AA, B)
        end
    end
end

function \{TA,Tb}(A::LQ{TA}, b::StridedVector{Tb})
    S = promote_type(TA,Tb)
    m = checksquare(A)
    m == length(b) || throw(DimensionMismatch("left hand side has $m rows, but right hand side has length $(length(b))"))
    AA = convert(Factorization{S}, A)
    x = A_ldiv_B!(AA, copy_oftype(b, S))
    return x
end
function \{TA,TB}(A::LQ{TA},B::StridedMatrix{TB})
    S = promote_type(TA,TB)
    m = checksquare(A)
    m == size(B,1) || throw(DimensionMismatch("left hand side has $m rows, but right hand side has $(size(B,1)) rows"))
    AA = convert(Factorization{S}, A)
    X = A_ldiv_B!(AA, copy_oftype(B, S))
    return X
end
# With a real lhs and complex rhs with the same precision, we can reinterpret
# the complex rhs as a real rhs with twice the number of columns
function (\){T<:BlasReal}(F::LQ{T}, B::VecOrMat{Complex{T}})
    c2r = reshape(transpose(reinterpret(T, B, (2, length(B)))), size(B, 1), 2*size(B, 2))
    x = A_ldiv_B!(F, c2r)
    return reinterpret(Complex{T}, transpose(reshape(x, div(length(x), 2), 2)),
        isa(B, AbstractVector) ? (size(F,2),) : (size(F,2), size(B,2)))
end

function A_ldiv_B!{T}(A::LQ{T}, B::StridedVecOrMat{T})
    Ac_mul_B!(A[:Q], A_ldiv_B!(LowerTriangular(A[:L]),B))
    return B
end

max_stars_repo_name: felixrehren/julia

max_issues_repo_name: felixrehren/julia

max_forks_repo_name: felixrehren/julia


text:
"""
function getDCResistivitySourcesAndRecAll

sets up dipole sources for DC resistivity survey. Sources and receivers
are spaced equidistantly on nodes of top x1-x2 surface of 3D cube with Mesh.n
cells.
Sources are dipoles of width 2*offset and are placed in x1 and x2 direction.
Receivers (with much shorter width) are placed on entire domain and measure
potential differences in x1 and x2 direction.

Input:

	Mesh - mesh for forward problem

Optional Inputs:

	srcSpacing - spacing between sources      (default=[2,2])
	srcPad     - padding from boundary        (default=[2,2])
	offset     - half the length of sources   (default=12)
	offsetRec  - half the length of receivers (default-1)

Outputs:

	Sources    - discretized sources
	Receivers  - discretized receivers
"""
function getDCResistivitySourcesAndRecAll(Mesh;srcSpacing::Vector=[2,2],srcPad::Vector=[4,4,4,4],
	                                                               offset=12,offsetRec=1)

	n    = Mesh.n
	i1   =  srcPad[1]+offset:srcSpacing[1]:n[1]-srcPad[2]-offset
	i2   =  srcPad[3]+offset:srcSpacing[2]:n[2]-srcPad[4]-offset
	nsrc = 2*length(i1)*length(i2)

	println("--- generate $nsrc sources ---")
	Sources   = spzeros(prod(Mesh.n+1),nsrc)
	cnt = 1
	# x1-direction source dipoles
	for i=i1
		for j=i2
			# source
			S                  = zeros(tuple(n+1...))
			S[i-offset,j,1]  = -1
			S[i+offset,j,1]  =  1
			Sources[:,cnt]     = sparse(vec(S))
			cnt +=1
		end
	end
	# y2 direction dipoles sources
	for i=i1
		for j=i2
			S               = zeros(tuple(n+1...))
			S[i,j-offset,1] = -1
			S[i,j+offset,1] =  1
			Sources[:,cnt]  = sparse(vec(S))
			cnt +=1
		end
	end

	i1   =  srcPad[1]+offsetRec:1:n[1]-srcPad[2]-offsetRec
	i2   =  srcPad[3]+offsetRec:1:n[2]-srcPad[4]-offsetRec
	nrec = 2*length(i1)*length(i2)

	println("--- generate $nrec receivers ---")
	# x1-direction receiver dipoles
	Receivers = Sources;
	return Sources,Receivers
end