text:
@testset "Agent" begin
    @testset "Default selection" begin
        agent = get_default_agent()
        @test agent !== nothing
        @test AMDGPU.device_type(agent) == :gpu

        agent_name = AMDGPU.get_name(agent)
        @test length(agent_name) > 0
        @test !occursin('\0', agent_name)

        if length(AMDGPU.get_agents()) > 1
            @testset "Multi-GPU" begin
                init_agent = AMDGPU.get_default_agent()
                init_dev = AMDGPU.device()
                @test init_dev == 1
                AMDGPU.device!(2)
                @test AMDGPU.device() == 2
                @test AMDGPU.get_default_agent() != init_agent
                AMDGPU.device!(1)
                @test AMDGPU.device() == 1
                @test AMDGPU.get_default_agent() == init_agent

                @test_throws BoundsError AMDGPU.device!(0)
                @test_throws BoundsError AMDGPU.device!(length(AMDGPU.get_agents())+1)
            end
        else
            @test_skip "Multi-GPU"
        end
    end

    @testset "ISAs" begin
        agent = get_default_agent()
        agent_isa = get_first_isa_string(agent)
        @test length(agent_isa) > 0
    end
end


text:
@doc """
  An axis aligned hyperbox - a set of intervals
  """ ->
type HyperBox{T} <: Domain{T}
  intervals::Vector{Interval{T}}
end

@doc doc"""A LazyBox is a hyperbox whose dimensions can be
  created dynamically and do not have to be contiguous""" ->
immutable LazyBox{T} <: Domain{T}
  intervals::Dict{Int64,Interval{T}}
end

@doc "Boxes are abstractions of Vector{T<:Real}" ->
typealias Boxes{T} Union{LazyBox{T},HyperBox{T}}

## Domain Operations
## =================
isrelational{T<:Boxes}(::Type{T}) = false
domaineq(x::Boxes,y::Boxes) = isequal(x,y)
function isequal(x::Boxes,y::Boxes)
  for dim in union(dims(x),dims(y))
    if !(hasdim(x,dim) && hasdim(y,dim) && isequal(x[dim],y[dim]))
      return false
    end
  end
  true
end

## Box Domain General functions
## ============================
function product_boxes{T}(b::Boxes, split_point::Dict{Int, T})
  @assert(length(keys(split_point)) <= ndims(b))

  splits = Vector{Interval{T}}[]
  for dim in dims(b)
    if haskey(split_point, dim)
      push!(splits, split_box(b[dim],split_point[dim]))
    else # Don't split along this dimension
      push!(splits, [b[dim]])
    end
  end
  prodsubboxes(splits)
end

# Find cartesian product of a bunch of subboxes (vector of intervals)
function prodsubboxes{T}(splits::Vector{Vector{Interval{T}}})
  @compat boxes = Tuple{Vararg{Interval}}[]
  for subbox in product(splits...)
    push!(boxes, subbox) # product returns tuple
  end
  return boxes
end

mid{T}(b::Boxes{T}) = T[mid(b[dim]) for dim in dims(b)]

@doc "Split box into 2^d equally sized boxes by cutting down middle of each axis" ->
mid_split(b::Boxes) = split_box(b, mid(b))

@doc "Do a partial split at the midpoints of dimensions `dims`" ->
mid_partial_split(b::Boxes, partial_dims::Vector{Int}) =
  @compat partial_split_box(b,Dict([dim => mid(b[dim]) for dim in partial_dims]))

## Sampling
## ========
rand(b::Boxes) = [rand(b[dim]) for dim in dims(b)]

include("boxes/hyperbox.jl")
include("boxes/lazybox.jl")

ext: jl


text:
using DifferentialEquations
using PyCall
# pygui(:qt)
using PyPlot
# pygui(true)

include("leap_frog_definitions.jl")
include("PSS_Definitions_Dimer_X.jl")

function plot_orbits_leapfrog_aref(H, max_time, q2, p2)
    h = 1 / (2 * H)
    a = (2 + h - 2 * sqrt(h + 1)) / h

    q1 = 0
    p1 = P1_find_dimer(q2, p2, H)
    if isempty(p1)
        p1 = P1_find_dimer_second(q2, p2, H)
    end

    Q1 = (q1 + q2) / sqrt(2)
    Q2 = (q1 - q2) / sqrt(2)
    P1 = (p1[1] + p2) / sqrt(2)
    P2 = (p1[1] - p2) / sqrt(2);

    u0 = [Q1;P2;Q2;P1;0;0] # X P Q Y
    Q0 = [Q1,Q2]; P0 = [P1,P2];
    tspan = (0.0, max_time);

    prob = ODEProblem(Eq_of_M_LAB_FRAME, u0, tspan)
    t, A = solve(prob, Vern9(), reltol = 1e-13, abstol = 1e-15, maxiters = 1e15);
    N = length(t.t)
    Q1 = A[1,:]
    Q2 = A[3,:]
    P1 = A[4,:]
    P2 = A[2,:]

    q1 = (Q1 + Q2) / sqrt(2)
    q2 = (Q1 - Q2) / sqrt(2)
    p1 = (P1 + P2) / sqrt(2)
    p2 = (P1 - P2) / sqrt(2);

    figure(dpi = 250)
    plot(Q1, P1, "b", linewidth = .1)
    plot(Q2, P2, "r", linewidth = .1)
    title("Aref coordinates");
end

function plot_orbits_leapfrog_dimer(H, max_time, q2, p2)
    h = 1 / (2 * H)
    a = (2 + h - 2 * sqrt(h + 1)) / h

    q1 = 0
    p1 = P1_find_dimer(q2, p2, H)
    if isempty(p1)
        p1 = P1_find_dimer_second(q2, p2, H)
    end

    Q1 = (q1 + q2) / sqrt(2)
    Q2 = (q1 - q2) / sqrt(2)
    P1 = (p1[1] + p2) / sqrt(2)
    P2 = (p1[1] - p2) / sqrt(2);

    u0 = [Q1;P2;Q2;P1;0;0] # X P Q Y
    Q0 = [Q1,Q2]; P0 = [P1,P2];
    tspan = (0.0, max_time);

    prob = ODEProblem(Eq_of_M_LAB_FRAME, u0, tspan)
    t, A = solve(prob, Vern9(), reltol = 1e-13, abstol = 1e-15, maxiters = 1e15);
    N = length(t.t)
    Q1 = A[1,:]
    Q2 = A[3,:]
    P1 = A[4,:]
    P2 = A[2,:]

    q1 = (Q1 + Q2) / sqrt(2)
    q2 = (Q1 - Q2) / sqrt(2)
    p1 = (P1 + P2) / sqrt(2)
    p2 = (P1 - P2) / sqrt(2);

    figure(dpi = 250)
    plot(q1, p1, "b", linewidth = .1)
    plot(q2, p2, "r", linewidth = .1)
    title("Dimer coordinates");
end

function plot_orbits_leapfrog_lab(H, max_time, q2, p2)
    h = 1 / (2 * H)
    a = (2 + h - 2 * sqrt(h + 1)) / h

    q1 = 0
    p1 = P1_find_dimer(q2, p2, H)
    if isempty(p1)
        p1 = P1_find_dimer_second(q2, p2, H)
    end

    Q1 = (q1 + q2) / sqrt(2)
    Q2 = (q1 - q2) / sqrt(2)
    P1 = (p1[1] + p2) / sqrt(2)
    P2 = (p1[1] - p2) / sqrt(2);

    u0 = [Q1;P2;Q2;P1;0;0] # X P Q Y
    Q0 = [Q1,Q2]; P0 = [P1,P2];
    tspan = (0.0, max_time);

    prob = ODEProblem(Eq_of_M_LAB_FRAME, u0, tspan)
    t, A = solve(prob, Vern9(), reltol = 1e-13, abstol = 1e-15, maxiters = 1e15);
    N = length(t.t)
    Centeroid = A[5,:] + im * A[6,:];
    Z = A[1,:] + im * A[2,:];
    W = A[3,:] + im * A[4,:];
    lin_impulse = im * (1 + a ) * ones(N);
    z_1pos = .5 * (Centeroid + Z + lin_impulse + W);
    z_1neg = .5 * (Centeroid + Z - lin_impulse - W);
    z_2pos = .5 * (Centeroid - Z + lin_impulse - W);
    z_2neg = .5 * (Centeroid - Z - lin_impulse + W);

    z_1posX = real(z_1pos)
    z_1posY = imag(z_1pos)

    z_2posX = real(z_2pos)
    z_2posY = imag(z_2pos)

    z_1negX = real(z_1neg)
    z_1negY = imag(z_1neg)

    z_2negX = real(z_2neg)
    z_2negY = imag(z_2neg)

    figure(dpi = 250)
    plot(z_1posX, z_1posY, c = "b", linewidth = 1)
    plot(z_1negX, z_1negY, c = "r", linewidth = 1)
    plot(z_2posX, z_2posY, c = "b", linewidth = .5, linestyle = "--")
    plot(z_2negX, z_2negY, c = "r", linewidth = .5, linestyle = "--")
    ylim(-3.2, 3.2)
    axis("equal")
    title("Lab coordinates");
end


text:
module Octrees

using FMMTrees
using StaticArrays

struct Data{T}
    sector::Int
    values::Vector{T}
end

"""
An octree behaves like a pointer based tree, specifically designed for the
storage of data residing in Euclidian space. The routing and insertion methods
maintain a Hilbert space filling curve order.

This makes the Octree ideal as the basis for ordering algorithms aimed at
preserving data locality (in the context of e.g. parallellism).
"""
struct Octree{T}
    pbtree::FMMTrees.PointerBasedTrees.PointerBasedTree{FMMTrees.PointerBasedTrees.Node{Data{T}}}
    function Octree(value::Vector{T}) where {T}
        N = FMMTrees.PointerBasedTrees.Node{Data{T}}
        new{T}(
            FMMTrees.PointerBasedTrees.PointerBasedTree(
                N[N(Data(0,value),0,0,0,0)], 1))
    end
end

Base.length(itr::FMMTrees.DepthFirstIterator{T}) where {T<:Octree} = length(itr.tree.pbtree.nodes)

# Tree API
FMMTrees.data(tree::Octree, node) = data(tree.pbtree, node)
FMMTrees.children(tree::Octree, node) = children(tree.pbtree, node)
FMMTrees.insert!(tree::Octree, data; next, parent, prev) = FMMTrees.insert!(tree.pbtree, data, next=next, parent=parent, prev=prev)
FMMTrees.root(tree::Octree) = root(tree.pbtree)

# PointerBasedTrees API
FMMTrees.PointerBasedTrees.getnode(tree::Octree, node)  = FMMTrees.PointerBasedTrees.getnode(tree.pbtree, node)

const hilbert_states = [
    [1, 2, 3, 2, 4, 5, 3, 5],
    [2, 6, 0, 7, 8, 8, 0, 7],
    [0, 9,10, 9, 1, 1,11,11],
    [6, 0, 6,11, 9, 0, 9, 8],
    [11,11, 0, 7, 5, 9, 0, 7],
    [4, 4, 8, 8, 0, 6,10, 6],
    [5, 7, 5, 3, 1, 1,11,11],
    [6, 1, 6,10, 9, 4, 9,10],
    [10, 3, 1, 1,10, 3, 5, 9],
    [4, 4, 8, 8, 2, 7, 2, 3],
    [7, 2,11, 2, 7, 5, 8, 5],
    [10, 3, 2, 6,10, 3, 4, 4]]

const hilbert_positions = [
    [0,1,3,2,7,6,4,5],
    [0,7,1,6,3,4,2,5],
    [0,3,7,4,1,2,6,5],
    [2,3,1,0,5,4,6,7],
    [4,3,5,2,7,0,6,1],
    [6,5,1,2,7,4,0,3],
    [4,7,3,0,5,6,2,1],
    [6,7,5,4,1,0,2,3],
    [2,5,3,4,1,6,0,7],
    [2,1,5,6,3,0,4,7],
    [4,5,7,6,3,2,0,1],
    [6,1,7,0,5,2,4,3]]

"""
Return the center of the child box and its sector based on the relative
positioning of the parent center `ct` and the target point `pt`. Descending
a level will half the *half size* `hs`.
"""
function sector_center_size(pt, ct, hs)
    hs = hs / 2
    bl = pt .> ct
    ct = ifelse.(bl, ct.+hs, ct.-hs)
    sc = sum(b ? 2^(i-1) : 0 for (i,b) in enumerate(bl))
    return sc, ct, hs
end

struct Router{T,P}
    smallest_box_size::T
    target_point::P
end

function route!(tree::Octree, state, router)

    point = router.target_point
    smallest_box_size = router.smallest_box_size

    node_idx, center, size, sfc_state = state
    size <= smallest_box_size && return state
    target_sector, target_center, target_size = sector_center_size(point, center, size)
    target_pos = hilbert_positions[sfc_state][target_sector+1] + 1
    target_sfc_state = hilbert_states[sfc_state][target_sector+1] + 1
    prev_child, next_child = 0, 0
    for child in FMMTrees.children(tree, node_idx)
        child_sector = FMMTrees.data(tree,child).sector
        child_pos = hilbert_positions[sfc_state][child_sector+1]+1
        target_pos < child_pos  && (next_child = child; break)
        if child_sector == target_sector
            return (child, target_center, target_size, target_sfc_state)
        end
        prev_child = child
    end
    data = Data(target_sector, Int[])
    new_node_idx = FMMTrees.insert!(tree, data, next=next_child, prev=prev_child, parent=node_idx)
    return new_node_idx, target_center, target_size, target_sfc_state
end

# function updater!(tree, (node, center, size), data)
#     push!(FMMTrees.data(tree, node).values, data)
# end

end # module


text:
module GAlgebra

using PyCall

import Base: @pure

export galgebra

const galgebra = PyCall.PyNULL()
const metric = PyCall.PyNULL()
const ga = PyCall.PyNULL()
const mv = PyCall.PyNULL()
const lt = PyCall.PyNULL()
const printer = PyCall.PyNULL()

include("macros.jl")
include("ga.jl")
include("mv.jl")

function __init__()
    copy!(galgebra, PyCall.pyimport_conda("galgebra", "galgebra"))
    copy!(metric, PyCall.pyimport_conda("galgebra.metric", "galgebra"))
    copy!(ga, PyCall.pyimport_conda("galgebra.ga", "galgebra"))
    copy!(mv, PyCall.pyimport_conda("galgebra.mv", "galgebra"))
    copy!(lt, PyCall.pyimport_conda("galgebra.lt", "galgebra"))
    copy!(printer, PyCall.pyimport_conda("galgebra.printer", "galgebra"))

    pytype_mapping(galgebra.mv.Mv, Mv)
end

end # module


text:
module TableDirectories
using Glob
using DataFrames

PatternElement = Union{Glob.GlobMatch,Glob.FilenameMatch,Regex,Base.String}
PatternArray = Array{PatternElement}
PatternType = Union{Glob.GlobMatch,Glob.FilenameMatch,Regex,Base.String,Array{PatternElement,1}}

"""
TableRoot(path, ::PatternType)

Defines the top directory and the file filter that will contain the rows of the table.
Example: TableRoot("/sys/bus/usb/devices", r"^[-0-9.]+") on a Linux system will generate a large table of USB device attributes.
The default regex is r".+" (i.e. all files in top directory) which for usb devices will also include end points, which have quite different attributes.
"""
struct TableRoot
    path::String
    filter::PatternType
    interesting_columns::Vector{String} # default to unique?
    uninteresting_columns::Vector{String} # default to constant?
end

TableRoot(path) = TableRoot(path, r".+", [], [])
TableRoot(path, filter) = TableRoot(path, filter, [], [])

"""
rows(table_root)

Returns the rows of the table.
Example: rows(TableRoot("/sys/bus/usb/devices", r"^[-0-9.]+")) on a Linux system will return a Vector of RowRoot objects one for each USB device.
"""
function rows(table_root::TableRoot)::Vector{RowRoot}
    map(glob([table_root.filter], table_root.path)) do row_filename
        RowRoot(table_root, row_filename[length(table_root.path)+2:end])
        #RowRoot(table_root, row_filename[1:end])
    end
end # rows

struct RowRoot
    table_root::TableRoot
    row_filename::String
end

function path(row_root::RowRoot)::String
    joinpath(row_root.table_root.path, row_root.row_filename)
end # path

function paths(row_root::RowRoot)::Vector{String}
    readdir(path(row_root); join = true)
end

struct Cell
    row_root::RowRoot
    subpath::String
end

function path(cell::Cell)::String
    joinpath(path(cell.row_root), cell.subpath)
end # Cell.path

function paths(cell::Cell)::Vector{String}
    readdir(path(cell); join = true)
end

"""
cell_value(cell)

Given a cell object returns its value
If file is a normal data file the value is the contents of the file.
If the file is not readable, the value is the string "not readable".
If the file is a link the value is the absolute path linked to.
If the file does not exist missing is returned. This was perferred over nothing since describe(df, :nmissing) will count them.
Generally the file does not exist because it is not meaningful or not available such as the power consumption of a virtual device.
"""
function cell_value(cell::Cell)
    cell_path = path(cell)
    if islink(cell_path)
        realpath(cell_path)
    else
        if isfile(cell_path)
            try
                chomp(read(cell_path, String))
            catch error
                "not readable"
            end
        else
            if isdir(cell_path)
                sub_directory_files = readdir(cell_path)
                if sub_directory_files == []
                    missing
                else
                    cell_path * " is a directory containing: " * join(sub_directory_files, ", ")
                end
            else
                missing
            end
        end
    end
end # cell_value

function subpath(parent_cell::Cell, added_filename::String)::String
    parent_column.subpath * "/" * added_filename
end # subpath

function column_name_from_path(cell_path::String, row_root::RowRoot)::Cell
    sub_path = cell_path[length(path(row_root))+1:end]
    Cell(row_root, sub_path)
end # column_name_from_path

function find_all(root::String, filename::String)::Vector{Union{Missing,String}}
    map(readdir(root)) do path
        if filename == basename(path)
            path
        elseif isdir(path)
            find_all(path, filename)
        end
    end
end

"""
column_names(::RowRoot, column_path)

Returns Vector of column names of row.
Example: column_names(TableRoot("/sys/bus/usb/devices", r"^[-0-9.]+")) on a Linux system will return a Vector of column names one for each USB device.
If you want the union of all row's column_names use: names(data_frame())
"""

function column_names(row_root::RowRoot, path = "")::Vector{String}
    columns_dir = joinpath(TableDirectories.path(row_root), path)
    ret = map(readdir(columns_dir)) do sub_column
        sub_column_path = joinpath(columns_dir, sub_column)
        if islink(sub_column_path)
            joinpath(path, sub_column)
        else
            if isfile(sub_column_path)
                joinpath(path, sub_column)
            else
                if isdir(sub_column_path)
                    column_names(row_root, joinpath(path, sub_column))
                else
                    missing
                end
            end
        end
    end
    vcat(ret...)
end # column_names

function dict(row_root::RowRoot, column_names::Vector{String})
    ret = Dict()
    ret["row_path"] = row_root.row_filename
    for column_subpath in column_names
        cell = Cell(row_root, column_subpath)
        ret[column_subpath] = cell_value(cell)
    end
    ret
end # dict
"""
data_frame(table_root, [column_names])

Return a DataFrame
"""
function data_frame(table_root::TableRoot, column_names::Vector{String})
    data_frame = DataFrame()
    for row_filename in rows(table_root)
        append!(data_frame, dict(row_filename, column_names), cols = :union)
    end
    data_frame
end # data_frame

data_frame(table_root) = data_frame(table_root, column_names(rows(table_root)[1], ""))

end # TableDirectories


text:
using Documenter
using SpectralClustering

makedocs(
    modules = [SpectralClustering],
    format = Documenter.HTML(prettyurls = true),
    source = "src",
    clean = false,

    sitename = "SpectralClustering.jl",
    pages = Any[
        "Home" => "index.md",
        "Getting Started"=>"start.md",
         "Main Modules" => Any[
                       "Graph Creation" => "man/graphcreation.md",
                       "Embedding"      => "man/embedding.md",
                       "Approximate Embedding" => "man/approximate.md",
                       "Eigenvector Clustering" => "man/clusterize.md",
                       "Co-Regularized" => "man/multiview.md",
                       "Incremental"    => "man/incremental.md"
         ],
        "Utility Modules" => Any[
                       "Data Access" => "man/data_access.md",
                       "Graph" => "man/graph.md",
                       "Landmarks Selection" => "man/landmark_selection.md"
        ],
   ],
  doctest = false

)
notebook_output_dir =   joinpath(dirname(@__FILE__), "build","notebooks")
mkpath(notebook_output_dir)
using IJulia
jupyter_path = first(IJulia.find_jupyter_subcommand("nbconvert"))
for file in readdir(joinpath(dirname(@__FILE__), "notebooks"))
   full_path = joinpath(dirname(@__FILE__), "notebooks", file)
    if (endswith(file,".ipynb"))
	    run(`$(jupyter_path) nbconvert --to html $full_path --output-dir=$notebook_output_dir`)
    elseif (file != ".ipynb_checkpoints")
        cp(full_path, joinpath(notebook_output_dir,file),  force=true)
    end
end

deploydocs(
    repo = "github.com/lucianolorenti/SpectralClustering.jl.git",
    deps = nothing,
    make = nothing,
    target = "build"
)


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.8.0-DEV
#     language: julia
#     name: julia-1.8
# ---

# %%
using Plots
using Random
Random.seed!(4649373)

n = 10
x, y = randn(n+1), rand(n+1)
u, v = diff(x), diff(y)
quiver(x[1:n], y[1:n]; quiver=(u, v))

# %%
plot(x, y; arrow=true)

# %%
xs = ys = -1:0.2:1
X, Y = [x for y in ys, x in xs], [y for y in ys, x in xs]
U, V = -0.3Y, 0.3X
quiver(vec(X-U/2), vec(Y-V/2); quiver=(vec(U), vec(V)), size=(420, 400))

# %%

max_stars_repo_path: 0007/quiver example.jl

max_stars_repo_name: genkuroki/public

max_issues_repo_path: 0007/quiver example.jl

max_issues_repo_name: genkuroki/public

max_forks_repo_path: 0007/quiver example.jl

max_forks_repo_name: genkuroki/public


text:
# This file is part of IntegerSequences.
# Copyright Peter Luschny. License is MIT.

(@__DIR__) ∉ LOAD_PATH && push!(LOAD_PATH, (@__DIR__))

module SeqUtils

using Nemo

export ModuleSeqUtils
export SeqShow, SeqPrint, SeqName, SeqNum, Println, ZArray

"""

Nemo is a library designed, developed and maintained by William Hart with the help of others. Many functions in our project are based on Nemo.
"""
const ModuleSeqUtils = ""

"""

The basic integer type of our project is 'fmpz'. Think of a 'fmpz' as a 'BigInt'. 'fmpz' stands for 'fast multiple precision zahl' (zahl=integer). The mathematical symbol for the ring of integers is the blackboard (double-struck) Z, also written fmpz. In reference to this fmpz(n) defines the integer n as a fmpz.
"""
function Nemofmpz(n::Int)
    return fmpz(n)
end

"""

Return an array of type `fmpz` of length ``n`` preset with ``0``.
"""
ZArray(len::Int) = fill(fmpz(0), len)

"""

Return an array of type fmpz and of size ``n`` filled by the values of the function ``f`` in the range `offset:n`.
"""
function ZArray(n::Int, f::Function, offset = 1)
    n ≤ 0 && return fmpz[]
    fmpz[f(k) for k ∈ offset:offset+n-1] # !
end

"""

Return the size of the array ``A``.
"""
SeqSize(A) = Base.length(A)

"""

Return the range of a SeqArray.
"""
SeqRange(A) = 1:SeqSize(A)

"""

Print the array ``A`` in the format ``n ↦ A[n]`` for n in the given range.
"""
function SeqShow(A::Array, R::AbstractRange, offset = 0)
    for k ∈ R
        if isassigned(A, k)
            println(k + offset, " ↦ ", A[k])
        else
            println(k + offset, " ↦ ", "undef")
        end
    end
end

"""

Print the array ``A`` in the format ``n ↦ A[n]`` for n in the range first:last.
"""
SeqShow(A, first::Int, last::Int) = SeqShow(A, first:last)

"""

Print the array ``A`` in the format ``n ↦ A[n]``.
"""
SeqShow(A::Array, offset = 0) = SeqShow(A, 1:length(A), offset - 1)

function print_without_type(io, v::AbstractVector)
    print(io, "[")
    for (i, el) ∈ enumerate(v)
        i > 1 && print(io, ", ")
        print(io, el)
    end
    println(io, "]")
end

"""

Print the array without typeinfo.
"""
Println(v::AbstractVector) = print_without_type(IOContext(stdout), v)

"""

Print the array without typeinfo.
"""
function Println(V::Array{Array{fmpz,1},1})
    for v ∈ V  v |> Println end
end

"""

Print the array with or without typeinfo.
"""
function SeqPrint(v::AbstractVector, typeinfo = false)
    typeinfo ? println(v) : Println(v)
end

"""

Valid prefixes to the numerical part of the OEIS A-numbers.

* C  => Coroutine (channel)
* F  => Filter (all below n)
* G  => Generating function
* I  => Iteration
* L  => List (array based)
* M  => Matrix (2-dim square)
* P  => Polynomial (over fmpz or QQ)
* R  => RealFunction (Float64)
* S  => Staircase (list iteration)
* T  => Triangle (list iteration)
* V  => Value (single term)
* is => is a (predicate), boolean
"""
const ValidPrefixes = ['C', 'F', 'G', 'I', 'L', 'M', 'R', 'S', 'T', 'V']

"""

Return the name of a OEIS sequence given a similar named function as a string.
"""
function SeqName(fun)
    aname = string(fun)
    for X ∈ ValidPrefixes
        aname = replace(aname, X => 'A')
    end

    if !occursin(r"^A[0-9]{6}$", aname)
        fullname = split(aname, ".")
        aname = String(fullname[end])
        if !occursin(r"^A[0-9]{6}$", aname)
            @error(aname, " is not a valid A-name!")
            return
        end
    end
    aname
end

"""

Return the A-number of a OEIS sequence given a similar named function as an integer.
"""
function SeqNum(seq)
    aname = SeqName(seq)
    parse(Int, aname[2:end])
end

#START-TEST-########################################################

using Test

function test() end

function demo()
    A = fmpz[1, 2, 3, 4, 5]

    println()
    SeqShow(A)
    println()

    SeqShow(A, -4)
    println()

    A |> Println

    println()
    V = Array{fmpz,1}[[1], [0, 1], [0, 1, 70], [0, 1, 990, 34650]]
    for v ∈ V  v |> Println end
end

function perf() end

function main()
    test()
    demo()
    perf()
end

main()

end # module

# a(f) = applicable(f, 1) ? f(1) : f(1,2)
# m(f) = first(methods(h)).nargs == 1 ? f(1) : f(1,2)

max_stars_repo_name: OpenLibMathSeq/Sequences

max_issues_repo_name: OpenLibMathSeq/Sequences


text:
# common utilities

## convenient type alias
#
#  These types signficantly reduces the need of using
#  type parameters in functions (which are often just
#  for the purpose of restricting the arrays to real)
#

typealias RealArray{T<:Real,N} AbstractArray{T,N}
typealias RealVector{T<:Real} AbstractArray{T,1}
typealias RealMatrix{T<:Real} AbstractArray{T,2}

typealias IntegerArray{T<:Integer,N} AbstractArray{T,N}
typealias IntegerVector{T<:Integer} AbstractArray{T,1}
typealias IntegerMatrix{T<:Integer} AbstractArray{T,2}

typealias RealFP Union(Float32, Float64)

## conversion from real to fp types

fptype(::Type{Float32}) = Float32
fptype(::Type{Float64}) = Float64
fptype(::Type{Bool}) = Float32
fptype(::Type{Int8}) = Float32
fptype(::Type{Uint8}) = Float32
fptype(::Type{Int16}) = Float32
fptype(::Type{Uint16}) = Float32
fptype(::Type{Int32}) = Float64
fptype(::Type{Uint32}) = Float64
fptype(::Type{Int64}) = Float64
fptype(::Type{Uint64}) = Float64
fptype(::Type{Int128}) = Float64
fptype(::Type{Uint128}) = Float64

## consistent way to raise error
#
# only change this, if we decide to change to way that
# dimension errors are raised
#

raise_dimerror(msg::String) = throw(DimensionMismatch(msg))
raise_dimerror() = raise_dimerror("")

# weight vector:
#
#    a wrapper that indicator a vector representing a sequence of weights
#
immutable WeightVec{W<:Real,Vec<:RealVector}
	values::Vec
	sum::W
end

WeightVec{Vec<:RealVector,W<:Real}(vs::Vec,wsum::W) = WeightVec{W,Vec}(vs, wsum)
WeightVec(vs::RealVector) = WeightVec(vs, sum(vs))

weights(vs::RealVector) = WeightVec(vs)
weights(vs::RealArray) = WeightVec(vec(vs))

eltype(wv::WeightVec) = eltype(wv.values)
length(wv::WeightVec) = length(wv.values)
values(wv::WeightVec) = wv.values
sum(wv::WeightVec) = wv.sum
isempty(wv::WeightVec) = isempty(wv.values)


text:
type Formula
    lhs::Union{Symbol, Expr, Void}
    rhs::Union{Symbol, Expr, Integer}
end

macro ~(lhs, rhs)
    ex = Expr(:call,
              :Formula,
              Base.Meta.quot(lhs),
              Base.Meta.quot(rhs))
    return ex
end

function Base.show(io::IO, f::Formula)
    print(io,
          string("Formula: ",
                 f.lhs == nothing ? "" : f.lhs, " ~ ", f.rhs))
end


text:
function _pointlayer(df::DataFrame, x::Symbol, y::Symbol; groupBy::Symbol=Symbol())
    groupBy==Symbol() ? layer(df,x=x,y=y,Geom.point) : layer(df,x=x,y=y,color=groupBy,Geom.point)
end
function _dataframelinearregression(df::DataFrame,x::Symbol,y::Symbol,group=[])
    dataX = convert(Array, df[x])
    dataY = convert(Array, df[y])
    xVals = collect(extrema(dataX))

    A = xVals[2]-xVals[1] > 1e-6 ?
        hcat(ones(dataX),dataX)\dataY : # linear regression: y = A[1] + A[2]x
        (0.0,0.0)

    isempty(group) ? DataFrame(x=xVals, y=A[1]+A[2]*xVals) :
                     DataFrame(x=xVals, y=A[1]+A[2]*xVals, group=group)
end
function _linearregressionlayer(df::DataFrame, x::Symbol, y::Symbol; groupBy::Symbol=Symbol())
    if groupBy==Symbol()
        dfLinear = _dataframelinearregression(df,x,y)
        layer(dfLinear,x=x,y=y,Geom.line)
    else
        dfLinear = DataFrame()
        for group in unique(df[groupBy])
            dfMasked = df[df[groupBy].==group, :]
            dfLinear = vcat(dfLinear, _dataframelinearregression(dfMasked,x,y,group))
        end
        layer(dfLinear,x=:x,y=:y,color=:group,Geom.line)
    end
end

function linearregressionplot(df::DataFrame, x::Symbol, y::Symbol, args...; groupBy::Symbol=Symbol(), kwargs...)
    plot( _pointlayer(df,x,y,groupBy=groupBy),
          _linearregressionlayer(df,x,y,groupBy=groupBy),
          args...;
          kwargs...)
end


text:
using Constructs
using Intervals
using Test

@testset "Constructs.jl" begin
    @testset "primitive io" begin
        @testset "primitive type $type" for type in (Bool, UInt8, UInt16, UInt32, UInt64, UInt128, Int8, Int16, Int32, Int64, Int128, Float16, Float32, Float64)
            @test estimatesize(type) == sizeof(type)
            @test deserialize(type, zeros(UInt8, sizeof(type))) == zero(type)
            @test serialize(zero(type)) == zeros(UInt8, sizeof(type))
        end
        @testset "char $c" for c in ['\x00', 'A', 'α', '啊', '🆗']
            codes = transcode(UInt8, [codepoint(c)])
            @test length(codes) in estimatesize(Char)
            @test deserialize(Char, codes) == c
            @test serialize(c) == codes
        end
    end
    @testset "singleton" begin
        @testset "auto singleton for $instance" for instance in (nothing, missing)
            type = typeof(instance)
            @test estimatesize(type) == 0
            @test deserialize(type, UInt8[]) === instance
            @test serialize(instance) == UInt8[]
        end
        @testset "non-singleton type $type" for type in (Bool, Union{}, Char, DataType)
            @test_throws ArgumentError Singleton(type)
        end
    end
    @testset "byte order" begin
        be = (
            (0x0102, [0x01, 0x02]),
            (0x01020304, [0x01, 0x02, 0x03, 0x04]),
            (0x0102030405060708, [0x01, 0x02, 0x03, 0x04, 0x05, 0x06, 0x07, 0x08]),
            (Int32(-2140118960), [0x80, 0x70, 0x60, 0x50]),
            (Inf16, [0x7c, 0x00]),
            (Float32(-1.1), [0xbf, 0x8c, 0xcc, 0xcd]),
            )
        @testset "little endian for $n" for (n, bs) in be
            type = typeof(n)
            lbs = reverse(bs)
            @test estimatesize(LittleEndian(type)) == sizeof(type)
            @test deserialize(LittleEndian(type), lbs) == n
            @test serialize(LittleEndian(type), n) == lbs
        end
        @testset "big endian for $n" for (n, bs) in be
            type = typeof(n)
            @test estimatesize(BigEndian(type)) == sizeof(type)
            @test deserialize(BigEndian(type), bs) == n
            @test serialize(BigEndian(type), n) == bs
        end
    end
    @testset "enum" begin
        @enum Fruit::Int8 begin
            apple
            banna
            orange
        end
        @testset "auto type" begin
            @test estimatesize(Fruit) == sizeof(Int8)
            @test deserialize(Fruit, [0x01]) == banna
            @test serialize(orange) == [0x02]
        end
        @testset "override base type" begin
            @test estimatesize(IntEnum(UInt8, Fruit)) == sizeof(UInt8)
            @test deserialize(IntEnum(UInt8, Fruit), [0x01]) == banna
            @test serialize(IntEnum(UInt8, Fruit), orange) == [0x02]
        end
        @testset "override base type construct" begin
            @test estimatesize(IntEnum(BigEndian(UInt16), Fruit)) == sizeof(UInt16)
            @test deserialize(IntEnum(BigEndian(UInt16), Fruit), [0x00, 0x01]) == banna
            @test serialize(IntEnum(BigEndian(UInt16), Fruit), orange) == [0x00, 0x02]
        end
    end
    @testset "const" begin
        @test estimatesize(Const(0x0102)) == sizeof(0x0102)
        @test estimatesize(Const(b"BMP")) == sizeof(b"BMP")
        @test estimatesize(Const(Int32, 0x0102)) == sizeof(Int32)
        @test estimatesize(Const(BigEndian(Int32), 0x0102)) == sizeof(Int32)
        @test deserialize(Const(BigEndian(UInt16), 0x0102), [0x01, 0x02]) == 0x0102
        @test_throws ValidationError deserialize(Const(LittleEndian(UInt16), 0x0102), [0x01, 0x02])
        @test serialize(Const(BigEndian(UInt16), 0x0102), 0x0102) == [0x01, 0x02]
        @test_throws ValidationError serialize(Const(0x0102), 0x0201)
    end
    @testset "collections" begin
        @testset "Array" begin
            @test estimatesize(@Array(Int64)) == sizeof(Int64)
            @test estimatesize(@Array(Int64, 10)) == 10*sizeof(Int64)
            @test estimatesize(@Array(Int64, 2, 3, 5)) == 2*3*5*sizeof(Int64)
            @test estimatesize(@Array(Int64, this.size)) == Interval(UInt(0), nothing)
            @test deserialize(@Array(Int8, 3), [0x01, 0xff, 0x00]) == Int8[1, -1, 0]
            @test serialize(@Array(Int8, 3), Int8[1, -1, 0]) == [0x01, 0xff, 0x00]
            @test deserialize(@Array(Int8, 2, 3), Vector{UInt8}(1:6)) == Int8[1 3 5; 2 4 6]
            @test serialize(@Array(Int8, 2, 3), Int8[1 2 3; 4 5 6]) == [0x01, 0x04, 0x02, 0x05, 0x03, 0x06]
        end
        @testset "GreedyVector" begin
            @test estimatesize(GreedyVector(Int8)) == Interval(UInt(0), nothing)
            @test deserialize(GreedyVector(Int8), [0x01, 0xff, 0x00]) == Int8[1, -1, 0]
            @test serialize(GreedyVector(Int8), Int8[1, -1, 0]) == [0x01, 0xff, 0x00]
            @test deserialize(GreedyVector(BigEndian(UInt16)), [0x01, 0xff, 0x02, 0xab, 0xcc]) == [0x01ff, 0x02ab]
            @test serialize(GreedyVector(BigEndian(UInt16)), [0x01ff, 0xcc0a]) == [0x01, 0xff, 0xcc, 0x0a]
        end
    end
    @testset "internal" begin
        @testset "construct type" begin
            @test Constructs.constructtype(Int32) == Int32
            @test Constructs.constructtype(JuliaSerializer()) == Any
            @test Constructs.constructtype(Padding(4)) == Nothing
            @test Constructs.constructtype(BigEndian(UInt)) == UInt
            @test Constructs.constructtype(Const(0x0102)) == UInt16
            @test Constructs.constructtype(Const(b"BMP")) == Vector{UInt8}
            @test Constructs.constructtype(@Array(Int)) == Array{Int, 0}
            @test Constructs.constructtype(@Array(Float64, 10)) == Array{Float64, 1}
            @test Constructs.constructtype(@Array(BigEndian(UInt16), 5, 17)) == Array{UInt16, 2}
            @test Constructs.constructtype(@Array(Char, this.size)) == Array{Char, 1}
            @test Constructs.constructtype(@Array(LittleEndian(UInt128), this.width, this.height)) == Array{UInt128, 2}
        end
    end
    @testset "macro expand" begin # just check if there is any error raised during macro expanding
        structonly = quote
            @construct struct Bitmap
                signature::Const(b"BMP")
                width::UInt32
                height::UInt32
                ::Padding(8)
                pixel::@Array(UInt8, (this.width, this.height))
                rest
            end
        end
        withname = quote
            @construct BitmapConstruct struct Bitmap
                signature::Const(b"BMP")
                width::UInt32
                height::UInt32
                ::Padding(8)
                pixel::@Array(UInt8, (this.width, this.height))
                rest
            end
        end
        @testset "@construct" for ex in [structonly, withname]
            @inferred Expr macroexpand(@__MODULE__, ex)
        end
    end
end


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

mutable struct StorageAccountProperties <: SwaggerModel
    provisioningState::Any # spec type: Union{ Nothing, String } # spec name: provisioningState
    primaryEndpoints::Any # spec type: Union{ Nothing, Endpoints } # spec name: primaryEndpoints
    primaryLocation::Any # spec type: Union{ Nothing, String } # spec name: primaryLocation
    statusOfPrimary::Any # spec type: Union{ Nothing, String } # spec name: statusOfPrimary
    lastGeoFailoverTime::Any # spec type: Union{ Nothing, DateTime } # spec name: lastGeoFailoverTime
    secondaryLocation::Any # spec type: Union{ Nothing, String } # spec name: secondaryLocation
    statusOfSecondary::Any # spec type: Union{ Nothing, String } # spec name: statusOfSecondary
    creationTime::Any # spec type: Union{ Nothing, DateTime } # spec name: creationTime
    customDomain::Any # spec type: Union{ Nothing, CustomDomain } # spec name: customDomain
    secondaryEndpoints::Any # spec type: Union{ Nothing, Endpoints } # spec name: secondaryEndpoints
    encryption::Any # spec type: Union{ Nothing, Encryption } # spec name: encryption
    accessTier::Any # spec type: Union{ Nothing, String } # spec name: accessTier
    supportsHttpsTrafficOnly::Any # spec type: Union{ Nothing, Bool } # spec name: supportsHttpsTrafficOnly
    networkAcls::Any # spec type: Union{ Nothing, NetworkRuleSet } # spec name: networkAcls

    function StorageAccountProperties(;provisioningState=nothing, primaryEndpoints=nothing, primaryLocation=nothing, statusOfPrimary=nothing, lastGeoFailoverTime=nothing, secondaryLocation=nothing, statusOfSecondary=nothing, creationTime=nothing, customDomain=nothing, secondaryEndpoints=nothing, encryption=nothing, accessTier=nothing, supportsHttpsTrafficOnly=false, networkAcls=nothing)
        o = new()
        validate_property(StorageAccountProperties, Symbol("provisioningState"), provisioningState)
        setfield!(o, Symbol("provisioningState"), provisioningState)
        validate_property(StorageAccountProperties, Symbol("primaryEndpoints"), primaryEndpoints)
        setfield!(o, Symbol("primaryEndpoints"), primaryEndpoints)
        validate_property(StorageAccountProperties, Symbol("primaryLocation"), primaryLocation)
        setfield!(o, Symbol("primaryLocation"), primaryLocation)
        validate_property(StorageAccountProperties, Symbol("statusOfPrimary"), statusOfPrimary)
        setfield!(o, Symbol("statusOfPrimary"), statusOfPrimary)
        validate_property(StorageAccountProperties, Symbol("lastGeoFailoverTime"), lastGeoFailoverTime)
        setfield!(o, Symbol("lastGeoFailoverTime"), lastGeoFailoverTime)
        validate_property(StorageAccountProperties, Symbol("secondaryLocation"), secondaryLocation)
        setfield!(o, Symbol("secondaryLocation"), secondaryLocation)
        validate_property(StorageAccountProperties, Symbol("statusOfSecondary"), statusOfSecondary)
        setfield!(o, Symbol("statusOfSecondary"), statusOfSecondary)
        validate_property(StorageAccountProperties, Symbol("creationTime"), creationTime)
        setfield!(o, Symbol("creationTime"), creationTime)
        validate_property(StorageAccountProperties, Symbol("customDomain"), customDomain)
        setfield!(o, Symbol("customDomain"), customDomain)
        validate_property(StorageAccountProperties, Symbol("secondaryEndpoints"), secondaryEndpoints)
        setfield!(o, Symbol("secondaryEndpoints"), secondaryEndpoints)
        validate_property(StorageAccountProperties, Symbol("encryption"), encryption)
        setfield!(o, Symbol("encryption"), encryption)
        validate_property(StorageAccountProperties, Symbol("accessTier"), accessTier)
        setfield!(o, Symbol("accessTier"), accessTier)
        validate_property(StorageAccountProperties, Symbol("supportsHttpsTrafficOnly"), supportsHttpsTrafficOnly)
        setfield!(o, Symbol("supportsHttpsTrafficOnly"), supportsHttpsTrafficOnly)
        validate_property(StorageAccountProperties, Symbol("networkAcls"), networkAcls)
        setfield!(o, Symbol("networkAcls"), networkAcls)
        o
    end
end # type StorageAccountProperties

const _property_map_StorageAccountProperties = Dict{Symbol,Symbol}(Symbol("provisioningState")=>Symbol("provisioningState"), Symbol("primaryEndpoints")=>Symbol("primaryEndpoints"), Symbol("primaryLocation")=>Symbol("primaryLocation"), Symbol("statusOfPrimary")=>Symbol("statusOfPrimary"), Symbol("lastGeoFailoverTime")=>Symbol("lastGeoFailoverTime"), Symbol("secondaryLocation")=>Symbol("secondaryLocation"), Symbol("statusOfSecondary")=>Symbol("statusOfSecondary"), Symbol("creationTime")=>Symbol("creationTime"), Symbol("customDomain")=>Symbol("customDomain"), Symbol("secondaryEndpoints")=>Symbol("secondaryEndpoints"), Symbol("encryption")=>Symbol("encryption"), Symbol("accessTier")=>Symbol("accessTier"), Symbol("supportsHttpsTrafficOnly")=>Symbol("supportsHttpsTrafficOnly"), Symbol("networkAcls")=>Symbol("networkAcls"))
const _property_types_StorageAccountProperties = Dict{Symbol,String}(Symbol("provisioningState")=>"String", Symbol("primaryEndpoints")=>"Endpoints", Symbol("primaryLocation")=>"String", Symbol("statusOfPrimary")=>"String", Symbol("lastGeoFailoverTime")=>"DateTime", Symbol("secondaryLocation")=>"String", Symbol("statusOfSecondary")=>"String", Symbol("creationTime")=>"DateTime", Symbol("customDomain")=>"CustomDomain", Symbol("secondaryEndpoints")=>"Endpoints", Symbol("encryption")=>"Encryption", Symbol("accessTier")=>"String", Symbol("supportsHttpsTrafficOnly")=>"Bool", Symbol("networkAcls")=>"NetworkRuleSet")
Base.propertynames(::Type{ StorageAccountProperties }) = collect(keys(_property_map_StorageAccountProperties))
Swagger.property_type(::Type{ StorageAccountProperties }, name::Symbol) = Union{Nothing,eval(Meta.parse(_property_types_StorageAccountProperties[name]))}
Swagger.field_name(::Type{ StorageAccountProperties }, property_name::Symbol) =  _property_map_StorageAccountProperties[property_name]

const _allowed_StorageAccountProperties_provisioningState = ["Creating", "ResolvingDNS", "Succeeded"]

const _allowed_StorageAccountProperties_statusOfPrimary = ["available", "unavailable"]

const _allowed_StorageAccountProperties_statusOfSecondary = ["available", "unavailable"]

const _allowed_StorageAccountProperties_accessTier = ["Hot", "Cool"]

function check_required(o::StorageAccountProperties)
    true
end

function validate_property(::Type{ StorageAccountProperties }, name::Symbol, val)
    if name === Symbol("provisioningState")
        Swagger.validate_param(name, "StorageAccountProperties", :enum, val, _allowed_StorageAccountProperties_provisioningState)
    end
    if name === Symbol("statusOfPrimary")
        Swagger.validate_param(name, "StorageAccountProperties", :enum, val, _allowed_StorageAccountProperties_statusOfPrimary)
    end
    if name === Symbol("statusOfSecondary")
        Swagger.validate_param(name, "StorageAccountProperties", :enum, val, _allowed_StorageAccountProperties_statusOfSecondary)
    end
    if name === Symbol("accessTier")
        Swagger.validate_param(name, "StorageAccountProperties", :enum, val, _allowed_StorageAccountProperties_accessTier)
    end
end


text:
"""
    KeypointTensor{N, T, M} <: Block

Block for encoded [`Keypoints`](#)`{N, T, M}` returned by
[`KeypointPreprocessing`](#).
"""
struct KeypointTensor{N, T, M} <: Block
    sz::NTuple{M, Int}
end

mockblock(block::KeypointTensor{N}) where N = rand(SVector{N, Float32}, block.sz)
function checkblock(block::KeypointTensor{N, T}, obs::AbstractArray{T}) where {N, T}
    return length(obs) == (prod(block.sz) * N)
end

"""
    KeypointPreprocessing(bounds) <: Encoding

Scale a `Keypoints` block falling in a rectangle of `bounds` so
that they lie between -1 and 1.
"""
struct KeypointPreprocessing{N, T<:Number} <: Encoding
    bounds::NTuple{N, Int}
end
KeypointPreprocessing(bounds::NTuple{N, Int}) where N = KeypointPreprocessing{N, Float32}(bounds)

function encode(enc::KeypointPreprocessing{N, T}, context, block::Keypoints{N}, obs) where {N, T}
    ks = map(k -> (SVector{N, T}(k) .* (convert(T, 2) ./ enc.bounds)) .- one(T), obs)
    return reinterpret(T, ks)
end

function decode(enc::KeypointPreprocessing{N, T}, context, block::KeypointTensor{N}, obs) where {N, T}
    ks = reshape(reinterpret(SVector{N, T}, obs), block.sz)
    return map(k -> ((k) .+ one(T)) ./ (convert(T, 2) ./ SVector{N, T}(enc.bounds)), ks)
end

encodedblock(::KeypointPreprocessing{N, T}, block::Keypoints{N, M}) where {N, T, M} = KeypointTensor{N, T, M}(block.sz)
decodedblock(::KeypointPreprocessing{N}, block::KeypointTensor{N}) where N = Keypoints{N}(block.sz)

# ## Optional interfaces

# The default loss function to compare encoded keypoints is Mean Squared Error:

function blocklossfn(outblock::KeypointTensor{N}, yblock::KeypointTensor{N}) where {N}
    outblock.sz == yblock.sz || error("Sizes of $outblock and $yblock differ!")
    return Flux.Losses.mse
end

# ## Tests

@testset "KeypointPreprocessing [encoding]" begin

    ks = [
        SVector{2, Float32}(10, 10),
        SVector{2, Float32}(50, 80),
    ]
    sz = (100, 100)

    block = Keypoints{2}(2)
    enc = KeypointPreprocessing(sz)
    ctx = Training()

    testencoding(enc, block, ks)
    y = encode(enc, ctx, block, ks)
    ks_ = decode(enc, ctx, encodedblock(enc, block), y)
    @test ks ≈ ks_

end


text:
module CategoricalAlgebra

include("Permutations.jl")

end


text:
using ActorInterfacesTests
using Test


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

# Collection of sources required to build MPCBuilder
name = "MPC"
version = v"1.1.0"
sources = [
    "https://ftp.gnu.org/gnu/mpc/mpc-$version.tar.gz" =>
    "6985c538143c1208dcb1ac42cedad6ff52e267b47e5f970183a3e75125b43c2e",
]

# Bash recipe for building across all platforms

script = """cd mpc-$version
"""*raw"""
./configure --prefix=$prefix --host=$target --enable-shared --disable-static --with-gmp=$prefix --with-mpfr=$prefix
make -j
make install
# On Windows, make sure non-versioned filename exists...
if [[ ${target} == *mingw* ]]; then
    cp -v ${prefix}/bin/libmpc-*.dll ${prefix}/bin/libmpc.dll
fi
"""
# These are the platforms we will build for by default, unless further
# platforms are passed in on the command line
platforms = supported_platforms()

# The products that we will ensure are always built
products(prefix) = [
    LibraryProduct(prefix, "libmpc", :libmpc)
]

# Dependencies that must be installed before this package can be built
dependencies = [
    "https://github.com/JuliaPackaging/Yggdrasil/releases/download/GMP-v6.1.2-1/build_GMP.v6.1.2.jl"
    "https://github.com/JuliaPackaging/Yggdrasil/releases/download/MPFR-v4.0.2-0/build_MPFR.v4.0.2.jl"
]

# Build the tarballs, and possibly a `build.jl` as well.
build_tarballs(ARGS, name, version, sources, script, platforms, products, dependencies)

max_stars_repo_name: UnofficialJuliaMirror/isuruf-MPCBuilder

max_issues_repo_name: UnofficialJuliaMirror/isuruf-MPCBuilder

max_forks_repo_name: UnofficialJuliaMirror/isuruf-MPCBuilder


text:
using Base.Threads
import DataFrames.AbstractDataFrame
import DataFrames.DataFrame
import IndexedTables.NDSparse

"""
  column(df, :colname)

Extract a column from an AbstractDataFrame
"""
function column(dt::AbstractDataFrame, col::Symbol)
  i = dt.colindex.lookup[col]
  dt.columns[i]
end

"""
  select(:col)

Return a funciton that obtains a column with the named symbol from an AbstractDataFrame or NDSparse
"""
function select(col::Symbol)
  return function(df::Union{AbstractDataFrame,NDSparse})
    column(df, col)
  end
end

# from https://discourse.julialang.org/t/whats-the-fastest-way-to-generate-1-2-n/7564/15?u=xiaodai
using Base.Threads
function fcollect(N::Integer, T = Int)
    nt = nthreads()
    n,r = divrem(N,nt)
    a = Vector{T}(undef, N)
    @threads for i=1:nt
        ioff = (i-1)*n
        nn = ifelse(i == nt, n+r, n)
        @inbounds for j=1:nn
            a[ioff+j] = ioff+j
        end
    end
    a
end

# a general isgrouped algorithm
function isgrouped(grps::AbstractVector)
  # find where the change happens
  a = BitArray(undef, 2^(sizeof(UInt32)*8))
  a .= false

  hindex = hash(grps[1])

  a[Base.trunc_int(UInt32,hindex) + 1] = true
  a[Base.trunc_int(UInt32,hindex >> 32) + 1] = true

  for i = 2:length(grps)
      if grps[i-1] != grps[i]
          hindex = hash(grps[i])
          hindex1 = Base.trunc_int(UInt32, hindex) + 1
          hindex2 = Base.trunc_int(UInt32, hindex >> 32) + 1
          if a[hindex1] && a[hindex2]
              return false
          else
              a[hindex1] = true
              a[hindex2] = true
          end
      end
  end
  return true
end

"Generate CategoricalArrays"
function genca(refs::Vector{U}, pool::Vector{T}) where {U<:Unsigned, T}
    CategoricalArray{T, 1}(rand(U(1):U(length(pool)), length(refs)), CategoricalPool(pools, false));
end

"Simple data structure for carrying a string vector and its index; this allows
`sorttwo!` to sort the radix of the string vector and reorder the string and its
index at the same time opening the way for faster sort_perm"
# struct StringIndexVector
#     svec::Vector{String}
#     index::Vector{Int}
# end

# function setindex!(siv::StringIndexVector, X::StringIndexVector, inds)
#     siv.svec[inds] = X.svec
#     siv.index[inds] = X.index
# end

# function setindex!(siv::StringIndexVector, X, inds)
#     siv.svec[inds] = X[1]
#     siv.index[inds] = X[2]
# end

# getindex(siv::StringIndexVector, inds::Integer) = siv.svec[inds], siv.index[inds]
# getindex(siv::StringIndexVector, inds...) = StringIndexVector(siv.svec[inds...], siv.index[inds...])
# similar(siv::StringIndexVector) = StringIndexVector(similar(siv.svec), similar(siv.index))


text:
# text
# the 'toy' API... "Any serious application should avoid them." :)

"""
    text(str)
    text(str, pos)
    text(str, x, y)
    text(str, pos, halign=:left)
    text(str, valign=:baseline)
    text(str, valign=:baseline, halign=:left)
    text(str, pos, valign=:baseline, halign=:left)

Draw the text in the string `str` at `x`/`y` or `pt`, placing the start of
the string at the point. If you omit the point, it's placed at the current `0/0`. In Luxor,
placing text doesn't affect the current point.

`:halign` can be `:left`, `:center`, or `:right`. `:valign` can be `:baseline`, `:top`,
`:middle`, or `:bottom`.

However, the `:valign` doesn't work properly because we're using
Cairo's so-called "toy" interface... :(

"""

function text(t, x=0, y=0; halign=:left, valign=:baseline)
  wdfactor = 0
  htfactor = 0
  textents = textextents(t)
  halign == :left         && (wdfactor = 0  )
  halign == :center       && (wdfactor = 0.5)
  halign == :right        && (wdfactor = 1  )
  valign == :baseline     && (htfactor = 0  )
  valign == :top          && (htfactor = -1)
  valign == :middle       && (htfactor = 0.5)
  valign == :bottom       && (htfactor = 1  )
  textpointx = x - (wdfactor * (textents[3]))
  textpointy = y - (htfactor * textents[4])
  gsave()
  Cairo.move_to(currentdrawing.cr, textpointx, textpointy)
  Cairo.show_text(currentdrawing.cr, t)
  grestore()
end

text(t, pt::Point; kwargs...) = text(t, pt.x, pt.y; kwargs...)

"""
    textcentred(str)
    textcentred(str, x, y)
    textcentred(str, pt)

Draw text in the string `str` centered at `x`/`y` or `pt`.
If you omit the point, it's placed at 0/0.

Text doesn't affect the current point!
"""

function textcentred(t, x=0, y=0)
  text(t, x, y, halign=:center)
#  textwidth = textextents(t)[5]
#  text(t, x - textwidth/2., y)
end

textcentred(t, pt::Point) = textcentred(t, pt.x, pt.y)

"""
Do I spell textcentred wrong...?
"""
textcentered = textcentred

"""
    textright(str)
    textright(str, x, y)
    textright(str, pt)

Draw text in the string `str` right-aligned at `x`/`y` or `pt`.
If you omit the point, it's placed at 0/0.

Text doesn't affect the current point!
"""

function textright(t, x=0, y=0)
  text(t, x, y, halign=:right)
  # textwidth = textextents(t)[5]
  # text(t, x - textwidth, y)
end

textright(t, pt::Point) = textright(t, pt.x, pt.y)

"""

    fontface(fontname)

Select a font to use. If the font is unavailable, it defaults to Helvetica/San Francisco
(on macOS).
"""
fontface(f) = Cairo.select_font_face(currentdrawing.cr, f, Cairo.FONT_SLANT_NORMAL, Cairo.FONT_WEIGHT_NORMAL)

"""

    fontsize(n)

Set the font size to `n` points. Default is 10pt.

"""
fontsize(n) = Cairo.set_font_size(currentdrawing.cr, n)

"""
    textextents(str)

Return the measurements of the string `str` when set using the current font settings:

- x_bearing
- y_bearing
- width
- height
- x_advance
- y_advance

The bearing is the displacement from the reference point to the upper-left corner of the
bounding box. It is often zero or a small positive value for x displacement, but can be
negative x for characters like j; it's almost always a negative value for y displacement.

The width and height then describe the size of the bounding box. The advance takes you to
the suggested reference point for the next letter. Note that bounding boxes for subsequent
blocks of text can overlap if the bearing is negative, or the advance is smaller than the
width would suggest.

Example:

    textextents("R")

returns

    [1.18652; -9.68335; 8.04199; 9.68335; 9.74927; 0.0]
"""

textextents(str) = Cairo.text_extents(currentdrawing.cr, str)

"""
    textpath(t)

Convert the text in string `t` to a new path, for subsequent filling/stroking etc...
"""
function textpath(t)
    Cairo.text_path(currentdrawing.cr, t)
end

"""
Place a string of text on a curve. It can spiral in or out.

```
textcurve(the_text,
          start_angle,
          start_radius,
          x_pos = 0,
          y_pos = 0;
          # optional keyword arguments:
          spiral_ring_step = 0,    # step out or in by this amount
          letter_spacing = 0,      # tracking/space between chars, tighter is (-), looser is (+)
          spiral_in_out_shift = 0, # + values go outwards, - values spiral inwards
          clockwise = true         #
          )
```

`start_angle` is relative to +ve x-axis, arc/circle is centred on `(x_pos,y_pos)` with
radius `start_radius`.
"""

function textcurve(the_text, start_angle, start_radius, x_pos=0, y_pos=0;
  # keyword optional arguments
  spiral_ring_step = 0,
  letter_spacing = 0, # tracking/space between chars, tighter is (-), looser is (+)
  spiral_in_out_shift = 0, # makes spiral outwards (+) or inwards (-)
  clockwise = true
  )
  refangle = start_angle
  current_radius = start_radius
  spiral_space_step = 0
  xx = 0
  yy = 0
  angle_step = 0
  radius_step = 0
  for i in the_text
    glyph = string(i)
    glyph_x_bearing, glyph_y_bearing, glyph_width,
      glyph_height, glyph_x_advance, glyph_y_advance = textextents(glyph)
    spiral_space_step = glyph_x_advance + letter_spacing
    cnter = (2pi * current_radius) / spiral_space_step
    radius_step = (spiral_ring_step + spiral_in_out_shift) / cnter
    current_radius += radius_step
    angle_step += (glyph_x_advance / 2.) + letter_spacing/2.
    if clockwise
      refangle += angle_step / current_radius
    else
      refangle -= angle_step / current_radius
    end
    angle_step = (glyph_x_advance / 2.) + letter_spacing/2.
    xx = cos(refangle) * current_radius + x_pos
    yy = sin(refangle) * current_radius + y_pos
    gsave()
    translate(xx, yy)
    if clockwise
      rotate(pi/2. + refangle)
    else
      rotate(-pi/2. + refangle)
    end
    textcentred(glyph, 0., 0.)
    grestore()
    current_radius < 10. && break
  end
end

textcurve(the_text, start_angle, start_radius, centre::Point; kwargs...) =
  textcurve(the_text, start_angle, start_radius, centre.x, centre.y; kwargs...)


text:
#=
    This code is used to generate the integer power lookup tables.

    result of ipow(ipow_base::T, ipow_power::T) where {T}
      - either an overflow exception is generated
      - orelse the result is a positive integer of type `T`
          - result = ipow_base^ipow_power
          - result == BigInt(ipow_base)^BigInt(ipow_power)

    Each result is determined through exaustive search for overflow
    using Base.Int128s and DecFP.Dec128s.

    The results are independently validated by comparison with
      the results obtained using a different algorithmic approach
      and computed using BigInts exclusively.
=#
#=
    (maxint^10 that does not overflow, maxint == 78)
    (maxint^10 that does not overflow, maxint == 78)

=#

function pospow_safemax(base::T) where {T}
    pows = T(1):T(4096)
    for p in pows
        pospow_is_unsafe(base, p) && return p-1
    end
end

function pospow_unsafemin(base::T) where {T<:Signed}
    pows = T(4096):T(-1):T(1)
    for p in pows
        pospow_is_safe(base, p) && return p+1
    end
    throw(ErrorException("Should not be reached"))
end

function pospow_is_safe(base::T, power::T) where {T<:Signed}
    issafe = true
    tresult = base^power
    bibase, bipower = BigInt(base), BigInt(power)
    biresult = bibase^bipower
    converted_result = zero(T)
    try
        converted_result = T(biresult)
    catch
        converted_result = biresult
        issafe = false
    end
    issafe = issafe  & iszero(cmp(biresult,tresult)) & iszero(cmp(converted_result, biresult))
    return issafe
end

function pospow_unsafemin(base::T) where {T<:Unsigned}
    pows = T(4096):T(-1):T(1)
    for p in pows
        pospow_is_safe(base, p) && return p+1
    end
    throw(ErrorException("Should not be reached"))
end

function pospow_is_safe(base::T, power::T) where {T<:Unsigned}
    issafe = true
    tresult = base^power
    bibase, bipower = BigInt(base), BigInt(power)
    biresult = bibase^bipower
    converted_result = zero(T)
    try
        converted_result = T(biresult)
    catch
        converted_result = biresult
        issafe = false
    end
    issafe = issafe  & iszero(cmp(biresult,tresult)) & iszero(cmp(converted_result, biresult))
    return issafe
end

pospow_is_unsafe(base, power) = !pospow_is_safe(base, power)

using Base.Checked
#=
julia> names(Base.Checked)
14-element Vector{Symbol}:
 :Checked
 :add_with_overflow
 :checked_abs
 :checked_add
 :checked_cld
 :checked_div
 :checked_fld
 :checked_mod
 :checked_mul
 :checked_neg
 :checked_rem
 :checked_sub
 :mul_with_overflow
 :sub_with_overflow
=#

import DecFP, Quadmath
const D128 = DecFP.Dec128
const F128 = Quadmath.Float128

#=
   maxintfloat(Float128) [1.0e34] has 113 significant bits
      bit113 is set, bits112..1 are cleared
   maxintfloat(Float128) ~= 1.0e34 D1
      10384,593717,069655,257060,992658,440192

   We choose BigFloat precision to be large enough to resolve any outlier roundings,
     setprecision(BigFloat, 5*128)
=#

setprecision(BigFloat, 5*128)
Base.BigFloat(x::F128) = BigFloat(string(x))
Base.BigFloat(x::D128) = BigFloat(string(x))
F128(x::BigFloat) = parse(F128, string(x))
D128(x::BigFloat) = parse(D128, string(x))
Base.BigInt(x::F128) = BigInt(BigFloat(x))
Base.BigInt(x::D128) = BigInt(BigFloat(x))
F128(x::BigInt) = F128(BigFloat(x))
D128(x::BigInt) = D128(BigFloat(x))

# Parametric Constraints' Typwa
const ParamT = Int128
const PartsT = DecFP.Dec128

# Parametric Constraints
const ParamT = Int128
const PartsT = DecFP.Dec128

# Parametric Intervals
basemin = ParamT(2)
basemax = ParamT(128)
basetop = ParamT(1024)
powmin  = ParamT(2)
powmax  = ParamT(128)
powtop  = ParamT(1024)

# parametric domains
const allowed_bases   = collect(basemin:basemax);
const allowed_powers  = collect(powmin:powmax);
const possible_bases  = collect(basemin:basetop);
const possible_powers = collect(powmin:powtop);

#=
julia> pospowsafemax = (pospow_safemax.(1:128)...,)
(255, 62, 39, 31, 27, 24, 22, 20, 19, 18, 18, 17, 17, 16, 16, 15, 15, 15, 14, 14, 14, 14, 13, 13, 13, 13, 13, 13, 12, 12, 12, 12, 12, 12, 12, 12, 12, 12, 11, 11, 11, 11, 11, 11, 11, 11, 11, 11, 11, 11, 11, 11, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 8)
(nothing, 30, 19, 15, 13, 11, 11, 10, 9, 9, 8, 8, 8, 8, 7, 7, 7, 7, 7, 7, 7, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4

safemaxInt8pows = Int8[0, 6, 4, 3, 3, 2, 2, 2, 2, 2, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1]

safemaxInt16pows = Int8[0, 14, 9, 7, 6, 5, 5, 4, 4, 4, 4, 4, 4, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2]

safemaxInt32pows = Int8[0, 30, 19, 15, 13, 11, 11, 10, 9, 9, 8, 8, 8, 8, 7, 7, 7, 7, 7, 7, 7, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4]

safemaxInt64pows = Int8[0, 62, 39, 31, 27, 24, 22, 20, 19, 18, 18, 17, 17, 16, 16, 15, 15, 15, 14, 14, 14, 14, 13, 13, 13, 13, 13, 13, 12, 12, 12, 12, 12, 12, 12, 12, 12, 12, 11, 11, 11, 11, 11, 11, 11, 11, 11, 11, 11, 11, 11, 11, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 8]

safemaxInt128pows = Int8[0, 126, 80, 63, 54, 49, 45, 42, 40, 38, 36, 35, 34, 33, 32, 31, 31, 30, 29, 29, 28, 28, 28, 27, 27, 27, 26, 26, 26, 25, 25, 25, 25, 24, 24, 24, 24, 24, 24, 23, 23, 23, 23, 23, 23, 22, 22, 22, 22, 22, 22, 22, 22, 22, 21, 21, 21, 21, 21, 21, 21, 21, 21, 21, 21, 21, 20, 20, 20, 20, 20, 20, 20, 20, 20, 20, 20, 20, 20, 20, 20, 19, 19, 19, 19, 19, 19, 19, 19, 19, 19, 19, 19, 19, 19, 19, 19, 19, 19, 19, 19, 19, 18, 18, 18, 18, 18, 18, 18, 18, 18, 18, 18, 18, 18, 18, 18, 18, 18, 18, 18, 18, 18, 18, 18, 18, 18, 18]

safemaxUInt128pows = Int8[0, 127, 80, 63, 55, 49, 45, 42, 40, 38, 37, 35, 34, 33, 32, 31, 31, 30, 30, 29, 29, 28, 28, 27, 27, 27, 26, 26, 26, 26, 25, 25, 25, 25, 24, 24, 24, 24, 24, 24, 23, 23, 23, 23, 23, 23, 23, 22, 22, 22, 22, 22, 22, 22, 22, 22, 21, 21, 21, 21, 21, 21, 21, 21, 21, 21, 21, 21, 20, 20, 20, 20, 20, 20, 20, 20, 20, 20, 20, 20, 20, 20, 20, 20, 19, 19, 19, 19, 19, 19, 19, 19, 19, 19, 19, 19, 19, 19, 19, 19, 19, 19, 19, 19, 19, 19, 18, 18, 18, 18, 18, 18, 18, 18, 18, 18, 18, 18, 18, 18, 18, 18, 18, 18, 18, 18, 18, 18]

safemaxUInt64pows = Int8[0, 63, 40, 31, 27, 24, 22, 21, 20, 19, 18, 17, 17, 16, 16, 15, 15, 15, 15, 14, 14, 14, 14, 13, 13, 13, 13, 13, 13, 13, 12, 12, 12, 12, 12, 12, 12, 12, 12, 12, 11, 11, 11, 11, 11, 11, 11, 11, 11, 11, 11, 11, 11, 11, 11, 11, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9]

safemaxUInt32pows = Int8[0, 31, 20, 15, 13, 12, 11, 10, 10, 9, 9, 8, 8, 8, 8, 7, 7, 7, 7, 7, 7, 7, 7, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4]

safemaxUInt16pows = Int8[0, 15, 10, 7, 6, 6, 5, 5, 5, 4, 4, 4, 4, 4, 4, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2]

safemaxUInt8pows = Int8[0, 7, 5, 3, 3, 3, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1]


text:
module ClosedLoopReachability

include("init.jl")
include("activation.jl")
include("network.jl")
include("problem.jl")
include("nnops.jl")
include("setops.jl")
include("split.jl")
include("utils.jl")
include("simulate.jl")
include("solve.jl")

# problem types
export Network,
       Layer,
       Id, Relu, Sigmoid, Tanh,
       ControlledPlant,
       BlackBoxController

# splitters
export BoxSplitter, ZonotopeSplitter,
       IndexedSplitter,
       SignSplitter

# solvers
export solve, forward, forward_network, simulate,
       DeepZ, SampledApprox, VertexSolver, BoxSolver, SplitSolver, BlackBoxSolver

# utility functions
export @modelpath, read_nnet, read_nnet_mat, read_nnet_yaml, read_nnet_sherlock,
       write_nnet_sherlock,
       print_timed

end


text:
using DataFrames
using DatasToolbox

const SEQ_LENGTH = 3

df = DataFrame(A=1:10, B=101:110, C=201:210)

dh = TimeSeriesHandler{Int64}(df, :A, SEQ_LENGTH, input_cols=[:B, :C], output_cols=[:B])

split!(dh, 3)
assign!(dh)

function predict(X)
    println(size(X))
    println(X)
    [X[1, end, 1] + 201]
end

function f(X)
    println("$(size(X)) being called!!!")
    Dict(2=>(X[2] + 1001))
end

info("generating sequence")
# seq = generateSequence(predict, dh, 3, on_matrix=false)
seq = generateSequence(predict, dh, 3, f, on_matrix=false)
info("done generating sequence")

max_stars_repo_name: JuliaTagBot/DatasToolbox

max_issues_repo_name: JuliaTagBot/DatasToolbox

max_forks_repo_name: JuliaTagBot/DatasToolbox


text:
module zfit
"""
Hello world greeting.

# Arguments
    None

# Examples
```julia
greet()
```
"""
function greet()
    print("Hello World!")
end

end # module

max_stars_repo_name: jonas-eschle/zfit-julia

max_issues_repo_name: jonas-eschle/zfit-julia

max_forks_repo_name: jonas-eschle/zfit-julia


text:
import Bio.Align.BAM

function convert( a::Array{UInt8,1})
    s = ""
    for i in a
        s = string(s, Char(i +  0x21))
        ##  println(Int(i))
    end
    return s
end

function get_AlignLen(file_bam)
    count = 0
    reader = open(BAM.Reader, file_bam)

    for record in reader
        if BAM.ismapped(record)

            QualString = convert( BAM.quality(record))
            println(string(BAM.tempname(record), "\t", string(BAM.alignlength(record)), "\t", BAM.refname(record), "\t", BAM.position(record), "\t", BAM.cigar(record), "\t", BAM.seqlength(record), "\t", BAM.mappingquality(record), "\t", BAM.seqlength(record), "\t", BAM.sequence(record), "\t", QualString, "\t", BAM.auxdata(record)["nM"] ))
        end
    end
    close(reader)
    return count
end

AlignNum =  get_AlignLen("C:\\Users\\admin\\Desktop\\s201904405_chrX.bam")
println(AlignNum)

## Ref : https://biojulia.net/Bio.jl/latest/man/reading/#supported-file-formats
## Bio.Align.BAM.flag
## Bio.Align.BAM.ismapped
## Bio.Align.BAM.isprimary
## Bio.Align.BAM.refid
## Bio.Align.BAM.refname
## Bio.Align.BAM.position
## Bio.Align.BAM.rightposition
## Bio.Align.BAM.isnextmapped
## Bio.Align.BAM.nextrefid
## Bio.Align.BAM.nextrefname
## Bio.Align.BAM.nextposition
## Bio.Align.BAM.mappingquality
## Bio.Align.BAM.cigar
## Bio.Align.BAM.cigar_rle
## Bio.Align.BAM.alignment
## Bio.Align.BAM.alignlength
## Bio.Align.BAM.tempname
## Bio.Align.BAM.templength
## Bio.Align.BAM.sequence
## Bio.Align.BAM.seqlength
## Bio.Align.BAM.quality
## Bio.Align.BAM.auxdata


text:
function setupTravelTimeTomography(filenamePrefix::String, Minv,Iact,sback,modFun)
###########################################################
### Read receivers and sources files
RCVfile = string(filenamePrefix,"_rcvMap.dat");
SRCfile = string(filenamePrefix,"_srcMap.dat");

srcNodeMap = readSrcRcvLocationFile(SRCfile,Minv);
rcvNodeMap = readSrcRcvLocationFile(RCVfile,Minv);

Q = generateSrcRcvProjOperators(Minv.n+1,srcNodeMap);
Q = Q.*1/(norm(Minv.h)^2);
P = generateSrcRcvProjOperators(Minv.n+1,rcvNodeMap);

println("Travel time tomography: ",size(Q,2)," sources.");
#############################################################################################

HO = false;

println("Reading data:");

(DobsEik,WdEik) =  readDataFileToDataMat(string(filenamePrefix,"_travelTime.dat"),srcNodeMap,rcvNodeMap);

N = prod(Minv.n+1);

########################################################################################################
##### Set up remote workers ############################################################################
########################################################################################################

EikMPIWorkers = nworkers(); # this just set the maximal MPI workers. To activate parallelism, run addprocs()

(pFor,contDiv,SourcesSubInd) = getEikonalInvParam(Minv,Q,P,HO,EikMPIWorkers);

misfun = SSDFun

Wd 		= Array{Array{Float64}}(length(pFor));
dobs 	= Array{Array{Float64}}(length(pFor));
for i=1:length(pFor)
	I_i = SourcesSubInd[i];
	Wd[i]   = 0.5*WdEik[:,I_i];
	dobs[i] = DobsEik[:,I_i];
end
WdEik = 0;
DobsEik = 0;

pMisRFs = getMisfitParam(pFor, Wd, dobs, misfun, Iact,sback,ones(length(pFor)),modFun);

return (Q,P,pMisRFs,SourcesSubInd,contDiv);
end


text:
"
    onesite dmrg
"
@with_kw struct Dmrg <: Algorithm
    tol::Float64 = Defaults.tol
    maxiter::Int = Defaults.maxiter
    verbose::Bool = Defaults.verbose
    finalize::Function = (iter,state,ham,pars) -> (state,pars);
end

function find_groundstate(state::Union{FiniteMPS,MPSComoving}, H::Hamiltonian,alg::Dmrg,parameters = params(state,H))
    tol=alg.tol;maxiter=alg.maxiter
    iter = 0; delta = 2*tol

    while iter < maxiter && delta > tol
        delta=0.0

        #finalize
        (state,parameters) = alg.finalize(iter,state,H,parameters);

        for pos = [1:(length(state)-1);length(state):-1:2]
            (eigvals,vecs) =  let state = state,parameters = parameters
                eigsolve(state.AC[pos],1,:SR,Lanczos()) do x
                    ac_prime(x,pos,state,parameters)
                end
            end
            delta = max(delta,1-abs(dot(state.AC[pos],vecs[1])))

            state.AC[pos] = vecs[1]
        end

        alg.verbose && @show (iter,delta)
        flush(stdout)

        iter += 1
    end

    return state,parameters,delta
end

"twosite dmrg"
@with_kw struct Dmrg2 <: Algorithm
    tol = Defaults.tol;
    maxiter = Defaults.maxiter;
    trscheme = truncerr(1e-6);
    verbose = Defaults.verbose
end

function find_groundstate(state::Union{FiniteMPS,MPSComoving}, H::Hamiltonian,alg::Dmrg2,parameters = params(state,H))
    tol=alg.tol;maxiter=alg.maxiter
    iter = 0; delta = 2*tol

    while iter < maxiter && delta > tol
        delta=0.0

        ealg = Lanczos()

        for pos=[1:(length(state)-1);length(state)-2:-1:1]
            @tensor ac2[-1 -2; -3 -4]:=state.AC[pos][-1,-2,1]*state.AR[pos+1][1,-3,-4]
            (eigvals,vecs) =  let state=state,parameters=parameters
                eigsolve(x->ac2_prime(x,pos,state,parameters),ac2,1,:SR,ealg)
            end
            newA2center = vecs[1]

            (al,c,ar) = tsvd(newA2center,trunc=alg.trscheme)

            #yeah, we need a different convergence criterium
            @tensor ov[-1,-2,-3,-4]:=al[-1,-2,1]*c[1,2]*ar[2,-3,-4]-al[1,2,3]*c[3,4]*ar[4,5,6]*conj(state.AC[pos][1,2,7])*conj(state.AR[pos+1][7,5,6])*state.AC[pos][-1,-2,9]*state.AR[pos+1][9,-3,-4]
            delta = max(delta,norm(ov))

            state.AC[pos] = (al,complex(c))
            state.AC[pos+1] = (complex(c),_permute_front(ar))
        end

        alg.verbose && @show (iter,delta)
        flush(stdout)
        #finalize
        iter += 1
    end

    return state,parameters,delta
end


text:
module Covid19Inc

using Colors, CSV, DataFrames, Dates, PGFPlotsX
mainpath = normpath(@__DIR__, "..")
datapath = normpath(mainpath, "data")
plotdicts_all = [
    Dict(:num => :CasesFemale, :den => :PopFemale, :lab => "incidens (fall) kvinnor"),
    Dict(:num => :CasesMale, :den => :PopMale, :lab => "incidens (fall) män"),
    Dict(:num => :CasesFemale1, :den => :PopFemale, :lab => "incidens (fall) kvinnor fas 1"),
    Dict(:num => :CasesMale1, :den => :PopMale, :lab => "incidens (fall) män fas 1"),
    Dict(:num => :CasesFemale2, :den => :PopFemale, :lab => "incidens (fall) kvinnor fas 2"),
    Dict(:num => :CasesMale2, :den => :PopMale, :lab => "incidens (fall) män fas 2"),
    Dict(:num => :HospFemale, :den => :PopFemale, :lab => "incidens (sjukhus) kvinnor"),
    Dict(:num => :HospMale, :den => :PopMale, :lab => "incidens (sjukhus) män"),
    Dict(:num => :ICUFemale, :den => :PopFemale, :lab => "incidens (IVA) kvinnor"),
    Dict(:num => :ICUMale, :den => :PopMale, :lab => "incidens (IVA) män"),
    Dict(:num => :DeathsFemale, :den => :PopFemale, :lab => "mortalitet kvinnor"),
    Dict(:num => :DeathsMale, :den => :PopMale, :lab => "mortalitet män")
]
locs = Dict(
    56 => "Belgien", 208 => "Danmark", 380 => "Italien", 578 => "Norge",
    620 => "Portugal", 724 => "Spanien", 752 => "Sverige"
)

wpp2019 = CSV.File(normpath(datapath, "WPP2019_PopulationByAgeSex_Medium.csv")) |> DataFrame

ageaggrs = Dict(
    :i10o90 => DataFrame(AgeGrpStart = wpp2019[:AgeGrpStart][1:21],
        AgeGrpAStart = vcat(map(x->fld(x,10)*10, wpp2019[:AgeGrpStart][1:20]),90)),
    :i10o80 => DataFrame(AgeGrpStart = wpp2019[:AgeGrpStart][1:21],
        AgeGrpAStart = vcat(map(x->fld(x,10)*10, wpp2019[:AgeGrpStart][1:18]),fill(80,3))),
    :i40i10o90 => DataFrame(AgeGrpStart = wpp2019[:AgeGrpStart][1:21],
        AgeGrpAStart = vcat(fill(0, 8), map(x->fld(x,10)*10, wpp2019[:AgeGrpStart][9:20]),90)),
    :i50i10i70i5o90 => DataFrame(AgeGrpStart = wpp2019[:AgeGrpStart][1:21],
        AgeGrpAStart = vcat(fill(0, 10), fill(50, 2), fill(60, 2), 70, 75, 80, 85, fill(90, 3)))
)

function ages10_loct(loc, ageaggr, t)
    inframe = ageaggrs[ageaggr]
    wpp_sub = wpp2019[(wpp2019[:LocID].==loc) .& (wpp2019[:Time].==t),
        [:AgeGrpStart, :PopMale, :PopFemale]]
    by(join(wpp_sub, inframe, on = :AgeGrpStart),
        :AgeGrpAStart,
        PopFemale = :PopFemale => x->(sum(x.*1000)),
        PopMale = :PopMale => x->(sum(x.*1000)))
end

function covid_pop(loc, ageaggr, enddate)
    covid = CSV.File(normpath(datapath, "$(loc)_covid19_$(enddate).csv")) |> DataFrame
    yr = Dates.year(Date(enddate))
    df = join(ages10_loct(loc, ageaggr, yr), covid, on = :AgeGrpAStart)
    df[:AgeGrpAMean] = vcat(map((a,b) -> Int(a+(b-a)*.8), df[:AgeGrpAStart][1:end-1],
        df[:AgeGrpAStart][2:end]), df[:AgeGrpAStart][end]+5)
    Dict(:loc => loc, :enddate => enddate, :df => df)
end

lcd_save(lcd) =
    CSV.write(normpath(datapath, "$(lcd[:loc])_covid19_pop_$(lcd[:enddate]).csv"))

function lcd_load(loc, enddate)
    df = CSV.File("$(loc)_covid19_pop_$(enddate).csv") |> DataFrame
    Dict(:loc => loc, :enddate => enddate, :df => df)
end

function covidpl(lcd)
    df = lcd[:df]
    plotdicts = filter(pd->pd[:num] in names(df), plotdicts_all)
    p = @pgf Axis({ymode="log", xlabel="Ålder",
            title="COVID-19 $(locs[lcd[:loc]]) till $(lcd[:enddate])",
            legend_pos="outer north east", xmajorgrids, ymajorgrids, yminorgrids})
    ages = df[:AgeGrpAMean]
    plotcolors = distinguishable_colors(length(plotdicts)+1, [RGB(1,1,1)])[2:end]
    for (i, pd) in enumerate(plotdicts)
        rate = df[pd[:num]] ./ df[pd[:den]]
        @pgf push!(p, PlotInc({mark_options="solid", color=plotcolors[i]},
            Table([ages, rate])), LegendEntry(pd[:lab]))
    end
    p
end

end # module


text:
include("synapse.jl")
include("compartment.jl")
include("dendrite.jl")
include("axon.jl")
include("soma.jl")

# The Cell is the implementation of a Neuron.
# The cell generates streams of spikes and thus is a bit stream too.
# The stream is written to disc in spans.

mutable struct Cell <: AbstractBitStream
    # Contains a Soma.
    soma::AbstractSoma
    model::Model.ModelData
    samples::Model.Samples

    # Interfaces with InterNeurons (IN)
    # Dendrites Compartments also interact with INs

    function Cell(soma::AbstractSoma, model::Model.ModelData, samples::Model.Samples)
        o = new()
        o.soma = soma
        o.model = model
        o.samples = samples
        o
    end
end

function initialize!(cell::Cell)
    initialize!(cell.soma)
end

function reset!(cell::Cell)
    reset!(cell.soma)
end

function integrate!(cell::Cell, span_t::Int64, t::Int64)
    integrate!(cell.soma, span_t, t)
end

function output(cell::Cell)
    output(cell.soma)
end


text:
using Statistics
using PyCall
pygui(:tk)
import PyPlot; const plt = PyPlot
import InteratomicPotentials as Potentials

n=50
#############################################################################################################################################
#############################################################################################################################################
using JLD
lj_samples = load("examples/Argon13/fit_lj/samples_lj_TEMP_65.jld")["samples"]
lj_samples = hcat(lj_samples...)
snap_samples = load("examples/Argon13/fit_snap/samples_snap_TEMP_65_date_9_30_2021.jld")["samples"]
snap_samples = hcat(snap_samples...)

## Histogram 2d
fig = plt.figure(figsize = (10, 10))
ax = fig.add_subplot(1,1,1)
counts, xedges, yedges, im = ax.hist2d(lj_samples[1, :], lj_samples[2, :], bins = 100, range = [[0.9, 1.1],[0.9, 1.1]], density = true, cmap="Blues")
ax.set_xlabel("σ")
ax.set_yllabel("ϵ")
ax.set_title("Distribution of Fitted LJ Parameters")
plt.gcf()

plt.savefig("distribution_of_lj_parameters_hexbin.png", dpi = 600, bbox_inches = "tight")

## Scatter with Contours
fig = plt.figure(figsize = (10,10))
ax = fig.add_subplot(1,1,1)
ax.hexbin(lj_samples[1, :], lj_samples[2, :], gridsize = 20, cmap="Blues")
ax.set_xlabel("σ")
ax.set_ylabel("ϵ")
ax.set_title("Distribution of LJ Parameters")
plt.gcf()

fig = plt.figure(figsize = (10, 10))
ax = fig.add_subplot(1,1,1)
ax.boxplot(snap_samples[1:30, :]', notch = true)

ax.set_xlabel("Parameter")
ax.set_ylabel("SNAP Parameter Values")
ax.set_title("Distribution of Fitted SNAP Parameters")
plt.gcf()
plt.savefig("boxplot_snap_parameters.png", dpi = 600, bbox_inches = "tight")

fig = plt.figure(figsize = (10, 10))
ax = fig.add_subplot(1,1,1)
ax.plot(1:31, std(snap_samples, dims = 2), linestyle = " ", marker = "*")
ax.set_ylabel("Standard Deviation of SNAP Parameters")
ax.set_xlabel("SNAP Parameters")
ax.set_yscale("log")

plt.gcf()
plt.savefig("standard_deviation_snap_parameters.png", dpi = 600, bbox_inches = "tight")


text:
using DiffEqParamEstim, Test

@time @testset "Tests on ODEs" begin
  include("tests_on_odes/test_problems.jl")
  include("tests_on_odes/l2loss_test.jl")
  include("tests_on_odes/optim_test.jl")
  include("tests_on_odes/lm_test.jl")
  include("tests_on_odes/lsoptim_test.jl")
  include("tests_on_odes/nlopt_test.jl")
  include("tests_on_odes/two_stage_method_test.jl")
  include("tests_on_odes/regularization_test.jl")
  include("tests_on_odes/blackboxoptim_test.jl")
  include("tests_on_odes/weighted_loss_test.jl")
  include("tests_on_odes/l2_colloc_grad_test.jl")
  #include("tests_on_odes/genetic_algorithm_test.jl") # Not updated to v0.6
end

@time @testset "Multiple Shooting Objective" begin include("multiple_shooting_objective_test.jl") end
@time @testset "Likelihood Loss" begin include("likelihood.jl") end
@time @testset "Out-of-place ODE Tests" begin include("out_of_place_odes.jl") end
@time @testset "Steady State Tests" begin include("steady_state_tests.jl") end
@time @testset "DDE Tests" begin include("dde_tests.jl") end
@time @testset "Test on Monte" begin include("test_on_monte.jl") end


text:
"""
    cut_generating_cone(optimizer, A, xfrac, b, π, π0, lower_bounds = false)

Given the primal problem:
```
max c ⋅ x
A x ≥ b
```
and some integrality constraints, a branching rule can be deduced,
of the form `π ⋅ x ≤ π0 ⋁ π ⋅ x ≥ π0 + 1`.
Build the Cut Generating Linear Problem (CGLP) and returns the problem,
dual multipliers and cutting coefficients in a named tuple.
If `lower_bounds`, dual multipliers are added for `x ≥ 0`.
"""
function cut_generating_cone(optimizer, A, xfrac, b, π, π0, lower_bounds = false)
    cglp = Model(optimizer)
    n = length(xfrac)
    size(A, 2) == n == length(π) || throw(DimensionError("A and x must be of fitting dimensions"))
    m = size(A, 1)
    @variables(cglp, begin
        u[1:m] >= 0
        u0 >= 0
        v[1:m] >= 0
        v0 >= 0
        γ[1:n]
        γ0
    end)

    @constraints(cglp, begin
        γ .== A' * u .- u0 .* π
        γ .== A' * v .+ v0 .* π
        γ0 == dot(u, b) - u0 * π0
        γ0 == dot(v, b) + v0 * (π0 + 1)
    end)
    @objective(cglp, Min, dot(γ, xfrac) - γ0)
    return (cglp, (γ = γ, γ0 = γ0, u = u, u0 = u0, v = v, v0 = v0))
end


text:
include("zeppelin.jl")
include("translate.jl")


text:
@testset "DWDataReader.File basics" begin
    @test_throws ArgumentError DWDataReader.File("")
end