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

using Test

# code_native / code_llvm (issue #8239)
# It's hard to really test these, but just running them should be
# sufficient to catch segfault bugs.

module ReflectionTest
using Test, Random

function test_ir_reflection(freflect, f, types)
    @test !isempty(freflect(f, types))
    nothing
end

function test_bin_reflection(freflect, f, types)
    iob = IOBuffer()
    freflect(iob, f, types)
    str = String(take!(iob))
    @test !isempty(str)
    nothing
end

function test_code_reflection(freflect, f, types, tester)
    tester(freflect, f, types)
    tester(freflect, f, (types.parameters...,))
    nothing
end

function test_code_reflections(tester, freflect)
    test_code_reflection(freflect, occursin,
                         Tuple{Regex, AbstractString}, tester) # abstract type
    test_code_reflection(freflect, +, Tuple{Int, Int}, tester) # leaftype signature
    test_code_reflection(freflect, +,
                         Tuple{Array{Float32}, Array{Float32}}, tester) # incomplete types
    test_code_reflection(freflect, Module, Tuple{}, tester) # Module() constructor (transforms to call)
    test_code_reflection(freflect, Array{Int64}, Tuple{Array{Int32}}, tester) # with incomplete types
    test_code_reflection(freflect, muladd, Tuple{Float64, Float64, Float64}, tester)
end

test_code_reflections(test_ir_reflection, code_lowered)
test_code_reflections(test_ir_reflection, code_typed)

io = IOBuffer()
Base.print_statement_costs(io, map, (typeof(sqrt), Tuple{Int}))
str = String(take!(io))
@test occursin("map(f, t::Tuple{Any})", str)
@test occursin("sitofp", str)
@test occursin(r"20 .*sqrt_llvm.*::Float64", str)

end # module ReflectionTest

# isbits, isbitstype

@test !isbitstype(Array{Int})
@test isbitstype(Float32)
@test isbitstype(Int)
@test !isbitstype(AbstractString)
@test isbitstype(Tuple{Int, Vararg{Int, 2}})
@test !isbitstype(Tuple{Int, Vararg{Int}})
@test !isbitstype(Tuple{Integer, Vararg{Int, 2}})
@test isbitstype(Tuple{Int, Vararg{Any, 0}})
@test isbitstype(Tuple{Vararg{Any, 0}})
@test isbits(1)
@test isbits((1,2))
@test !isbits([1])
@test isbits(nothing)

# issue #16670
@test isconcretetype(Int)
@test isconcretetype(Vector{Int})
@test isconcretetype(Tuple{Int, Vararg{Int, 2}})
@test !isconcretetype(Tuple{Any})
@test !isconcretetype(Tuple{Integer, Vararg{Int, 2}})
@test !isconcretetype(Tuple{Int, Vararg{Int}})
@test !isconcretetype(Type{Tuple{Integer, Vararg{Int}}})
@test !isconcretetype(Type{Vector})
@test !isconcretetype(Type{Int})
@test !isconcretetype(Tuple{Type{Int}})
@test isconcretetype(DataType)
@test isconcretetype(Union)
@test !isconcretetype(Union{})
@test isconcretetype(Tuple{Union{}})
@test !isconcretetype(Complex)
@test !isconcretetype(Complex.body)
@test !isconcretetype(AbstractArray{Int,1})
struct AlwaysHasLayout{T}
    x
end
@test !isconcretetype(AlwaysHasLayout) && !isconcretetype(AlwaysHasLayout.body)
@test isconcretetype(AlwaysHasLayout{Any})
@test isconcretetype(Ptr{Cvoid})
@test !isconcretetype(Ptr) && !isconcretetype(Ptr.body)

# issue #10165
i10165(::Type) = 0
i10165(::Type{AbstractArray{T,n}}) where {T,n} = 1
@test i10165(AbstractArray{Int,n} where n) == 0
@test which(i10165, Tuple{Type{AbstractArray{Int,n} where n},}).sig == Tuple{typeof(i10165),Type}

# fullname
@test fullname(Base) == (:Base,)
@test fullname(Base.Iterators) == (:Base, :Iterators)

const a_const = 1
not_const = 1
@test isconst(@__MODULE__, :a_const) == true
@test isconst(Base, :pi) == true
@test isconst(@__MODULE__, :pi) == true
@test isconst(@__MODULE__, :not_const) == false
@test isconst(@__MODULE__, :is_not_defined) == false

@test ismutable(1) == false
@test ismutable([]) == true
@test ismutabletype(Int) == false
@test ismutabletype(Vector{Any}) == true
@test ismutabletype(Union{Int, Vector{Any}}) == false

## find bindings tests
@test ccall(:jl_get_module_of_binding, Any, (Any, Any), Base, :sin)==Base

# For curmod_*
include("testenv.jl")

module TestMod7648
using Test
import Base.convert
import ..curmod_name, ..curmod
export a9475, foo9475, c7648, foo7648, foo7648_nomethods, Foo7648

const c7648 = 8
d7648 = 9
const f7648 = 10
foo7648(x) = x
function foo7648_nomethods end
mutable struct Foo7648 end

module TestModSub9475
    using Test
    using ..TestMod7648
    import ..curmod_name
    export a9475, foo9475
    a9475 = 5
    b9475 = 7
    foo9475(x) = x
    let
        @test Base.binding_module(@__MODULE__, :a9475) == @__MODULE__
        @test Base.binding_module(@__MODULE__, :c7648) == TestMod7648
        @test Base.nameof(@__MODULE__) == :TestModSub9475
        @test Base.fullname(@__MODULE__) == (curmod_name..., :TestMod7648, :TestModSub9475)
        @test Base.parentmodule(@__MODULE__) == TestMod7648
    end
end # module TestModSub9475

using .TestModSub9475

let
    @test Base.binding_module(@__MODULE__, :d7648) == @__MODULE__
    @test Base.binding_module(@__MODULE__, :a9475) == TestModSub9475
    @test Base.nameof(@__MODULE__) == :TestMod7648
    @test Base.parentmodule(@__MODULE__) == curmod
end
end # module TestMod7648

let
    @test Base.binding_module(TestMod7648, :d7648) == TestMod7648
    @test Base.binding_module(TestMod7648, :a9475) == TestMod7648.TestModSub9475
    @test Base.binding_module(TestMod7648.TestModSub9475, :b9475) == TestMod7648.TestModSub9475
    @test Set(names(TestMod7648))==Set([:TestMod7648, :a9475, :foo9475, :c7648, :foo7648, :foo7648_nomethods, :Foo7648])
    @test Set(names(TestMod7648, all = true)) == Set([:TestMod7648, :TestModSub9475, :a9475, :foo9475, :c7648, :d7648, :f7648,
                                                :foo7648, Symbol("#foo7648"), :foo7648_nomethods, Symbol("#foo7648_nomethods"),
                                                :Foo7648, :eval, Symbol("#eval"), :include, Symbol("#include")])
    @test Set(names(TestMod7648, all = true, imported = true)) == Set([:TestMod7648, :TestModSub9475, :a9475, :foo9475, :c7648, :d7648, :f7648,
                                                      :foo7648, Symbol("#foo7648"), :foo7648_nomethods, Symbol("#foo7648_nomethods"),
                                                      :Foo7648, :eval, Symbol("#eval"), :include, Symbol("#include"),
                                                      :convert, :curmod_name, :curmod])
    @test isconst(TestMod7648, :c7648)
    @test !isconst(TestMod7648, :d7648)
end

let
    using .TestMod7648
    @test Base.binding_module(@__MODULE__, :a9475) == TestMod7648.TestModSub9475
    @test Base.binding_module(@__MODULE__, :c7648) == TestMod7648
    @test nameof(foo7648) == :foo7648
    @test parentmodule(foo7648, (Any,)) == TestMod7648
    @test parentmodule(foo7648) == TestMod7648
    @test parentmodule(foo7648_nomethods) == TestMod7648
    @test parentmodule(foo9475, (Any,)) == TestMod7648.TestModSub9475
    @test parentmodule(foo9475) == TestMod7648.TestModSub9475
    @test parentmodule(Foo7648) == TestMod7648
    @test nameof(Foo7648) == :Foo7648
    @test basename(functionloc(foo7648, (Any,))[1]) == "reflection.jl"
    @test first(methods(TestMod7648.TestModSub9475.foo7648)) == which(foo7648, (Int,))
    @test TestMod7648 == which(@__MODULE__, :foo7648)
    @test TestMod7648.TestModSub9475 == which(@__MODULE__, :a9475)
end

@test which(===, Tuple{Int, Int}) isa Method
@test length(code_typed(===, Tuple{Int, Int})) === 1
@test only(Base.return_types(===, Tuple{Int, Int})) === Any

module TestingExported
using Test
include("testenv.jl") # for curmod_str
import Base.isexported
global this_is_not_defined
export this_is_not_defined
@test_throws ErrorException("\"this_is_not_defined\" is not defined in module Main") which(Main, :this_is_not_defined)
@test_throws ErrorException("\"this_is_not_exported\" is not defined in module Main") which(Main, :this_is_not_exported)
@test isexported(@__MODULE__, :this_is_not_defined)
@test !isexported(@__MODULE__, :this_is_not_exported)
const a_value = 1
@test which(@__MODULE__, :a_value) === @__MODULE__
@test_throws ErrorException("\"a_value\" is not defined in module Main") which(Main, :a_value)
@test which(Main, :Core) === Main
@test !isexported(@__MODULE__, :a_value)
end

# PR 13825
let ex = :(a + b)
    @test string(ex) == "a + b"
end
foo13825(::Array{T, N}, ::Array, ::Vector) where {T, N} = nothing
@test startswith(string(first(methods(foo13825))),
                 "foo13825(::Array{T, N}, ::Array, ::Vector) where {T, N} in")

mutable struct TLayout
    x::Int8
    y::Int16
    z::Int32
end
tlayout = TLayout(5,7,11)
@test fieldnames(TLayout) == (:x, :y, :z) == Base.propertynames(tlayout)
@test hasfield(TLayout, :y)
@test !hasfield(TLayout, :a)
@test hasfield(Complex, :re)
@test !hasfield(Complex, :qxq)
@test hasproperty(tlayout, :x)
@test !hasproperty(tlayout, :p)
@test [(fieldoffset(TLayout,i), fieldname(TLayout,i), fieldtype(TLayout,i)) for i = 1:fieldcount(TLayout)] ==
    [(0, :x, Int8), (2, :y, Int16), (4, :z, Int32)]
@test fieldnames(Complex) === (:re, :im)
@test_throws BoundsError fieldtype(TLayout, 0)
@test_throws ArgumentError fieldname(TLayout, 0)
@test_throws BoundsError fieldoffset(TLayout, 0)
@test_throws BoundsError fieldtype(TLayout, 4)
@test_throws ArgumentError fieldname(TLayout, 4)
@test_throws BoundsError fieldoffset(TLayout, 4)

@test fieldtype(Tuple{Vararg{Int8}}, 1) === Int8
@test fieldtype(Tuple{Vararg{Int8}}, 10) === Int8
@test_throws BoundsError fieldtype(Tuple{Vararg{Int8}}, 0)
# issue #30505
@test fieldtype(Union{Tuple{Char},Tuple{Char,Char}},2) === Char
@test_throws BoundsError fieldtype(Union{Tuple{Char},Tuple{Char,Char}},3)

@test fieldnames(NTuple{3, Int}) == ntuple(i -> fieldname(NTuple{3, Int}, i), 3) == (1, 2, 3)
@test_throws ArgumentError fieldnames(Union{})
@test_throws BoundsError fieldname(NTuple{3, Int}, 0)
@test_throws BoundsError fieldname(NTuple{3, Int}, 4)

@test fieldnames(NamedTuple{(:z,:a)}) === (:z,:a)
@test fieldname(NamedTuple{(:z,:a)}, 1) === :z
@test fieldname(NamedTuple{(:z,:a)}, 2) === :a
@test_throws ArgumentError fieldname(NamedTuple{(:z,:a)}, 3)
@test_throws ArgumentError fieldnames(NamedTuple)
@test_throws ArgumentError fieldnames(NamedTuple{T,Tuple{Int,Int}} where T)
@test_throws ArgumentError fieldnames(Real)
@test_throws ArgumentError fieldnames(AbstractArray)

@test fieldtype((NamedTuple{T,Tuple{Int,String}} where T), 1) === Int
@test fieldtype((NamedTuple{T,Tuple{Int,String}} where T), 2) === String
@test_throws BoundsError fieldtype((NamedTuple{T,Tuple{Int,String}} where T), 3)

@test fieldtype(NamedTuple, 42) === Any
@test_throws BoundsError fieldtype(NamedTuple, 0)
@test_throws BoundsError fieldtype(NamedTuple, -1)

@test fieldtype(NamedTuple{(:a,:b)}, 1) === Any
@test fieldtype(NamedTuple{(:a,:b)}, 2) === Any
@test fieldtype((NamedTuple{(:a,:b),T} where T<:Tuple{Vararg{Integer}}), 2) === Integer
@test_throws BoundsError fieldtype(NamedTuple{(:a,:b)}, 3)

# issue #32697
@test fieldtype(NamedTuple{(:x,:y), T} where T <: Tuple{Int, Union{Float64, Missing}}, :x) == Int
@test fieldtype(NamedTuple{(:x,:y), T} where T <: Tuple{Int, Union{Float64, Missing}}, :y) == Union{Float64, Missing}

@test fieldtypes(NamedTuple{(:a,:b)}) == (Any, Any)
@test fieldtypes((NamedTuple{T,Tuple{Int,String}} where T)) === (Int, String)
@test fieldtypes(TLayout) === (Int8, Int16, Int32)

import Base: datatype_alignment, return_types
@test datatype_alignment(UInt16) == 2
@test datatype_alignment(TLayout) == 4
let rts = return_types(TLayout)
    @test length(rts) == 2 # general constructor and specific constructor
    @test all(rts .== TLayout)
end

# issue #15447
f15447_line = @__LINE__() + 1
@noinline function f15447(s, a)
    if s
        return a
    else
        nb = 0
        return nb
    end
end
@test functionloc(f15447)[2] == f15447_line

# issue #14346
@noinline function f14346(id, mask, limit)
    if id <= limit && mask[id]
        return true
    end
end
@test functionloc(f14346)[2] == @__LINE__() - 5

# issue #15714
# show variable names for slots and suppress spurious type warnings
function f15714(array_var15714)
    for index_var15714 in eachindex(array_var15714)
        array_var15714[index_var15714] += 0
    end
end

function g15714(array_var15714)
    for index_var15714 in eachindex(array_var15714)
        array_var15714[index_var15714] += 0
    end
    let index_var15714
        for outer index_var15714 in eachindex(array_var15714)
            array_var15714[index_var15714] += 0
        end
        index_var15714
    end
    let index_var15714
        for outer index_var15714 in eachindex(array_var15714)
            array_var15714[index_var15714] += 0
        end
        index_var15714
    end
end

import InteractiveUtils.code_warntype

used_dup_var_tested15714 = false
used_unique_var_tested15714 = false
function test_typed_ir_printing(Base.@nospecialize(f), Base.@nospecialize(types), must_used_vars)
    src, rettype = code_typed(f, types, optimize=false)[1]
    dupnames = Set()
    slotnames = Set()
    for name in src.slotnames
        if name in slotnames || name === Symbol("")
            push!(dupnames, name)
        else
            push!(slotnames, name)
        end
    end
    # Make sure must_used_vars are in slotnames
    for name in must_used_vars
        @test name in slotnames
    end
    must_used_checked = Dict{Symbol,Bool}()
    for sym in must_used_vars
        must_used_checked[sym] = false
    end
    for str in (sprint(io -> code_warntype(io, f, types, optimize=false)),
                repr("text/plain", src))
        for var in must_used_vars
            @test occursin(string(var), str)
        end
        # Check that we are not printing the bare slot numbers
        for i in 1:length(src.slotnames)
            name = src.slotnames[i]
            if name in dupnames
                if name in must_used_vars && occursin(Regex("_$i\\b"), str)
                    must_used_checked[name] = true
                    global used_dup_var_tested15714 = true
                end
            else
                @test !occursin(Regex("_$i\\b"), str)
                if name in must_used_vars
                    global used_unique_var_tested15714 = true
                end
            end
        end
    end
    for sym in must_used_vars
        if sym in dupnames
            @test must_used_checked[sym]
        end
        must_used_checked[sym] = false
    end
    # Make sure printing an AST outside CodeInfo still works.
    str = sprint(show, src.code)
    # Check that we are printing the slot numbers when we don't have the context
    # Use the variable names that we know should be present in the optimized AST
    for i in 2:length(src.slotnames)
        name = src.slotnames[i]
        if name in must_used_vars && occursin(Regex("_$i\\b"), str)
            must_used_checked[name] = true
        end
    end
    for sym in must_used_vars
        @test must_used_checked[sym]
    end
end
test_typed_ir_printing(f15714, Tuple{Vector{Float32}},
                       [:array_var15714,  :index_var15714])
test_typed_ir_printing(g15714, Tuple{Vector{Float32}},
                       [:array_var15714,  :index_var15714])
#This test doesn't work with the new optimizer because we drop slotnames
#We may want to test it against debug info eventually
#@test used_dup_var_tested15715
@test used_unique_var_tested15714

let li = typeof(fieldtype).name.mt.cache.func::Core.MethodInstance,
    lrepr = string(li),
    mrepr = string(li.def),
    lmime = repr("text/plain", li),
    mmime = repr("text/plain", li.def)

    @test lrepr == lmime == "MethodInstance for fieldtype(...)"
    @test mrepr == mmime == "fieldtype(...) in Core"
end

# Linfo Tracing test
function tracefoo end
# Method Tracing test
methtracer(x::Ptr{Cvoid}) = (@test isa(unsafe_pointer_to_objref(x), Method); global didtrace = true; nothing)
let cmethtracer = @cfunction(methtracer, Cvoid, (Ptr{Cvoid},))
    ccall(:jl_register_newmeth_tracer, Cvoid, (Ptr{Cvoid},), cmethtracer)
end
didtrace = false
tracefoo2(x, y) = x*y
@test didtrace
didtrace = false
tracefoo(x::Int64, y::Int64) = x*y
@test didtrace
didtrace = false
ccall(:jl_register_newmeth_tracer, Cvoid, (Ptr{Cvoid},), C_NULL)

# test for reflection over large method tables
for i = 1:100; @eval fLargeTable(::Val{$i}, ::Any) = 1; end
for i = 1:100; @eval fLargeTable(::Any, ::Val{$i}) = 2; end
fLargeTable(::Any...) = 3
@test length(methods(fLargeTable, Tuple{})) == 1
fLargeTable(::Complex, ::Complex) = 4
fLargeTable(::Union{ComplexF32, ComplexF64}...) = 5
@test length(methods(fLargeTable, Tuple{})) == 1
fLargeTable() = 4
@test length(methods(fLargeTable)) == 204
@test length(methods(fLargeTable, Tuple{})) == 1
@test fLargeTable(1im, 2im) == 4
@test fLargeTable(1.0im, 2.0im) == 5
@test_throws MethodError fLargeTable(Val(1), Val(1))
@test fLargeTable(Val(1), 1) == 1
@test fLargeTable(1, Val(1)) == 2
fLargeTable(::Union, ::Union) = "a"
@test fLargeTable(Union{Int, Missing}, Union{Int, Missing}) == "a"
fLargeTable(::Union, ::Union) = "b"
@test length(methods(fLargeTable)) == 205
@test fLargeTable(Union{Int, Missing}, Union{Int, Missing}) == "b"

# issue #15280
function f15280(x) end
@test functionloc(f15280)[2] > 0

# bug found in #16850, Base.url with backslashes on Windows
function module_depth(from::Module, to::Module)
    if from === to || parentmodule(to) === to
        return 0
    else
        return 1 + module_depth(from, parentmodule(to))
    end
end
function has_backslashes(mod::Module)
    for n in names(mod, all = true, imported = true)
        isdefined(mod, n) || continue
        Base.isdeprecated(mod, n) && continue
        f = getfield(mod, n)
        if isa(f, Module) && module_depth(Main, f) <= module_depth(Main, mod)
            continue
        end
        h = has_backslashes(f)
        h === nothing || return h
    end
    return nothing
end
function has_backslashes(f::Function)
    for m in methods(f)
        h = has_backslashes(m)
        h === nothing || return h
    end
    return nothing
end
function has_backslashes(meth::Method)
    if '\\' in string(meth.file)
        return meth
    else
        return nothing
    end
end
has_backslashes(x) = nothing
h16850 = has_backslashes(Base)
if Sys.iswindows()
    if h16850 === nothing
        @warn """No methods found in Base with backslashes in file name,
                 skipping test for `Base.url`"""
    else
        @test !('\\' in Base.url(h16850))
    end
else
    @test h16850 === nothing
end

# PR #18888: code_typed shouldn't cache, return_types should
f18888() = nothing
let
    world = Core.Compiler.get_world_counter()
    m = first(methods(f18888, Tuple{}))
    ft = typeof(f18888)

    code_typed(f18888, Tuple{}; optimize=false)
    @test !isempty(m.specializations) # uncached, but creates the specializations entry
    mi = Core.Compiler.specialize_method(m, Tuple{ft}, Core.svec())
    interp = Core.Compiler.NativeInterpreter(world)
    @test !Core.Compiler.haskey(Core.Compiler.code_cache(interp), mi)
    @test !isdefined(mi, :cache)

    code_typed(f18888, Tuple{}; optimize=true)
    @test !isdefined(mi, :cache)

    Base.return_types(f18888, Tuple{})
    @test Core.Compiler.getindex(Core.Compiler.code_cache(interp), mi) === mi.cache
    @test mi.cache isa Core.CodeInstance
    @test !isdefined(mi.cache, :next)
end

# code_typed_by_type
@test Base.code_typed_by_type(Tuple{Type{<:Val}})[1][2] == Val
@test Base.code_typed_by_type(Tuple{typeof(sin), Float64})[1][2] === Float64

# New reflection methods in 0.6
struct ReflectionExample{T<:AbstractFloat, N}
    x::Tuple{T, N}
end

@test !isabstracttype(Union{})
@test !isabstracttype(Union{Int,Float64})
@test isabstracttype(AbstractArray)
@test isabstracttype(AbstractSet{Int})
@test !isabstracttype(ReflectionExample)
@test !isabstracttype(Int)
@test !isabstracttype(TLayout)

@test !isprimitivetype(Union{})
@test !isprimitivetype(Union{Int,Float64})
@test !isprimitivetype(AbstractArray)
@test !isprimitivetype(AbstractSet{Int})
@test !isprimitivetype(ReflectionExample)
@test isprimitivetype(Int)
@test !isprimitivetype(TLayout)

@test !isstructtype(Union{})
@test !isstructtype(Union{Int,Float64})
@test !isstructtype(AbstractArray)
@test !isstructtype(AbstractSet{Int})
@test isstructtype(ReflectionExample)
@test !isstructtype(Int)
@test isstructtype(TLayout)

let
    wrapperT(T) = Base.typename(T).wrapper
    @test @inferred wrapperT(ReflectionExample{Float64, Int64}) == ReflectionExample
    @test @inferred wrapperT(ReflectionExample{Float64, N} where N) == ReflectionExample
    @test @inferred wrapperT(ReflectionExample{T, Int64} where T) == ReflectionExample
    @test @inferred wrapperT(ReflectionExample) == ReflectionExample
    @test @inferred wrapperT(Union{ReflectionExample{Union{},1},ReflectionExample{Float64,1}}) == ReflectionExample
    @test_throws(ErrorException("typename does not apply to unions whose components have different typenames"),
                 Base.typename(Union{Int, Float64}))
end

# sizeof and nfields
@test sizeof(Int16) == 2
@test sizeof(ComplexF64) == 16
primitive type ParameterizedByte__{A,B} 8 end
@test sizeof(ParameterizedByte__) == 1
@test sizeof(nothing) == 0
@test sizeof(()) == 0
struct TypeWithIrrelevantParameter{T}
    x::Int32
end
@test sizeof(TypeWithIrrelevantParameter) == sizeof(Int32)
@test sizeof(TypeWithIrrelevantParameter{Int8}) == sizeof(Int32)
@test sizeof(:abc) == 3
@test sizeof(Symbol("")) == 0
@test_throws(ErrorException("Abstract type Real does not have a definite size."),
             sizeof(Real))
@test sizeof(Union{ComplexF32,ComplexF64}) == 16
@test sizeof(Union{Int8,UInt8}) == 1
@test_throws ErrorException sizeof(AbstractArray)
@test_throws ErrorException sizeof(Tuple)
@test_throws ErrorException sizeof(Tuple{Any,Any})
@test_throws ErrorException sizeof(String)
@test_throws ErrorException sizeof(Vector{Int})
@test_throws ErrorException sizeof(Symbol)
@test_throws ErrorException sizeof(Core.SimpleVector)
@test_throws ErrorException sizeof(Union{})

@test nfields((1,2)) == 2
@test nfields(()) == 0
@test nfields(nothing) == fieldcount(Nothing) == 0
@test nfields(1) == 0
@test_throws ArgumentError fieldcount(Union{})
@test fieldcount(Tuple{Any,Any,T} where T) == 3
@test fieldcount(Complex) == fieldcount(ComplexF32) == 2
@test fieldcount(Union{ComplexF32,ComplexF64}) == 2
@test fieldcount(Int) == 0
@test_throws(ArgumentError("type does not have a definite number of fields"),
             fieldcount(Union{Complex,Pair}))
@test_throws ArgumentError fieldcount(Real)
@test_throws ArgumentError fieldcount(AbstractArray)
@test_throws ArgumentError fieldcount(Tuple{Any,Vararg{Any}})

# PR #22979

function test_similar_codeinfo(a, b)
    @test a.code == b.code
    @test a.slotnames == b.slotnames
    @test a.slotflags == b.slotflags
end

@generated f22979(x...) = (y = 1; :(x[1] + x[2]))
let
    x22979 = (1, 2.0, 3.0 + im)
    T22979 = Tuple{typeof(f22979), typeof.(x22979)...}
    world = Core.Compiler.get_world_counter()
    match = Base._methods_by_ftype(T22979, -1, world)[1]
    instance = Core.Compiler.specialize_method(match)
    cinfo_generated = Core.Compiler.get_staged(instance)
    @test_throws ErrorException Base.uncompressed_ir(match.method)

    test_similar_codeinfo(code_lowered(f22979, typeof(x22979))[1], cinfo_generated)

    cinfos = code_lowered(f22979, typeof.(x22979), generated=true)
    @test length(cinfos) == 1
    cinfo = cinfos[1]
    test_similar_codeinfo(cinfo, cinfo_generated)
    @test_throws ErrorException code_lowered(f22979, typeof.(x22979), generated=false)
end

module MethodDeletion
using Test, Random

# Deletion after compiling top-level call
bar1(x) = 1
bar1(x::Int) = 2
foo1(x) = bar1(x)
faz1(x) = foo1(x)
@test faz1(1) == 2
@test faz1(1.0) == 1
m = first(methods(bar1, Tuple{Int}))
Base.delete_method(m)
@test bar1(1) == 1
@test bar1(1.0) == 1
@test foo1(1) == 1
@test foo1(1.0) == 1
@test faz1(1) == 1
@test faz1(1.0) == 1

# Deletion after compiling middle-level call
bar2(x) = 1
bar2(x::Int) = 2
foo2(x) = bar2(x)
faz2(x) = foo2(x)
@test foo2(1) == 2
@test foo2(1.0) == 1
m = first(methods(bar2, Tuple{Int}))
Base.delete_method(m)
@test bar2(1.0) == 1
@test bar2(1) == 1
@test foo2(1) == 1
@test foo2(1.0) == 1
@test faz2(1) == 1
@test faz2(1.0) == 1

# Deletion after compiling low-level call
bar3(x) = 1
bar3(x::Int) = 2
foo3(x) = bar3(x)
faz3(x) = foo3(x)
@test bar3(1) == 2
@test bar3(1.0) == 1
m = first(methods(bar3, Tuple{Int}))
Base.delete_method(m)
@test bar3(1) == 1
@test bar3(1.0) == 1
@test foo3(1) == 1
@test foo3(1.0) == 1
@test faz3(1) == 1
@test faz3(1.0) == 1

# Deletion before any compilation
bar4(x) = 1
bar4(x::Int) = 2
foo4(x) = bar4(x)
faz4(x) = foo4(x)
m = first(methods(bar4, Tuple{Int}))
Base.delete_method(m)
@test bar4(1) == 1
@test bar4(1.0) == 1
@test foo4(1) == 1
@test foo4(1.0) == 1
@test faz4(1) == 1
@test faz4(1.0) == 1

# Methods with keyword arguments
fookw(x; direction=:up) = direction
fookw(y::Int) = 2
@test fookw("string") == :up
@test fookw(1) == 2
m = collect(methods(fookw))[2]
Base.delete_method(m)
@test fookw(1) == 2
@test_throws MethodError fookw("string")

# functions with many methods
types = (Float64, Int32, String)
for T1 in types, T2 in types, T3 in types
    @eval foomany(x::$T1, y::$T2, z::$T3) = y
end
@test foomany(Int32(5), "hello", 3.2) == "hello"
m = first(methods(foomany, Tuple{Int32, String, Float64}))
Base.delete_method(m)
@test_throws MethodError foomany(Int32(5), "hello", 3.2)

struct EmptyType end
Base.convert(::Type{EmptyType}, x::Integer) = EmptyType()
m = first(methods(convert, Tuple{Type{EmptyType}, Integer}))
Base.delete_method(m)
@test_throws MethodError convert(EmptyType, 1)

# parametric methods
parametric(A::Array{T,N}, i::Vararg{Int,N}) where {T,N} = N
@test parametric(rand(2,2), 1, 1) == 2
m = first(methods(parametric))
Base.delete_method(m)
@test_throws MethodError parametric(rand(2,2), 1, 1)

# Deletion and ambiguity detection
foo(::Int, ::Int) = 1
foo(::Real, ::Int) = 2
foo(::Int, ::Real) = 3
Base.delete_method(first(methods(foo)))
@test_throws MethodError foo(1, 1)
foo(::Int, ::Int) = 1
foo(1, 1)
Base.delete_method(first(methods(foo)))
@test_throws MethodError foo(1, 1)

# multiple deletions and ambiguities
typeparam(::Type{T}, a::Array{T}) where T<:AbstractFloat = 1
typeparam(::Type{T}, a::Array{T}) where T = 2
for mth in collect(methods(typeparam))
    Base.delete_method(mth)
end
typeparam(::Type{T}, a::AbstractArray{T}) where T<:AbstractFloat = 1
typeparam(::Type{T}, a::AbstractArray{T}) where T = 2
@test typeparam(Float64, rand(2))  == 1
@test typeparam(Int, rand(Int, 2)) == 2

# prior ambiguities (issue #28899)
uambig(::Union{Int,Nothing}) = 1
uambig(::Union{Float64,Nothing}) = 2
@test uambig(1) == 1
@test uambig(1.0) == 2
@test_throws MethodError uambig(nothing)
m = which(uambig, Tuple{Int})
Base.delete_method(m)
@test_throws MethodError uambig(1)
@test uambig(1.0) == 2
@test uambig(nothing) == 2

end

module HasmethodKwargs
using Test
f(x::Int; y=3) = x + y
@test hasmethod(f, Tuple{Int})
@test hasmethod(f, Tuple{Int}, ())
@test hasmethod(f, Tuple{Int}, (:y,))
@test !hasmethod(f, Tuple{Int}, (:jeff,))
@test !hasmethod(f, Tuple{Int}, (:y,), world=typemin(UInt))
g(; b, c, a) = a + b + c
h(; kwargs...) = 4
for gh = (g, h)
    @test hasmethod(gh, Tuple{})
    @test hasmethod(gh, Tuple{}, ())
    @test hasmethod(gh, Tuple{}, (:a,))
    @test hasmethod(gh, Tuple{}, (:a, :b))
    @test hasmethod(gh, Tuple{}, (:a, :b, :c))
end
@test !hasmethod(g, Tuple{}, (:a, :b, :c, :d))
@test hasmethod(h, Tuple{}, (:a, :b, :c, :d))
end

# issue #31353
function f31353(f, x::Array{<:Dict})
end
@test  hasmethod(f31353, Tuple{Any, Array{D}} where D<:Dict)
@test !hasmethod(f31353, Tuple{Any, Array{D}} where D<:AbstractDict)

# issue #26267
module M26267
import Test
foo(x) = x
end
@test !(:Test in names(M26267, all=true, imported=false))
@test :Test in names(M26267, all=true, imported=true)
@test :Test in names(M26267, all=false, imported=true)

# issue #20872
f20872(::Val{N}, ::Val{N}) where {N} = true
f20872(::Val, ::Val) = false
@test which(f20872, Tuple{Val{N},Val{N}} where N).sig == Tuple{typeof(f20872), Val{N}, Val{N}} where N
@test which(f20872, Tuple{Val,Val}).sig == Tuple{typeof(f20872), Val, Val}
@test which(f20872, Tuple{Val,Val{N}} where N).sig == Tuple{typeof(f20872), Val, Val}
@test_throws ErrorException which(f20872, Tuple{Any,Val{N}} where N)
@test which(Tuple{typeof(f20872), Val{1}, Val{2}}).sig == Tuple{typeof(f20872), Val, Val}

module M29962 end
# make sure checking if a binding is deprecated does not resolve it
@test !Base.isdeprecated(M29962, :sin) && !Base.isbindingresolved(M29962, :sin)

# @locals
using Base: @locals
let
    local x, y
    global z
    @test isempty(keys(@locals))
    x = 1
    @test @locals() == Dict{Symbol,Any}(:x=>1)
    y = ""
    @test @locals() == Dict{Symbol,Any}(:x=>1,:y=>"")
    for i = 8:8
        @test @locals() == Dict{Symbol,Any}(:x=>1,:y=>"",:i=>8)
    end
    for i = 42:42
        local x
        @test @locals() == Dict{Symbol,Any}(:y=>"",:i=>42)
    end
    @test @locals() == Dict{Symbol,Any}(:x=>1,:y=>"")
    x = (y,)
    @test @locals() == Dict{Symbol,Any}(:x=>("",),:y=>"")
end

function _test_at_locals1(::Any, ::Any)
    x = 1
    @test @locals() == Dict{Symbol,Any}(:x=>1)
end
_test_at_locals1(1,1)
function _test_at_locals2(a::Any, ::Any, c::T) where T
    x = 2
    @test @locals() == Dict{Symbol,Any}(:x=>2,:a=>a,:c=>c,:T=>typeof(c))
end
_test_at_locals2(1,1,"")
_test_at_locals2(1,1,0.5f0)

@testset "issue #31687" begin
    import InteractiveUtils._dump_function

    @noinline f31687_child(i) = f31687_nonexistent(i)
    f31687_parent() = f31687_child(0)
    params = Base.CodegenParams()
    _dump_function(f31687_parent, Tuple{},
                   #=native=#false, #=wrapper=#false, #=strip=#false,
                   #=dump_module=#true, #=syntax=#:att, #=optimize=#false, :none,
                   #=binary=#false,
                   params)
end

@test nameof(Any) === :Any
@test nameof(:) === :Colon
@test nameof(Core.Intrinsics.mul_int) === :mul_int
@test nameof(Core.Intrinsics.arraylen) === :arraylen

module TestMod33403
f(x) = 1
f(x::Int) = 2
g() = 3

module Sub
import ..TestMod33403: f
f(x::Char) = 3
end
end

@testset "methods with module" begin
    using .TestMod33403: f, g
    @test length(methods(f)) == 3
    @test length(methods(f, (Int,))) == 1

    @test length(methods(f, TestMod33403)) == 2
    @test length(methods(f, [TestMod33403])) == 2
    @test length(methods(f, (Int,), TestMod33403)) == 1
    @test length(methods(f, (Int,), [TestMod33403])) == 1

    @test length(methods(f, TestMod33403.Sub)) == 1
    @test length(methods(f, [TestMod33403.Sub])) == 1
    @test length(methods(f, (Char,), TestMod33403.Sub)) == 1
    @test length(methods(f, (Int,), TestMod33403.Sub)) == 0

    @test length(methods(g, ())) == 1
end

module BodyFunctionLookup
f1(x, y; a=1) = error("oops")
f2(f::Function, args...; kwargs...) = f1(args...; kwargs...)
end

@testset "bodyfunction" begin
    m = first(methods(BodyFunctionLookup.f1))
    f = Base.bodyfunction(m)
    @test occursin("f1#", String(nameof(f)))
    m = first(methods(BodyFunctionLookup.f2))
    f = Base.bodyfunction(m)
    @test f !== Core._apply_iterate
    @test f !== Core._apply
    @test occursin("f2#", String(nameof(f)))
end

@testset "code_typed(; world)" begin
    mod = @eval module $(gensym()) end

    @eval mod foo() = 1
    world1 = Base.get_world_counter()
    @test only(code_typed(mod.foo, ())).second == Int
    @test only(code_typed(mod.foo, (); world=world1)).second == Int

    @eval mod foo() = 2.
    world2 = Base.get_world_counter()
    @test only(code_typed(mod.foo, ())).second == Float64
    @test only(code_typed(mod.foo, (); world=world1)).second == Int
    @test only(code_typed(mod.foo, (); world=world2)).second == Float64
end

@testset "default_tt" begin
    m = Module()
    @eval m f1() = return
    @test Base.default_tt(m.f1) == Tuple{}
    @eval m f2(a) = return
    @test Base.default_tt(m.f2) == Tuple{Any}
    @eval m f3(a::Integer) = return
    @test Base.default_tt(m.f3) == Tuple{Integer}
    @eval m f4() = return
    @eval m f4(a) = return
    @test Base.default_tt(m.f4) == Tuple
end

max_stars_repo_name: JuliaLabs/julia

max_issues_repo_name: JuliaLabs/julia

max_forks_repo_name: JuliaLabs/julia


text:
# ---
# jupyter:
#   jupytext:
#     formats: ipynb,jl:hydrogen
#     text_representation:
#       extension: .jl
#       format_name: hydrogen
#       format_version: '1.3'
#       jupytext_version: 1.11.2
#   kernelspec:
#     display_name: Julia 1.6.2
#     language: julia
#     name: julia-1.6
# ---

# %%
VERSION

# %%
f(x) = 3.0 * x
g(x) = x * 3.0

# %%
@code_llvm debuginfo=:none f(4.0)

# %%
@code_llvm debuginfo=:none g(4.0)

# %%
@code_native debuginfo=:none f(4.0)

# %%
@code_native debuginfo=:none g(4.0)

# %%

ext: jl

max_stars_repo_path: 0021/commutativity of multiplication.jl

max_stars_repo_name: genkuroki/public

max_issues_repo_path: 0021/commutativity of multiplication.jl

max_issues_repo_name: genkuroki/public

max_forks_repo_path: 0021/commutativity of multiplication.jl

max_forks_repo_name: genkuroki/public


text:
# Copyright (c) 2016, Jernej Kovacic
# Licensed under the MIT "Expat" License.
#
# See LICENSE.md for more details.

module IncGammaBeta

export inc_gamma_upper, inc_gamma_lower, inc_gamma_upper_reg, inc_gamma_lower_reg
export inc_gamma_lower_inv, inc_gamma_upper_inv, inc_gamma_lower_reg_inv, inc_gamma_upper_reg_inv
export inc_beta_lower, inc_beta_upper, inc_beta_lower_reg, inc_beta_upper_reg
export inc_beta_lower_inv, inc_beta_upper_inv, inc_beta_lower_reg_inv, inc_beta_upper_reg_inv
export SpecfunUndefinedError, NoConvergenceError

using SpecialFunctions

include("specfun.jl")
include("exception.jl")

"""
A Julia package with functions that implement incomplete gamma and beta functions
of all kinds (upper and lower, regularized and unegularized). All functions except
inverse ones also support complex arguments.

The package is written entirely in Julia and does not require any 3rd party libraries
or additional Julia packages.
"""
IncGammaBeta

end  # module


text:
module ChunkIterators
using ..BigArrays
using ..BigArrays.Indexes

import ..BigArrays: AbstractBigArray

export ChunkIterator

struct ChunkIterator{N}
    globalRange     ::CartesianIndices{N}
    chunkSize       ::NTuple{N}
    chunkIDRange    ::CartesianIndices{N}
    # this offset really means the starting coordinate of the real data
    offset          ::CartesianIndex{N}
end

function ChunkIterator( globalRange::CartesianIndices{N},
                           chunkSize::NTuple{N};
                           offset::CartesianIndex{N} = CartesianIndex{N}() - 1 ) where N
    chunkIDStart = index2chunkid( first(globalRange), chunkSize; offset=offset )
    chunkIDStop  = index2chunkid( last(globalRange),  chunkSize; offset=offset )
    ranges = map((x,y)->x:y, chunkIDStart, chunkIDStop)
    chunkIDRange = CartesianIndices(ranges)
    ChunkIterator( globalRange, chunkSize, chunkIDRange, offset )
end

function ChunkIterator(   idxes::Tuple,
                     chunkSize::NTuple{N};
                     offset::CartesianIndex{N} = CartesianIndex{N}(0)) where N
    # the offset in neuroglancer really means the starting coordinate of valid data
    # since bigarray assumes infinite data range,
    # here we only need to use it to adjust the alignment of chunks
    # so we only use the mod to make offset
    offset = CartesianIndex( map((o,c) -> mod(o,c), offset.I, chunkSize) )
    #idxes = map(index2unit_range, idxes)
    globalRange = CartesianIndices(idxes)
    ChunkIterator( globalRange, chunkSize; offset=offset )
end

function ChunkIterator( ba::AbstractBigArray )
    ChunkIterator( ba.globalRange, ba.chunkSize )
end

function Base.length( iter::ChunkIterator )
    length( iter.globalRange )
end

function Base.eltype( iter::ChunkIterator )
    eltype( iter.globalRange )
end

"""
    Base.iterate(iter::ChunkIterator, state=true)
the initial state is true, meaning it is the start of iteration
"""
function Base.iterate(iter::ChunkIterator{N},
                      state::Union{Nothing, Bool, CartesianIndex{N}}=true) where N
    if isa(state, Bool) && state == true
        # this is the first iteration
        chunkID, nextState = iterate(iter.chunkIDRange)
    elseif state == last(iter.chunkIDRange)
        return nothing
    else
        chunkID, nextState = iterate(iter.chunkIDRange, state)
    end
    # get current global range in this chunk
    start = map((x,y,z,o)->max((x-1)*y+1+o, z), chunkID.I,
                            iter.chunkSize, first(iter.globalRange).I,
                            iter.offset.I )
    stop  = map((x,y,z,o)->min(x*y+o, z),       chunkID.I,
                            iter.chunkSize, last(iter.globalRange).I,
                            iter.offset.I )
    # the global range of the cutout in this chunk
    cutoutGlobalRange = CartesianIndices(map((x,y)->x:y, start, stop))
    # the range inside this chunk
    rangeInChunk  = global_range2chunk_range( cutoutGlobalRange, iter.chunkSize; offset=iter.offset)
    # the range inside the buffer
    rangeInBuffer = global_range2buffer_range(cutoutGlobalRange, iter.globalRange)
    # the global range of this chunk
    chunkGlobalRange = chunkid2global_range( chunkID, iter.chunkSize; offset=iter.offset )

    return (chunkID, chunkGlobalRange, cutoutGlobalRange, rangeInChunk, rangeInBuffer), nextState
end

end # end of module


text:
using PyCall
using Plots; pyplot()
using Statistics
using LinearAlgebra

vtk = pyimport("vtk")
np = pyimport("numpy")

reader = vtk.vtkUnstructuredGridReader()
reader.SetFileName("../data/h06.vtk")
reader.Update()
output = reader.GetOutput()
points = output.GetPoints()
r = [collect(points.GetPoint(i-1)) for i in 1:points.GetNumberOfPoints()]
length(r)

z = [r[i][3] for i=1:262144]

println(sum( abs.(z.-mean(z)))/length(z))
println(std(z))

z = reshape(z,512,512)
heatmap(z)

include("dfa.jl")
dfa2d(z)

function fbmArea(z)
    nx, ny = size(z)
    area = 0.0
    for i=1:nx-1
        for j=1:ny-1
            v1 = [0; 1; z[i,j+1] - z[i,j]]
            v2 = [1; 0; z[i+1,j] - z[i,j]]
            v3 = [0; -1; z[i+1,j] - z[i+1,j+1]]
            v4 = [-1; 0; z[i,j+1] - z[i+1,j+1]]
            area += (norm(v2 × v1)+norm(v4 × v3))/2
        end
    end

    Lp  = 0.0
    for i=1:nx-1 # eval Lp for each line

        for j=1:ny-1
            v1 = [1; 0; z[i+1,j] - z[i,j]]
            Lp += norm(v1)
        end

    end

    return area/(nx*ny), (Lp/ny)/nx
end

fbmArea(z)

2^9


text:
immutable TightenedSITB <: Cut{2}; groups; end
export TightenedSITB

function getcut(prob::Problem, vars::Variables, ::TightenedSITB, i, j)
	N = prob.N
	Lˣ,  Lʸ  = prob.W,   prob.H
	lbˣ, lbʸ = prob.wlb, prob.hlb
	(cˣ,cʸ,ℓˣ,ℓʸ,_,_,_,zˣ,zʸ) = unpack(vars)
	@LinearConstraints(begin
		cˣ[i] ≥      0.5ℓˣ[i] + lbˣ[j]*zˣ[j,i]
		cˣ[i] ≤ Lˣ - 0.5ℓˣ[i] - lbˣ[j]*zˣ[i,j]
		cʸ[i] ≥      0.5ℓʸ[i] + lbʸ[j]*zʸ[j,i]
		cʸ[i] ≤ Lʸ - 0.5ℓʸ[i] - lbʸ[j]*zʸ[i,j]
	end)
end


text:
function rankordinal{T<:Real}(a::Array{T,1})
    collect(1:length(a))
end

max_stars_repo_name: ethansaxenian/RosettaDecode

max_issues_repo_name: ethansaxenian/RosettaDecode

max_forks_repo_name: ethansaxenian/RosettaDecode


text:
@enum QuonConst begin
    Zero
    Pi
    HalfPi
    NegHalfPi
    InfZero
    InfPi
end

const QuonComplex = Union{QuonConst, ComplexF64}

function Base.Complex(p::QuonConst)
    p === Zero && return 0.0im
    p === Pi && return π*im
    p === HalfPi && return π/2*im
    p === NegHalfPi && return -π/2*im
    @warn "Representing infinity greatness with a Complex number will cause numerical unstable"
    p === InfZero && return -Inf + 0.0im
    p === InfPi && return Inf + 0.0im
end

function Base.:(-)(p::QuonConst)
    p === HalfPi && return NegHalfPi
    p === NegHalfPi && return HalfPi
    p === InfZero && return InfPi
    p === InfPi && return InfZero
    return p
end

function Base.:(+)(p1::QuonConst, p2::QuonConst)
    p1 === Zero && return p2
    p2 === Zero && return p1
    if p1 === Pi
        p2 === Pi && return Zero
        p2 === HalfPi && return NegHalfPi
        p2 === NegHalfPi && return HalfPi
    elseif p1 === HalfPi
        p2 === Pi && return NegHalfPi
        p2 === HalfPi && return Pi
        p2 === NegHalfPi && return Zero
    elseif p1 === NegHalfPi
        p2 === Pi && return HalfPi
        p2 === HalfPi && return Zero
        p2 === NegHalfPi && return Pi
    end

    (p1 in (InfZero, InfPi)) && return p1
    return p2
end
function Base.:(+)(p1::Complex, p2::QuonConst)
    p = p1
    p2 === Zero && return to_quon_const(p)
    p2 === Pi && return to_quon_const(p + π*im)
    p2 === HalfPi && return to_quon_const(p + π/2*im)
    p2 === NegHalfPi && return to_quon_const(p - π/2*im)
    return p2
end
Base.:(+)(p1::QuonConst, p2::Complex) = p2 + p1

Base.zero(::QuonConst) = Zero
Base.zero(::Type{QuonConst}) = Zero
Base.zero(::Type{QuonComplex}) = Zero
function Base.imag(p::QuonConst)
    p === Pi && return pi
    p === HalfPi && return pi/2
    p === NegHalfPi && return -pi/2
    return 0.0
end
Base.real(p::QuonConst) = real(Complex(p))

function change_direction(p::QuonConst)
    p === Zero && return InfZero
    p === InfZero && return Zero
    p === Pi && return InfPi
    p === InfPi && return Pi
    p === HalfPi && return NegHalfPi
    p === NegHalfPi && return HalfPi
end

is_pi(p::QuonConst) = (p === Pi)
is_zero(p::QuonConst) = (p === Zero)
is_half_pi(p::QuonConst) = (p === HalfPi)
is_neg_half_pi(p::QuonConst) = (p === NegHalfPi)
is_inf_pi(p::QuonConst) = (p === InfPi)
is_inf_zero(p::QuonConst) = (p === InfZero)

to_quon_const(p::QuonConst) = p
function to_quon_const(p::Complex)
    is_pi(p) && return Pi
    is_zero(p) && return Zero
    is_half_pi(p) && return HalfPi
    is_neg_half_pi(p) && return NegHalfPi
    is_inf_pi(p) && return InfPi
    is_inf_zero(p) && return InfZero

    return p
end

# one const
function yang_baxter_param(α::QuonConst, β::QuonComplex, γ::QuonComplex)
    if α === Zero
        return to_quon_const(β), to_quon_const(γ), Zero
    elseif α === Pi
        return to_quon_const(-β), γ+Pi, Zero
    elseif α === InfZero
        return Zero, InfZero, change_direction(change_direction(β) + γ)
    elseif α === InfPi
        return Pi, InfZero, change_direction(-change_direction(β) + γ)
    end
    return yang_baxter_param(Complex(α), β, γ)
end

function yang_baxter_param(α::QuonComplex, β::QuonConst, γ::QuonComplex)
    if β === Zero
        return Zero, α+γ, Zero
    elseif β === Pi
        return Pi, -α+γ, Zero
    elseif β === InfZero
        return change_direction(α), InfZero, change_direction(γ)
    elseif β === InfPi
        return -change_direction(α), InfZero, -change_direction(γ)
    end
    return yang_baxter_param(α, Complex(β), γ)
end

function yang_baxter_param(α::QuonComplex, β::QuonComplex, γ::QuonConst)
    c, b, a = yang_baxter_param(γ, β, α)
    return a, b, c
end

# two consts
function yang_baxter_param(α::QuonConst, β::QuonConst, γ::QuonComplex)
    if α === Zero
        return to_quon_const(β), to_quon_const(γ), Zero
    elseif α === Pi
        return -β, γ+Pi, Zero
    elseif α === InfZero
        return Zero, InfZero, change_direction(change_direction(β) + γ)
    elseif α === InfPi
        return Pi, InfZero, change_direction(-change_direction(β) + γ)
    end
    if β === Zero
        return Zero, α+γ, Zero
    elseif β === Pi
        return Pi, -α+γ, Zero
    elseif β === InfZero
        return change_direction(α), InfZero, change_direction(γ)
    elseif β === InfPi
        return -change_direction(α), InfZero, -change_direction(γ)
    end
    return yang_baxter_param(Complex(α), Complex(β), γ)
end

function yang_baxter_param(α::QuonComplex, β::QuonConst, γ::QuonConst)
    c, b, a = yang_baxter_param(γ, β, α)
    return a, b, c
end

function yang_baxter_param(α::QuonConst, β::QuonComplex, γ::QuonConst)
    if α === Zero
        return to_quon_const(β), γ, Zero
    elseif α === Pi
        return to_quon_const(-β), γ+Pi, Zero
    elseif α === InfZero
        return Zero, InfZero, change_direction(change_direction(β) + γ)
    elseif α === InfPi
        return Pi, InfZero, change_direction(-change_direction(β) + γ)
    end
    if γ === Zero
        return Zero, α, to_quon_const(β)
    elseif γ === Pi
        return Zero, α+Pi, to_quon_const(-β)
    elseif γ === InfZero
        return change_direction(change_direction(β) + α), InfZero, Zero
    elseif γ === InfPi
        return change_direction(-change_direction(β) + α), InfZero, Pi
    end
    return yang_baxter_param(Complex(α), β, Complex(γ))
end

# three consts
function yang_baxter_param(α::QuonConst, β::QuonConst, γ::QuonConst)
    if (α, β, γ) ⊆ (HalfPi, NegHalfPi)
        if α === HalfPi && β === HalfPi && γ === HalfPi
            return HalfPi, HalfPi, HalfPi
        elseif α === HalfPi && β === HalfPi && γ === NegHalfPi
            return NegHalfPi, HalfPi, HalfPi
        elseif α === HalfPi && β === NegHalfPi && γ === HalfPi
            return NegHalfPi, HalfPi, NegHalfPi
        elseif α === HalfPi && β === NegHalfPi && γ === NegHalfPi
            return NegHalfPi, NegHalfPi, HalfPi
        elseif α === NegHalfPi && β === HalfPi && γ === HalfPi
            return HalfPi, HalfPi, NegHalfPi
        elseif α === NegHalfPi && β === HalfPi && γ === NegHalfPi
            return HalfPi, NegHalfPi, HalfPi
        elseif α === NegHalfPi && β === NegHalfPi && γ === HalfPi
            return HalfPi, NegHalfPi, NegHalfPi
        elseif α === NegHalfPi && β === NegHalfPi && γ === NegHalfPi
            return NegHalfPi, NegHalfPi, NegHalfPi
        end
    end
    if α === Zero
        return β, γ, Zero
    elseif α === Pi
        return -β, γ+Pi, Zero
    elseif α === InfZero
        return Zero, InfZero, change_direction(change_direction(β) + γ)
    elseif α === InfPi
        return Pi, InfZero, change_direction(-change_direction(β) + γ)
    end
    if β === Zero
        return Zero, α+γ, Zero
    elseif β === Pi
        return Pi, -α+γ, Zero
    elseif β === InfZero
        return change_direction(α), InfZero, change_direction(γ)
    elseif β === InfPi
        return -change_direction(α), InfZero, -change_direction(γ)
    end
    if γ === Zero
        return Zero, α, β
    elseif γ === Pi
        return Zero, α+Pi, -β
    elseif γ === InfZero
        return change_direction(change_direction(β) + α), InfZero, Zero
    elseif γ === InfPi
        return change_direction(-change_direction(β) + α), InfZero, Pi
    end
end


text:
# Helper functions
sameSense(pra::Int, ra::Int) = mod(pra,2)==mod(ra,2)
strengthen(pra::Int, ra::Int) = ra>pra
weaken(pra::Int,ra::Int) = ra < pra
strongAlert(ra::Int) = ra > 2

# Reward function for HCAS MDP
function POMDPs.reward(mdp::HCAS_MDP, s::stateType, ra::actType)
    r = s[1]; t = s[2]; p = s[3]; vown = s[4]; vint = s[5]; pra = s[6]
    tau = mdp.currentTau
    rew = 0.0

    relx = r*cos(t); rely = r*sin(t); dx = vint*cos(p)-vown; dy = vint*sin(p)
    dv2 = dx*dx+dy*dy
    dCPA = r
    tCPA = 0.0
    if dv2>0.0
        tCPA = (-relx*dx - rely*dy)/dv2
        xT = relx + dx*tCPA
        yT = rely + dy*tCPA
        if tCPA>0.0
            dCPA = sqrt(xT*xT + yT*yT)
        else
            tCPA=0.0
        end
    end

    if tau==0.0
        if r<=500.0
            rew-=1.0
        else
            rew -= 1.0*exp(-(r-500.0)/500.0) # /300.0 for v5
        end
    end
    factor=1.0
    if pra!=COC
        factor=0.1 # factor always equals 1.0 for v5
    end
    if ra != COC
        rew-=5e-4 *factor#1e-9
        if strongAlert(ra)
            rew-=2e-3*factor
        end
        if (pra!=COC) && (!sameSense(pra,ra))
            rew-=5e-2
        elseif strengthen(pra,ra)
            rew-=1e-3
        elseif weaken(pra,ra)
            rew-=2e-4
        end

    # ra = COC
    else
        rew -= 1e-2*exp(-dCPA/500.0)*exp(-tCPA/10.0) /factor #dCPA/300.0 for v5
    end
    return rew
end

max_stars_repo_name: cherrywoods/HorizontalCAS

max_issues_repo_name: ahmed-irfan/HorizontalCAS

max_forks_repo_name: ahmed-irfan/HorizontalCAS


text:
abstract type AbstractTensorProductLayer <: Function end
"""
Constructs the Tensor Product Layer, which takes as input an array of n tensor
product basis, [B_1, B_2, ..., B_n] a data point x, computes
z[i] = W[i,:] ⨀ [B_1(x[1]) ⨂ B_2(x[2]) ⨂ ... ⨂ B_n(x[n])], where W is the layer's weight,
and returns [z[1], ..., z[out]].

```julia
TensorLayer(model,out,p=nothing)
```
Arguments:
- `model`: Array of TensorProductBasis [B_1(n_1), ..., B_k(n_k)], where k corresponds to the dimension of the input.
- `out`: Dimension of the output.
- `p`: Optional initialization of the layer's weight. Initialized to standard normal by default.
"""
struct TensorLayer{M<:Array{TensorProductBasis},P<:AbstractArray,Int} <: AbstractTensorProductLayer
    model::M
    p::P
    in::Int
    out::Int
    function TensorLayer(model,out,p=nothing)
        number_of_weights = 1
        for basis in model
            number_of_weights *= basis.n
        end
        if p === nothing
            p = randn(out*number_of_weights)
        end
        new{Array{TensorProductBasis},typeof(p),Int}(model,p,length(model),out)
    end
end

function (layer::TensorLayer)(x,p=layer.p)
    model,out = layer.model,layer.out
    W = reshape(p, out, Int(length(p)/out))
    tensor_prod = model[1](x[1])
    for i in 2:length(model)
        tensor_prod = kron(tensor_prod,model[i](x[i]))
    end
    z = W*tensor_prod
    return z
end


text:
# Find and add MDP's path to Julia search path
cdir = pwd(); ii = findlast("MDP", cdir);
if ii===nothing
    # MDPpath = /path/to/MDP;
    display("MDP's path is not found. Please uncomment the above line and set the path to MDP.");
else
    MDPpath = cdir[1:ii[end]] * "/";
    include(MDPpath * "Installation/setpath.jl");
end

# External packages
using Revise, DelimitedFiles, SymPy

# MDP packages
using Preprocessing, Gencode, Postprocessing

app = Preprocessing.initializeapp(cdir,MDPpath);
#app.buildexec = 1;
app.platform = "cpu";    # change it to app.platform = "gpu" => MDP runs on Nvidia GPU
app.unitstyle = "metal"; # unit system
app.pbc = [1 1 1];       # periodic boundary conditions
nlattice = 4;            # note that # of atoms = 2 * nlattice^3
app.appname = "snapInP" * string(nlattice);

#  neighbor list
app.neighskin = 1.0;
app.neighcell = 0;       # 0 -> O(N^2) algorithm, 1 -> Cell-linked list algorithm to form neighbor list
app.neighmax = 150;      # maximum number of neighbors allowed
app.neighevery = 1;      # perform neighbor rebuild list for "every" iterations
app.neighdelay = 0;      # delay neighbor rebuild
app.neighcheck = 1;      # 0 -> neighbor rebuild by delay and every settings,
                         # 1 -> neighbor rebuild occurs if some atom has moved more than half the skin distance
app.neighpair = 0;       # 0 -> full neighbor list, 1 -> half neighbor list for pair potentials

# atom
app.atommasses = [114.76 30.98];
app.atomcharges = [0 0];

# lattice, region, domain
basisatomtype = [1 1 1 1 2 2 2 2]; # 1 => In, 2 => P
app.lattice = setlattice("diamond", 5.83,[1.0 0.0 0.0],[0.0 1.0 0.0],[0.0 0.0 1.0],[1.0 1.0 1.0],basisatomtype);
app.region = setregion([nlattice nlattice nlattice]);

#app.lattice = Preprocessing.setlattice("bcc", 3.316);

# nonparametric potential descriptors
app.descriptor = "snap";    # (snap, shp)
app.snaprcutfac = 1.0;
app.snaptwojmax = 6;
app.snaprfac0 = 0.99363;
app.snaprmin0 = 0;
app.snapbzeroflag = 1;
app.snapquadraticflag = 0;
app.snapswitchflag = 1;
app.snapchemflag = 1;
app.snapbnormflag = 1;
app.snapwselfallflag = 1;
app.snapnelem = 2;    # # of elements
app.snapncoeff = 241;  # # of coeffients per element
app.snapelemradius = [3.81205 3.82945]; # radius per element
app.snapelemweight = [1 0.929316];   # weight per element
app.snapcoefffile = "snapInP.coeff";
app.rcutml = 2*app.snaprcutfac*maximum(app.snapelemradius[:]); # maximum cut-off radius for nonparametric potential

# parametric potential descriptors
app.potentialfile = "ZBLpotentialInP.jl";   # empirical potential file
include(app.potentialfile);  # include the potential file
zblcutinner = 4;
zblcutouter = 4.2;
zblz1 = 49;
zblz2 = 15;
app.pot2a = reshape([1],1,1);     # list of nonbonded pair potentials
app.mu2a = [zblcutinner zblcutouter zblz1 zblz2];    # parameters for ZBL potential
app.rcut2a = reshape([zblcutouter],1,1);  # cut-off radius

# time integration and outputs
app.ensemblemode = "nve";
app.time = 0.0;             # initial time
app.dt = 0.5e-3;            # time step size
app.ntimesteps = 100;       # # time steps
app.globalfreq = 10;

temp = 300;  seed0 = 4928459; distflag = 0; sumflag = 0; loopflag = 2; momflag = 1; rotflag = 0;
app.createvelocity = [temp seed0 distflag sumflag loopflag momflag rotflag];

# generate input files and run simulation
app = mdp(app);

max_stars_repo_name: cesmix-mit/MDP

max_issues_repo_name: cesmix-mit/MDP

max_forks_repo_name: cesmix-mit/MDP


text:
using Test, Sandbox, Scratch

function with_temp_scratch(f::Function)
    mktempdir() do scratch_dir
        with_scratch_directory(scratch_dir) do
            f()
        end
    end
end

if executor_available(DockerExecutor)
    @testset "Docker" begin
        with_temp_scratch() do
            # With a temporary scratch directory, let's start by testing load/save timestamps
            @test !isfile(Sandbox.timestamps_path())
            @test isempty(Sandbox.load_timestamps())
            Sandbox.save_timestamp("foo", 1.0)
            @test isfile(Sandbox.timestamps_path())
            timestamps = Sandbox.load_timestamps()
            @test timestamps["foo"] == 1.0

            # Next, let's actually create a docker image out of our alpine rootfs image
            mktempdir() do rootfs_path
                cp(Sandbox.alpine_rootfs(), rootfs_path; force=true)
                @test Sandbox.should_build_docker_image(rootfs_path)
                @test_logs (:info, r"Building docker image") match_mode=:any begin
                    Sandbox.build_docker_image(rootfs_path; verbose=true)
                end

                # Ensure that it doesn't try to build again since the content is unchanged
                @test !Sandbox.should_build_docker_image(rootfs_path)
                @test_logs begin
                    Sandbox.build_docker_image(rootfs_path; verbose=true)
                end

                # Change the content
                chmod(joinpath(rootfs_path, "bin", "busybox"), 0o775)
                @test Sandbox.should_build_docker_image(rootfs_path)
                @test_logs (:info, r"Building docker image") match_mode=:any begin
                    Sandbox.build_docker_image(rootfs_path; verbose=true)
                end

                # Ensure that it once again doesn't try to build
                @test !Sandbox.should_build_docker_image(rootfs_path)
                @test_logs begin
                    Sandbox.build_docker_image(rootfs_path; verbose=true)
                end
            end
        end

        @testset "probe_executor" begin
            @test_logs (:info, "Testing Docker Executor (read-only, read-write)") begin
                with_executor(DockerExecutor) do exe
                    @test probe_executor(exe; test_read_only_map=true, test_read_write_map=true, verbose=true)
                end
            end
        end
    end
else
    @error("Skipping Docker tests, as it does not seem to be available")
end


text:
import Base: convert, copyto!, Array, vec
import LinearAlgebra: norm
import Plots: plot, plot!

function get_cost_matrices(solver::Solver)
    solver.obj.cost.Q, solver.obj.cost.R, solver.obj.cost.Qf
end

function get_sizes(solver::Solver)
    return solver.model.n, solver.model.m, solver.N
end

function get_sizes(X::Vector{T}, U::Vector{T}) where {T<:MVector}
    N = length(X)
    n,m = length(X[1]), length(U[1])
end

function get_N(solver::Solver,method::Symbol)
    get_N(solver.N,method)
end

function get_N(N0::Int,method::Symbol)
    if method == :midpoint
        N,N_ = N0,N0
    elseif method == :trapezoid
        N,N_ = N0,N0
    elseif method == :hermite_simpson_separated
        N,N_ = 2N0-1,2N0-1
    elseif method == :hermite_simpson
        N,N_ = N0,2N0-1
    end
end

function get_sizes(X::AbstractArray,U::AbstractArray)
    n,N = size(X)
    m = size(U,1)
    return n,m,N
end

function norm(X::Vector{MVector{S,Float64}}) where {S}
    sqrt(norm2(X))
end

function norm2(X::Union{Vector{MVector{S,Float64} where S}, Vector{Vector{Float64}}})
    v = 0.
    for i = 1:size(X,1)
        v += X[i]'X[i]
    end
    v
end

function norm2(X::Union{MVector{S,Float64},Vector{Float64}}) where {S}
    norm(X)^2
end

function norm(X::Vector{MVector{S,Float64}},inds::UnitRange{Int64}) where {S}
    sqrt(norm2(X,inds))
end

function norm2(X::Union{Vector{MVector{S,Float64}} where S, Vector{Vector{Float64}}},inds::UnitRange{Int64})
    v = 0.
    for i = 1:size(X,1)
        v += X[i][inds]'X[i][inds]
    end
    v
end

function norm2(x::Vector{Real})
    v = 0.
    for val in x
        v += val^2
    end
    return v
end

function to_array(X::Vector{Vector{Float64}})
    N = length(X)
    n = length(X[1])
    Y = zeros(n,N)
    for i = 1:N
        Y[:,i] = X[i]
    end
    Y
end

function to_array(X::Vector{Matrix{Float64}})
    N = length(X)
    n,m = size(X[1])
    Y = zeros(n,m,N)
    for i = 1:N
        Y[:,:,i] = X[i]
    end
    Y
end

function to_array(X::Vector{T}) where {T<:MArray}
    N = length(X)
    Y = zeros(size(X[1])...,N)
    ax = axes(Y)
    for i = 1:N
        Y[ax[1:end-1]...,i] = X[i]
    end
    Y
end

function to_array(A::Vector{D} where {D<:Diagonal})
    n = size(A[1],1)
    N = length(A)
    B = zeros(n,n,N)
    for i = 1:N
        B[:,:,i] = A[i]
    end
    return B
end

function vec(A::Trajectory)
	vec(to_array(A))
end

function to_dvecs(X::AbstractArray)
    N = size(X)[end]
    ax = axes(X)[1:end-1]
    [X[ax...,i] for i = 1:N]
end

function to_svecs(X::AbstractArray)
    N = size(X)[end]
    s = size(X)[1:end-1]
    ax = axes(X)[1:end-1]
    [MArray{Tuple{s...}}(X[ax...,1]) for i = 1:N]
end

function copyto!(A::Vector{T}, B::AbstractArray{Float64}) where {T<:Union{MArray,VecOrMat}}
    N = size(B)[end]
    ax = axes(B)[1:end-1]
    for i = 1:N
        A[i] = B[ax...,i]
    end
    A
end

function copyto!(A::Vector{T}, B::Vector{T}) where {T<:VecOrMat{Float64}}
    A .= copy.(B)
end

function copyto!(A::Vector{T}, B::Vector{T}) where {T<:SArray{S,Float64,N,L} where {S,N,L}}
    A .= copy.(B)
end

function Array{Float64,3}(X::Vector{D}) where {D<:Diagonal}
    to_array(X)
end

function convert(::Type{Array{Float64,3}}, X::Vector{Diagonal{Float64,V}}) where {V}
    to_array(X)
end

function convert(::Type{Array}, X::Vector{T}) where {T<:MArray}
    to_array(X)
end

function println(level::Symbol, msg::String)::Nothing
    if level_priorities[level] ≥ level_priorities[debug_level::Symbol]
        println(msg)
    end
end

function root_dir()
    joinpath(dirname(pathof(TrajectoryOptimization)),"..")
end

print_info(msg) = println(:info,msg)
print_debug(msg) = println(:debug,msg)

"""
$(SIGNATURES)

    Convert ZYX Euler angles to quaternion (q =[v;s])
"""
function eul2quat(eul)
    ## Translate the given Euler angle with a specified axis rotation sequence into the corresponding quaternion:

    quat = zeros(4,1)
    # case 'ZYX'
	c_1 = cos(eul[1,1]*0.5)
	s_1 = sin(eul[1,1]*0.5)
	c_2 = cos(eul[2,1]*0.5)
	s_2 = sin(eul[2,1]*0.5)
	c_3 = cos(eul[3,1]*0.5)
	s_3 = sin(eul[3,1]*0.5)

    # https://ntrs.nasa.gov/archive/nasa/casi.ntrs.nasa.gov/19770024290.pdf, the original algorithm has q = [s;v]
	quat[4,1] = c_1*c_2*c_3 + s_1*s_2*s_3
	quat[1,1] = c_1*c_2*s_3 - s_1*s_2*c_3
	quat[2,1] = c_1*s_2*c_3 + s_1*c_2*s_3
	quat[3,1] = s_1*c_2*c_3 - c_1*s_2*s_3

    quat
end

function print_solver(solver::Solver,name::String,io::IO=STDOUT)
    println(io,"###  $name  ###")

    println(io,"\nModel Props")
    println(io,"-----------")
    println(io,"\t n: $(solver.model.n)")
    println(io,"\t m: $(solver.model.m)")
    println(io,"\t inplace dynamics?: $(is_inplace_dynamics(solver.model))")

    println(io,"\nObjective")
    println(io,"----------")
    println(io,"\t tf: $(obj.tf)")
    println(io,"\t x0: $(obj.x0)")
    println(io,"\t xf: $(obj.xf)")
    println(io,"\t Q: $(diag(obj.Q))")
    println(io,"\t R: $(diag(obj.R))")

    println(io,"\nSolver Settings")
    println(io,"-----------------")
    println(io,"\t dt: $(solver.dt)")
    println(io,"\t N: $(solver.N)")
    println(io,"\t integration: $(solver.integration)")

end

"""
$(SIGNATURES)
Generate random circle constraints
"""
function generate_random_circle_obstacle_field(n_circles::Int64,x_rand::Float64=10.0,y_rand::Float64=10.0,r_rand::Float64=0.5)
    x0 = x_rand*rand(n_circles)
    y0 = y_rand*rand(n_circles)
    r = r_rand*ones(n_circles)

    function constraints!(c,x,u)
        for i = 1:n_circles
            c[i] = circle_constraint(x,x0[i],y0[i],r[i])
        end
    end
    return constraints!, zip(x0, y0, r)
end

"""
$(SIGNATURES)
Generate random sphere constraints
"""
function generate_random_sphere_obstacle_field(n_spheres::Int64,x_rand::Float64=10.0,y_rand::Float64=10.0,z_rand::Float64=10.0,r_rand::Float64=0.5)
    x0 = x_rand*rand(n_sphere)
    y0 = y_rand*rand(n_spheres)
    z0 = z_rand*rand(n_spheres)
    r = r_rand*ones(n_spheres)

    function constraints(x,u)::Array
        c = zeros(typeof(x[1]),n_spheres)

        for i = 1:n_spheres
            c[i] = sphere_constraint(x,x0[i],y0[i],z0[i],r[i])
        end
        c
    end
    constraints, (x0, y0, z0, r)
end

"""
$(SIGNATURES)
Plots a circle given a center and radius
"""
function plot_circle!(center,radius;kwargs...)
    pts = Plots.partialcircle(0,2π,100,radius)
    x,y = Plots.unzip(pts)
    x .+= center[1]
    y .+= center[2]
    plot!(Shape(x,y);kwargs...)
end

function plot_obstacles(circles,clr=:red)
    for circle in circles
        x,y,r = circle
        plot_circle!((x,y),r,color=clr,linecolor=clr,label="")
    end
end

function plot_trajectory!(X::AbstractMatrix;kwargs...)
    plot!(X[1,:],X[2,:];kwargs...)
end

function plot_trajectory!(res::TrajectoryOptimization.SolverVectorResults; kwargs...)
	plot_trajectory!(to_array(res.X); kwargs...)
end

function plot_trajectory!(res::DircolVars;kwargs...)
    plot_trajectory!(res.X; kwargs...)
end

"""
$(SIGNATURES)
    Skew-symmetric cross-product matrix (a.k.a. hat map)
"""
function hat(x)
      S = [  0   -x[3] x[2];
	  	    x[3]   0  -x[1];
		   -x[2]  x[1]  0]
end

"""
$(SIGNATURES)
    Convert quaternion (q = [s;v]) to rotation matrix
"""
function quat2rot(q)
      q = q./norm(q)
      s = q[1]
	  v = q[2:4]

      R = Matrix{Float64}(I,3,3) + 2*hat(v)*(hat(v) + s.*Matrix{Float64}(I,3,3))
end

function ispossemidef(A)
	eigs = eigvals(A)
	if any(real(eigs) .< 0)
		return false
	else
		return true
	end
end

function convergence_rate(stats::Dict;tail::Float64=0.5,plot_fit=false)
    total_iters = stats["iterations"]
    start = Int(ceil((1-tail)*total_iters))
    iters = collect(start:total_iters)
    grad = stats["gradient_norm"][iters]
    coeffs = convergence_rate(iters,grad)
    if plot_fit
        x = log.(1:total_iters)
        p = plot(x,log.(stats["gradient_norm"]),xlabel="iterations (log)",ylabel="gradient (log)")
        line = @. coeffs[2]*x + coeffs[1]
        plot!(x,line)
        ylim = collect(ylims(p))
        plot_vertical_lines!(log.(stats["outer_updates"]),ylim,linecolor=:black,linestyle=:dash,label="")
        display(p)
    end
    return coeffs[2]
end

function convergence_rate(x::Vector,y::Vector)
    n = length(x)
    X = [ones(n) log.(x)]
    coeffs = X\log.(y)
    return coeffs
end

function plot_vertical_lines!(p::Plots.Plot,x::Vector; kwargs...)
	ylim = collect(ylims(p))
	plot_vertical_lines!(x, ylim; kwargs...)
end

function plot_vertical_lines!(x,ylim=[-100,100]; kwargs...)
    ys = [ylim for val in x]
    xs = [[val; val] for val in x]
    plot!(xs,ys, linestyle=:dash, color=:black, label=""; kwargs...)
end

plot(X::Trajectory; kwargs...) = plot(to_array(X)'; kwargs...)
plot!(X::Trajectory; kwargs...) = plot!(to_array(X)'; kwargs...)
plot(X::Trajectory, inds::UnitRange; kwargs...) = plot(to_array(X)[inds,:]'; kwargs...)


text:
using CatViews
import Base.LinAlg: axpy!

# data params
m,n = 1000,500  # dimensions
r = 5           # rank
ξ = 1e-4        # noise level

# True loadings and components
Wtr = randn(m,r)
Ctr = randn(r,n)

# noisy, low-rank data
A = Wtr*Ctr + randn(m,n)*ξ

# random, initial guess for optimization
W = randn(m,r)
C = randn(r,n)

# storage for gradients
∇W = similar(W)
∇C = similar(C)

## PCA model ##
abstract Model
immutable PCA{T<:AbstractFloat} <: Model
    A::AbstractMatrix{T}
    W::AbstractMatrix{T}
    C::AbstractMatrix{T}
    ∇W::AbstractMatrix{T}
    ∇C::AbstractMatrix{T}
end
model = PCA(A,W,C,∇W,∇C)

value{T}(model::PCA{T}) = convert(T,0.5)*sumabs2(model.A - model.W*model.C)

function value_grad!{T}(m::PCA{T})
    resid = m.W*m.C - m.A
    At_mul_B!(m.∇C,m.W,resid)
    A_mul_Bt!(m.∇W,resid,m.C)
    return convert(T,0.5)*sumabs2(resid) # objective val
end

## Gradient Descent Optimization Type ##
abstract Optimizer
type GradDescent{T<:AbstractFloat} <: Optimizer
    ρ::T    # stepsize
    β::T    # backtracking linesearch
    function GradDescent(ρ::T,β::T)
        0 < β < 1 || error("β needs to be between 0 and 1.")
        ρ > 0 || error("ρ must be positive")
        new(ρ,β)
    end
end
GradDescent{T<:AbstractFloat}(ρ::T=1.0) = GradDescent{T}(ρ,convert(T,0.5))
opt = GradDescent()

# collect optimization variables into a vector
x = CatView(view(W,:),view(C,:))
∇x = CatView(view(∇W,:),view(∇C,:))

function update!(x::CatView,∇x::CatView,model::Model,opt::GradDescent)
    f = value_grad!(model)      # update ∇x
    axpy!(-opt.ρ,∇x,x)          # updates model params
    f_next = value(model)
    while f_next > f
        s = opt.ρ*opt.β         # smaller stepsize
        axpy!(opt.ρ-s,∇x,x)     # backtrack
        opt.ρ = s               # update stepsize
        f_next = value(model)
        opt.ρ < 1e-10 && break
    end
    return f_next
end

iterations = 100
f_trace = [ update!(x,∇x,model,opt) for iter in 1:iterations ]


text:
using JLD

import JLD.@write

ftp = "ftp://ssd.jpl.nasa.gov/pub/eph/planets/ascii"
denums = [
    "431",
    "430",
    "424",
    "423",
    "422",
    "421",
    "418",
    "414",
    "413",
    "410",
    "406",
    "405",
    "403",
    # DE202 is currently incomplete
    #= "202", =#
    "200",
    "102"
]
headers = Dict(zip(denums[3:end],["header.$d" for d in denums[3:end]]))
merge!(headers, Dict("431"=>"header.431_572","430"=>"header.430_572"))

datafiles = Dict("430"=>["ascp$i.430" for i = 1550:100:2550],
"431"=>[[@sprintf("ascm%05d.431", i) for i = 1000:1000:13000]; [@sprintf("ascp%05d.431", i) for i = 0:1000:16000]],
"424"=>[[@sprintf("ascm%04d.424", i) for i = 100:100:3000]; [@sprintf("ascp%04d.424", i) for i = 0:100:2900]],
"423"=>["ascp$i.423" for i = 1800:50:2150],
"422"=>[[@sprintf("ascm%04d.422", i) for i = 100:100:3000]; [@sprintf("ascp%04d.422", i) for i = 0:100:2900]],
"421"=>["ascp1900.421", "ascp2050.421"],
"418"=>["ascp1900.418"],
"414"=>["ascp$i.414" for i = 1600:100:2100],
"413"=>["ascp$i.413" for i = 1900:25:2025],
"410"=>["ascp$i.410" for i = 1960:20:2000],
"406"=>[[@sprintf("ascm%04d.406", i) for i = 100:100:3000]; [@sprintf("ascp%04d.406", i) for i = 0:100:2900]],
"405"=>["ascp$i.405" for i = 1600:20:2200],
"403"=>["ascp$i.403" for i = 1600:100:2100],
"202"=>["ascp$i.202" for i = 1900:50:2150],
"200"=>["ascp$i.200" for i = 1600:20:2160],
"102"=>[[@sprintf("ascm%04d.102", i) for i = 200:300:1400]; [@sprintf("ascp%04d.102", i) for i = 100:300:2800]])

function download(ftp, file, out)
    errorstring = "\nPlease download all files from '$ftp' manually, place them in '$PATH' and re-run 'getephem'.\n"
    try
        if OS_NAME == :Windows
            run(`where curl`)
        else
            run(`which curl`)
        end
    catch
        error("Curl could not be found. $errorstring")
    end

    try
        run(`curl $ftp/$file -o $out`)
    catch
        error("Curl could not connect to the server. $errorstring")
    end
end

function getephem(denum; force=false, debug=false)
    denum = string(denum)
    if ~haskey(datafiles, denum)
        error("Unknown ephemeris 'DE$denum'.")
    end
    outfile = "$PATH/de$denum.jld"
    if isfile(outfile) && force
        rm(outfile)
    end
    if ~isfile(outfile)
        testfile = "testpo.$denum"
        testlocal = "$PATH/$testfile"
        header = headers[denum]
        headerlocal = "$PATH/$header"
        data = datafiles[denum]
        datalocal = ["$PATH/$d" for d in data]
        println("Building ephemeris DE$denum.")
        if ~isfile("$PATH/$header") || force
            download("$ftp/de$denum", header, headerlocal)
        end
        if ~isfile("$PATH/$testfile") || force
            download("$ftp/de$denum", testfile, testlocal)
        end
        for (f,d) in zip(data, datalocal)
            if ~isfile(d) || force
                download("$ftp/de$denum", f, d)
            end
        end
        read_ascii(headerlocal, datalocal, outfile)
        if ~debug
            println("Removing ASCII files.")
            rm(headerlocal)
            for f in datalocal
                rm(f)
            end
        end
    end
end

function rmephem(denum)
    rm("$PATH/de$denum.jld")
    rm("$PATH/testpo.$denum")
end

function read_ascii(header, datafiles, outfile)
    startepoch = 0.0
    finalepoch = 0.0
    ncoeff = 0
    constantnames = AbstractString[]
    constantvalues = Float64[]
    ind = zeros(Int64,3,13)
    intervall = zeros(13)
    dtable = Array(Vector{Float64},13)
    for i = 1:13
        dtable[i] = Float64[]
    end

    println("Parsing header file.")
    l = open(readlines, header)
    ncoeff = parse(Int,(split(l[1])[4]))
    for i = 1:length(l)
        if startswith(l[i], "GROUP   1030")
            startepoch, finalepoch = float(split(l[i+2])[1:2])
        elseif startswith(l[i], "GROUP   1040")
            n = parse(Int,(l[i+2]))
            firstline = i+3
            lastline = i+3+div(n,10)
            for j = firstline:lastline
                append!(constantnames, split(l[j]))
            end
        elseif startswith(l[i], "GROUP   1041")
            n = parse(Int,(l[i+2]))
            firstline = i+3
            lastline = i+3+div(n,3)
            for j = firstline:lastline
                append!(constantvalues, float(split(replace(l[j],"D","e"))))
            end
        elseif startswith(l[i], "GROUP   1050")
            firstline = i+2
            lastline = i+4
            for j = firstline:lastline
                # The header of DE431 refers to 15 bodies instead of 13 as in the
                # other ephemerides. The file contains no coefficients for these
                # additional bodies though. Therefore only the first 13 indices
                # are used.
                ind[j-firstline+1,:] = [parse(Int,s) for s=split(l[j])[1:13]]
            end
        end
    end
    constants = Dict(zip(constantnames, constantvalues))

    #= file = jldopen(outfile, "w") =#
    try
        jldopen(outfile, "w") do file
            @write file constants
            @write file startepoch
            @write file finalepoch
            @write file ind
        end
        for i = 1:length(datafiles)
            println("Parsing data file $i/$(length(datafiles)).")
            parsedatafile(datafiles[i], ncoeff, ind, dtable, outfile)
        end
        for i = 1:length(dtable)
            sort!(dtable[i])
            intervall[i] = dtable[i][2] - dtable[i][1]
        end
        jldopen(outfile, "r+") do file
            @write file dtable
            @write file intervall
        end
    catch
        #= close(file) =#
        rm(outfile)
        rethrow()
    #= finally =#
    #=     close(file) =#
    end
end

function fromfortran(l)
    return float(split(replace(l, "D", "e")))
end

function parsedatafile(datafile, ncoeff, ind, dtable, outfile)
    coeff = Float64[]
    header = r"^\s+[0-9]+\s+[0-9]+"
    f = open(datafile, "r")
    seekend(f)
    filesize = Base.position(f)
    seekstart(f)
    while ~eof(f)
        l = readline(f)
        if ismatch(header, l)
            date, finaldate, c = fromfortran(readline(f))
            #println("Processing coefficients for $date-$finaldate.")
            push!(coeff, c)
            while ~eof(f)
                mark(f)
                l = readline(f)
                if ~ismatch(header, l)
                    append!(coeff, fromfortran(l))
                else
                    reset(f)
                    break
                end
            end
            savecoeff!(coeff, date, finaldate, ind, dtable, outfile)
        end
    end
end

function savecoeff!(coeff, date, finaldate, ind, dtable, outfile)
    for j = 1:13
        dt = (finaldate - date)/ind[3,j]
        if ~in(date, dtable[j])
            push!(dtable[j], date)
        elseif in(date, dtable[j])
            continue
        end
        if j == 12
            n = 2
        else
            n = 3
        end
        offset = ind[2,j]*n
        i1 = ind[1,j] - 2
        i2 = i1 + offset - 1
        jldopen(outfile, "r+") do file
            write(file, "$j-$date", reshape(coeff[i1:i2], ind[2,j], n))
        end
        if ind[3,j] != 1
            dt = (finaldate - date)/ind[3,j]
            for k = 1:ind[3,j]-1
                date1 = date + dt*k
                push!(dtable[j], date1)
                i1k = i1 + k*offset
                i2k = i1 + (k+1)*offset - 1
                jldopen(outfile, "r+") do file
                    write(file, "$j-$date1",
                        reshape(coeff[i1k:i2k], ind[2,j], n))
                end
            end
        end
    end
    empty!(coeff)
end


text:
module UnitParser

export destructure_units, reduce_units_expr

using YAML

const aliases_data = YAML.load(open("src/aliases.yaml"))
const prefixes = aliases_data["prefixes"]
const aliases = aliases_data["units"]

const r_integer = r"(-?[0-9]+)"
const r_unit_name = r"([A-Za-zμΩ]+)"
const r_power = r" *(?:\*\*|\^)? *"
const r_divide = r" *(?:/|per|PER) *"
const r_multiply = r" *[\*\. ] *"

"""
    consume(string, m::RegexMatch)

Remove the matching part of the RegexMatch and return the shorter string.
"""
consume(string, m::RegexMatch) = string[length(m.match)+1:end]
consume(string, ::Nothing) = string

match_start(pattern, string) = match(r"^" * pattern, string)

function reduce_units_expr(string::AbstractString)
    factors = []
    for (prefix, core, exponent) in destructure_units(string)
        if exponent != 1
            factor = "$prefix$core^$exponent"
        else
            factor = "$prefix$core"
        end
        push!(factors, factor)
    end
    return join(factors, "*")
end

"""
    destructure_units(string)

Deconstruct a string representing units into a list of (prefix, name, exponent),
one for each factor.
"""
function destructure_units(string::AbstractString)
    string, factors = destructure_factors(string)
    string = consume(string, match_start(r_divide, string))
    string, divisors = destructure_factors(string)

    if length(string) > 0
        @error "The input string could not be fully parsed"
    end

    divisors = map(divisors) do (prefix, symbol, exponent)
        (prefix, symbol, -exponent)
    end

    append!(factors, divisors)
    return factors
end

function destructure_factors(string::AbstractString)
    factors = []
    while (match_start(r_divide, string) === nothing &&
          (m = match_start(r_unit_name, string)) !== nothing)
        unit_name = m.match
        string = consume(string, m)

        if (m = match_start(r_power * r_integer, string)) !== nothing
            string = consume(string, m)
            exponent = parse(Int, m.captures[1])
        else
            exponent = 1
        end

        push!(factors, (short_form(unit_name)..., exponent))

        string = consume(string, match_start(r_multiply, string))
    end

    return string, factors
end

"""
    short_form(units_name::AstractString)

Reduce the natural name of a units to its symbolic form for prefix and name.

If no match is found, return the original string.

Example
=======
julia> short_form("micrometers")
("μ", "m")

Return
======
unit_prefix: String
    Short form representing the prefix of the units, e.g. "k" for "kilo".
    Set to empty string if the units had no prefix.

core_unit: String
    Short form representing the units, e.g. "A" for "ampere".
"""
function short_form(units_name::AbstractString)
    units_name = lowercase(units_name)

    # Remove trailing 's' to account for pluralized units
    if units_name[end] == 's'
        units_name = units_name[1:end-1]
    end

    unit_prefix = ""

    for (prefix, symbol) in prefixes
        n = length(prefix)
        if length(units_name) >= n && units_name[1:n] == prefix
            unit_prefix = symbol
            units_name = units_name[n+1:end]
            break
        end
    end

    core_unit = units_name

    for (alias, symbol) in aliases
        if units_name == alias
            core_unit = symbol
            break
        end
    end

    return unit_prefix, core_unit
end

end # module


text:
using SignatureGB
using Singular
P, (x1, x2, x3, x4, x5, x6, x7, x8, x9, x10)  = Singular.PolynomialRing(Fp(65521), ["x1", "x2", "x3", "x4", "x5", "x6", "x7", "x8", "x9", "x10"]);
gb = SignatureGB.f5([x1, x2, x3, x4, x5, x6, x7, x8, x9, x10]);F = [8179*x1^4-1134*x1^3*x7-810*x1^2*x10^2-10710*x1^2*x10*x7+2714*x1^2*x10*x8-3618*x1^2*x2*x7+6804*x1^2*x2*x9-1242*x1^2*x3*x5+10206*x1^2*x3*x6-4104*x1^2*x3*x7+36*x1^2*x5^2+6678*x1^2*x6^2+4956*x1^2*x6*x8+441*x1^2*x7^2+630*x1*x10^2*x7+2814*x1*x2*x7^2+996*x1*x3*x5*x6+402*x1*x3*x5*x7+3192*x1*x3*x7^2+720*x1*x5*x9^2+225*x10^4+2010*x10^2*x2*x7+690*x10^2*x3*x5+2280*x10^2*x3*x7+7229*x10^2*x7^2+529*x10^2*x8^2+1369*x10^2*x9^2+4662*x10*x2*x4*x9+1332*x10*x2*x5*x9-9180*x10*x2*x7*x9-13770*x10*x3*x6*x7+44*x10*x3*x7^2+116*x10*x3*x7*x9-40*x10*x4*x7^2+4514*x10*x4*x8*x9+1406*x10*x5^2*x9-9010*x10*x6^2*x7+1932*x10*x6*x8^2+3969*x2^2*x4^2+2268*x2^2*x4*x5+324*x2^2*x5^2+4489*x2^2*x7^2+2916*x2^2*x9^2+3082*x2*x3*x5*x7+8748*x2*x3*x6*x9+10184*x2*x3*x7^2+7686*x2*x4^2*x8+2394*x2*x4*x5^2+2196*x2*x4*x5*x8+684*x2*x5^3+5724*x2*x6^2*x9+529*x3^2*x5^2+3496*x3^2*x5*x7+13450*x3^2*x6^2-7802*x3^2*x6*x7+8106*x3^2*x7^2+638*x3^2*x7*x9+841*x3^2*x9^2-220*x3*x4*x7^2-580*x3*x4*x7*x9+8586*x3*x6^3+9960*x3*x6*x9^2-5640*x3*x7*x9^2+100*x4^2*x7^2+3721*x4^2*x8^2+2318*x4*x5^2*x8+361*x5^4+2809*x6^4+1764*x6^2*x8^2+3600*x9^4-9828*x1^2*x10+10738*x1^2*x3+10856*x1^2*x4-576*x1^2*x5+9794*x1^2*x9-72*x1*x10*x7-7968*x1*x3*x6+4116*x1*x3*x7-1044*x1*x3*x9+360*x1*x4*x7+120*x1*x5*x7-5760*x1*x9^2+13260*x10^2*x7-8424*x10*x2*x9-12636*x10*x3*x6+4186*x10*x3*x8+4232*x10*x4*x8-8268*x10*x6^2+222*x10*x6*x9-136*x10*x7*x8+3818*x10*x8*x9+378*x2*x4*x6+108*x2*x5*x6+1660*x3*x6*x7+7644*x3*x6*x8-940*x3*x7^2-748*x3*x7*x8-1972*x3*x8*x9+8094*x4*x6*x8+680*x4*x7*x8+114*x5^2*x6+6972*x6*x8*x9+1200*x7*x9^2+2628*x1^2-960*x1*x7+1224*x1*x8+6084*x10^2+8281*x3^2+16744*x3*x4+15106*x3*x9+8464*x4^2+15272*x4*x9+9*x6^2+100*x7^2+1156*x8^2+6889*x9^2, -3618*x1^2*x7+6804*x1^2*x9+2814*x1*x7^2+2010*x10^2*x7+4662*x10*x4*x9+1332*x10*x5*x9-9180*x10*x7*x9+7938*x2*x4^2+4536*x2*x4*x5+648*x2*x5^2+8978*x2*x7^2+5832*x2*x9^2+3082*x3*x5*x7+8748*x3*x6*x9+10184*x3*x7^2+7686*x4^2*x8+2394*x4*x5^2+2196*x4*x5*x8+684*x5^3+5724*x6^2*x9-8424*x10*x9+378*x4*x6+108*x5*x6, -1242*x1^2*x5+10206*x1^2*x6-4104*x1^2*x7+996*x1*x5*x6+402*x1*x5*x7+3192*x1*x7^2+690*x10^2*x5+2280*x10^2*x7-13770*x10*x6*x7+44*x10*x7^2+116*x10*x7*x9+3082*x2*x5*x7+8748*x2*x6*x9+10184*x2*x7^2+1058*x3*x5^2+6992*x3*x5*x7+26900*x3*x6^2-15604*x3*x6*x7+16212*x3*x7^2+1276*x3*x7*x9+1682*x3*x9^2-220*x4*x7^2-580*x4*x7*x9+8586*x6^3+9960*x6*x9^2-5640*x7*x9^2+10738*x1^2-7968*x1*x6+4116*x1*x7-1044*x1*x9-12636*x10*x6+4186*x10*x8+1660*x6*x7+7644*x6*x8-940*x7^2-748*x7*x8-1972*x8*x9+16562*x3+16744*x4+15106*x9, 4662*x10*x2*x9-40*x10*x7^2+4514*x10*x8*x9+7938*x2^2*x4+2268*x2^2*x5+15372*x2*x4*x8+2394*x2*x5^2+2196*x2*x5*x8-220*x3*x7^2-580*x3*x7*x9+200*x4*x7^2+7442*x4*x8^2+2318*x5^2*x8+10856*x1^2+360*x1*x7+4232*x10*x8+378*x2*x6+8094*x6*x8+680*x7*x8+16744*x3+16928*x4+15272*x9, -1242*x1^2*x3+72*x1^2*x5+996*x1*x3*x6+402*x1*x3*x7+720*x1*x9^2+690*x10^2*x3+1332*x10*x2*x9+2812*x10*x5*x9+2268*x2^2*x4+648*x2^2*x5+3082*x2*x3*x7+4788*x2*x4*x5+2196*x2*x4*x8+2052*x2*x5^2+1058*x3^2*x5+3496*x3^2*x7+4636*x4*x5*x8+1444*x5^3-576*x1^2+120*x1*x7+108*x2*x6+228*x5*x6, 10206*x1^2*x3+13356*x1^2*x6+4956*x1^2*x8+996*x1*x3*x5-13770*x10*x3*x7-18020*x10*x6*x7+1932*x10*x8^2+8748*x2*x3*x9+11448*x2*x6*x9+26900*x3^2*x6-7802*x3^2*x7+25758*x3*x6^2+9960*x3*x9^2+11236*x6^3+3528*x6*x8^2-7968*x1*x3-12636*x10*x3-16536*x10*x6+222*x10*x9+378*x2*x4+108*x2*x5+1660*x3*x7+7644*x3*x8+8094*x4*x8+114*x5^2+6972*x8*x9+18*x6, -1134*x1^3-10710*x1^2*x10-3618*x1^2*x2-4104*x1^2*x3+882*x1^2*x7+630*x1*x10^2+5628*x1*x2*x7+402*x1*x3*x5+6384*x1*x3*x7+2010*x10^2*x2+2280*x10^2*x3+14458*x10^2*x7-9180*x10*x2*x9-13770*x10*x3*x6+88*x10*x3*x7+116*x10*x3*x9-80*x10*x4*x7-9010*x10*x6^2+8978*x2^2*x7+3082*x2*x3*x5+20368*x2*x3*x7+3496*x3^2*x5-7802*x3^2*x6+16212*x3^2*x7+638*x3^2*x9-440*x3*x4*x7-580*x3*x4*x9-5640*x3*x9^2+200*x4^2*x7-72*x1*x10+4116*x1*x3+360*x1*x4+120*x1*x5+13260*x10^2-136*x10*x8+1660*x3*x6-1880*x3*x7-748*x3*x8+680*x4*x8+1200*x9^2-960*x1+200*x7, 2714*x1^2*x10+4956*x1^2*x6+1058*x10^2*x8+4514*x10*x4*x9+3864*x10*x6*x8+7686*x2*x4^2+2196*x2*x4*x5+7442*x4^2*x8+2318*x4*x5^2+3528*x6^2*x8+4186*x10*x3+4232*x10*x4-136*x10*x7+3818*x10*x9+7644*x3*x6-748*x3*x7-1972*x3*x9+8094*x4*x6+680*x4*x7+6972*x6*x9+1224*x1+2312*x8, 6804*x1^2*x2+1440*x1*x5*x9+2738*x10^2*x9+4662*x10*x2*x4+1332*x10*x2*x5-9180*x10*x2*x7+116*x10*x3*x7+4514*x10*x4*x8+1406*x10*x5^2+5832*x2^2*x9+8748*x2*x3*x6+5724*x2*x6^2+638*x3^2*x7+1682*x3^2*x9-580*x3*x4*x7+19920*x3*x6*x9-11280*x3*x7*x9+14400*x9^3+9794*x1^2-1044*x1*x3-11520*x1*x9-8424*x10*x2+222*x10*x6+3818*x10*x8-1972*x3*x8+6972*x6*x8+2400*x7*x9+15106*x3+15272*x4+13778*x9, -1620*x1^2*x10-10710*x1^2*x7+2714*x1^2*x8+1260*x1*x10*x7+900*x10^3+4020*x10*x2*x7+1380*x10*x3*x5+4560*x10*x3*x7+14458*x10*x7^2+1058*x10*x8^2+2738*x10*x9^2+4662*x2*x4*x9+1332*x2*x5*x9-9180*x2*x7*x9-13770*x3*x6*x7+44*x3*x7^2+116*x3*x7*x9-40*x4*x7^2+4514*x4*x8*x9+1406*x5^2*x9-9010*x6^2*x7+1932*x6*x8^2-9828*x1^2-72*x1*x7+26520*x10*x7-8424*x2*x9-12636*x3*x6+4186*x3*x8+4232*x4*x8-8268*x6^2+222*x6*x9-136*x7*x8+3818*x8*x9+12168*x10];
sort!(F, by = p -> leading_monomial(p))
println(@elapsed SignatureGB.f5(F, verbose=true));exit()


text:
# Autogenerated wrapper script for samtools_jll for powerpc64le-linux-gnu
export samtools

using Ncurses_jll
using htslib_jll
JLLWrappers.@generate_wrapper_header("samtools")
JLLWrappers.@declare_executable_product(samtools)
function __init__()
    JLLWrappers.@generate_init_header(Ncurses_jll, htslib_jll)
    JLLWrappers.@init_executable_product(
        samtools,
        "bin/samtools",
    )

    JLLWrappers.@generate_init_footer()
end  # __init__()