text:
using LinearAlgebra
import Adapt

@testset "constructors" begin
  xs = CuArray{Int}(undef, 2, 3)
  @test collect(CuArray([1 2; 3 4])) == [1 2; 3 4]
  @test collect(cu[1, 2, 3]) == [1, 2, 3]
  @test collect(cu([1, 2, 3])) == [1, 2, 3]
  @test testf(vec, rand(5,3))
  @test cu(1:3) === 1:3
  @test Base.elsize(xs) == sizeof(Int)
  @test CuArray{Int, 2}(xs) === xs

  # test aggressive conversion to Float32, but only for floats, and only with `cu`
  @test cu([1]) isa AbstractArray{Int}
  @test cu(Float64[1]) isa AbstractArray{Float32}
  @test cu(ComplexF64[1+1im]) isa AbstractArray{ComplexF32}
  @test Adapt.adapt(CuArray, Float64[1]) isa AbstractArray{Float64}
  @test Adapt.adapt(CuArray, ComplexF64[1]) isa AbstractArray{ComplexF64}
  @test Adapt.adapt(CuArray{Float16}, Float64[1]) isa AbstractArray{Float16}

  @test_throws ArgumentError Base.unsafe_convert(Ptr{Int}, xs)
  @test_throws ArgumentError Base.unsafe_convert(Ptr{Float32}, xs)

  # unsafe_wrap
  @test Base.unsafe_wrap(CuArray, CU_NULL, 1; own=false).state == CUDA.ARRAY_UNMANAGED
  ## compare structs without using mapreduce
  structeq(a,b) = false
  structeq(a::T,b::T) where T = all(field->getfield(a,field) == getfield(b,field), fieldnames(T))
  @test structeq(Base.unsafe_wrap(CuArray, CU_NULL, 2),                CuArray{Nothing,1}(CU_NULL, (2,)))
  @test structeq(Base.unsafe_wrap(CuArray{Nothing}, CU_NULL, 2),       CuArray{Nothing,1}(CU_NULL, (2,)))
  @test structeq(Base.unsafe_wrap(CuArray{Nothing,1}, CU_NULL, 2),     CuArray{Nothing,1}(CU_NULL, (2,)))
  @test structeq(Base.unsafe_wrap(CuArray, CU_NULL, (1,2)),            CuArray{Nothing,2}(CU_NULL, (1,2)))
  @test structeq(Base.unsafe_wrap(CuArray{Nothing}, CU_NULL, (1,2)),   CuArray{Nothing,2}(CU_NULL, (1,2)))
  @test structeq(Base.unsafe_wrap(CuArray{Nothing,2}, CU_NULL, (1,2)), CuArray{Nothing,2}(CU_NULL, (1,2)))

  @test collect(CUDA.zeros(2, 2)) == zeros(Float32, 2, 2)
  @test collect(CUDA.ones(2, 2)) == ones(Float32, 2, 2)

  @test collect(CUDA.fill(0, 2, 2)) == zeros(Float32, 2, 2)
  @test collect(CUDA.fill(1, 2, 2)) == ones(Float32, 2, 2)
end

@testset "adapt" begin
  A = rand(Float32, 3, 3)
  dA = CuArray(A)
  @test Adapt.adapt(Array, dA) == A
  @test Adapt.adapt(CuArray, A) isa CuArray
  @test Array(Adapt.adapt(CuArray, A)) == A
end

@testset "view" begin
  @test testf(rand(5)) do x
    y = x[2:4]
    y .= 1
    x
  end

  @test testf(rand(5)) do x
    y = view(x, 2:4)
    y .= 1
    x
  end

  @test testf(x->view(x, :, 1:4, 3), rand(Float32, 5, 4, 3))

  let x = CUDA.rand(Float32, 5, 4, 3)
    @test_throws BoundsError view(x, :, :, 1:10)
  end

  # bug in parentindices conversion
  let x = CuArray{Int}(undef, 1, 1)
    x[1,:] .= 42
    @test Array(x)[1,1] == 42
  end

  # bug in conversion of indices (#506)
  show(devnull, cu(view(ones(1), [1])))

  # performance loss due to Array indices
  let x = CuArray{Int}(undef, 1)
    i = [1]
    y = view(x, i)
    @test parent(y) isa CuArray
    @test parentindices(y) isa Tuple{CuArray}
  end

  @testset "GPU array source" begin
      a = rand(3)
      i = rand(1:3, 2)
      @test testf(view, a, i)
      @test testf(view, a, view(i, 2:2))
  end

  @testset "CPU array source" begin
      a = rand(3)
      i = rand(1:3, 2)
      @test testf(view, a, i)
      @test testf(view, a, view(i, 2:2))
  end

  @testset "unmanaged view" begin
    a = CuArray([1,2,3])
    ptr = pointer(a, 2)

    b = unsafe_wrap(CuArray, ptr, 2)
    @test Array(b) == [2,3]

    c = view(b, 2:2)
    @test Array(c) == [3]
  end
end

@testset "reshape" begin
  A = [1 2 3 4
       5 6 7 8]
  gA = reshape(CuArray(A),1,8)
  _A = reshape(A,1,8)
  _gA = Array(gA)
  @test all(_A .== _gA)
  A = [1,2,3,4]
  gA = reshape(CuArray(A),4)

  @testset "unmanaged reshape" begin
    a = CuArray([1,2,3])
    ptr = pointer(a, 2)

    b = unsafe_wrap(CuArray, ptr, 2)
    @test Array(b) == [2,3]

    c = reshape(b, (1,2))
    @test Array(c) == [2 3]
  end
end

@testset "reinterpret" begin
  A = Int32[-1,-2,-3]
  dA = CuArray(A)
  dB = reinterpret(UInt32, dA)
  @test reinterpret(UInt32, A) == Array(dB)

  @testset "unmanaged reinterpret" begin
    a = CuArray(Int32[-1,-2,-3])
    ptr = pointer(a, 2)

    b = unsafe_wrap(CuArray, ptr, 2)
    @test Array(b) == Int32[-2,-3]

    c = reinterpret(UInt32, b)
    @test Array(c) == reinterpret(UInt32, Int32[-2,-3])
  end
end

@testset "Dense derivatives" begin
  a = CUDA.rand(Int64, 5, 4, 3)
  @test a isa CuArray

  # Contiguous views should return new CuArray
  @test view(a, :, 1, 2) isa CuVector{Int64}
  @test view(a, 1:4, 1, 2) isa CuVector{Int64}
  @test view(a, :, 1:4, 3) isa CuMatrix{Int64}
  @test view(a, :, :, 1) isa CuMatrix{Int64}
  @test view(a, :, :, :) isa CuArray{Int64,3}
  @test view(a, :) isa CuVector{Int64}
  @test view(a, 1:3) isa CuVector{Int64}
  @test view(a, 1, 1, 1) isa CuArray{Int64}

  # Non-contiguous views should fall back to base's SubArray
  @test view(a, 1:3, 1:3, 3) isa SubArray
  @test view(a, 1, :, 3) isa SubArray
  @test view(a, 1, 1:4, 3) isa SubArray
  @test view(a, :, 1, 1:3) isa SubArray
  @test view(a, :, 1:2:4, 1) isa SubArray
  @test view(a, 1:2:5, 1, 1) isa SubArray

  b = reshape(a, (6,10))
  @test b isa CuArray
  @test b isa StridedCuArray
  @test view(b, :, :, 1) isa DenseCuArray

  b = reshape(a, :)
  @test b isa CuArray

  b = reinterpret(Float64, a)
  @test b isa CuArray
  @test b isa StridedCuArray
  @test view(b, :, :, 1) isa DenseCuArray
end

@testset "StridedArray" begin
  a = CUDA.rand(Int64, 2,2,2)
  @test a isa StridedCuArray

  @test view(a, :, :, 1) isa StridedCuArray
  @test view(a, :, 1, :) isa StridedCuArray
  @test view(a, 1, :, :) isa StridedCuArray

  b = reshape(a, (2,4))
  @test b isa CuArray
  @test b isa StridedCuArray
  @test view(b, :, 1, :) isa StridedCuArray

  b = reinterpret(Float64, a)
  @test b isa CuArray
  @test b isa StridedCuArray
  @test view(b, :, 1, :) isa StridedCuArray
end

@testset "accumulate" begin
  for n in (0, 1, 2, 3, 10, 10_000, 16384, 16384+1) # small, large, odd & even, pow2 and not
    @test testf(x->accumulate(+, x), rand(n))
  end

  # multidimensional
  for (sizes, dims) in ((2,) => 2,
                        (3,4,5) => 2,
                        (1, 70, 50, 20) => 3)
    @test testf(x->accumulate(+, x; dims=dims), rand(Int, sizes))
  end

  # using initializer
  for (sizes, dims) in ((2,) => 2,
                        (3,4,5) => 2,
                        (1, 70, 50, 20) => 3)
    @test testf(x->accumulate(+, x; dims=dims, init=100.), rand(Int, sizes))
  end

  # in place
  @test testf(x->(accumulate!(+, x, copy(x)); x), rand(2))

  # specialized
  @test testf(cumsum, rand(2))
  @test testf(cumprod, rand(2))
end

@testset "logical indexing" begin
  @test CuArray{Int}(undef, 2)[CuArray{Bool}(undef, 2)] isa CuArray
  @test CuArray{Int}(undef, 2, 2)[CuArray{Bool}(undef, 2, 2)] isa CuArray
  @test CuArray{Int}(undef, 2, 2, 2)[CuArray{Bool}(undef, 2, 2, 2)] isa CuArray

  @test CuArray{Int}(undef, 2)[Array{Bool}(undef, 2)] isa CuArray
  @test CuArray{Int}(undef, 2, 2)[Array{Bool}(undef, 2, 2)] isa CuArray
  @test CuArray{Int}(undef, 2, 2, 2)[Array{Bool}(undef, 2, 2, 2)] isa CuArray

  @test testf((x,y)->x[y], rand(2), rand(Bool, 2))
  @test testf((x,y)->x[y], rand(2, 2), rand(Bool, 2, 2))
  @test testf((x,y)->x[y], rand(2, 2, 2), rand(Bool, 2, 2, 2))

  @test testf(x -> x[x .> 0.5], rand(2))
  @test testf(x -> x[x .> 0.5], rand(2,2))
  @test testf(x -> x[x .> 0.5], rand(2,2,2))

  @test testf(x -> filter(y->y .> 0.5, x), rand(2))
  @test testf(x -> filter(y->y .> 0.5, x), rand(2,2))
  @test testf(x -> filter(y->y .> 0.5, x), rand(2,2,2))
end

@testset "reverse" begin
    # 1-d out-of-place
    @test testf(x->reverse(x), rand(1000))
    @test testf(x->reverse(x, 10), rand(1000))
    @test testf(x->reverse(x, 10, 90), rand(1000))

    # 1-d in-place
    @test testf(x->reverse!(x), rand(1000))
    @test testf(x->reverse!(x, 10), rand(1000))
    @test testf(x->reverse!(x, 10, 90), rand(1000))

    # n-d out-of-place
    for shape in ([1, 2, 4, 3], [4, 2], [5], [2^5, 2^5, 2^5]),
        dim in 1:length(shape)
      @test testf(x->reverse(x; dims=dim), rand(shape...))

      cpu = rand(shape...)
      gpu = CuArray(cpu)
      reverse!(gpu; dims=dim)
      @test Array(gpu) == reverse(cpu; dims=dim)
    end

    # wrapped array
    @test testf(x->reverse(x), reshape(rand(2,2), 4))
end

@testset "findall" begin
    # 1D
    @test testf(x->findall(x), rand(Bool, 100))
    @test testf(x->findall(y->y>0.5, x), rand(100))

    # ND
    let x = rand(Bool, 10, 10)
      @test findall(x) == Array(findall(CuArray(x)))
    end
    let x = rand(10, 10)
      @test findall(y->y>0.5, x) == Array(findall(y->y>0.5, CuArray(x)))
    end
end

@testset "findfirst" begin
    # 1D
    @test testf(x->findfirst(x), rand(Bool, 100))
    @test testf(x->findfirst(y->y>0.5, x), rand(100))
    let x = fill(false, 10)
      @test findfirst(x) == findfirst(CuArray(x))
    end

    # ND
    let x = rand(Bool, 10, 10)
      @test findfirst(x) == findfirst(CuArray(x))
    end
    let x = rand(10, 10)
      @test findfirst(y->y>0.5, x) == findfirst(y->y>0.5, CuArray(x))
    end
end

@testset "findmax & findmin" begin
  let x = rand(Float32, 100)
      @test findmax(x) == findmax(CuArray(x))
      @test findmax(x; dims=1) == Array.(findmax(CuArray(x); dims=1))

      x[32] = x[33] = x[55] = x[66] = NaN32
      @test isequal(findmax(x), findmax(CuArray(x)))
      @test isequal(findmax(x; dims=1), Array.(findmax(CuArray(x); dims=1)))
  end
  let x = rand(Float32, 10, 10)
      @test findmax(x) == findmax(CuArray(x))
      @test findmax(x; dims=1) == Array.(findmax(CuArray(x); dims=1))
      @test findmax(x; dims=2) == Array.(findmax(CuArray(x); dims=2))

      x[rand(CartesianIndices((10, 10)), 10)] .= NaN
      @test isequal(findmax(x), findmax(CuArray(x)))
      @test isequal(findmax(x; dims=1), Array.(findmax(CuArray(x); dims=1)))
  end
  let x = rand(Float32, 10, 10, 10)
      @test findmax(x) == findmax(CuArray(x))
      @test findmax(x; dims=1) == Array.(findmax(CuArray(x); dims=1))
      @test findmax(x; dims=2) == Array.(findmax(CuArray(x); dims=2))
      @test findmax(x; dims=3) == Array.(findmax(CuArray(x); dims=3))

      x[rand(CartesianIndices((10, 10, 10)), 20)] .= NaN
      @test isequal(findmax(x), findmax(CuArray(x)))
      @test isequal(findmax(x; dims=1), Array.(findmax(CuArray(x); dims=1)))
      @test isequal(findmax(x; dims=2), Array.(findmax(CuArray(x); dims=2)))
      @test isequal(findmax(x; dims=3), Array.(findmax(CuArray(x); dims=3)))
  end

  let x = rand(Float32, 100)
      @test findmin(x) == findmin(CuArray(x))
      @test findmin(x; dims=1) == Array.(findmin(CuArray(x); dims=1))

      x[32] = x[33] = x[55] = x[66] = NaN32
      @test isequal(findmin(x), findmin(CuArray(x)))
      @test isequal(findmin(x; dims=1), Array.(findmin(CuArray(x); dims=1)))
  end
  let x = rand(Float32, 10, 10)
      @test findmin(x) == findmin(CuArray(x))
      @test findmin(x; dims=1) == Array.(findmin(CuArray(x); dims=1))
      @test findmin(x; dims=2) == Array.(findmin(CuArray(x); dims=2))

      x[rand(CartesianIndices((10, 10)), 10)] .= NaN
      @test isequal(findmin(x), findmin(CuArray(x)))
      @test isequal(findmin(x; dims=1), Array.(findmin(CuArray(x); dims=1)))
      @test isequal(findmin(x; dims=2), Array.(findmin(CuArray(x); dims=2)))
      @test isequal(findmin(x; dims=3), Array.(findmin(CuArray(x); dims=3)))
  end
  let x = rand(Float32, 10, 10, 10)
      @test findmin(x) == findmin(CuArray(x))
      @test findmin(x; dims=1) == Array.(findmin(CuArray(x); dims=1))
      @test findmin(x; dims=2) == Array.(findmin(CuArray(x); dims=2))
      @test findmin(x; dims=3) == Array.(findmin(CuArray(x); dims=3))

      x[rand(CartesianIndices((10, 10, 10)), 20)] .= NaN
      @test isequal(findmin(x), findmin(CuArray(x)))
      @test isequal(findmin(x; dims=1), Array.(findmin(CuArray(x); dims=1)))
      @test isequal(findmin(x; dims=2), Array.(findmin(CuArray(x); dims=2)))
      @test isequal(findmin(x; dims=3), Array.(findmin(CuArray(x); dims=3)))
  end
end

@testset "argmax & argmin" begin
    @test testf(argmax, rand(Int, 10))
    @test testf(argmax, -rand(Int, 10))

    @test testf(argmin, rand(Int, 10))
    @test testf(argmin, -rand(Int, 10))
end

@testset "issue #543" begin
  x = CUDA.rand(ComplexF32, 1)
  @test x isa CuArray{Complex{Float32}}

  y = exp.(x)
  @test y isa CuArray{Complex{Float32}}
end

@testset "resizing" begin
    a = CuArray([1,2,3])

    resize!(a, 3)
    @test length(a) == 3
    @test Array(a) == [1,2,3]

    resize!(a, 5)
    @test length(a) == 5
    @test Array(a)[1:3] == [1,2,3]

    resize!(a, 2)
    @test length(a) == 2
    @test Array(a)[1:2] == [1,2]

    GC.@preserve a begin
      b = unsafe_wrap(CuArray{Int}, pointer(a), 2)
      @test_throws ArgumentError resize!(b, 3)
    end
end

@testset "aliasing" begin
  x = CuArray([1,2])
  y = view(x, 2:2)
  @test Base.mightalias(x, x)
  @test Base.mightalias(x, y)
  z = view(x, 1:1)
  @test Base.mightalias(x, z)
  @test !Base.mightalias(y, z)

  a = copy(y)::typeof(x)
  @test !Base.mightalias(x, a)
  a .= 3
  @test Array(y) == [2]

  b = Base.unaliascopy(y)::typeof(y)
  @test !Base.mightalias(x, b)
  b .= 3
  @test Array(y) == [2]
end


text:
# How to build an executable

using PackageCompiler

PackageCompiler.create_app(
    ".",
    "plutotool";
    # app_name="plutotool",
    include_transitive_dependencies=true,
    force=true)

# The did somethinbg but I have no clue what to do with the result.

#=

I ran it like this:

plutotool\bin\plutotool.exe help

and it printed out the help.

I expect it needs all of the files that it put there.

=#

ext: jl


text:
using Documenter

push!(LOAD_PATH, "../src/")
using TwoDimensional
using TwoDimensional: WeightedPoint, WeightedPoint2D
using TwoDimensional.Suffixed

DEPLOYDOCS = (get(ENV, "CI", nothing) == "true")

makedocs(
    sitename = "TwoDimensional.jl Package",
    format = Documenter.HTML(
        prettyurls = DEPLOYDOCS,
    ),
    authors = "Éric Thiébaut and contributors",
    pages = ["index.md", "install.md", "AffineTransforms.md",
             "Points.md", "BoundingBoxes.md", "reference.md"]
)

if DEPLOYDOCS
    deploydocs(
        repo = "github.com/emmt/TwoDimensional.jl.git",
    )
end


text:
using Survey
using Test

@testset "Survey.jl" begin
    # Write your tests here.
end

max_stars_repo_name: samplics-org/survey-jl


text:
module Associates

# import InverseFunctions: inverse

import Base: show, inv, ∘, \, /, ∪, ∩
import Base.Iterators: product

#export types
export GMap

#export functions
export show, dom, codom, image, preimage, inverse, inv, ∘, compose, \, /, ∪, ∩, ×

struct GMap
    forward::Dict
    backward::Dict
end

function show(io::IO, m::GMap)
    dt = keytype(m.forward)
    ct = keytype(m.backward)

    println(io, "$dt -> $ct")

    compact = get(io, :compact, false)
    # if !compact
    #     show(io, m.forward)
    # end
end

GMap(forward::Dict) = GMap(forward, invert(forward))

GMap(fun::Function, d) = GMap(fun2dict(fun, d))

(m::GMap)(x) = applymap(m, x)
applymap(m::GMap, x) = m.forward[x]
dom(m::GMap) = keys(m.forward)
inverse(m::GMap) = GMap(m.backward, m.forward)
inv(m::GMap) = inverse(m)

(⋓)(a, b) = (∪)([a], [b])
(⋓)(a::T, b) where {T<:Union{AbstractVector,AbstractSet}} = (∪)(a, [b])
(⋓)(a, b::T) where {T<:Union{AbstractVector,AbstractSet}} = (∪)([a], b)
(⋓)(a::T, b::N) where {T<:Union{AbstractVector,AbstractSet},N<:Union{AbstractVector,AbstractSet}} = (∪)(a, b)
image(m::GMap, x::T) where {T<:Union{AbstractVector,AbstractSet}} = reduce(⋓, m.(x))
image(m::GMap, x) = reduce(⋓, m(x))
image(m::GMap) = codom(m)
preimage(m::GMap) = dom(m)
preimage(m::GMap, x) = image(inverse(m), x)
codom(m::GMap) = dom(inverse(m))

pull_back(m2, m1) = dom(m2) ∩ codom(m1) |> x -> preimage(m1, x)
push_forward(m2, m1) = m2.(m1.(pull_back(m2, m1)))
compose(m2, m1) = GMap(Dict(zip(pull_back(m2, m1), push_forward(m2, m1))))
# compose(m2, m1) = GMap(x->m2(x), collect(pull_back(m2,m1)))

(∘)(m2::GMap, m1::GMap) = compose(m2, m1)

(∘)(s::Set, m1::GMap) = Set(s ∩ codom(m1))

(∘)(m2::GMap, S::Set) = Set(image(m2, S ∩ dom(m2)))

(∘)(f::Function, m1::GMap) = compose(GMap(f, collect(codom(m1))), m1)

(∪)(m1::GMap, m2::GMap) = GMap(Dict(m1.forward ∪ m2.forward))

(∩)(m1::GMap, m2::GMap) = GMap(Dict(m1.forward ∩ m2.forward))

(/)(m1::GMap, m2::GMap) = m2 ∘ inverse(m1)

(\)(m1::GMap, m2::GMap) = inverse(m2) ∘ m1

(×)(m1::GMap, m2::GMap) = GMap(Dict(zip(product(keys(m2.forward), keys(m1.forward)), product(values(m2.forward), values(m1.forward)))))

#some utilities
hashkv(D) = Dict(hash(value) => hash(key) for (key, value) in D)

fun2dict(fun, d) = Dict(zip(d, map(fun, d)))

function unwrapsinglet(x)
    if isa(x, Vector)
        length(x) == 1 ? x[1] : x
    else
        x
    end
end

function invert(D::Dict)
    invD = Dict()
    for k in keys(D)
        if D[k] in keys(invD)
            push!(invD[D[k]], k)
        else
            invD[D[k]] = [k]
        end
    end
    for k in keys(invD)
        invD[k] = length(invD[k]) == 1 ? invD[k][1] : invD[k]
    end
    K = typejoin(typeof.(keys(invD))...)
    V = typejoin(typeof.(values(invD))...)
    Dict{K,V}(invD)
end

end


text:
using JSON
import JSON: json

export json

# JSON Export for reports, etc.
joreduce(x) = x
joreduce{S,T}(a::Array{S,T}) = map(joreduce, a)
joreduce(m::Monetary) = stringmime("text/plain", m)
joreduce(b::Basket) = joreduce(collect(b))

joreduce{S<:AbstractString,B<:Basket}(b::Pair{S,B}) = Dict(
    "name" => b.first,
    "value" => joreduce(b.second))

joreduce(s::BSSection) = Dict(
    "heading" => joreduce(s.heading),
    "entries" => joreduce(s.entries),
    "subsections" => joreduce(s.subsections),
    "total" => joreduce(total(s)))

joreduce(s::BalanceSheet) = Dict(
    "assets" => joreduce(s.assets),
    "liabilitiesandequity" => joreduce(s.le))

json(s::BalanceSheet) = json(joreduce(s))


text:
function dasherize(cur_string::AbstractString)
  cur_string = _decamelize(cur_string)

  kebabcase(cur_string)
end


text:
#
println("Ci) Test of the DecayYield module with ASF from an internally generated multiplet.")
#
setDefaults("print summary: open", "zzz-DecayYield.sum")
grid = Radial.Grid(Radial.Grid(false), rnt = 4.0e-6, h = 5.0e-2, hp = 1.0e-2, rbox = 10.0)
wa   = Atomic.Computation(Atomic.Computation(), name="xx", grid=grid, nuclearModel=Nuclear.Model(12.),
                          configs=[Configuration("1s 2s^2 2p^6")],
                          propertySettings=[ DecayYield.Settings("SCA", true, LevelSelection()) ] )

wb = perform(wa)
setDefaults("print summary: close", "")


text:
function gen_newton_al_functions(solver::Solver)
    n,m,N = get_sizes(solver)
    p,pI,pE = get_num_constraints(solver)
    p_N,pI_N,pE_N = get_num_terminal_constraints(solver)

    c_function!, c_jacobian!, c_labels, cI!, cE! = generate_constraint_functions(solver.obj)
    f = solver.fd
    s_x = reshape(1:(N*n),n,N)
    s_u = reshape(N*n .+ (1:(m*(N-1))),m,N-1)
    s_λ = (N*n + (N-1)*m) .+ (1:((N-1)*pE + pE_N))
    s_μ = (N*n + (N-1)*m + length(s_λ)) .+ (1:((N-1)*pI + pI_N))
    s_ν = (N*n + (N-1)*m + length(s_λ) + length(s_μ)) .+ (1:N*n)
    s_s = (N*n + (N-1)*m + length(s_λ) + length(s_μ) + N*n) .+ (1:length(s_μ))
    sze = (X=s_x, U=s_u, λ=s_λ, μ=s_μ, ν=s_ν, s=s_s)

    function newton_cost(V::Vector,ρ=1)
        X,U,λ,μ,ν,s = unpackV(V)
        e = cE(V)
        i = cI(V) + 0.5s.^2
        dd = d(V)
        J = cost(solver,X,U) + λ'e + μ'i + ν'dd + 0.5*ρ*(e'e + i'i + dd'dd)
        return J
    end

    function cI(V)
        X,U,λ,μ,ν,s = unpackV(V)
        C = zeros(eltype(V),pI,N-1)
        for k = 1:N-1
            cI!(view(C,1:pI,k),X[:,k],U[:,k])
        end
        CN = zeros(eltype(V),pI_N)
        cI!(CN,X[:,N])
        CI = [vec(C); CN]
        return CI
    end

    function cE(V)
        X,U,λ,μ,ν,s = unpackV(V)
        C = zeros(eltype(V),pE,N-1)
        for k = 1:N-1
            cE!(view(C,1:pE,k),X[:,k],U[:,k])
        end
        CN = zeros(eltype(V),pE_N)
        cE!(CN,X[:,N])
        CE = [vec(C); CN]
        return CE
    end

    function d(V)
        X,U,λ,μ,ν,s = unpackV(V)
        D = zeros(eltype(V),n,N)
        D[:,1] = X[:,1] - solver.obj.x0
        for k = 2:N
            f(view(D,1:n,k),X[:,k-1],U[:,k-1])
            D[:,k] -= X[:,k]
        end
        return vec(D)
    end

    function active_set(V,eps=1e-2)
        X,U,λ,μ,ν,s = unpackV(V)
        a = trues(length(V))
        ci = cI(V)
        c_inds = ci .>= -eps
        a[sze.μ] = c_inds
        a[sze.s] = c_inds
        return a
    end

    function packV(X,U,λ,μ,ν,s)
        V = [vec(X);vec(U);λ;μ;ν;s]
    end

    function unpackV(V)
        X = V[sze.X]
        U = V[sze.U]
        λ = V[sze.λ]
        μ = V[sze.μ]
        ν = V[sze.ν]
        s = V[sze.s]
        return X,U,λ,μ,ν,s
    end
    return newton_cost, packV, unpackV, cI, cE, d, active_set, sze
end

function init_newton_results(results::SolverIterResults,Solver::Solver)
    n,m,N = get_sizes(solver)
    p,pI,pE = get_num_constraints(solver)
    pN,pI_N,pE_N = get_num_terminal_constraints(solver)

    X = to_array(res.X)
    U = to_array(res.U)
    μ = vcat([k < N ? res.λ[k][1:pI] : res.λ[k][1:pI_N] for k = 1:N]...)
    λ = vcat([k < N ? res.λ[k][pI .+ (1:pE)] : res.λ[k][pI_N .+ (1:pE_N)] for k = 1:N]...)
    ν = vec(to_array(res.s))
    s = vcat([k < N ? res.C[k][1:pI] : res.C[k][1:pI_N] for k = 1:N]...)
    s = sqrt.(2*max.(0,-s))
    packV(X,U,λ,μ,ν,s)
end

function line_search(V,H,g,a=trues(length(V)))
    J0 = mycost(V)
    α = 1.
    V_ = copy(V)5.17, 6.2 & A5.1, 6.9 & A5.2, A3.18, A3.20, A11.4, and A13.2.
    delta = (H[a,a]\g[a])
    V_[a] = V[a] - α*delta
    J = mycost(V_)
    iter = 1
    while J > J0
        α /= 2
        V_[a] = V[a] - α*delta
        J = mycost(V_)
        iter += 1
        @show J, α
        iter > 10 ? break : nothing
    end
    @show J0-J
    return V_
end

# Block move
model, obj = TrajectoryOptimization.Dynamics.double_integrator
u_min = -0.2
u_max = 0.2
x_bnd = 2
obj_con = TrajectoryOptimization.ConstrainedObjective(obj, tf=5.0,use_xf_equality_constraint=true,u_min=u_min, u_max=u_max, x_max=x_bnd, x_min=-x_bnd)#, x_max=x_max)

solver = TrajectoryOptimization.Solver(model,obj_con,integration=:rk4,N=5,opts=opts)
Random.seed!(1)
U = rand(solver.model.m, solver.N)

res, stats = TrajectoryOptimization.solve(solver,U)
@assert max_violation(res) < opts.constraint_tolerance

newton_cost, packV, unpackV, cI, cE, d, active_set = gen_newton_al_functions(solver)
max_c(V) = max(maximum(cI(V)),norm(cE(V),Inf),norm(d(V),Inf))
V = init_newton_results(res,solver)

cost(solver,res)
newton_cost(V)
max_violation(res)
max_c(V)
g = ForwardDiff.gradient(newton_cost,V)
H = ForwardDiff.hessian(newton_cost,V)
a = active_set(V,1.0)
rank(H[a,a])
size(H[a,a])
V[a] = V[a] - (H[a,a] + 0I)\g[a]
# V = V - H\g
newton_cost(V)
max_c(V)

# dubinscar
model,obj = Dynamics.dubinscar_parallelpark
solver = Solver(model,obj,N=21)
solver.opts.cost_tolerance = 1e-4
solver.opts.constraint_tolerance = 1e-3
solver.opts.penalty_max = 1e4
n,m,N = get_sizes(solver)
Random.seed!(1)
U0 = rand(m,N-1)
res,stats = solve(solver,U0)
cost(solver,res)
λ_update_default!(res,solver);
update_constraints!(res,solver)
@assert max_violation(res) < solver.opts.constraint_tolerance
ρ = maximum(maximum.(res.μ))
maximum(vcat(res.C...) .* vcat(res.active_set...))

newton_cost, packV, unpackV, cI, cE, d, active_set = gen_newton_al_functions(solver)
mycost(V) = newton_cost(V,1e4)
max_c(V) = max(maximum(cI(V)),norm(cE(V),Inf),norm(d(V),Inf))
# c_all = [cI(V); abs.(cE(V)); abs.(d(V))]
# argmax(c_all)
V = init_newton_results(res,solver)

cost(solver,res)
mycost(V)
max_violation(res)
max_c(V)
g = ForwardDiff.gradient(mycost,V)
H = ForwardDiff.hessian(mycost,V)
a = active_set(V,1e-3)
rank(H)
rank(H[a,a])
size(H[a,a])
minimum(abs.(eigvals(H)))
# V[a] - H[a,a]\g[a]
V1 = line_search(V,H+1000I*0,g,a)
# V = V - H\g
mycost(V1)
max_c(V1)
ForwardDiff.gradient(mycost,V1)
copyto!(V,V1)
max_c(V1)

J = Float64[]
c_max = Float64[]

## quadrotor
model,obj = Dynamics.quadrotor
n = 13
m = 4

# Unconstrained
Q = (1e-3)*Matrix(I,n,n)
R = (1e-2)*Matrix(I,m,m)
# R = (1e-1)*Matrix(I,m,m)
Qf = (1000.0)*Matrix(I,n,n)
tf = 1.0
dt = 0.1

# -initial state
x0 = zeros(n)
x0[1:3] = [0.; 0.; 0.]
q0 = [1.;0.;0.;0.]
x0[4:7] = q0

# -final state
xf = copy(x0)
xf[1:3] = [0.;10.;0.] # xyz position
xf[4:7] = q0
u_max = 10.0
u_min = 0.0

obj_uncon = LQRObjective(Q, R, Qf, tf, x0, xf)
obj_con = TrajectoryOptimization.ConstrainedObjective(obj_uncon,u_min=u_min,u_max=u_max)#,u_min=u_min)
solver = Solver(model,obj_con,N=10)

solver.obj
solver.opts.cost_tolerance = 1e-5
solver.opts.constraint_tolerance = 1e-4
solver.opts.penalty_max = 1e5
n,m,N = get_sizes(solver)
Random.seed!(1)
U = 0.5*9.81/4.0*ones(solver.model.m, solver.N-1)
res,stats = solve(solver,U)
plot(res.U)
plot(res.X)
cost(solver,res)
λ_update_default!(res,solver);
update_constraints!(res,solver)
max_violation(res)
@assert max_violation(res) < solver.opts.constraint_tolerance
ρ = maximum(maximum.(res.μ))
maximum(vcat(res.C...) .* vcat(res.active_set...))

newton_cost, packV, unpackV, cI, cE, d, active_set = gen_newton_al_functions(solver)
mycost(V) = newton_cost(V,ρ)
max_c(V) = max(maximum(cI(V)),norm(cE(V),Inf),norm(d(V),Inf))
# max_c(V) = max(maximum(cE(V)),norm(d(V),Inf))

# c_all = [cI(V); abs.(cE(V)); abs.(d(V))]
# argmax(c_all)
V = init_newton_results(res,solver)

cost(solver,res)
mycost(V)
max_violation(res)
max_c(V)
g = ForwardDiff.gradient(mycost,V)
H = ForwardDiff.hessian(mycost,V)
a = active_set(V,1e-3)
X,U,λ,μ,ν,s = unpackV(V)
13*10 + 4*9 + length(λ) + length(ν)# + 36*2
sum(a)
rank(H)
size(H)
rank(H[a,a])
size(H[a,a])
rank(H[1:309,1:309])
# V[a] - H[a,a]\g[a]
V1 = line_search(V,H,g)
# V = V - H\g
mycost(V1)
max_c(V1)
ForwardDiff.gradient(mycost,V1)
copyto!(V,V1)
max_c(V1)
minimum(abs.(eigvals(H)))
J = Float64[]
c_max = Float64[]
for i = 1:10
    g = ForwardDiff.gradient(mycost,V)
    H = ForwardDiff.hessian(mycost,V)
    a = active_set(V,1e-3)
    V1 = line_search(V,H+0I,g,a)

    push!(J,mycost(V1))
    push!(c_max,max_c(V1))
    push!(grad, norm(ForwardDiff.gradient(mycost,V1)))
    copyto!(V,V1)
end

plot(-diff(J),yscale=:log10)
plot(c_max,yscale=:log10)
plot(grad,yscale=:log10)


text:
@testset "ch12_dp_kf" begin
    dt = 0.1
    # environment
    xlim = [-5.0, 5.0]
    ylim = [-5.0, 5.0]
    # id of landmark must start from 0 with 1 step
    landmarks =
        [Landmark([1.0, 4.0], 0), Landmark([4.0, 1.0], 1), Landmark([-4.0, -4.0], 2)]
    envmap = Map()
    push!(envmap, landmarks)
    world = PuddleWorld(xlim, ylim)
    push!(world, envmap)
    # goal
    goal = Goal(-3.0, -3.0)
    push!(world, goal)
    # robot side
    initial_pose = [2.0, 2.0, 0.0]
    estimator = KalmanFilter(envmap, initial_pose)
    reso = [0.1, 0.1, pi / 20]
    dp_agent = DpPolicyAgent(0.2, 10.0 * pi / 180, dt, estimator, goal, reso)
    init_policy(dp_agent, "ch12_policy.txt")
    robot = RealRobot(initial_pose, dp_agent, RealCamera(landmarks); color = "red")
    push!(world, robot)
    # puddles
    push!(world, Puddle([-2.0, 0.0], [0.0, 2.0], 0.1))
    push!(world, Puddle([-0.5, -2.0], [2.5, 1.0], 0.1))

    anim = @animate for i = 1:100
        t = dt * i
        annota = "t = $(round(t, sigdigits=3))[s]"
        # t
        update_status(world)
        z = observations(robot.sensor_, robot.pose_; noise = true, bias = true)
        p = draw(world, annota)

        # t+1
        v, ω = decision(dp_agent, z)
        state_transition(robot, v, ω, dt; move_noise = true, vel_bias_noise = true)
    end
    if GUI
        gif(anim, "ch12_dp_kf.gif", fps = 10)
    end
end

@testset "ch12_dp_mcl" begin
    dt = 0.1
    # environment
    xlim = [-5.0, 5.0]
    ylim = [-5.0, 5.0]
    # id of landmark must start from 0 with 1 step
    landmarks =
        [Landmark([1.0, 4.0], 0), Landmark([4.0, 1.0], 1), Landmark([-4.0, -4.0], 2)]
    envmap = Map()
    push!(envmap, landmarks)
    world = PuddleWorld(xlim, ylim)
    push!(world, envmap)
    # goal
    goal = Goal(-3.0, -3.0)
    push!(world, goal)
    # robot side
    initial_pose = [2.0, 2.0, 0.0]
    # estimator = KalmanFilter(envmap, initial_pose)
    estimator = Mcl(initial_pose, 100)
    reso = [0.1, 0.1, pi / 20]
    dp_agent = DpPolicyAgent(0.2, 10.0 * pi / 180, dt, estimator, goal, reso)
    init_policy(dp_agent, "ch12_policy.txt")
    robot = RealRobot(initial_pose, dp_agent, RealCamera(landmarks); color = "red")
    push!(world, robot)
    # puddles
    push!(world, Puddle([-2.0, 0.0], [0.0, 2.0], 0.1))
    push!(world, Puddle([-0.5, -2.0], [2.5, 1.0], 0.1))

    anim = @animate for i = 1:180
        t = dt * i
        annota = "t = $(round(t, sigdigits=3))[s]"
        # t
        update_status(world)
        z = observations(robot.sensor_, robot.pose_; noise = true, bias = true)
        p = draw(world, annota)

        # t+1
        v, ω = decision(dp_agent, z, envmap; resample = true)
        state_transition(robot, v, ω, dt; move_noise = true, vel_bias_noise = true)
    end
    if GUI
        gif(anim, "ch12_dp_mcl.gif", fps = 10)
    end
end

@testset "ch12_qmdp12" begin
    dt = 0.1
    # environment
    xlim = [-5.0, 5.0]
    ylim = [-5.0, 5.0]
    # id of landmark must start from 0 with 1 step
    landmarks =
        [Landmark([1.0, 4.0], 0), Landmark([4.0, 1.0], 1), Landmark([-4.0, -4.0], 2)]
    envmap = Map()
    push!(envmap, landmarks)
    world = PuddleWorld(xlim, ylim)
    push!(world, envmap)
    # goal
    goal = Goal(-3.0, -3.0)
    push!(world, goal)
    # robot side
    initial_pose = [2.0, 2.0, 0.0]
    # estimator = KalmanFilter(envmap, initial_pose)
    estimator = Mcl(initial_pose, 100)
    reso = [0.1, 0.1, pi / 20]
    dp_agent = QMDPAgent(0.2, 10.0 * pi / 180, dt, estimator, goal, reso)
    init_policy(dp_agent, "ch12_policy.txt")
    init_value(dp_agent, "ch12_value.txt")
    robot = RealRobot(initial_pose, dp_agent, RealCamera(landmarks); color = "red")
    push!(world, robot)
    # puddles
    push!(world, Puddle([-2.0, 0.0], [0.0, 2.0], 0.1))
    push!(world, Puddle([-0.5, -2.0], [2.5, 1.0], 0.1))

    anim = @animate for i = 1:10
        t = dt * i
        annota = "t = $(round(t, sigdigits=3))[s]"
        # t
        update_status(world)
        z = observations(robot.sensor_, robot.pose_; noise = true, bias = true)
        p = draw(world, annota)

        # t+1
        v, ω = decision(dp_agent, z, envmap; resample = true)
        state_transition(robot, v, ω, dt; move_noise = true, vel_bias_noise = true)
    end
    if GUI
        gif(anim, "ch12_qmdp1.gif", fps = 10)
    end
end

@testset "ch12_amdp2345" begin
    dt = 0.1
    # environment
    xlim = [-5.0, 5.0]
    ylim = [-5.0, 5.0]
    # id of landmark must start from 0 with 1 step
    landmarks =
        [Landmark([1.0, 4.0], 0), Landmark([4.0, 1.0], 1), Landmark([-4.0, -4.0], 2)]
    envmap = Map()
    push!(envmap, landmarks)
    world = PuddleWorld(xlim, ylim)
    push!(world, Puddle([-2.0, 0.0], [0.0, 2.0], 0.1))
    # goal
    goal = Goal(-3.0, -3.0)
    push!(world, goal)
    push!(world, Puddle([-0.5, -2.0], [2.5, 1.0], 0.1))
    push!(world, envmap)
    # robot side
    initial_pose = [2.0, 2.0, 0.0]
    # estimator = KalmanFilter(envmap, initial_pose)
    estimator = Mcl(initial_pose, 100)
    reso = [0.1, 0.1, pi / 20]
    dp_agent = BeliefDP([0.2, 0.2, pi / 18], Goal(-3.0, -3.0); dt = 0.1)
    sampling_num = 10
    init_value(dp_agent)
    init_policy(dp_agent)
    init_depth(dp_agent, world, sampling_num = sampling_num)
    init_state_transition_probs(dp_agent, sampling_num = sampling_num)

    value_iteration_sweep(dp_agent)
    value_iteration_sweep(dp_agent)
    v = dp_agent.value_function_[:, :, 18, 1]

    init_motion_sigma_transition_probs(dp_agent)
    landmarks = [
        Landmark([1.0, 4.0], 0),
        Landmark([4.0, 1.0], 1),
        Landmark([-4.0, 1.0], 2),
        Landmark([-2.0, 1.0], 3),
    ]
    camera = IdealCamera(landmarks)
    init_obs_sigma_transition_probs(dp_agent, camera)
    init_expected_depths(dp_agent, world)
    value_iteration_sweep(dp_agent)
end

@testset "ch12_amdp_policy_agent" begin
    dt = 0.1
    # environment
    xlim = [-5.0, 5.0]
    ylim = [-5.0, 5.0]
    # id of landmark must start from 0 with 1 step
    landmarks = [
        Landmark([-4.0, 2.0], 0),
        Landmark([2.0, -3.0], 1),
        Landmark([4.0, 4.0], 2),
        Landmark([-4.0, -4.0], 3)
    ]
    envmap = Map()
    push!(envmap, landmarks)
    world = PuddleWorld(xlim, ylim)
    push!(world, envmap)
    # goal
    goal = Goal(-3.0, -3.0)
    push!(world, goal)
    # robot side
    initial_pose = [2.0, 2.0, 0.0]
    estimator = KalmanFilter(envmap, initial_pose)
    # this resolution must be consistant with the resolution
    # with which policy_amdp.txt was generated in BeliefDP.
    reso = [0.2, 0.2, pi / 18]
    dp_agent = AMDPPolicyAgent(0.2, 10.0 * pi / 180, dt, estimator, goal, reso)
    sampling_num = 10
    init_policy(dp_agent, "policy_amdp.txt")

    robot = RealRobot(initial_pose, dp_agent, RealCamera(landmarks); color = "red")
    push!(world, robot)
    # puddles
    push!(world, Puddle([-2.0, 0.0], [0.0, 2.0], 0.1))
    push!(world, Puddle([-0.5, -2.0], [2.5, 1.0], 0.1))

    anim = @animate for i = 1:10
        t = dt * i
        annota = "t = $(round(t, sigdigits=3))[s]"
        # t
        update_status(world)
        z = observations(robot.sensor_, robot.pose_; noise = true, bias = true)
        p = draw(world, annota)

        # t+1
        v, ω = decision(dp_agent, z)
        # v, ω = policy(dp_agent, agent.estimator_.pose_)
        state_transition(robot, v, ω, dt; move_noise = true, vel_bias_noise = true)
    end
    if GUI
        gif(anim, "ch12_amdp_policy_agent.gif", fps = 10)
    end
end

max_stars_repo_name: soblin/JuliaProbo

max_issues_repo_name: soblin/JuliaProbo

max_forks_repo_name: soblin/JuliaProbo


text: Dict{String,Any}("multiline" => Dict{String,Any}("value" => "This string\nhas ' a quote character\nand more than\none newline\nin it.","type" => "string"),"firstnl" => Dict{String,Any}("value" => "This string has a ' quote character.","type" => "string"),"oneline" => Dict{String,Any}("value" => "This string has a ' quote character.","type" => "string"))

max_stars_repo_name: vanillajonathan/julia

max_issues_repo_name: vanillajonathan/julia

max_forks_repo_name: vanillajonathan/julia


text:
module QuantumComputing
import Base: convert, show, kron, replace, length, sub, split

export N,E,Id,X,Y,Z,H,CNOT,SWAP,Q
export Ket,Bra
export GGate,GCNot,CNot,Not,Swap,Hadamard,VentureMartinez,M
export append,appendhead,extract

# |0> NULL
# |1> EINS
const N = sparse(complex([1.0,0.0]))
const E = sparse(complex([0.0,1.0]))

# Unit
const Id = complex(speye(2))

# Pauli matrices
const X = sparse(complex([0.0 1.0;1.0 0.0]))
const Y = sparse([0.0 -im; im 0.0])
const Z = sparse(complex([1.0 0.0; 0.0 -1.0]))

# Hadamard matrix
const H = sparse(1.0/sqrt(2)*complex([1.0 1.0; 1.0 -1.0]))

# 2 q-bits operator
const CNOT = blkdiag(Id,X)
const SWAP = hcat(kron(N,N),kron(E,N),kron(N,E),kron(E,E))

const Q = H*N

const sbra = "\u27E8" # <
const svbar = "\u007C" # |
const sket = "\u27E9" # >

abstract QGate

# 0->N, 1->E
function parse01(x)
    if x==0
       return N
    elseif x == 1
       return E
    else
       return x
    end
end

immutable Ket
    v::SparseMatrixCSC
    function Ket(x...)
        v = map(parse01, x)
        new(kron(v...))
    end
    function Ket(v::SparseMatrixCSC)
        new(v)
    end
end

immutable GGate <: QGate
    v::SparseMatrixCSC
    function GGate(v)
        new(v)
    end
end

immutable CNot <: QGate
    v::SparseMatrixCSC
    t::Integer
    c::Integer
    function CNot(n,t,c)
        M1 = Swap(n,1,c).v
        M2 = Swap(n,2,t==1?2:t).v
        v = M1*M2*kron(CNOT,[Id for i=1:n-2]...)*M2*M1
        new(v,t,c)
    end
end

immutable GCNot <: QGate
    v::SparseMatrixCSC
    t::Integer
    c
    function GCNot(n,t,c...)
        nc = length(c)
        v = blkdiag([Id for i=1:2^nc]...)
        v[end-1:end,end-1:end] = X

        Mf = SparseMatrixCSC[]
        c = sort([c...])
        idx=1
        tp = t
        for i in filter(x->x>nc,c)
            while (idx in c)
                idx+=1
                @assert(idx<=nc)
            end
            push!(Mf,Swap(n,idx,i).v)
            if idx==t
                tp = i
            end
            idx+=1
        end
        push!(Mf,Swap(n,nc+1,tp).v)
        v = reduce(*,Mf)*kron(v,[Id for i=1:n-nc-1]...)*reduce(*,reverse(Mf))

        new(v,t,c)
    end
end

immutable Not <: QGate
    v::SparseMatrixCSC
    t::Integer
    function Not(n,t)
        v = [Id for i=1:n]
        v[t] = X
        new(kron(v...),t)
    end
end

immutable Swap <: QGate
    v::SparseMatrixCSC
    t1::Integer
    t2::Integer
    function Swap(n,t1,t2)
        mint = min(t1,t2)
        maxt = max(t1,t2)
        v = kron([Id for i=1:n]...)
        for t=mint:maxt-1
            tv = kron([Id for i=1:t-1]...,SWAP,[Id for i=t+1:n-1]...)
            v *= tv
        end
        for t=maxt-2:-1:mint
            tv = kron([Id for i=1:t-1]...,SWAP,[Id for i=t+1:n-1]...)
            v *= tv
        end
        new(v,t1,t2)
    end
end

immutable Hadamard <: QGate
    v::SparseMatrixCSC
    t::Integer
    function Hadamard(n,t)
        v = [Id for i=1:n]
        v[t] = H
        new(kron(v...),t)
    end
end

immutable VentureMartinez <: QGate
    v::SparseMatrixCSC
    t::Integer
    c::Integer
    j::Integer
    s::Integer
    function VentureMartinez(n,t,c,j,s=-1)
        S = complex(blkdiag(
        speye(2),
        sparse(
        [sqrt((j-1.0)/j) -s/sqrt(j);s/sqrt(j) sqrt((j-1.0)/j)])
        ))
        M1 = Swap(n,1,c).v
        M2 = Swap(n,2,t==1?2:t).v
        v = M1*M2*kron(S, [Id for i=1:n-2]...)*M2*M1
        new(v,t,c,j,s)
    end
end

immutable M <: QGate
    v::SparseMatrixCSC
    m::Complex
    function M(n,m)
        v = complex(spzeros(2^n,2^n))
        v[1,1] = zero(eltype(v))
        for i=2:2^n
            v[i,i] = sqrt(1/(1-m))
        end
        new(v,m)
    end
end

kron(ket1::Ket, ket2::Ket) = Ket(kron(ket1.v, ket2.v))
kron(v) = v
length(ket::Ket) = int(log(size(ket.v,1))/log(2))
show(io::IO, ket::Ket) = print(io,convert(String,ket))

*(val::Number, ket::Ket) = Ket(ket.v*val)
*(ket::Ket, val::Number) = val*ket.v
*(op1::QGate, op2::QGate) = GGate(op1.v*op2.v)
*(ket::Ket, op::QGate) = op*ket
*(op::QGate, ket::Ket) = Ket(op.v*ket.v)
==(ket1::Ket, ket2::Ket) = ket1.v==ket2.v
+(ket1::Ket, ket2::Ket) = Ket(ket1.v+ket2.v)

function convert(::Type{String}, ket::Ket)
    vp = prob(ket.v)
    str=map(find(vp)) do x
        string(@sprintf("%1.13f",vp[x]),"\t", svbar,bin(x-1,length(ket)),sket," \n")
    end
    join(sort(str))
end

function prob(v::SparseMatrixCSC)
    prob(x) = real(x*x)
    map(prob,v)
end

# (|0000>,3) -> |000>
function sub(ket::Ket,n)
    nk = length(ket)
    @assert(n<=nk)
    v = prob(ket.v)
    for i=1:nk-n
        v = reshape(v,2,length(v)>>1)'
        v = v[:,1] + v[:,2]
    end
    kets=map(find(v)) do x
        sqrt(v[x]) * Ket(map(parseint,split(bin(x-1,n),""))...)
    end
    reduce(+,kets)
end

# (|abcd>,3) -> |c>
function extract(ket::Ket,n)
    sub(Swap(length(ket),1,n)*ket,1)
end

#  (|0000>,a,b,c) -> |0000abc>
function append(ket::Ket, x...)
    av = map(parse01, x)
    Ket(kron(ket.v, av...))
end

#  (|0000>,a,b,c) -> |abc0000>
function appendhead(ket::Ket, x...)
    nk = length(ket)
    nx = length(x)
    av = map(parse01, x)

    v = ket.v
    for i=1:nx
        v = kron(v, av[i])
        nk = nk+1
        v = Swap(nk,i,nk).v*v
        v = reduce(*, [Swap(nk,j,j+1).v for j=i+1:nk-1])*v
    end

    Ket(v)
end

end # module


text:
module CallableCollections

export ∘

typealias IndexedCollection Union{AbstractArray, Associative}
typealias SetCollection Union{Set, IntSet}

Base.call(arr::IndexedCollection, ind...) = getindex(arr, ind...)
Base.call(set::SetCollection, ind) = ind ∈ set

Base.union(fn1, fn2) = x -> fn1(x) || fn2(x)
Base.intersect(fn1, fn2) = x -> fn1(x) && fn2(x)

Base.map(fn, fns...) = function(xs...)
    fn([f(xs...) for f in fns]...)
end

∘(f, g) = x -> f(g(x))

end  # module


text:
function greet()
    println("Welcome to JHmat...")
    return 1
end


text:
using Revise
using westbrick

obs_dict, ax, fig = visualize_environment!([])
bobby = load_robot()
traj, obj_traj = init_visualization()

pose = [5.0, 4., 0.0]
move!(bobby, pose, traj, obj_traj)

pose = [1.0, 5.0, 0.0]
move!(bobby, pose, traj, obj_traj)

pose = [3., 5.0, 0.0]
move!(bobby, pose, traj, obj_traj)

pose = [0, 5.0, 0.0]
move!(bobby, pose, traj, obj_traj)

visualize_trajectory!(bobby, traj, obj_traj, obs_dict, ax, fig;name="media/move_test.gif")


text:
"""
    AbstractHistory{L,M}

Abstract supertype for event histories with locations of type `L` and marks of type `M`.
"""
abstract type AbstractHistory{L,M} end


text:
module TestDyingPhilosopher

sleep(x) = SLEEP_IMPL[](x)
const SLEEP_IMPL = Ref{Any}(Base.sleep)

println(x) = PRINTLN_IMPL[](x)
const PRINTLN_IMPL = Ref{Any}(Base.println)

include("../../../examples/dying_philosopher.jl")

using Test

function patch_with(f, ref, impl)
    old = ref[]
    ref[] = impl
    try
        f()
    finally
        ref[] = old
    end
end

function run_dying_philosopher(nphils = 5)
    outputs = Channel(Inf)  # not using Julio's to avoid cancellation
    function println_impl(x)
        put!(outputs, x)
    end

    nsleeps = Threads.Atomic{Int}(0)
    alllocked = Julio.Promise{Nothing}()
    function sleep_impl(time)
        if time <= 1
            # sleep(rand()) call
            if Threads.atomic_add!(nsleeps, 1) + 1 == nphils
                alllocked[] = nothing
            else
                alllocked[]
            end
        else
            # sleep(10) call
            alllocked[]
        end
    end

    patch_with(SLEEP_IMPL, sleep_impl) do
        patch_with(PRINTLN_IMPL, println_impl) do
            dying_philosopher(nphils)
        end
    end

    close(outputs)
    return collect(outputs)
end

function test_dying_philosopher(nphils = 5)
    outputs = run_dying_philosopher(nphils)
    @test startswith(outputs[1], "Killing")
    @test startswith(outputs[2], "But what is death?")
    @test all(startswith.(outputs[3:end][1:nphils-1], "eating"))
end

end  # module


text:
int main() {
   return f(3);
}

int f(int x) {
    if (x<0)
       return x;
}

max_stars_repo_name: carlostome/Javalette-Compiler

max_issues_repo_name: carlostome/Javalette-Compiler

max_forks_repo_name: carlostome/Javalette-Compiler


text:
#!/usr/bin/env julia
using StatsBase

# The number of syllables per word you want to generate: #
sylnum = 2
# The number of words you want to generate: #
outputlines = 3000

# Function for random choice: #
function random(chars, weights)
	normWeights = ProbabilityWeights(weights)
	choice = sample(chars, normWeights)
	return choice
end

# Writes a wordlist: #

function make_words(outputlines)
	words = fill("", outputlines)
	# Do this as many times as the number of words you want in the output: #
	for i in 1:outputlines
		# Do this for each syllable: #
		oldsyl = ""
		for j in 1:sylnum
			syl = ""
			# Generate nucleus: #
			vowels = ["a", "e", "i", "o", "u"]
			vowelweights = [38.56, 13.96, 19.13, 13.11, 15.24]
			nuc = random(vowels, vowelweights)
			# Generate onset: #
			# First syllable: #
			if oldsyl == ""
				onsc1 = ["", "m", "n", "ny", "p", "t", "k", "`", "b", "d", "g", "f", "s", "z", "x", "c", "l", "r", "y", "w", "h"]
				onsc1weights = [17.69, 9.09, 5.25, 0.87, 4.67, 7.44, 10.00, 0.82, 5.69, 4.35, 3.50, 1.90, 4.48, 1.03, 1.14, 2.77, 4.43, 2.07, 2.78, 2.49, 3.68]
				onsc1 = random(onsc1, onsc1weights)
			end
			# After first syllable: #
			if oldsyl != ""
				if oldsyl[1] == "p" || oldsyl[1] == "b" || oldsyl[1] == "f"
					onsc1= ["", "m", "n", "ny", "p", "t", "k", "`", "b", "d", "g", "f", "s", "z", "x", "c", "l", "r", "y", "w", "h"]
					onsc1weights = [10.29, 6.55, 8.62, 1.05, 0.3, 9.25, 7.36, 1.64, 0.48, 2.89, 3.72, 0.17, 4.20, 0.64, 0.97, 1.33, 9.27, 9.65, 2.25, 2.18, 3.93]
					onsc1 = random(onsc1, onsc1weights)
				elseif oldsyl[1] == "t" || oldsyl[1] == "d" || oldsyl[1] == "s"
					onsc1 = ["", "m", "n", "ny", "p", "t", "k", "`", "b", "d", "g", "f", "s", "z", "x", "c", "l", "r", "y", "w", "h"]
					onsc1weights = [9.47, 6.55, 8.62, 1.05, 3.05, 0.93, 7.36, 1.64, 4.81, 0.29, 3.72, 1.67, 0.42, 0.64, 0.97, 1.33, 9.27, 9.65, 2.25, 2.18, 3.93]
					onsc1 = random(onsc1, onsc1weights)
				elseif oldsyl[1] == "k" || oldsyl[1] == "g"
					onsc1 = ["", "m", "n", "ny", "p", "t", "k", "`", "b", "d", "g", "f", "s", "z", "x", "c", "l", "r", "y", "w", "h"]
					onsc1weights = [10.11, 6.55, 8.62, 1.05, 3.05, 9.25, 0.74, 1.64, 4.81, 2.89, 0.37, 1.67, 4.20, 0.64, 0.97, 1.33, 9.27, 9.65, 2.25, 2.18, 3.93]
					onsc1 = random(onsc1, onsc1weights)
				elseif oldsyl[1] == "m" || oldsyl[1] == "n"
					onsc1 = ["", "m", "n", "ny", "p", "t", "k", "`", "b", "d", "g", "f", "s", "z", "x", "c", "l", "r", "y", "w", "h"]
					onsc1weights = [9.49, 0.66, 0.86, 0.11, 3.05, 9.25, 7.36, 1.64, 4.81, 2.89, 3.72, 1.67, 4.20, 0.64, 0.97, 1.33, 9.27, 9.65, 2.25, 2.18, 3.93]
					onsc1 = random(onsc1, onsc1weights)
				elseif oldsyl[1] == "l" || oldsyl[1] == "r"
					onsc1 = ["", "m", "n", "ny", "p", "t", "k", "`", "b", "d", "g", "f", "s", "z", "x", "c", "l", "r", "y", "w", "h"]
					onsc1weights = [9.16, 6.55, 8.62, 1.05, 3.05, 9.25, 7.36, 1.64, 4.81, 2.89, 3.72, 1.67, 4.20, 0.64, 0.97, 1.33, 0.93, 0.97, 2.25, 2.18, 3.93]
					onsc1 = random(onsc1, onsc1weights)
				elseif oldsyl[1] == "z" || oldsyl[1] == "x" || oldsyl[1] == "c"
					onsc1 = ["", "m", "n", "ny", "p", "t", "k", "`", "b", "d", "g", "f", "s", "z", "x", "c", "l", "r", "y", "w", "h"]
					onsc1weights = [11.09, 6.55, 8.62, 1.05, 3.05, 9.25, 7.36, 1.64, 4.81, 2.89, 3.72, 1.67, 4.20, 0.06, 0.09, 0.13, 9.27, 9.65, 2.25, 2.18, 3.93]
					onsc1 = random(onsc1, onsc1weights)
				elseif oldsyl[1] == "h" || oldsyl[1] == "`"
					onsc1 = ["", "m", "n", "ny", "p", "t", "k", "`", "b", "d", "g", "f", "s", "z", "x", "c", "l", "r", "y", "w", "h"]
					onsc1weights = [10.78, 6.55, 8.62, 1.05, 3.05, 9.25, 7.36, 0.16, 4.81, 2.89, 3.72, 1.67, 4.20, 0.64, 0.97, 1.33, 9.27, 9.65, 2.25, 2.18, 0.39]
					onsc1 = random(onsc1, onsc1weights)
				elseif oldsyl[1] == "w" || oldsyl[1] == "y"
					onsc1 = ["", "m", "n", "ny", "p", "t", "k", "`", "b", "d", "g", "f", "s", "z", "x", "c", "l", "r", "y", "w", "h"]
					onsc1weights = [10.92, 6.55, 8.62, 1.05, 3.05, 9.25, 7.36, 1.64, 4.81, 2.89, 3.72, 1.67, 4.20, 0.64, 0.97, 1.33, 9.27, 9.65, 0.22, 0.22, 3.93]
					onsc1 = random(onsc1, onsc1weights)
				else
					onsc1 = ["", "m", "n", "ny", "p", "t", "k", "`", "b", "d", "g", "f", "s", "z", "x", "c", "l", "r", "y", "w", "h"]
					onsc1weights = [11.45, 6.55, 8.62, 1.05, 3.05, 9.25, 7.36, 1.64, 4.81, 2.89, 3.72, 1.67, 4.20, 0.64, 0.97, 1.33, 9.27, 9.65, 2.25, 2.18, 3.93]
					onsc1 = random(onsc1, onsc1weights)
				end
			end
			# Generate coda: #
			# Before penultimate syllable: #
			if j < sylnum-1
				if length(oldsyl) > 2
					if oldsyl[3] == "m" || oldsyl[3] == "n" || oldsyl[3] == "ng"
						coda = ["", "m", "n", "ng", "p", "t", "k", "s", "l", "r", "y", "w"]
						codaweights = [201.68, 0.66, 2.12, 0.74, 3.03, 8.53, 8.41, 3.72, 7.03, 7.17, 6.46, 2.55]
						coda = random(coda, codaweights)
						if !(nuc[1] in ("a", "e")) && coda in ("y", "w")
							coda = ""
						end
					elseif oldsyl[3] == "p"
						coda = ["", "m", "n", "ng", "p", "t", "k", "s", "l", "r", "y", "w"]
						codaweights = [317.44, 6.63, 21.18, 7.38, 0.3, 8.53, 8.41, 3.72, 7.03, 7.17, 6.46, 2.55]
						coda = random(coda, codaweights)
						if !(nuc[1] in ("a", "e")) && coda in ("y", "w")
							coda = ""
						end
					elseif oldsyl[3] == "t"
						coda = ["", "m", "n", "ng", "p", "t", "k", "s", "l", "r", "y", "w"]
						codaweights = [297.64, 6.63, 21.18, 7.38, 3.03, 0.85, 8.41, 3.72, 7.03, 7.17, 6.46, 2.55]
						coda = random(coda, codaweights)
						if !(nuc[1] in ("a", "e")) && coda in ("y", "w")
							coda = ""
						end
					elseif oldsyl[3] == "k"
						coda = ["", "m", "n", "ng", "p", "t", "k", "s", "l", "r", "y", "w"]
						codaweights = [298.08, 6.63, 21.18, 7.38, 3.03, 8.53, 0.84, 3.72, 7.03, 7.17, 6.46, 2.55]
						coda = random(coda, codaweights)
						if !(nuc[1] in ("a", "e")) && coda in ("y", "w")
							coda = ""
						end
					elseif oldsyl[3] == "s"
						coda = ["", "m", "n", "ng", "p", "t", "k", "s", "l", "r", "y", "w"]
						codaweights = [314.96, 6.63, 21.18, 7.38, 3.03, 8.53, 8.41, 0.37, 7.03, 7.17, 6.46, 2.55]
						coda = random(coda, codaweights)
						if !(nuc[1] in ("a", "e")) && coda in ("y", "w")
							coda = ""
						end
					elseif oldsyl[3] == "l" || oldsyl[3] == "r"
						coda = ["", "m", "n", "ng", "p", "t", "k", "s", "l", "r", "y", "w"]
						codaweights = [277.24, 6.63, 21.18, 7.38, 3.03, 8.53, 8.41, 3.72, 0.7, 0.72, 6.46, 2.55]
						coda = random(coda, codaweights)
						if !(nuc[1] in ("a", "e")) && coda in ("y", "w")
							coda = ""
						end
					elseif oldsyl[3] == "y" || oldsyl[3] == "w"
						coda = ["", "m", "n", "ng", "p", "t", "k", "s", "l", "r", "y", "w"]
						codaweights = [295.96, 6.63, 21.18, 7.38, 3.03, 8.53, 8.41, 3.72, 7.03, 7.17, 0.65, 0.26]
						coda = random(coda, codaweights)
						if !(nuc[1] in ("a", "e")) && coda in ("y", "w")
							coda = ""
						end
					else
						coda = ["", "m", "n", "ng", "p", "t", "k", "s", "l", "r", "y", "w"]
						codaweights = [328.36, 6.63, 21.18, 7.38, 3.03, 8.53, 8.41, 3.72, 7.03, 7.17, 6.46, 2.55]
						coda = random(coda, codaweights)
						if !(nuc[1] in ("a", "e")) && coda in ("y", "w")
							coda = ""
						end
					end
				else
					coda = ["", "m", "n", "ng", "p", "t", "k", "s", "l", "r", "y", "w"]
					codaweights = [328.36, 6.63, 21.18, 7.38, 3.03, 8.53, 8.41, 3.72, 7.03, 7.17, 6.46, 2.55]
					coda = random(coda, codaweights)
					if !(nuc[1] in ("a", "e")) && coda in ("y", "w")
						coda = ""
					end
				end
			end
			# Penultimate syllable: #
			if j == sylnum-1
				if length(oldsyl) > 2
					if oldsyl[3] == "m" || oldsyl[3] == "n" || oldsyl[3] == "ng"
						coda = ["", "m", "n", "ng", "p", "t", "k", "s", "l", "r", "y", "w"]
						codaweights = [117.65, 0.66, 2.12, 0.74, 3.03, 8.53, 8.41, 3.72, 7.03, 7.17, 6.46, 2.55]
						coda = random(coda, codaweights)
						if !(nuc[1] in ("a", "e")) && coda in ("y", "w")
							coda = ""
						end
					elseif oldsyl[3] == "p"
						coda = ["", "m", "n", "ng", "p", "t", "k", "s", "l", "r", "y", "w"]
						codaweights = [185.17, 6.63, 21.18, 7.38, 0.3, 8.53, 8.41, 3.72, 7.03, 7.17, 6.46, 2.55]
						coda = random(coda, codaweights)
						if !(nuc[1] in ("a", "e")) && coda in ("y", "w")
							coda = ""
						end
					elseif oldsyl[3] == "t"
						coda = ["", "m", "n", "ng", "p", "t", "k", "s", "l", "r", "y", "w"]
						codaweights = [173.62, 6.63, 21.18, 7.38, 3.03, 0.85, 8.41, 3.72, 7.03, 7.17, 6.46, 2.55]
						coda = random(coda, codaweights)
						if !(nuc[1] in ("a", "e")) && coda in ("y", "w")
							coda = ""
						end
					elseif oldsyl[3] == "k"
						coda = ["", "m", "n", "ng", "p", "t", "k", "s", "l", "r", "y", "w"]
						codaweights = [173.88, 6.63, 21.18, 7.38, 3.03, 8.53, 0.84, 3.72, 7.03, 7.17, 6.46, 2.55]
						coda = random(coda, codaweights)
						if !(nuc[1] in ("a", "e")) && coda in ("y", "w")
							coda = ""
						end
					elseif oldsyl[3] == "s"
						coda = ["", "m", "n", "ng", "p", "t", "k", "s", "l", "r", "y", "w"]
						codaweights = [183.73, 6.63, 21.18, 7.38, 3.03, 8.53, 8.41, 0.37, 7.03, 7.17, 6.46, 2.55]
						coda = random(coda, codaweights)
						if !(nuc[1] in ("a", "e")) && coda in ("y", "w")
							coda = ""
						end
					elseif oldsyl[3] == "l" || oldsyl[3] == "r"
						coda = ["", "m", "n", "ng", "p", "t", "k", "s", "l", "r", "y", "w"]
						codaweights = [161.72, 6.63, 21.18, 7.38, 3.03, 8.53, 8.41, 3.72, 0.7, 0.72, 6.46, 2.55]
						coda = random(coda, codaweights)
						if !(nuc[1] in ("a", "e")) && coda in ("y", "w")
							coda = ""
						end
					elseif oldsyl[3] == "y" || oldsyl[3] == "w"
						coda = ["", "m", "n", "ng", "p", "t", "k", "s", "l", "r", "y", "w"]
						codaweights = [172.64, 6.63, 21.18, 7.38, 3.03, 8.53, 8.41, 3.72, 7.03, 7.17, 0.65, 0.26]
						coda = random(coda, codaweights)
						if !(nuc[1] in ("a", "e")) && coda in ("y", "w")
							coda = ""
						end
					else
						coda = ["", "m", "n", "ng", "p", "t", "k", "s", "l", "r", "y", "w"]
						codaweights = [191.54, 6.63, 21.18, 7.38, 3.03, 8.53, 8.41, 3.72, 7.03, 7.17, 6.46, 2.55]
						coda = random(coda, codaweights)
						if !(nuc[1] in ("a", "e")) && coda in ("y", "w")
							coda = ""
						end
					end
				else
					coda = ["", "m", "n", "ng", "p", "t", "k", "s", "l", "r", "y", "w"]
					codaweights = [191.54, 6.63, 21.18, 7.38, 3.03, 8.53, 8.41, 3.72, 7.03, 7.17, 6.46, 2.55]
					coda = random(coda, codaweights)
					if !(nuc[1] in ("a", "e")) && coda in ("y", "w")
						coda = ""
					end
				end
			end
			# Final syllable: #
			if j == sylnum
				if length(oldsyl) > 2
					if oldsyl[3] == "m" || oldsyl[3] == "n" || oldsyl[3] == "ng"
						coda = ["", "m", "n", "ng", "p", "t", "k", "s", "l", "r", "y", "w"]
						codaweights = [89.25, 0.66, 2.12, 0.74, 3.03, 8.53, 8.41, 3.72, 7.03, 7.17, 6.46, 2.55]
						coda = random(coda, codaweights)
						if !(nuc[1] in ("a", "e")) && coda in ("y", "w")
							coda = ""
						end
					elseif oldsyl[3] == "p"
						coda = ["", "m", "n", "ng", "p", "t", "k", "s", "l", "r", "y", "w"]
						codaweights = [140.45, 6.63, 21.18, 7.38, 0.3, 8.53, 8.41, 3.72, 7.03, 7.17, 6.46, 2.55]
						coda = random(coda, codaweights)
						if !(nuc[1] in ("a", "e")) && coda in ("y", "w")
							coda = ""
						end
					elseif oldsyl[3] == "t"
						coda = ["", "m", "n", "ng", "p", "t", "k", "s", "l", "r", "y", "w"]
						codaweights = [131.71, 6.63, 21.18, 7.38, 3.03, 0.85, 8.41, 3.72, 7.03, 7.17, 6.46, 2.55]
						coda = random(coda, codaweights)
						if !(nuc[1] in ("a", "e")) && coda in ("y", "w")
							coda = ""
						end
					elseif oldsyl[3] == "k"
						coda = ["", "m", "n", "ng", "p", "t", "k", "s", "l", "r", "y", "w"]
						codaweights = [131.9, 6.63, 21.18, 7.38, 3.03, 8.53, 0.84, 3.72, 7.03, 7.17, 6.46, 2.55]
						coda = random(coda, codaweights)
						if !(nuc[1] in ("a", "e")) && coda in ("y", "w")
							coda = ""
						end
					elseif oldsyl[3] == "s"
						coda = ["", "m", "n", "ng", "p", "t", "k", "s", "l", "r", "y", "w"]
						codaweights = [139.38, 6.63, 21.18, 7.38, 3.03, 8.53, 8.41, 0.37, 7.03, 7.17, 6.46, 2.55]
						coda = random(coda, codaweights)
						if !(nuc[1] in ("a", "e")) && coda in ("y", "w")
							coda = ""
						end
					elseif oldsyl[3] == "l" || oldsyl[3] == "r"
						coda = ["", "m", "n", "ng", "p", "t", "k", "s", "l", "r", "y", "w"]
						codaweights = [122.68, 6.63, 21.18, 7.38, 3.03, 8.53, 8.41, 3.72, 0.7, 0.72, 6.46, 2.55]
						coda = random(coda, codaweights)
						if !(nuc[1] in ("a", "e")) && coda in ("y", "w")
							coda = ""
						end
					elseif oldsyl[3] == "y" || oldsyl[3] == "w"
						coda = ["", "m", "n", "ng", "p", "t", "k", "s", "l", "r", "y", "w"]
						codaweights = [130.97, 6.63, 21.18, 7.38, 3.03, 8.53, 8.41, 3.72, 7.03, 7.17, 0.65, 0.26]
						coda = random(coda, codaweights)
						if !(nuc[1] in ("a", "e")) && coda in ("y", "w")
							coda = ""
						end
					else
						coda = ["", "m", "n", "ng", "p", "t", "k", "s", "l", "r", "y", "w"]
						codaweights = [145.3, 6.63, 21.18, 7.38, 3.03, 8.53, 8.41, 3.72, 7.03, 7.17, 6.46, 2.55]
						coda = random(coda, codaweights)
						if !(nuc[1] in ("a", "e")) && coda in ("y", "w")
							coda = ""
						end
					end
				else
					coda = ["", "m", "n", "ng", "p", "t", "k", "s", "l", "r", "y", "w"]
					codaweights = [145.3, 6.63, 21.18, 7.38, 3.03, 8.53, 8.41, 3.72, 7.03, 7.17, 6.46, 2.55]
					coda = random(coda, codaweights)
					if !(nuc[1] in ("a", "e")) && coda in ("y", "w")
						coda = ""
					end
				end
			end
			syl = syl * onsc1 * nuc * coda * "."
			oldsyl = syl
			words[i] = words[i] * syl
		end
		# Write to output file: #

	end
	open("output.txt", "w") do f
		for w in words
			write(f, w)
			write(f, "\n")
		end
	end
end

make_words(outputlines)


text:
module TestCorr
using PlmDCA, Test

function corrtests_simple()
    W,Z,N,M,q = PlmDCA.ReadFasta("test.fasta",1.0, 0.0 , true)
    @test q == 2
    @test N == 3
    @test M == 2
    @test Z == [1 1; 2 2; 1 2]
    @test W == [0.5, 0.5]
end
corrtests_simple()
printstyled("All TestCorr passed!\n",color=:light_green,bold=true)
end # end module TestGauge

max_stars_repo_name: giovact/PlmDCA

max_issues_repo_name: giovact/PlmDCA

max_forks_repo_name: giovact/PlmDCA


text:
using Random, LinearAlgebra

function strictlyuppertriangular(n)
    M = (rand(n, n) .+ 0.5)
    for c in CartesianIndices(M)
        if c[1] ≤ c[2]
            M[c] = 0
        end
    end
    return M
end

Random.seed!(1234)

n = 5
ratio = (1, 1)
σ = identity
x₀ = rand(n) .+ 0.5
W = strictlyuppertriangular(n)

expensive_sequence = [x₀]
let xcurrent = x₀
    for i in 1:(n-1)
        xnew = W * σ.(xcurrent) + x₀
        push!(expensive_sequence, (xnew - xcurrent))
        xcurrent = xnew
    end
end

optimized_sequence = [x₀]
let xcurrent = x₀
    for i in 1:(n-1)
        xnew = copy(xcurrent)
        xnew[i+1] += W[i+1, :] ⋅ σ.(xcurrent)
        push!(optimized_sequence, (xnew - xcurrent))
        xcurrent = xnew
    end
end

using AlgebraOfGraphics, CairoMakie
set_aog_theme!()

y = hcat(optimized_sequence, expensive_sequence)
x = map(_ -> 0:n-1, y)

plt = mapping(
    x => "node",
    y => "activation value",
    stack=0:n-1 => nonnumeric => "step",
    color=0:n-1 => nonnumeric => "step",
    row=[1 2] => renamer(1 => "optimized sequence", 2 => "expensive sequence")
)

##

using GraphMakie, CairoMakie, Graphs
using AlgebraOfGraphics
set_aog_theme!()

font = AlgebraOfGraphics.firasans("Light")

fig = Figure(resolution=(600, 600))

ax1 = Axis(fig[1, 1:ratio[1]], aspect=1)
ax2 = Axis(fig[2, 1:ratio[1]], aspect=1)

colormap = collect(cgrad(:sienna, n, categorical=true, rev=true))

g = complete_graph(n)
colors = map(edges(g)) do edge
    colormap[dst(edge)]
end

arrow_size = 25
edge_width_base = 3
node_size=45

function layout(g::AbstractGraph)
    N = nv(g)
    return map(1:N) do i
        θ = 2π*(i - 1)/N
        return Point(cos(π - θ), sin(π - θ))
    end
end

graphplot!(
    ax1, g, arrow_show=true, arrow_size=arrow_size, edge_width=edge_width_base, edge_color=colors,
    nlabels=map(string, 0:n-1), nlabels_align=(:center, :center), nlabels_attr=(; font, textsize=16f0),
    node_color=:white, node_size=node_size, layout=layout
)

edge_width = map(edges(g)) do edge
    edge_width_base * src(edge)
end

graphplot!(
    ax2, g, edge_width=edge_width, edge_color=colors,
    nlabels=map(string, 0:n-1), nlabels_align=(:center, :center), nlabels_attr=(; font, textsize=16f0),
    node_color=:white, node_size=node_size, layout=layout
)

foreach(hidedecorations!, (ax1, ax2))
foreach(hidespines!, (ax1, ax2))

ag = draw!(
    fig[:, ratio[1] .+ (1:ratio[2])],
    plt * visual(BarPlot),
    palettes=(color=colormap,),
    axis = (xticks=0:n-1,),
)

legend!(fig[0, 1:sum(ratio)], ag, orientation=:horizontal)

display(fig)


text:
module EarCut

using GeometryBasics

using EarCut_jll

include("cwrapper.jl")

export triangulate

end # module


text:
using Base.Test
using QuantumOptics

N = 3
Ncutoff = 2
T = [0.:0.1:1.;]

Ω = [0 2 3;
     2 0 1;
     3 1 0]

ω = [1., 1.2, 1.5]

basis_fock = FockBasis(Ncutoff)
basis = tensor([basis_fock for i=1:N]...)

a = destroy(basis_fock)
at = create(basis_fock)
I = identity(basis_fock)

rho0 = cumulantexpansion.ProductDensityOperator([coherentstate(basis_fock, i%Ncutoff) for i=1:N]...)

J = LazyTensor[]
Jdagger = LazyTensor[]

# Interaction picture
Hrot = LazyTensor[]
for i=1:N, j=1:N
    if i==j
        continue
    end
    h = LazyTensor(basis, [i,j], [a, at], Ω[i,j])
    push!(Hrot, h)
end
Hrot = LazySum(Hrot...)

# Schroedinger picture
function f(t, rho)
    H = LazyTensor[LazyTensor(basis, i, at*a, ω[i]) for i=1:N]
    for i=1:N, j=1:N
        if i==j
            continue
        end
        h = LazyTensor(basis, [i,j], [a, at], exp(1im*(ω[i]-ω[j])*t)*Ω[i,j])
        push!(H, h)
    end
    LazySum(H...), J, Jdagger
end

tout, rho_rot_t = cumulantexpansion.master(T, rho0, Hrot, J)
tout, rho_t = cumulantexpansion.master_timedependent(T, rho0, f)

for (i, t) in enumerate(tout)
    R = prod([embed(basis, i, expm(1im*ω[i]*t*full(at*a))) for i=1:N])
    rho_rot = rho_rot_t[i]
    rho = rho_t[i]
    @test tracedistance(full(rho_rot), full(R)*full(rho)*dagger(full(R))) < 1e-5
end


text:
"""
    auto_r0(f, g, a...)

Create a K_ object and call f on it before finalizing it with r0.

  Usage:

  ```julia
  @test auto_r0(kj, 42) do x
    # work with x
  end
  ```
"""
function auto_r0(f, g, a...)
  x = g(a...)
  try
    return f(x)
  finally
    r0(x)
  end
end

NUMBER_TYPES = [UInt8, Int16, Int32, Int64, Float32, Float64]

function roundtrip_scalar(jk, kj, x)
    k = kj(x)
    r = jk(k)
    r0(k)
    return x == r
end

function empty_vector(t, typ)
  x = asarray(ktn(t, 0))
  res = (eltype(x) === typ
    && length(x) == 0
    && collect(x) == typ[])
  return res
end
# Adds r0(x) to the end and runs @test
macro xtest(e::Expr)
  quote
    @test $(esc(Expr(:block, Expr(:(=), :r, e), :(r0(x)), :r)))
  end
end

@testset "Low level (_k)" begin
  @testset "Type names" begin
    @test string(K_) == "K_"
    # @test string(K) == "K"
    @test string(K_short) == "K_short"
  end
  @testset "Reference counts" begin
    @xtest begin
      x = kj(0)
      n = [xr(x)]
      push!(n, xr(r1(x)))
      r0(x)
      push!(n, xr(x))
      n == [0, 1, 0]
    end
    @xtest begin
      x = kj(0)
      y = K_new(x)
      n = xr(y)
      r0(y)
      x === y && n == 1
    end
    @test begin
      x = K_Ref(kj(666))
      y = K_Ref(r1(x.x))
      finalize(x)
      xr(y.x) == 0
    end
  end
  @testset "Scalar constructors" begin
    @xtest (x = kb(1); xt(x) == -KB && xg(x) === G_(1))
    @xtest (x = kg(8); xt(x) == -KG && xg(x) === G_(8))
    @xtest (x = kh(100); xt(x) == -KH && xh(x) === H_(100))
    @xtest (x = ki(100); xt(x) == -KI && xi(x) === I_(100))
    @xtest (x = kj(100); xt(x) == -KJ && xj(x) === J_(100))
    @xtest (x = ke(1.5); xt(x) == -KE && xe(x) === E_(1.5))
    @xtest (x = kf(1.5); xt(x) == -KF && xf(x) === F_(1.5))
    @xtest (x = kc(10); xt(x) == -KC && xg(x) == G_(10))
    @xtest (x = ks("a"); xt(x) == -KS && xs(x) == "a")
    @test auto_r0(ktj, I_(101), I_(0)) do x
      xt(x) == 101 && xj(x) == 0
    end
  end
  @testset "Interning" begin
    @test (x = sn("abcd", 3); unsafe_string(x) == "abc")
    @test (x = sn("abcd", 2); x === ss("ab"))
  end
  @testset "Date conversions" begin
    @test ymd(2000, 1, 1) == 0
    @test ymd(1970, 1, 1) == -10957
    @test ymd(2017, 1, 1) == 6210
    @test dj(0) == 20000101
    @test dj(-10957) == 19700101
    @test dj(6210) == 20170101
  end
  @testset "Scalar roundtrip" begin
    @test roundtrip_scalar(xg, kg, I_(42))
    @test roundtrip_scalar(xh, kh, I_(101))
    @test roundtrip_scalar(xi, ki, I_(10101))
    @test roundtrip_scalar(xj, kj, J_(1010101))
    @test roundtrip_scalar(xf, kf, F_(1e10))
    @test roundtrip_scalar(xs, ks, "abc")
  end
  @testset "Vector types" begin
    @test empty_vector(KB, B_)
    @test empty_vector(KG, G_)
    @test empty_vector(KH, H_)
    @test empty_vector(KI, I_)
    @test empty_vector(KJ, J_)
    @test empty_vector(KE, E_)
    @test empty_vector(KF, F_)
  end
  @testset "Vector accessors S" begin
    auto_r0(ktn, KS, 5) do x
      @test x|>n == 5
      @test x|>t == KS
      @test (a = kS(x); a[:] = s = ss("a"); a[2:4] == [s, s, s])
    end
  end
  @testset "Vector accessors K" begin
    auto_r0(ktn, KK, 5) do x
      @test x|>n == 5
      @test x|>t == KK
      @test begin
        a = kK(x)
        a[:] = k = ktj(101, 0)
        map(r1, a[2:end])
        a[2:4] == [k, k, k]
      end
    end
  end
  @testset "Vector accessors 2" for T in "GHIJEFC"
    t_, f_ = map(eval, [Symbol("K", T), Symbol("k", T)])
    auto_r0(ktn, t_, 5) do x
      @test x|>n == 5
      @test x|>t == t_
      @test (a = f_(x); a[:] = 1:5; a[2:4] == [2, 3, 4])
    end
  end
  @testset "Mixed types list" begin
    @xtest begin
      x = knk(3, kb(I_(0)), ku(U_(10)), ku(0))
      xn(x) == 3
    end
  end
  @testset "Table and dict" begin
    @xtest begin
      a = ktn(KI, 0)
      b = ktn(KJ, 0)
      x = xD(a, b)
      xx(x) === a && xy(x) == b
    end
    @xtest begin
      a = ktn(KS, 1)
      fill!(asarray(a, false), ss("a"))
      b = knk(1, ktn(KJ, 0))
      d = xD(a, b)
      x = xT(d)
      xk(x) === d
    end
  end
  @testset "Vector ops" begin
    let x = K(ktn(KJ, 3))
      @test eltype(x) === J_
      @test length(x) == 3
      @test fill!(x, 42) == [42, 42, 42]
      @test copy!(x, [1, 2, 3]) == [1, 2, 3]
    end
  end
  @testset "Vector extend" begin
    @xtest begin
      x = ktn(KS, 0)
      x = js(Ref{K_}(x), ss("a"))
      x= js(Ref{K_}(x), ss("b"))
      map(unsafe_string, asarray(x, false)) == ["a", "b"]
    end
    @xtest begin
      x = ktn(KK, 0)
      x = jk(Ref{K_}(x), ktn(0, 0))
      x = jk(Ref{K_}(x), ktn(0, 0))
      xn(x) == 2
    end
  end
  @testset "Serializarion/deserialization" begin
    @test auto_r0(kj, 42) do x
      s = b9(0, x)
      y = d9(s)
      try
        return okx(s) == 1 && xj(y)== xj(x)
      finally
        r0(s)
        r0(y)
      end
    end
  end
  @testset "communications" begin
    server() do port
      let h = khp("", port)
        @test h > 0
        # XXX: On the server, the first call to k(h, ..) gets
        # 'rcv. OS reports: Resource temporarily unavailable.
        # TODO: Figure out a cause and write a robust "hget".
        # Q.GOT_Q && r0(ee(k(h, "666")))
        # @test 666 == auto_r0(k, h, "666") do x
        #  xj(x)
        # end
      end
    end
  end
end  # "Low level"

@testset "asarray" begin
  @test begin
    a = asarray(kj(42))
    a[] == 42
  end
  @test begin
    x = ktn(KJ, 5)
    a = asarray(x)
    a[:] = 1:5
    kJ(x) == a
  end
  @test begin
    n = 0x0102030405060708090a0b0c0d0e0f10
    x = ku(n)
    a = asarray(x)
    a[] == n
  end
  @test begin
    d = xD(ktn(KJ, 0), ktn(KS, 0))
    a = asarray(d)
    x, y = a
    asarray(x, false) == asarray(y, false) == []
  end
  @test begin
    d = xD(ktn(KS, 2), knk(2, ktn(KJ, 3), ktn(KF, 3)))
    a = asarray(d)
    k, v = a
    asarray(k, false)[:] = map(ss, ["a", "b"])
    c1 = asarray(v, false)[1]
    c2 = asarray(v, false)[2]
    asarray(c1, false)[:] = 1:3
    asarray(c2, false)[:] = 3.14
    table = xT(r1(d))
    b = asarray(table)
    b[] == d
  end
  @test (a = asarray(K_new(nothing)); a[] == 0)
end


text:
"""
    greedy(g)

Create a tour which starts at node 1 and goes to the nearest unvisited neighbor in each step.
Return the tour as a list of nodes and the cost of the tour.
If no tour is possible return `nothing, Inf`. If no tour was found but there isn't a proof that no tour exists
return `nothing, NaN`.
"""
function greedy(g)
    N = nv(g)
    sorted_neighbors = Vector{Vector{Int}}(undef, N)
    for i in 1:N
        vec_neighbors = collect(neighbors(g,i))
        NN = length(vec_neighbors)
        weights = zeros(NN)
        for j in 1:NN
            weights[j] = get_prop(g, i, vec_neighbors[j], :weight)
        end
        order = sortperm(weights)
        sorted_neighbors[i] = vec_neighbors[order]
        p = findfirst(>(0),weights[order])
        # if three edges are fixed
        if p > 3
            return nothing, Inf
        end
        # if both end points are fixed we can fix them here
        if p == 3
            sorted_neighbors[i] = sorted_neighbors[i][1:2]
        end
    end

    visited = zeros(Bool, N)
    node_next = ones(Int, N)
    tour = zeros(Int, N)
    tour[1] = 1
    visited[1] = true
    cp = 1
    backtrack_c = 0
    while 1 <= cp < N
        # simply give up searching :D
        if backtrack_c == 1000000
            return nothing, NaN
        end
        neighbor = next_neighbor!(g, tour, cp, node_next, visited, sorted_neighbors[tour[cp]])
        if neighbor === nothing
            backtrack_c += 1
            node_next[tour[cp]] = 1
            visited[tour[cp]] = false
            cp -= 1
            if cp == 0
                break
            end
        else
            visited[neighbor] = true
            cp += 1
            tour[cp] = neighbor
        end
    end
    return tour, get_tour_cost(g, tour)
end

function next_neighbor!(g, tour, cp, node_next, visited, sorted_neighbors)
    for neighbor_id in node_next[tour[cp]]:length(sorted_neighbors)
        node_next[tour[cp]] += 1
        neighbor = sorted_neighbors[neighbor_id]
        visited[neighbor] && continue
        if (cp == length(tour)-1 && !has_edge(g, neighbor, 1))
            continue
        end

        return sorted_neighbors[neighbor_id]
    end
    return nothing
end

function two_opt(tour, tour_cost, cost_matrix)
    swapped = true
    while swapped
        swapped = false
        for i in 1:length(tour)
            for j in i+1:length(tour)
                swap_cost = get_swap_cost(tour, i, j, cost_matrix)
                if swap_cost < -1e-6
                    tour_cost += swap_cost
                    swap!(tour, i, j)
                    swapped = true
                    break
                end
            end
            swapped && break
        end
    end
    return tour, tour_cost
end

function get_swap_cost(tour, i, j, cost_matrix)
    cost = -cost_matrix[tour[i], tour[mod1(i+1,end)]]
    cost -= cost_matrix[tour[j], tour[mod1(j+1,end)]]

    cost += cost_matrix[tour[i], tour[j]]
    cost += cost_matrix[tour[mod1(i+1,end)], tour[mod1(j+1,end)]]
    return cost
end

function swap!(tour, i, j)
    reverse!(tour, i+1, j)
end


text:
module Filtering

using IPFitting: LsqDB
export filter_basis, filter_configs,
       anyf, allf

function filter_basis(db::LsqDB, args...)
   @assert length(args) > 0
   choose_b = fill(false, length(db.basis))
   for (i, b) in enumerate(db.basis)
      choose_b[i] = all( a(b) for a in args )
   end
   return findall(choose_b)
end

function filter_configs(db::LsqDB, args...)
   @assert length(args) > 0
   choose_cfg = fill(false, length(db.configs))
   for (i, cfg) in enumerate(db.configs)
      choose_cfg[i] = all( a(cfg) for a in args )
   end
   return findall(choose_cfg)
end

anyf(args...) = x -> any( f(x) for f in args )
allf(args...) = x -> all( f(x) for f in args )

end

# TODO: create a filter that picks a random subset of the configurations,
# this used to be controlled by the following:
# """
# `struct ConfigFitInfo` : stores `configtype`-dependent
# weights and also proportion of configurations to fit to
# """
# struct ConfigFitInfo
#    weight::Float64
#    proportion::Float64
#    order::Symbol
# end
#
# ConfigFitInfo(cfi::ConfigFitInfo) = deepcopy(cfi)
# ConfigFitInfo(cfi::Tuple) = ConfigFitInfo(cfi...)
# ConfigFitInfo(w::Real) = ConfigFitInfo(w, 1.0, :ignore)
# ConfigFitInfo(w::Real, p::Real) = ConfigFitInfo(w, p, :ignore)
# ConfigFitInfo(w::Real, p::Real, o::Any) = ConfigFitInfo(w, p, Symbol(o))


text:
#   This file is part of VerTeX.jl. It is licensed under the MIT license
#   Copyright (C) 2019 Michael Reed

#---------#
# Options #
#---------#

struct OptionSpec
    name::String
    short_name::Union{Nothing,String}
    api::Pair{Symbol, Any}
    takes_arg::Bool
end

# TODO assert names matching lex regex
# assert now so that you don't fail at user time
# see function `REPLMode.APIOptions`
function OptionSpec(;name::String,
                    short_name::Union{Nothing,String}=nothing,
                    takes_arg::Bool=false,
                    api::Pair{Symbol,<:Any})::OptionSpec
    takes_arg && @assert hasmethod(api.second, Tuple{String})
    return OptionSpec(name, short_name, api, takes_arg)
end

function OptionSpecs(decs::Vector{OptionDeclaration})::Dict{String, OptionSpec}
    specs = Dict()
    for x in decs
        opt_spec = OptionSpec(;x...)
        @assert get(specs, opt_spec.name, nothing) === nothing # don't overwrite
        specs[opt_spec.name] = opt_spec
        if opt_spec.short_name !== nothing
            @assert get(specs, opt_spec.short_name, nothing) === nothing # don't overwrite
            specs[opt_spec.short_name] = opt_spec
        end
    end
    return specs
end

#-----------#
# Arguments #
#-----------#
struct ArgSpec
    count::Pair
    parser::Function
end

struct CommandSpec
    kind::CommandKind
    canonical_name::String
    short_name::Union{Nothing,String}
    handler::Union{Nothing,Function}
    argument_spec::ArgSpec
    option_specs::Dict{String,OptionSpec}
    completions::Union{Nothing,Function}
    description::String
    help::Union{Nothing,Markdown.MD}
end

function CommandSpec(;kind::Union{Nothing,CommandKind}      = nothing,
                     name::Union{Nothing,String}            = nothing,
                     short_name::Union{Nothing,String}      = nothing,
                     handler::Union{Nothing,Function}       = nothing,
                     option_spec::Vector{OptionDeclaration} = OptionDeclaration[],
                     help::Union{Nothing,Markdown.MD}       = nothing,
                     description::Union{Nothing,String}     = nothing,
                     completions::Union{Nothing,Function}   = nothing,
                     arg_count::Pair                        = (0=>0),
                     arg_parser::Function                   = unwrap,
                     )::CommandSpec
    @assert kind !== nothing "Register and specify a `CommandKind`"
    @assert name !== nothing "Supply a canonical name"
    @assert description !== nothing "Supply a description"
    # TODO assert isapplicable completions dict, string
    return CommandSpec(kind, name, short_name, handler, ArgSpec(arg_count, arg_parser),
                       OptionSpecs(option_spec), completions, description, help)
end

function CommandSpecs(declarations::Vector{CommandDeclaration})::Dict{String,CommandSpec}
    specs = Dict()
    for dec in declarations
        spec = CommandSpec(;dec...)
        @assert !haskey(specs, spec.canonical_name) "duplicate spec entry"
        specs[spec.canonical_name] = spec
        if spec.short_name !== nothing
            @assert !haskey(specs, spec.short_name) "duplicate spec entry"
            specs[spec.short_name] = spec
        end
    end
    return specs
end

function SuperSpecs(compound_commands)::Dict{String,Dict{String,CommandSpec}}
    super_specs = Dict()
    for x in compound_commands
        name = x.first
        spec = CommandSpecs(x.second)
        @assert !haskey(super_specs, name) "duplicate super spec entry"
        super_specs[name] = spec
    end
    return super_specs
end

###########
# Parsing #
###########

# QString: helper struct for retaining quote information
struct QString
    raw::String
    isquoted::Bool
end
unwrap(xs::Vector{QString}) = map(x -> x.raw, xs)

#---------#
# Options #
#---------#
struct Option
    val::String
    argument::Union{Nothing,String}
    Option(val::AbstractString) = new(val, nothing)
    Option(val::AbstractString, arg::Union{Nothing,String}) = new(val, arg)
end
Base.show(io::IO, opt::Option) = print(io, "--$(opt.val)", opt.argument == nothing ? "" : "=$(opt.argument)")
wrap_option(option::String)  = length(option) == 1 ? "-$option" : "--$option"
is_opt(word::AbstractString) = first(word) == '-' && word != "-"

function parse_option(word::AbstractString)::Option
    m = match(r"^(?: -([a-z]) | --([a-z]{2,})(?:\s*=\s*(\S*))? )$"ix, word)
    m === nothing && vtxerror("malformed option: ", repr(word))
    option_name = m.captures[1] !== nothing ? m.captures[1] : m.captures[2]
    option_arg  = m.captures[3] === nothing ? nothing : String(m.captures[3])
    return Option(option_name, option_arg)
end

#-----------#
# Statement #
#-----------#
# Statement: text-based representation of a command
Base.@kwdef mutable struct Statement
    super::Union{Nothing,String}                  = nothing
    spec::Union{Nothing,CommandSpec}              = nothing
    options::Union{Vector{Option},Vector{String}} = String[]
    arguments::Vector{QString}                    = QString[]
    preview::Bool                                 = false
end

function lex(cmd::String)::Vector{QString}
    in_doublequote = false
    in_singlequote = false
    qstrings = QString[]
    token_in_progress = Char[]

    push_token!(is_quoted) = begin
        push!(qstrings, QString(String(token_in_progress), is_quoted))
        empty!(token_in_progress)
    end

    for c in cmd
        if c == '"'
            if in_singlequote # raw char
                push!(token_in_progress, c)
            else # delimiter
                in_doublequote ? push_token!(true) : push_token!(false)
                in_doublequote = !in_doublequote
            end
        elseif c == '\''
            if in_doublequote # raw char
                push!(token_in_progress, c)
            else # delimiter
                in_singlequote ? push_token!(true) : push_token!(false)
                in_singlequote = !in_singlequote
            end
        elseif c == ' '
            if in_doublequote || in_singlequote # raw char
                push!(token_in_progress, c)
            else # delimiter
                push_token!(false)
            end
        elseif c == ';'
            if in_doublequote || in_singlequote # raw char
                push!(token_in_progress, c)
            else # special delimiter
                push_token!(false)
                push!(qstrings, QString(";", false))
            end
        else
            push!(token_in_progress, c)
        end
    end
    (in_doublequote || in_singlequote) ? vtxerror("unterminated quote") : push_token!(false)
    # to avoid complexity in the main loop, empty tokens are allowed above and
    # filtered out before returning
    return filter(x->!isempty(x.raw), qstrings)
end

function tokenize(cmd::String)
    cmd = replace(replace(cmd, "\r\n" => "; "), "\n" => "; ") # for multiline commands
    qstrings = lex(cmd)
    statements = foldl(qstrings; init=[QString[]]) do collection, next
        (next.raw == ";" && !next.isquoted) ?
            push!(collection, QString[]) :
            push!(collection[end], next)
        return collection
    end
    return statements
end

function core_parse(words::Vector{QString}; only_cmd=false)
    statement = Statement()
    word = nothing
    function next_word!()
        isempty(words) && return false
        word = popfirst!(words)
        return true
    end

    # begin parsing
    next_word!() || return statement, ((word === nothing) ? nothing : word.raw)
    if word.raw == "preview"
        statement.preview = true
        next_word!() || return statement, word.raw
    end
    # handle `?` alias for help
    # It is special in that it requires no space between command and args
    if word.raw[1]=='?' && !word.isquoted
        length(word.raw) > 1 && pushfirst!(words, QString(word.raw[2:end],false))
        word = QString("?", false)
    end
    # determine command
    super = get(super_specs, word.raw, nothing)
    if super !== nothing # explicit
        statement.super = word.raw
        next_word!() || return statement, word.raw
        command = get(super, word.raw, nothing)
        command !== nothing || return statement, word.raw
    else # try implicit package
        super = super_specs["package"]
        command = get(super, word.raw, nothing)
        command !== nothing || return statement, word.raw
    end
    statement.spec = command

    only_cmd && return statement, word.raw # hack to hook in `help` command

    next_word!() || return statement, word.raw

    # full option parsing is delayed so that the completions parser can use the raw string
    while is_opt(word.raw)
        push!(statement.options, word.raw)
        next_word!() || return statement, word.raw
    end

    pushfirst!(words, word)
    statement.arguments = words
    return statement, words[end].raw
end

parse(input::String) =
    map(Base.Iterators.filter(!isempty, tokenize(input))) do words
        statement, _ = core_parse(words)
        statement.spec === nothing && vtxerror("Could not determine command")
        statement.options = map(parse_option, statement.options)
        statement
    end

#------------#
# APIOptions #
#------------#
const APIOptions = Dict{Symbol, Any}
function APIOptions(options::Vector{Option},
                    specs::Dict{String, OptionSpec},
                    )::APIOptions
    api_options = Dict{Symbol, Any}()
    for option in options
        spec = specs[option.val]
        api_options[spec.api.first] = spec.takes_arg ?
            spec.api.second(option.argument) :
            spec.api.second
    end
    return api_options
end
Context!(ctx::APIOptions)::Context = Types.Context!(collect(ctx))

#---------#
# Command #
#---------#
Base.@kwdef struct Command
    spec::Union{Nothing,CommandSpec} = nothing
    options::APIOptions              = APIOptions()
    arguments::Vector                = []
    preview::Bool                    = false
end

function enforce_option(option::Option, specs::Dict{String,OptionSpec})
    spec = get(specs, option.val, nothing)
    spec !== nothing || vtxerror("option '$(option.val)' is not a valid option")
    if spec.takes_arg
        option.argument !== nothing ||
            vtxerror("option '$(option.val)' expects an argument, but no argument given")
    else # option is a switch
        option.argument === nothing ||
            vtxerror("option '$(option.val)' does not take an argument, but '$(option.argument)' given")
    end
end

"""
checks:
- options are understood by the given command
- options do not conflict (e.g. `rm --project --manifest`)
- options which take an argument are given arguments
- options which do not take arguments are not given arguments
"""
function enforce_option(options::Vector{Option}, specs::Dict{String,OptionSpec})
    unique_keys = Symbol[]
    get_key(opt::Option) = specs[opt.val].api.first

    # per option checking
    foreach(x->enforce_option(x,specs), options)
    # checking for compatible options
    for opt in options
        key = get_key(opt)
        if key in unique_keys
            conflicting = filter(opt->get_key(opt) == key, options)
            vtxerror("Conflicting options: $conflicting")
        else
            push!(unique_keys, key)
        end
    end
end

"""
Final parsing (and checking) step.
This step is distinct from `parse` in that it relies on the command specifications.
"""
function Command(statement::Statement)::Command
    # arguments
    arg_spec = statement.spec.argument_spec
    arguments = arg_spec.parser(statement.arguments)
    if !(arg_spec.count.first <= length(arguments) <= arg_spec.count.second)
        vtxerror("Wrong number of arguments")
    end
    # options
    opt_spec = statement.spec.option_specs
    enforce_option(statement.options, opt_spec)
    options = APIOptions(statement.options, opt_spec)
    return Command(statement.spec, options, arguments, statement.preview)
end

######################
# REPL mode creation #
######################

# Provide a string macro pkg"cmd" that can be used in the same way
# as the REPLMode `pkg> cmd`. Useful for testing and in environments
# where we do not have a REPL, e.g. IJulia.
struct MiniREPL <: REPL.AbstractREPL
    display::TextDisplay
    t::REPL.Terminals.TTYTerminal
end
function MiniREPL()
    MiniREPL(TextDisplay(stdout), REPL.Terminals.TTYTerminal(get(ENV, "TERM", Sys.iswindows() ? "" : "dumb"), stdin, stdout, stderr))
end
REPL.REPLDisplay(repl::MiniREPL) = repl.display


text:
"""
    TexasHoldem

A no-limit Texas Holdem simulator.

# Terminology
 - `game` a single "game", where players are dealt hands,
   winner(s) are declared once.
 - `state` a point or process in the game, including
   `PreFlop`, `Flop`, `Turn`, `River`.
 - `round` the process of each player deciding which
   actions to take, until no further actions are taking.
"""
module TexasHoldem

using PlayingCards
using PokerHandEvaluator
using Printf

export AbstractGameState, PreFlop, Flop, Turn, River

abstract type AbstractGameState end
struct PreFlop <: AbstractGameState end
struct Flop <: AbstractGameState end
struct Turn <: AbstractGameState end
struct River <: AbstractGameState end

include("player_types.jl")
include("transactions.jl")
include("table.jl")
include("game.jl")
include("player_actions.jl")
include("player_options.jl")
include("config_game.jl")

end # module


text:
using Plots
# pyplot()
do_viz = true

@views function heat_1D_2procs()
    # Physics
    Tl  = 10.0   # left Temperature
    Tr  = 1.0    # right Temperature
    λ   = 1.0
    ρCp = 1.0
    nt  = 200    # number of time steps
    # Numerics
    nx  = 32     # number of local grid points
    dx  = 1.0    # cell size
    # Derived numerics
    dt  = dx^2/ρCp/λ/2.1
    # Initial condition
    TL  = Tl*ones(nx)
    TR  = Tr*ones(nx)
    T   = [TL[1:end-1]; TR[2:end]]
    Tg  = T
    # Time loop
    for it = 1:nt
        # Compute physics
        # TODO: Implement ∆TL/dt = λ/ρCp ∂^2(TL)/dx^2
        # TODO: Implement ∆TR/dt = λ/ρCp ∂^2(TR)/dx^2

        # Update boundaries (MPI)
        # TODO: Implement the boundary exchange - the core of distributed memory parallelisation

        # Global picture
        T .= [TL[1:end-1]; TR[2:end]]
        # Compute physics locally
        Tg[2:end-1] .= # update Tg (∆Tg/dt = λ/ρCp ∂^2(Tg)/dx^2) to cross check your fake parallelisation implementation
        # Visualise
        if do_viz
            plot(Tg, legend=false, linewidth=0, markershape=:circle, markersize=5)
            display(plot!(T, legend=false, linewidth=2, framestyle=:box, xlabel="lx", ylabel="heat", title="diffusion it=$(it)"))
        end
    end
    return
end

@time heat_1D_2procs()

max_stars_repo_name: luraess/geo-hpc-course

max_issues_repo_name: luraess/geo-hpc-course

max_forks_repo_name: luraess/geo-hpc-course


text:
"""
Default implementation of mass matrix entries. Keeps the default values in the identity matrix
"""
mass_matrix_entries(::PSY.DynamicInjection) = nothing


text:
# Mahalanobis distances

struct Mahalanobis{M<:AbstractMatrix} <: Metric
    qmat::M
end

struct SqMahalanobis{M<:AbstractMatrix} <: SemiMetric
    qmat::M
end

function result_type(d::Mahalanobis, ::Type{T1}, ::Type{T2}) where {T1,T2}
    z = zero(T1) - zero(T2)
    return typeof(sqrt(z * zero(eltype(d.qmat)) * z))
end

function result_type(d::SqMahalanobis, ::Type{T1}, ::Type{T2}) where {T1,T2}
    z = zero(T1) - zero(T2)
    return typeof(z * zero(eltype(d.qmat)) * z)
end

# SqMahalanobis

function (dist::SqMahalanobis)(a::AbstractVector, b::AbstractVector)
    if length(a) != length(b)
        throw(DimensionMismatch("first array has length $(length(a)) which does not match the length of the second, $(length(b))."))
    end

    Q = dist.qmat
    z = a - b
    return dot(z, Q * z)
end

sqmahalanobis(a::AbstractVector, b::AbstractVector, Q::AbstractMatrix) = SqMahalanobis(Q)(a, b)

function colwise!(r::AbstractArray, dist::SqMahalanobis, a::AbstractMatrix, b::AbstractMatrix)
    Q = dist.qmat
    get_colwise_dims(size(Q, 1), r, a, b)
    z = a - b
    dot_percol!(r, Q * z, z)
end

function colwise!(r::AbstractArray, dist::SqMahalanobis, a::AbstractVector, b::AbstractMatrix)
    Q = dist.qmat
    get_colwise_dims(size(Q, 1), r, a, b)
    z = a .- b
    dot_percol!(r, Q * z, z)
end

function _pairwise!(r::AbstractMatrix, dist::SqMahalanobis,
                    a::AbstractMatrix, b::AbstractMatrix)
    Q = dist.qmat
    m, na, nb = get_pairwise_dims(size(Q, 1), r, a, b)

    Qa = Q * a
    Qb = Q * b
    sa2 = dot_percol(a, Qa)
    sb2 = dot_percol(b, Qb)
    mul!(r, a', Qb)

    for j = 1:nb
        @simd for i = 1:na
            @inbounds r[i, j] = max(sa2[i] + sb2[j] - 2 * r[i, j], 0)
        end
    end
    r
end

function _pairwise!(r::AbstractMatrix, dist::SqMahalanobis,
                    a::AbstractMatrix)
    Q = dist.qmat
    m, n = get_pairwise_dims(size(Q, 1), r, a)

    Qa = Q * a
    sa2 = dot_percol(a, Qa)
    mul!(r, a', Qa)

    for j = 1:n
        for i = 1:(j - 1)
            @inbounds r[i, j] = r[j, i]
        end
        r[j, j] = 0
        for i = (j + 1):n
            @inbounds r[i, j] = max(sa2[i] + sa2[j] - 2 * r[i, j], 0)
        end
    end
    r
end

# Mahalanobis

function (dist::Mahalanobis)(a::AbstractVector, b::AbstractVector)
    sqrt(SqMahalanobis(dist.qmat)(a, b))
end

mahalanobis(a::AbstractVector, b::AbstractVector, Q::AbstractMatrix) = Mahalanobis(Q)(a, b)

function colwise!(r::AbstractArray, dist::Mahalanobis, a::AbstractMatrix, b::AbstractMatrix)
    sqrt!(colwise!(r, SqMahalanobis(dist.qmat), a, b))
end

function colwise!(r::AbstractArray, dist::Mahalanobis, a::AbstractVector, b::AbstractMatrix)
    sqrt!(colwise!(r, SqMahalanobis(dist.qmat), a, b))
end

function _pairwise!(r::AbstractMatrix, dist::Mahalanobis,
                    a::AbstractMatrix, b::AbstractMatrix)
    sqrt!(_pairwise!(r, SqMahalanobis(dist.qmat), a, b))
end

function _pairwise!(r::AbstractMatrix, dist::Mahalanobis,
                    a::AbstractMatrix)
    sqrt!(_pairwise!(r, SqMahalanobis(dist.qmat), a))
end


text:
function initialization_performed(obj::IterationListener, arg0::IterationEvent)
    return jcall(obj, "initializationPerformed", void, (IterationEvent,), arg0)
end

function iteration_performed(obj::IterationListener, arg0::IterationEvent)
    return jcall(obj, "iterationPerformed", void, (IterationEvent,), arg0)
end

function iteration_started(obj::IterationListener, arg0::IterationEvent)
    return jcall(obj, "iterationStarted", void, (IterationEvent,), arg0)
end

function termination_performed(obj::IterationListener, arg0::IterationEvent)
    return jcall(obj, "terminationPerformed", void, (IterationEvent,), arg0)
end


text:
using REPL
using Base.Docs: catdoc, modules, DocStr, Binding, MultiDoc, isfield, namify, bindingexpr,
    defined, resolve, getdoc, meta, aliasof, signature, isexpr
import Base.Docs: doc, formatdoc, parsedoc, apropos,
    DocStr, lazy_iterpolate, quot, metadata, docexpr, keyworddoc

function REPL.lookup_doc(ex)
    locale = get(ENV, "REPL_LOCALE", "en_US")

    if locale == ""
        locale == "en_US"
    end

    kw_dict = locale == "zh_CN" ? keywords : Docs.keywords

    if haskey(kw_dict, ex)
        Docs.parsedoc(kw_dict[ex])
    elseif isa(ex, Union{Expr, Symbol})
        binding = esc(bindingexpr(namify(ex)))
        if isexpr(ex, :call) || isexpr(ex, :macrocall)
            sig = esc(signature(ex))
            :($(doc)($binding, $sig))
        else
            :($(doc)($binding))
        end
    else
        :($(doc)($(typeof)($(esc(ex)))))
    end
end

function Base.Docs.keyworddoc(__source__, __module__, str, def::Base.BaseDocs.Keyword)
    @nospecialize str
    docstr = esc(docexpr(__source__, __module__, lazy_iterpolate(str), metadata(__source__, __module__, def, false)))
    return :($setindex!($(keywords), $docstr, $(esc(quot(def.name)))); nothing)
end


text:
using Distributed, AzManagers, AzSessions

#
# create a session-bundle that will be sent to the detached job.
#
sessionbundle!(
    management = AzSession())

#
# create a detached job, and persist the server
#
job1 = @detach vm(;session=sessionbundle(:management),persist=true) begin
    write(stdout, "some stdout stuff, 1\n")
    write(stderr, "some stderr stuff, 1\n")
    @info "write a message using Logging, 1"
end

status(job1)
read(job1)
read(job1;stdio=stderr)

#
# send a new job to the server started above
#
job2 = @detachat job1.vm["ip"] begin
    write(stdout, "some stdout stuff, 2\n")
    write(stderr, "some stderr stuff, 2\n")
    @info "write a message using Logging, 2"
end

#
# get job status
#
status(job1)
status(job2)

#
# wait for jobs to finish
#
wait(job1)
wait(job2)

#
# get stdout from each job
#
write(stdout, read(job1))
write(stdout, read(job2))

#
# shut-down the detached server
#
rmproc(job1.vm; session=sessionbundle(:management))

#
# create a new job on a new detached server that auto-destructs upon completion of its work
#
job3 = @detach vm(;session=sessionbundle(:management),persist=false) begin
    @info "running detached stuff, 3"
end

#
# more involved example with storage and scale-sets
#
sessionbundle!(
    management = AzSession(),
    storage = AzSession(;scope="openid+offline_access+https://storage.azure.com/user_impersonation"))

job4 = @detach vm(;session=sessionbundle(:management),persist=true) begin
    using Distributed, AzManagers, AzStorage, AzSessions
    addprocs("cbox04", 4; session=sessionbundle(:management))

    for pid in workers()
        hostname = remotecall_fetch(gethostname, pid)
        @info "pid=$pid, hostname=$hostname"
    end

    using Random
    x = randstring('a':'z',4)
    container = AzContainer("foobar-$x"; storageaccount="samkaplan0", session=sessionbundle(:storage))
    mkpath(container)
    write(container, "foobar-$x", "hello")
    @info "contents=$(read(container, "foobar-$x", String))"
    rm(container)

    rmprocs(workers())
end

status(job4)
wait(job4)
read(job4)

rmproc(job4.vm; session=sessionbundle(:management))


text:
# Copyright (c) 2017: Miles Lubin and contributors
# Copyright (c) 2017: Google Inc.
#
# Use of this source code is governed by an MIT-style license that can be found
# in the LICENSE.md file or at https://opensource.org/licenses/MIT.

module TestConstraintSquare

using Test

using MathOptInterface
const MOI = MathOptInterface

function runtests()
    for name in names(@__MODULE__; all = true)
        if startswith("$(name)", "test_")
            @testset "$(name)" begin
                getfield(@__MODULE__, name)()
            end
        end
    end
    return
end

include("../utilities.jl")

function test_Square()
    mock = MOI.Utilities.MockOptimizer(MOI.Utilities.Model{Float64}())
    config = MOI.Test.Config()
    bridged_mock = MOI.Bridges.Constraint.Square{Float64}(mock)
    MOI.Test.runtests(
        bridged_mock,
        config,
        include = [
            "test_basic_VectorOfVariables_PositiveSemidefiniteConeSquare",
            "test_basic_VectorAffineFunction_PositiveSemidefiniteConeSquare",
            "test_basic_VectorQuadraticFunction_PositiveSemidefiniteConeSquare",
        ],
    )
    MOI.empty!(bridged_mock)
    mock.optimize! =
        (mock::MOI.Utilities.MockOptimizer) -> MOI.Utilities.mock_optimize!(
            mock,
            ones(4),
            (MOI.ScalarAffineFunction{Float64}, MOI.EqualTo{Float64}) =>
                [2, 2],
        )
    MOI.Test.test_conic_PositiveSemidefiniteConeSquare_VectorOfVariables(
        bridged_mock,
        config,
    )
    MOI.empty!(bridged_mock)
    mock.optimize! =
        (mock::MOI.Utilities.MockOptimizer) -> MOI.Utilities.mock_optimize!(
            mock,
            ones(4),
            (
                MOI.VectorAffineFunction{Float64},
                MOI.PositiveSemidefiniteConeTriangle,
            ) => [[1, -1, 1]],
            (MOI.ScalarAffineFunction{Float64}, MOI.EqualTo{Float64}) =>
                [2, 2],
        )
    MOI.Test.test_conic_PositiveSemidefiniteConeSquare_VectorAffineFunction(
        bridged_mock,
        config,
    )
    ci = first(
        MOI.get(
            bridged_mock,
            MOI.ListOfConstraintIndices{
                MOI.VectorAffineFunction{Float64},
                MOI.PositiveSemidefiniteConeSquare,
            }(),
        ),
    )
    _test_delete_bridge(
        bridged_mock,
        ci,
        4,
        (
            (
                MOI.VectorAffineFunction{Float64},
                MOI.PositiveSemidefiniteConeTriangle,
                0,
            ),
            (MOI.ScalarAffineFunction{Float64}, MOI.EqualTo{Float64}, 1),
        ),
    )
    return
end

end  # module

TestConstraintSquare.runtests()

max_stars_repo_name: JuliaOpt/MathOptInterface

max_issues_repo_name: JuliaOpt/MathOptInterface

max_forks_repo_name: JuliaOpt/MathOptInterface


text:
AbstractCalculated = Union{Real, Void}
AbstractSymbol = Union{Real, SymEngine.Basic}

AbstractKey = Union{Symbol, Void}

abstract type AbstractEquationSet end
abstract type AbstractEquation end
abstract type AbstractSampling end
abstract type AbstractReactor end
abstract type AbstractStudy end
abstract type AbstractScan end

AbstractNullableReactor = Union{Void, AbstractReactor}