text:
# This file is based on LinearAlgebra/test/lu.jl, a part of Julia. License is MIT: https://julialang.org/license

using Test, LinearAlgebra, Random, MatrixFactorizations
using LinearAlgebra: ldiv!, BlasReal, BlasInt, BlasFloat, rdiv!

@testset "UL" begin
    n = 10

    # Split n into 2 parts for tests needing two matrices
    n1 = div(n, 2)
    n2 = 2*n1

    Random.seed!(1234321)

    areal = randn(n,n)/2
    aimg  = randn(n,n)/2
    bimg  = randn(n,2)/2
    breal = randn(n,2)/2
    creal = randn(n)/2
    cimg  = randn(n)/2
    dureal = randn(n-1)/2
    duimg  = randn(n-1)/2
    dlreal = randn(n-1)/2
    dlimg  = randn(n-1)/2
    dreal = randn(n)/2
    dimg  = randn(n)/2

    @testset "UL simple test" begin
        A = randn(n,n)
        U,L = ul(A,Val(false))
        @test U*L ≈ A
        U,L,p = ul(A)
        @test (U*L) ≈ A[p,:]
        P = ul(A).P
        @test P'*U*L ≈ A
        L,U,p = lu(A)
        @test (L*U) ≈ A[p,:] # for context
        P = lu(A).P
        @test P'*L*U ≈ A
    end

    @testset for eltya in (Float32, Float64, ComplexF32, ComplexF64, BigFloat, Int)
        a = eltya == Int ? rand(1:7, n, n) :
            convert(Matrix{eltya}, eltya <: Complex ? complex.(areal, aimg) : areal)
        d = if eltya == Int
            Tridiagonal(rand(1:7, n-1), rand(1:7, n), rand(1:7, n-1))
        elseif eltya <: Complex
            convert(Tridiagonal{eltya}, Tridiagonal(
                complex.(dlreal, dlimg), complex.(dreal, dimg), complex.(dureal, duimg)))
        else
            convert(Tridiagonal{eltya}, Tridiagonal(dlreal, dreal, dureal))
        end
        ε = εa = eps(abs(float(one(eltya))))

        if eltya <: BlasFloat
            @testset "UL factorization for Number" begin
                num = rand(eltya)
                @test (ul(num)...,) == (hcat(one(eltya)), hcat(num), [1])
                @test convert(Array, ul(num)) ≈ eltya[num]
            end
        end
        κ  = cond(a,1)
        κd    = cond(Array(d),1)
        # @testset "Tridiagonal UL" begin
        #     lud   = ul(d)
        #     @test LinearAlgebra.issuccess(lud)
        #     @test ul(lud) == lud
        #     @test_throws ErrorException lud.Z
        #     @test lud.L*lud.U ≈ lud.P*Array(d)
        #     @test lud.L*lud.U ≈ Array(d)[lud.p,:]
        #     @test AbstractArray(lud) ≈ d
        #     @test Array(lud) ≈ d
        #     if eltya != Int
        #         dlu = convert.(eltya, [1, 1])
        #         dia = convert.(eltya, [-2, -2, -2])
        #         tri = Tridiagonal(dlu, dia, dlu)
        #         @test_throws ArgumentError ul!(tri)
        #     end
        # end
        @testset for eltyb in (Float32, Float64, ComplexF32, ComplexF64, Int)
            b  = eltyb == Int ? rand(1:5, n, 2) :
                convert(Matrix{eltyb}, eltyb <: Complex ? complex.(breal, bimg) : breal)
            c  = eltyb == Int ? rand(1:5, n) :
                convert(Vector{eltyb}, eltyb <: Complex ? complex.(creal, cimg) : creal)
            εb = eps(abs(float(one(eltyb))))
            ε  = max(εa,εb)
            @testset "(Automatic) Square UL decomposition" begin
                lua   = factorize(a)
                let Bs = copy(b), Cs = copy(c)
                    for (bb, cc) in ((Bs, Cs), (view(Bs, 1:n, 1), view(Cs, 1:n)))
                        @test norm(a*(lua\bb) - bb, 1) < ε*κ*n*2 # Two because the right hand side has two columns
                        @test norm(a'*(lua'\bb) - bb, 1) < ε*κ*n*2 # Two because the right hand side has two columns
                        @test norm(a'*(lua'\a') - a', 1) < ε*κ*n^2
                        @test norm(a*(lua\cc) - cc, 1) < ε*κ*n # cc is a vector
                        @test norm(a'*(lua'\cc) - cc, 1) < ε*κ*n # cc is a vector
                        @test AbstractArray(lua) ≈ a
                        @test norm(transpose(a)*(transpose(lua)\bb) - bb,1) < ε*κ*n*2 # Two because the right hand side has two columns
                        @test norm(transpose(a)*(transpose(lua)\cc) - cc,1) < ε*κ*n
                    end

                    # Test whether Ax_ldiv_B!(y, UL, x) indeed overwrites y
                    resultT = typeof(oneunit(eltyb) / oneunit(eltya))

                    b_dest = similar(b, resultT)
                    c_dest = similar(c, resultT)

                    ldiv!(b_dest, lua, b)
                    ldiv!(c_dest, lua, c)
                    @test norm(b_dest - lua \ b, 1) < ε*κ*2n
                    @test norm(c_dest - lua \ c, 1) < ε*κ*n

                    ldiv!(b_dest, transpose(lua), b)
                    ldiv!(c_dest, transpose(lua), c)
                    @test norm(b_dest - transpose(lua) \ b, 1) < ε*κ*2n
                    @test norm(c_dest - transpose(lua) \ c, 1) < ε*κ*n

                    ldiv!(b_dest, adjoint(lua), b)
                    ldiv!(c_dest, adjoint(lua), c)
                    @test norm(b_dest - lua' \ b, 1) < ε*κ*2n
                    @test norm(c_dest - lua' \ c, 1) < ε*κ*n

                    if eltyb != Int && !(eltya <: Complex) || eltya <: Complex && eltyb <: Complex
                        p = Matrix(b')
                        q = Matrix(c')
                        p_dest = copy(p)
                        q_dest = copy(q)
                        rdiv!(p_dest, lua)
                        rdiv!(q_dest, lua)
                        @test norm(p_dest - p / lua, 1) < ε*κ*2n
                        @test norm(q_dest - q / lua, 1) < ε*κ*n
                    end
                end
                if eltya <: BlasFloat && eltyb <: BlasFloat
                    e = rand(eltyb,n,n)
                    @test norm(e/lua - e/a,1) < ε*κ*n^2
                end
            end
            # @testset "Tridiagonal UL" begin
            #     lud   = factorize(d)
            #     f = zeros(eltyb, n+1)
            #     @test_throws DimensionMismatch lud\f
            #     @test_throws DimensionMismatch transpose(lud)\f
            #     @test_throws DimensionMismatch lud'\f
            #     @test_throws DimensionMismatch LinearAlgebra.ldiv!(transpose(lud), f)
            #     let Bs = copy(b)
            #         for bb in (Bs, view(Bs, 1:n, 1))
            #             @test norm(d*(lud\bb) - bb, 1) < ε*κd*n*2 # Two because the right hand side has two columns
            #             if eltya <: Real
            #                 @test norm((transpose(lud)\bb) - Array(transpose(d))\bb, 1) < ε*κd*n*2 # Two because the right hand side has two columns
            #                 if eltya != Int && eltyb != Int
            #                     @test norm(LinearAlgebra.ldiv!(transpose(lud), copy(bb)) - Array(transpose(d))\bb, 1) < ε*κd*n*2
            #                 end
            #             end
            #             if eltya <: Complex
            #                 @test norm((lud'\bb) - Array(d')\bb, 1) < ε*κd*n*2 # Two because the right hand side has two columns
            #             end
            #         end
            #     end
            #     if eltya <: BlasFloat && eltyb <: BlasFloat
            #         e = rand(eltyb,n,n)
            #         @test norm(e/lud - e/d,1) < ε*κ*n^2
            #         @test norm((transpose(lud)\e') - Array(transpose(d))\e',1) < ε*κd*n^2
            #         #test singular
            #         du = rand(eltya,n-1)
            #         dl = rand(eltya,n-1)
            #         dd = rand(eltya,n)
            #         dd[1] = zero(eltya)
            #         du[1] = zero(eltya)
            #         dl[1] = zero(eltya)
            #         zT = Tridiagonal(dl,dd,du)
            #         @test !LinearAlgebra.issuccess(ul(zT; check = false))
            #     end
            # end
            # @testset "Thin UL" begin
            #     lua   = @inferred ul(a[:,1:n1])
            #     @test lua.U*lua.L ≈ lua.P'*a[:,1:n1]
            # end
            # @testset "Fat UL" begin
            #     lua   = ul(a[1:n1,:])
            #     @test lua.L*lua.U ≈ lua.P*a[1:n1,:]
            # end
        end

        # @testset "UL of Symmetric/Hermitian" begin
        #     for HS in (Hermitian(a'a), Symmetric(a'a))
        #         luhs = ul(HS)
        #         @test luhs.L*luhs.U ≈ luhs.P*Matrix(HS)
        #     end
        # end
    end

    @testset "Singular matrices" for T in (Float64, ComplexF64)
        A = T[1 2; 0 0]
        @test_throws SingularException ul(A)
        @test_throws SingularException ul!(copy(A))
        @test_throws SingularException ul(A; check = true)
        @test_throws SingularException ul!(copy(A); check = true)
        @test !issuccess(ul(A; check = false))
        @test !issuccess(ul!(copy(A); check = false))
        @test_throws ZeroPivotException ul(A, Val(false))
        @test_throws ZeroPivotException ul!(copy(A), Val(false))
        @test_throws ZeroPivotException ul(A, Val(false); check = true)
        @test_throws ZeroPivotException ul!(copy(A), Val(false); check = true)
        @test !issuccess(ul(A, Val(false); check = false))
        @test !issuccess(ul!(copy(A), Val(false); check = false))
        F = ul(A; check = false)
        @test sprint((io, x) -> show(io, "text/plain", x), F) ==
            "Failed factorization of type $(typeof(F))"
    end

    # @testset "conversion" begin
    #     Random.seed!(3)
    #     a = Tridiagonal(rand(9),rand(10),rand(9))
    #     fa = Array(a)
    #     falu = ul(fa)
    #     alu = ul(a)
    #     falu = convert(typeof(falu),alu)
    #     @test Array(alu) == fa
    #     @test AbstractArray(alu) == fa
    # end

    @testset "Rational Matrices" begin
        ## Integrate in general tests when more linear algebra is implemented in julia
        a = convert(Matrix{Rational{BigInt}}, rand(1:10//1,n,n))/n
        b = rand(1:10,n,2)
        @inferred ul(a)
        lua   = factorize(a)
        l,u,p = lua.L, lua.U, lua.p
        @test l*u ≈ a[p,:]
        @test l[invperm(p),:]*u ≈ a
        @test a*inv(lua) ≈ Matrix(I, n, n)
        let Bs = b
            for b in (Bs, view(Bs, 1:n, 1))
                @test a*(lua\b) ≈ b
            end
        end
        @test @inferred(det(a)) ≈ det(Array{Float64}(a))
    end

    @testset "Rational{BigInt} and BigFloat Hilbert Matrix" begin
        ## Hilbert Matrix (very ill conditioned)
        ## Testing Rational{BigInt} and BigFloat version
        nHilbert = 50
        H = Rational{BigInt}[1//(i+j-1) for i = 1:nHilbert,j = 1:nHilbert]
        Hinv = Rational{BigInt}[(-1)^(i+j)*(i+j-1)*binomial(nHilbert+i-1,nHilbert-j)*
            binomial(nHilbert+j-1,nHilbert-i)*binomial(i+j-2,i-1)^2
            for i = big(1):nHilbert,j=big(1):nHilbert]
        @test inv(H) == Hinv
        setprecision(2^10) do
            @test norm(Array{Float64}(inv(float(H)) - float(Hinv))) < 1e-100
        end
    end

    @testset "logdet" begin
        @test @inferred(logdet(ComplexF32[1.0f0 0.5f0; 0.5f0 -1.0f0])) === 0.22314355f0 + 3.1415927f0im
        @test_throws DomainError logdet([1 1; 1 -1])
    end

    @testset "REPL printing" begin
            bf = IOBuffer()
            show(bf, "text/plain", ul(Matrix(I, 4, 4)))
            seekstart(bf)
            @test String(take!(bf)) == """
UL{Float64,Array{Float64,2},Array{Int64,1}}
U factor:
4×4 Array{Float64,2}:
 1.0  0.0  0.0  0.0
 0.0  1.0  0.0  0.0
 0.0  0.0  1.0  0.0
 0.0  0.0  0.0  1.0
L factor:
4×4 Array{Float64,2}:
 1.0  0.0  0.0  0.0
 0.0  1.0  0.0  0.0
 0.0  0.0  1.0  0.0
 0.0  0.0  0.0  1.0"""
    end

    @testset "propertynames" begin
        names = sort!(collect(string.(Base.propertynames(ul(rand(3,3))))))
        @test names == ["L", "P", "U", "p"]
        allnames = sort!(collect(string.(Base.propertynames(ul(rand(3,3)), true))))
        @test allnames == ["L", "P", "U", "factors", "info", "ipiv", "p"]
    end

    @testset "Issue #30917. Determinant of integer matrix" begin
        @test det([1 1 0 0 1 0 0 0
                1 0 1 0 0 1 0 0
                1 0 0 1 0 0 1 0
                0 1 1 1 0 0 0 0
                0 1 0 0 0 0 1 1
                0 0 1 0 1 0 0 1
                0 0 0 1 1 1 0 0
                0 0 0 0 1 1 0 1]) ≈ 6
    end

    @testset "Issue #33177. No ldiv!(UL, Adjoint)" begin
        A = [1 0; 1 1]
        B = [1 2; 2 8]
        F = ul(B)
        @test (A  / F') * B == A
        @test (A' / F') * B == A'

        a = complex.(randn(2), randn(2))
        @test (a' / F') * B ≈ a'
        @test (transpose(a) / F') * B ≈ transpose(a)

        A = complex.(randn(2, 2), randn(2, 2))
        @test (A' / F') * B ≈ A'
        @test (transpose(A) / F') * B ≈ transpose(A)
    end

    @testset "0x0 matrix" begin
        A = ones(0, 0)
        F = ul(A)
        @test F.U == ones(0, 0)
        @test F.L == ones(0, 0)
        @test F.P == ones(0, 0)
        @test F.p == []
    end

    @testset "more ldiv! methods" begin
        for elty in (Float16, Float64, ComplexF64), transform in (transpose, adjoint)
            A = randn(elty, 5, 5)
            B = randn(elty, 5, 5)
            @test ldiv!(transform(ul(A)), copy(B)) ≈ transform(ul(A)) \ B ≈ transform(A) \ B
            @test ldiv!(transform(ul(A)), transform(copy(B))) ≈ transform(ul(A)) \ transform(B) ≈ transform(A) \ transform(B)
        end
    end

    @testset "more rdiv! methods" begin
        A = randn(5,5)
        B = randn(5,5)
        a = randn(5)
        @test rdiv!(copy(A), ul(B)) ≈ A / B

        for elty in (Float16, Float64, ComplexF64), transform in (transpose, adjoint)
            A = randn(elty, 5, 5)
            C = copy(A)
            B = randn(elty, 5, 5)
            @test rdiv!(transform(A), transform(ul(B))) ≈ transform(C) / transform(ul(B)) ≈ transform(C) / transform(B)
        end
    end

    @testset "conversions" begin
        A = randn(5,5)
        F = ul(A)
        @test ul(F) ≡ F ≡ Factorization{Float64}(F)
        @test UL{Float32}(F).factors ≈ UL{Float32,Matrix{Float32}}(F).factors ≈
                Factorization{Float32}(F).factors ≈ ul(Float32.(A)).factors
        @test copy(F).factors == F.factors
    end

    @testset "det" begin
        A = randn(5,5)
        @test det(A) ≈ det(ul(A))
        @test all(logabsdet(A) .≈ logabsdet(ul(A)))
    end

    @testset "inv" begin
        A = randn(5,5)
        @test inv(ul(A)) ≈ inv(A)
    end

    @testset "MyVector/MyMatrix" begin
        A = randn(5,5)
        b = randn(5)
        B = randn(5,3)
        b̃ = MyVector(copy(b))
        B̃ = MyMatrix(copy(B))
        @test ldiv!(ul(A), b̃) === b̃
        @test b̃ ≈ A \ b
        @test ldiv!(ul(A), B̃) === B̃
        @test B̃ ≈ A \ B
    end

    @testset "Tridiagonal" begin
        n = 10
        A = Tridiagonal(randn(n-1), randn(n), randn(n-1))+10I
        U,L = ul(A, Val(false))
        @test U isa Bidiagonal
        @test L isa Bidiagonal
        Ũ,L̃ = ul(Matrix(A), Val(false))
        @test U ≈ Ũ
        @test L ≈ L̃
        @test U*L ≈ A
    end
end


text:
@inline function _cumsum(v::Vector{Int}, endidx::Int)
    s = 0
    @inbounds for i in 1:endidx
        s += v[i]
    end
    return s
end

@inline function _sumiter(v::Vector{Int}, endidx::Int)
    s = 0
    @inbounds for i in 1:endidx
        s += v[i]
    end
    return s
end

ext: jl


text:
function print_block(b::blockrec)
    ccall((:printb, CSDP.csdp), Nothing, (blockrec,), b)
end

function printm(A::blockmatrix)
    ccall((:printm, CSDP.csdp), Nothing, (blockmatrix,), A)
end
printm(A::BlockMatrix) = printm(A.csdp)
export printm

print_sparseblock(A::sparseblock) = print_sparseblock(pointer_from_objref(A))
function print_sparseblock(a::Ptr{sparseblock})
    ccall((:print_sparse_block, CSDP.csdp), Nothing, (Ptr{sparseblock},), a)
end

function print_sizeof()
    ccall((:print_sizeof, CSDP.csdp), Nothing, ())
end

function print_constraints(C::Vector{constraintmatrix})
    ccall((:print_constraints, CSDP.csdp), Nothing, (CSDP_INT, Ptr{constraintmatrix}), length(C), fptr(C))
end

function Base.show(io::IO, b::sparseblock)
    println(io, "\nsparseblock(", pointer_from_objref(b))
    println(io, " next          :  ",  b.next         )
    println(io, " nextbyblock   :  ",  b.nextbyblock  )
    println(io, " entries       :  ",  b.entries      )
    println(io, " iindices      :  ",  b.iindices     )
    println(io, " jindices      :  ",  b.jindices     )
    println(io, " numentries    :  ",  b.numentries   )
    println(io, " blocknum      :  ",  b.blocknum     )
    println(io, " blocksize     :  ",  b.blocksize    )
    println(io, " constraintnum :  ",  b.constraintnum)
    println(io, " issparse      :  ",  b.issparse     )
    println(io, ")")
end


text:
"""
    AnovaEffect

Supertype for all ANOVA data items: `AnovaValue`, `AnovaFactor`, `AnovaResult`
"""
abstract type AnovaEffect end

import Base.-, Base.+
-(x::AnovaEffect, y::AnovaEffect) = AnovaValue(x.name, x.ss - y.ss, x.df - y.df)
+(x::AnovaEffect, y::AnovaEffect) = AnovaValue(x.name, x.ss + y.ss, x.df + y.df)

Broadcast.broadcastable(a::T) where {T <: AnovaEffect} = (a,) # workaround for current behavior


text:
# Sitemap: $(joinpath(path(:site), "sitemap.xml"))
#
# User-agent: *
# Disallow:
#

const DISALLOW = Vector{String}()

"""
$SIGNATURES

Add an entry to `DISALLOW`.
"""
function add_disallow_item()
    loc = url_curpage()
    loc in DISALLOW || push!(DISALLOW, loc)
    return nothing
end

"""
$SIGNATURES

Generate a `robots.txt`, if one already exists, it will be replaced.
"""
function robots_generator()
    dst = joinpath(path(:site), "robots.txt")
    isfile(dst) && rm(dst)
    io = IOBuffer()
    globvar(:generate_sitemap)::Bool && println(io, """
        Sitemap: $(joinpath(globvar(:website_url)::String, "sitemap.xml"))
        """)
    print(io, """
        User-agent: *
        """)
    if !(all(isempty, (DISALLOW, globvar(:robots_disallow)::Vector{String})))
        for page in DISALLOW
            print(io, """
                Disallow: $page
                """)
        end
        for dir in globvar(:robots_disallow)::Vector{String}
            if dir == "/"
                dir = ""
            end
            print(io, """
                Disallow: /$dir
                """)
        end
    else
        print(io, """
            Disallow:
            """)
    end
    write(dst, take!(io))
    return nothing
end


text:
using BinaryProvider

# This is where all binaries will get installed
const prefix = Prefix(!isempty(ARGS) ? ARGS[1] : joinpath(@__DIR__,"usr"))

libsundials_nvecserial       = LibraryProduct(prefix, String["libsundials_nvecserial"], :libsundials_nvecserial)
libsundials_kinsol           = LibraryProduct(prefix, String["libsundials_kinsol"], :libsundials_kinsol)
libsundials_cvodes           = LibraryProduct(prefix, String["libsundials_cvodes"], :libsundials_cvodes)
libsundials_idas             = LibraryProduct(prefix, String["libsundials_idas"], :libsundials_idas)
libsundials_sunlinsolspfgmr  = LibraryProduct(prefix, String["libsundials_sunlinsolspfgmr"], :libsundials_sunlinsolspfgmr)
libsundials_sunmatrixdense   = LibraryProduct(prefix, String["libsundials_sunmatrixdense"], :libsundials_sunmatrixdense)
libsundials_sunlinsolspbcgs  = LibraryProduct(prefix, String["libsundials_sunlinsolspbcgs"], :libsundials_sunlinsolspbcgs)
libsundials_sunlinsoldense   = LibraryProduct(prefix, String["libsundials_sunlinsoldense"], :libsundials_sunlinsoldense)
libsundials_sunlinsolspgmr   = LibraryProduct(prefix, String["libsundials_sunlinsolspgmr"], :libsundials_sunlinsolspgmr)
libsundials_sunlinsolpcg     = LibraryProduct(prefix, String["libsundials_sunlinsolpcg"], :libsundials_sunlinsolpcg)
libsundials_sunlinsolsptfqmr = LibraryProduct(prefix, String["libsundials_sunlinsolsptfqmr"], :libsundials_sunlinsolsptfqmr)
libsundials_sunmatrixsparse  = LibraryProduct(prefix, String["libsundials_sunmatrixsparse"], :libsundials_sunmatrixsparse)
libsundials_sunlinsolband    = LibraryProduct(prefix, String["libsundials_sunlinsolband"], :libsundials_sunlinsolband)
libsundials_sunmatrixband    = LibraryProduct(prefix, String["libsundials_sunmatrixband"], :libsundials_sunmatrixband)
libsundials_arkode           = LibraryProduct(prefix, String["libsundials_arkode"], :libsundials_arkode)
libsundials_sunlinsolklu     = LibraryProduct(prefix, String["libsundials_sunlinsolklu"], :libsundials_sunlinsolklu)

products = [
    libsundials_nvecserial,
    libsundials_kinsol,
    libsundials_cvodes,
    libsundials_idas,
    libsundials_sunlinsolspfgmr,
    libsundials_sunmatrixdense,
    libsundials_sunlinsolspbcgs,
    libsundials_sunlinsoldense,
    libsundials_sunlinsolspgmr,
    libsundials_sunlinsolpcg,
    libsundials_sunlinsolsptfqmr,
    libsundials_sunmatrixsparse,
    libsundials_sunlinsolband,
    libsundials_sunmatrixband,
    libsundials_arkode,
    libsundials_sunlinsolklu,
]

# Download binaries from hosted location
bin_prefix = "https://github.com/JuliaDiffEq/SundialsBuilder/releases/download/sundials-3.1.0-3"

# Listing of files generated by BinaryBuilder:
download_info = Dict(
    BinaryProvider.Linux(:aarch64, :glibc) => ("$bin_prefix/Sundials.aarch64-linux-gnu.tar.gz", "ded75e83b80710015bd5e06b950e47bf5e58cf13c8c7203fbfdc6d20652c8610"),
    BinaryProvider.Linux(:armv7l, :glibc) => ("$bin_prefix/Sundials.arm-linux-gnueabihf.tar.gz", "989f8535dc03940964bbef79918056ddf3922e490c905cf7482752bc8c99d604"),
    BinaryProvider.Linux(:i686, :glibc) => ("$bin_prefix/Sundials.i686-linux-gnu.tar.gz", "6fbc1769dba5534f3c53b15164ef64132027bb44f7a8b33ea06a1a863d7c3503"),
    BinaryProvider.Windows(:i686) => ("$bin_prefix/Sundials.i686-w64-mingw32.tar.gz", "939106914a1bd438ec6cedab3644f09560ab192fa5cabfe3fde9aca04fb26049"),
    BinaryProvider.Linux(:powerpc64le, :glibc) => ("$bin_prefix/Sundials.powerpc64le-linux-gnu.tar.gz", "1e94174b63a6edb052fe2b5374a07a76937833fbcb0236d4b319f19871de267d"),
    BinaryProvider.MacOS() => ("$bin_prefix/Sundials.x86_64-apple-darwin14.tar.gz", "7c9b0a73c334c9e7cd3b7b43a6df8c9275feb7fd02c7240e3804ea29a50dbedb"),
    BinaryProvider.Linux(:x86_64, :glibc) => ("$bin_prefix/Sundials.x86_64-linux-gnu.tar.gz", "10cf64bc564a45459359b07a5a9708b73b0808daff47ef6fc86a530ad1e49b27"),
    BinaryProvider.Windows(:x86_64) => ("$bin_prefix/Sundials.x86_64-w64-mingw32.tar.gz", "4e3e1b1e60b0ec59922b64bab70a2eac66f287aa2f0ac0826481953cd3f81664"),
)
if platform_key() in keys(download_info)
    # First, check to see if we're all satisfied
    if any(!satisfied(p; verbose=true) for p in products)
        # Download and install binaries
        url, tarball_hash = download_info[platform_key()]
        install(url, tarball_hash; prefix=prefix, force=true, verbose=true)
    end

    # Finally write out a deps.jl file
    write_deps_file(joinpath(@__DIR__, "deps.jl"), products)
else
    error("Your platform $(triplet(platform_key())) is not supported by this package!")
end


text:
module Classes

using DataStructures
using MacroTools
using MacroTools:combinedef, combinestructdef
using InteractiveUtils: subtypes

export @class, Class, AbstractClass, isclass, classof, superclass, superclasses, issubclass, subclasses, absclass

abstract type AbstractClass end            # supertype of all shadow class types
abstract type Class <: AbstractClass end   # superclass of all concrete classes

include("utilities.jl")

_cache = nothing

# Return info about a class in a named tuple
function _class_info(::Type{T}) where {T <: AbstractClass}
    global _cache
    _cache === nothing && (_cache = Dict())
    haskey(_cache, T) && return _cache[T]

    # @info "_class_info($T)"

    typ = (typeof(T) === UnionAll ? Base.unwrap_unionall(T) : T)

    # note: must extract symbol from type to create required expression
    ivars = (isabstracttype(typ) ? Expr[] : [:($vname::$(typename(vtype))) for (vname, vtype) in zip(fieldnames(typ), typ.types)])
    wheres = typ.parameters

    d = Dict(t.name=>regensym(t.name) for t in wheres)    # create mapping of type params to gen'd symbols
    ivars = [_translate_ivar(d, iv) for iv in ivars]      # translate types to use gensyms
    wheres = [_translate_where(d, w) for w in wheres]

    result = (wheres=wheres, ivars=ivars, super=superclass(typ))
    _cache[T] = result
    return result
end

"""
    superclass(t::Type{Class})

Returns the type of the concrete superclass of the given class, or `nothing`
for `Class`, which is the root of the class hierarchy.
"""
superclass(::Type{Class}) = nothing

"""
    superclasses(::Type{T}) where {T <: AbstractClass}

Returns a vector of superclasses from the superclass of the current class
to `Class`, in order.
"""
superclasses(::Type{Class}) = []

function superclasses(::Type{T}) where {T <: AbstractClass}
    super = superclass(T)
    [super, superclasses(super)...]
end

"""
    isclass(X)

Return `true` if `X` is a concrete subclass of `AbstractClass`, or is `Class`, which is abstract.
"""
isclass(any) = false

# Note that !isabstracttype(T) != isconcretetype(T): parameterized types return false for both
isclass(::Type{T}) where {T <: AbstractClass} = (T === Class || !isabstracttype(T))

"""
    issubclass(t1::DataType, t2::DataType)

Returns `true` if `t1` is a subclass of `t2`, else false.
"""
# identity
issubclass(::Type{T}, ::Type{T}) where {T <: AbstractClass} = true

issubclass(::Type{T1}, ::Type{T2}) where {T1 <: AbstractClass, T2 <: AbstractClass} = T1 in Set(superclasses(T2))

"""
    classof(::Type{T}) where {T <: AbstractClass}

Compute the concrete class associated with abstract class `T`, which must
be a subclass of `AbstractClass`.
"""
function classof(::Type{T}) where {T <: AbstractClass}
    if isclass(T)
        return T
    end

    # Abstract types should have only one concrete subtype
    concrete = filter(isconcretetype, subtypes(T))

    if length(concrete) == 1
        return concrete[1]
    end

    # Should never happen unless user manually creates errant subtypes
    error("Abstract class supertype $T has multiple concrete subtypes: $concrete")
end

"""
    subclasses(::Type{T}) where {T <: AbstractClass}

Compute the vector of subclasses for a given class.
"""
function subclasses(::Type{T}) where {T <: AbstractClass}
    # immediate supertype is "our" entry in the type hierarchy
    super = supertype(T)

    # collect immediate subclasses
    subs = [classof(t) for t in subtypes(super) if isabstracttype(t)]

    # recurse on subclasses
    return [subs; [subclasses(t) for t in subs]...]
end

"""
    absclass(::Type{T}) where {T <: AbstractClass}

Returns the abstract type associated with the concrete class `T`.
"""
function absclass(::Type{T}) where {T <: AbstractClass}
    isclass(T) ? supertype(T) : error("absclass(T) must be called on concrete classes; $T is abstract.")
end

include("@class.jl")

end # module


text:
using Test,PiecewiseSplicedArrays,PiecewiseIncreasingRanges

@testset "initialization and indexing" begin
	ax = PiecewiseUnitRange([1:3,6:7]);
	a = zeros(5,5);
	psa = PiecewiseSplicedArray(a,ax,ax);
	@test parent(psa) === a
	for i in eachindex(psa)
       psa[i]=i
       @test psa[i]==i
    end
    @test PiecewiseSplicedArrays.parentinds(psa,(2,6)) == (2,4)

    b=@view psa[1:3,6:7];
    @. b= 99
    @test psa[3,6] == 99
end

@testset "size and axes" begin
	ax = PiecewiseUnitRange([1:3,6:7]);
	psa = PiecewiseSplicedArray(zeros(5,5),ax,ax);
	for d in 1:ndims(psa)
		@test axes(psa,d)==ax
	end
	@test length(psa)==25
	@test size(psa)==(5,5)
end


text:
function c2x(c::Array{Complex{Float64},1},tinfo::ProjectedGPE.Tinfo)
    return tinfo.Tx*c
end

function c2x(c::Array{Complex{Float64},2},tinfo::ProjectedGPE.Tinfo)
    return tinfo.Tx*c*tinfo.Ty'
end

function c2x(c::Array{Complex{Float64},3},tinfo::ProjectedGPE.Tinfo)
    sa = size(c);sx = size(tinfo.Tx);sy = size(tinfo.Ty);sz = size(tinfo.Tz)
    Ax = reshape(c,(sa[1],sa[2]*sa[3]))
    Ax = reshape(tinfo.Tx*Ax,(sx[1],sa[2],sa[3]))
    Ax = permutedims(Ax,(2,3,1));sa = size(Ax)

    Ay = reshape(Ax,(sa[1],sa[2]*sa[3]))
    Ay = reshape(tinfo.Ty*Ay,(sy[1],sa[2],sa[3]))
    Ay = permutedims(Ay,(2,3,1));sa = size(Ay)

    Az = reshape(Ay,(sa[1],sa[2]*sa[3]))
    Az = reshape(tinfo.Tz*Az,(sz[1],sa[2],sa[3]))
    return permutedims(Az,(2,3,1))
end


text:
"""
    Spectrum!(spec::AbstractVector{Tv}, f::Function, size::Ti) where {Tv<:Complex, Ti<:Integer}

Set the spectrum from user provided function.

## Examples
```jldoctest; setup = :(using SMG2S; using LinearAlgebra)
julia> function f(idx::Integer, size = 10)
           return 10 * cos(((idx-1) / 2) * 2 * π / size) + sin(((idx-1) / 2) * 2 * π / size) * im
       end
f (generic function with 2 methods)

julia> vec=zeros(ComplexF64, 10);

julia> Spectrum!(vec, f, 10);

julia> vec
10-element Vector{ComplexF64}:
                  10.0 + 0.0im
     9.510565162951535 + 0.3090169943749474im
     8.090169943749475 + 0.5877852522924731im
     5.877852522924732 + 0.8090169943749475im
    3.0901699437494745 + 0.9510565162951535im
 6.123233995736766e-16 + 1.0im
   -3.0901699437494736 + 0.9510565162951536im
     -5.87785252292473 + 0.8090169943749475im
    -8.090169943749473 + 0.5877852522924732im
    -9.510565162951535 + 0.3090169943749475im

```

"""
function Spectrum!(spec::AbstractVector{Tv}, f::Function, size::Ti) where {Tv<:Complex, Ti<:Integer}
    for i = 1:size
        spec[i] = Tv(f(i))
    end
end

"""
    Spectrum!(spec::AbstractVector{Tv}, vec::AbstractVector{Tv}, size::Ti) where {Tv<:Complex, Ti<:Integer}

Set the spectrum from user provided vector. In fact, this function only check if the size of given spectrum vector is large enough.
And they copy the first `size` elements of `vec` into `spec`.

"""
function Spectrum!(spec::AbstractVector{Tv}, vec::AbstractVector{Tv}, size::Ti) where {Tv<:Complex, Ti<:Integer}

    if size(vec, 1) < size
        error("the provided spectrum vector is too small compared to the matrix to be generated.")
    end

    for i = 1:size
        spec[i] = Tv(vec[i])
    end
end

##only for non-symmetric cases
function checkSpectrum(spectrum::AbstractVector{Tv}, size::Ti) where {Tv <: Complex, Ti <: Integer}
    idx = 1

    while idx < size
        if imag(spectrum[idx]) == 0
            step = 1
        else
            if (imag(spectrum[idx]) != -imag(spectrum[idx+1])) && (real(spectrum[idx]) != real(spectrum[idx+1]))
                @info spectrum[idx], spectrum[idx+1],conj(spectrum[idx+1])
                error(
                    "for initialisationg of non-symmetric matrix, the given spectrum is invalid, please follow the instruction",
                    )
            end
            step = 2
        end

        idx += step
    end
end


text:
using Overdot # Bring in the top-level sim environment.
# import HDF5 # For loading the log file (Pkg.add("HDF5");)

# If the out directory doesn't exist, make it.
if !isdir("out")
    mkdir("out")
end

# Run the sim on the simplest little scenario.
scenario = simulate("test-scenario.yaml", @__MODULE__)

# Open the logs and plot some things.
if !isempty(scenario.sim.log_file)

    # # Read some data from the logs. This do block takes care of
    # # error-handling, so we can just read whatever we want and return it.
    # t, q_BI, ω_BI_B = HDF5.h5open(scenario.sim.log_file, "r") do logs
    #     (read(logs, "/cube1/body/state/time"),
    #      read(logs, "/cube1/body/state/data/q_BI"),
    #      read(logs, "/cube1/body/state/data/ω_BI_B"))
    # end
    #
    # # Scalars are logged as 1-by-n (don't know why). Convert to just n for
    # # the sake of plotting.
    # t = dropdims(t, dims=1)

end # plotting


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

JLLWrappers.@generate_main_file_header("GoogleTest")
JLLWrappers.@generate_main_file("GoogleTest", UUID("60d6b632-e33d-5ba5-b345-b8acd80c8739"))
end  # module GoogleTest_jll


text:
seats = reduce(vcat, permutedims.(collect.(readlines("day11/input.txt"))))

function step_seats!(seats)
    copy_seats = copy(seats)
    function execute_seat_rules(i, j)
        hi_i, hi_j = size(seats)
        adjacent_seats = copy_seats[max(1, i-1):min(i+1, hi_i), max(1, j-1):min(j+1, hi_j)]
        if copy_seats[i, j] == 'L' && count(s -> s != '#', adjacent_seats) == length(adjacent_seats)
            seats[i,j] = '#'
        elseif copy_seats[i, j] == '#' && count(s -> s == '#', adjacent_seats) >= 5
            seats[i,j] = 'L'
        end
    end

    for i in 1:size(seats, 1), j in 1:size(seats, 2)
        execute_seat_rules(i, j)
    end
end

old_seats = nothing
while old_seats != seats
    global old_seats = copy(seats)
    step_seats!(seats)
end

println(count(s -> s == '#', seats))


text:
######################## Mesh Generation function with GMSH ####################
function MeshGenerator(L,h1,h2,h3,hd,dpml,l0,ld,lpml)
    gmsh.initialize()
    gmsh.option.setNumber("General.Terminal", 1)
    gmsh.option.setNumber("Mesh.Algorithm", 6)
    gmsh.clear()
    gmsh.model.add("geometry") # name it whatever you want

    # Add points
    gmsh.model.geo.addPoint(-L/2-dpml, -hd-h2-h3-dpml, 0, lpml,  1)
    gmsh.model.geo.addPoint(-L/2     , -hd-h2-h3-dpml, 0, lpml,  2)
    gmsh.model.geo.addPoint(-L/2-dpml, -hd-h2-h3     , 0, lpml,  3)
    gmsh.model.geo.addPoint(-L/2     , -hd-h2-h3     , 0,   l0,  4)
    gmsh.model.geo.addPoint( L/2     , -hd-h2-h3-dpml, 0, lpml,  5)
    gmsh.model.geo.addPoint( L/2     , -hd-h2-h3     , 0,   l0,  6)
    gmsh.model.geo.addPoint( L/2+dpml, -hd-h2-h3-dpml, 0, lpml,  7)
    gmsh.model.geo.addPoint( L/2+dpml, -hd-h2-h3     , 0, lpml,  8)
    gmsh.model.geo.addPoint(-L/2-dpml, -hd-h2        , 0,   l0,  9)
    gmsh.model.geo.addPoint(-L/2     , -hd-h2        , 0,   l0, 10)
    gmsh.model.geo.addPoint( L/2     , -hd-h2        , 0,   l0, 11)
    gmsh.model.geo.addPoint( L/2+dpml, -hd-h2        , 0,   l0, 12)
    gmsh.model.geo.addPoint(-L/2-dpml, -hd           , 0,   ld, 13)
    gmsh.model.geo.addPoint(-L/2     , -hd           , 0,   ld, 14)
    gmsh.model.geo.addPoint( L/2     , -hd           , 0,   ld, 15)
    gmsh.model.geo.addPoint( L/2+dpml, -hd           , 0,   ld, 16)
    gmsh.model.geo.addPoint(-L/2-dpml, 0.            , 0,   ld, 17)
    gmsh.model.geo.addPoint(-L/2     , 0.            , 0,   ld, 18)
    gmsh.model.geo.addPoint( L/2     , 0.            , 0,   ld, 19)
    gmsh.model.geo.addPoint( L/2+dpml, 0.            , 0,   ld, 20)
    gmsh.model.geo.addPoint(-L/2-dpml, h1            , 0, lpml, 21)
    gmsh.model.geo.addPoint(-L/2     , h1            , 0,   l0, 22)
    gmsh.model.geo.addPoint( L/2     , h1            , 0,   l0, 23)
    gmsh.model.geo.addPoint( L/2+dpml, h1            , 0, lpml, 24)
    gmsh.model.geo.addPoint(-L/2-dpml, h1+dpml       , 0, lpml, 25)
    gmsh.model.geo.addPoint(-L/2     , h1+dpml       , 0, lpml, 26)
    gmsh.model.geo.addPoint( L/2     , h1+dpml       , 0, lpml, 27)
    gmsh.model.geo.addPoint( L/2+dpml, h1+dpml       , 0, lpml, 28)
    #gmsh.model.geo.addPoint(-Lt/2, h1/2-ht/2, 0, l0, 29)
    #gmsh.model.geo.addPoint( Lt/2, h1/2-ht/2, 0, l0, 30)
    #gmsh.model.geo.addPoint( Lt/2, h1/2+ht/2, 0, l0, 31)
    #gmsh.model.geo.addPoint(-Lt/2, h1/2+ht/2, 0, l0, 32)
    #gmsh.model.geo.addPoint(x_s,y_s, 0, l_s, 17)
    # Add lines
    gmsh.model.geo.addLine(  1,  2,  1)
    gmsh.model.geo.addLine(  2,  4,  2)
    gmsh.model.geo.addLine(  4,  3,  3)
    gmsh.model.geo.addLine(  3,  1,  4)
    gmsh.model.geo.addLine(  2,  5,  5)
    gmsh.model.geo.addLine(  5,  6,  6)
    gmsh.model.geo.addLine(  6,  4,  7)
    gmsh.model.geo.addLine(  5,  7,  8)
    gmsh.model.geo.addLine(  7,  8,  9)
    gmsh.model.geo.addLine(  8,  6, 10)
    gmsh.model.geo.addLine(  4, 10, 11)
    gmsh.model.geo.addLine( 10,  9, 12)
    gmsh.model.geo.addLine(  9,  3, 13)
    gmsh.model.geo.addLine(  6, 11, 14)
    gmsh.model.geo.addLine( 11, 10, 15)
    gmsh.model.geo.addLine(  8, 12, 16)
    gmsh.model.geo.addLine( 12, 11, 17)
    gmsh.model.geo.addLine( 10, 14, 18)
    gmsh.model.geo.addLine( 14, 13, 19)
    gmsh.model.geo.addLine( 13,  9, 20)
    gmsh.model.geo.addLine( 11, 15, 21)
    gmsh.model.geo.addLine( 15, 14, 22)
    gmsh.model.geo.addLine( 12, 16, 23)
    gmsh.model.geo.addLine( 16, 15, 24)
    gmsh.model.geo.addLine( 14, 18, 25)
    gmsh.model.geo.addLine( 18, 17, 26)
    gmsh.model.geo.addLine( 17, 13, 27)
    gmsh.model.geo.addLine( 15, 19, 28)
    gmsh.model.geo.addLine( 19, 18, 29)
    gmsh.model.geo.addLine( 16, 20, 30)
    gmsh.model.geo.addLine( 20, 19, 31)
    gmsh.model.geo.addLine( 18, 22, 32)
    gmsh.model.geo.addLine( 22, 21, 33)
    gmsh.model.geo.addLine( 21, 17, 34)
    gmsh.model.geo.addLine( 19, 23, 35)
    gmsh.model.geo.addLine( 23, 22, 36)
    gmsh.model.geo.addLine( 20, 24, 37)
    gmsh.model.geo.addLine( 24, 23, 38)
    gmsh.model.geo.addLine( 22, 26, 39)
    gmsh.model.geo.addLine( 26, 25, 40)
    gmsh.model.geo.addLine( 25, 21, 41)
    gmsh.model.geo.addLine( 23, 27, 42)
    gmsh.model.geo.addLine( 27, 26, 43)
    gmsh.model.geo.addLine( 24, 28, 44)
    gmsh.model.geo.addLine( 28, 27, 45)
    #gmsh.model.geo.addLine( 29, 30, 46)
    #gmsh.model.geo.addLine( 30, 31, 47)
    #gmsh.model.geo.addLine( 31, 32, 48)
    #gmsh.model.geo.addLine( 32, 29, 49)
    # Construct curve loops and surfaces
    gmsh.model.geo.addCurveLoop([  1,  2,  3,  4], 1)
    gmsh.model.geo.addCurveLoop([  5,  6,  7, -2], 2)
    gmsh.model.geo.addCurveLoop([  8,  9, 10, -6], 3)
    gmsh.model.geo.addCurveLoop([ 11, 12, 13, -3], 4)
    gmsh.model.geo.addCurveLoop([ -7, 14, 15,-11], 5)
    gmsh.model.geo.addCurveLoop([ 16, 17,-14,-10], 6)
    gmsh.model.geo.addCurveLoop([ 18, 19, 20,-12], 7)
    gmsh.model.geo.addCurveLoop([-15, 21, 22,-18], 8)
    gmsh.model.geo.addCurveLoop([-17, 23, 24,-21], 9)
    gmsh.model.geo.addCurveLoop([-19, 25, 26, 27],10)
    gmsh.model.geo.addCurveLoop([-22, 28, 29,-25],11)
    gmsh.model.geo.addCurveLoop([-24, 30, 31,-28],12)
    gmsh.model.geo.addCurveLoop([-26, 32, 33, 34],13)
    gmsh.model.geo.addCurveLoop([-29, 35, 36,-32],14)
    gmsh.model.geo.addCurveLoop([-31, 37, 38,-35],15)
    gmsh.model.geo.addCurveLoop([-33, 39, 40, 41],16)
    gmsh.model.geo.addCurveLoop([-36, 42, 43,-39],17)
    gmsh.model.geo.addCurveLoop([-38, 44, 45,-42],18)
    #gmsh.model.geo.addCurveLoop([ 46, 47, 48, 49],19)
    #gmsh.model.geo.addCurveLoop([29, 30], 11)
    gmsh.model.geo.addPlaneSurface([ 1],  1)
    gmsh.model.geo.addPlaneSurface([ 2],  2)
    gmsh.model.geo.addPlaneSurface([ 3],  3)
    gmsh.model.geo.addPlaneSurface([ 4],  4)
    gmsh.model.geo.addPlaneSurface([ 5],  5)
    gmsh.model.geo.addPlaneSurface([ 6],  6)
    gmsh.model.geo.addPlaneSurface([ 7],  7)
    gmsh.model.geo.addPlaneSurface([ 8],  8)
    gmsh.model.geo.addPlaneSurface([ 9],  9)
    gmsh.model.geo.addPlaneSurface([10], 10)
    gmsh.model.geo.addPlaneSurface([11], 11)
    gmsh.model.geo.addPlaneSurface([12], 12)
    gmsh.model.geo.addPlaneSurface([13], 13)
    gmsh.model.geo.addPlaneSurface([14], 14)
    gmsh.model.geo.addPlaneSurface([15], 15)
    gmsh.model.geo.addPlaneSurface([16], 16)
    gmsh.model.geo.addPlaneSurface([17], 17)
    gmsh.model.geo.addPlaneSurface([18], 18)
    #gmsh.model.geo.addPlaneSurface([19], 19)
    # Physical groups
    gmsh.model.addPhysicalGroup(0, [1,3,21,25], 1)
    gmsh.model.setPhysicalName(0, 1, "LeftNodes")
    gmsh.model.addPhysicalGroup(0, [7,8,24,28], 2)
    gmsh.model.setPhysicalName(0, 2, "RightNodes")
    gmsh.model.addPhysicalGroup(0, [2,5], 3)
    gmsh.model.setPhysicalName(0, 3, "BottomNodes")
    gmsh.model.addPhysicalGroup(0, [26,27], 4)
    gmsh.model.setPhysicalName(0, 4, "TopNodes")
    gmsh.model.addPhysicalGroup(0, [13,17], 5)
    gmsh.model.setPhysicalName(0, 5, "LeftDNodes")
    gmsh.model.addPhysicalGroup(0, [16,20], 6)
    gmsh.model.setPhysicalName(0, 6, "RightDNodes")
    gmsh.model.addPhysicalGroup(0, [14,15,18,19], 7)
    gmsh.model.setPhysicalName(0, 7, "TBDNodes")
    gmsh.model.addPhysicalGroup(0, [4,6,22,23], 8)
    gmsh.model.setPhysicalName(0, 8, "InnerNodes")
    gmsh.model.addPhysicalGroup(0, [9,10,11,12], 9)
    gmsh.model.setPhysicalName(0, 9, "SourceNodes")

    gmsh.model.addPhysicalGroup(1, [4,13,20,34,41], 10)
    gmsh.model.setPhysicalName(1, 10, "LeftEdge")
    gmsh.model.addPhysicalGroup(1, [9,16,23,37,44], 11)
    gmsh.model.setPhysicalName(1, 11, "RightEdge")
    gmsh.model.addPhysicalGroup(1, [1,5,8], 12)
    gmsh.model.setPhysicalName(1, 12, "BottomEdge")
    gmsh.model.addPhysicalGroup(1, [40,43,45], 13)
    gmsh.model.setPhysicalName(1, 13, "TopEdge")
    gmsh.model.addPhysicalGroup(1, [27], 14)
    gmsh.model.setPhysicalName(1, 14, "LeftDEdge")
    gmsh.model.addPhysicalGroup(1, [30], 15)
    gmsh.model.setPhysicalName(1, 15, "RightDEdge")
    gmsh.model.addPhysicalGroup(1, [19,22,24,26,29,31], 16)
    gmsh.model.setPhysicalName(1, 16, "TBDEdge")
    gmsh.model.addPhysicalGroup(1, [2,3,6,7,10,11,14,18,21,25,28,32,33,35,36,38,39,42], 17)
    gmsh.model.setPhysicalName(1, 17, "InnerEdge")
    gmsh.model.addPhysicalGroup(1, [46,47,48,49], 18)
    gmsh.model.setPhysicalName(1, 18, "TargetEdge")
    gmsh.model.addPhysicalGroup(1, [12,15,17], 19)
    gmsh.model.setPhysicalName(1, 19, "SourceEdge")

    gmsh.model.addPhysicalGroup(2, [10,11,12], 20)
    gmsh.model.setPhysicalName(2, 20, "Design")
    gmsh.model.addPhysicalGroup(2, [16,17,18], 21)
    gmsh.model.setPhysicalName(2, 21, "TopPML")
    gmsh.model.addPhysicalGroup(2, [1,2,3], 22)
    gmsh.model.setPhysicalName(2, 22, "BottomPML")
    gmsh.model.addPhysicalGroup(2, [14], 23)
    gmsh.model.setPhysicalName(2, 23, "Target")
    gmsh.model.addPhysicalGroup(2, [4,5,6,7,8,9], 24)
    gmsh.model.setPhysicalName(2, 24, "Incident")
    gmsh.model.addPhysicalGroup(2, [13,15], 25)
    gmsh.model.setPhysicalName(2, 25, "LRTarget")
    gmsh.model.geo.synchronize()
    gmsh.model.mesh.generate(2)

    # ... and save it to disk
    gmsh.write("geometry.msh")
    gmsh.finalize()
end

max_stars_repo_name: Tiamat-Tech/hpinn

max_issues_repo_name: louhz/hpinn

max_forks_repo_name: louhz/hpinn


text:
rlxex1 = ex1()
rlx1 = MPCCSolver.RlxMPCC(rlxex1, 1.0, 1.0, 0.0, 1.0) #Relaxation KDB avec t=1.0

#Test 1: meta check
@test rlx1.mod.meta.ncc == 1
@test rlx1.meta.nvar == 2
@test rlx1.meta.ncon == 4

#Test 3: check obj et grad
#rlxex1 has a linear objective function, so the gradient is constant and hessian vanishes
@test MPCCSolver.obj(rlx1,rlx1.meta.x0) == 0.0
@test MPCCSolver.grad(rlx1,rlx1.meta.x0) == MPCCSolver.grad(rlxex1,[Inf,Inf])
@test MPCCSolver.hess(rlx1,rlx1.meta.x0) == sparse(zeros(2,2))
#Htest = MPCCSolver.hess(rlx1, [1.,1.], obj_weight = 1.0, y = ones(4))

@test MPCCSolver.cons(rlx1, [1.,1.])             == [-1.0, 2.0, 2.0, 0.0]
@test norm(MPCCSolver.viol(rlx1, [1.,1.]),Inf)   == 0.0
@test norm(MPCCSolver.viol(rlx1, [0.,1.]),Inf)   == 0.0
@test norm(MPCCSolver.viol(rlx1, [0.,-1.5]),Inf) == 2.5

#Test 4: Stationarité
Spoint = [0.0, 1.0]
Srlx, yS = [-1.,1.], [-1.,-1.,0.,0.]

@test norm(MPCCSolver.jac(rlx1, Srlx)'*yS+MPCCSolver.grad(rlx1, Srlx),Inf) == 0.
@test norm(MPCCSolver.viol(rlx1,Srlx),Inf) == 0.0

#Test 5: on va maintenant résoudre ce problème avec IpOpt:
x0 = [1.0, 1.0]
stats = ipopt(rlx1, print_level = 0, x0 = x0)
#Use y0 (general), zL (lower bound), zU (upper bound)
#for initial guess of Lagrange multipliers.
@show stats.solution, stats.status

nlp_at_x = NLPAtX(x0)
stop = NLPStopping(rlx1, nlp_at_x, optimality_check = KKT)

solveIpopt(stop)
@show stop.current_state.x, status(stop)


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

mutable struct VirtualMachineScaleSetUpdateOSProfile <: SwaggerModel
    customData::Any # spec type: Union{ Nothing, String } # spec name: customData
    windowsConfiguration::Any # spec type: Union{ Nothing, WindowsConfiguration } # spec name: windowsConfiguration
    linuxConfiguration::Any # spec type: Union{ Nothing, LinuxConfiguration } # spec name: linuxConfiguration
    secrets::Any # spec type: Union{ Nothing, Vector{VaultSecretGroup} } # spec name: secrets

    function VirtualMachineScaleSetUpdateOSProfile(;customData=nothing, windowsConfiguration=nothing, linuxConfiguration=nothing, secrets=nothing)
        o = new()
        validate_property(VirtualMachineScaleSetUpdateOSProfile, Symbol("customData"), customData)
        setfield!(o, Symbol("customData"), customData)
        validate_property(VirtualMachineScaleSetUpdateOSProfile, Symbol("windowsConfiguration"), windowsConfiguration)
        setfield!(o, Symbol("windowsConfiguration"), windowsConfiguration)
        validate_property(VirtualMachineScaleSetUpdateOSProfile, Symbol("linuxConfiguration"), linuxConfiguration)
        setfield!(o, Symbol("linuxConfiguration"), linuxConfiguration)
        validate_property(VirtualMachineScaleSetUpdateOSProfile, Symbol("secrets"), secrets)
        setfield!(o, Symbol("secrets"), secrets)
        o
    end
end # type VirtualMachineScaleSetUpdateOSProfile

const _property_map_VirtualMachineScaleSetUpdateOSProfile = Dict{Symbol,Symbol}(Symbol("customData")=>Symbol("customData"), Symbol("windowsConfiguration")=>Symbol("windowsConfiguration"), Symbol("linuxConfiguration")=>Symbol("linuxConfiguration"), Symbol("secrets")=>Symbol("secrets"))
const _property_types_VirtualMachineScaleSetUpdateOSProfile = Dict{Symbol,String}(Symbol("customData")=>"String", Symbol("windowsConfiguration")=>"WindowsConfiguration", Symbol("linuxConfiguration")=>"LinuxConfiguration", Symbol("secrets")=>"Vector{VaultSecretGroup}")
Base.propertynames(::Type{ VirtualMachineScaleSetUpdateOSProfile }) = collect(keys(_property_map_VirtualMachineScaleSetUpdateOSProfile))
Swagger.property_type(::Type{ VirtualMachineScaleSetUpdateOSProfile }, name::Symbol) = Union{Nothing,eval(Meta.parse(_property_types_VirtualMachineScaleSetUpdateOSProfile[name]))}
Swagger.field_name(::Type{ VirtualMachineScaleSetUpdateOSProfile }, property_name::Symbol) =  _property_map_VirtualMachineScaleSetUpdateOSProfile[property_name]

function check_required(o::VirtualMachineScaleSetUpdateOSProfile)
    true
end

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


text:
using Arbiter.ArbiterTasks
using Arbiter.ArbiterSchedulers
import Arbiter.ArbiterGraphs: NodeSet, ImmutableNodeSet

function create_task(name, dependencies=())
    ArbiterTask(name, () -> nothing, dependencies)
end

facts("empty") do
    scheduler = Scheduler()

    @fact completed(scheduler) --> ImmutableNodeSet()
    @fact failed(scheduler) --> ImmutableNodeSet()
    @fact running(scheduler) --> ImmutableNodeSet()
    @fact runnable(scheduler) --> ImmutableNodeSet()
    @fact scheduler --> isfinished
    @fact start_task!(scheduler) --> isnull

    remove_unrunnable!(scheduler)

    @fact completed(scheduler) --> ImmutableNodeSet()
    @fact failed(scheduler) --> ImmutableNodeSet()
    @fact running(scheduler) --> ImmutableNodeSet()
    @fact runnable(scheduler) --> ImmutableNodeSet()
    @fact scheduler --> isfinished
    @fact start_task!(scheduler) --> isnull

    fail_remaining!(scheduler)

    @fact completed(scheduler) --> ImmutableNodeSet()
    @fact failed(scheduler) --> ImmutableNodeSet()
    @fact running(scheduler) --> ImmutableNodeSet()
    @fact runnable(scheduler) --> ImmutableNodeSet()
    @fact scheduler --> isfinished
    @fact start_task!(scheduler) --> isnull
end

facts("add_task") do
    # add a task to Scheduler

    scheduler = Scheduler()

    # no dependencies
    add_task!(scheduler, create_task(:foo))

    @fact completed(scheduler) --> ImmutableNodeSet()
    @fact failed(scheduler) --> ImmutableNodeSet()
    @fact running(scheduler) --> ImmutableNodeSet()
    @fact runnable(scheduler) --> ImmutableNodeSet((:foo,))
    @fact scheduler --> not(isfinished)

    # 1 dependency
    add_task!(scheduler, create_task(:bar, (:foo,)))

    @fact completed(scheduler) --> ImmutableNodeSet()
    @fact failed(scheduler) --> ImmutableNodeSet()
    @fact running(scheduler) --> ImmutableNodeSet()
    @fact runnable(scheduler) --> ImmutableNodeSet((:foo,))
    @fact scheduler --> not(isfinished)

    # invalid tasks

    # circular dependencies
    add_task!(scheduler, create_task(:ouroboros, (:ouroboros,)))

    @fact completed(scheduler) --> ImmutableNodeSet()
    @fact failed(scheduler) --> ImmutableNodeSet((:ouroboros,))
    @fact running(scheduler) --> ImmutableNodeSet()
    @fact runnable(scheduler) --> ImmutableNodeSet((:foo,))
    @fact scheduler --> not(isfinished)

    # dependency made circular
    add_task!(scheduler, create_task(:tick, (:tock,)))

    @fact completed(scheduler) --> ImmutableNodeSet()
    @fact failed(scheduler) --> ImmutableNodeSet((:ouroboros,))
    @fact running(scheduler) --> ImmutableNodeSet()
    @fact runnable(scheduler) --> ImmutableNodeSet((:foo,))
    @fact scheduler --> not(isfinished)

    add_task!(scheduler, create_task(:tock, (:tick,)))

    @fact completed(scheduler) --> ImmutableNodeSet()
    @fact failed(scheduler) --> ImmutableNodeSet((:ouroboros, :tick, :tock))
    @fact runnable(scheduler) --> ImmutableNodeSet((:foo,))
    @fact scheduler --> not(isfinished)

    at_init = Scheduler(
        (
            create_task(:foo),
            create_task(:bar, (:foo,)),
            create_task(:ipsum, (:lorem,)),
            create_task(:ouroboros, (:ouroboros,)),
            create_task(:tick, (:tock,)),
            create_task(:tock, (:tick,)),
        )
    )

    @fact completed(at_init) --> ImmutableNodeSet()
    @fact failed(at_init) --> ImmutableNodeSet((:ouroboros, :tick, :tock))
    @fact runnable(at_init) --> ImmutableNodeSet((:foo,))
    @fact at_init --> not(isfinished)
end

facts("remove_unrunnable!") do
    # remove unrunnable Scheduler tasks

    scheduler = Scheduler(
        (
            create_task(:foo),
            create_task(:bar, (:foo,)),
            create_task(:baz, (:bar,)),
            create_task(:ipsum, (:lorem,)),
            create_task(:dolor, (:ipsum,)),
            create_task(:sit, (:dolor, :stand)),
            create_task(:stand),
        )
    )

    @fact completed(scheduler) --> ImmutableNodeSet()
    @fact failed(scheduler) --> ImmutableNodeSet()
    @fact running(scheduler) --> ImmutableNodeSet()
    @fact runnable(scheduler) --> ImmutableNodeSet((:foo, :stand))
    @fact scheduler --> not(isfinished)

    remove_unrunnable!(scheduler)

    @fact completed(scheduler) --> ImmutableNodeSet()
    @fact failed(scheduler) --> ImmutableNodeSet((:ipsum, :dolor, :sit))
    @fact running(scheduler) --> ImmutableNodeSet()
    @fact runnable(scheduler) --> ImmutableNodeSet((:foo, :stand))
    @fact scheduler --> not(isfinished)
end

facts("start_task!") do
    # start a task

    scheduler = Scheduler(
        (
            create_task(:foo),
            create_task(:fighters, (:foo,)),
            create_task(:bar, (:foo,)),
            create_task(:baz, (:bar,)),
            create_task(:bell, (:bar,)),
            create_task(:node),
        )
    )

    @fact completed(scheduler) --> ImmutableNodeSet()
    @fact failed(scheduler) --> ImmutableNodeSet()
    @fact running(scheduler) --> ImmutableNodeSet()
    @fact runnable(scheduler) --> ImmutableNodeSet((:foo, :node))
    @fact scheduler --> not(isfinished)

    # start a specific task
    task = start_task!(scheduler, :node)
    @fact task --> not(isnull)
    @fact get(task).name --> :node

    @fact completed(scheduler) --> ImmutableNodeSet()
    @fact failed(scheduler) --> ImmutableNodeSet()
    @fact running(scheduler) --> ImmutableNodeSet((:node,))
    @fact runnable(scheduler) --> ImmutableNodeSet((:foo,))
    @fact scheduler --> not(isfinished)

    # start tasks invalidly
    @fact_throws ErrorException start_task!(scheduler, :node)
    @fact_throws ErrorException start_task!(scheduler, :bar)
    @fact_throws ErrorException start_task!(scheduler, :fake)

    @fact completed(scheduler) --> ImmutableNodeSet()
    @fact failed(scheduler) --> ImmutableNodeSet()
    @fact running(scheduler) --> ImmutableNodeSet((:node,))
    @fact runnable(scheduler) --> ImmutableNodeSet((:foo,))
    @fact scheduler --> not(isfinished)

    # is node still stoppable
    end_task!(scheduler, :node)

    @fact completed(scheduler) --> ImmutableNodeSet((:node,))
    @fact failed(scheduler) --> ImmutableNodeSet()
    @fact running(scheduler) --> ImmutableNodeSet()
    @fact runnable(scheduler) --> ImmutableNodeSet((:foo,))
    @fact scheduler --> not(isfinished)

    # start an arbitrary task
    task = start_task!(scheduler)
    @fact task --> not(isnull)
    @fact get(task).name --> :foo

    @fact completed(scheduler) --> ImmutableNodeSet((:node,))
    @fact failed(scheduler) --> ImmutableNodeSet()
    @fact running(scheduler) --> ImmutableNodeSet((:foo,))
    @fact runnable(scheduler) --> ImmutableNodeSet()
    @fact scheduler --> not(isfinished)

    # no startable tasks
    @fact start_task!(scheduler) --> isnull

    @fact completed(scheduler) --> ImmutableNodeSet((:node,))
    @fact failed(scheduler) --> ImmutableNodeSet()
    @fact running(scheduler) --> ImmutableNodeSet((:foo,))
    @fact runnable(scheduler) --> ImmutableNodeSet()
    @fact scheduler --> not(isfinished)

    # start an arbitrary task
    end_task!(scheduler, :foo)

    started = start_task!(scheduler)
    @fact started --> not(isnull)
    name = get(started).name
    @fact name --> anyof(:bar, :fighters)

    @fact completed(scheduler) --> ImmutableNodeSet((:node, :foo))
    @fact failed(scheduler) --> ImmutableNodeSet()

    if name === :bar
        @fact running(scheduler) --> ImmutableNodeSet((:bar,))
        @fact runnable(scheduler) --> ImmutableNodeSet((:fighters,))
    else
        @fact running(scheduler) --> ImmutableNodeSet((:fighters,))
        @fact runnable(scheduler) --> ImmutableNodeSet((:bar,))
    end

    @fact scheduler --> not(isfinished)
end

facts("end_task!") do
    # end a task

    scheduler = Scheduler(
        (
            create_task(:foo),
            create_task(:fighters, (:foo,)),
            create_task(:bar, (:foo,)),
            create_task(:baz, (:bar,)),
            create_task(:qux, (:baz,)),
            create_task(:bell, (:bar,)),
        )
    )

    @fact completed(scheduler) --> ImmutableNodeSet()
    @fact failed(scheduler) --> ImmutableNodeSet()
    @fact running(scheduler) --> ImmutableNodeSet()
    @fact runnable(scheduler) --> ImmutableNodeSet((:foo,))
    @fact scheduler --> not(isfinished)

    # end a task
    start_task!(scheduler, :foo)
    end_task!(scheduler, :foo)

    @fact completed(scheduler) --> ImmutableNodeSet((:foo,))
    @fact failed(scheduler) --> ImmutableNodeSet()
    @fact running(scheduler) --> ImmutableNodeSet()
    @fact runnable(scheduler) --> ImmutableNodeSet((:bar, :fighters))
    @fact scheduler --> not(isfinished)

    # invalid ends
    @fact_throws KeyError end_task!(scheduler, :foo)
    @fact_throws KeyError end_task!(scheduler, :bar)
    @fact_throws KeyError end_task!(scheduler, :baz)

    @fact completed(scheduler) --> ImmutableNodeSet((:foo,))
    @fact failed(scheduler) --> ImmutableNodeSet()
    @fact running(scheduler) --> ImmutableNodeSet()
    @fact runnable(scheduler) --> ImmutableNodeSet((:bar, :fighters))
    @fact scheduler --> not(isfinished)

    # fail a task
    start_task!(scheduler, :bar)
    end_task!(scheduler, :bar, false)

    @fact completed(scheduler) --> ImmutableNodeSet((:foo,))
    @fact failed(scheduler) --> ImmutableNodeSet((:bar, :baz, :qux, :bell))
    @fact running(scheduler) --> ImmutableNodeSet()
    @fact runnable(scheduler) --> ImmutableNodeSet((:fighters,))
    @fact scheduler --> not(isfinished)
end

facts("fail_remaining!") do
    # stop the scheduler

    scheduler = Scheduler(
        (
            create_task(:foo),
            create_task(:fighters, (:foo,)),
            create_task(:bar, (:foo,)),
            create_task(:baz, (:bar,)),
            create_task(:qux, (:baz,)),
            create_task(:bell, (:bar,)),
            create_task(:node),
        )
    )

    start_task!(scheduler, :foo)
    end_task!(scheduler, :foo)
    start_task!(scheduler, :bar)

    fail_remaining!(scheduler)

    @fact completed(scheduler) --> ImmutableNodeSet((:foo,))
    @fact failed(scheduler) --> ImmutableNodeSet((:bar, :baz, :qux, :bell, :fighters, :node))
    @fact running(scheduler) --> ImmutableNodeSet()
    @fact runnable(scheduler) --> ImmutableNodeSet()
    @fact scheduler --> isfinished

    # did that break adding tasks
    add_task!(scheduler, create_task(:restart))

    @fact completed(scheduler) --> ImmutableNodeSet((:foo,))
    @fact failed(scheduler) --> ImmutableNodeSet((:bar, :baz, :qux, :bell, :fighters, :node))
    @fact running(scheduler) --> ImmutableNodeSet()
    @fact runnable(scheduler) --> ImmutableNodeSet((:restart,))
    @fact scheduler --> not(isfinished)
end

facts("do notation") do
    # use do notation with the scheduler

    completed_nodes = NodeSet()
    failed_nodes = NodeSet()
    tasks = (
        create_task(:foo),
        create_task(:bar, (:foo,)),
        create_task(:baz, (:bar,)),
        create_task(:bell, (:bar,)),
        create_task(:lorem, ()),
        create_task(:ipsum, (:lorem,)),
        create_task(:node),
        create_task(:failed, (:fake,)),
    )

    Scheduler(tasks, completed_nodes, failed_nodes) do scheduler
        @fact completed_nodes --> NodeSet()
        @fact failed_nodes --> NodeSet((:failed,))

        start_task!(scheduler, :foo)
        end_task!(scheduler, :foo)

        @fact completed_nodes --> NodeSet((:foo,))
        @fact failed_nodes --> NodeSet((:failed,))

        start_task!(scheduler, :lorem)
        end_task!(scheduler, :lorem, false)

        @fact completed_nodes --> NodeSet((:foo,))
        @fact failed_nodes --> NodeSet((:failed, :lorem, :ipsum))

        start_task!(scheduler, :bar)
    end

    @fact completed_nodes --> NodeSet((:foo,))
    @fact failed_nodes --> NodeSet((:failed, :lorem, :ipsum, :bar, :baz, :bell, :node))
end


text:
# ---
# jupyter:
#   jupytext:
#     text_representation:
#       extension: .jl
#       format_name: light
#       format_version: '1.5'
#       jupytext_version: 1.5.0
#   kernelspec:
#     display_name: Julia 1.4.2
#     language: julia
#     name: julia-1.4
# ---

using LazySets
using LazySets.Approximations
using Plots

]st LazySets

hs1 = HalfSpace([-1.0, 1.0], 0.0)
hs2 = HalfSpace([1.,1.],0.)
hs3 = HalfSpace([0.,-1.],0.5)
p=plot(aspect_ratio=:equal,xlim=[-1,1])
plot!(p, hs1, color=:red)
plot!(p, hs2, color=:blue)
plot!(p, hs3, color=:green)

plot([hs1, hs2, hs3], aspect_ratio=:equal)

# +
p=plot(xlim=[-1,1],ylim=[-1,1],aspect_ratio=:equal)

ph = HPolyhedron([hs1, hs2, hs3])
P_as_polytope = convert(HPolytope, ph);
P_as_polytope |> plot!


text:
module Rabbit
    # Подключаем библиотеку, связанную с временем
    include("../utils/Utils.jl")
    import .Utils

    # Подключаем библиотеку для работы с обменником
    include("./src/Exchanges.jl")
    import .Exchanges
    # Подключаем библиотеку для работы с очередью
    include("./src/Queues.jl")
    import .Queues
    # Подключаем библиотеку для работы с сообщениями
    include("./src/Messages.jl")
    import .Messages

    using AMQPClient

    export shotdown, get_channel, get_auth_params, get_connection

    """Функция закрывает соединение.

    connection: активное соединение с сервисом очередей
    """
    function shotdown(conection::AMQPClient.MessageChannel)
        # Выводим информационное сообщение
        @info Utils.Chronometer.message_with_time("Вызывана функция остановки сервера...")
        # Логер для режима отладки
        @debug "Параметры вызова функции:"
        @debug "> соединение с сервером очередей RabbitMQ" conection
        # Проверяем открыто ли соединение на данный момент
        if AMQPClient.isopen(conection)
            # Закрываем соединение
            AMQPClient.close(conection)
            # Закрытие соединеня асинхронная операция. Ждём завершения для всех каналов
            AMQPClient.wait_for_state(conection, AMQPClient.CONN_STATE_CLOSED)
        end
    end

    """Функция создания канала для соединения RabbitMQ.

    connection: активное соединение с сервисом очередей,
    chanid: идентификатор канала,
    create: опция, которая показывает требуется ли создать канал, если его не существует
    """
    function get_channel(connection::AMQPClient.MessageChannel, chanid::Union{String, Nothing}, create::Bool)::AMQPClient.MessageChannel
        #
        @info Utils.Chronometer.message_with_time("Вызвана функция создания канала для соединения в очереди на сервисе очередей...")
        # Логерд для режима отладки
        @debug "Параметры вызова функции:"
        @debug "> соединение с сервисом очередей" connection
        @debug "> идентификатор канала" chanid
        @debug "> необходимость создания канала (если его нет)" create
        # Инициализируем идентификатор канала.
        # Если он передан в параметрах функции - то присваиваем его
        # Если он не передан в параметрах функции - то генерим начение
        channel_id = chanid == nothing ? AMQPClient.UNUSED_CHANNEL : parse(Int64, chanid)
        # Создаем новый канал для соединения
        try
            # Выполняем операцию создания/поиска канала
            result::AMQPClient.MessageChannel = AMQPClient.channel(connection, channel_id, create)
            # Выводим информационное сообщение
            @info "=> Канал успешно создан!"
            # Возвращаем результат
            return result
        catch error
            # Обрабатываем ошибки: выводим в терминал и кидаем исключение
            @error Utils.Chronometer.message_with_time("Ошибка при создании канала") error
            throw(error)
        end
    end

    """Функция для получения параметров авторизации на RabbitMQ.

    login: логин пользователя на сервисе очередей,
    password: пароль от пользователя на сервисе очередей
    """
    function get_auth_params(login::String, password::String, mechanism::String)::Dict{String, Any}
        # Выводим информационное сообщение
        @info Utils.Chronometer.message_with_time("Вызвана функция сборки параметров соединеня с сервисом очередей...")
        # Логирование для режима отладки
        @debug "Параметры вызова функции:"
        @debug "> имя пользователя на сервисе очередей" login
        @debug "> пароль пользователя на сервисе очередей" password
        @debug "> механизм соединеня на сервисе очередей" mechanism
        # Создаём словарь с параметрами авторизации на сервисе очередей
        # Важно! Тип данного словаря должен быть именно {String, Any}!
        # При смене типа словаря параметров соединене с сервером отказывается
        # устаналиваться.
        #Такой тип задан в файле:
        # * https://github.com/JuliaComputing/AMQPClient.jl/blob/master/src/auth.jl
        result::Dict{String,Any} = Dict{String, Any}(
            "MECHANISM" => mechanism,
            "LOGIN" => login,
            "PASSWORD" => password
        )
        @debug "=> Результат сборки параметров авторизации" result
        # Возвращаем результат (словарь настроек)
        return result
    end

    """Функция для получения установки соединения с сервером RabbitMQ.

    vhost: виртуальный хост на сервисе очередей,
    host: хост для подключенияя к сервису очередей,
    port: порт для подключения к сервису очередей,
    auth_params: параметры авторизации на сервисе очередей
    """
    function get_connection(vhost::String, host::String, port::Int64, auth_params::Dict{String, Any})::AMQPClient.MessageChannel
        #
        @info Utils.Chronometer.message_with_time("Вызвана функция соединения с сервером очередей RabbitMQ...")
        # Логе для режима отладки
        @debug "Параметры вызова функции:"
        @debug "> виртуальный хост" vhost
        @debug "> хост сервера RabbitMQ" host
        @debug "> порт сервера RabbitMQ" port
        @debug "> параметры авторизации на сервере RabbitMQ" auth_params
        # Получаем соединение с сервером очередей RabbitMQ
        try
            # Выполняем операцию соединения
            result::AMQPClient.MessageChannel = AMQPClient.connection(
                ;virtualhost=vhost,
                host=host,
                port=port,
                auth_params=auth_params
            )
            # Выводим информационное сообщение
            @info "=> Соединене успешно установлено!"
            # Возвращаем результат
            return result
        catch error
            # Обрабатываем ошибки: выводим в терминал и кидаем исключение
            @error Utils.Chronometer.message_with_time("Не удалось установить соединене с RabbitMQ. Ошибка") error
            throw(error)
        end
    end
end


text:
using Base.Filesystem
lnname ="oho.ln"
ln=Base.Filesystem.readlink(lnname)
println(ln)
lnisdir = Base.Filesystem.isdir( lnname )
tgrisdir = Base.Filesystem.isdir( ln)
println( lnisdir,tgrisdir)

println( "-----",Base.Filesystem.readlink(lnname) )

lnisdir = Base.Filesystem.islink( lnname )
tgrisdir = Base.Filesystem.islink( ln)
println( lnisdir,tgrisdir)

lnisdir = Base.Filesystem.ispath( lnname )
tgrisdir = Base.Filesystem.ispath( ln)
println( lnisdir,tgrisdir)

println(Base.Filesystem.homedir() )

mtime= Base.Filesystem.mtime(lnname)

println(mtime)

using Dates
a=Dates.unix2datetime( mtime)
println(a)


text:
# provided by https://github.com/JuliaLang/julia/pull/35094
function register_hints()
    if isdefined(Base, :Experimental) && isdefined(Base.Experimental, :register_error_hint)
        Base.Experimental.register_error_hint(MethodError) do io, exc, argtypes, kwargs
            if exc.f in (zero, one) && argtypes[1] <: Union{Type{<:AbstractRGB}, AbstractRGB}
                print(io, "\nYou may need to `using ColorVectorSpace`.")
            end
            if exc.f in (zeros, ones) && argtypes[1] <: Type{<:AbstractRGB}
                print(io, "\nYou may need to `using ColorVectorSpace`.")
            end
        end
    end
end


text:
println(1)

function add(lv, rv)
    return lv + rv
end

println("add(1, 3) = ", add(1, 3))

max_stars_repo_name: dannypsnl/languages-learn

max_issues_repo_name: dannypsnl/languages-learn

max_forks_repo_name: dannypsnl/languages-learn


text:
__precompile__()

module ContMechTensors

import Base.@pure
using StaticArrays

export AbstractTensor, SymmetricTensor, Tensor, Vec, FourthOrderTensor, SecondOrderTensor

export otimes, ⊗, ⊡, dcontract, dev, vol, symmetric, skew, minorsymmetric, majorsymmetric
export minortranspose, majortranspose, isminorsymmetric, ismajorsymmetric
export tdot, dott, dotdot
export hessian#, gradient
export basevec, eᵢ
export rotate

@deprecate extract_components(tensor) Array(tensor)

#########
# Types #
#########
abstract AbstractTensor{order, dim, T <: Real} <: AbstractArray{T, order}

immutable SymmetricTensor{order, dim, T <: Real, M} <: AbstractTensor{order, dim, T}
   data::SVector{M, T}
end

immutable Tensor{order, dim, T <: Real, M} <: AbstractTensor{order, dim, T}
   data::SVector{M, T}
end

###############
# Typealiases #
###############
typealias Vec{dim, T, M} Tensor{1, dim, T, dim}

typealias AllTensors{dim, T} Union{SymmetricTensor{2, dim, T}, Tensor{2, dim, T},
                                   SymmetricTensor{4, dim, T}, Tensor{4, dim, T},
                                   Vec{dim, T}}

typealias SecondOrderTensor{dim, T} Union{SymmetricTensor{2, dim, T}, Tensor{2, dim, T}}
typealias FourthOrderTensor{dim, T} Union{SymmetricTensor{4, dim, T}, Tensor{4, dim, T}}
typealias SymmetricTensors{dim, T} Union{SymmetricTensor{2, dim, T}, SymmetricTensor{4, dim, T}}
typealias Tensors{dim, T} Union{Tensor{2, dim, T}, Tensor{4, dim, T}, Vec{dim, T}}

include("utilities.jl")
include("indexing.jl")
include("promotion_conversion.jl")
include("tensor_products.jl")
include("transpose.jl")
include("symmetric.jl")
include("math_ops.jl")
include("special_ops.jl")
include("tensor_ops_errors.jl")
include("automatic_differentiation.jl")

##############################
# Utility/Accessor Functions #
##############################
get_data(t::AbstractTensor) = t.data.data
tovector(t::AbstractTensor) = t.data

function tomatrix{dim}(t::Tensor{4, dim})
    N = n_components(Tensor{2,dim})
    return SMatrix{N, N}(tovector(t))
end

function tomatrix{dim}(t::Tensor{2, dim})
    N = n_components(Tensor{1,dim})
    return SMatrix{N, N}(tovector(t))
end

@pure n_components{dim}(::Type{SymmetricTensor{2, dim}}) = dim*dim - div((dim-1)*dim, 2)
@pure function n_components{dim}(::Type{SymmetricTensor{4, dim}})
    n = n_components(SymmetricTensor{2, dim})
    return n*n
end
@pure n_components{order, dim}(::Type{Tensor{order, dim}}) = dim^order

@pure get_base{order, dim, T, M}(::Type{SymmetricTensor{order, dim, T, M}}) = SymmetricTensor{order, dim}
@pure get_base{order, dim, T}(::Type{SymmetricTensor{order, dim, T}}) = SymmetricTensor{order, dim}
@pure get_base{order, dim}(::Type{SymmetricTensor{order, dim}}) = SymmetricTensor{order, dim}
@pure get_base{order, dim, T, M}(::Type{Tensor{order, dim, T, M}}) = Tensor{order, dim}
@pure get_base{order, dim, T}(::Type{Tensor{order, dim, T}}) = Tensor{order, dim}
@pure get_base{order, dim}(::Type{Tensor{order, dim}}) = Tensor{order, dim}

############################
# Abstract Array interface #
############################
Base.linearindexing{T <: SymmetricTensor}(::Type{T}) = Base.LinearSlow()
Base.linearindexing{T <: Tensor}(::Type{T}) = Base.LinearFast()

get_type{X}(::Type{Type{X}}) = X

########
# Size #
########
Base.size{dim}(::Vec{dim}) = (dim,)
Base.size{dim}(::SecondOrderTensor{dim}) = (dim, dim)
Base.size{dim}(::FourthOrderTensor{dim}) = (dim, dim, dim, dim)

#########################
# Internal constructors #
#########################
function _constructor_check{order, dim}(Tt::Union{Type{Tensor{order, dim}}, Type{SymmetricTensor{order, dim}}})
    !(order in (1,2,4)) && throw(ArgumentError("tensors only supported for order 1, 2 and 4"))
    !(dim in (1,2,3)) && throw(ArgumentError("tensors only supported for dim 1, 2 and 3"))
    order == 1 && Tt <: SymmetricTensor && throw(ArgumentError("symmetric tensors only supported for order 2 and 4"))
    return nothing
end

for TensorType in (SymmetricTensor, Tensor)
    @eval begin
        @inline (::Type{$TensorType{order, dim}}){order, dim, T <: Real, M}(t::NTuple{M, T}) = $TensorType{order, dim}(SVector{M,T}(t))
        @inline (::Type{$TensorType{order, dim, T1}}){order, dim, T1 <: Real, T2 <: Real, M}(t::NTuple{M, T2}) = $TensorType{order, dim}(SVector{M,T1}(t))
        @inline (::Type{$TensorType{order, dim, T1, M}}){order, dim, T1 <: Real, T2 <: Real, M}(t::NTuple{M, T2}) = $TensorType{order, dim}(SVector{M,T1}(t))

        @inline (::Type{$TensorType{order, dim}}){order, dim, T <: Real, M}(t::SVector{M, T}) = _constructor_check($TensorType{order, dim})
    end
    for order in (2,4), dim in (1,2,3)
        M = n_components(TensorType{order,dim})
        @eval @inline (::Type{$TensorType{$order, $dim}}){T <: Real}(t::SVector{$M, T}) = $TensorType{$order, $dim, T, $M}(t)
        @eval @inline (::Type{$TensorType{$order, $dim}}){T <: Real, Q}(t::SMatrix{Q, Q, T, $M}) = $TensorType{$order, $dim, T, $M}(t.data)
    end
    if TensorType == Tensor
        for dim in (1,2,3)
            @eval @inline (::Type{Tensor{1, $dim}}){T <: Real}(t::SVector{$dim, T}) = Tensor{1, $dim, T, $dim}(t)
        end
    end
end
# Special for Vec
@inline (Tt::Type{Vec{dim}}){dim}(data) = Tensor{1, dim}(data)

# General fallbacks
@inline          (Tt::Type{Tensor{order, dim, T}}){order, dim, T}(data::Union{AbstractArray, Tuple, Function}) = convert(Tensor{order, dim, T}, Tensor{order, dim}(data))
@inline (Tt::Type{SymmetricTensor{order, dim, T}}){order, dim, T}(data::Union{AbstractArray, Tuple, Function}) = convert(SymmetricTensor{order, dim, T}, SymmetricTensor{order, dim}(data))
#@inline          (Tt::Type{Tensor{order, dim, T, M}}){order, dim, T, M}(data::Union{AbstractArray, Tuple, Function}) = Tensor{order, dim, T}(data)
#@inline (Tt::Type{SymmetricTensor{order, dim, T, M}}){order, dim, T, M}(data::Union{AbstractArray, Tuple, Function}) = SymmetricTensor{order, dim, T}(data)

###############
# Simple Math #
###############

# op(Number, Tensor): *, /
for TensorType in (SymmetricTensor, Tensor)
    @eval begin
        @inline Base.:*{order, dim, T, N}(n::Number, t::$TensorType{order, dim, T, N}) = $TensorType{order, dim}(n * tovector(t))
        @inline Base.:*{order, dim, T, N}(t::$TensorType{order, dim, T, N}, n::Number) = $TensorType{order, dim}(tovector(t) * n)
        @inline Base.:/{order, dim, T, N}(t::$TensorType{order, dim, T, N}, n::Number) = $TensorType{order, dim}(tovector(t) / n)

        # Unary -, +
        @inline Base.:-{order, dim, T, N}(t::$TensorType{order, dim, T, N}) = $TensorType{order, dim}(-tovector(t))
        @inline Base.:+{order, dim, T, N}(t::$TensorType{order, dim, T, N}) = $TensorType{order, dim}(+tovector(t))
    end
end

# Binary op(Tensor, Tensor): +, -, .+, .-, .*, ./
for op in (:+, :-, :.+, :.-, :.*, :./)
    for TensorType in (SymmetricTensor, Tensor)
        @eval begin
            @inline function Base.$op{order, dim, T1, T2, N}(t1::$TensorType{order, dim, T1, N}, t2::$TensorType{order, dim, T2, N})
                $TensorType{order, dim}($op(tovector(t1), tovector(t2)))
            end
        end
    end
    @eval begin
        Base.$op{order, dim}(t1::AbstractTensor{order, dim}, t2::AbstractTensor{order, dim}) = $op(promote(t1, t2)...)
    end
end

######################
# Basic constructors #
######################

# zero, rand, ones
for op in (:zero, :rand, :ones)
    for TensorType in (SymmetricTensor, Tensor)
        @eval begin
            @inline Base.$op{order, dim}(Tt::Type{$TensorType{order, dim}}) = $op($TensorType{order, dim, Float64})
            @inline Base.$op{order, dim, T, M}(Tt::Type{$TensorType{order, dim, T, M}}) = $op($TensorType{order, dim, T})
            @inline function Base.$op{order, dim, T}(Tt::Type{$TensorType{order, dim, T}})
                N = n_components($TensorType{order, dim})
                return $TensorType{order, dim}($op(SVector{N, T}))
            end
        end
    end
    # Special case for Vec
    @eval @inline Base.$op{dim}(Tt::Type{Vec{dim}}) = $op(Vec{dim, Float64})

    # zero, rand or ones of a tensor
    @eval @inline Base.$op(t::AllTensors) = $op(typeof(t))
end

# zeros, ones
for (op, el) in ((:zeros, :zero), (:ones, :one))
    for TensorType in (SymmetricTensor, Tensor)
        @eval begin
            @inline Base.$op{order, dim}(Tt::Type{$TensorType{order, dim}}, dims...) = $op($TensorType{order, dim, Float64}, dims...)
            @inline function Base.$op{order, dim, T}(Tt::Type{$TensorType{order, dim, T}}, dims...)
                N = n_components($TensorType{order, dim})
                return $op($TensorType{order, dim, T, N}, dims...)
            end
            @inline Base.$op{order, dim, T, M}(Tt::Type{$TensorType{order, dim, T, M}}, dims...) =
                fill!(Array{$TensorType{order, dim, T, M}}(dims...), $el($TensorType{order, dim, T}))
        end
    end
    @eval @inline Base.$op{dim}(::Type{Vec{dim}}, dims...) = $op(Vec{dim, Float64}, dims...)
end

# diagm
for TensorType in (SymmetricTensor, Tensor)
    @eval begin
        @generated function Base.diagm{dim}(Tt::Type{$(TensorType){2, dim}}, v::Union{AbstractVector, Tuple})
            f = (i,j) -> i == j ? :(v[$i]) : :($(zero(eltype(v))))
            exp = tensor_create(get_type(Tt), f)
            return quote
                $(Expr(:meta, :inline))
                @inbounds t = $exp
                $($TensorType){2, dim}(t)
            end
        end

        @inline Base.diagm{dim, T}(Tt::Type{$(TensorType){2, dim}}, v::T) = v * one($(TensorType){2, dim, T})
    end
end

# one (identity tensor)
for TensorType in (SymmetricTensor, Tensor)
    @eval begin
        @inline Base.one{order, dim}(Tt::Type{$(TensorType){order, dim}}) = one($TensorType{order, dim, Float64})
        @inline Base.one{order, dim, T, M}(Tt::Type{$(TensorType){order, dim, T, M}}) = one($TensorType{order, dim, T})
        @inline Base.one{order, dim, T}(Tt::$TensorType{order, dim, T}) = one($TensorType{order, dim, T})

        @generated function Base.one{order, dim, T}(Tt::Type{$(TensorType){order, dim, T}})
            !(order in (2,4)) && throw(ArgumentError("`one` only defined for order 2 and 4"))
            δ = (i,j) -> i == j ? :($(one(T))) : :($(zero(T)))
            if order == 2
                f = (i,j) -> δ(i,j)
            elseif order == 4 && $TensorType == Tensor
                f = (i,j,k,l) -> δ(i,k) * δ(j,l)
            else # order == 4 && TensorType == SymmetricTensor
                f = (i,j,k,l) -> (δ(i,k) * δ(j,l) + δ(i,l) * δ(j,k)) / 2
            end
            exp = tensor_create(get_base(get_type(Tt)), f)
            return quote
                $(Expr(:meta, :inline))
                $($TensorType){order, dim}($exp)
            end
        end
    end
end

# Tensor from function
@generated function (Tt::Union{Type{Tensor{order, dim}}, Type{SymmetricTensor{order, dim}}}){order, dim}(f::Function)
    _constructor_check(get_base(get_type(Tt)))
    if order == 1
        exp = tensor_create(get_base(get_type(Tt)), (i) -> :(f($i)))
    elseif order == 2
        exp = tensor_create(get_base(get_type(Tt)), (i,j) -> :(f($i, $j)))
    elseif order == 4
        exp = tensor_create(get_base(get_type(Tt)), (i,j,k,l) -> :(f($i, $j, $k, $l)))
    end

    return :(get_base(Tt)($exp))
end
Base.convert{dim, T}(::Type{SVector{dim, T}}, f::Function) = SVector{dim, T}(ntuple(f, Val{dim}))

# Tensor from AbstractArray
@generated function (Tt::Type{Tensor{order, dim}}){order, dim}(data::AbstractArray)
    _constructor_check(get_base(get_type(Tt)))
    n = n_components(Tensor{order,dim})
    return quote
        if length(data) != $n
            throw(ArgumentError("wrong number of vector elements, expected $($n), got $(length(data))"))
        end
        Tensor{order, dim}(SVector{$n}(data))
    end
end

# SymmetricTensor from AbstractArray
@generated function (Tt::Type{SymmetricTensor{order, dim}}){order, dim}(data::AbstractArray)
    _constructor_check(get_base(get_type(Tt)))
    n = n_components(Tensor{order,dim})
    m = n_components(SymmetricTensor{order,dim})
    return quote
        if length(data) != $n && length(data) != $m
            throw(ArgumentError("wrong number of vector elements, expected $($n) or $($m), got $(length(data))"))
        end
        if length(data) == $m
            return SymmetricTensor{order, dim}(SVector{$m}(data))
        end
        S = Tensor{order, dim}(SVector{$n}(data))
        return convert(SymmetricTensor{order, dim}, S)
    end
end

"""
```julia
basevec(::Type{Vec{dim, T}})
basevec(::Type{Vec{dim, T}}, i)
basevec(::Vec{dim, T})
basevec(::Vec{dim, T}, i)
```
Return a tuple with the base vectors corresponding to the dimension `dim` and type
`T`. An optional integer `i` can be used to extract the i:th base vector.
The alias `eᵢ` can also be used, written `e\\_i<TAB>`.

**Example:**

```jldoctest
julia> eᵢ(Vec{2, Float64})
([1.0,0.0],[0.0,1.0])

julia> eᵢ(Vec{2, Float64}, 2)
2-element ContMechTensors.Tensor{1,2,Float64,2}:
 0.0
 1.0
```
"""
@inline function basevec{T}(::Type{Vec{1, T}})
    o = one(T)
    return (Vec{1, T}((o,)), )
end
@inline function basevec{T}(::Type{Vec{2, T}})
    o = one(T)
    z = zero(T)
    return (Vec{2, T}((o, z)),
            Vec{2, T}((z, o)))
end
@inline function basevec{T}(::Type{Vec{3, T}})
    o = one(T)
    z = zero(T)
    return (Vec{3, T}((o, z, z)),
            Vec{3, T}((z, o, z)),
            Vec{3, T}((z, z, o)))
end

@inline basevec{dim}(::Type{Vec{dim}}) = basevec(Vec{dim, Float64})
@inline basevec{dim, T}(::Type{Vec{dim, T}}, i::Int) = basevec(Vec{dim, T})[i]
@inline basevec{dim}(::Type{Vec{dim}}, i::Int) = basevec(Vec{dim, Float64})[i]
@inline basevec{dim, T}(v::Vec{dim, T}) = basevec(typeof(v))
@inline basevec{dim, T}(v::Vec{dim, T}, i::Int) = basevec(typeof(v), i)

const eᵢ = basevec

end # module


text:
get_arch_os_libc(target::AbstractString) = get_arch(target), get_os(target), get_libc(target)

get_arch(target::AbstractString) = first(split(target, '-'))

function get_os(target::AbstractString)
    _, vendor, _ = split(target, '-')
    if vendor == "apple"
        os = "macos"
    elseif vendor == "w64"
        os = "windows"
    elseif vendor == "unknown"
        os = "freebsd"
    elseif vendor == "linux"
        os = "linux"
    else
        error("Unknown OS: $target")
    end
    return os
end

function get_libc(target::AbstractString)
    _, _, env = split(target, '-')
    if startswith(env, "gnu")
        return "glibc"
    elseif startswith(env, "musl")
        return "musl"
    else
        return ""
    end
end


text:
module Keywords

export Keyword

import Base: ==, hash
using ..Lists
import ..Lists.unparse

struct Keyword
    sym::String
end
x::Keyword == y::Keyword = x.sym == y.sym
hash(x::Keyword, n::UInt) = hash(Keyword, hash(x.sym, n))

unparse(kw::Keyword) = string(kw)
Base.show(io::IO, kw::Keyword) = print(io, "#:", kw.sym)

end


text:
using BenchmarkTools, Random
using Mill, Flux

using StatsPlots, LaTeXStrings; pgfplotsx()

Random.seed!(42)
params = Base.Iterators.product([10, 100, 1000], [100, 200, 300]) |> collect |> vec
for (sl,sn) in params
    @show sn, sl
    # A = rand(100, 2053)
    A = PostImputingMatrix(rand(100, 2053))
    B = NGramMatrix([rand() < 0.5 ? missing : randstring(sl) for _ in 1:sn], 3, 256, 2053)
    # B = NGramMatrix(fill(missing, sn), 3, 256, 2053)
    # B = NGramMatrix([randstring(sl) for _ in 1:sn], 3, 256, 2053)
    @btime *($A, $B);
end

xticklabels = [latexstring("$a|$b") for (a,b) in params]

t1 = [56.427, 425.149, 4.043*1000, 112.646, 850.431, 8.274*1000, 169.531, 1.274*1000, 12.801*1000]
l1 = "Mill 1.5.2: Matrix * NGramMatrix{String, …}"
t2 = [55.109, 440.126, 4*1000, 116.665, 860.89, 8.186*1000, 183.248, 1.294*1000, 13.093*1000]
l2 = "Mill 2.0: Matrix * NGramMatrix{String, …}"
t4 = [38.523, 38.552, 38.452, 38.420, 38.650, 38.565, 37.717, 38.652, 38.789]
l4 = "Mill 2.0: PostImputingMatrix * NGramMatrix{Missing, …}"
t3 = [37.686, 256.430, 2179, 87.778, 450.185, 4497, 110.967, 773.528, 6761]
l3 = "Mill 2.0: PostImputingMatrix * NGramMatrix{Union{Missing, String}, …}"
groupedbar(xticklabels, hcat(t1, t2, t3, t4); ylabel=L"\mu s", xlabel="n. of strings | string length", labels=[l1 l2 l3 l4], title="PostImputing NGramMatrix multiplication")

savefig("ngram_mult.png")

Random.seed!(42)
params = Base.Iterators.product([10, 100, 1000], [100, 200, 300]) |> collect |> vec
for (sl,sn) in params
    @show sn, sl
    # A = rand(100, 2053)
    A = PostImputingMatrix(rand(100, 2053))
    # B = NGramMatrix([rand() < 0.5 ? missing : randstring(sl) for _ in 1:sn], 3, 256, 2053)
    B = NGramMatrix(fill(missing, sn), 3, 256, 2053)
    # B = NGramMatrix([randstring(sl) for _ in 1:sn], 3, 256, 2053)
    @btime gradient((A,B) -> sum(A*B), $A, $B)
end

xticklabels = [latexstring("$a|$b") for (a,b) in params]

t1 = [280.452, 1104, 8696, 442.584, 1848, 16047, 541.064, 2754, 25244]
l1 = "Mill 1.5.2: Matrix * NGramMatrix{String, …}"
t2 = [287.431, 1108, 9185, 421.466, 2061, 15871, 560.260, 2944, 28016]
l2 = "Mill 2.0: Matrix * NGramMatrix{String, …}"
t4 = [28.615, 28.226, 29.704, 52.057, 49.217, 54.538, 82.527, 77.738, 72.033]
l4 = "Mill 2.0: PostImputingMatrix * NGramMatrix{Missing, …}"
t3 = [201.775, 843.369, 4098, 370.220, 1302, 10952, 506.579, 1668, 16352]
l3 = "Mill 2.0: PostImputingMatrix * NGramMatrix{Union{Missing, String}, …}"
groupedbar(xticklabels, hcat(t1, t2, t3, t4); ylabel=L"\mu s", xlabel="n. of strings | string length", labels=[l1 l2 l3 l4], title="PostImputing NGramMatrix multiplication (grad)")

savefig("ngram_mult_grad.png")


text:
# Static Discrete Choice exercise
function sdc(N, numSch;
                     γ = 1.2,
                     σ = sqrt((pi^2)/6.0),
                     μ = 0,
                     seed = 0)

# Generating the data
# Main parameters:
Random.seed!(seed+1);
δ    = randn(numSch,1) .+ μ; #generating the parameter deltas (mean utility of each school)
dist = abs.(randn(N,numSch)); #generating the observable distance

# Draw the utilities
Random.seed!(seed+2);
ϵ=σ*randn(N,numSch);
util = repeat(δ',N,1) .- γ*dist .+ ϵ;

#Let's introduce a variation in the menu
menuSize = zeros(Int64,N,1);
choice   = zeros(Int64,N,1);
menu     = zeros(Bool, N,numSch);

#Loop
Random.seed!(seed+3);
for i=1:N
    #Number of schools in the menu
    menuSize[i]=rand(DiscreteUniform(1,numSch-1)) + 1;
    # Generate a random integer uniformaly distributed between 1 and (numSch-1)
    # I add 1 to assure that they will have at least 2 options

    #Identity of schools in the menu
    temp=randperm(numSch);
    # return a random permutation of integers between 1 and NumSch

    for j=1:numSch
        if any(temp[1:menuSize[i]].==j)
            # any() returns equal 1 if any element of the array is non-zero
            # 1:menuSize[i] is the menu of individual i
            # Thus I will put one in the matrix menu when the i have the
            # option j available
            menu[i,j]=true;
        end
    end

    #Choose the school that provides highest utility
    choice[i]=findall(util[i,:].== findmax(util[i,menu[i,:]])[1])[1]; # I am assuming no ties!
    # I am returning the position of the option with the highest utility
    # that is available in the menu

end

# Let's compute the likelihood funtion
# We will assume ϵ∼T1EV(0,1)
function dcobj_t1ev(x,choice,dist,N,numSch,menu)
    δ = [x[1:numSch-1];0.0];
    γ = x[numSch];
    v = exp.(δ'.-γ*dist).*menu;
    PredProb = v ./ repeat(sum(v,dims=2),1,numSch);
    chocProb = [PredProb[i,choice[i]] for i in 1:N];
    objVal   =-sum(log, chocProb);

    return objVal
end

#Now lets assume that  ϵ∼N(0,π^2/6)
#In this case we have to simulate the choices in order to get the likelihood function.
function dcobj_normal(x,choice,dist,N,numSch, menu,seed)
    δ = [x[1:numSch-1];0.0];
    γ = x[numSch];
    σ = sqrt((pi^2)/6);

    # Draw the epsilon to simulate the choice probabilities
    Random.seed!(seed+4);
    numDraws=numSch*20;
    ϵsim = σ*randn(numDraws,numSch, N);

    # Compute the simulated utilities and the choice for each draws
    # Note that here we need only the probability of the best option, the other
    # probabilities do not matter here
    simChoice = zeros(Int64,numDraws,1);
    simProb   = zeros(N,1);

    for i=1:N
        simUtil = repeat(δ',numDraws,1) - γ*repeat(dist[i,:]',numDraws,1) + ϵsim[:,:,i];
        for k=1:numDraws
            simChoice[k]=findall(simUtil[k,:].== findmax(simUtil[k,menu[i,:]])[1])[1];
        end
        simProb[i] = sum(simChoice.==choice[i])/numDraws;
    end

    objVal   =-sum(log, simProb);
    return objVal
end

# Let's call the optimizer
xtrue =     [δ[1:numSch-1].-δ[numSch];γ];
xinit = 0.9*[δ[1:numSch-1].-δ[numSch];γ];
#dcobj_t1ev(xinit,choice,dist,N,numSch,menu)
#dcobj_normal(xinit,choice,dist,N,numSch,menu,seed)

out_t1ev   = optimize(x -> dcobj_t1ev(x,choice,dist,N,numSch, menu) , xinit)
out_normal = optimize(x -> dcobj_normal(x,choice,dist,N,numSch,menu,seed), xinit)

println("ϵ-T1EV:    minimum = $(out_t1ev.minimum) with argmin = $(out_t1ev.minimizer) in $(out_t1ev.iterations) iterations")
println("ϵ-Normal:  minimum = $(out_normal.minimum) with argmin = $(out_normal.minimizer) in $(out_normal.iterations) iterations")
println("θtrue = $(xtrue), θhat_t1ev = $(out_t1ev.minimizer), θhat_normal = $(out_normal.minimizer)")

# Organizing the outputs
comparison = (θtrue = xtrue, θhat_t1ev = out_t1ev.minimizer, θhat_normal = out_normal.minimizer);
data = ( δ=δ, dist=dist, ϵ=ϵ, choice=choice, menu=menu);
params = (γ = γ, σ = σ);
output = (comparison = comparison, data = data, params=params);

return output;
end


text:
include("components.jl")
using DataFrames

type SPNdata
  unit_size::Int
  num_unit::Int
  plaintexts::Array
  ciphertexts::Array
  traces::Array
  res_dir::String
end

type TA_TRAINdata
  unit_size::Int
  num_unit::Int
  plaintexts::Array
  ciphertexts::Array
  keys::Array
  traces::Array
  res_dir::String
end

type TA_TESTdata
  unit_size::Int
  num_unit::Int
  plaintexts::Array
  ciphertexts::Array
  traces::Array
  res_dir::String
end

function loadSPNtrace(path_to_csvtraces, unit_size)
  #creating directory to store the results
  part_names = split(path_to_csvtraces,"/")
  path_to_resdir = ""
  for a in part_names[1:end-1]
    path_to_resdir = join([path_to_resdir, a, "/"])
  end
  path_to_resdir = join([path_to_resdir, part_names[end], "results"])
  mkdir(path_to_resdir)

  pct = readtable(path_to_csvtraces,header=false)

  plaintexts = pct[1]; #all the plaintexts are suppose to be of the same length and in hex string
  ciphertexts = pct[2]; #all the ciphertexts are suppose to be of the same length and in hex string
  traces = Array(pct[:,3:end]);
  num_unit = div(length(plaintexts[1])*4, unit_size)

  #the HW creation
  for i in 0:2^unit_size-1
    HW[i] = hw_look_up_maker(bin(i,unit_size))
  end

  x = SPNdata(unit_size,num_unit,plaintexts,ciphertexts,traces,path_to_resdir)
end
#add more loading function for other types of cipher if needed

#this function will return the info (stored in array types) of the total number
#of templates that have to be built before classification of any new trace
function loadTA_TRAINtrace(path_to_csvtraces, unit_size)
  #creating the director to store the trained templates for later testing
  part_names = split(path_to_csvtraces,"/")
  path_to_resdir = ""
  for a in part_names[1:end-1]
    path_to_resdir = join([path_to_resdir, a, "/"])
  end
  path_to_resdir = join([path_to_resdir, part_names[end], "templates"])
  mkdir(path_to_resdir)

  pct = readtable(path_to_csvtraces,header=false)

  #HW creation may be needed in some kind of templates
  for i in 0:2^unit_size-1
    HW[i] = hw_look_up_maker(bin(i,unit_size))
  end

  plaintexts = pct[1]
  ciphertexts = pct[2]
  keys = pct[3]
  traces = Array(pct[:,4:end])
  println("dbg: plaintexts[1] is: ",plaintexts[1])
  println("dbg: length(plaintexts[1)] is: ",length(plaintexts[1]))
  num_unit = div(length(plaintexts[1])*4, unit_size)

  return TA_TRAINdata(unit_size, num_unit, plaintexts, ciphertexts, keys, traces, path_to_resdir)
end

function loadTA_TESTtrace(path_to_csvtraces, unit_size)
  #for the testing case the traces need not have plaintexts and ciphertexts
  #for some template models, In such casses user have to ensure that there are
  #plaintexts and ciphertexts cols are filled with 000...0 with apporpriate number of
  #zeros in so that num_unit can be calculated

  #creating the directory to store results
  part_names = split(path_to_csvtraces,"/")
  path_to_resdir = ""
  for a in part_names[1:end-1]
    path_to_resdir = join([path_to_resdir, a, "/"])
  end
  path_to_resdir = join([path_to_resdir, part_names[end], "results"])
  mkdir(path_to_resdir)

  pct = readtable(path_to_csvtraces,header=false)

  plaintexts = pct[1]; #all the plaintexts are suppose to be of the same length and in hex string
  ciphertexts = pct[2]; #all the ciphertexts are suppose to be of the same length and in hex string
  traces = Array(pct[:,3:end]);
  num_unit = div(length(plaintexts[1])*4, unit_size)

  #the HW creation
  for i in 0:2^unit_size-1
    HW[i] = hw_look_up_maker(bin(i,unit_size))
  end

  return TA_TESTdata(unit_size,num_unit,plaintexts,ciphertexts,traces,path_to_resdir)
end

max_stars_repo_name: udion/zypher

max_issues_repo_name: udion/zypher

max_forks_repo_name: udion/zypher


text:
st_affine(ga::GeoArray) = ga.f

function get_affine_map(ds::ArchGDAL.Dataset)
    # ArchGDAL fails hard on datasets without
    # an affinemap. GDAL documents that on fail
    # a default affinemap should be returned.
    try
        global gt = ArchGDAL.getgeotransform(ds)
    catch y
        @warn y.msg
        global gt = [0.,1.,0.,0.,0.,1.]
    end
    geotransform_to_affine(SVector{6, Float64}(gt))
end

function geotransform_to_affine(gt::SVector{6, Float64})
    # See https://lists.osgeo.org/pipermail/gdal-dev/2011-July/029449.html
    # for an explanation of the geotransform format
    AffineMap(SMatrix{2,2}([gt[2] gt[3]; gt[5] gt[6]]), SVector{2}([gt[1], gt[4]]))
end
geotransform_to_affine(A::Vector{Float64}) = geotransform_to_affine(SVector{6}(A))

function affine_to_geotransform(am::AffineMap{SArray{Tuple{2,2},Float64,2,4},SArray{Tuple{2},Float64,1,2}})
    l = am.linear
    t = am.translation
    [t[1], l[1], l[3], t[2], l[2], l[4]]
end

function affine_to_geotransform(am::AffineMap{Array{Float64,2},Array{Float64,1}})
    l = am.linear
    t = am.translation
    (length(l) != 4 || length(t) != 2) || error("AffineMap has wrong dimensions.")
    [t[1], l[1], l[3], t[2], l[2], l[4]]
end

"""Check wether the AffineMap of a GeoArray contains rotations."""
function is_rotated(ga::GeoArray)
    ga.f.linear[2] != 0. || ga.f.linear[3] != 0.
end

function unitrange_to_affine(x::StepRangeLen, y::StepRangeLen)
    δx, δy = step(x), step(y)
    AffineMap(
        SMatrix{2,2}(δx, 0, 0, δy),
        SVector(x[1] - δx/2, y[1] - δy/2)
    )
end

# Extend CoordinateTransformations
CoordinateTransformations.compose(ga::GeoArray, t2::AffineMap) = CoordinateTransformations.compose(ga.f, t2)
CoordinateTransformations.compose(ga::GeoArray, t2::LinearMap) = CoordinateTransformations.compose(ga.f, t2)
CoordinateTransformations.compose(ga::GeoArray, t2::Translation) = CoordinateTransformations.compose(ga.f, t2)

"""Transform an GeoArray by applying a Transformation."""
function compose!(ga::GeoArray, t2::AffineMap)
    ga.f = compose(ga, t2)
end
function compose!(ga::GeoArray, t2::LinearMap)
    ga.f = compose(ga, t2)
end
function compose!(ga::GeoArray, t2::Translation)
    ga.f = compose(ga, t2)
end

export st_affine


text:
export getfirstnotes, filterpitches, separatepitches, firstnotes, removepitches
export estimate_delay, estimate_delay_recursive

"""
    firstnotes(notes, grid)

Return the notes that first appear in each grid point, *without quantizing them*.

This function *does not* consider the notes modulo the quarter note! Different
quarter notes have different grid points.
"""
function firstnotes(notes::Notes{N}, grid) where {N<:AbstractNote}
    isgrid(grid)

    clas = classify(notes, grid)

    # Notes with dif == 0 belong to the same grid bin
    dif = clas[2:end] - clas[1:end-1]

    posdif = [Int(notes[j].position - notes[j-1].position) for j in 2:length(notes)]
    # minpos is the minimum position difference notes can have while being in the
    # same quarter note
    tpq = notes.tpq
    minpos = tpq - tpq*(grid[end] - grid[end-1] + grid[2] - grid[1])/2

    # The first notes are notes that have dif == 0 OR posdif > minpos
    toadd = Vector{Note}()
    for j in 1:length(notes)-1
        if (dif[j] == 0 && posdif[j] > minpos) || dif[j] != 0
            push!(toadd, notes[j+1])
        end
    end

    # Take care of first note in notes:
    if clas[1] != clas[2] || notes[2].position - notes[1].position > minpos
        pushfirst!(toadd, notes[1])
    elseif notes[2].position < notes[1].position
        toadd[1] = notes[1]
    end
    return Notes(toadd, notes.tpq)
end

"""
    filterpitches(notes::Notes, filters) -> newnotes

Only keep the notes that have a pitch specified in `filters` (one or many pitches).
"""
function filterpitches(notes::Notes{N}, filters) where {N<:AbstractNote}
    n = N[]
    for i ∈ 1:length(notes)
        notes[i].pitch ∈ filters && push!(n, copy(notes[i]))
    end
    return Notes(n, notes.tpq)
end

"""
    removepitches(notes::Notes, remove) -> newnotes

Remove notes that have a pitch specified in `remove`.
"""
function removepitches(notes::Notes{N}, remove) where {N}
    n = N[]
    for i ∈ 1:length(notes)
        notes[i].pitch ∉ remove && push!(n, copy(notes[i]))
    end
    return Notes(n, notes.tpq)
end

struct Everything end
Base.:∈(x, ::Everything) = true

"""
    separatepitches(notes::Notes [, allowed])

Get a dictionary \"pitch\"=>\"notes of that pitch\".
Optionally only keep pitches that are contained in `allowed`.
"""
function separatepitches(notes::Notes{N}, allowed = Everything()) where {N}
    separated = Dict{UInt8, Notes{N}}()
    for note in notes
        (note.pitch ∈ allowed) && _add_note_to_dict!(separated, note, notes.tpq)
    end
    return separated
end

function _add_note_to_dict!(separated, note::N, tpq) where {N<:AbstractNote}
    if haskey(separated, note.pitch)
        push!(separated[note.pitch], copy(note))
    else
        separated[note.pitch] = Notes{N}(Vector{N}[], tpq)
        push!(separated[note.pitch], copy(note))
    end
end

estimate_delay(notes, sub::Int) = estimate_delay(notes,  0:(1/sub):1)

"""
    estimate_delay(notes, grid)
Estimate the average temporal deviation of the given `notes` from the
quarter note grid point. The notes are classified according to the `grid`
and only notes in the first and last grid bins are used. Their position
is subtracted from the nearby quarter note and the returned value
is the average of this operation.
"""
function estimate_delay(notes::Notes, grid::AbstractVector)
    # TODO: this can be optimized by looping over the notes directly
    # and classifying one by one, and adding to `d` one by one
    clas = classify(notes, grid)
    base = notes[findall(x -> x == 1 || x == length(grid), clas)]
    d = 0.0
    for n in base
        pos = Int(n.position % notes.tpq)
        res = pos ≤ grid[2]*notes.tpq ? pos : pos - notes.tpq
        d += res
    end
    return d / length(base)
end

"""
    estimate_delay_recursive(notes, grid, m)
Do the same as [`estimate_delay`](@ref) but for `m` times, while in each step
shifting the notes by the previously found delay. This improves the accuracy
of the algorithm, because the distribution of the quarter notes is estimated
better and better each time. The function should typically converge
after a couple of `m`.

The returned result is the estimated delay, in integer (ticks), as only integers
can be used to actually shift the notes around.
"""
function estimate_delay_recursive(notes, grid, m)
    totaldelay = 0
    xnotes = notes
    for i in 1:m
        delay = round(Int, estimate_delay(xnotes, grid))
        totaldelay += delay
        xnotes = translate(xnotes, -delay)
    end
    return totaldelay
end


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

# # Generating local matrices for problems

using JuliaFEM

# Plane stress Quad4 element with linear material model:

X = Dict(1 => [0.0, 0.0], 2 => [2.0, 0.0], 3 => [2.0, 2.0], 4 => [0.0, 2.0])
element = Element(Quad4, [1, 2, 3, 4])
update!(element, "geometry", X)
update!(element, "youngs modulus", 288.0)
update!(element, "poissons ratio", 1/3)
problem = Problem(Elasticity, "test problem", 2)
problem.properties.formulation = :plane_stress
add_elements!(problem, [element])
assemble!(problem, 0.0)
K = round.(Matrix(problem.assembly.K), 5)