text:
"""
    solve_with_ipopt(problem, robot;
                     initial_guess=Float64[],
                     use_inv_dyn=false,
                     minimise_τ=false,
                     user_options=Dict())

Solve the nonlinear optimization problem with Ipopt.

Further options can be set using the keyword arguments. See [Solver Interfaces](@ref).

# Keyword arguments
- `initial_guess::Array{Float64}=Float64[]`: the starting point for the solver.
- `use_inv_dyn::Bool=false`: if true, enables the use of inverse dynamics instead of forward dynamics.
- `minimise_τ::Bool=false`: if true, activates a cost function to minimize the joint torques.
- `user_options::Dict=Dict()`: the user options for Ipopt.

See also: [`solve_with_knitro`](@ref)
"""
function solve_with_ipopt(problem::Problem, robot::Robot;
                          initial_guess::Array{Float64}=Float64[],
                          use_inv_dyn::Bool=false,
                          minimise_τ::Bool=false,
                          user_options::Dict=Dict())
    # # # # # # # # # # # # # # # #
    # Variables and their bounds  #
    # # # # # # # # # # # # # # # #

    n1 = robot.n_q *  problem.num_knots       # positions
    n2 = robot.n_v *  problem.num_knots       # velocities
    n3 = robot.n_τ * (problem.num_knots - 1)  # torques

    n = n1 + n2 + n3  # total number of decision variables

    x_L = Array{Float64,2}(undef, robot.n_q, 0)
    x_U = Array{Float64,2}(undef, robot.n_q, 0)

    for k = 1:problem.num_knots
        if k ∈ keys(problem.fixed_q)
            x_L = [x_L problem.fixed_q[k]]
            x_U = [x_U problem.fixed_q[k]]
        else
            x_L = [x_L robot.q_lo]
            x_U = [x_U robot.q_hi]
        end

        if k ∈ keys(problem.fixed_v)
            x_L = [x_L problem.fixed_v[k]]
            x_U = [x_U problem.fixed_v[k]]
        else
            x_L = [x_L robot.v_lo]
            x_U = [x_U robot.v_hi]
        end

        if k <= problem.num_knots - 1
            if k ∈ keys(problem.fixed_τ)
                x_L = [x_L problem.fixed_τ[k]]
                x_U = [x_U problem.fixed_τ[k]]
            else
                x_L = [x_L robot.τ_lo]
                x_U = [x_U robot.τ_hi]
            end
        end
    end

    # # # # # # # #
    # Constraints #
    # # # # # # # #

    use_m₁ = true  # nonlinear equality: dynamics
    use_m₂ = true  # nonlinear equality: ee xyz-position

    m₁ = !use_m₁ ? 0 : (robot.n_q + robot.n_v) * (problem.num_knots - 1)  # equations of motion
    m₂ = !use_m₂ ? 0 : 3 * length(problem.ee_pos)  # xyz-position (3 NL equality)

    m = m₁ + m₂  # total number of constraints

    ind_con_dyn    = (1:m₁)
    ind_con_ee_pos = (1:m₂) .+ (m₁)

    g_L, g_U = Float64[], Float64[]

    if use_m₁
        append!(g_L, zeros(m₁))
        append!(g_U, zeros(m₁))
    end

    if use_m₂
        knots_con_ee = sort(collect(keys(problem.ee_pos)))
        con_ee = hcat([[problem.ee_pos[k].data...] for k in knots_con_ee]...)
        append!(g_L, vec(con_ee))
        append!(g_U, vec(con_ee))
    end

    @assert length(g_L) == length(g_U)
    @assert eltype(g_L) == eltype(g_U) == Float64

    # dimension of each mesh point
    nₓ = robot.n_q + robot.n_v + robot.n_τ

    if use_inv_dyn
        jacdata_dyn = problem.jacdata_inv_dyn
        dynamics_defects! = inverse_dynamics_defects!
    else
        jacdata_dyn = problem.jacdata_fwd_dyn
        dynamics_defects! = forward_dynamics_defects!
    end

    function eval_g(x, g)
        if use_m₁
            length_c = size(jacdata_dyn.jac, 1)  # length of constraint per pair of knots
            for i = 0:problem.num_knots - 2
                # Calculate the indices of the appropriate decision variables
                ind_vars = range(1 + i * nₓ, length=size(jacdata_dyn.jac, 2))

                # Evaluate constraints
                offset_con = i * length_c
                ind_cons = (1:length_c) .+ offset_con
                @views dynamics_defects!(g[ind_cons], robot, x[ind_vars], problem.dt)
            end
        end

        if use_m₂
            length_c = 3
            for (i, k) = enumerate(sort(collect(keys(problem.ee_pos))))
                ind_vars = range(1 + (k - 1) * nₓ, length=robot.n_q)  # indices of the decision variables
                offset_con = (i - 1) * length_c + (m₁)
                ind_cons = (1:length_c) .+ offset_con
                @views ee_position!(g[ind_cons], robot, x[ind_vars])  # constraint evaluation at that point
            end
        end
    end

    jacIndexConsCB = Cint[]
    jacIndexVarsCB = Cint[]

    if use_m₁
        for i = 0:problem.num_knots - 2
            offset_con = i * size(jacdata_dyn.jac, 1)
            ind_cons = rowvals(jacdata_dyn.jac) .+ offset_con
            append!(jacIndexConsCB, ind_cons)

            for j = 1:size(jacdata_dyn.jac, 2)
                ind_vars = fill(j + i * nₓ, length(nzrange(jacdata_dyn.jac, j)))
                append!(jacIndexVarsCB, ind_vars)
            end
        end

        @assert length(jacIndexConsCB) == length(jacIndexVarsCB)
    end

    if use_m₂
        ind_offset = m₁
        for k = knots_con_ee
            inds = rowvals(problem.jacdata_ee_position.jac) .+ ind_offset
            ind_offset += size(problem.jacdata_ee_position.jac, 1)
            append!(jacIndexConsCB, inds)
        end

        for k = knots_con_ee
            vals = vcat([fill(j + (k - 1) * nₓ, length(nzrange(problem.jacdata_ee_position.jac, j)))
                         for j = 1:size(problem.jacdata_ee_position.jac, 2)]...)
            append!(jacIndexVarsCB, vals)
        end

        @assert length(jacIndexConsCB) == length(jacIndexVarsCB)
    end

    function eval_jac_g(x, mode, rows, cols, values)
        if mode == :Structure
            for (i, r, c) in zip(1:length(jacIndexConsCB), jacIndexConsCB, jacIndexVarsCB)
                rows[i] = r
                cols[i] = c
            end
        else
            offset_prev = 0

            if use_m₁
                # Dynamics constraints
                for i = 0:problem.num_knots - 2
                    # Calculate the indices of the appropriate decision variables
                    ind_vars = range(1 + i * nₓ, length=size(jacdata_dyn.jac, 2))

                    # Evaluate the Jacobian at that point
                    jacdata_dyn(x[ind_vars])

                    # Pass the Jacobian to Knitro
                    offset_jac = i * jacdata_dyn.length_jac
                    ind_jac = (1:jacdata_dyn.length_jac) .+ offset_jac
                    values[ind_jac] = nonzeros(jacdata_dyn.jac)
                end
                offset_prev += (problem.num_knots - 1) * jacdata_dyn.length_jac
            end

            if use_m₂
                # End-effector constraints
                for (i, k) = enumerate(sort(collect(keys(problem.ee_pos))))
                    ind_qᵢ = range(1 + (k - 1) * nₓ, length=robot.n_q)  # indices of the decision variables
                    problem.jacdata_ee_position(x[ind_qᵢ])              # jacobian evaluation at that point
                    offset_jac = (i - 1) * problem.jacdata_ee_position.length_jac + offset_prev
                    ind_jac = (1:problem.jacdata_ee_position.length_jac) .+ offset_jac
                    values[ind_jac] = nonzeros(problem.jacdata_ee_position.jac)  # pass jacobian to Knitro
                end
            end
        end
    end

    # # # # # # #
    # Objective #
    # # # # # # #

    function eval_f(x)
        return 0
    end

    function eval_grad_f(x, grad_f)
        grad_f[:] = zeros(n)
    end

    if minimise_τ
        ind_τ = hcat([range(1 + (i * nₓ) + robot.n_q + robot.n_v, length=robot.n_τ)
                      for i = (1:problem.num_knots - 1) .- 1]...)
        indexVars, coefs = vec(ind_τ), fill(1 / (problem.num_knots - 1), n3)
        @assert length(indexVars) == length(coefs)

        eval_f = function(x)
            # Minimize τ, i.e., the necessary joint torques.
            return sum(coefs .* x[indexVars] .* x[indexVars])
        end

        eval_grad_f = function(x, grad_f)
            grad_f[:] = zeros(n)
            grad_f[indexVars] .= coefs .* 2x[indexVars]
        end
    end

    # # # # # # # # # #
    # Create Problem  #
    # # # # # # # # # #

    # number of nonzeros in the Jacobian of the constraints
    nele_jac = 0
    nele_jac += !use_m₁ ? 0 : jacdata_dyn.length_jac * (problem.num_knots - 1)
    nele_jac += !use_m₂ ? 0 : problem.jacdata_ee_position.length_jac * length(problem.ee_pos)

    prob = Ipopt.createProblem(n, vec(x_L), vec(x_U), m, g_L, g_U, nele_jac, 0,
                               eval_f, eval_g, eval_grad_f, eval_jac_g)

    # # # # # # # # #
    # Initial guess #
    # # # # # # # # #

    if initial_guess |> isempty
        initial_guess = zeros(n)
    end

    @assert length(initial_guess) == n

    # Set starting solution
    prob.x = copy(initial_guess)

    # # # # # # # # #
    # User Options  #
    # # # # # # # # #

    addOption(prob, "hessian_approximation", "limited-memory")
    addOption(prob, "mu_strategy", "adaptive")
    # addOption(prob, "linear_solver", "ma57")
    # addOption(prob, "ma57_pre_alloc", 2.0)
    addOption(prob, "tol", 1.0e-3)
    addOption(prob, "max_cpu_time", 10.0)

    foreach(x -> addOption(prob, x...), user_options)

    solver_log = SolverLog(n)

    function intermediate(alg_mod::Int, iter_count::Int, obj_value::Float64, inf_pr::Float64, inf_du::Float64, mu::Float64,
                          d_norm::Float64, regularization_size::Float64, alpha_du::Float64, alpha_pr::Float64, ls_trials::Int)
        print("")  # Flush the output to the Jupyter notebook cell

        update!(solver_log, abs_feas_error=inf_pr, obj_value=obj_value)

        return true
    end

    setIntermediateCallback(prob, intermediate)

    # # Perform a derivative check.
    # addOption(prob, "derivative_test", "first-order")

    # # # # #
    # Solve #
    # # # # #

    cpu_time = @elapsed status = Ipopt.solveProblem(prob)

    # println(Ipopt.ApplicationReturnStatus[status])
    # println(prob.x)
    # println(prob.obj_val)

    return cpu_time, prob.x, solver_log
end

export solve_with_ipopt


text:
using Base.Test

# Tutorials
for f in [
        "passthrough",
        "global_hypothesis_verification",
        "voxel_grid",
        "correspondence_grouping",
        "planar_segmentation",
        "statistical_removal",
        "region_growing_rgb_segmentation",
        "extract_indices",
        "tilt_compensation",
        "cluster_extraction",
        ]
    @testset "$f" begin
        include(joinpath(string(f, ".jl")))
    end
end

if parse(Int, get(ENV, "PCLJL_RUN_VISUALIZATION_TESTS", "1")) != 0
    @testset "offscreen_rendering" begin
        include("offscreen_rendering.jl")
    end
end

ext: jl


text:
# This file is auto-generated by AWSMetadata.jl
using AWS
using AWS.AWSServices: apigatewaymanagementapi

using AWS.Compat
using AWS.UUIDs
"""
    DeleteConnection()

Delete the connection with the provided id.

# Required Parameters
- `connectionId`:

"""

delete_connection(connectionId; aws_config::AWSConfig=global_aws_config()) = apigatewaymanagementapi("DELETE", "/@connections/$(connectionId)"; aws_config=aws_config)
delete_connection(connectionId, args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = apigatewaymanagementapi("DELETE", "/@connections/$(connectionId)", args; aws_config=aws_config)

"""
    GetConnection()

Get information about the connection with the provided id.

# Required Parameters
- `connectionId`:

"""

get_connection(connectionId; aws_config::AWSConfig=global_aws_config()) = apigatewaymanagementapi("GET", "/@connections/$(connectionId)"; aws_config=aws_config)
get_connection(connectionId, args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = apigatewaymanagementapi("GET", "/@connections/$(connectionId)", args; aws_config=aws_config)

"""
    PostToConnection()

Sends the provided data to the specified connection.

# Required Parameters
- `Data`: The data to be sent to the client specified by its connection id.
- `connectionId`: The identifier of the connection that a specific client is using.

"""

post_to_connection(Data, connectionId; aws_config::AWSConfig=global_aws_config()) = apigatewaymanagementapi("POST", "/@connections/$(connectionId)", Dict{String, Any}("Data"=>Data); aws_config=aws_config)
post_to_connection(Data, connectionId, args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = apigatewaymanagementapi("POST", "/@connections/$(connectionId)", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("Data"=>Data), args)); aws_config=aws_config)


text:
"""
    ZAdjoin

  Let α be an element satisfying a monic polynomial of degree n > 1 over the integers, Z. That is,

  ```
  α^n = p[1] + p[2]α + ... + p[n]α^(n-1),  (I)
  ```

  where the p[i]ϵZ. Let Z[α] be the ring formed by adjoining α to the Z. This module provides code for basic representation and manipulation of such rings.

  Elements of Z[α] can be represented as polynomials, q[1] + q[2]α + ... + q[n]α^(n-1), of degree n-1 in α, since higher powers of α can be reduced modulo the defining relationship, α^n=p(α). Thus, under addition, Z[α] is a module of dimension n over Z. Since, under multiplication, elements of Z[α] induce linear transformations of this module, elements may also be represented as matrices over Z. Specifically, α may be represented as the companion matrix of the defining polynomial, α^n-p(α), powers of α up to n-1 as corresponding powers of this companion matrix, and general elements as integer combinations of these powers.

  A defining relationship, essentially n integer coefficients p[1],...,p[n] as in equation (I) above, is represented by one of the types Modulus, KummerModulus, or QIModulus (QI for Quadratic Integer,) all concrete subtypes of AbstractModulus. Ring elements are represented as integer vectors and a modulus (instance of a concrete subtype of AbstractModulus) combined in the Element type. Addition, subtraction, multiplication, and exponentiation by integer powers of Elements are supported. Conversion of Elements to matrices is also supported.

  # Examples

  ```
  # A Kummer Ring of degree 4 is the Gaussian (complex) Integers
  julia> g = KummerModulus(4)
  ZAdjoin.KummerModulus{Int64}
  α^2 = -1

  # Function generator(modulus) will return the adjoined element.
  julia> α = generator(g)

  # Other elements can be formed from the adjoined element.
  julia> β = 3+4α # Another element
  ZAdjoin.Element{Int64}
  3+4α

  julia> α + β
  ZAdjoin.Element{Int64}
  3+5α

  julia> α*β
  ZAdjoin.Element{Int64}
  -4+3α

  # Compare with complex multiplication
  julia> im*(3+4im)
  -4 + 3im

  julia> α^2
  ZAdjoin.Element{Int64}
  -1

  # Matrix representation
  julia> Matrix(β)
  2×2 Array{Int64,2}:
   3  -4
   4   3

  # Representations are arithmetcally compatible (isomorphic)
  julia> Matrix(α+β^5)
  2×2 Array{Int64,2}:
    -237  3115
   -3115  -237

  julia> Matrix(α)+Matrix(β)^5
  2×2 Array{Int64,2}:
    -237  3115
   -3115  -237

  # The determinant of an element is the determinant of its Matrix
  # and acts as a multiplicative homomorphism into the integers.
  julia> det(β)
  25
  julia> det(1+α)
  2
  julia> det((1+α)*β)
  50

  ```

"""
module ZAdjoin

import Base: +, -, *, one, det, ==, string, show, Matrix
import Primes
export AbstractModulus, Modulus, KummerModulus, Element, generator

include("types.jl")
include("operations.jl")
include("utilities.jl")

end

max_stars_repo_name: WilCrofter/ZAdjoin

max_issues_repo_name: WilCrofter/ZAdjoin

max_forks_repo_name: WilCrofter/ZAdjoin


text:
using BenchmarkTools, GeometryTypes
import SearchBall

@benchmark SearchBall.rotate2d_left90([1, 1]) # realistic example

a = rand(1_000_000)
b = rand(1_000_000)
@benchmark SearchBall.vector(a, b) # maximum example
@benchmark SearchBall.vector([1.0, 2.0], [2.3, 4.5]) # realistic example

@benchmark SearchBall.intersection(LineSegment(Point2(2.0, 2.0), Point2(3.0, 3.0)), LineSegment(Point2(3.0, 2.0), Point2(2.0, 3.0))) == (true, Point2(2.5, 2.5)) # realistic example

@benchmark SearchBall.intersection(LineSegment(Point2(2.0, 2.0), Point2(3.0, 3.0)), HyperRectangle{2,Float64}([0.0, 0.0], [7.0, 4.0])) == (true, Point2(4.0, 4.0)) # realistic example

@benchmark SearchBall.does_intersect(LineSegment(Point2(2.0, 2.0), Point2(3.0, 3.0)), Circle(Point2(0.0, 0.0), 1.0)) # realistic example

@benchmark SearchBall.polygon_intersection([Point2(0.0, 0.0), Point2(6.0, 0.0), Point2(6.0, 2.0), Point2(3.0, 2.0), Point2(0.0, 2.0)], [Point2(1.0, 1.0), Point2(3.0, 3.0), Point2(5.0, 1.0), Point2(3.0, 4.0)]) # realistic example


text:
using Revise
	using DiffEqOperators, Setfield, Parameters
	using PseudoArcLengthContinuation, LinearAlgebra, Plots, SparseArrays
	const PALC = PseudoArcLengthContinuation

heatmapsol(x) = heatmap(reshape(Array(x), Nx, Ny)', color=:viridis)

Nx = 151*1
	Ny = 100*1
	lx = 8pi
	ly = 2*2pi/sqrt(3)

function Laplacian2D(Nx, Ny, lx, ly, bc = :Neumann)
	hx = 2lx/Nx
	hy = 2ly/Ny
	D2x = CenteredDifference(2, 2, hx, Nx)
	D2y = CenteredDifference(2, 2, hy, Ny)
	if bc == :Dirichlet
		Qx = Dirichlet0BC(typeof(hx))
		Qy = Dirichlet0BC(typeof(hy))
	elseif bc == :Neumann
		Qx = Neumann0BC(hx)
		Qy = Neumann0BC(hy)
	elseif bc == :Periodic
		Qx = PeriodicBC(hx)
		Qy = PeriodicBC(hy)
	end
	# @show norm(D2x - Circulant(D2x[1,:]))
	A = kron(sparse(I, Ny, Ny), sparse(D2x * Qx)[1]) + kron(sparse(D2y * Qy)[1], sparse(I, Nx, Nx))
	return A, D2x
end

function F_sh(u, p)
	@unpack l, ν, L1 = p
	return -L1 * u .+ (l .* u .+ ν .* u.^2 .- u.^3)
end

function dF_sh(u, p)
	@unpack l, ν, L1 = p
	return -L1 .+ spdiagm(0 => l .+ 2 .* ν .* u .- 3 .* u.^2)
end

X = -lx .+ 2lx/(Nx) * collect(0:Nx-1)
Y = -ly .+ 2ly/(Ny) * collect(0:Ny-1)

sol0 = [(cos(x) .+ cos(x/2) * cos(sqrt(3) * y/2) ) for x in X, y in Y]
		sol0 .= sol0 .- minimum(vec(sol0))
		sol0 ./= maximum(vec(sol0))
		sol0 = sol0 .- 0.25
		sol0 .*= 1.7
		heatmap(sol0', color=:viridis)

Δ, D2x = Laplacian2D(Nx, Ny, lx, ly, :Neumann)
L1 = (I + Δ)^2
par = (l = -0.1, ν = 1.3, L1 = L1)

optnew = PALC.NewtonPar(verbose = true, tol = 1e-8, maxIter = 20)
# optnew = PALC.NewtonPar(verbose = true, tol = 1e-8, maxIter = 20, eigsolver = EigArpack(0.5, :LM))
	sol_hexa, hist, flag = @time PALC.newton(
				x ->  F_sh(x, par),
				u -> dF_sh(u, par),
				vec(sol0),
				optnew)
	println("--> norm(sol) = ", norm(sol_hexa, Inf64))
	heatmapsol(sol_hexa)

heatmapsol(0.2vec(sol_hexa) .* vec([exp(-(x+0lx)^2/25) for x in X, y in Y]))

###################################################################################################
# recherche de solutions
deflationOp = DeflationOperator(2.0, (x, y)->dot(x, y), 1.0, [sol_hexa])

optnew = @set optnew.maxIter = 250
outdef, _, flag, _ = @time PALC.newton(
				x -> F_sh(x, par),
				u -> dF_sh(u, par),
				# 0.4vec(sol_hexa) .* vec([exp(-1(x+1lx)^2/25) for x in X, y in Y]),
				0.4vec(sol_hexa) .* vec([1 .- exp(-1(x+0lx)^2/55) for x in X, y in Y]),

				optnew, deflationOp, normN = x -> norm(x, Inf))
	println("--> norm(sol) = ", norm(outdef))
	heatmapsol(outdef) |> display
	flag && push!(deflationOp, outdef)

heatmapsol(deflationOp[4])

heatmapsol(0.4vec(sol_hexa) .* vec([1 .- exp(-1(x+0lx)^2/55) for x in X, y in Y]))
###################################################################################################
normC = x -> norm(x, Inf)
optcont = ContinuationPar(dsmin = 0.0001, dsmax = 0.005, ds= -0.001, pMax = -0.095, pMin = -1.0, newtonOptions = setproperties(optnew; tol = 1e-9, maxIter = 15), maxSteps = 145, detectBifurcation = 2, nev = 40, detectFold = false, dsminBisection =1e-7, saveSolEveryNsteps = 4)
	optcont = @set optcont.newtonOptions.eigsolver = EigArpack(0.1, :LM)

	br, u1 = @time PALC.continuation(
		(x, p) ->  F_sh(x, @set par.l = p),
		(x, p) -> dF_sh(x, @set par.l = p),
		deflationOp[1], -0.1, optcont;
		plot = true, verbosity = 3,
		tangentAlgo = BorderedPred(),
		# linearAlgo = MatrixBLS(),
		plotSolution = (x; kwargs...) -> (heatmap!(X, Y, reshape(x, Nx, Ny)'; color=:viridis, label="", kwargs...);ylims!(-1,1,subplot=4);xlims!(-.5,.3,subplot=4)),
		printSolution = (x,p) -> norm(x),
		finaliseSolution = (z, tau, step, contResult) -> (Base.display(contResult.eig[end].eigenvals) ;true) ,
		normC = x -> norm(x, Inf))
###################################################################################################
using IncompleteLU
prec = ilu(L1 + I,τ=0.05)
prec = lu(L1 + I)
ls = GMRESIterativeSolvers(tol = 1e-5, N = Nx*Ny, Pl = prec)

sol_hexa, _, flag = @time PALC.newton(
			x ->  F_sh(x, par),
			u -> dF_sh(u, par),
			vec(sol0),
			@set optnew.linsolver = ls)
	println("--> norm(sol) = ", norm(sol_hexa, Inf64))
	heatmapsol(sol_hexa)


text:
#   Copyright (c) 2020-2022 Velexi Corporation
#
#   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.

"""
Unit tests for isapprox(x, y) and !isapprox(x, y)
"""

# --- Imports

# Standard library
import InteractiveUtils.subtypes
using Test

# GeometricAlgebra.jl
using GeometricAlgebra

# --- File inclusions

# Test utilities
include("test_utils.jl")

# --- Tests

# ------ B::Blade

@testset "isapprox(B::Blade, C::Blade)" begin
    # Preparations
    vectors = [3 3; 4 4; 0 1]
    B = Blade(vectors)

    # dim(B) == dim(C), grade(B) == grade(C), volume(B) ≈ volume(C)
    # basis(B) ≈ basis(C)
    for precision_type1 in subtypes(AbstractFloat)
        for precision_type2 in subtypes(AbstractFloat)
            B = Blade(convert(Array{precision_type1}, vectors))
            C = Blade(convert(Array{precision_type2}, vectors))
            @test B ≈ C
        end
    end

    # dim(B) != dim(C)
    for precision_type1 in subtypes(AbstractFloat)
        vectors2 = [3 3; 4 4; 0 1; 0 1]
        for precision_type2 in subtypes(AbstractFloat)
            B = Blade(convert(Array{precision_type1}, vectors))
            C = Blade(convert(Array{precision_type2}, vectors2))
            @test B ≉ C
        end
    end

    # grade(B) != grade(C)
    for precision_type1 in subtypes(AbstractFloat)
        vectors2 = [3 3 3; 4 4 4; 0 1 0; 0 0 1]
        for precision_type2 in subtypes(AbstractFloat)
            B = Blade(convert(Array{precision_type1}, vectors))
            C = Blade(convert(Array{precision_type2}, vectors2))
            @test B ≉ C
        end
    end

    # volume(B) ≉ volume(C)
    for precision_type1 in subtypes(AbstractFloat)
        for precision_type2 in subtypes(AbstractFloat)
            B = Blade(convert(Array{precision_type1}, vectors))
            C = Blade(convert(Array{precision_type2}, vectors),
                      volume=2*volume(B))
            @test B ≉ C
        end
    end

    # basis(B) ≉ basis(C)
    test_volume = 5.0
    for precision_type1 in subtypes(AbstractFloat)
        vectors2 = [3 4; 4 10; 0 1]
        for precision_type2 in subtypes(AbstractFloat)
            B = Blade(convert(Array{precision_type1}, vectors),
                      volume=test_volume)
            C = Blade(convert(Array{precision_type2}, vectors2),
                      volume=test_volume)
            @test B ≉ C
        end
    end

    # B and C have opposite orientations
    test_volume = 5.0
    for precision_type1 in subtypes(AbstractFloat)
        for precision_type2 in subtypes(AbstractFloat)
            B = Blade(convert(Array{precision_type1}, vectors),
                      volume=test_volume)
            C = Blade(convert(Array{precision_type2}, vectors),
                      volume=-test_volume)
            @test B ≉ C
        end
    end
end

@testset "!isapprox(B::Blade, C::Pseudoscalar)" begin
    test_vectors = [3 3; 4 4; 0 1]

    test_dim = 3
    test_value = 5

    for precision_type in subtypes(AbstractFloat)
        B = Blade{precision_type}(test_vectors)
        C = Pseudoscalar{precision_type}(test_dim, test_value)
        @test B ≉ C
    end
end

@testset "!isapprox(B::Blade, C::Scalar)" begin
    test_vectors = [3 3; 4 4; 0 1]
    test_value = 5

    for precision_type in subtypes(AbstractFloat)
        B = Blade{precision_type}(test_vectors)
        C = Scalar{precision_type}(test_value)
        @test B ≉ C
    end
end

@testset "!isapprox(B::Blade, C::One)" begin
    test_vectors = [3 3; 4 4; 0 1]
    for precision_type in subtypes(AbstractFloat)
        B = Blade{precision_type}(test_vectors)
        C = One{precision_type}()
        @test B ≉ C
    end
end

@testset "!isapprox(B::Blade, C::Zero)" begin
    test_vectors = [3 3; 4 4; 0 1]
    for precision_type in subtypes(AbstractFloat)
        B = Blade{precision_type}(test_vectors)
        C = Zero{precision_type}()
        @test B ≉ C
    end
end

@testset "!isapprox(B::Blade, C::Real)" begin
    test_vectors = [3 3; 4 4; 0 1]
    test_value = 5
    for precision_type in subtypes(AbstractFloat)
        B = Blade{precision_type}(test_vectors)
        C = precision_type(test_value)
        @test B ≉ C
    end
end

@testset "isapprox(B::Blade, C::Vector)" begin
    # --- Preparations

    test_vector = [3; 4; 1]

    # B ≉ C
    test_vectors = [3 3; 4 4; 0 1]
    for precision_type in subtypes(AbstractFloat)
        B = Blade{precision_type}(test_vectors)
        C = Vector{precision_type}(test_vector)
        @test B ≉ C
    end

    test_vectors = [3; 4; 10]
    for precision_type in subtypes(AbstractFloat)
        B = Blade{precision_type}(test_vectors)
        C = Vector{precision_type}(test_vector)
        @test B ≉ C
    end

    # B ≈ C
    for precision_type1 in subtypes(AbstractFloat)
        for precision_type2 in subtypes(AbstractFloat)
            B = Blade{precision_type1}(test_vector)
            C = Vector{precision_type2}(test_vector)
            @test B ≈ C
        end
    end
end

# ------ B::Pseudoscalar

@testset "!isapprox(B::Pseudoscalar, C::Blade)" begin
    test_dim = 3
    test_value = 5

    test_vectors = [3 3; 4 4; 0 1]

    for precision_type in subtypes(AbstractFloat)
        B = Pseudoscalar{precision_type}(test_dim, test_value)
        C = Blade{precision_type}(test_vectors)
        @test B ≉ C
    end
end

@testset "isapprox(B::Pseudoscalar, C::Pseudoscalar)" begin
    # Preparations
    dim = 10

    test_value = get_random_value(1)  # add 1 to keep value away from 0

    # dim(B) == dim(C), value(B) == value(C)
    for precision_type1 in subtypes(AbstractFloat)
        for precision_type2 in subtypes(AbstractFloat)
            B = Pseudoscalar(dim, precision_type1(test_value))
            C = Pseudoscalar(dim, precision_type2(test_value))
            @test B ≈ C
        end
    end

    # dim(B) != dim(C), value(B) == value(C)
    for precision_type1 in subtypes(AbstractFloat)
        for precision_type2 in subtypes(AbstractFloat)
            B = Pseudoscalar(dim, precision_type1(test_value))
            C = Pseudoscalar(dim + 1, precision_type2(test_value))
            @test B ≉ C
        end
    end

    # dim(B) == dim(C), value(B) != value(C)
    for precision_type1 in subtypes(AbstractFloat)
        for precision_type2 in subtypes(AbstractFloat)
            B = Pseudoscalar(dim, precision_type1(test_value))
            C = Pseudoscalar(dim, precision_type2(test_value) + 1)
            @test B ≉ C
        end
    end
end

@testset "isapprox(B::Pseudoscalar, C::Scalar)" begin
    test_dim = 10
    test_value = 5

    for precision_type in subtypes(AbstractFloat)
        B = Pseudoscalar{precision_type}(test_dim, test_value)
        C = Scalar{precision_type}(test_value)
        @test B ≉ C
    end
end

@testset "isapprox(B::Pseudoscalar, C::One)" begin
    test_dim = 10
    test_value = 5

    for precision_type in subtypes(AbstractFloat)
        B = Pseudoscalar{precision_type}(test_dim, test_value)
        C = One{precision_type}()
        @test B ≉ C
    end
end

@testset "isapprox(B::Pseudoscalar, C::Zero)" begin
    test_dim = 10
    test_value = 5

    for precision_type in subtypes(AbstractFloat)
        B = Pseudoscalar{precision_type}(test_dim, test_value)
        C = Zero{precision_type}()
        @test B ≉ C
    end
end

@testset "isapprox(B::Pseudoscalar, C::Real)" begin
    test_dim = 10
    test_value = 5

    for precision_type in subtypes(AbstractFloat)
        B = Pseudoscalar{precision_type}(test_dim, test_value)
        C = precision_type(test_value)
        @test B ≉ C
    end
end

@testset "isapprox(B::Pseudoscalar, C::Vector)" begin
    test_value = 5

    test_dim_1 = 3
    test_vector_1 = [3; 4; 1]

    test_dim_2 = 1
    test_vector_2 = [test_value]

    for precision_type1 in subtypes(AbstractFloat)
        # B ≉ C
        B = Pseudoscalar{precision_type1}(test_dim_1, test_value)
        C = Vector{precision_type1}(test_vector_1)
        @test B ≉ C

        # B ≈ C
        B = Pseudoscalar{precision_type1}(test_dim_2, test_value)
        for precision_type2 in subtypes(AbstractFloat)
            C = Vector{precision_type2}(test_vector_2)
            @test B ≈ C
        end
    end
end

# ------ B::Scalar

@testset "!isapprox(B::Scalar, C::Blade)" begin
    test_value = 5
    test_vectors = [3 3; 4 4; 0 1]
    for precision_type in subtypes(AbstractFloat)
        B = Scalar{precision_type}(test_value)
        C = Blade{precision_type}(test_vectors)
        @test B ≉ C
    end
end

@testset "isapprox(B::Scalar, C::Pseudoscalar)" begin
    test_value = 5
    test_dim = 10
    for precision_type in subtypes(AbstractFloat)
        B = Scalar{precision_type}(test_value)
        C = Pseudoscalar{precision_type}(test_dim, test_value)
        @test B ≉ C
    end
end

@testset "isapprox(B::Scalar, C::Scalar)" begin
    # Preparations

    test_value = get_random_value(2)  # add 2 to keep value away from 0 and 1

    # B::Scalar, C::Scalar
    for precision_type1 in subtypes(AbstractFloat)
        for precision_type2 in subtypes(AbstractFloat)
            # B ≈ C
            B = Scalar(precision_type1(test_value))
            C = Scalar(precision_type2(test_value))
            @test B ≈ C

            # B ≉ C
            B = Scalar(precision_type1(test_value))
            C = Scalar(precision_type2(test_value + 1))
            @test B ≉ C
        end
    end
end

@testset "isapprox(B::Scalar, C::One)" begin
    for precision_type in subtypes(AbstractFloat)
        C = One{precision_type}()

        # B ≈ C
        test_value = 1 - (√eps(precision_type)) / 2
        B = Scalar{precision_type}(test_value)
        @test B ≈ C

        # B ≉ C
        test_value = 5
        B = Scalar{precision_type}(test_value)
        @test B ≉ C
    end
end

@testset "isapprox(B::Scalar, C::Zero)" begin
    test_value = 5
    for precision_type in subtypes(AbstractFloat)
        B = Scalar{precision_type}(test_value)
        C = Zero{precision_type}()
        @test B ≉ C

        # B ≈ C is not being tested as isapprox between a Scalar and a
        # Zero can never return true with the default atol and rtol value
    end
end

@testset "isapprox(B::Scalar, C::Real)" begin
    test_value = 5
    for precision_type1 in subtypes(AbstractFloat)
        B = Scalar{precision_type1}(test_value)

        # B ≉ C
        C = precision_type1(test_value + 1)
        @test B ≉ C

        # B ≈ C
        for precision_type2 in subtypes(AbstractFloat)
            C = precision_type2(test_value)
            @test B ≈ C
        end
    end
end

@testset "isapprox(B::Scalar, C::Vector)" begin
    test_value = 5
    test_vector = [3; 4; 1]
    for precision_type in subtypes(AbstractFloat)
        B = Scalar{precision_type}(test_value)
        C = Vector{precision_type}(test_vector)
        @test B ≉ C
    end
end

# ------ B::One

@testset "!isapprox(B::One, C::Blade)" begin
    test_vectors = [3 3; 4 4; 0 1]
    for precision_type in subtypes(AbstractFloat)
        B = One{precision_type}()
        C = Blade{precision_type}(test_vectors)
        @test B ≉ C
    end
end

@testset "isapprox(B::One, C::Pseudoscalar)" begin
    test_dim = 10
    test_value = 5
    for precision_type in subtypes(AbstractFloat)
        B = One{precision_type}()
        C = Pseudoscalar{precision_type}(test_dim, test_value)
        @test B ≉ C
    end
end

@testset "isapprox(B::One, C::Scalar)" begin
    for precision_type in subtypes(AbstractFloat)
        B = One{precision_type}()

        # B ≈ C
        test_value = 1 - (√eps(precision_type)) / 2
        C = Scalar{precision_type}(test_value)
        @test B ≈ C

        # B ≉ C
        test_value = 5
        C = Scalar{precision_type}(test_value)
        @test B ≉ C
    end
end

@testset "isapprox(B::One, C::One)" begin
    for precision_type1 in subtypes(AbstractFloat)
        for precision_type2 in subtypes(AbstractFloat)
            B = One{precision_type1}()
            C = One{precision_type2}()
            @test B ≈ C
        end
    end
end

@testset "isapprox(B::One, C::Zero)" begin
    for precision_type in subtypes(AbstractFloat)
        B = One{precision_type}()
        C = Zero{precision_type}()
        @test B ≉ C
    end
end

@testset "isapprox(B::One, C::Real)" begin
    test_value = 5
    for precision_type1 in subtypes(AbstractFloat)
        B = One{precision_type1}()

        # B ≉ C
        C = precision_type1(test_value)
        @test B ≉ C

        # B ≈ C
        for precision_type2 in subtypes(AbstractFloat)
            C = precision_type2(1)
            @test B ≈ C
        end
    end
end

@testset "isapprox(B::One, C::Vector)" begin
    test_vector = [3; 4; 1]
    for precision_type in subtypes(AbstractFloat)
        B = One{precision_type}()
        C = Vector{precision_type}(test_vector)
        @test B ≉ C
    end
end

# ------ B::Zero

@testset "!isapprox(B::Zero, C::Blade)" begin
    test_vectors = [3 3; 4 4; 0 1]
    for precision_type in subtypes(AbstractFloat)
        B = Zero{precision_type}()
        C = Blade{precision_type}(test_vectors)
        @test B ≉ C
    end
end

@testset "isapprox(B::Zero, C::Pseudoscalar)" begin
    test_dim = 10
    test_value = 5
    for precision_type in subtypes(AbstractFloat)
        B = Zero{precision_type}()
        C = Pseudoscalar{precision_type}(test_dim, test_value)
        @test B ≉ C
    end
end

@testset "isapprox(B::Zero, C::Scalar)" begin
    test_value = 5
    for precision_type in subtypes(AbstractFloat)
        B = Zero{precision_type}()
        C = Scalar{precision_type}(test_value)
        @test B ≉ C

        # B ≈ C is not being tested as isapprox between a Zero and a
        # Scalar can never return true with the default atol and rtol value
    end
end

@testset "isapprox(B::Zero, C::One)" begin
    for precision_type in subtypes(AbstractFloat)
        B = Zero{precision_type}()
        C = One{precision_type}()
        @test B ≉ C
    end
end

@testset "isapprox(B::Zero, C::Zero)" begin
    for precision_type1 in subtypes(AbstractFloat)
        for precision_type2 in subtypes(AbstractFloat)
            B = Zero{precision_type1}()
            C = Zero{precision_type2}()
            @test B ≈ C
        end
    end
end

@testset "isapprox(B::Zero, C::Real)" begin
    test_value = 5
    for precision_type1 in subtypes(AbstractFloat)
        B = Zero{precision_type1}()

        # B ≉ C
        C = precision_type1(test_value)
        @test B ≉ C

        # B ≈ C
        for precision_type2 in subtypes(AbstractFloat)
            C = precision_type2(0)
            @test B ≈ C
        end
    end
end

@testset "isapprox(B::Zero, C::Vector)" begin
    # B ≉ C
    test_vector = [3; 4; 1]
    for precision_type in subtypes(AbstractFloat)
        B = Zero{precision_type}()
        C = Vector{precision_type}(test_vector)
        @test B ≉ C
    end

    test_vector = [0.001; 0.002; 0.003]
    for precision_type in subtypes(AbstractFloat)
        B = Zero{precision_type}()
        C = Vector{precision_type}(test_vector)
        @test !isapprox(B, C, atol=0.001)
    end

    # B ≈ C
    test_vector = [0; 0; 0]
    for precision_type in subtypes(AbstractFloat)
        B = Zero{precision_type}()
        C = Vector{precision_type}(test_vector)
        @test B ≈ C
    end

    test_vector = [0.001; 0.002; 0.003]
    for precision_type in subtypes(AbstractFloat)
        B = Zero{precision_type}()
        C = Vector{precision_type}(test_vector)
        @test isapprox(B, C, atol=0.01)
    end
end

# ------ B::Real

@testset "!isapprox(B::Real, C::Blade)" begin
    test_value = 5
    test_vectors = [3 3; 4 4; 0 1]
    for precision_type in subtypes(AbstractFloat)
        B = precision_type(test_value)
        C = Blade{precision_type}(test_vectors)
        @test B ≉ C
    end
end

@testset "isapprox(B::Real, C::Pseudoscalar)" begin
    test_value = 5
    test_dim = 10
    for precision_type in subtypes(AbstractFloat)
        B = precision_type(test_value)
        C = Pseudoscalar{precision_type}(test_dim, test_value)
        @test B ≉ C
    end
end

@testset "isapprox(B::Real, C::Scalar)" begin
    test_value = 5
    for precision_type1 in subtypes(AbstractFloat)
        B = precision_type1(test_value)

        # B ≉ C
        C = Scalar{precision_type1}(test_value + 1)
        @test B ≉ C

        # B ≈ C
        for precision_type2 in subtypes(AbstractFloat)
            C = Scalar{precision_type2}(test_value)
            @test B ≈ C
        end
    end
end

@testset "isapprox(B::Real, C::One)" begin
    test_value = 5
    for precision_type1 in subtypes(AbstractFloat)
        # B ≉ C
        B = precision_type1(test_value)
        C = One{precision_type1}()
        @test B ≉ C

        # B ≈ C
        B = precision_type1(1)
        for precision_type2 in subtypes(AbstractFloat)
            C = One{precision_type2}()
            @test B ≈ C
        end
    end
end

@testset "isapprox(B::Real, C::Zero)" begin
    test_value = 5
    for precision_type1 in subtypes(AbstractFloat)
        # B ≉ C
        B = precision_type1(test_value)
        C = Zero{precision_type1}()
        @test B ≉ C

        # B ≈ C
        B = precision_type1(0)
        for precision_type2 in subtypes(AbstractFloat)
            C = Zero{precision_type2}()
            @test B ≈ C
        end
    end
end

# ------ B::Vector

@testset "isapprox(B::Vector, C::Blade)" begin
    # --- Preparations

    test_vector = [3; 4; 1]

    # B ≉ C
    test_vectors = [3 3; 4 4; 0 1]
    for precision_type in subtypes(AbstractFloat)
        B = Vector{precision_type}(test_vector)
        C = Blade{precision_type}(test_vectors)
        @test B ≉ C
    end

    test_vectors = [3; 4; 10]
    for precision_type in subtypes(AbstractFloat)
        B = Vector{precision_type}(test_vector)
        C = Blade{precision_type}(test_vectors)
        @test B ≉ C
    end

    # B ≈ C
    for precision_type1 in subtypes(AbstractFloat)
        for precision_type2 in subtypes(AbstractFloat)
            B = Vector{precision_type1}(test_vector)
            C = Blade{precision_type2}(test_vector)
            @test B ≈ C
        end
    end
end

@testset "isapprox(B::Vector, C::Pseudoscalar)" begin
    test_value = 5

    test_dim_1 = 3
    test_vector_1 = [3; 4; 1]

    test_dim_2 = 1
    test_vector_2 = [test_value]

    for precision_type1 in subtypes(AbstractFloat)
        # B ≉ C
        B = Vector{precision_type1}(test_vector_1)
        C = Pseudoscalar{precision_type1}(test_dim_1, test_value)
        @test B ≉ C

        # B ≈ C
        B = Vector{precision_type1}(test_vector_2)
        for precision_type2 in subtypes(AbstractFloat)
            C = Pseudoscalar{precision_type2}(test_dim_2, test_value)
            @test B ≈ C
        end
    end
end

@testset "isapprox(B::Vector, C::Scalar)" begin
    test_vector = [3; 4; 1]
    test_value = 5
    for precision_type in subtypes(AbstractFloat)
        B = Vector{precision_type}(test_vector)
        C = Scalar{precision_type}(test_value)
        @test B ≉ C
    end
end

@testset "isapprox(B::Vector, C::One)" begin
    test_vector = [3; 4; 1]
    for precision_type in subtypes(AbstractFloat)
        B = Vector{precision_type}(test_vector)
        C = One{precision_type}()
        @test B ≉ C
    end
end

@testset "isapprox(B::Vector, C::Zero)" begin
    # B ≉ C
    test_vector = [3; 4; 1]
    for precision_type in subtypes(AbstractFloat)
        B = Vector{precision_type}(test_vector)
        C = Zero{precision_type}()
        @test B ≉ C
    end

    test_vector = [0.001; 0.002; 0.003]
    for precision_type in subtypes(AbstractFloat)
        B = Vector{precision_type}(test_vector)
        C = Zero{precision_type}()
        @test !isapprox(B, C, atol=0.001)
    end

    # B ≈ C
    test_vector = [0; 0; 0]
    for precision_type in subtypes(AbstractFloat)
        B = Vector{precision_type}(test_vector)
        C = Zero{precision_type}()
        @test B ≈ C
    end

    test_vector = [0.001; 0.002; 0.003]
    for precision_type in subtypes(AbstractFloat)
        B = Vector{precision_type}(test_vector)
        C = Zero{precision_type}()
        @test isapprox(B, C, atol=0.01)
    end
end


text:
export AbstractPolyhedron

"""
    AbstractPolyhedron{N} <: LazySet{N}

Abstract type for compact convex polyhedral sets.

### Notes

Every concrete `AbstractPolyhedron` must define the following functions:
- `constraints_list(::AbstractPolyhedron{N})` -- return a list of all facet
    constraints

```jldoctest; setup = :(using LazySets: subtypes)
julia> subtypes(AbstractPolyhedron)
8-element Vector{Any}:
 AbstractPolytope
 HPolyhedron
 HalfSpace
 Hyperplane
 Line
 Line2D
 Star
 Universe
```

Polyhedra are defined as the intersection of a finite number of closed
half-spaces.
As such, polyhedra are closed and convex but not necessarily bounded.
Bounded polyhedra are called *polytopes* (see [`AbstractPolytope`](@ref)).
"""
abstract type AbstractPolyhedron{N} <: LazySet{N} end

isconvextype(::Type{<:AbstractPolyhedron}) = true

# --- common AbstractPolyhedron functions ---

# supertype for all linear map algorithms
abstract type AbstractLinearMapAlgorithm end

# To account for the compilation order, functions are defined in the file
# AbstractPolyhedron_functions.jl


text:
using FASTQDemultiplexer
using Bio.Seq.FASTQ
using Base.Test

import FASTQDemultiplexer: idtoname, gen_id, Barcodes, Barcode,
    Protocol, InterpretedRecord, getgroup, interpret!

# utils
@testset "utils" begin
    ain = UInt8['A','C','T','A','C','A','T','G','A','C']
    id = 0x00000000000184b2
    @test gen_id(ain) == id
    @test idtoname(id,length(ain))==String(ain)
end

# barcodes
@testset "barcodes" begin
    ain1 = UInt8['A','C','T','A','C','A','T','G','A','C']
    ain2 = reverse(ain1)
    ain3 = UInt8['C','A','A','A','A','A','T','G','A','C']
    ain4 = reverse(ain3)
    ainnone = UInt8['A','A','A','A','A','A','A','A','A','A']

    group1 = Set(map(gen_id,[ain1, ain2]))
    group2 = Set(map(gen_id,[ain3, ain4]))

    barcodes = Barcodes((group1,group2),("group1","group2"))

    @test typeof(Barcode(ain1)) == Barcode{length(ain1)}

    @test getgroup(barcodes,Barcode(ain1)) == (1, "group1")
    @test getgroup(barcodes,Barcode(ain4)) == (2, "group2")
    @test getgroup(barcodes,Barcode(ainnone)) == (-1, "unassigned")

    barcodes2 = Barcodes(["data/barcodes.txt"],["some-barcodes"])
    @test getgroup(barcodes2,Barcode(Vector{UInt8}("AAGG"))) == (1,"some-barcodes")
end

@testset "protocol" begin

    mars = Protocol(:marsseq)
    drop = Protocol(:dropseq)
    x10 = Protocol(:x10)

    @test typeof(mars) == Protocol{2,10,6}
    @test typeof(drop) == Protocol{2,12,8}
    @test typeof(x10) == Protocol{3,16,10}

    @test mars.insertread == 1
    @test mars.insertpos == 8:60
    @test mars.cellpos == (4:7,1:6)
    @test mars.umipos == (1:0,7:12)

    @test drop.insertread == 2
    @test drop.insertpos == 1:39
    @test drop.cellpos == (1:12,1:0)
    @test drop.umipos == (13:20,1:0)

    @test x10.insertread  == 2
    @test x10.insertpos == 1:98
    @test x10.cellpos == (1:16,1:0,1:0)
    @test x10.umipos == (17:26,1:0,1:0)

end

@testset "record" begin
    N = 2; C = 4; U = 2;
    inter = Protocol{N,C,U}(1,5:12,(2:3,3:4),(1:0,1:2),("R1","R2"))
    ir = InterpretedRecord(inter)

    ir.records = map(("data/r1.fastq","data/r2.fastq")) do file
        open(FASTQ.Reader,file) do f
            read(f)
        end
    end

    @test typeof(ir) == InterpretedRecord{N,C,U}
    @test length(ir.umi) == U
    @test length(ir.cell) == C

    cell = zeros(UInt8,C)
    FASTQDemultiplexer.extract!(cell,ir.records,inter.cellpos)
    @test cell == Vector{UInt8}("AAGG")

    umi = zeros(UInt8,U)
    FASTQDemultiplexer.extract!(umi,ir.records,inter.umipos)
    @test umi == Vector{UInt8}("TT")

    bcodes = Barcodes(["data/barcodes.txt"],["some-barcodes"])

    interpret!(ir,inter,bcodes)

    @test ir.cell == cell
    @test ir.umi == umi
    @test ir.cellid == Barcode(cell)
    @test ir.umiid == Barcode(umi)
    @test ir.groupid == 1
    @test ir.groupname == "some-barcodes"

    # read the second record
    ir.records = map(("data/r1.fastq","data/r2.fastq")) do file
        open(FASTQ.Reader,file) do f
            read(f)
            read(f)
        end
    end

    cell = zeros(UInt8,C)
    FASTQDemultiplexer.extract!(cell,ir.records,inter.cellpos)
    @test cell == Vector{UInt8}("CCTT")

    umi = zeros(UInt8,U)
    FASTQDemultiplexer.extract!(umi,ir.records,inter.umipos)
    @test umi == Vector{UInt8}("GG")

    interpret!(ir,inter,bcodes)

    @test ir.cell == cell
    @test ir.umi == umi
    @test ir.cellid == Barcode(cell)
    @test ir.umiid == Barcode(umi)
    @test ir.groupid == -1
    @test ir.groupname == "unassigned"
end

max_stars_repo_name: pwl/FASTQDemultiplexer

max_issues_repo_name: pwl/FASTQDemultiplexer

max_forks_repo_name: pwl/FASTQDemultiplexer


text:
module Gondola

using ..Ahorn, Maple

const placements = Ahorn.PlacementDict(
    "Gondola" => Ahorn.EntityPlacement(
        Maple.Gondola,
        "line",
    )
)

frontSprite = "objects/gondola/front"
backSprite = "objects/gondola/back"
topSprite = "objects/gondola/top"
leverSprite = "objects/gondola/lever01"
leftSprite = "objects/gondola/cliffsideLeft"
rightSprite = "objects/gondola/cliffsideRight"

renderOffsetY = -64

Ahorn.nodeLimits(entity::Maple.Gondola) = 1, 1

function Ahorn.selection(entity::Maple.Gondola)
    x, y = Ahorn.position(entity)
    nodes = get(entity.data, "nodes", ((x + 128, y + 128)))
    destX, destY = Int.(nodes[1])

    return Ahorn.Rectangle[
        Ahorn.getSpriteRectangle(frontSprite, x, y + renderOffsetY, jx=0.5, jy=0.0),
        Ahorn.getSpriteRectangle(rightSprite, destX + 144, destY - 104, jx=0.0, jy=0.5, sx=-1)
    ]
end

function Ahorn.renderAbs(ctx::Ahorn.Cairo.CairoContext, entity::Maple.Gondola, room::Maple.Room)
    x, y = Ahorn.position(entity)

    nodes = get(entity.data, "nodes", ())

    if !isempty(nodes)
        top = Ahorn.getSprite(topSprite, "Gameplay")

        destX, destY = Int.(nodes[1])

        wireLeftX = x - 124 + 40
        wireLeftY = y - 12

        wireRightX = destX + 144 - 40
        wireRightY = destY - 104 - 4

        topX = x
        topY = y + renderOffsetY + top.height

        Ahorn.Cairo.save(ctx)

        Ahorn.setSourceColor(ctx, (0.0, 0.0, 0.0, 1.0))
        Ahorn.set_antialias(ctx, 1)
        Ahorn.set_line_width(ctx, 1);

        for i in 0:1
            Ahorn.move_to(ctx, wireLeftX, wireLeftY + i)
            Ahorn.line_to(ctx, topX, topY + i)

            Ahorn.stroke(ctx)

            Ahorn.move_to(ctx, topX, topY + i)
            Ahorn.line_to(ctx, wireRightX, wireRightY + i)

            Ahorn.stroke(ctx)
        end

        # Don't rotate the top sprite, it's super broken in game
        Ahorn.drawSprite(ctx, frontSprite, x, y + renderOffsetY, jx=0.5, jy=0.0)
        Ahorn.drawSprite(ctx, topSprite, x, y + renderOffsetY, jx=0.5, jy=0.0)
        Ahorn.drawSprite(ctx, leverSprite, x, y + renderOffsetY, jx=0.5, jy=0.0)
        Ahorn.drawSprite(ctx, backSprite, x, y + renderOffsetY, jx=0.5, jy=0.0)

        Ahorn.drawSprite(ctx, leftSprite, x - 124, y, jx=0.0, jy=1.0)
        Ahorn.drawSprite(ctx, rightSprite, destX + 144, destY - 104, jx=0.0, jy=0.5, sx=-1)

        Ahorn.Cairo.restore(ctx)
    end
end

end


text:
#=
julia --color=yes guruguru.jl
=#

using Crayons
using Crayons.Box
using ColorTypes

const SPACE = " " # one white space
BG_COLORS=[
    #BLACK_BG         ,
    RED_BG           ,
    GREEN_BG         ,
    YELLOW_BG        ,
    BLUE_BG          ,
    MAGENTA_BG       ,
    CYAN_BG          ,
    LIGHT_GRAY_BG    ,
    DEFAULT_BG       ,
    DARK_GRAY_BG     ,
    LIGHT_RED_BG     ,
    LIGHT_GREEN_BG   ,
    LIGHT_YELLOW_BG  ,
    LIGHT_BLUE_BG    ,
    LIGHT_MAGENTA_BG ,
    LIGHT_CYAN_BG    ,
    WHITE_BG         ,
]

function draw(c,H,W; init::Bool=false)
    H,W = 5,10
    buf = IOBuffer()
    if !init
        print(buf, "\x1b[999D\x1b[$(H)A") #rollback H-times
    end
    for i in 1:H
        msg = repeat("$(rand(1:9))",W)
        println(buf, c, msg, Crayon(reset=true))
    end
    print(buf |> take! |> String)
end

function clean(H)
    buf = IOBuffer()
    for i in 1:H+1
        print(buf, "\x1b[2K") # clear line
        print(buf, "\x1b[999D\x1b[$(1)A") # rollback
    end
    print(buf |> take! |> String)
end

function partypeople()
    H,W=5,10
    draw(BLACK_BG,H,W,init=true)
    for c in BG_COLORS
        draw(c,H,W)
        sleep(0.1)
    end
    clean(H)
end

partypeople()

max_stars_repo_name: terasakisatoshi/juliaExer

max_issues_repo_name: terasakisatoshi/juliaExer

max_forks_repo_name: terasakisatoshi/juliaExer


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

JLLWrappers.@generate_main_file_header("spglib")
JLLWrappers.@generate_main_file("spglib", UUID("ac4a9f1e-bdb2-5204-990c-47c8b2f70d4e"))
end  # module spglib_jll


text:
function glcm(img::AbstractArray{T, 2}, distance::Integer, angle::Real, mat_size::Integer = 16) where T<:Real
    max_img = maximum(img)
    img_rescaled = map(i -> max(1, Int(ceil(i * mat_size / max_img))), img)
    _glcm(img_rescaled, distance, angle, mat_size)
end

function glcm(img::AbstractArray{T, 2}, distances::Array{Int, 1}, angles::Array{A, 1}, mat_size::Integer = 16) where {T<:Real, A<:Real}
    max_img = maximum(img)
    img_rescaled = map(i -> max(1, Int(ceil(i * mat_size / max_img))), img)
    glcm_matrices = [_glcm(img_rescaled, d, a, mat_size) for d in distances, a in angles]
    glcm_matrices
end

function glcm(img::AbstractArray{T, 2}, distance::Integer, angle::Real, mat_size::Integer = 16) where T<:Colorant
    img_rescaled = map(i -> max(1, Int(ceil((reinterpret(gray(convert(Gray{N0f8}, i)))) * mat_size / 256))), img)
    _glcm(img_rescaled, distance, angle, mat_size)
end

function glcm(img::AbstractArray{T, 2}, distances::Array{Int, 1}, angles::Array{A, 1}, mat_size::Integer = 16) where {T<:Colorant, A<:Real}
    img_rescaled = map(i -> max(1, Int(ceil((reinterpret(gray(convert(Gray{N0f8}, i)))) * mat_size / 256))), img)
    glcm_matrices = [_glcm(img_rescaled, d, a, mat_size) for d in distances, a in angles]
    glcm_matrices
end

glcm(img::AbstractArray{T, 2}, distances::Int, angles::Array{A, 1}, mat_size::Integer = 16) where {T<:Union{Colorant, Real}, A<:Real} = glcm(img, [distances], angles, mat_size)
glcm(img::AbstractArray{T, 2}, distances::Array{Int, 1}, angles::Real, mat_size::Integer = 16) where {T<:Union{Colorant, Real}} = glcm(img, distances, [angles], mat_size)

function _glcm(img::AbstractArray{T, 2}, distance::Integer, angle::Number, mat_size::Integer) where T
    co_oc_matrix = zeros(Int, mat_size, mat_size)
    sin_angle = sin(angle)
    cos_angle = cos(angle)
    for j = 1:size(img, 1)
        j_off = j + Int(round(sin_angle * distance))
        checkbounds(Bool, img, j_off, 1) || continue
        for i = 1:size(img, 2)
            int_one = img[j, i]
            i_off = i + Int(round(cos_angle * distance))

            if checkbounds(Bool, img, j_off, i_off)
                int_two = img[j_off, i_off]
                co_oc_matrix[int_one, int_two] += 1
            end
        end
    end
    co_oc_matrix
end

function glcm_symmetric(img::AbstractArray, distance::Integer, angle::Real, mat_size)
    co_oc_matrix = glcm(img, distance, angle, mat_size)
    co_oc_matrix_trans = co_oc_matrix'
    co_oc_matrix_symm = co_oc_matrix + co_oc_matrix_trans
    co_oc_matrix_symm
end

function glcm_symmetric(img::AbstractArray, distances, angles, mat_size)
    co_oc_matrices = glcm(img, distances, angles, mat_size)
    co_oc_matrices_sym = map(gmat -> gmat + gmat', co_oc_matrices)
    co_oc_matrices_sym
end

function glcm_norm(img::AbstractArray, distance::Integer, angle::Real, mat_size)
    co_oc_matrix = glcm(img, distance, angle, mat_size)
    co_oc_matrix_norm = co_oc_matrix / Float64(sum(co_oc_matrix))
    co_oc_matrix_norm
end

function glcm_norm(img::AbstractArray, distances, angles, mat_size)
    co_oc_matrices = glcm(img, distances, angles, mat_size)
    co_oc_matrices_norm = map(gmat -> gmat /= Float64(sum(gmat)), co_oc_matrices)
    co_oc_matrices_norm
end

function glcm_prop(gmat::Array{T, 2}, window_height::Integer, window_width::Integer, property::Function) where T<:Real
    k_h = Int(floor(window_height / 2))
    k_w = Int(floor(window_width / 2))
    glcm_size = size(gmat)
    R = CartesianIndices(glcm_size)
    prop_mat = zeros(Float64, glcm_size)
    for I in R
        prop_mat[I] = property(gmat[max(1, I[1] - k_h) : min(glcm_size[1], I[1] + k_h), max(1, I[2] - k_w) : min(glcm_size[2], I[2] + k_w)])
    end
    prop_mat
end

glcm_prop(gmat::Array{T, 2}, window_size::Integer, property::Function) where {T<:Real} = glcm_prop(gmat, window_size, window_size, property)

function glcm_prop(gmat::Array{T, 2}, property::Function) where T<:Real
    property(gmat)
end

function contrast(glcm_window::Array{T, 2}) where T<:Real
    sum([(id[1] - id[2]) ^ 2 * glcm_window[id] for id in CartesianIndices(size(glcm_window))])
end

function dissimilarity(glcm_window::Array{T, 2}) where T<:Real
    sum([glcm_window[id] * abs(id[1] - id[2]) for id in CartesianIndices(size(glcm_window))])
end

function glcm_entropy(glcm_window::Array{T, 2}) where T<:Real
    -sum(map(i -> i * log(i), glcm_window))
end

function ASM(glcm_window::Array{T, 2}) where T<:Real
    sum(map(i -> i ^ 2, glcm_window))
end

function IDM(glcm_window::Array{T, 2}) where T<:Real
    sum([glcm_window[id] / (1 + (id[1] - id[2]) ^ 2) for id in CartesianIndices(size(glcm_window))])
end

function glcm_mean_ref(glcm_window::Array{T, 2}) where T<:Real
    sumref = sum(glcm_window, dims=2)
    meanref = sum([id * sumref[id] for id = 1:size(glcm_window)[1]])
    meanref
end

function glcm_mean_neighbour(glcm_window::Array{T, 2}) where T<:Real
    sumneighbour = sum(glcm_window, dims=1)
    meanneighbour = sum([id * sumneighbour[id] for id = 1:size(glcm_window)[2]])
    meanneighbour
end

function glcm_var_ref(glcm_window::Array{T, 2}) where T<:Real
    mean_ref = glcm_mean_ref(glcm_window)
    sumref = sum(glcm_window, dims=2)
    var_ref = sum([((id - mean_ref) ^ 2) * sumref[id] for id = 1:size(glcm_window)[1]])
    var_ref ^ 0.5
end

function glcm_var_neighbour(glcm_window::Array{T, 2}) where T<:Real
    mean_neighbour = glcm_mean_neighbour(glcm_window)
    sumneighbour = sum(glcm_window, dims=1)
    var_neighbour = sum([((id - mean_neighbour) ^ 2) * sumneighbour[id] for id = 1:size(glcm_window)[2]])
    var_neighbour ^ 0.5
end

function correlation(glcm_window::Array{T, 2}) where T<:Real
    mean_ref = glcm_mean_ref(glcm_window)
    var_ref = glcm_var_ref(glcm_window)
    mean_neighbour = glcm_mean_neighbour(glcm_window)
    var_neighbour = glcm_var_neighbour(glcm_window)
    if var_ref == 0 || var_neighbour == 0
        return 1
    end
    sum([glcm_window[i, j] * (i - mean_ref) * (j - mean_neighbour) for i = 1:size(glcm_window)[1], j = 1:size(glcm_window)[2]]) / (var_ref * var_neighbour)
end

function max_prob(glcm_window::Array{T, 2}) where T<:Real
    maxfinite(glcm_window)
end

function energy(glcm_window::Array{T, 2}) where T<:Real
    ASM(glcm_window) ^ 0.5
end


text:
# This file is a part of BAT.jl, licensed under the MIT License (MIT).

# ToDo: This is type piracy, find a cleaner solution!
@recipe function f(h::Histogram,
                parsel::NTuple{2,Integer};
                intervals = standard_confidence_vals,
                colors = standard_colors,
                diagonal = Dict(),
                upper = Dict(),
                right = Dict())

    pi_x, pi_y = parsel

    seriestype = get(plotattributes, :seriestype, :histogram2d)

    if seriestype == :histogram2d || seriestype == :histogram || seriestype == :hist
        _plots_module() != nothing || throw(ErrorException("Package Plots not available, but required for this operation"))

        @series begin
            seriestype := :bins2d
            colorbar --> true
            h.edges[1], h.edges[2], _plots_module().Surface(h.weights)
        end

    elseif seriestype == :smallest_intervals_contour || seriestype == :smallest_intervals_contourf
        _plots_module() != nothing || throw(ErrorException("Package Plots not available, but required for this operation"))

        colors = colors[sortperm(intervals, rev=true)]

        if seriestype == :smallest_intervals_contour
            plotstyle = :contour
        else
            plotstyle = :contourf
        end

        lev = calculate_levels(h.weights, intervals)
        x, y, m = bin_centers(h)

        # quick fix: needed when plotting contour on top of histogram
        # otherwise scaling of histogram colorbar would change scaling
        lev = lev/10000
        m = m/10000

        colorbar --> false

       if _plots_module().backend() == _plots_module().PyPlotBackend()
            @series begin
                seriestype := plotstyle
                levels --> lev
                linewidth --> 2
                color --> colors # only works with pyplot
                (x, y, m')
            end
        else
             @series begin
                seriestype := plotstyle
                levels --> lev
                linewidth --> 2
                (x, y, m')
            end
        end

    elseif seriestype == :smallest_intervals
        _plots_module() != nothing || throw(ErrorException("Package Plots not available, but required for this operation"))

        colors = colors[sortperm(intervals, rev=true)]

        hists, orig_hist, realintervals = split_smallest(h, intervals)

        for (i, int) in enumerate(realintervals)
            @series begin
                seriestype := :bins2d
                color --> _plots_module().cgrad([colors[i], colors[i]])
                label --> "smallest $(@sprintf("%.2f", realintervals[i]*100))% interval(s)"
                hists[i].edges[1], hists[i].edges[2], _plots_module().Surface(hists[i].weights)
            end
            # fake a legend
            @series begin
                seriestype := :shape
                fillcolor --> colors[i]
                linewidth --> 0
                label --> "smallest $(@sprintf("%.2f", realintervals[i]*100))% interval(s)"
                colorbar --> false
                [hists[i].edges[1][1], hists[i].edges[1][1]], [hists[i].edges[2][1], hists[i].edges[2][1]]
            end
        end

    elseif seriestype == :marginal
        _plots_module() != nothing || throw(ErrorException("Package Plots not available, but required for this operation"))

        size --> (900, 600)

        layout --> _plots_module().grid(2,2, widths=(0.8, 0.2), heights=(0.2, 0.8))
        link --> :both

        if get(diagonal, "seriestype", :histogram) != :histogram
            colorbar --> false
        end

        @series begin
            subplot := 1
            xguide --> "\$p(\\theta_$(pi_x))\$"
            seriestype := get(upper, "seriestype", :histogram)
            bins --> get(upper, "nbins", 200)
            normalize --> get(upper, "normalize", true)
            colors --> get(upper, "colors", standard_colors)
            intervals --> get(upper, "intervals", standard_confidence_vals)
            #mean --> get(upper, "mean", false)
            #std_dev --> get(upper, "std_dev", false)
            #globalmode --> get(upper, "globalmode", false)
            #localmode --> get(upper, "localmode", false)
            h, pi_x
        end

        # empty plot (needed since @layout macro not available)
        @series begin
            seriestype := :scatter
            subplot := 2
            grid := false
            xaxis := false
            yaxis := false
            markersize := 0.1
            markerstrokewidth := 0
            markeralpha := 1
            markerstrokealpha := 1
            legend := false
            label := ""
            xguide := ""
            yguide := ""
            [(0,0)]
        end

        @series begin
            subplot := 3
            seriestype := get(diagonal, "seriestype", :histogram)
            legend --> false

            normalize --> get(diagonal, "normalize", true)
            bins --> get(diagonal, "nbins", 200)
            colors --> get(diagonal, "colors", standard_colors)
            intervals --> get(diagonal, "intervals", standard_confidence_vals)
            #mean --> get(diagonal, "mean", false)
            #std_dev --> get(diagonal, "std_dev", false)
            #globalmode --> get(diagonal, "globalmode", false)
            #localmode --> get(diagonal, "localmode", false)

            h, (pi_x, pi_y)
        end

        @series begin
            subplot := 4
            seriestype := get(right, "seriestype", :histogram)
            orientation := :horizontal
            xguide --> "\$p(\\theta_$(pi_y))\$"
            normalize --> get(right, "normalize", true)
            bins --> get(right, "nbins", 200)
            colors --> get(right, "colors", standard_colors)
            intervals --> get(right, "intervals", standard_confidence_vals)
            #mean --> get(right, "mean", false)
            #std_dev --> get(right, "std_dev", false)
            #globalmode --> get(right, "globalmode", false)
            #localmode --> get(right, "localmode", false)

            h, pi_y
        end

    else
        error("seriestype $seriestype not supported")
    end

end

#--- rectangle bounds ------------------------
@recipe function f(bounds::HyperRectBounds, parsel::NTuple{2,Integer})
    pi_x, pi_y = parsel

    vol = spatialvolume(bounds)
    vhi = vol.hi[[pi_x, pi_y]]; vlo = vol.lo[[pi_x, pi_y]]
    rect_xy = rectangle_path(vlo, vhi)
    bext = 0.1 * (vhi - vlo)
    xlims = (vlo[1] - bext[1], vhi[1] + bext[1])
    ylims = (vlo[2] - bext[2], vhi[2] + bext[2])

    @series begin
        seriestype := :path
        label --> "bounds"
        color --> :darkred
        alpha --> 0.75
        linewidth --> 2
        xlims --> xlims
        ylims --> ylims
        (rect_xy[:,1], rect_xy[:,2])
    end

    nothing
end


text:
include("matrix.jl")

function normalize!(vec::AbstractVector)
    return vec[:]=vec/norm(vec)
end

max_stars_repo_name: Roger-luo/QuCmp

max_issues_repo_name: Roger-luo/QuCmp

max_forks_repo_name: Roger-luo/QuCmp


text:
module ElectronTests

using Electron, JSServe, URIParser
using JSServe.HTTP: Request
using JSServe: Session, Dependency, @js_str, JSCode
import JSServe: start, evaljs, evaljs_value
using JSServe.Hyperscript: Node, HTMLSVG
using JSServe.DOM
using Base: RefValue
using Test

"""
Our Javascript library, defining a few JS helper functions!
"""
const JSTest = Dependency(:Test, [joinpath(@__DIR__, "Test.js")])

"""
    TestSession(handler)

Main construct, which will lunch an electron browser session, serving the application
created by `handler(testsession)::DOM.div`.
Can be used via the testsession function:
```julia
testsession(handler; url="0.0.0.0", port=8081, timeout=300)) do testsession
    # test code using testsession
end
```
This will ensure proper setup and teardown once done with the tests.
The testsession object allows to communicate with the browser session, run javascript
and get values from the html dom!
"""
mutable struct TestSession
    url::URI
    initialized::Bool
    error_in_handler::Any
    server::JSServe.Server
    window::Electron.Window
    session::Session
    dom::Node{HTMLSVG}
    request::Request

    function TestSession(url::URI)
        return new(url, false, nothing)
    end

    function TestSession(url::URI, server::JSServe.Server, window::Electron.Window, session::Session)
        testsession = new(url, true, nothing, server, window, session)
        return testsession
    end
end

JSServe.session(testsession::TestSession) = testsession.session

function Base.show(io::IO, testsession::TestSession)
    print(io, "TestSession")
end

function check_and_close_display()
    # For some reason, when running code in Atom, it happens very easily,
    # That JSServe display server gets started!
    # Maybe better to PR an option in JSServe to prohibit starting it in the first place
    if isassigned(JSServe.GLOBAL_SERVER) && JSServe.isrunning(JSServe.GLOBAL_SERVER[])
        @warn "closing JSServe display server, which interfers with testing!"
        close(JSServe.GLOBAL_SERVER[])
    end
end

function TestSession(handler; url="0.0.0.0", port=8081, timeout=300)
    check_and_close_display()
    testsession = TestSession(URI(string("http://localhost:", port)))
    testsession.server = JSServe.Server(url, port) do session, request
        try
            dom = handler(session, request)
            testsession.dom = dom
            testsession.session = session
            testsession.request = request
            return DOM.div(JSTest, dom)
        catch e
            testsession.error_in_handler = (e, Base.catch_backtrace())
        end
    end
    try
        start(testsession; timeout=timeout)
        return testsession
    catch e
        close(testsession)
        rethrow(e)
    end
end

"""
```julia
    testsession(f, handler; url="0.0.0.0", port=8081, timeout=300)

testsession(handler; url="0.0.0.0", port=8081, timeout=300)) do testsession
    # test code using testsession
end
```
This function will ensure proper setup and teardown once done with the tests or whenever an error occurs.
The testsession object passed to `f` allows to communicate with the browser session, run javascript
and get values from the html dom!
"""
function testsession(f, handler; kw...)
    testsession = TestSession(handler; kw...)
    try
        f(testsession)
    catch e
        rethrow(e)
    finally
        close(testsession)
    end
end

function testsession(handler; kw...)
    return TestSession(handler; kw...)
end

"""
    wait(testsession::TestSession; timeout=300)

Wait for testsession to be fully loaded!
Note, if you call wait on a fully loaded test
"""
function wait(testsession::TestSession; timeout=300)
    testsession.initialized && return true
    if !testsession.window.exists
        error("Window isn't open, can't wait for testsession to be initialized")
    end
    Electron.toggle_devtools(testsession.window)
    while testsession.window.exists
        # We done!
        isdefined(testsession, :session) && isopen(testsession.session) && break
        if testsession.error_in_handler !== nothing
            e, backtrace = testsession.error_in_handler
            Base.show_backtrace(stderr, backtrace)
            throw(e)
        end
        # Again, we need to sleep instead of just waiting on `take!`
        # But if we don't do this, on an error in serving, we'd wait indefinitely
        # even if the window gets closed...And since Julia can't deal with interrupting
        # Wait, that'd mean killing Julia completely
        sleep(0.01)
    end
    if !isopen(testsession.session)
        error("Window closed before getting a message from serving request")
    end
    on_timeout = "Timed out when waiting for JS to being loaded! Likely an error happend on the JS side, or your testsession is taking longer than $(timeout) seconds. If no error in console, try increasing timeout!"
    tstart = time()
    while time() - tstart < timeout
        if isready(testsession.session.js_fully_loaded)
            # Error on js during init! We can't continue like this :'(
            if testsession.session.init_error[] !== nothing
                throw(testsession.session.init_error[])
            end
            break
        end
        sleep(0.01)
    end
    testsession.initialized = true
    return true
end

"""
    reload!(testsession::TestSession)

Reloads the served application and waits untill all state is initialized.
"""
function reload!(testsession::TestSession; timeout=300)
    check_and_close_display()
    testsession.initialized = true # we need to put it to true, otherwise handler will block!
    # Make 100% sure we're serving something, since otherwise, well block forever
    @assert JSServe.isrunning(testsession.server)
    # Extra long time out for compilation!
    response = JSServe.HTTP.get(string(testsession.url), readtimeout=500)
    @assert response.status == 200
    testsession.initialized = false
    testsession.error_in_handler = nothing
    Electron.load(testsession.window, testsession.url)
    wait(testsession; timeout=timeout)
    @assert testsession.initialized
    return true
end

"""
    start(testsession::TestSession)

Start the testsession and make sure everything is loaded correctly.
Will close all connections, if any error occurs!
"""
function JSServe.start(testsession::TestSession; timeout=300)
    check_and_close_display()
    try
        if !JSServe.isrunning(testsession.server)
            start(testsession.server)
        end
        if !isdefined(testsession, :window) || !testsession.window.exists
            testsession.window = Window()
        end
        reload!(testsession; timeout=timeout)
    catch e
        close(testsession)
        rethrow(e)
    end
    return true
end

"""
    close(testsession::TestSession)

Close the testsession and clean up the state!
"""
function Base.close(testsession::TestSession)
    if isdefined(testsession, :server)
        Electron.close(testsession.server)
    end
    if isdefined(testsession, :window)
        # testsession.window.app.exists && close(testsession.window.app)
        testsession.window.exists && close(testsession.window)
    end
    testsession.initialized = false
end

"""
    evaljs(testsession::TestSession, js::JSCode)
Runs javascript code `js` in testsession.
Will return the return value of `js`. Might return garbage data, if return value
isn't json serializable.

Example:
```julia
evaljs(testsession, js"document.getElementById('the-id')")
```
"""
function evaljs(testsession::TestSession, js::JSCode)
    JSServe.evaljs_value(testsession.session, js)
end

"""
    wait_test(condition)
Waits for condition expression to become true and then tests it!
"""
macro wait_for(condition, timeout=5)
    return quote
        tstart = time()
        while !$(esc(condition)) && (time() - tstart) < $(timeout)
            sleep(0.001)
        end
        @test $(esc(condition))
    end
end

"""
    trigger_keyboard_press(testsession::TestSession, code::String, element=nothing)
Triggers a keyboard press on `element`! If element is `nothing`, the event will be
triggered for the whole `document`!
Find out the key code string to pass at `http://keycode.info/`
"""
function trigger_keyboard_press(testsession::TestSession, code::String, element=nothing)
    JSTest.trigger_keyboard_press(testsession.session, code, element)
end

"""
    trigger_mouse_move(testsession::TestSession, code::String, position::Tuple{Int, Int}, element=nothing)
Triggers a MouseMove event! If element == nothing, it will try to trigger on any canvas element
found in the DOM.
"""
function trigger_mouse_move(testsession::TestSession, position::Tuple{Int, Int}, element=nothing)
    JSTest.trigger_mouse_move(testsession.session, position, element)
end

"""
    query_testid(id::String)
Returns a js string, that queries for `id`.
"""
function query_testid(id::String)
    js"document.querySelector('[data-test-id=$(id)]')"
end

"""
    query_testid(testsession::TestSession, id::String)
Returns a JSObject representing the object found for `id`.
"""
function query_testid(testsession::TestSession, id::String)
    return evaljs_value(testsession, query_testid(id))
end

export @wait_for, evaljs, testsession, trigger_keyboard_press, trigger_mouse_move, query_testid

end # module


text:
@test SpecificV(0.0035, 300.0) ≈ 0.394_913_866E2
@test SpecificV(0.0035u"MPa", 300.0u"K") ≈ 0.394_913_866E2u"m^3/kg"

@test SpecificV(0.0035, 700.0) ≈ 0.923_015_898E2
@test SpecificV(0.0035u"MPa", 700.0u"K") ≈ 0.923_015_898E2u"m^3/kg"

@test SpecificV(30.0, 700.0) ≈ 0.542_946_619E-2
@test SpecificV(30.0u"MPa", 700.0u"K") ≈ 0.542_946_619E-2u"m^3/kg"

@test SpecificH(0.0035, 300.0) ≈ 0.254_991_145E4
@test SpecificH(0.0035u"MPa", 300.0u"K") ≈ 0.254_991_145E4u"kJ/kg"

@test SpecificH(0.0035, 700.0) ≈ 0.333_568_375E4
@test SpecificH(0.0035u"MPa", 700.0u"K") ≈ 0.333_568_375E4u"kJ/kg"

@test SpecificH(30.0, 700.0) ≈ 0.263_149_474E4
@test SpecificH(30.0u"MPa", 700.0u"K") ≈ 0.263_149_474E4u"kJ/kg"

@test SpecificU(0.0035, 300.0) ≈ 0.241_169_160E4
@test SpecificU(0.0035u"MPa", 300.0u"K") ≈ 0.241_169_160E4u"kJ/kg"

@test SpecificU(0.0035, 700.0) ≈ 0.301_262_819E4
@test SpecificU(0.0035u"MPa", 700.0u"K") ≈ 0.301_262_819E4u"kJ/kg"

@test SpecificU(30.0, 700.0) ≈ 0.246_861_076E4
@test SpecificU(30.0u"MPa", 700.0u"K") ≈ 0.246_861_076E4u"kJ/kg"

@test SpecificS(0.0035, 300.0) ≈ 0.852_238_967E1
@test SpecificS(0.0035u"MPa", 300.0u"K") ≈ 0.852_238_967E1u"kJ/kg/K"

@test SpecificS(0.0035, 700.0) ≈ 0.101_749_996E2
@test SpecificS(0.0035u"MPa", 700.0u"K") ≈ 0.101_749_996E2u"kJ/kg/K"

@test SpecificS(30.0, 700.0) ≈ 0.517_540_298E1
@test SpecificS(30.0u"MPa", 700.0u"K") ≈ 0.517_540_298E1u"kJ/kg/K"

@test SpecificCP(0.0035, 300.0) ≈ 0.191_300_162E1
@test SpecificCP(0.0035u"MPa", 300.0u"K") ≈ 0.191_300_162E1u"kJ/kg/K"

@test SpecificCP(0.0035, 700.0) ≈ 0.208_141_274E1
@test SpecificCP(0.0035u"MPa", 700.0u"K") ≈ 0.208_141_274E1u"kJ/kg/K"

@test SpecificCP(30.0, 700.0) ≈ 0.103_505_092E2
@test SpecificCP(30.0u"MPa", 700.0u"K") ≈ 0.103_505_092E2u"kJ/kg/K"

@test SpeedOfSound(0.0035, 300.0) ≈ 0.427_920_172E3
@test SpeedOfSound(0.0035u"MPa", 300.0u"K") ≈ 0.427_920_172E3u"m/s"

@test SpeedOfSound(0.0035, 700.0) ≈ 0.644_289_068E3
@test SpeedOfSound(0.0035u"MPa", 700.0u"K") ≈ 0.644_289_068E3u"m/s"

@test SpeedOfSound(30.0, 700.0) ≈ 0.480_386_523E3
@test SpeedOfSound(30.0u"MPa", 700.0u"K") ≈ 0.480_386_523E3u"m/s"


text:
using Compat
using Test
using MPI

MPI.Init()

#MPI.mpitype_dict[Boundary] = MPI.mpitype_dict[Int]
comm_size = MPI.Comm_size(MPI.COMM_WORLD)
comm_rank = MPI.Comm_rank(MPI.COMM_WORLD) + 1

# send to next higher process, with wraparound
dest = (comm_rank % comm_size) + 1
if comm_rank > 1
  src = comm_rank - 1
else
  src = comm_size
end

# test simple type

mutable struct NotABits
  a::Any
end

@test_throws ArgumentError MPI.type_create(NotABits)

struct Boundary
  c::UInt16  # force some padding to be inserted
  a::Int
  b::UInt8
end

MPI.mpitype(Boundary)

arr = Array{Boundary}(undef, 3)
for i=1:3
  arr[i] = Boundary( (comm_rank + i) % 127, i + comm_rank, i % 64)
end

req_send = MPI.Isend(arr, dest - 1, 1, MPI.COMM_WORLD)

# receive the message
arr_recv = Array{Boundary}(undef, 3)
req_recv = MPI.Irecv!(arr_recv, src - 1, 1, MPI.COMM_WORLD)

MPI.Wait!(req_send)
MPI.Wait!(req_recv)

# check received array
for i=1:3
  bndry_i = arr_recv[i]
  @test bndry_i.a == (src + i)
  @test bndry_i.b == i % 64
  @test bndry_i.c == (src + i) % 127
end

# test nested types
struct Boundary2
  a::UInt32
  b::Tuple{Int, UInt8}
end

MPI.mpitype(Boundary2)

arr = Array{Boundary2}(undef,3)
arr_recv = Array{Boundary2}(undef,3)

for i=1:3
  arr[i] = Boundary2( (comm_rank + i) % 127, ( Int(i + comm_rank), UInt8(i % 64) ) )
end

req_send = MPI.Isend(arr, dest - 1, 1, MPI.COMM_WORLD)
req_recv = MPI.Irecv!(arr_recv, src - 1, 1, MPI.COMM_WORLD)

MPI.Wait!(req_send)
MPI.Wait!(req_recv)

# check received array
for i=1:3
  bndry_i = arr_recv[i]
  @test bndry_i.a == (src + i) % 127
  @test bndry_i.b[1] == (src + i)
  @test bndry_i.b[2] == (i % 64)
end

MPI.Barrier(MPI.COMM_WORLD)
MPI.Finalize()


text:
using Test
using JOLI
using Zygote
using JOLI4Flux

x = rand(Float32, 3)

W = randn(Float32, 3, 3)
W_JOLI = joMatrix(W)

g = gradient(x -> sum(W*x), x)[1]
g_JOLI = gradient(x -> sum(W_JOLI*x), x)[1]

@test isapprox(sum((g - g_JOLI).^2)/sum(g_JOLI.^2), 0f0; atol=1f-8)


text:
module SkeletonTriangulationsTests

using Test
using Gridap.TensorValues
using Gridap.Arrays
using Gridap.Fields
using Gridap.ReferenceFEs
using Gridap.Geometry
using Gridap.Geometry: IN, OUT
using Gridap.Geometry: DiscreteModelMock

model = DiscreteModelMock()

strian = SkeletonTriangulation(model)

trian = get_volume_triangulation(strian)

fun(x) = sin(pi*x[1])*sin(pi*x[2])

q2x = get_cell_map(trian)

q2fun = compose(fun,q2x)

s2fun = restrict(q2fun,strian)

s = CompressedArray([Point{1,Float64}[(0.25,),(0.75,)]],get_cell_type(strian))

fs = evaluate(s2fun.left-s2fun.right,s)

r = fill(zeros(2),length(fs))

test_array(fs,r)

#function polar(q)
#  r, t, z = q
#  x = r*cos(t)
#  y = r*sin(t)
#  Point(x,y,z)
#end
#domain = (1,2,0,pi,0,0.5)
#partition = (10,30,4)
#model = CartesianDiscreteModel(domain,partition,polar)

domain = (0,1,0,1,0,1)
partition = (3,3,3)
model = CartesianDiscreteModel(domain,partition)

strian = SkeletonTriangulation(model)
test_triangulation(strian)

trian = get_volume_triangulation(strian)

s = CompressedArray([Point{2,Float64}[(0.25,0.25),(0.75,0.75)]],get_cell_type(strian))

s2x = get_cell_map(strian)
x = evaluate(s2x,s)

nvec = get_normal_vector(strian)
nx = evaluate(nvec,s)
collect(nx)

fun(x) = sin(pi*x[1])*cos(pi*x[2])

q2x = get_cell_map(trian)

funq = compose(fun,q2x)

fun_gamma = restrict(funq,strian)

@test isa(fun_gamma, SkeletonPair)

cellids = collect(1:num_cells(trian))

cellids_gamma = reindex(cellids,strian)
@test isa(fun_gamma, SkeletonPair)
@test cellids_gamma.left == get_face_to_cell(strian.left)
@test cellids_gamma.right == get_face_to_cell(strian.right)

ids = get_cell_id(strian)
@test isa(ids,SkeletonPair)

#using Gridap.Visualization
#
#writevtk(trian,"trian")
#writevtk(strian,"strian")
#writevtk(x,"x",nodaldata=["nvec" => nx])

model = DiscreteModelMock()

trian = Triangulation(model)

cell_to_is_left = [true,true,false,true,false]

itrian = InterfaceTriangulation(model,cell_to_is_left)

ltrian = get_left_boundary(itrian)
rtrian = get_right_boundary(itrian)

ni = get_normal_vector(itrian)
nl = get_normal_vector(ltrian)
nr = get_normal_vector(rtrian)

#using Gridap.Visualization
#
#writevtk(trian,"trian")
#writevtk(itrian,"itrian",nsubcells=10,cellfields=["ni"=>ni,"nl"=>nl,"nr"=>nr])

domain = (0,1,0,1)
partition = (10,10)
model = CartesianDiscreteModel(domain,partition)

cell_to_inout = zeros(Int8,num_cells(model))
cell_to_inout[1:13] .= OUT
cell_to_inout[14:34] .= IN

itrian = InterfaceTriangulation(model,cell_to_inout)
@test num_cells(itrian) == 11

itrian = InterfaceTriangulation(model,1:13,14:34)
@test num_cells(itrian) == 11

#ltrian = get_left_boundary(itrian)
#rtrian = get_right_boundary(itrian)
#ni = get_normal_vector(itrian)
#nl = get_normal_vector(ltrian)
#nr = get_normal_vector(rtrian)
#trian = Triangulation(model)
#
#using Gridap.Visualization
#writevtk(trian,"trian",celldata=["inout"=>cell_to_inout])
#writevtk(itrian,"itrian",nsubcells=10,cellfields=["ni"=>ni,"nl"=>nl,"nr"=>nr])

reffe = LagrangianRefFE(Float64,QUAD,(2,2))
conf = CDConformity((CONT,DISC))
face_own_dofs = get_face_own_dofs(reffe,conf)
strian = SkeletonTriangulation(model,reffe,face_own_dofs)
ns = get_normal_vector(strian)

#using Gridap.Visualization
#writevtk(strian,"strian",cellfields=["normal"=>ns])

end # module


text:
format_spec =
[
    (:_type, AbstractString, nothing),
    (:parse, Dict{Any, Any}, nothing),
    (:property, AbstractString, nothing),
    (:feature, AbstractString, nothing),
    (:mesh, AbstractString, nothing),
    (:children, Dict{Any, Any}, nothing)
]

primitivefactory(:VegaFormat, format_spec)

max_stars_repo_name: jpthor/vega


text:
struct LinearObservation{T}
    H::Array{T}
    R::Array{T}
end

function observe(s, lo::LinearObservation)
    return lo.H * s + Gaussian(lo.R)
end


text:
using CategoricalArrays

struct OrderedRuleSet <: ParticleClassifier
    name::String
    rules::Vector{Tuple{String,Function}}

    OrderedRuleSet(name::String, rules::Tuple{String,Function}...) = new(name, collect(rules))

    function OrderedRuleSet(orss::OrderedRuleSet...)
        rules = Vector{Tuple{String,Function}}()
        name = String[]
        for ors in orss
            append!(rules, ors.rules)
            push!(name, ors.name)
        end
        return new(join(name, "+"), rules)
    end
end

Base.show(io::IO, ors::OrderedRuleSet) = print(io, "JORS[$(ors.name)]")
classnames(sr::OrderedRuleSet) = collect(map(r -> r[1], sr.rules))

"""
    classify(zep::Zeppelin, ruleset::AbstractString, classnames::AbstractVector, clsidx::AbstractVector{Int})

Constructs a new Zeppelin item and replaces the "CLASS" column with the classes named in `classnames` and indexed into by `clsidx`.
"""
function classify(zep::Zeppelin, ruleset::AbstractString, classnames::AbstractVector{<:AbstractString}, clsidx::AbstractVector{<:Integer})
    @assert length(clsidx)==nrow(zep.data) "The number of assigned classes must match the number of rows."
    @assert minimum(clsidx)>=-1 "The minimum class index must be -1 or greater."
    @assert maximum(clsidx)<length(classnames) "The maximum class index must less than the number of class names."
    clscol = findfirst(isequal("CLASS"), names(zep.data))
    # Remove "CLASS" and "VERIFIED_CLASS" columns
    rd = copy(zep.data)[:, Cols(x -> !(x in ("CLASS", "VERIFIEDCLASS")))]
    result = Zeppelin("$(zep.headerfile[1:end-4])[$ruleset].hdz", copy(zep.header), rd, classnames)
    if isnothing(clscol)
        tmp = findfirst(isequal("FIRSTELM"), rd)
        clscol = isnothing(tmp) ? ncol(rd) + 1 : tmp
    end
    insertcols!(result.data, clscol, "CLASS" => map(idx->ZepClass(result, idx-1), clsidx))
    insertcols!(result.data, clscol+1, "VERIFIEDCLASS" => fill(ZepClass(result, -1), nrow(rd)))
    return result
end

"""
    classify(zep::Zeppelin, ruleset::AbstractString, catcls::AbstractVector{String})

Create a new Zeppelin file and assign the contents of this CategoricalArray to the "CLASS" column.
"""
classify(zep::Zeppelin, ruleset::AbstractString, catcls::CategoricalVector{String}) =
    classify(zep, ruleset, levels(catcls), catcls.refs)

struct _RowAdapter{T <: DataFrames.DataFrameRow}
    row::T
    default::Dict{Symbol,Any}

    function _RowAdapter(r::DataFrames.DataFrameRow, def::Dict{Symbol,<:Any})
        new{typeof(r)}(r, def)
    end
end

Base.getindex(rw::_RowAdapter, key::Symbol) = #
    get(rw.row, key) do
        getindex(rw.default, key)
    end

Base.get(rw::_RowAdapter, key::Symbol, def) = #
    get(rw.row, key) do
        get(rw.default, key, def)
    end

function Base.getproperty(rw::_RowAdapter, key::Symbol)
    if key in ( :row, :default )
        getfield(rw, key)
    else
        get(getfield(rw, :row), key) do
            getindex(getfield(rw, :default), key)
        end
    end
end

Base.sum(rw::_RowAdapter, keys::Symbol...) = sum(key->rw[key], keys)

between(val::AbstractFloat, min::AbstractFloat, max::AbstractFloat) = (val>=min) && (val<=max)
between(rw::_RowAdapter, key::Symbol, min::AbstractFloat, max::AbstractFloat) = between(rw[key], min, max)

"""
    classify(zep::Zeppelin, sr::OrderedRuleSet)::Zeppelin

Classify the rows in `zep`, row-by-row by calling the applying the `OrderedRuleSet` to each row. Returns a new
Zeppelin object.
"""
function classify(zep::Zeppelin, sr::OrderedRuleSet)::Zeppelin
    cat = empty!(categorical(append!(classnames(sr), ["Other", "Missing", "Error"]), compress=false, ordered = true)) # set up the levels
    errcx = 0
    defs = Dict(Symbol(uppercase(e.symbol))=>0.0 for e in elements)
    cres = ThreadsX.map(eachrow(zep.data)) do row
        try
            input = _RowAdapter(row, defs)
            i = findfirst(rule -> rule[2](input), sr.rules)
            isnothing(i) ? "Other" : sr.rules[i][1]
        catch err
            if (errcx += 1) <= 4
                @info "ERROR[$errcx]: " exception=(err, catch_backtrace())
            end
            "Error"
        end
    end
    append!(cat, cres)
    res = classify(zep, repr(sr), cat)
    res.header["RULE_FILE"]=repr(sr)
    return res
end

max_stars_repo_name: NicholasWMRitchie/NeXLParticle

max_issues_repo_name: NicholasWMRitchie/NeXLParticle

max_forks_repo_name: NicholasWMRitchie/NeXLParticle


text:
@testset "Testing `pmemobj_open` class functions" begin
    using Persistence.Lib

    testfile = joinpath(@__DIR__, "test.pool")
    ispath(testfile) && rm(testfile)

    # Just try some simple opening and closing.
    pool = Lib.create(testfile, 2^24)
    @test Lib.isnull(pool) == false

    Lib.close(pool)
    @test Lib.check(testfile) == true
    pool = Lib.open(testfile)
    @test Lib.isnull(pool) == false
    Lib.close(pool)

    # Since the file already exists, if we try to make another one with the same name,
    # we should get a null pool in return.
    pool = Lib.create(testfile, 2^24)
    @test Lib.isnull(pool) == true

    # Now, do some checks that the "layout" argument works
    rm(testfile)
    pool = Lib.create(testfile, 2^24; layout = "layout_1")
    Lib.close(pool)

    pool = Lib.open(testfile; layout = "billy_bob")
    @test Lib.isnull(pool) == true

    pool = Lib.open(testfile; layout = "layout_1")
    @test Lib.isnull(pool) == false

    # Cleanup
    rm(testfile)
end

@testset "Testing `pmemobj_root` class functions" begin
    using Persistence.Lib

    testfile = joinpath(@__DIR__, "test.pool")
    ispath(testfile) && rm(testfile)

    pool = Lib.create(testfile, 2^24)
    # Allocate 4 Int64's for the root.
    oid = Lib.root(pool, 4 * sizeof(Int64), Int64)
    @test !Lib.isnull(oid)

    # Memory should be zeroed out - make sure this is the case using the scary unsafe loads.
    ptr = Lib.direct(oid)
    @test isa(ptr, Ptr{Int64})
    @test unsafe_load(ptr, 1) == 0
    @test unsafe_load(ptr, 2) == 0
    @test unsafe_load(ptr, 3) == 0
    @test unsafe_load(ptr, 4) == 0

    # Now do some unsafe stores - don't worry about performing transactions just yet.
    unsafe_store!(ptr, 1, 1)
    unsafe_store!(ptr, 2, 2)
    unsafe_store!(ptr, 3, 3)
    unsafe_store!(ptr, 4, 4)

    @test unsafe_load(ptr, 1) == 1
    @test unsafe_load(ptr, 2) == 2
    @test unsafe_load(ptr, 3) == 3
    @test unsafe_load(ptr, 4) == 4

    # Close the file and reopen it - make sure these values persist.
    Lib.close(pool)
    pool = Lib.open(testfile)
    this_oid = Lib.root(pool, 4 * sizeof(Int64), Int64)

    # The persistent OID should be the same across loads
    @test this_oid == oid
    ptr = Lib.direct(this_oid)

    @test unsafe_load(ptr, 1) == 1
    @test unsafe_load(ptr, 2) == 2
    @test unsafe_load(ptr, 3) == 3
    @test unsafe_load(ptr, 4) == 4

    # Test resizing - make sure upper memory is zeroed
    oid = Lib.root(pool, 6 * sizeof(Int64), Int64)
    ptr = Lib.direct(oid)
    @test unsafe_load(ptr, 1) == 1
    @test unsafe_load(ptr, 2) == 2
    @test unsafe_load(ptr, 3) == 3
    @test unsafe_load(ptr, 4) == 4
    @test unsafe_load(ptr, 5) == 0
    @test unsafe_load(ptr, 6) == 0

    # Cleanup
    Lib.close(pool)
    rm(testfile)
end

@testset "Testing `oid_is_null` class functions" begin
    using Persistence.Lib

    testfile = joinpath(@__DIR__, "test.pool")
    ispath(testfile) && rm(testfile)

    pool = Lib.create(testfile, 2^24)

    # Create a root and make sure all the conversion functions work like we expect.
    root_oid = Lib.root(pool, sizeof(Int64))
    root_ptr = Lib.direct(root_oid)

    @test Lib.pool(root_oid) == pool
    @test Lib.pool(root_ptr) == pool

    @test Lib.oid(root_ptr) == root_oid

    # Cleanup
    Lib.close(pool)
    rm(testfile)
end


text:
import YAML: load

type DirectoryStructure
    root::AbstractString
    tmp::AbstractString
    src::AbstractString
    bin::AbstractString
    store::AbstractString

    function DirectoryStructure(root)
        new(
            root,
            joinpath(root, "tmp"),
            joinpath(root, "src"),
            joinpath(root, "bin"),
            joinpath(root, "share"),
        )
    end
end

type Config
    dir::DirectoryStructure
    default_playground_path::AbstractString
    default_prompt::AbstractString
    default_shell::AbstractString
    default_git_address::AbstractString
    default_git_revision::AbstractString
    isolated_shell_history::Bool
    isolated_julia_history::Bool

    function Config(kwargs::Dict)
        new(
            DirectoryStructure(kwargs["root"]),
            abspath(kwargs["default_playground_path"]),
            kwargs["default_prompt"],
            haskey(kwargs, "default_shell") ? kwargs["default_shell"] : "",
            kwargs["default_git_address"],
            kwargs["default_git_revision"],
            kwargs["isolated_shell_history"],
            kwargs["isolated_julia_history"],
        )
    end
end

type PlaygroundConfig
    name::AbstractString
    root_path::AbstractString
    log_path::AbstractString
    bin_path::AbstractString
    pkg_path::AbstractString
    julia_path::AbstractString
    default_shell::AbstractString
    isolated_shell_history::Bool
    isolated_julia_history::Bool

    function PlaygroundConfig(config::Config, dir::AbstractString, name)
        root_path = get_playground_dir(config, dir, name)
        new(
            name,
            root_path,
            joinpath(root_path, "log"),
            joinpath(root_path, "bin"),
            joinpath(root_path, "packages"),
            joinpath(root_path, "bin", "julia"),
            config.default_shell,
            config.isolated_shell_history,
            config.isolated_julia_history
        )
    end
end

function init(config::Config)
    mkpath(config.dir.root)
    mkpath(config.dir.tmp)
    mkpath(config.dir.src)
    mkpath(config.dir.bin)
    mkpath(config.dir.store)
end

function load_config(config::AbstractString, root::AbstractString)
    is_path = false

    try
        is_path = ispath(config)
    end

    if is_path
        config_dict = load(open(config))
    else
        config_dict = load(config)
    end

    config_dict["root"] = root

    return Config(config_dict)
end

function create_paths(pg::PlaygroundConfig)
    mkpath(pg.root_path)
    mkpath(pg.bin_path)
    mkpath(pg.log_path)
    mkpath(pg.pkg_path)
end

function set_envs(pg::PlaygroundConfig)
    ENV["PATH"] = "$(pg.bin_path):" * ENV["PATH"]
    ENV["PLAYGROUND_ENV"] = pg.name == "" ? basename(pg.root_path) : pg.name
    ENV["JULIA_PKGDIR"] = pg.pkg_path

    if pg.default_shell != ""
        ENV["SHELL"] = pg.default_shell
    end

    if pg.isolated_julia_history
        ENV["JULIA_HISTORY"] = joinpath(pg.root_path, ".julia_history")
    end

    if pg.isolated_shell_history
        ENV["HISTFILE"] = joinpath(pg.root_path, ".shell_history")
    end
end

config_path() = joinpath(homedir(), ".playground")


text:
const CUBLAS_STATUS_SUCCESS = 0
const CUBLAS_STATUS_NOT_INITIALIZED = 1
const CUBLAS_STATUS_ALLOC_FAILED = 3
const CUBLAS_STATUS_INVALID_VALUE = 7
const CUBLAS_STATUS_ARCH_MISMATCH = 8
const CUBLAS_STATUS_MAPPING_ERROR = 11
const CUBLAS_STATUS_EXECUTION_FAILED = 13
const CUBLAS_STATUS_INTERNAL_ERROR = 14
const CUBLAS_STATUS_NOT_SUPPORTED = 15
const CUBLAS_STATUS_LICENSE_ERROR = 16

const ERROR_MESSAGE = Dict(
    CUBLAS_STATUS_SUCCESS => "SUCCESS",
    CUBLAS_STATUS_NOT_INITIALIZED => "NOT_INITIALIZED",
    CUBLAS_STATUS_ALLOC_FAILED => "ALLOC_FAILED",
    CUBLAS_STATUS_INVALID_VALUE => "INVALID_VALUE",
    CUBLAS_STATUS_ARCH_MISMATCH => "ARCH_MISMATCH",
    CUBLAS_STATUS_MAPPING_ERROR => "MAPPING_ERROR",
    CUBLAS_STATUS_EXECUTION_FAILED => "EXECUTION_FAILED",
    CUBLAS_STATUS_INTERNAL_ERROR => "INTERNAL_ERROR",
    CUBLAS_STATUS_NOT_SUPPORTED => "NOT_SUPPORTED",
    CUBLAS_STATUS_LICENSE_ERROR => "LICENSE_ERROR"
)

const DEFINE = Dict{Symbol,Symbol}()

DEFINE[:cublasCreate] = :cublasCreate_v2
DEFINE[:cublasDestroy] = :cublasDestroy_v2
DEFINE[:cublasGetVersion] = :cublasGetVersion_v2
DEFINE[:cublasSetStream] = :cublasSetStream_v2
DEFINE[:cublasGetStream] = :cublasGetStream_v2
DEFINE[:cublasGetPointerMode] = :cublasGetPointerMode_v2
DEFINE[:cublasSetPointerMode] = :cublasSetPointerMode_v2

# Blas3 RoutinesS

DEFINE[:cublasSnrm2] = :cublasSnrm2_v2
DEFINE[:cublasDnrm2] = :cublasDnrm2_v2
DEFINE[:cublasScnrm2] = :cublasScnrm2_v2
DEFINE[:cublasDznrm2] = :cublasDznrm2_v2

DEFINE[:cublasSdot] = :cublasSdot_v2
DEFINE[:cublasDdot] = :cublasDdot_v2
DEFINE[:cublasCdotu] = :cublasCdotu_v2
DEFINE[:cublasCdotc] = :cublasCdotc_v2
DEFINE[:cublasZdotu] = :cublasZdotu_v2
DEFINE[:cublasZdotc] = :cublasZdotc_v2

DEFINE[:cublasSscal] = :cublasSscal_v2
DEFINE[:cublasDscal] = :cublasDscal_v2
DEFINE[:cublasCscal] = :cublasCscal_v2
DEFINE[:cublasCsscal] = :cublasCsscal_v2
DEFINE[:cublasZscal] = :cublasZscal_v2
DEFINE[:cublasZdscal] = :cublasZdscal_v2

DEFINE[:cublasSaxpy] = :cublasSaxpy_v2
DEFINE[:cublasDaxpy] = :cublasDaxpy_v2
DEFINE[:cublasCaxpy] = :cublasCaxpy_v2
DEFINE[:cublasZaxpy] = :cublasZaxpy_v2

DEFINE[:cublasScopy] = :cublasScopy_v2
DEFINE[:cublasDcopy] = :cublasDcopy_v2
DEFINE[:cublasCcopy] = :cublasCcopy_v2
DEFINE[:cublasZcopy] = :cublasZcopy_v2

DEFINE[:cublasSswap] = :cublasSswap_v2
DEFINE[:cublasDswap] = :cublasDswap_v2
DEFINE[:cublasCswap] = :cublasCswap_v2
DEFINE[:cublasZswap] = :cublasZswap_v2

DEFINE[:cublasIsamax] = :cublasIsamax_v2
DEFINE[:cublasIdamax] = :cublasIdamax_v2
DEFINE[:cublasIcamax] = :cublasIcamax_v2
DEFINE[:cublasIzamax] = :cublasIzamax_v2

DEFINE[:cublasIsamin] = :cublasIsamin_v2
DEFINE[:cublasIdamin] = :cublasIdamin_v2
DEFINE[:cublasIcamin] = :cublasIcamin_v2
DEFINE[:cublasIzamin] = :cublasIzamin_v2

DEFINE[:cublasSasum] = :cublasSasum_v2
DEFINE[:cublasDasum] = :cublasDasum_v2
DEFINE[:cublasScasum] = :cublasScasum_v2
DEFINE[:cublasDzasum] = :cublasDzasum_v2

DEFINE[:cublasSrot] = :cublasSrot_v2
DEFINE[:cublasDrot] = :cublasDrot_v2
DEFINE[:cublasCrot] = :cublasCrot_v2
DEFINE[:cublasCsrot] = :cublasCsrot_v2
DEFINE[:cublasZrot] = :cublasZrot_v2
DEFINE[:cublasZdrot] = :cublasZdrot_v2

DEFINE[:cublasSrotg] = :cublasSrotg_v2
DEFINE[:cublasDrotg] = :cublasDrotg_v2
DEFINE[:cublasCrotg] = :cublasCrotg_v2
DEFINE[:cublasZrotg] = :cublasZrotg_v2

DEFINE[:cublasSrotm] = :cublasSrotm_v2
DEFINE[:cublasDrotm] = :cublasDrotm_v2

DEFINE[:cublasSrotmg] = :cublasSrotmg_v2
DEFINE[:cublasDrotmg] = :cublasDrotmg_v2

# Blas2 Routines

DEFINE[:cublasSgemv] = :cublasSgemv_v2
DEFINE[:cublasDgemv] = :cublasDgemv_v2
DEFINE[:cublasCgemv] = :cublasCgemv_v2
DEFINE[:cublasZgemv] = :cublasZgemv_v2

DEFINE[:cublasSgbmv] = :cublasSgbmv_v2
DEFINE[:cublasDgbmv] = :cublasDgbmv_v2
DEFINE[:cublasCgbmv] = :cublasCgbmv_v2
DEFINE[:cublasZgbmv] = :cublasZgbmv_v2

DEFINE[:cublasStrmv] = :cublasStrmv_v2
DEFINE[:cublasDtrmv] = :cublasDtrmv_v2
DEFINE[:cublasCtrmv] = :cublasCtrmv_v2
DEFINE[:cublasZtrmv] = :cublasZtrmv_v2

DEFINE[:cublasStbmv] = :cublasStbmv_v2
DEFINE[:cublasDtbmv] = :cublasDtbmv_v2
DEFINE[:cublasCtbmv] = :cublasCtbmv_v2
DEFINE[:cublasZtbmv] = :cublasZtbmv_v2

DEFINE[:cublasStpmv] = :cublasStpmv_v2
DEFINE[:cublasDtpmv] = :cublasDtpmv_v2
DEFINE[:cublasCtpmv] = :cublasCtpmv_v2
DEFINE[:cublasZtpmv] = :cublasZtpmv_v2

DEFINE[:cublasStrsv] = :cublasStrsv_v2
DEFINE[:cublasDtrsv] = :cublasDtrsv_v2
DEFINE[:cublasCtrsv] = :cublasCtrsv_v2
DEFINE[:cublasZtrsv] = :cublasZtrsv_v2

DEFINE[:cublasStpsv] = :cublasStpsv_v2
DEFINE[:cublasDtpsv] = :cublasDtpsv_v2
DEFINE[:cublasCtpsv] = :cublasCtpsv_v2
DEFINE[:cublasZtpsv] = :cublasZtpsv_v2

DEFINE[:cublasStbsv] = :cublasStbsv_v2
DEFINE[:cublasDtbsv] = :cublasDtbsv_v2
DEFINE[:cublasCtbsv] = :cublasCtbsv_v2
DEFINE[:cublasZtbsv] = :cublasZtbsv_v2

DEFINE[:cublasSsymv] = :cublasSsymv_v2
DEFINE[:cublasDsymv] = :cublasDsymv_v2
DEFINE[:cublasCsymv] = :cublasCsymv_v2
DEFINE[:cublasZsymv] = :cublasZsymv_v2
DEFINE[:cublasChemv] = :cublasChemv_v2
DEFINE[:cublasZhemv] = :cublasZhemv_v2

DEFINE[:cublasSsbmv] = :cublasSsbmv_v2
DEFINE[:cublasDsbmv] = :cublasDsbmv_v2
DEFINE[:cublasChbmv] = :cublasChbmv_v2
DEFINE[:cublasZhbmv] = :cublasZhbmv_v2

DEFINE[:cublasSspmv] = :cublasSspmv_v2
DEFINE[:cublasDspmv] = :cublasDspmv_v2
DEFINE[:cublasChpmv] = :cublasChpmv_v2
DEFINE[:cublasZhpmv] = :cublasZhpmv_v2

DEFINE[:cublasSger] = :cublasSger_v2
DEFINE[:cublasDger] = :cublasDger_v2
DEFINE[:cublasCgeru] = :cublasCgeru_v2
DEFINE[:cublasCgerc] = :cublasCgerc_v2
DEFINE[:cublasZgeru] = :cublasZgeru_v2
DEFINE[:cublasZgerc] = :cublasZgerc_v2

DEFINE[:cublasSsyr] = :cublasSsyr_v2
DEFINE[:cublasDsyr] = :cublasDsyr_v2
DEFINE[:cublasCsyr] = :cublasCsyr_v2
DEFINE[:cublasZsyr] = :cublasZsyr_v2
DEFINE[:cublasCher] = :cublasCher_v2
DEFINE[:cublasZher] = :cublasZher_v2

DEFINE[:cublasSspr] = :cublasSspr_v2
DEFINE[:cublasDspr] = :cublasDspr_v2
DEFINE[:cublasChpr] = :cublasChpr_v2
DEFINE[:cublasZhpr] = :cublasZhpr_v2

DEFINE[:cublasSsyr2] = :cublasSsyr2_v2
DEFINE[:cublasDsyr2] = :cublasDsyr2_v2
DEFINE[:cublasCsyr2] = :cublasCsyr2_v2
DEFINE[:cublasZsyr2] = :cublasZsyr2_v2
DEFINE[:cublasCher2] = :cublasCher2_v2
DEFINE[:cublasZher2] = :cublasZher2_v2

DEFINE[:cublasSspr2] = :cublasSspr2_v2
DEFINE[:cublasDspr2] = :cublasDspr2_v2
DEFINE[:cublasChpr2] = :cublasChpr2_v2
DEFINE[:cublasZhpr2] = :cublasZhpr2_v2

# Blas3 Routines

DEFINE[:cublasSgemm] = :cublasSgemm_v2
DEFINE[:cublasDgemm] = :cublasDgemm_v2
DEFINE[:cublasCgemm] = :cublasCgemm_v2
DEFINE[:cublasZgemm] = :cublasZgemm_v2

DEFINE[:cublasSsyrk] = :cublasSsyrk_v2
DEFINE[:cublasDsyrk] = :cublasDsyrk_v2
DEFINE[:cublasCsyrk] = :cublasCsyrk_v2
DEFINE[:cublasZsyrk] = :cublasZsyrk_v2
DEFINE[:cublasCherk] = :cublasCherk_v2
DEFINE[:cublasZherk] = :cublasZherk_v2

DEFINE[:cublasSsyr2k] = :cublasSsyr2k_v2
DEFINE[:cublasDsyr2k] = :cublasDsyr2k_v2
DEFINE[:cublasCsyr2k] = :cublasCsyr2k_v2
DEFINE[:cublasZsyr2k] = :cublasZsyr2k_v2
DEFINE[:cublasCher2k] = :cublasCher2k_v2
DEFINE[:cublasZher2k] = :cublasZher2k_v2

DEFINE[:cublasSsymm] = :cublasSsymm_v2
DEFINE[:cublasDsymm] = :cublasDsymm_v2
DEFINE[:cublasCsymm] = :cublasCsymm_v2
DEFINE[:cublasZsymm] = :cublasZsymm_v2
DEFINE[:cublasChemm] = :cublasChemm_v2
DEFINE[:cublasZhemm] = :cublasZhemm_v2

DEFINE[:cublasStrsm] = :cublasStrsm_v2
DEFINE[:cublasDtrsm] = :cublasDtrsm_v2
DEFINE[:cublasCtrsm] = :cublasCtrsm_v2
DEFINE[:cublasZtrsm] = :cublasZtrsm_v2

DEFINE[:cublasStrmm] = :cublasStrmm_v2
DEFINE[:cublasDtrmm] = :cublasDtrmm_v2
DEFINE[:cublasCtrmm] = :cublasCtrmm_v2
DEFINE[:cublasZtrmm] = :cublasZtrmm_v2


text:
"By default every transition from a state to each neighbour costs 1"
defaultcost(s1, s2) = one(Int64)

"By default, the herustic returns 0 (Breadth First Search)"
defaultheuristic(state, goal) = zero(Int64)

"By default, the goal is determined using =="
defaultisgoal(state, goal) = state == goal

"reconstruct the path of states up to the found final node"
function reconstructpath(n)
  res = [n.data]
  while !isnothing(n.parent)
    n = n.parent
    push!(res, n.data)
  end
  return reverse!(res)
end


text:
include("../src/tree.jl")
include("tree.jl")
include("../src/pvtree.jl")
include("pvtree.jl")
import ProfileView
include("test.jl")