text:
## wrappers for multi-gpu functionality in cusolverMg

# auxiliary functionality

# NOTE: in the cublasMg preview, which also relies on this functionality, a separate library
#       called 'cudalibmg' is introduced. factor this out when we actually ship that.

mutable struct MatrixDescriptor
    desc::cudaLibMgMatrixDesc_t

    function MatrixDescriptor(a, grid; rowblocks = size(a, 1), colblocks = size(a, 2), elta=eltype(a) )
        desc = Ref{cudaLibMgMatrixDesc_t}()
        cusolverMgCreateMatrixDesc(desc, size(a, 1), size(a, 2), rowblocks, colblocks, elta, grid)
        return new(desc[])
    end
end

Base.unsafe_convert(::Type{cudaLibMgMatrixDesc_t}, obj::MatrixDescriptor) = obj.desc

mutable struct DeviceGrid
    desc::cudaLibMgGrid_t

    function DeviceGrid(num_row_devs, num_col_devs, deviceIds, mapping)
        @assert num_row_devs == 1 "Only 1-D column block cyclic is supported, so numRowDevices must be equal to 1."
        desc = Ref{cudaLibMgGrid_t}()
        cusolverMgCreateDeviceGrid(desc, num_row_devs, num_col_devs, deviceIds, mapping)
        return new(desc[])
    end
end

Base.unsafe_convert(::Type{cudaLibMgGrid_t}, obj::DeviceGrid) = obj.desc

function allocateBuffers(n_row_devs, n_col_devs, mat::Matrix)
    mat_row_block_size = div(size(mat, 1), n_row_devs)
    mat_col_block_size = div(size(mat, 2), n_col_devs)
    mat_buffers  = Vector{CuMatrix{eltype(mat)}}(undef, ndevices())
    mat_numRows  = Vector{Int64}(undef, ndevices())
    mat_numCols  = Vector{Int64}(undef, ndevices())
    typesize = sizeof(eltype(mat))
    ldas = Vector{Int64}(undef, ndevices())
    current_dev = device()
    mat_cpu_bufs = Vector{Matrix{eltype(mat)}}(undef, ndevices())
    for (di, dev) in enumerate(devices())
        ldas[di]    = mat_col_block_size
        dev_row     = mod(di - 1, n_row_devs) + 1
        dev_col     = div(di - 1, n_row_devs) + 1

        mat_row_inds     = ((dev_row-1)*mat_row_block_size+1):min(dev_row*mat_row_block_size, size(mat, 1))
        mat_col_inds     = ((dev_col-1)*mat_col_block_size+1):min(dev_col*mat_col_block_size, size(mat, 2))
        mat_cpu_bufs[di] = Array(mat[mat_row_inds, mat_col_inds])
    end
    for (di, dev) in enumerate(devices())
        device!(dev)
        mat_gpu_buf = CuMatrix{eltype(mat)}(undef, size(mat_cpu_bufs[di]))
        copyto!(mat_gpu_buf, mat_cpu_bufs[di])
        mat_buffers[di] = mat_gpu_buf
    end
    for (di, dev) in enumerate(devices())
        device!(dev)
        synchronize()
    end
    device!(current_dev)
    return mat_buffers
end

function returnBuffers(n_row_devs, n_col_devs, row_block_size, col_block_size, dDs, D)
    row_block_size = div(size(D, 1), n_row_devs)
    col_block_size = div(size(D, 2), n_col_devs)
    numRows  = [row_block_size for dev in 1:ndevices()]
    numCols  = [col_block_size for dev in 1:ndevices()]
    typesize = sizeof(eltype(D))
    current_dev = device()
    cpu_bufs = Vector{Matrix{eltype(D)}}(undef, ndevices())
    for (di, dev) in enumerate(devices())
        device!(dev)
        dev_row = mod(di - 1, n_row_devs) + 1
        dev_col = div(di - 1, n_row_devs) + 1
        row_inds = ((dev_row-1)*row_block_size+1):min(dev_row*row_block_size, size(D, 1))
        col_inds = ((dev_col-1)*col_block_size+1):min(dev_col*col_block_size, size(D, 2))
        cpu_bufs[di] = Matrix{eltype(D)}(undef, length(row_inds), length(col_inds))
        copyto!(cpu_bufs[di], dDs[di])
    end
    for (di, dev) in enumerate(devices())
        device!(dev)
        synchronize()
        dev_row = mod(di - 1, n_row_devs) + 1
        dev_col = div(di - 1, n_row_devs) + 1
        row_inds = ((dev_row-1)*row_block_size+1):min(dev_row*row_block_size, size(D, 1))
        col_inds = ((dev_col-1)*col_block_size+1):min(dev_col*col_block_size, size(D, 2))
        D[row_inds, col_inds] = cpu_bufs[di]
    end
    device!(current_dev)
    return D
end

## wrappers

function mg_syevd!(jobz::Char, uplo::Char, A; dev_rows=1, dev_cols=ndevices()) # one host-side array A
    dev = device()
    grid = DeviceGrid(1, ndevices(), devices(), CUDALIBMG_GRID_MAPPING_COL_MAJOR)
    if uplo != 'L'
        throw(ArgumentError("only lower fill mode (uplo='L') supported"))
    end
    m, n    = size(A)
    N       = div(size(A, 2), ndevices()) # dimension of the sub-matrix
    lwork         = Ref{Int64}(0)
    workspace     = Vector{CuArray}(undef, ndevices())
    W             = Vector{real(eltype(A))}(undef, n)
    desc = MatrixDescriptor(A, grid; colblocks=N)
    A_arr     = allocateBuffers(dev_rows, dev_cols, A)
    IA            = 1 # for now
    JA            = 1
    GC.@preserve A_arr begin
        cusolverMgSyevd_bufferSize(mg_handle(), jobz, uplo, n, pointer.(A_arr), IA, JA, desc, W, real(eltype(A)), eltype(A), lwork)
    end
    for (di, dev) in enumerate(devices())
        device!(dev)
        workspace[di] = CUDA.zeros(eltype(A), lwork[])
        synchronize()
    end
    device!(dev)
    info = Ref{Cint}(C_NULL)
    GC.@preserve A_arr workspace begin
        cusolverMgSyevd(mg_handle(), jobz, uplo, n, pointer.(A_arr), IA, JA, desc, W, real(eltype(A)), eltype(A), pointer.(workspace), lwork[], info)
    end
    chkargsok(BlasInt(info[]))
    A = returnBuffers(dev_rows, dev_cols, div(size(A, 1), dev_rows), div(size(A, 2), dev_cols), A_arr, A)
    if jobz == 'N'
        return W
    elseif jobz == 'V'
        return W, A
    end
end

function mg_potrf!(uplo::Char, A; dev_rows=1, dev_cols=ndevices()) # one host-side array A
    dev = device()
    grid = DeviceGrid(1, ndevices(), devices(), CUDALIBMG_GRID_MAPPING_COL_MAJOR)
    if uplo != 'L'
        throw(ArgumentError("only lower fill mode (uplo='L') supported"))
    end
    m, n    = size(A)
    N       = div(size(A, 2), ndevices()) # dimension of the sub-matrix
    lwork         = Ref{Int64}(0)
    workspace     = Vector{CuArray}(undef, ndevices())
    desc = MatrixDescriptor(A, grid; colblocks=N)
    A_arr     = allocateBuffers(dev_rows, dev_cols, A)
    IA      = 1 # for now
    JA      = 1
    GC.@preserve A_arr begin
        cusolverMgPotrf_bufferSize(mg_handle(), uplo, n, pointer.(A_arr), IA, JA, desc, eltype(A), lwork)
    end
    for (di, dev) in enumerate(devices())
        device!(dev)
        workspace[di]     = CUDA.zeros(eltype(A), lwork[])
        synchronize()
    end
    device!(dev)
    info = Ref{Cint}(C_NULL)
    GC.@preserve A_arr workspace begin
        cusolverMgPotrf(mg_handle(), uplo, n, pointer.(A_arr), IA, JA, desc, eltype(A), pointer.(workspace), lwork[], info)
    end
    chkargsok(BlasInt(info[]))

    returnBuffers(dev_rows, dev_cols, div(size(A, 1), dev_rows), div(size(A, 2), dev_cols), A_arr, A)
end

function mg_potri!(uplo::Char, A; dev_rows=1, dev_cols=ndevices()) # one host-side array A
    dev = device()
    grid = DeviceGrid(1, ndevices(), devices(), CUDALIBMG_GRID_MAPPING_COL_MAJOR)
    if uplo != 'L'
        throw(ArgumentError("only lower fill mode (uplo='L') supported"))
    end
    m, n    = size(A)
    N       = div(size(A, 2), ndevices()) # dimension of the sub-matrix
    lwork         = Ref{Int64}(0)
    workspace     = Vector{CuArray}(undef, ndevices())
    desc = MatrixDescriptor(A, grid; colblocks=N)
    A_arr     = allocateBuffers(dev_rows, dev_cols, A)
    IA      = 1 # for now
    JA      = 1
    GC.@preserve A_arr begin
        cusolverMgPotri_bufferSize(mg_handle(), uplo, n, pointer.(A_arr), IA, JA, desc, eltype(A), lwork)
    end
    for (di, dev) in enumerate(devices())
        device!(dev)
        workspace[di]     = CUDA.zeros(eltype(A), lwork[])
        synchronize()
    end
    device!(dev)
    info = Ref{Cint}(C_NULL)
    GC.@preserve A_arr workspace begin
        cusolverMgPotri(mg_handle(), uplo, n, pointer.(A_arr), IA, JA, desc, eltype(A), pointer.(workspace), lwork[], info)
    end
    chkargsok(BlasInt(info[]))

    returnBuffers(dev_rows, dev_cols, div(size(A, 1), dev_rows), div(size(A, 2), dev_cols), A_arr, A)
end

function mg_potrs!(uplo::Char, A, B; dev_rows=1, dev_cols=ndevices()) # one host-side array A
    dev = device()
    grid = DeviceGrid(1, ndevices(), devices(), CUDALIBMG_GRID_MAPPING_COL_MAJOR)
    if uplo != 'L'
        throw(ArgumentError("only lower fill mode (uplo='L') supported"))
    end
    ma, na   = size(A)
    mb, nb   = size(B)
    NA       = div(size(A, 2), ndevices()) # dimension of the sub-matrix
    NB       = div(size(B, 2), ndevices()) # dimension of the sub-matrix
    lwork         = Ref{Int64}(0)
    workspace     = Vector{CuArray}(undef, ndevices())
    descA = MatrixDescriptor(A, grid; colblocks=NA)
    descB = MatrixDescriptor(A, grid; colblocks=NB)
    A_arr     = allocateBuffers(dev_rows, dev_cols, A)
    B_arr     = allocateBuffers(dev_rows, dev_cols, B)
    IA      = 1 # for now
    JA      = 1
    IB      = 1 # for now
    JB      = 1
    GC.@preserve A_arr B_arr begin
        cusolverMgPotrs_bufferSize(mg_handle(), uplo, na, nb, pointer.(A_arr), IA, JA, descA, pointer.(B_arr), IB, JB, descB, eltype(A), lwork)
    end
    for (di, dev) in enumerate(devices())
        device!(dev)
        workspace[di]     = CUDA.zeros(eltype(A), lwork[])
        synchronize()
    end
    device!(dev)
    info = Ref{Cint}(C_NULL)
    GC.@preserve A_arr B_arr workspace begin
        cusolverMgPotrs(mg_handle(), uplo, na, nb, pointer.(A_arr), IA, JA, descA, pointer.(B_arr), IB, JB, descB, eltype(A), pointer.(workspace), lwork[], info)
    end
    chkargsok(BlasInt(info[]))

    returnBuffers(dev_rows, dev_cols, div(size(B, 1), dev_rows), div(size(B, 2), dev_cols), B_arr, B)
end

function mg_getrf!(A; dev_rows=1, dev_cols=ndevices()) # one host-side array A
    dev = device()
    grid = DeviceGrid(1, ndevices(), devices(), CUDALIBMG_GRID_MAPPING_COL_MAJOR)
    m, n    = size(A)
    N       = div(size(A, 2), ndevices()) # dimension of the sub-matrix
    lwork         = Ref{Int64}(0)
    ipivs         = Vector{CuVector{Cint}}(undef, ndevices())
    workspace     = Vector{CuArray}(undef, ndevices())
    desc = MatrixDescriptor(A, grid; colblocks=N)
    A_arr     = allocateBuffers(dev_rows, dev_cols, A)
    IA      = 1 # for now
    JA      = 1
    for (di, dev) in enumerate(devices())
        device!(dev)
        ipivs[di]     = CUDA.zeros(Cint, N)
        synchronize()
    end
    device!(dev)
    GC.@preserve A_arr ipivs begin
        cusolverMgGetrf_bufferSize(mg_handle(), m, n, pointer.(A_arr), IA, JA, desc, pointer.(ipivs), eltype(A), lwork)
    end
    device_synchronize()
    for (di, dev) in enumerate(devices())
        device!(dev)
        workspace[di]     = CUDA.zeros(eltype(A), lwork[])
        synchronize()
    end
    device!(dev)
    info = Ref{Cint}(C_NULL)
    GC.@preserve A_arr ipivs workspace begin
        cusolverMgGetrf(mg_handle(), m, n, pointer.(A_arr), IA, JA, desc, pointer.(ipivs), eltype(A), pointer.(workspace), lwork[], info)
    end
    device_synchronize()
    chkargsok(BlasInt(info[]))

    A = returnBuffers(dev_rows, dev_cols, div(size(A, 1), dev_rows), div(size(A, 2), dev_cols), A_arr, A)
    ipiv = Vector{Int}(undef, n)
    for (di, dev) in enumerate(devices())
        device!(dev)
        ipiv[((di-1)*N + 1):min((di*N), n)] = collect(ipivs[di])
    end
    device!(dev)
    return A, ipiv
end

function mg_getrs!(trans, A, ipiv, B; dev_rows=1, dev_cols=ndevices()) # one host-side array A
    dev = device()
    grid = DeviceGrid(1, ndevices(), devices(), CUDALIBMG_GRID_MAPPING_COL_MAJOR)
    ma, na  = size(A)
    mb, nb  = size(B)
    NA      = div(size(A, 2), ndevices()) # dimension of the sub-matrix
    NB      = div(size(B, 2), ndevices()) # dimension of the sub-matrix
    lwork         = Ref{Int64}(0)
    ipivs         = Vector{CuVector{Cint}}(undef, ndevices())
    workspace     = Vector{CuArray}(undef, ndevices())
    descA = MatrixDescriptor(A, grid; colblocks=NA)
    descB = MatrixDescriptor(A, grid; colblocks=NB)
    A_arr     = allocateBuffers(dev_rows, dev_cols, A)
    B_arr     = allocateBuffers(dev_rows, dev_cols, B)
    IA      = 1 # for now
    JA      = 1
    IB      = 1 # for now
    JB      = 1
    for (di, dev) in enumerate(devices())
        device!(dev)
        local_ipiv    = Cint.(ipiv[(di-1)*NA+1:min(di*NA,length(ipiv))])
        ipivs[di]     = CuArray(local_ipiv)
        synchronize()
    end
    device!(dev)
    GC.@preserve A_arr B_arr ipivs begin
        cusolverMgGetrs_bufferSize(mg_handle(), trans, na, nb, pointer.(A_arr), IA, JA, descA, pointer.(ipivs), pointer.(B_arr), IB, JB, descB, eltype(A), lwork)
    end
    for (di, dev) in enumerate(devices())
        device!(dev)
        workspace[di]     = CUDA.zeros(eltype(A), lwork[])
        synchronize()
    end
    device!(dev)
    info = Ref{Cint}(C_NULL)
    GC.@preserve A_arr B_arr ipivs workspace begin
        cusolverMgGetrs(mg_handle(), trans, na, nb, pointer.(A_arr), IA, JA, descA, pointer.(ipivs), pointer.(B_arr), IB, JB, descB, eltype(A), pointer.(workspace), lwork[], info)
    end
    chkargsok(BlasInt(info[]))

    returnBuffers(dev_rows, dev_cols, div(size(B, 1), dev_rows), div(size(B, 2), dev_cols), B_arr, B)
end


text:
# Benchmark Catlab.Graphs against Graphs.jl and MetaGraphs.jl.
using BenchmarkTools
const SUITE = BenchmarkGroup()

using Random
import Graphs as SimpleGraphs
import MetaGraphs
const LG, MG = SimpleGraphs, MetaGraphs

using Catlab, Catlab.CategoricalAlgebra, Catlab.Graphs
using Catlab.Graphs.BasicGraphs: TheoryGraph, TheoryLabeledGraph
using Catlab.WiringDiagrams: query
using Catlab.Programs: @relation

testdatadir = joinpath(dirname(@__FILE__), "..", "test", "testdata")

# Example Graphs
#
################

# Stolen from the Lightgraphs benchmark suite

dg1fn = joinpath(testdatadir, "graph-50-500.jgz")

LG_GRAPHS = Dict{String,LG.DiGraph}(
    "complete100"   => LG.complete_digraph(100),
    # "5000-50000"    => LG.loadgraph(dg1fn)["graph-5000-50000"],
    "path500"       => LG.path_digraph(500)
)

GRAPHS = Dict(k => Graph(g) for (k,g) in LG_GRAPHS)

LG_SYMGRAPHS = Dict{String,LG.Graph}(
    "complete100"   => LG.complete_graph(100),
    "tutte"         => LG.smallgraph(:tutte),
    "path500"       => LG.path_graph(500),
    # "5000-49947"    => LG.SimpleGraph(DIGRAPHS["5000-50000"])
)

SYMGRAPHS = Dict(k => SymmetricGraph(g) for (k,g) in LG_SYMGRAPHS)

# Helpers
#########

# `bench_iter_edges` and `bench_has_edge` adapted from Graphs:
# https://github.com/JuliaGraphs/Graphs.jl/blob/master/benchmark/core.jl

function bench_iter_edges(g)
  count = 0
  for e in edges(g)
    s, t = src(g,e), tgt(g,e)
    count += 1
  end
  count
end
function bench_iter_edges_vectorized(g)
  count = 0
  for (s,t) in zip(src(g), tgt(g))
    count += 1
  end
  count
end

function bench_has_edge(g)
  Random.seed!(1)
  n = nv(g)
  srcs, tgts = rand(1:n, n÷4), rand(1:n, n÷4)
  count = 0
  for (s, t) in zip(srcs, tgts)
    if has_edge(g, s, t)
      count += 1
    end
  end
  count
end

function bench_iter_neighbors(g)
  count = 0
  for v in vertices(g)
    count += length(neighbors(g,v))
  end
  count
end

@inline Graphs.nv(g::LG.AbstractGraph) = LG.nv(g)
@inline Graphs.vertices(g::LG.AbstractGraph) = LG.vertices(g)
@inline Graphs.edges(g::LG.AbstractGraph) = LG.edges(g)
@inline Graphs.src(g::LG.AbstractGraph, e::LG.AbstractEdge) = LG.src(e)
@inline Graphs.tgt(g::LG.AbstractGraph, e::LG.AbstractEdge) = LG.dst(e)
@inline Graphs.has_edge(g::LG.AbstractGraph, args...) = LG.has_edge(g, args...)
@inline Graphs.neighbors(g::LG.AbstractGraph, args...) = LG.neighbors(g, args...)

function Graphs.connected_component_projection(g::LG.AbstractGraph)
  label = Vector{Int}(undef, LG.nv(g))
  LG.connected_components!(label, g)
  normalized = searchsortedfirst.(Ref(unique!(sort(label))), label)
  FinFunction(normalized)
end

abstract type FindTrianglesAlgorithm end
struct TriangleBacktrackingSearch <: FindTrianglesAlgorithm end
struct TriangleQuery <: FindTrianglesAlgorithm end

""" Number of triangles in a graph.
"""
function ntriangles(g::T, ::TriangleBacktrackingSearch) where T
  triangle = T(3)
  add_edges!(triangle, [1,2,1], [2,3,3])
  length(homomorphisms(triangle, g, alg=BacktrackingSearch()))
end
function ntriangles(g, ::TriangleQuery)
  length(query(g, ntriangles_query))
end
const ntriangles_query = @relation (;) begin
  E(src=v1, tgt=v2)
  E(src=v2, tgt=v3)
  E(src=v1, tgt=v3)
end

# Graphs
########

bench = SUITE["Graph"] = BenchmarkGroup()
clbench = bench["Catlab"] = BenchmarkGroup()
clvecbench = bench["Catlab-vectorized"] = BenchmarkGroup()
lgbench = bench["LightGraphs"] = BenchmarkGroup()

n = 10000
clbench["make-path"] = @benchmarkable path_graph(Graph,n)
lgbench["make-path"] = @benchmarkable begin
  g = LG.DiGraph()
  LG.add_vertices!(g, n)
  for v in 1:(n-1)
    LG.add_edge!(g, v, v+1)
  end
end

g = path_graph(Graph, n)
lg = LG.DiGraph(g)

clbench["iter-edges"] = @benchmarkable bench_iter_edges($g)
clvecbench["iter-edges"] = @benchmarkable bench_iter_edges_vectorized($g)
lgbench["iter-edges"] = @benchmarkable bench_iter_edges($lg)
clbench["has-edge"] = @benchmarkable bench_has_edge($g)
lgbench["has-edge"] = @benchmarkable bench_has_edge($lg)
clbench["iter-neighbors"] = @benchmarkable bench_iter_neighbors($g)
lgbench["iter-neighbors"] = @benchmarkable bench_iter_neighbors($lg)

bench = SUITE["GraphConnComponents"] = BenchmarkGroup()
clbench = bench["Catlab"] = BenchmarkGroup()
lgbench = bench["LightGraphs"] = BenchmarkGroup()

n₀ = 2000
g₀ = path_graph(Graph, n₀)
g = ob(coproduct(fill(g₀, 5)))
lg = LG.DiGraph(g)
clbench["path-graph"] =
  @benchmarkable connected_component_projection_bfs($g)
lgbench["path-graph"] =
  @benchmarkable connected_component_projection($lg)

g₀ = star_graph(Graph, n₀)
g = ob(coproduct(fill(g₀, 5)))
lg = LG.DiGraph(g)
clbench["star-graph"] =
  @benchmarkable connected_component_projection_bfs($g)
lgbench["star-graph"] =
  @benchmarkable connected_component_projection($lg)

for gn in keys(GRAPHS)
  clbench[gn] = @benchmarkable connected_component_projection_bfs($(GRAPHS[gn]))
  lgbench[gn] = @benchmarkable connected_component_projection($(LG_GRAPHS[gn]))
end

bench = SUITE["GraphTriangles"] = BenchmarkGroup()
clbench = bench["Catlab"] = BenchmarkGroup()
lgbench = bench["LightGraphs"] = BenchmarkGroup()

n = 100
g = wheel_graph(Graph, n)
lg = LG.DiGraph(g)
clbench["wheel-graph-triangles-hom"] =
  @benchmarkable ntriangles($g, TriangleBacktrackingSearch())
clbench["wheel-graph-triangles"] =
  @benchmarkable ntriangles($g, TriangleQuery())

# Symmetric graphs
##################

bench = SUITE["SymmetricGraph"] = BenchmarkGroup()
clbench = bench["Catlab"] = BenchmarkGroup()
clvecbench = bench["Catlab-vectorized"] = BenchmarkGroup()
lgbench = bench["LightGraphs"] = BenchmarkGroup()

n = 10000
clbench["make-path"] = @benchmarkable path_graph(SymmetricGraph, n)
lgbench["make-path"] = @benchmarkable begin
  g = LG.Graph()
  LG.add_vertices!(g, n)
  for v in 1:(n-1)
    LG.add_edge!(g, v, v+1)
  end
end

g = path_graph(SymmetricGraph, n)
lg = LG.Graph(g)

clbench["iter-edges"] = @benchmarkable bench_iter_edges($g)
clvecbench["iter-edges"] = @benchmarkable bench_iter_edges_vectorized($g)
lgbench["iter-edges"] = @benchmarkable bench_iter_edges($lg)
clbench["has-edge"] = @benchmarkable bench_has_edge($g)
lgbench["has-edge"] = @benchmarkable bench_has_edge($lg)
clbench["iter-neighbors"] = @benchmarkable bench_iter_neighbors($g)
lgbench["iter-neighbors"] = @benchmarkable bench_iter_neighbors($lg)

bench = SUITE["SymmetricGraphConnComponent"] = BenchmarkGroup()
clbench = bench["Catlab"] = BenchmarkGroup()
lgbench = bench["LightGraphs"] = BenchmarkGroup()

n₀ = 2000
g₀ = path_graph(SymmetricGraph, n₀)
g = ob(coproduct(fill(g₀, 5)))
lg = LG.Graph(g)
clbench["path-graph-components"] =
  @benchmarkable connected_component_projection($g)
lgbench["path-graph-components"] =
  @benchmarkable connected_component_projection($lg)

g₀ = star_graph(SymmetricGraph, n₀)
g = ob(coproduct(fill(g₀, 5)))
lg = LG.Graph(g)
clbench["star-graph-components"] =
  @benchmarkable connected_component_projection($g)
lgbench["star-graph-components"] =
  @benchmarkable connected_component_projection($lg)

for gn in keys(SYMGRAPHS)
  clbench[gn] = @benchmarkable connected_component_projection_bfs($(SYMGRAPHS[gn]))
  lgbench[gn] = @benchmarkable connected_component_projection($(LG_SYMGRAPHS[gn]))
end

bench = SUITE["SymmetricGraphTriangles"] = BenchmarkGroup()
clbench = bench["Catlab"] = BenchmarkGroup()
lgbench = bench["LightGraphs"] = BenchmarkGroup()

n = 100
g = wheel_graph(SymmetricGraph, n)
lg = LG.Graph(g)
clbench["wheel-graph-triangles-hom"] =
  @benchmarkable ntriangles($g, TriangleBacktrackingSearch())
# clbench["wheel-graph-triangles-query"] =
clbench["wheel-graph-triangles"] =
  @benchmarkable ntriangles($g, TriangleQuery())
lgbench["wheel-graph-triangles"] = @benchmarkable sum(LG.triangles($lg))

# Weighted graphs
#################

bench = SUITE["WeightedGraph"] = BenchmarkGroup()
clbench = bench["Catlab"] = BenchmarkGroup()
clvecbench = bench["Catlab-vectorized"] = BenchmarkGroup()
lgbench = bench["MetaGraphs"] = BenchmarkGroup()

n = 10000
g = path_graph(WeightedGraph{Float64}, n; E=(weight=range(0,1,length=n-1),))
mg = MG.MetaDiGraph(g)

clvecbench["sum-weights"] = @benchmarkable sum(weight($g))
clbench["sum-weights"] = @benchmarkable begin
  total = 0.0
  for e in edges($g)
    total += weight($g, e)
  end
  total
end
lgbench["sum-weights"] = @benchmarkable begin
  total = 0.0
  for e in MG.edges($mg)
    total += MG.get_prop($mg, e, :weight)
  end
  total
end

clvecbench["increment-weights"] = @benchmarkable begin
  $g[:weight] .= $g[:weight] .+ 1.0
end

clbench["increment-weights"] = @benchmarkable begin
  for e in edges($g)
    $g[e,:weight] += 1.0
  end
end
lgbench["increment-weights"] = @benchmarkable begin
  for e in MG.edges($mg)
    MG.set_prop!($mg, e, :weight, MG.get_prop($mg, e, :weight) + 1.0)
  end
end

# Labeled graphs
################

bench = SUITE["LabeledGraph"] = BenchmarkGroup()
clbench = bench["Catlab"] = BenchmarkGroup()
lgbench = bench["LightGraphs"] = BenchmarkGroup()

@acset_type IndexedLabeledGraph(TheoryLabeledGraph, index=[:src,:tgt],
                                unique_index=[:label]) <: AbstractLabeledGraph

function discrete_labeled_graph(n::Int; indexed::Bool=false)
  g = (indexed ? IndexedLabeledGraph{String} : LabeledGraph{String})()
  add_vertices!(g, n, label=("v$i" for i in 1:n))
  g
end

function discrete_labeled_metagraph(n::Int; indexed::Bool=false)
  mg = MG.MetaDiGraph()
  for i in 1:n
    MG.add_vertex!(mg, :label, "v$i")
  end
  if indexed; MG.set_indexing_prop!(mg, :label) end
  mg
end

n = 5000
clbench["make-discrete"] = @benchmarkable discrete_labeled_graph($n)
lgbench["make-discrete"] = @benchmarkable discrete_labeled_metagraph($n)
clbench["make-discrete-indexed"] =
  @benchmarkable discrete_labeled_graph($n, indexed=true)
lgbench["make-discrete-indexed"] =
  @benchmarkable discrete_labeled_metagraph($n, indexed=true)

n = 10000
g = discrete_labeled_graph(n)
mg = discrete_labeled_metagraph(n)
clbench["iter-labels"] = @benchmarkable begin
  for v in vertices($g)
    label = $g[v,:label]
  end
end
lgbench["iter-labels"] = @benchmarkable begin
  for v in MG.vertices($mg)
    label = MG.get_prop($mg, v, :label)
  end
end

g = discrete_labeled_graph(n, indexed=true)
mg = discrete_labeled_metagraph(n, indexed=true)
Random.seed!(1)
σ = randperm(n)
clbench["indexed-lookup"] = @benchmarkable begin
  for i in $σ
    @assert incident($g, "v$i", :label) == i
  end
end
lgbench["indexed-lookup"] = @benchmarkable begin
  for i in $σ
    @assert $mg["v$i", :label] == i
  end
end

# Random Graphs
###############

bench = SUITE["RandomGraph"] = BenchmarkGroup()
clbench = bench["Catlab"] = BenchmarkGroup()
lgbench = bench["LightGraphs"] = BenchmarkGroup()

sizes = [10000]
ps = [0.001]
for size in sizes, p in ps
  clbench["erdos_renyi-$size-$p"] =
    @benchmarkable erdos_renyi($Graph, $size, $(p/2))
  lgbench["erdos_renyi-$size-$p"] =
    @benchmarkable LG.erdos_renyi($size, $p)
end

ks = [10]

for size in sizes, k in ks
  clbench["expected_degree_graph-$size-$k"] =
    @benchmarkable expected_degree_graph($Graph, $([min(k,size-1) for _ in 1:size]))
  lgbench["expected_degree_graph-$size-$k"] =
    @benchmarkable LG.expected_degree_graph($([min(k,size-1) for _ in 1:size]))
end

for size in sizes, k in ks
  clbench["watts_strogatz-$size-$k"] =
    @benchmarkable watts_strogatz($Graph, $size, $(min(k,size-1)), 0.5)
  lgbench["watts_strogatz-$size-$k"] =
    @benchmarkable LG.watts_strogatz($size, $(min(k,size-1)), 0.5)
end

# Searching
###########

bench = SUITE["Searching"] = BenchmarkGroup()
clbench = bench["Catlab"] = BenchmarkGroup()
lgbench = bench["LightGraphs"] = BenchmarkGroup()

for size in sizes, p in ps
  local g = erdos_renyi(Graph, size, p)
  local lg = LG.SimpleDiGraph(g)
  clbench["bfs_erdos_renyi-$size-$p"] = @benchmarkable bfs_parents($g,1)
  lgbench["bfs_erdos_renyi-$size-$p"] = @benchmarkable LG.bfs_parents($lg,1)
  clbench["dfs_erdos_renyi-$size-$p"] = @benchmarkable dfs_parents($g,1)
  lgbench["dfs_erdos_renyi-$size-$p"] = @benchmarkable LG.dfs_parents($lg,1)
end

ext: jl

max_forks_repo_name: epatters/CompCat


text:
#!/usr/bin/env julia
#=
   Author: Tim Sterne-Weiler, 3/16/2015
   e-mail: tim.sterne.weiler@utoronto.ca

=#
using GZip

# slow performance of HDF5 has promted [savedict, loaddict]:
function savedict{K,V}( file::AbstractString, dict::Dict{K,V} )
  GZip.open(file, "w") do fh
    println(fh, "@Dict{" * string(K) * "," * string(V) * "}:" * string(length(dict)))
    for k in keys(dict)
      println(fh, string(k) * "˧" * string(dict[k]))
    end
  end
end #--> nothing

function loaddict(file::AbstractString)
   myconvert(x::Type{Char}, y) = string(y)[1]
   myconvert{T <: AbstractString}(x::Type{T}, y) = convert(T, string(y))
   myconvert{T <: Number, S <: AbstractString}(x::Type{T}, y::S) = parse(T, y)
   myconvert(x, y) = convert(x, y)
   #hint = parse(chomp(readall(pipe(`zcat $file`, `wc -l`)))) # --deprecated 6/10/15 tsw
   dict = nothing
   GZip.open(file, "r") do fh
    heads = split(chomp( readline(fh) ), r"@Dict{|,|}:")
    @assert( length(heads) == 4 )
    @assert( ismatch(r"^[A-Z|a-z|0-9]+$", heads[2]) ) # must look like a type
    @assert( ismatch(r"^[A-Z|a-z|0-9]+$", heads[3]) ) # before we allow parse()
    ktype = eval(parse(heads[2]))
    vtype = eval(parse(heads[3]))
    dict = Dict{ktype,vtype}()
    hint  = parse(Int, heads[4])
    sizehint!(dict, hint)
    for l in eachline(fh)
      sub = split(chomp(l), '˧')
      key = myconvert(ktype, sub[1])
      val = myconvert(vtype, sub[2])
      dict[key] = val
    end
  end
  dict
end #--> Dict{K,V}

# dictionary value increment
# this should be standard imo.
function dinc!{K,V}(dict::Dict{K,V}, key::K, val::V=one(V))
  indx = Base.ht_keyindex(dict, key)
  if indx <= 0
    dict[key] = val
  else
    @assert( isa(dict.vals[indx], Number) )
    dict.vals[indx] += val
  end
end #--> nothing

# dictionary [] create or push if exists
# this also should be standard imo
function dush!{K,V}(dict::Dict{K,V}, key::K, val; arrType = Any)
  indx = Base.ht_keyindex(dict, key)
  @assert( isa(val, arrType) )
  if indx <= 0
    dict[key] = arrType[ val ]
  else
    push!(dict.vals[indx], val)
  end
end #--> nothing

# normalize numeric dictionary values
function dnorm!{K, V<:Number}(dict::Dict{K,V}, n = sum(collect(values(dict))) )
  @assert( n > 0 )
  for i in keys(dict)
    dict[i] /= n
  end #normalize
end #--> nothing

# normalize numeric dictionary values
function dnorm!(dict::Dict, n = sum(collect(values(dict))) )
  @assert( n > 0 )
  for i in keys(dict)
    @assert( isa(dict[i], Number) )
    dict[i] /= n
  end #normalize
end #--> nothing

max_forks_repo_name: timbitz/aligater


text:
#! /usr/bin/env julia --color=yes

using Base.Filesystem
using Pkg

scriptsDir = abspath(dirname(@__FILE__))
baseDir = abspath(joinpath(scriptsDir, "../"))
srcDir = abspath(joinpath(baseDir, "src/"))
docsDir = abspath(joinpath(scriptsDir, "../docs/"))
docSrcDir = abspath(joinpath(docsDir, "src/"))
buildDir = abspath(joinpath(scriptsDir, "../docs/build/"))

Pkg.activate(docsDir)

using FileWatching
using Logging
using HTTP
using Sockets

function watchDir(changes,dir)
  @async begin
    while true
      @info ">>>  Waiting for changes in $dir..."
      change = watch_folder(dir)
      put!(changes, change)
    end
  end
end

function contentType(filePath)
  if filePath[end-4:end] == ".html"
    return "text/html"
  elseif filePath[end-2:end] == ".js"
    return "text/javascript"
  elseif filePath[end-3:end] == ".css"
    return "text/css"
  elseif (filePath[end-3:end] == ".jpg" || filePath[end-4:end] == ".jpeg")
    return "image/jpeg"
  elseif filePath[end-3:end] == ".png"
    return "image/png"
  elseif filePath[end-4:end] == ".json"
    return "application/json"
  elseif filePath[end-3:end] == ".ico"
    return "image/vnd.microsoft.icon"
  else
    return "application/octet-stream"
  end
end

function serveDir(dir)
  @async begin
    HTTP.serve(Sockets.localhost,8080) do req::HTTP.Request
        path = req.target
        @info "Serving $path"
        if path[begin] == '/'
          path = path[begin+1:end]
        end
        filePath = abspath(joinpath(buildDir,path))
        if filePath[end] == '/'
          filePath *= "index.html"
        end
        if ispath(filePath)
          response = HTTP.Response(200, read(filePath))
          HTTP.setheader(response,"Content-Type" => contentType(filePath))
          return response
        else
          return HTTP.Response(404,"Not found")
        end
    end
  end
end

@sync begin
  changes = Channel()
  watchDir(changes, docSrcDir)
  watchDir(changes, srcDir)
  serveDir(buildDir)
  @async begin
    @info ">>>  Triggering initial generation of documentation"
    put!(changes, :first)
  end
  for change in changes
    cd(docsDir) do
      try
        run(`julia --color=yes make.jl`)
        @info ">>>  Generation complete"
      catch ex
        @error "!!!  Error generating docs: $ex"
      end
    end
  end
end


text:
@inline function to_ecode(x::Int32)
    if x === Int32(2)
        return :DICOM
    elseif x === Int32(4)
        return :AFNI
    elseif x === Int32(6)
        return :Comment
    elseif x === Int32(8)
        return :XCEDE
    elseif x === Int32(10)
        return :JimDimInfo
    elseif x === Int32(12)
        return :WorkflowFWDS
    elseif x === Int32(14)
        return :Freesurfer
    elseif x === Int32(16)
        return :PyPickly
    elseif x === Int32(18)
        return :MindIdent
    elseif x === Int32(20)
        return :BValue
    elseif x === Int32(22)
        return :SphericalDirection
    elseif x === Int32(24)
        return :DTComponent
    elseif x === Int32(26)
        return :SHCDegreeOrder
    elseif x === Int32(28)
        return :Voxbo
    elseif x === Int32(30)
        return :Caret
    elseif x === Int32(32)
        return :CIfTI
    elseif x === Int32(34)
        return :VariableFrameTiming
    elseif x === Int32(38)
        return :Eval
    elseif x === Int32(40)
        return :Matlab
    else
        return :Ignore
    end
end

@inline function ecode_to_int(x::Symbol)
    if x === :DICOM
        return Int32(2)
    elseif x === :AFNI
        return Int32(4)
    elseif x === :Comment
        return Int32(6)
    elseif x === :XCEDE
        return Int32(8)
    elseif x === :JimDimInfo
        return Int32(10)
    elseif x === :WorkflowFWDS
        return Int32(12)
    elseif x === :Freesurfer
        return Int32(14)
    elseif x === :PyPickly
        return Int32(16)
    elseif x === :MindIdent
        return Int32(18)
    elseif x === :BValue
        return Int32(20)
    elseif x === :SphericalDirection
        return Int32(22)
    elseif x === :DTComponent
        return Int32(24)
    elseif x === :SHCDegreeOrder
        return Int32(26)
    elseif x === :Voxbo
        return Int32(28)
    elseif x === :Caret
        return Int32(30)
    elseif x === :CIfTI
        return Int32(32)
    elseif x === :VariableFrameTiming
        return Int32(34)
    elseif x === :Eval
        return Int32(38)
    elseif x === :Matlab
        return Int32(40)
    else
        return Int32(0)
    end
end

"""
    NIfTIExtension

* `ecode::Int32`: extension code that can be interpreted used `ecode(::NIfTIExtension)`
* `edata::Vector{UInt8}`: raw data with no byte-swapping
"""
struct NIfTIExtension
    ecode::Int32
    edata::Vector{UInt8}
end

"""
    NIfTI Extension Codes

* Ignore
* DICOM: intended for raw DICOM attributes
* AFNI: Robert W Cox: rwcox@nih.gov https://afni.nimh.nih.gov/afni
* Comment: plain ASCII text only
* XCEDE: http://www.nbirn.net/Resources/Users/Applications/xcede/index.html
* JIMDimInfo: Dimensionalinformation for the JIM software (XML format)
* WorkflowFWDS:
* Freesurfer:
* PyPickle: embedded Python objects
* MiNDIdent: LONI MiND codes: http://www.loni.ucla.edu/twiki/bin/view/Main/MiND
* BValue
* SphericalDirection
* DTComponent
* SHCDegreeOrder
* Voxbo: www.voxbo.org
* Caret: http://brainvis.wustl.edu/wiki/index.php/Caret:Documentation
* CIfTI
* VariableFrameTiming
* AgilentProcpar
* Eval
* Matlab
"""
ecode(x::NIfTIExtension) = to_ecode(x.ecode)
ecode(x::Vector{NIfTIExtension}) = map(ecode, x)

# Calculates the size of a NIfTI extension
"""
    esize(ex::NIfTIExtension)

NIfTIExtensions should be of a byte size that is a mulitple of 16. This includes
raw encoding of the the `ecode` (as an Int32) and the esize itself (also as an
Int32). Therefore, `8 + sizeof(ex.edata)` should be divisible by 16.
"""
function esize(ex::NIfTIExtension)
    ret = 8 + length(ex.edata)
    @assert ret%16 != 0 "NIfTIExtension has innapropriate size. See docstrings for more details."
    return ret
end

function read_extensions(io, n)
    ret = NIfTIExtension[]
    if eof(io)
        return ret
    else
        b1 = read(io, UInt8)
        # GZIP doesn't skip so we read and throw away
        read(io, UInt8)
        read(io, UInt8)
        read(io, UInt8)
        if b1 === zero(UInt8)
            return ret
        else
            counter = 0
            while counter < (n - 1)
                esize = read(io, Int32)
                ec = read(io, Int32)
                push!(ret, NIfTIExtension(ec, read!(io, Array{UInt8}(undef, esize - 8))))
                counter += esize
            end
            return ret
        end
    end
end

function write(io::IO, x::Vector{NIfTIExtension})
    if isempty(x)
        write(io, fill(UInt8(0), 4))
    else
        write(io, UInt8[1, 0, 0, 0])
        for ex in x
            write(io, Int32(esize(ex)))
            write(io, ex.ecode)
            write(io, ex.edata)
#            write(io, zeros(UInt8, esize(ex) - length(ex.edata)))
        end
    end
end


text:
abstract type Mutability end

struct Mutable <: Mutability
end

struct Immutable <: Mutability
end

"""
    AbstractValue{M<:Mutability}

A value can be marked as either Mutable or Immutable.

An immutable value can be assigned only once.
A mutable value can be assigned more than once.
"""
abstract type AbstractValue{M<:Mutability} end

const MutableValue = AbstractValue{Mutable}
const ImmutableValue = AbstractValue{Immutable}

@generated function is_mutable(::Type{T}) :: Bool where {M<:Mutability, T<:AbstractValue{M}}
    M == Mutable
end

@generated function is_immutable(::Type{T}) :: Bool where {M<:Mutability, T<:AbstractValue{M}}
    M == Immutable
end

@generated function is_mutable(v::AbstractValue) :: Bool
    is_mutable(v)
end

@generated function is_immutable(v::AbstractValue) :: Bool
    is_immutable(v)
end

struct NullPointer <: ImmutableValue
end

"""
A pointer to an instruction
that computes a value.
"""
struct SSAValue <: ImmutableValue
    address::Int # instruction location
end

"""
Constant value to be encoded directly into the IR.
The address is the value itself.
"""
struct Const{T} <: ImmutableValue
    val::T
end

"""
    Variable{M<:Mutability}

Creates a variable identified by a symbol
that can be either mutable or immutable.

# Alias

For convenience, the following constants are defined in the package:

```julia
const MutableVariable = Variable{Mutable}
const ImmutableVariable = Variable{Immutable}
```

# Examples

```julia
m = OptimizingIR.MutableVariable(:varmut) # a mutable variable
im = OptimizingIR.ImmutableVariable(:varimut) # an immutable variable
```
"""
struct Variable{M} <: AbstractValue{M}
    symbol::Symbol
end

const MutableVariable = Variable{Mutable}
const ImmutableVariable = Variable{Immutable}

"""
Sets which optimizations are
allowed to an `Op`.
"""
struct OptimizationRule{T, M<:Union{Integer, Tuple, NTuple, Nothing}}
    pure::Bool              # whether the Op itself is pure or impure
    commutative::Bool
    hasleftidentity::Bool   # [left=element] op      right
    hasrightidentity::Bool  #     left       op [right=element]
    identity_element::T
    mutable_arg::M          # which function arguments are mutated (check if mutable)

    function OptimizationRule(pure::Bool, commutative::Bool,
            hasleftidentity::Bool, hasrightidentity::Bool, identity_element::T,
            mutable_arg::M) where {T, M}

        if commutative || hasleftidentity || hasrightidentity
            @assert pure "Can't apply commutative or identity optimization on impure op."
        end

        return new{T,M}(pure, commutative, hasleftidentity, hasrightidentity, identity_element, mutable_arg)
    end
end

# typeof(O) == OptimizationRule
struct Op{F<:Function, O}
    op::F

    function Op(f::F, o::OptimizationRule) where {F<:Function}
        return new{F, o}(f)
    end
end

"""
    Op(f::Function;
            pure::Bool=false,
            commutative::Bool=false,
            hasleftidentity::Bool=false,
            hasrightidentity::Bool=false,
            identity_element::T=NULL_IDENTITY_ELEMENT,
            mutable_arg=nothing)

Defines a basic instruction with optimization annotations.

# Arguments

* `f` is a Julia function to be executed by the `Op`.

* `pure`: marks the function as pure (`true`) or impure (`false`) .

* `commutative`: marks the `Op` as commutative.

* `hasleftidentity`: marks the `Op` as having an identity when operating from the left, which means that `f(I, v) = v`, where `I` is the `identity_element`.

* `hasrightidentity`: marks the `Op` as having an identity when operating from the right, which means that `f(v, I) = v`, where `I` is the `identity_element`.

* `mutable_arg`: either `nothing`, or the index or tuple of indexes of the arguments that need to be mutable.

# Purity

A function is considered [pure](https://en.wikipedia.org/wiki/Pure_function)
if its return value is the same for the same arguments, and has no side-effects.

Operations marked as **pure** are suitable for [Value-Numbering](https://en.wikipedia.org/wiki/Value_numbering)
optimization.

When marked as impure, all optimization passes are disabled.

# Examples

```julia
# this Op allows using `+` as a pure, commutative function. Sets 0 as identity element.
const OP_SUM = OIR.Op(+, pure=true, commutative=true, hasleftidentity=true, hasrightidentity=true, identity_element=0)

# `-` as a pure function. Identity is zero but only `x - 0 = x` case is checked.
const OP_SUB = OIR.Op(-, pure=true, hasrightidentity=true, identity_element=0)

# `*` as pure commutative function. Sets 1 as identity element.
const OP_MUL = OIR.Op(*, pure=true, commutative=true, hasleftidentity=true, hasrightidentity=true, identity_element=1)

# `/` as a pure function. Identity is checked to the right: `a / 1 = a`.
const OP_DIV = OIR.Op(/, pure=true, hasrightidentity=true, identity_element=1)

# power function
const OP_POW = OIR.Op(^, pure=true, hasrightidentity=true, identity_element=1)

# an Op that uses an arbitrary Julia function
foreign_fun(a, b, c) = a^3 + b^2 + c
const OP_FOREIGN_FUN = OIR.Op(foreign_fun, pure=true)

# An Op that is impure: every time we run `zeros` a different Array is returned.
const OP_ZEROS = OIR.Op(zeros)
```
"""
function Op(f::F;
            pure::Bool=false,
            commutative::Bool=false,
            hasleftidentity::Bool=false,
            hasrightidentity::Bool=false,
            identity_element::T=NULL_IDENTITY_ELEMENT,
            mutable_arg=nothing) where {F<:Function, T}

    return Op(f, OptimizationRule(pure, commutative, hasleftidentity, hasrightidentity, identity_element, mutable_arg))
end

# defines functor for Op so that op(arg1, ...) will call op.op(arg1, ...)
# see https://docs.julialang.org/en/v1/manual/methods/#Function-like-objects-1
(op::Op)(args...) = op.op(args...)

abstract type AbstractCall end
abstract type AbstractOpCall{OP<:Op} <: AbstractCall end
abstract type Instruction end
abstract type LinearInstruction{T<:AbstractCall} <: Instruction end
abstract type BranchInstruction <: Instruction end

# lhs = rhs
struct Assignment{V1<:Variable, V2<:ImmutableValue} <: AbstractCall
    lhs::V1
    rhs::V2
end

"""
A `PureInstruction` is a call to an operation
that always returns the same value
if the same arguments are passed to the instruction.
It is suitable for memoization, in the sense that
it can be optimized in the Value-Number algorithm
inside a Basic Block.

An instruction is considered pure if its `call`
has an pure `Op` and all of its arguments
are immutable.
"""
struct PureInstruction{T} <: LinearInstruction{T}
    call::T

    function PureInstruction(call::G) where {G<:AbstractOpCall}
        # An AbstractOpCall has on Op that is either pure or impure
        @assert is_pure(call) "Call is not pure ($call)."
        return new{G}(call)
    end

    function PureInstruction(call::G) where {G<:AbstractCall}
        # catch other AbstractCalls (GetIndex/SetIndex)
        return new{G}(call)
    end
end

"""
An `ImpureInstruction` is a call to an operation
that not always returns the same value
if the same arguments are passed to the instruction.
It is not suitable for memoization,
and the Value-Number optimization must be
disabled for this call.

An instruction is considered impure if its `call`
has an impure `Op`, or if onde of its arguments
is mutable.

Marking as mutable avoids Value-Numbering for this call.
"""
mutable struct ImpureInstruction{T} <: LinearInstruction{T}
    call::T

    function ImpureInstruction(call::G) where {G<:AbstractOpCall}
        # An AbstractOpCall has on Op that is either pure or impure
        @assert is_impure(call) "Call is not impure ($call)."
        return new{G}(call)
    end

    function ImpureInstruction(call::G) where {G<:AbstractCall}
        # catch other AbstractCalls (GetIndex/SetIndex)
        return new{G}(call)
    end

    function ImpureInstruction(call::Assignment)
        return new{Assignment}(call)
    end
end

struct CallUnary{OP, A<:AbstractValue} <: AbstractOpCall{OP}
    op::OP
    arg::A
end

struct CallBinary{OP, A<:AbstractValue, B<:AbstractValue} <: AbstractOpCall{OP}
    op::OP
    arg1::A
    arg2::B
end

struct Call3Args{OP, A<:AbstractValue, B<:AbstractValue, C<:AbstractValue} <: AbstractOpCall{OP}
    op::OP
    arg1::A
    arg2::B
    arg3::C
end

struct Call4Args{OP, A<:AbstractValue, B<:AbstractValue, C<:AbstractValue, D<:AbstractValue} <: AbstractOpCall{OP}
    op::OP
    arg1::A
    arg2::B
    arg3::C
    arg4::D
end

struct Call5Args{OP, A<:AbstractValue, B<:AbstractValue, C<:AbstractValue, D<:AbstractValue, E<:AbstractValue} <: AbstractOpCall{OP}
    op::OP
    arg1::A
    arg2::B
    arg3::C
    arg4::D
    arg5::E
end

struct Call6Args{OP, A<:AbstractValue, B<:AbstractValue, C<:AbstractValue, D<:AbstractValue, E<:AbstractValue, F<:AbstractValue} <: AbstractOpCall{OP}
    op::OP
    arg1::A
    arg2::B
    arg3::C
    arg4::D
    arg5::E
    arg6::F
end

abstract type Program end

mutable struct BasicBlock <: Program
    instructions::LookupTable{LinearInstruction}
    inputs::LookupTable{ImmutableVariable}
    mutable_locals::LookupTable{MutableVariable}
    immutable_locals::Dict{ImmutableVariable, ImmutableValue}
    outputs::LookupTable{Variable}
#    branch::Union{Nothing, BranchInstruction}
#    next::Union{Nothing, BasicBlock}
#    cfg::Union{Nothing, Program}
end

"""
A CompiledBasicBlock mirrors BasicBlock
but the instructions LookupTable is replaced
with Vector to save memory.
"""
struct CompiledBasicBlock <: Program
    instructions::Vector{LinearInstruction}
    inputs::LookupTable{ImmutableVariable}
    mutable_locals::LookupTable{MutableVariable}
    immutable_locals::Dict{ImmutableVariable, ImmutableValue}
    outputs::LookupTable{Variable}
end

#struct Goto <: BranchInstruction
#    target::BasicBlock
#end

#struct GotoIf{A<:StaticAddress} <: BranchInstruction
#    cond::A
#    target::BasicBlock
#end

#struct GotoIfNot{A<:StaticAddress} <: BranchInstruction
#    cond::A
#    target::BasicBlock
#end

#mutable struct CFG <: Program
#    start::BasicBlock
#    globals::Dict{Variable, StaticAddress}
#end

abstract type AbstractMachine end

"""
    compile(::Type{T}, program::Program) where {T<:AbstractMachine}

Compiles the IR to a Julia function.
It returns a function or a callable object (functor).

```julia
const OIR = OptimizingIR
bb = OIR.BasicBlock()
# (...) add instructions to basic block
finterpreter = OIR.compile(OIR.BasicBlockInterpreter, bb)
fnative = OIR.compile(OIR.Native, bb)
```
"""
function compile end


text:
#!/usr/bin/env julia

using Base.Test
using Drs.DrsMathProgSolverInterface

# Test case
A = [1 0; 0 2; 3 2]
b = [180, 150, 300]
c = [-3, -5]

lb = [-Inf, 150, 300]
ub = [180, 150, Inf]

mA = [1 0 1 0 0 0;
      0 2 0 0 1 0;
      3 2 0 -1 0 1]
mc = [-3, -5, 0, 0, 0, 0]

m = DrsMathProgModel(A, b, c)
Drs.DrsMathProgSolverInterface.DrsTransformToStandardForm!(m, lb, ub, :Min)

@test m.A == mA
@test m.b == b
@test m.c == mc

m = DrsMathProgModel(A, b, c)
Drs.DrsMathProgSolverInterface.DrsTransformToStandardForm!(m, lb, ub, :Max)

@test m.A == mA
@test m.b == b
@test m.c == -mc

# Test case
A = [1 0; 0 2; 3 2]
b = [180, 150, 300]
c = [-3, -5]

lb = [180, 150, 300]
ub = [Inf, Inf, Inf]

mA = [1 0 -1 0 0 1 0 0;
      0 2 0 -1 0 0 1 0;
      3 2 0 0 -1 0 0 1]
mc = [-3, -5, 0, 0, 0, 0, 0, 0]

m = DrsMathProgModel(A, b, c)
Drs.DrsMathProgSolverInterface.DrsTransformToStandardForm!(m, lb, ub, :Min)

@test m.A == mA
@test m.b == b
@test m.c == mc

m = DrsMathProgModel(A, b, c)
Drs.DrsMathProgSolverInterface.DrsTransformToStandardForm!(m, lb, ub, :Max)

@test m.A == mA
@test m.b == b
@test m.c == -mc

# Test case
A = [1 0; 0 2; 3 2]
b = [-180, -150, -300]
c = [-3, -5]

lb = [-Inf, -150, -300]
ub = [-180, -150, Inf]

mA = [-1 0 -1 0 0 1;
      0 -2 0 0 1 0;
      -3 -2 0 1 0 0]
mb = [180, 150, 300]
mc = [-3, -5, 0, 0, 0, 0]

m = DrsMathProgModel(A, b, c)
Drs.DrsMathProgSolverInterface.DrsTransformToStandardForm!(m, lb, ub, :Min)

@test m.A == mA
@test m.b == mb
@test m.c == mc

# Test case
A = [1 0; 0 2; 3 2]
b = [0, 0, 0]
c = [-3, -5]

lb = [-Inf, 0, 0]
ub = [0, 0, Inf]

mA = [1 0 1 0 0 0;
      0 2 0 0 1 0;
      3 2 0 -1 0 1]
mb = [0, 0, 0]
mc = [-3, -5, 0, 0, 0, 0]

m = DrsMathProgModel(A, b, c)
Drs.DrsMathProgSolverInterface.DrsTransformToStandardForm!(m, lb, ub, :Min)

@test m.A == mA
@test m.b == mb
@test m.c == mc

# Test case
A = [3 5 2; 4 4 4; 2 4 5]
b = [60, 72, 100]
c = [5, 10, 8]

lb = [60, 72, -Inf]
ub = [60, Inf, 100]

mA = [3 5 2 0 0 1 0; 4 4 4 -1 0 0 1; 2 4 5 0 1 0 0]
mb = [60, 72, 100]
mc = [5, 10, 8, 0, 0, 0, 0]

m = DrsMathProgModel(A, b, c)
Drs.DrsMathProgSolverInterface.DrsTransformToStandardForm!(m, lb, ub, :Min)

@test m.A == mA
@test m.b == mb
@test m.c == mc

m = DrsMathProgModel(A, b, c)
Drs.DrsMathProgSolverInterface.DrsTransformToStandardForm!(m, lb, ub, :Max)

@test m.A == mA
@test m.b == mb
@test m.c == -mc

max_stars_repo_name: abogushov/lplib

max_issues_repo_name: abogushov/lplib

max_forks_repo_name: abogushov/lplib


text:
module Test06

import Automa
import Automa.RegExp: @re_str
using Test

@testset "Test06" begin
    re = Automa.RegExp

    foo = re.cat("foo")
    foos = re.rep(re.cat(foo, re" *"))
    foo.actions[:exit]  = [:foo]
    actions = Dict(:foo => :(push!(ret, state.p:p-1); @escape))

    machine = Automa.compile(foos)

    @eval mutable struct MachineState
        p::Int
        cs::Int
        function MachineState()
            $(Automa.generate_init_code(Automa.CodeGenContext(), machine))
            return new(p, cs)
        end
    end

    for generator in (:table, :inline, :goto), checkbounds in (true, false), clean in (true, false)
        ctx = Automa.CodeGenContext(generator=generator, checkbounds=checkbounds, clean=clean)
        run! = @eval function (state, data)
            ret = []
            p = state.p
            cs = state.cs
            p_end = p_eof = lastindex(data)
            $(Automa.generate_exec_code(ctx, machine, actions))
            state.p = p
            state.cs = cs
            return ret
        end
        state = MachineState()
        data = b"foo foofoo   foo"
        @test run!(state, data) == [1:3]
        @test run!(state, data) == [5:7]
        @test run!(state, data) == [9:10]
        @test run!(state, data) == [12:16]
        @test run!(state, data) == []
        @test run!(state, data) == []
    end
end

end


text:
golfscorecalculator(parlist, scorelist) =
    sum( map( parscore -> parscore[2]-parscore[1], zip(parlist, scorelist)) )

max_stars_repo_name: UlrichBerntien/Codewars-Katas

max_issues_repo_name: UlrichBerntien/Codewars-Katas

max_forks_repo_name: UlrichBerntien/Codewars-Katas


text:
"""
    Scale <: Augmentor.AffineOperation

Description
--------------

Multiplies the image height and image width by the specified
`factors`. This means that the size of the output image depends
on the size of the input image.

The provided `factors` can either be numbers or vectors of
numbers.

- If numbers are provided, then the operation is deterministic
  and will always scale the input image with the same factors.

- In the case vectors are provided, then each time the operation
  is applied a valid index is sampled and the elements
  corresponding to that index are used as scaling factors.

The scaling is performed relative to the image center, which can
be useful when following the operation with [`CropNative`](@ref).

Usage
--------------

    Scale(factors)

    Scale(factors...)

Arguments
--------------

- **`factors`** : `NTuple` or `Vararg` of `Real` or
    `AbstractVector` that denote the scale factor(s) for each
    array dimension. If only one variable is specified it is
    assumed that height and width should be scaled by the same
    factor(s).

See also
--------------

[`Zoom`](@ref), [`Resize`](@ref), [`augment`](@ref)

Examples
--------------

```julia
using Augmentor
img = testpattern()

# half the image size
augment(img, Scale(0.5))

# uniformly scale by a random factor from 1.2, 1.3, or 1.4
augment(img, Scale([1.2, 1.3, 1.4]))

# scale by either 0.5x0.7 or by 0.6x0.8
augment(img, Scale([0.5, 0.6], [0.7, 0.8]))
```
"""
struct Scale{N,T<:AbstractVector} <: AffineOperation
    factors::NTuple{N,T}

    function Scale{N}(factors::NTuple{N,T}) where {N,T<:AbstractVector}
        eltype(T) <: Real || throw(ArgumentError("The specified factors must be vectors of Real. Actual: $T"))
        n = length(factors[1])
        n > 0 || throw(ArgumentError("The specified factors must all have a length greater than 0"))
        all(f->length(f)==n, factors) || throw(ArgumentError("The specified factors must all be of the same length for each dimension. Actual: $(factors)"))
        new{N,T}(factors)
    end
end
Scale{N}() where N = throw(MethodError(Scale{N}, ()))
Scale() = throw(MethodError(Scale, ()))
Scale{N}(::Tuple{}) where N = throw(MethodError(Scale{N}, ((), )))
Scale(::Tuple{}) = throw(MethodError(Scale, ((),)))
Scale(factors...) = Scale(factors)
Scale(factor::Union{AbstractVector,Real}) = Scale((factor, factor))
Scale(factors::NTuple{N,Any}) where {N} = Scale(map(vectorize, factors))
Scale(factors::NTuple{N,AbstractRange}) where {N} = Scale{N}(promote(factors...))
function Scale(factors::NTuple{N,AbstractVector}) where N
    Scale{N}(map(Vector{Float64}, factors))
end
function (::Type{Scale{N}})(factors::NTuple{N,Any}) where N
    Scale(map(vectorize, factors))
end

@inline supports_eager(::Type{<:Scale}) = false

randparam(op::Scale, imgs::Tuple) = randparam(op, imgs[1])

function randparam(op::Scale, img::AbstractArray{T,N}) where {T,N}
    i = safe_rand(1:length(op.factors[1]))
    ntuple(j -> Float64(op.factors[j][i]), Val(N))
end

function toaffinemap(op::Scale{2}, img::AbstractMatrix, idx)
    @inbounds tfm = recenter(@SMatrix([Float64(idx[1]) 0.; 0. Float64(idx[2])]), center(img))
    tfm
end

function showconstruction(io::IO, op::Scale)
    fct = length(op.factors[1]) == 1 ? map(first,op.factors) : op.factors
    print(io, typeof(op).name.name, '(', join(map(string, fct),", "), ')')
end

function Base.show(io::IO, op::Scale{N}) where N
    if get(io, :compact, false)
        str = join(map(t->join(t,"×"), collect(zip(op.factors...))), ", ")
        if length(op.factors[1]) == 1
            print(io, "Scale by $(str)")
        else
            print(io, "Scale by I ∈ {$(str)}")
        end
    else
        print(io, "Augmentor.")
        fct = length(op.factors[1]) == 1 ? map(first,op.factors) : op.factors
        print(io, typeof(op).name, "{$N}($(fct))")
    end
end


text:
using MPI
using CLIMA
using CLIMA.Mesh.Topologies
using CLIMA.Mesh.Grids
using CLIMA.MPIStateArrays
using CLIMA.DGBalanceLawDiscretizations
using Printf
using LinearAlgebra
using Logging

const ArrayTypes = (CLIMA.array_type(),)

const _nauxstate = 7
@inline function auxiliary_state_initialization!(aux, x, y, z, dim)
  @inbounds begin
    aux[1] = x
    aux[2] = y
    aux[3] = z
    if dim == 2
      aux[4] = x*y + z*y
      aux[5] = 2*x*y + sin(x)*y^2/2 - (z-1)^2*y^3/3
    else
      aux[4] = x*z + z^2/2
      aux[5] = 2*x*z + sin(x)*y*z - (1+(z-1)^3)*y^2/3
    end
  end
end
@inline function integral_knl(val, Q, aux)
  @inbounds begin
    x, y, z = aux[1], aux[2], aux[3]
    val[1] = x + z
    val[2] = 2*x + sin(x)*y - (z-1)^2*y^2
  end
end

using Test
function run(mpicomm, dim, ArrayType, Ne, N, FT)

  brickrange = ntuple(j->range(FT(0); length=Ne[j]+1, stop=3), dim)
  topl = StackedBrickTopology(mpicomm, brickrange,
                              periodicity=ntuple(j->true, dim))

  grid = DiscontinuousSpectralElementGrid(topl,
                                          FloatType = FT,
                                          DeviceArray = ArrayType,
                                          polynomialorder = N,
                                         )

  spacedisc = DGBalanceLaw(grid = grid,
                           length_state_vector = 0,
                           flux! = (x...) -> (),
                           numerical_flux! = (x...) -> (),
                           auxiliary_state_length = _nauxstate,
                           auxiliary_state_initialization! = (x...) ->
                           auxiliary_state_initialization!(x..., dim))

  Q = MPIStateArray(spacedisc)

  DGBalanceLawDiscretizations.indefinite_stack_integral!(spacedisc,
                                                         integral_knl, Q,
                                                         (6, 7))

  # Wrapping in Array ensure both GPU and CPU code use same approx
  @test Array(spacedisc.auxstate.data[:, 4, :]) ≈ Array(spacedisc.auxstate.data[:, 6, :])
  @test Array(spacedisc.auxstate.data[:, 5, :]) ≈ Array(spacedisc.auxstate.data[:, 7, :])
end

let
  CLIMA.init()
  mpicomm = MPI.COMM_WORLD
  ll = uppercase(get(ENV, "JULIA_LOG_LEVEL", "INFO"))
  loglevel = ll == "DEBUG" ? Logging.Debug :
  ll == "WARN"  ? Logging.Warn  :
  ll == "ERROR" ? Logging.Error : Logging.Info
  logger_stream = MPI.Comm_rank(mpicomm) == 0 ? stderr : devnull
  global_logger(ConsoleLogger(logger_stream, loglevel))

  numelem = (5, 5, 5)
  lvls = 1

  polynomialorder = 4

  @testset "$(@__FILE__)" for ArrayType in ArrayTypes
    for FT in (Float64,) #Float32)
      for dim = 2:3
        err = zeros(FT, lvls)
        for l = 1:lvls
          @info (ArrayType, FT, dim)
          run(mpicomm, dim, ArrayType, ntuple(j->2^(l-1) * numelem[j], dim),
              polynomialorder, FT)
        end
      end
    end
  end
end

nothing

max_stars_repo_name: emassoud/CLIMA

max_issues_repo_name: emassoud/CLIMA

max_forks_repo_name: emassoud/CLIMA


text:
### A Pluto.jl notebook ###
# v0.14.5

using Markdown
using InteractiveUtils

# ╔═╡ b99805ba-3931-11eb-0e18-6ffb5497630d
function ingredients(path::String)
    # this is from the Julia source code (evalfile in base/loading.jl)
    # but with the modification that it returns the module instead of the last object
    name = Symbol(basename(path))
    m = Module(name)
    Core.eval(m,
        Expr(:toplevel,
            :(eval(x) = $(Expr(:core, :eval))($name, x)),
            :(include(x) = $(Expr(:top, :include))($name, x)),
            :(include(mapexpr::Function, x) = $(Expr(:top, :include))(mapexpr, $name, x)),
            :(include($path))))
    m
end;

# ╔═╡ d840fcce-3931-11eb-0d08-59ed4b60b2f9
begin
    import Pkg
    #Pkg.activate("Mage")

    using GameOne
    using Random
    using Serialization
    using PlaymatSimulator

    Image = ImageActor
    GIF = GIFActor

    gd = "C:/Users/dusty/Documents/PlaymatProjects/PlaymatGames/Mage"
    US = deserialize("$gd/../../PlaymatSimulator/tmp/user_selection.jls")

    md"""
    ## EDH GAME STATE
    Provide an initial game state function for the game engine to load upon launch by defining a `gs` Dict{Symbol,Any}(...) object below with all of the game assets.
    """
end

# ╔═╡ ebd61216-4da3-11eb-2e77-5900d3980b8e
STAGE = game_include("$gd/Base/notebooks/game_stage.jl")

# ╔═╡ f1e8642e-4da3-11eb-26bc-234df79ed8a6
GS = game_include("$gd/Base/notebooks/game_settings.jl")

# ╔═╡ 2e97f2d2-4c1d-11eb-23d5-3fc9bcd4bc47
#const SCREEN_WIDTH = Int32(1080)

# ╔═╡ 93af0788-4c12-11eb-2c17-8ddd05dfc94d
#const SCREEN_BORDER = Int32(10)

# ╔═╡ 988ca2cc-4c12-11eb-284d-f51814cbf3a3
#const SCREEN_HEIGHT = Int32(1920)

# ╔═╡ 22163bfc-3932-11eb-37cc-0701ba61483e
gs = Dict{Symbol,Any}(
    :GAME_NAME => "EDH",
    :MOUSE_POS => Int32[0, 0],
    :MOUSE_OFFSETS => [Int32[0, 0]],
    :terminal => TextActor(">", "$gd/Base/fonts/OpenSans-Regular.ttf", alpha = 0),
    :music => readdir("$gd/Base/music"),
    :ui => OrderedDict(
        :cursor_icon => Image("mouse_cursor_icon",
            load("$gd/Base/ui/icons/RavenmoreIconPack/64/swordWood.png")),
        :cursor => Image("mouse_cursor",
            load("$gd/Base/ui/zones/area_wht.png"), alpha = 0, w = 1, h = 1),
    ),
    :stage => STAGE,
    :resources => Dict{Symbol,Any}(
        :mana => Dict{Symbol,Int32}(
            :white => 0,
            :blue => 0,
            :black => 0,
            :red => 0,
            :green => 0,
            :colorless => 0,
        ),
        :energy => 0,
        :life => 40,  # EDH life adjustment
    ),
    :resource_spinners => OrderedDict{Symbol,Actor}(
        :life => Text("40:L ",
            "$gd/Base/fonts/OpenSans-Semibold.ttf",
        ),
        :white_mana => Text(" 0:W ",
            "$gd/Base/fonts/OpenSans-Semibold.ttf", font_color = [255, 255, 255, 255]
        ),
        :blue_mana => Text(" 0:U ",
            "$gd/Base/fonts/OpenSans-Semibold.ttf", font_color = [0, 0, 255, 255]
        ),
        :black_mana => Text(" 0:B ",
            "$gd/Base/fonts/OpenSans-Semibold.ttf", font_color = [0, 0, 0, 255]
        ),
        :red_mana => Text(" 0:R ",
            "$gd/Base/fonts/OpenSans-Semibold.ttf", font_color = [255, 0, 0, 255]
        ),
        :green_mana => Text(" 0:G ",
            "$gd/Base/fonts/OpenSans-Semibold.ttf", font_color = [0, 255, 0, 255]
        ),
        :colorless_mana => Text(" 0:C ",
            "$gd/Base/fonts/OpenSans-Semibold.ttf", font_color = [125, 125, 125, 255]
        ),
        :energy => Text(" 0:E ",
            "$gd/Base/fonts/OpenSans-Semibold.ttf",
        ),
        :poison => Text(" 0:P ",
            "$gd/Base/fonts/OpenSans-Semibold.ttf", font_color = [255, 0, 255, 255]
        ),
    ),
    :sfx => OrderedDict(
        :shade_wht => Image("shade_wht",
            load("$gd/Base/ui/zones/area_wht.png")),
        :shade_blk => Image("shade_blk",
            load("$gd/Base/ui/zones/area_blk.png")),
        :sel_box => Image("sel_box",
            load("$gd/Base/ui/zones/area_wht.png"), alpha = 35, w = 0, h = 0),
    ),
    :overlay => OrderedDict(
        :texts => Actor[],
        :effects => Actor[],
        :shades => Actor[],
    ),
    :zone => OrderedDict(
        "Command" => Card[],  # EDH-specific zone
        "Graveyard" => Card[],
        "Battlefield" => Card[],
        "Library" => Card[],
        "Hand" => Card[],
        #"Exile" => Card[],  # replaced with Command
        #"Stack" => Union{Card,Spell}[]  # not implemented
    ),
    :group => OrderedDict(
        :clickables => Actor[],
        :selected => Actor[],
    ),
    :stack => OrderedDict(),
)

# ╔═╡ 4339c7e0-4c1f-11eb-1372-5b42b11fa88e
begin
    merge!(gs, GS)
    gs
end

# ╔═╡ Cell order:
# ╠═b99805ba-3931-11eb-0e18-6ffb5497630d
# ╟─d840fcce-3931-11eb-0d08-59ed4b60b2f9
# ╠═ebd61216-4da3-11eb-2e77-5900d3980b8e
# ╠═f1e8642e-4da3-11eb-26bc-234df79ed8a6
# ╠═2e97f2d2-4c1d-11eb-23d5-3fc9bcd4bc47
# ╠═93af0788-4c12-11eb-2c17-8ddd05dfc94d
# ╟─988ca2cc-4c12-11eb-284d-f51814cbf3a3
# ╠═22163bfc-3932-11eb-37cc-0701ba61483e
# ╠═4339c7e0-4c1f-11eb-1372-5b42b11fa88e


text:
function test_isa(go, term1, term2)
    @test is_a(go, term1, term2)
    @test !is_a(go, term2, term1)
end

let
    go = GOGraph("$testdir/data/go_mini.obo")
    term1 = Term("GO:0000001", "one")
    term2 = Term("GO:0000002", "two")
    term4 = Term("GO:0000004", "four")
    term5 = Term("GO:0000005", "five")

    test_isa(go, term1, term2)
    test_isa(go, term4, term2)
    test_isa(go, term5, term4)
    test_isa(go, term5, term2)

    @test !is_a(go, term1, term5)
    @test !is_a(go, term5, term1)

    @test parents(go, term1) == [term2]
    @test isempty(parents(go, term2))
    @test parents(go, term4) == [term2]
    @test parents(go, term5) == [term4]

    @test isempty(children(go, term1))
    @test children(go, term2) == [term1, term4]
    @test children(go, term4) == [term5]
    @test isempty(children(go, term5))
end


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

using Tables
using DataFrames

###############################################################################
# Config.

on_stock = true
# PCB, PSU...
fixed_cost = 0

###############################################################################

families = read_families()

#Read offers with params from XLS.
if !isfile("tmp/offers.xls")
	include("get_offers.jl")
end

offers = read_table("tmp/offers.xls", "offers")

###############################################################################

# Super table with all kind of params needed for check bellow.
uber_table = offers

# Hack for having some FPGAs in pinout zips.
a = families["Artix-7"]["dev_pack_banks"]
a[("XC7A12T", "CPG236")] = a[("XC7A15T", "CPG236")]
a[("XC7A12T", "CSG325")] = a[("XC7A15T", "CSG325")]
a[("XC7A25T", "CPG236")] = a[("XC7A35T", "CPG236")]
a[("XC7A25T", "CSG325")] = a[("XC7A35T", "CSG325")]

N_rows = size(uber_table)[1]

#TODO Should make metric for PCB price depending on package.
# Then little larger FPGAs will be better.
cost = uber_table[:price] .+ fixed_cost

###############################################################################

uber_table[:pins] = Vector{Float64}(N_rows)
for r in eachrow(uber_table)
	dpp_t = families[r[:family]]["dev_pack_pins"]
	# Make device-package pair.
	dp = (r[:device], r[:package])

	# Search for pins info with such device-package pair.
	i = findfirst(dpp_t[:dev_pack], dp)

	# No such device-package pair.
	if i == 0
		return 0
	end

	# Extract row for device-package pair, to harvest pins data.
	dpp_r = dpp_t[i, :]

	c = 0
	if haskey(dpp_r, :HR)
		c += dpp_r[:HR][1]
	end
	if haskey(dpp_r, :HP)
		c += dpp_r[:HP][1]
	end

	r[:pins] = c
end

uber_table[:cost_per_pin] = cost./uber_table[:pins]

###############################################################################

uber_table[:HR_DDR_MTps] = Vector{Int}(N_rows)
uber_table[:HP_DDR_MTps] = Vector{Int}(N_rows)
for r in eachrow(uber_table)
	ds = families[r[:family]]["ddr_speeds"]
	r[:HR_DDR_MTps] = ds["HR"][r[:speed_grade]]
	if haskey(ds, "HP")
		r[:HP_DDR_MTps] = ds["HP"][r[:speed_grade]]
	else
		r[:HP_DDR_MTps] = 0
	end
end

###############################################################################

pins_for_byte_group = 12 #TODO Could it be lower?

function DDR_bytes_groups_per_cols(pin_type, col)
	counts = Int[]
	for r in eachrow(uber_table)
		dpb = families[r[:family]]["dev_pack_banks"]
		# Make device-package pair.
		dp = (r[:device], r[:package])

		count = 0

		# No such device-package pair.
		if !haskey(dpb, dp)
			warn("No bank data for ", dp)
		else
			# Bank table.
			bt = dpb[dp]

			pin_type = string(pin_type)

			for r in eachrow(bt)
				bank = r[:bank]
				if bank != "NA" && length(bank) == 2 &&
					r[:pin_type] == pin_type &&	parse(Int, bank[1]) == col

					# All 4 byte groups needed.
					if r[:byte_group_0_pin_num] >= pins_for_byte_group &&
						r[:byte_group_1_pin_num] >= pins_for_byte_group &&
						r[:byte_group_2_pin_num] >= pins_for_byte_group &&
						r[:byte_group_3_pin_num] >= pins_for_byte_group
						count += 4
					end
				end
			end
		end

		push!(counts, count)
	end

	return counts
end
uber_table[:HR_DDR_BG_c1] = DDR_bytes_groups_per_cols(:HR, 1)
uber_table[:HR_DDR_BG_c2] = DDR_bytes_groups_per_cols(:HR, 2)
uber_table[:HR_DDR_BG_c3] = DDR_bytes_groups_per_cols(:HR, 3)
uber_table[:HR_DDR_BG_c4] = DDR_bytes_groups_per_cols(:HR, 4)
uber_table[:HP_DDR_BG_c1] = DDR_bytes_groups_per_cols(:HP, 1)
uber_table[:HP_DDR_BG_c2] = DDR_bytes_groups_per_cols(:HP, 2)
uber_table[:HP_DDR_BG_c3] = DDR_bytes_groups_per_cols(:HP, 3)
uber_table[:HP_DDR_BG_c4] = DDR_bytes_groups_per_cols(:HP, 4)
@assert sum(DDR_bytes_groups_per_cols(:HP, 5)) == 0

###############################################################################

function BG_to_B(BG)
	#TODO What if could make 6 bytes?
	if BG >= 11*4
		warn("Could set 4 DIMMs on one column?")
		return 8*4
	elseif BG >= 11*3
		return 8*3
	elseif BG >= 11*2
		return 8*2
	elseif BG >= 11
		return 8
	elseif BG >= 7
		return 4
	elseif BG >= 5
		return 2
	elseif BG >= 4
		return 1
	else
		return 0
	end
end

uber_table[:HR_DDR_B] = Vector{Int}(N_rows)
uber_table[:HP_DDR_B] = Vector{Int}(N_rows)
for r in eachrow(uber_table)
	r[:HR_DDR_B] =
		BG_to_B(r[:HR_DDR_BG_c1]) +
		BG_to_B(r[:HR_DDR_BG_c2]) +
		BG_to_B(r[:HR_DDR_BG_c3]) +
		BG_to_B(r[:HR_DDR_BG_c4])
	r[:HP_DDR_B] =
		BG_to_B(r[:HP_DDR_BG_c1]) +
		BG_to_B(r[:HP_DDR_BG_c2]) +
		BG_to_B(r[:HP_DDR_BG_c3]) +
		BG_to_B(r[:HP_DDR_BG_c4])
end
uber_table[:DDR_B] = uber_table[:HR_DDR_B] + uber_table[:HP_DDR_B]

###############################################################################

function BG_to_DIMM(BG)
	if BG >= 11*4
		warn("Could set 4 DIMMs on one column?")
		return 4
	elseif BG >= 11*3
		return 3
	elseif BG >= 11*2
		return 2
	elseif BG >= 11
		return 1
	else
		return 0
	end
end

N = size(uber_table)[1]
uber_table[:HR_DIMM] = Vector{Int}(N)
uber_table[:HP_DIMM] = Vector{Int}(N)
for r in eachrow(uber_table)
	r[:HR_DIMM] =
		BG_to_DIMM(r[:HR_DDR_BG_c1]) +
		BG_to_DIMM(r[:HR_DDR_BG_c2]) +
		BG_to_DIMM(r[:HR_DDR_BG_c3]) +
		BG_to_DIMM(r[:HR_DDR_BG_c4])
	r[:HP_DIMM] =
		BG_to_DIMM(r[:HP_DDR_BG_c1]) +
		BG_to_DIMM(r[:HP_DDR_BG_c2]) +
		BG_to_DIMM(r[:HP_DDR_BG_c3]) +
		BG_to_DIMM(r[:HP_DDR_BG_c4])
end
uber_table[:DIMM] = uber_table[:HR_DIMM] + uber_table[:HP_DIMM]

###############################################################################

uber_table[:cost_per_DDR_B] = cost./uber_table[:DDR_B]
uber_table[:cost_per_DIMM] = cost./uber_table[:DIMM]

###############################################################################

uber_table[:HR_DDR_MBps] = uber_table[:HR_DDR_B] .* uber_table[:HR_DDR_MTps]
uber_table[:HP_DDR_MBps] = uber_table[:HP_DDR_B] .* uber_table[:HP_DDR_MTps]
uber_table[:HR_DIMM_MBps] = 8*uber_table[:HR_DIMM] .* uber_table[:HR_DDR_MTps]
uber_table[:HP_DIMM_MBps] = 8*uber_table[:HP_DIMM] .* uber_table[:HP_DDR_MTps]

uber_table[:DDR_MBps] = uber_table[:HR_DDR_MBps] + uber_table[:HP_DDR_MBps]
uber_table[:DIMM_MBps] = uber_table[:HR_DIMM_MBps] + uber_table[:HP_DIMM_MBps]

###############################################################################

uber_table[:cost_per_DDR_MBps] = cost./uber_table[:DDR_MBps]
uber_table[:cost_per_DIMM_MBps] = cost./uber_table[:DIMM_MBps]

###############################################################################

macro save(table)
	quote
		n = $(string(table))
		write_table("tmp/" * n * ".xlsx", n, $table)
	end
end

# Save it just for documentation.
@save uber_table

###############################################################################

# Select those who have needed stock.
if on_stock
	uber_table = uber_table[uber_table[:stock_vs_need] .>= 0, :]
end

###############################################################################
# Cost.

price = deepcopy(uber_table)
sort!(price, cols = :price)

###############################################################################
# Cost per pin.

cost_per_pin = deepcopy(uber_table)
sort!(cost_per_pin, cols = :cost_per_pin)

if false
	# Cheapest from family:
	println("Cheapest Artix-7:")
	println(
		sort(
			cost_per_pin[cost_per_pin[:family] .== "Artix-7", :],
			cols = :price
		)[1, :]
	)
	println("Cheapest Kintex-7:")
	println(
		sort(
			cost_per_pin[cost_per_pin[:family] .== "Kintex-7", :],
			cols = :price
		)[1, :]
	)
	println("Cheapest Virtex-7:")
	println(
		sort(
			cost_per_pin[cost_per_pin[:family] .== "Virtex-7", :],
			cols = :price
		)[1, :]
	)

	@save cost_per_pin
	println("Cheapers per pin:")
	best = cost_per_pin[1, :]
	println(best)
end

###############################################################################
# Cost per memory width.

cost_per_DDR_B = deepcopy(uber_table)
sort!(cost_per_DDR_B, cols = :cost_per_DDR_B)

cost_per_DIMM = deepcopy(uber_table[uber_table[:DIMM] .!= 0, :])
sort!(cost_per_DIMM, cols = :cost_per_DIMM)

###############################################################################
# Cost per memory bandwidth.

cost_per_DDR_MBps = deepcopy(uber_table)
sort!(cost_per_DDR_MBps, cols = :cost_per_DDR_MBps)

cost_per_DIMM_MBps = deepcopy(uber_table[uber_table[:DIMM] .!= 0, :])
sort!(cost_per_DIMM_MBps, cols = :cost_per_DIMM_MBps)
@save cost_per_DIMM_MBps

best_for_DDR_is_best_for_DIMM =
	cost_per_DDR_MBps[1, :] == cost_per_DIMM_MBps[1, :]

@show best_for_DDR_is_best_for_DIMM

println("Best (cost per bandwidth):")
println(cost_per_DIMM_MBps[1, :])

###############################################################################

#TODO Check DSP performance and bandwidth with Wang ratio.

#TODO Check:
# + cost per pin
# + cost per pin bandwidth
# + cost per SO-DIMM badwidth (sharing common pins or not)
# - compare price with and without stock
# - cost with and without PCIe
# - without PCIe and user rest of SO-DIMM pins for parallel bus.

###############################################################################

println("End")

###############################################################################


text:
#= For any given array there are 2^N subsets. Thus generating each and every subset using BitMasking.
We follow the binary representation of whole numbers and use them as boolean value if a character in the Set
will occur or not. We will remove the duplicate subsets and display the non null unique subsets.
=#

## Function

function subsets(array, n)
    println("The subsets are: ")
    powerSet = Set()
    for i = 1:(1<<n)
        s = ""
        for j = 1:n
            if (i & (1 << j) > 0)
                temp = string(array[j])
                s *= temp
            end
        end
        if (length(s) > 0)
            println(s)
        end
    end
end

## Input

print("Enter the length of array: ")
n = readline()
n = parse(Int64, n)
print("Enter the elements of array: ")
array = Int64[]
for i = 1:n
    temp = readline()
    temp = parse(Int64, temp)
    push!(array, temp)
end

## Calling the Function

subsets(array, n)

#=
Sample Test Case:
Input:
    Enter the length of array: 5
    Enter the elements of array: 1 2 3 4 5
Output:
    The subsets are:
    124
    234
    2
    1234
    134
    34
    24
    4
    123
    3
    5
    23
    13
    14
    12
    1
Time complexity: O(2^N)
Space complexity: O(1)
=#

max_stars_repo_name: TechSpiritSS/NeoAlgo


text:
module mt3gen1
using StaticArrays
using MeshCore: nshapes
using MeshSteward: vtkwrite
using MeshSteward: T3blockx
using Test
function test()
    connectivity = T3blockx([0.0, 1.0, 3.0], [0.0, 1.0, 3.0], :a)
    @test (nshapes(connectivity.right), nshapes(connectivity.left)) == (9, 8)
    vtkwrite("mt3gen1a", connectivity)
    try rm("mt3gen1a.vtu"); catch end
    true
end
end
using .mt3gen1
mt3gen1.test()

module mt3gen1b
using StaticArrays
using MeshCore: nshapes
using MeshSteward: vtkwrite
using MeshSteward: T3blockx
using Test
function test()
    connectivity = T3blockx([0.0, 1.0, 3.0], [0.0, 1.0, 3.0], :b)
    @test (nshapes(connectivity.right), nshapes(connectivity.left)) == (9, 8)
    vtkwrite("mt3gen1b", connectivity)
    try rm("mt3gen1b.vtu"); catch end
    true
end
end
using .mt3gen1b
mt3gen1b.test()

module mt3gen2
using StaticArrays
using MeshCore: nshapes
using MeshSteward: vtkwrite
using MeshSteward: T3block
using Test
function test()
    connectivity = T3block(1.0, 2.0, 7, 9, :a)
    # @show (nshapes(connectivity.right), nshapes(connectivity.left))
    @test (nshapes(connectivity.right), nshapes(connectivity.left)) == (80, 63*2)
    vtkwrite("mt3gen2", connectivity)
    try rm("mt3gen2.vtu"); catch end
    true
end
end
using .mt3gen2
mt3gen2.test()

module mt3gen3
using StaticArrays
using MeshCore: nshapes
using MeshSteward: vtkwrite
using MeshSteward: export_MESH
using MeshSteward: T3blockx, gradedspace
using Test
function test()
    connectivity = T3blockx([1.0, 2.0, 3.0], gradedspace(0.0, 5.0, 7, 2), :a)
    # @show (nshapes(connectivity.right), nshapes(connectivity.left))
    @test (nshapes(connectivity.right), nshapes(connectivity.left)) == (21, 24)
    vtkwrite("mt3gen3", connectivity)
    try rm("mt3gen3.vtu"); catch end
    export_MESH("mt3gen3", connectivity)
    try rm("mt3gen3.mesh"); catch end
    try rm("mt3gen3-xyz.dat"); catch end
    try rm("mt3gen3-conn.dat"); catch end
    true
end
end
using .mt3gen3
mt3gen3.test()

module mt3gen5
using StaticArrays
using MeshCore: nshapes
using MeshSteward: vtkwrite
using MeshSteward: T3blockx
using Test
function test()
    connectivity = T3blockx([0.0, 1.0, 3.0], [0.0, 1.0, 3.0], :a)
    @test (nshapes(connectivity.right), nshapes(connectivity.left)) == (9, 8)
    vtkwrite("mt3gen5a", connectivity)
    try rm("mt3gen5a.vtu"); catch end
    true
end
end
using .mt3gen5
mt3gen5.test()

module mt6gen1
using StaticArrays
using MeshCore: nshapes
using MeshSteward: vtkwrite
using MeshSteward: T6blockx
using Test
function test()
    connectivity = T6blockx([0.0, 1.0, 3.0], [0.0, 1.0, 3.0], :b)
    @test (nshapes(connectivity.right), nshapes(connectivity.left)) == (25, 8)
    vtkwrite("mt6gen1", connectivity)
    try rm("mt6gen1.vtu"); catch end
    true
end
end
using .mt6gen1
mt6gen1.test()

module mt6gen2
using StaticArrays
using MeshCore: nshapes
using MeshSteward: vtkwrite
using MeshSteward: T6block
using Test
function test()
    connectivity = T6block(2.0, 0.75*pi, 6, 5, :b)
    locs = connectivity.right.attributes["geom"]
    for i in 1:length(locs)
        r, a = locs[i][1]+2.7, locs[i][2]
        locs[i] = (cos(a)*r, sin(a)*r)
    end
    @test nshapes(connectivity.left) == 60
    vtkwrite("mt6gen2", connectivity)
    try rm("mt6gen2.vtu"); catch end
    true
end
end
using .mt6gen2
mt6gen2.test()

module mt6gen3
using StaticArrays
using MeshCore: nshapes
using MeshSteward: vtkwrite
using MeshSteward: T3block, T6block, T3toT6
using Test
function test()
    connectivity = T3block(2.0, 0.75*pi, 6, 5, :b)
    locs = connectivity.right.attributes["geom"]
    for i in 1:length(locs)
        r, a = locs[i][1]+2.7, locs[i][2]
        locs[i] = (cos(a)*r, sin(a)*r)
    end
    @test nshapes(connectivity.left) == 60
    connectivity = T3toT6(connectivity)
    vtkwrite("mt6gen3", connectivity)
    try rm("mt6gen3.vtu"); catch end
    true
end
end
using .mt6gen3
mt6gen3.test()

module mt6gen4
using StaticArrays
using MeshCore: nshapes
using MeshSteward: vtkwrite
using MeshSteward: T3block, T6block, T3toT6, T6toT3
using Test
function test()
    connectivity = T6block(2.0, 0.75*pi, 6, 5, :b)
    locs = connectivity.right.attributes["geom"]
    for i in 1:length(locs)
        r, a = locs[i][1]+2.7, locs[i][2]
        locs[i] = (cos(a)*r, sin(a)*r)
    end
    @test nshapes(connectivity.left) == 60
    connectivity3 = T6toT3(connectivity)
    @test nshapes(connectivity3.left) == nshapes(connectivity.left)
    @test nshapes(connectivity3.right) != nshapes(connectivity.right)
    @test nshapes(connectivity3.right) == 42
    vtkwrite("mt6gen4", connectivity3)
    try rm("mt6gen4.vtu"); catch end
    true
end
end
using .mt6gen4
mt6gen4.test()


text:
using Pkg, HTTP

export downloadOpenMPIData, downloadCalibrationDataHighRes, downloadCalibrationDataLowRes,
       download2DData, download3DData, basedir

basedir() =  joinpath(dirname(pathof(OpenMPIData)),"..")

function makeFolders()
  mkpath(joinpath(basedir(),"data","reconstructions"))
  mkpath(joinpath(basedir(),"data","calibrations"))
  mkpath(joinpath(basedir(),"data","measurements","resolutionPhantom"))
  mkpath(joinpath(basedir(),"data","measurements","shapePhantom"))
  mkpath(joinpath(basedir(),"data","measurements","concentrationPhantom"))
end

function downloadOpenMPIData()
  download1DData()
  download2DData()
  download3DData()
  downloadCalibrationDataLowRes()
  #downloadCalibrationDataHighRes()  # This should be done explicitely
end

function download1DData()
  makeFolders()

  files = ["data/measurements/resolutionPhantom/1.mdf",
           "data/measurements/shapePhantom/1.mdf",
           "data/measurements/concentrationPhantom/1.mdf"]

  for file in files
    HTTP.open("GET", "http://media.tuhh.de/ibi/openMPIData/"*file) do http
      open(joinpath(basedir(),file), "w") do file
        write(file, http)
      end
    end
  end
end
function download2DData()
  makeFolders()

  files = ["data/measurements/resolutionPhantom/2.mdf",
           "data/measurements/shapePhantom/2.mdf",
           "data/measurements/concentrationPhantom/2.mdf"]

  for file in files
    HTTP.open("GET", "http://media.tuhh.de/ibi/openMPIData/"*file) do http
      open(joinpath(basedir(),file), "w") do file
        write(file, http)
      end
    end
  end
end

function download3DData()
  makeFolders()
  files = ["data/measurements/resolutionPhantom/3.mdf",
           "data/measurements/shapePhantom/3.mdf",
           "data/measurements/concentrationPhantom/3.mdf"]

  for file in files
    HTTP.open("GET", "http://media.tuhh.de/ibi/openMPIData/"*file) do http
      open(joinpath(basedir(),file), "w") do file
        write(file, http)
      end
    end
  end
end

function downloadCalibrationDataLowRes()
  makeFolders()

  files = ["data/calibrations/1.mdf",
           "data/calibrations/2.mdf",
           "data/calibrations/3.mdf",]

  for file in files
    HTTP.open("GET", "http://media.tuhh.de/ibi/openMPIData/"*file) do http
      open(joinpath(basedir(),file), "w") do file
        write(file, http)
      end
    end
  end
end

function downloadCalibrationDataHighRes()
  makeFolders()

  files = ["data/calibrations/4.mdf",
           "data/calibrations/5.mdf",
           "data/calibrations/6.mdf",]

  for file in files
    HTTP.open("GET", "http://media.tuhh.de/ibi/openMPIData/"*file) do http
      open(joinpath(basedir(),file), "w") do file
        write(file, http)
      end
    end
  end
end


text:
using Random
Random.seed!(1)

N = 10^5
prob1 = 0.7
eps0, eps1 = 0.1, 0.05

flipWithProb(bit,prob) = rand() < prob ? xor(bit,1) : bit

TxData = rand(N) .< prob1
RxData = [x == 0 ? flipWithProb(x,eps0) : flipWithProb(x,eps1) for x in TxData]

numTx1 = 0
totalRx1 = 0
for i in 1:N
   if RxData[i] == 1
        global totalRx1 += 1
        global numTx1 += TxData[i]
    end
end

monteCarlo = numTx1/totalRx1
analytic = ((1-eps1)*0.7)/((1-eps1)*0.7+0.3*eps0)

println("Monte Carlo: ", monteCarlo, "\t\tAnalytic: ", analytic)

max_stars_repo_name: Yoshinobu-Ishizaki/StatsWithJuliaBook

max_issues_repo_name: Yoshinobu-Ishizaki/StatsWithJuliaBook

max_forks_repo_name: Yoshinobu-Ishizaki/StatsWithJuliaBook


text:
using Documenter, Noise

using Images, TestImages, ImageIO, Random
 # set seed fixed for documentation
Random.seed!(42)

DocMeta.setdocmeta!(Noise, :DocTestSetup, :(using Noise, Images, TestImages, ImageIO); recursive=true)
makedocs(modules = [Noise],
         sitename = "Noise.jl",
         pages = ["index.md",
                  "man/additive_white_gaussian.md",
                  "man/mult_gauss.md",
                  "man/salt_pepper.md",
                  "man/poisson.md",
                  "man/quantization.md",
                  "man/function_references.md",
                 ]
        )

deploydocs(repo = "github.com/roflmaostc/Noise.jl.git")

max_issues_repo_name: roflmaostc/Noise

max_forks_repo_name: roflmaostc/Noise


text:
# Copyright 2018 IBM Corp.
#
# Licensed under the Apache License, Version 2.0 (the "License");
# you may not use this file except in compliance with the License.
# You may obtain a copy of the License at
#
#     http://www.apache.org/licenses/LICENSE-2.0
#
# Unless required by applicable law or agreed to in writing, software
# distributed under the License is distributed on an "AS IS" BASIS,
# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
# See the License for the specific language governing permissions and
# limitations under the License.

module TestSemanticEnrichment
using Test

using Catlab.WiringDiagrams
using SemanticFlowGraphs

const db_filename = joinpath(@__DIR__, "ontology", "data", "employee.json")
db = OntologyDB()
load_ontology_file(db, db_filename)

# Convenience methods to make raw boxes and ports with or without annotations.
const prefix = "opendisc/employee"
add_raw_box!(f::WiringDiagram, args...; kw...) =
  add_box!(f, raw_box(args...; kw...))
raw_box(inputs, outputs; kw...) =
  Box(RawNode(; kw...), raw_ports(inputs), raw_ports(outputs))
raw_box(name::String, inputs, outputs) =
  Box(RawNode(annotation="$prefix/$name"), raw_ports(inputs), raw_ports(outputs))
raw_ports(n::Int) = RawPort[ RawPort() for i in 1:n ]
raw_ports(xs::Vector) = RawPort[ raw_port(x) for x in xs ]
raw_port(::Nothing) = RawPort()
raw_port(name::String) = RawPort(annotation="$prefix/$name")
raw_port(name::String, index::Int) =
  RawPort(annotation="$prefix/$name", annotation_index=index)
raw_port(args::Tuple) = raw_port(args...)

# Expand single annotated box.
f = WiringDiagram(raw_ports(["employee"]), raw_ports(["employee"]))
v = add_raw_box!(f, "manager", [("employee",1)], [("employee",1)])
add_wires!(f, [
  (input_id(f),1) => (v,1),
  (v,1) => (output_id(f),1)
])
actual = to_semantic_graph(db, f)
target = WiringDiagram(concepts(db, ["employee"]), concepts(db, ["employee"]))
reports_to = add_box!(target, Box(concept(db, "reports-to")))
add_wires!(target, [
  (input_id(target), 1) => (reports_to, 1),
  (reports_to, 1) => (output_id(target), 1),
])
@test actual == target

# Expand slots.
f = WiringDiagram(raw_ports(["employee"]), raw_ports(["department","str","str"]))
manager = add_raw_box!(f, "manager", [("employee",1)], [("employee",1)])
dept = add_raw_box!(f, "employee-department", [("employee",1)], [("department",1)])
first_name = add_raw_box!(f, [("employee",1)], [("str",1)],
                          annotation="$prefix/employee",
                          annotation_index=1, annotation_kind=SlotAnnotation)
last_name = add_raw_box!(f, [("employee",1)], [("str",1)],
                         annotation="$prefix/employee",
                         annotation_index=2, annotation_kind=SlotAnnotation)
add_wires!(f, [
  (input_id(f), 1) => (manager, 1),
  (manager, 1) => (dept, 1),
  (manager, 1) => (first_name, 1),
  (manager, 1) => (last_name, 1),
  (dept, 1) => (output_id(f), 1),
  (first_name, 1) => (output_id(f), 2),
  (last_name, 1) => (output_id(f), 3),
])
actual = to_semantic_graph(db, f)
target = WiringDiagram(concepts(db, ["employee"]),
                       concepts(db, ["department", "string", "string"]))
reports_to = add_box!(target, Box(concept(db, "reports-to")))
works_in = add_box!(target, Box(concept(db, "works-in")))
first_name = add_box!(target, Box(concept(db, "person-first-name")))
last_name = add_box!(target, Box(concept(db, "person-last-name")))
add_wires!(target, [
  (input_id(target), 1) => (reports_to, 1),
  (reports_to, 1) => (works_in, 1),
  (reports_to, 1) => (first_name, 1),
  (reports_to, 1) => (last_name, 1),
  (works_in, 1) => (output_id(target), 1),
  (first_name, 1) => (output_id(target), 2),
  (last_name, 1) => (output_id(target), 3),
])
@test actual == target

# Collapse adjacent unannotated boxes.
f = WiringDiagram(raw_ports(1), raw_ports(1))
u = add_raw_box!(f, 1, 1)
v = add_raw_box!(f, 1, 1)
add_wires!(f, [
  (input_id(f),1) => (u,1),
  (u,1) => (v,1),
  (v,1) => (output_id(f),1)
])
actual = to_semantic_graph(db, f)
target = WiringDiagram([nothing], [nothing])
v = add_box!(target, Box(nothing, [nothing], [nothing]))
add_wires!(target, [
  (input_id(target), 1) => (v,1),
  (v,1) => (output_id(target), 1),
])
@test actual == target

# Don't collapse adjacent unannotated boxes with an intermediate annotated box.
f = WiringDiagram(raw_ports(0), raw_ports(0))
u = add_raw_box!(f, [], [nothing, "employee"])
v = add_raw_box!(f, [nothing, "employee"], [])
manager = add_raw_box!(f, "manager", [("employee",1)], [("employee",1)])
add_wires!(f, [
  (u,1) => (v,1),
  (u,2) => (manager,1),
  (manager,1) => (v,2),
])
actual = to_semantic_graph(db, f)
target = WiringDiagram([], [])
u = add_box!(target, Box(nothing, [], [nothing, concept(db,"employee")]))
v = add_box!(target, Box(nothing, [nothing, concept(db,"employee")], []))
reports_to = add_box!(target, Box(concept(db,"reports-to")))
add_wires!(target, [
  (u,1) => (v,1),
  (u,2) => (reports_to,1),
  (reports_to,1) => (v,2),
])
@test actual == target

# Don't assume that collapsibility is a transitive relation.
f = WiringDiagram(raw_ports(0), raw_ports(0))
u = add_raw_box!(f, [], [nothing, "employee"])
v = add_raw_box!(f, 1, 1)
w = add_raw_box!(f, [nothing, "employee"], [])
manager = add_raw_box!(f, "manager", [("employee",1)], [("employee",1)])
add_wires!(f, [
  (u,1) => (v,1),
  (v,1) => (w,1),
  (u,2) => (manager,1),
  (manager,1) => (w,2),
])
actual = to_semantic_graph(db, f)
target = WiringDiagram([], [])
u = add_box!(target, Box(nothing, [], [concept(db,"employee"), nothing]))
v = add_box!(target, Box(nothing, [nothing, concept(db,"employee")], []))
reports_to = add_box!(target, Box(concept(db,"reports-to")))
add_wires!(target, [
  (u,2) => (v,1),
  (u,1) => (reports_to,1),
  (reports_to,1) => (v,2),
])
@test actual == target

end

max_stars_repo_name: IBM/semanticflowgraph

max_issues_repo_name: epatters/semanticflowgraph

max_forks_repo_name: IBM/semanticflowgraph


text:
# # Poisson Equation
#
# We want to solve the [poisson equation](https://en.wikipedia.org/wiki/Poisson_equation) on the unit interval. It is given by
# `Δu = f` with boundary conditions `u(0) = a` and `u(1) = b`
# First of all let us choose some values for the parameters and remark, that there is an exact solution:
f = 1.0
a = -1.0
b = 2.0

u_analytic(x) = f/2*x^2 + (b-a-f/2) * x + a

# We would like to recompute this solution numerically
using DiffEqOperators

nknots = 10
h = 1.0/(nknots+1)
ord_deriv = 2
ord_approx = 2

Δ = CenteredDifference(ord_deriv, ord_approx, h, nknots)
bc = DirichletBC(a, b)

# Before solving the equation, lets take a look at Δ and bc:
# display(Array(Δ))
# display(bc*zeros(nknots))
# We see that `Δ` is a (lazy) matrix with the laplace stencil extended over the boundaries.
# And `bc` acts by padding the values just outside the boundaries.

u = (Δ*bc) \ fill(f, nknots)
knots = range(h, step=h, length=nknots)

# Since we used a second order approximation and the analytic solution itself was a second order
# polynomial, we expect that they are equal up to rounding errors:
using Test
@test u ≈ u_analytic.(knots)


text:
import Pkg
Pkg.add("Plots")
Pkg.add("PyPlot")
import Plots as plt
plt.pyplot()

include("./solver.jl")
import .solver as sol
import LinearAlgebra as lina

# function f(t::Float64, y::AbstractVector{Float64})::AbstractVector{Float64}
#     return y
# end

# a = 0.0
# b = 10.0
# h = 0.1
# y_euler = sol.euler([1.0], f, a, b, h)
# y_midpoint = sol.midpoint([1.0], f, a, b, h)

# k = ceil(Integer, (b - a) / h)
# t_range = range(a, stop=b, length=k)
# y_true = exp.(t_range)
# p = plt.plot(t_range, y_true, label="true")
# plt.plot!(t_range, lina.transpose(y_euler), label="Euler")
# plt.plot!(t_range, lina.transpose(y_midpoint), label="midpoint")
# plt.png(p, "solver/exp_plot.png")

function f_lorenz(t::Float64, y::AbstractVector{Float64})::AbstractVector{Float64}
    return [
   -y[2] - y[3]
   y[1] + 0.2 * y[2]
   0.2 + y[1] * y[3] - 5.7 * y[3]
 ]
end

a = 0.0
b = 200.0
h = 0.01
println("starting")
y_midpoint = sol.midpoint([-1.0; 0; 0], f_lorenz, a, b, h)
println("finished")
y_trans = lina.transpose(y_midpoint)
p = plt.plot(y_trans[:,1], y_trans[:,2], y_trans[:,3], label="")
plt.png(p, "solver/lorenz.png")

max_stars_repo_name: applied-math-coding/article-snippets

max_issues_repo_name: applied-math-coding/article-snippets

max_forks_repo_name: applied-math-coding/article-snippets


text:
# core formatting functions

### auxiliary functions

### print char n times

function _repprint(out::IO, c::AbstractChar, n::Int)
    while n > 0
        print(out, c)
        n -= 1
    end
end

### print string or char

function _pfmt_s(out::IO, fs::FormatSpec, s::Union{AbstractString,AbstractChar})
    wid = fs.width
    slen = length(s)
    if wid <= slen
        print(out, s)
    else
        a = fs.align
        if a == '<'
            print(out, s)
            _repprint(out, fs.fill, wid-slen)
        else
            _repprint(out, fs.fill, wid-slen)
            print(out, s)
        end
    end
end

### print integers

_mul(x::Integer, ::_Dec) = x * 10
_mul(x::Integer, ::_Bin) = x << 1
_mul(x::Integer, ::_Oct) = x << 3
_mul(x::Integer, ::Union{_Hex, _HEX}) = x << 4

_div(x::Integer, ::_Dec) = div(x, 10)
_div(x::Integer, ::_Bin) = x >> 1
_div(x::Integer, ::_Oct) = x >> 3
_div(x::Integer, ::Union{_Hex, _HEX}) = x >> 4

function _ndigits(x::Integer, op)  # suppose x is non-negative
    m = 1
    q = _div(x, op)
    while q > 0
        m += 1
        q = _div(q, op)
    end
    return m
end

_ipre(op) = ""
_ipre(::Union{_Hex, _HEX}) = "0x"
_ipre(::_Oct) = "0o"
_ipre(::_Bin) = "0b"

_digitchar(x::Integer, ::_Bin) = x == 0 ? '0' : '1'
_digitchar(x::Integer, ::_Dec) = Char(Int('0') + x)
_digitchar(x::Integer, ::_Oct) = Char(Int('0') + x)
_digitchar(x::Integer, ::_Hex) = Char(x < 10 ? '0' + x : 'a' + (x - 10))
_digitchar(x::Integer, ::_HEX) = Char(x < 10 ? '0' + x : 'A' + (x - 10))

_signchar(x::Real, s::AbstractChar) = signbit(x) ? '-' :
                                s == '+' ? '+' :
                                s == ' ' ? ' ' : '\0'

function _pfmt_int(out::IO, sch::AbstractChar, ip::ASCIIStr, zs::Integer, ax::Integer,
                   op::Op) where {Op}
    # print sign
    sch != '\0' && print(out, sch)
    # print prefix
    !isempty(ip) && print(out, ip)
    # print padding zeros
    zs > 0 && _repprint(out, '0', zs)
    # print actual digits
    ax == 0 ? print(out, '0') : _pfmt_intdigits(out, ax, op)
    nothing
end

function _pfmt_intdigits(out::IO, ax::T, op::Op) where {Op, T<:Integer}
    b_lb = _div(ax, op)
    b = one(T)
    while b <= b_lb
        b = _mul(b, op)
    end
    r = ax
    while b > 0
        (q, r) = divrem(r, b)
        print(out, _digitchar(q, op))
        b = _div(b, op)
    end
end

function _pfmt_i(out::IO, fs::FormatSpec, x::Integer, op::Op) where {Op}
    # calculate actual length
    ax = abs(x)
    xlen = _ndigits(abs(x), op)
    # sign char
    sch = _signchar(x, fs.sign)
    if sch != '\0'
        xlen += 1
    end
    # prefix (e.g. 0x, 0b, 0o)
    ip = ""
    if fs.ipre
        ip = _ipre(op)
        xlen += length(ip)
    end

    # printing
    wid = fs.width
    if wid <= xlen
        _pfmt_int(out, sch, ip, 0, ax, op)
    elseif fs.zpad
        _pfmt_int(out, sch, ip, wid-xlen, ax, op)
    else
        a = fs.align
        if a == '<'
            _pfmt_int(out, sch, ip, 0, ax, op)
            _repprint(out, fs.fill, wid-xlen)
        else
            _repprint(out, fs.fill, wid-xlen)
            _pfmt_int(out, sch, ip, 0, ax, op)
        end
    end
end

### print floating point numbers

function _pfmt_float(out::IO, sch::AbstractChar, zs::Integer, intv::Real, decv::Real, prec::Int)
    # print sign
    sch != '\0' && print(out, sch)

    # print padding zeros
    zs > 0 && _repprint(out, '0', zs)

    idecv = round(Integer, decv * exp10(prec))
    # print integer part
    if intv == 0
        print(out, '0')
    else
        _pfmt_intdigits(out, intv, _Dec())
    end
    # print decimal part
    if prec > 0
        print(out, '.')
        nd = _ndigits(idecv, _Dec())
        nd < prec && _repprint(out, '0', prec - nd)
        _pfmt_intdigits(out, idecv, _Dec())
    end
end

function _pfmt_f(out::IO, fs::FormatSpec, x::AbstractFloat)
    # separate sign, integer, and decimal part
    rax = round(abs(x); digits = fs.prec)
    sch = _signchar(x, fs.sign)
    intv = trunc(Integer, rax)
    decv = rax - intv

    # calculate length
    xlen = _ndigits(intv, _Dec()) + ifelse(fs.prec > 0, fs.prec + 1, 0)
    sch != '\0' && (xlen += 1)

    # print
    wid = fs.width
    if wid <= xlen
        _pfmt_float(out, sch, 0, intv, decv, fs.prec)
    elseif fs.zpad
        _pfmt_float(out, sch, wid-xlen, intv, decv, fs.prec)
    else
        a = fs.align
        if a == '<'
            _pfmt_float(out, sch, 0, intv, decv, fs.prec)
            _repprint(out, fs.fill, wid-xlen)
        else
            _repprint(out, fs.fill, wid-xlen)
            _pfmt_float(out, sch, 0, intv, decv, fs.prec)
        end
    end
end

function _pfmt_floate(out::IO, sch::AbstractChar, zs::Integer, u::Real, prec::Int, e::Integer,
                      ec::AbstractChar)
    intv = trunc(Integer, u)
    decv = u - intv
    if intv == 0 && decv != 0
        intv = 1
        decv -= 1
    end
    _pfmt_float(out, sch, zs, intv, decv, prec)
    print(out, ec)
    if e >= 0
        print(out, '+')
    else
        print(out, '-')
        e = -e
    end
    e < 10 && print(out, '0')
    _pfmt_intdigits(out, e, _Dec())
end

function _pfmt_e(out::IO, fs::FormatSpec, x::AbstractFloat)
    # extract sign, significand, and exponent
    ax = abs(x)
    sch = _signchar(x, fs.sign)
    if ax == 0.0
        e = 0
        u = zero(x)
    else
        rax = round(ax; sigdigits = fs.prec + 1)
        e = floor(Integer, log10(rax))  # exponent
        u = rax * exp10(-e)  # significand
        i = 0
        v10 = 1
        while isinf(u)
            i += 1
            i > 18 && (u = 0.0; e = 0; break)
            v10 *= 10
            u = v10 * rax * exp(-e - i)
        end
    end

    # calculate length
    xlen = 6 + fs.prec
    abs(e) > 99 && (xlen += _ndigits(abs(e), _Dec()) - 2)
    sch != '\0' && (xlen += 1)

    # print
    ec = isuppercase(fs.typ) ? 'E' : 'e'
    wid = fs.width
    if wid <= xlen
        _pfmt_floate(out, sch, 0, u, fs.prec, e, ec)
    elseif fs.zpad
        _pfmt_floate(out, sch, wid-xlen, u, fs.prec, e, ec)
    else
        a = fs.align
        if a == '<'
            _pfmt_floate(out, sch, 0, u, fs.prec, e, ec)
            _repprint(out, fs.fill, wid-xlen)
        else
            _repprint(out, fs.fill, wid-xlen)
            _pfmt_floate(out, sch, 0, u, fs.prec, e, ec)
        end
    end
end

function _pfmt_g(out::IO, fs::FormatSpec, x::AbstractFloat)
    # number decomposition
    ax = abs(x)
    if 1.0e-4 <= ax < 1.0e6
        _pfmt_f(out, fs, x)
    else
        _pfmt_e(out, fs, x)
    end
end

function _pfmt_specialf(out::IO, fs::FormatSpec, x::AbstractFloat)
    if isinf(x)
        _pfmt_s(out, fs, x > 0 ? "Inf" : "-Inf")
    else
        @assert isnan(x)
        _pfmt_s(out, fs, "NaN")
    end
end


text:
"""
MT19937 (32 bits)
"""

module MersenneTwister

export MtParam, MtRandom, mt19937_param, extract_number!

# Further reading (why Mersenne Twister sucks):
#
# https://cs.stackexchange.com/questions/50059/why-is-the-mersenne-twister-regarded-as-good

# Warning: Bad design. The type of parameters should not be fixed so we can easily extend it to the
# 64-bit version.
"Parameters for Mersenne Twister"
struct MtParam
    "Word size"
    w::UInt32  # This field makes no sense

    "Degree of recurrence"
    n::UInt32

    "Middle word"
    m::UInt32

    "Separation point of one word"
    r::UInt32

    "Coefficients of the rational normal form twist matrix"
    a::UInt32

    "TGFSR(R) tempering bitmasks `b`"
    b::UInt32

    "TGFSR(R) tempering bitmasks `c`"
    c::UInt32

    "TGFSR(R) tempering bit shifts `s`"
    s::UInt32

    "TGFSR(R) tempering bit shifts `t`"
    t::UInt32

    "Additional MT tempering bit shifts/masks `u`"
    u::UInt32

    "Additional MT tempering bit shifts/masks `d`"
    d::UInt32

    "Additional MT tempering bit shifts/masks `l`"
    l::UInt32

    "(Not part of the algorithm)"
    f::UInt32
end

function mt19937_param()::MtParam
    (w, n, m, r) = (UInt32(32), UInt32(624), UInt32(397), UInt32(31))
    a = 0x9908_B0DF
    (u, d) = (UInt32(11), 0xFFFF_FFFF)
    (s, b) = (UInt32(7), 0x9D2C_5680)
    (t, c) = (UInt32(15), 0xEFC6_0000)
    l = UInt32(18)
    f = 0x6c07_8965

    MtParam(
        w, n, m, r,
        a, b, c,
        s, t,
        u, d, l,
        f,
    )
end

mutable struct MtRandom
    param::MtParam
    state::Vector{UInt32}
    index::UInt32
    lower_mask::UInt32
    upper_mask::UInt32

    function MtRandom(seed::UInt32 = UInt32(5489), param::MtParam = mt19937_param())
        index = UInt32(0)  # Rust has got `NonNull`, but inferior PLs don't.
        lower_mask = (UInt32(1) << param.r) - UInt32(1)
        upper_mask = ~lower_mask
        @assert UInt32((1 << param.w) - 1) == upper_mask | lower_mask

        state = zeros(UInt32, param.n)
        state[1] = seed
        for i in UInt32(2) : param.n
            state[i] = param.f * (state[i - 1] ⊻ (state[i - 1] >> (param.w - 2))) + (i - UInt32(1))
        end

        new(param, state, index, lower_mask, upper_mask)
    end
end

function extract_number!(mt::MtRandom)::UInt32
    if iszero(mt.index)
        twist!(mt)
    end

    y = mt.state[mt.index]
    y ⊻= (y >> mt.param.u) & mt.param.d
    y ⊻= (y << mt.param.s) & mt.param.b
    y ⊻= (y << mt.param.t) & mt.param.c
    y ⊻= (y >> mt.param.l)

    mt.index = (mt.index + UInt32(1)) % (mt.param.n + UInt32(1))
    y
end

function twist!(mt::MtRandom)
    @inbounds for i in 1:mt.param.n
        x = (mt.state[i] & mt.upper_mask) |
            (mt.state[mod1(i + 1, mt.param.n)] & mt.lower_mask)
        xA = x >> 1
        if isodd(x)
            xA ⊻= mt.param.a
        end
        mt.state[i] = mt.state[mod1(i + mt.param.m, mt.param.n)] ⊻ xA
    end
    mt.index = UInt32(1)
end

end  # module

max_stars_repo_name: gyk/TrivialSolutions

max_issues_repo_name: gyk/TrivialSolutions

max_forks_repo_name: gyk/TrivialSolutions