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

module Libdl
@doc """
Interface to libdl. Provides dynamic linking support.
""" Libdl

import Base.DL_LOAD_PATH

export DL_LOAD_PATH, RTLD_DEEPBIND, RTLD_FIRST, RTLD_GLOBAL, RTLD_LAZY, RTLD_LOCAL,
    RTLD_NODELETE, RTLD_NOLOAD, RTLD_NOW, dlclose, dlopen, dlopen_e, dlsym, dlsym_e,
    dlpath, find_library, dlext, dllist

"""
    DL_LOAD_PATH

When calling [`dlopen`](@ref), the paths in this list will be searched first, in
order, before searching the system locations for a valid library handle.
"""
DL_LOAD_PATH

# note: constants to match JL_RTLD_* in src/julia.h, translated
#       to system-specific values by JL_RTLD macro in src/dlload.c
const RTLD_LOCAL     = 0x00000001
const RTLD_GLOBAL    = 0x00000002
const RTLD_LAZY      = 0x00000004
const RTLD_NOW       = 0x00000008
const RTLD_NODELETE  = 0x00000010
const RTLD_NOLOAD    = 0x00000020
const RTLD_DEEPBIND  = 0x00000040
const RTLD_FIRST     = 0x00000080

"""
    RTLD_DEEPBIND
    RTLD_FIRST
    RTLD_GLOBAL
    RTLD_LAZY
    RTLD_LOCAL
    RTLD_NODELETE
    RTLD_NOLOAD
    RTLD_NOW

Enum constant for [`dlopen`](@ref). See your platform man page for details, if
applicable.
"""
(RTLD_DEEPBIND, RTLD_FIRST, RTLD_GLOBAL, RTLD_LAZY, RTLD_LOCAL, RTLD_NODELETE, RTLD_NOLOAD, RTLD_NOW)

"""
    dlsym(handle, sym)

Look up a symbol from a shared library handle, return callable function pointer on success.
"""
function dlsym(hnd::Ptr, s::Union{Symbol,AbstractString}; throw_error::Bool = true)
    hnd == C_NULL && throw(ArgumentError("NULL library handle"))
    val = Ref(Ptr{Cvoid}(0))
    symbol_found = ccall(:jl_dlsym, Cint,
        (Ptr{Cvoid}, Cstring, Ref{Ptr{Cvoid}}, Cint),
        hnd, s, val, Int64(throw_error)
    )
    if symbol_found == 0
        return nothing
    end
    return val[]
end

"""
    dlsym_e(handle, sym)

Look up a symbol from a shared library handle, silently return `NULL` pointer on lookup failure.  It is preferred to use dlsym(handle, sym; throw_error=false).
"""
function dlsym_e(hnd::Ptr, s::Union{Symbol,AbstractString})
    return dlsym(hnd, s; throw_error=false)
end

"""
    dlopen(libfile::AbstractString [, flags::Integer])

Load a shared library, returning an opaque handle.

The extension given by the constant `dlext` (`.so`, `.dll`, or `.dylib`)
can be omitted from the `libfile` string, as it is automatically appended
if needed.   If `libfile` is not an absolute path name, then the paths
in the array `DL_LOAD_PATH` are searched for `libfile`, followed by the
system load path.

The optional flags argument is a bitwise-or of zero or more of `RTLD_LOCAL`, `RTLD_GLOBAL`,
`RTLD_LAZY`, `RTLD_NOW`, `RTLD_NODELETE`, `RTLD_NOLOAD`, `RTLD_DEEPBIND`, and `RTLD_FIRST`.
These are converted to the corresponding flags of the POSIX (and/or GNU libc and/or MacOS)
dlopen command, if possible, or are ignored if the specified functionality is not available
on the current platform. The default flags are platform specific. On MacOS the default
`dlopen` flags are `RTLD_LAZY|RTLD_DEEPBIND|RTLD_GLOBAL` while on other platforms the
defaults are `RTLD_LAZY|RTLD_DEEPBIND|RTLD_LOCAL`. An important usage of these flags is to
specify non default behavior for when the dynamic library loader binds library references to
exported symbols and if the bound references are put into process local or global scope. For
instance `RTLD_LAZY|RTLD_DEEPBIND|RTLD_GLOBAL` allows the library's symbols to be available
for usage in other shared libraries, addressing situations where there are dependencies
between shared libraries.

If the library cannot be found, this method returns `nothing`.
"""
function dlopen end

dlopen(s::Symbol, flags::Integer = RTLD_LAZY | RTLD_DEEPBIND; kwargs...) =
    dlopen(string(s), flags; kwargs...)

function dlopen(s::AbstractString, flags::Integer = RTLD_LAZY | RTLD_DEEPBIND; throw_error::Bool = true)
    ret = ccall(:jl_load_dynamic_library, Ptr{Cvoid}, (Cstring,UInt32,Cint), s, flags, Cint(throw_error))
    if ret == C_NULL
        return nothing
    end
    return ret
end

"""
    dlopen_e(libfile::AbstractString [, flags::Integer])

Similar to [`dlopen`](@ref), except returns a `NULL` pointer instead of raising errors.
It is preferred to directly call dlopen(libfile, flags; throw_error=false)`
"""
dlopen_e(args...) = dlopen(args...; throw_error=false)

"""
    dlclose(handle)

Close shared library referenced by handle.
"""
function dlclose(p::Ptr)
    0 == ccall(:jl_dlclose, Cint, (Ptr{Cvoid},), p)
end

"""
    dlclose(::Nothing)

For the very common pattern usage pattern of

    try
        hdl = dlopen(library_name)
        ... do something
    finally
        dlclose(hdl)
    end

We define a `dlclose()` method that accepts a parameter of type `Nothing`, so
that user code does not have to change its behavior for the case that `library_name`
was not found.
"""
function dlclose(p::Nothing)
end

"""
    find_library(names, locations)

Searches for the first library in `names` in the paths in the `locations` list,
`DL_LOAD_PATH`, or system library paths (in that order) which can successfully be dlopen'd.
On success, the return value will be one of the names (potentially prefixed by one of the
paths in locations). This string can be assigned to a `global const` and used as the library
name in future `ccall`'s. On failure, it returns the empty string.
"""
function find_library(libnames, extrapaths=String[])
    for lib in libnames
        for path in extrapaths
            l = joinpath(path, lib)
            p = dlopen(l, RTLD_LAZY; throw_error=false)
            if p !== nothing
                dlclose(p)
                return l
            end
        end
        p = dlopen(lib, RTLD_LAZY; throw_error=false)
        if p !== nothing
            dlclose(p)
            return lib
        end
    end
    return ""
end
find_library(libname::Union{Symbol,AbstractString}, extrapaths=String[]) =
    find_library([string(libname)], extrapaths)

function dlpath(handle::Ptr{Cvoid})
    p = ccall(:jl_pathname_for_handle, Cstring, (Ptr{Cvoid},), handle)
    s = unsafe_string(p)
    Sys.iswindows() && Libc.free(p)
    return s
end

function dlpath(libname::Union{AbstractString, Symbol})
    handle = dlopen(libname)
    path = dlpath(handle)
    dlclose(handle)
    return path
end

if Sys.isapple()
    const dlext = "dylib"
elseif Sys.iswindows()
    const dlext = "dll"
else
    #assume Sys.islinux, or similar
    const dlext = "so"
end

"""
    dlext

File extension for dynamic libraries (e.g. dll, dylib, so) on the current platform.
"""
dlext

if Sys.islinux()
    struct dl_phdr_info
        # Base address of object
        addr::Cuint

        # Null-terminated name of object
        name::Ptr{UInt8}

        # Pointer to array of ELF program headers for this object
        phdr::Ptr{Cvoid}

        # Number of program headers for this object
        phnum::Cshort
    end

    # This callback function called by dl_iterate_phdr() on Linux
    function dl_phdr_info_callback(di::dl_phdr_info, size::Csize_t, dynamic_libraries::Array{AbstractString,1})
        # Skip over objects without a path (as they represent this own object)
        name = unsafe_string(di.name)
        if !isempty(name)
            push!(dynamic_libraries, name)
        end
        return convert(Cint, 0)::Cint
    end
end # linux-only

if Sys.isbsd() && !Sys.isapple()
    # DL_ITERATE_PHDR(3) on freebsd
    struct dl_phdr_info
        # Base address of object
        addr::Cuint

        # Null-terminated name of object
        name::Ptr{UInt8}

        # Pointer to array of ELF program headers for this object
        phdr::Ptr{Cvoid}

        # Number of program headers for this object
        phnum::Cshort
    end

    function dl_phdr_info_callback(di::dl_phdr_info, size::Csize_t, dy_libs::Vector{AbstractString})
        name = unsafe_string(di.name)
        if !isempty(name)
            push!(dy_libs, name)
        end
        return convert(Cint, 0)::Cint
    end
end # bsd family

function dllist()
    dynamic_libraries = Vector{AbstractString}()

    @static if Sys.islinux()
        callback = @cfunction(dl_phdr_info_callback, Cint,
                              (Ref{dl_phdr_info}, Csize_t, Ref{Vector{AbstractString}}))
        ccall(:dl_iterate_phdr, Cint, (Ptr{Cvoid}, Ref{Vector{AbstractString}}), callback, dynamic_libraries)
    end

    @static if Sys.isapple()
        numImages = ccall(:_dyld_image_count, Cint, ())

        # start at 1 instead of 0 to skip self
        for i in 1:numImages-1
            name = unsafe_string(ccall(:_dyld_get_image_name, Cstring, (UInt32,), i))
            push!(dynamic_libraries, name)
        end
    end

    @static if Sys.iswindows()
        ccall(:jl_dllist, Cint, (Any,), dynamic_libraries)
    end

    @static if Sys.isbsd() && !Sys.isapple()
        callback = @cfunction(dl_phdr_info_callback, Cint,
                              (Ref{dl_phdr_info}, Csize_t, Ref{Vector{AbstractString}}))
        ccall(:dl_iterate_phdr, Cint, (Ptr{Cvoid}, Ref{Vector{AbstractString}}), callback, dynamic_libraries)
        popfirst!(dynamic_libraries)
    end

    return dynamic_libraries
end

end # module

max_stars_repo_name: ylxdzsw/julia

max_issues_repo_name: ylxdzsw/julia

max_forks_repo_name: ylxdzsw/julia


text:
# Previous versions assumed that bytesavailable could accurately check for an empty CSV, but
# this doesn't work reliably for streams because bytesavailable only checks buffered bytes
# (see issue #77). This test verifies that even when bytesavailable would return 0 on a stream
# the full stream is still read.

mutable struct MultiStream{S<:IO} <: IO
    streams::Array{S}
    index::Int
end

function MultiStream(streams::AbstractArray{S}) where {S <: IO}
    MultiStream(streams, 1)
end

function refill(s::MultiStream)
    while eof(s.streams[s.index]) && s.index < length(s.streams)
        close(s.streams[s.index])
        s.index += 1
    end
end

function Base.close(s::MultiStream)
    for i in s.index:length(s.streams)
        close(s.streams[i])
    end
    s.index = length(s.streams)
    nothing
end

function Base.eof(s::MultiStream)
    eof(s.streams[s.index]) && s.index == length(s.streams)
end

function Base.read(s::MultiStream, ::Type{UInt8})
    refill(s)
    read(s.streams[s.index], UInt8)::UInt8
end

function Base.bytesavailable(s::MultiStream)
    bytesavailable(s.streams[s.index])
end

stream = MultiStream(
    [IOBuffer(""), IOBuffer("a,b,c\n1,2,3\n"), IOBuffer(""), IOBuffer("4,5,6")]
)

@test bytesavailable(stream) == 0
@test CSV.read(stream) == CSV.read(IOBuffer("a,b,c\n1,2,3\n4,5,6"))

ext: jl


text:
@doc "T_A_1" ->
abstract T_A_1

@doc """T_A_2""" ->
abstract T_A_2

@doc "T_M_1" ->
type T_M_1
end

@doc """T_M_2""" ->
type T_M_2
end

@doc "T_I_1" ->
immutable T_I_1
end

@doc """T_I_2""" ->
immutable T_I_2
end

@doc "BT_1" ->
bitstype 8 BT_1

@doc """BT_2""" ->
bitstype 16 BT_2

@doc "TA_1" ->
typealias TA_1 Int

@doc """TA_2""" ->
typealias TA_2 Int


text:
################################################################################
#  Copyright 2021, Marta Vanin                                                 #
################################################################################
function run_case_study_D(path_to_csv_result::String, nlsolver; power_base::Float64=1e5, rsc::Int64=100)

    # Input data
    rm_transfo = true
    rd_lines = true

    season = "summer"
    time_step = 144
    elm = ["load", "pv"]
    pfs = [0.95, 0.90]

    msr_path = joinpath(mktempdir(),"temp.csv")
    # Set solve
    solver = _PMD.optimizer_with_attributes(nlsolver...)

    df = _DF.DataFrame(solve_time=Float64[], n_bus=Int64[],
                termination_status=String[], objective=Float64[], criterion=String[], rescaler = Float64[], eq_model = String[],
                     err_max_1=Float64[], err_max_2=Float64[], err_max_3=Float64[], err_avg_1 = Float64[], err_avg_2 = Float64[], err_avg_3 = Float64[], n_meas = Float64[], pbase=Float64[])

    # Load the data
    data = _PMD.parse_file(_PMDSE.get_enwl_dss_path(1, 1),data_model=_PMD.ENGINEERING);
    if rm_transfo _PMDSE.rm_enwl_transformer!(data) end
    if rd_lines _PMDSE.reduce_enwl_lines_eng!(data) end
    data["settings"]["sbase_default"] = power_base

    # Insert the load profiles
    _PMDSE.insert_profiles!(data, season, elm, pfs, t = time_step)

    # Transform data model
    data = _PMD.transform_data_model(data);

    # Solve the power flow
    pf_results = _PMD.solve_mc_pf(data, _PMD.ACPUPowerModel, solver)

    # Write measurements based on power flow
    v_pu = data["settings"]["vbases_default"]["1"]* data["settings"]["voltage_scale_factor"] # divider [V] to get the voltage in per units.
    v_max_err = 1.15 # maximum error of voltage measurement = 0.5% or 1.15 V
    σ_v = 1/3*v_max_err/v_pu

    p_pu = data["settings"]["sbase"] # divider [kW] to get the power in per units.
    p_max_err = 0.01 # maximum error of power measurement = 10W, or 0.01 kW
    σ_p = 1/3*p_max_err/p_pu

    # Write measurements based on power flow
    σ_dict = Dict("load" => Dict("load" => σ_p,
                    "bus"  => σ_v),
                    "gen"  => Dict("gen" => σ_p,
                    "bus"  => σ_v)
                            )

    _PMDSE.write_measurements!(_PMD.ACPUPowerModel, data, pf_results, msr_path, σ = σ_dict)

    # Read-in measurement data and set initial values
    _PMDSE.add_measurements!(data, msr_path, actual_meas = false, seed = 2)

    # Set se settings
    data["se_settings"] = Dict{String,Any}("criterion" => "rwlav",
                            "rescaler" => rsc)

    load_ids = []
    bus_ids = []
    for (l, load) in data["load"]
        push!(load_ids, parse(Int64, l))
        push!(bus_ids, load["load_bus"])
    end

    meas_ids = []
    for i in 1:length(load_ids)
        for (m, meas) in data["meas"]
            if meas["var"] == :vm && meas["cmp_id"] == bus_ids[i]
                push!(meas_ids, m)
            elseif meas["var"] != :vm && meas["cmp_id"] == load_ids[i]
                push!(meas_ids, m)
            end
        end
    end
    delete_bus = []
    for i in 1:length(load_ids)
        se_results = _PMDSE.solve_ivr_red_mc_se(data, solver)
        se_res = deepcopy(se_results)
        pf_res = deepcopy(pf_results)

        if i != 1
            for db in delete_bus
                if db == 1
                    pf_res["solution"]["bus"]["1"] = Dict{String, Any}("vm" => [1.0, 1.0, 1.0])
                    se_res["solution"]["bus"]["1"] = Dict{String, Any}("vr" => [1.0, 1.0, 1.0])
                    se_res["solution"]["bus"]["1"]["vi"] = [0.0, 0.0, 0.0]
                else
                    delete!(se_res["solution"]["bus"], "$db")
                    delete!(pf_res["solution"]["bus"], "$db")
                end
            end
        end

        delta_1, delta_2, delta_3, max_1, max_2, max_3, mean_1, mean_2, mean_3 = _PMDSE.calculate_voltage_magnitude_error_perphase(se_results, pf_results)

        # store result
        push!(df, [se_results["solve_time"], length(data["bus"]),
                string(se_results["termination_status"]),
                se_results["objective"], "rwlav", rsc, "rIVR", max_1, max_2, max_3, mean_1, mean_2, mean_3, length(data["meas"]), power_base])

        push!(delete_bus, data["meas"][meas_ids[3*(i-1)+1]]["cmp_id"])
        delete!(data["meas"], meas_ids[3*(i-1)+1])
        delete!(data["meas"], meas_ids[3*(i-1)+2])
        delete!(data["meas"], meas_ids[3*(i-1)+3])
    end
    CSV.write(path_to_csv_result, df)
end


text:
using Test
using BCVUMPS
using BCVUMPS:bcvumps_env,magnetisation,magofβ,energy,eneofβ,Z,Zofβ
using Random

@testset "$(Ni)x$(Nj) ising" for Ni = [1,2,3], Nj = [1,2,3]
    Random.seed!(100)
    model = Ising(Ni,Nj)
    for β = 0.2:0.2:0.8
        env = bcvumps_env(model, β, 2; tol=1e-10, maxiter=20, verbose = false)
        @test isapprox(magnetisation(env,model,β), magofβ(model,β), atol=1e-5)
        @test isapprox(energy(env,model,β), eneofβ(model,β), atol=1e-2)
        @test isapprox(Z(env), Zofβ(model,β), atol=1e-3)
    end
end

@testset "J1-J2-2x2-Ising" begin
    Random.seed!(100)
    model = Ising22(1)
    for β = 0.2:0.2:0.8
        env = bcvumps_env(model, β, 2; tol=1e-10, maxiter=20, verbose = false)
        @test isapprox(magnetisation(env, model, β), magofβ(Ising(),β), atol=1e-5)
        @test isapprox(energy(env, model, β), eneofβ(Ising(),β), atol=1e-2)
        @test isapprox(Z(env), Zofβ(Ising(),β), atol=1e-3)
    end

    Random.seed!(100)
    model = Ising22(2)
    for β = 0.4:0.2:0.8
        env = bcvumps_env(model, β, 10; tol=1e-10, maxiter=20, verbose = false)
        mag22,ene22 = MCMC(model,16,β,10000,100000)
        @test isapprox(magnetisation(env,model,β), mag22, atol=1e-3)
        @test isapprox(energy(env,model,β), ene22, atol=1e-3)
    end
end


text:
"""
for a given type of analysis
    Get element stiffness matrix formulation
        Get B matrix

"""

function get_strain_disp_matrix_3D(
    shfg_xyz::Union{SMatrix,Matrix},
)::Matrix{Float64}
    num_shape_func::Int64 = size(shfg_xyz)[2]
    Bmatrix::Matrix{Float64} = zeros(Float64, 6, 3 * num_shape_func)
    for i in 1:num_shape_func
        dNdX, dNdY, dNdZ = shfg_xyz[:, i]
        Bmatrix[:, ((3*i)-2):(3*i)] = [
            dNdZ 0.0 0.0
            0.0 dNdX 0.0
            0.0 0.0 dNdY
            0.0 dNdY dNdX
            dNdY 0.0 dNdZ
            dNdX dNdZ 0.0
        ]
    end
    return Bmatrix
end

function get_piezo_Bmat(
    shfg_xyz::Union{SMatrix,Matrix},
)::Matrix{Float64}
    num_shape_func::Int64 = size(shfg_xyz)[2]
    Bmatrix::Matrix{Float64} = zeros(Float64, 9, 4 * num_shape_func)
    for i in 1:num_shape_func
        dNdX, dNdY, dNdZ = shfg_xyz[:, i]
        Bmatrix[:, ((4*i)-3):(4*i)] = [
            dNdZ 0.0 0.0 0.0
            0.0 dNdX 0.0 0.0
            0.0 0.0 dNdY 0.0
            0.0 dNdY dNdX 0.0
            dNdY 0.0 dNdZ 0.0
            dNdX dNdZ 0.0 0.0
            0.0 0.0 0.0 -dNdZ
            0.0 0.0 0.0 -dNdX
            0.0 0.0 0.0 -dNdY
        ]
    end
    return Bmatrix
end

function get_MEE_Bmat(
    shfg_xyz::Union{SMatrix,Matrix},
)::Matrix{Float64}
    num_shape_func::Int64 = size(shfg_xyz)[2]
    Bmatrix::Matrix{Float64} = zeros(Float64, 12, 5 * num_shape_func)
    for i in 1:num_shape_func
        dNdX, dNdY, dNdZ = shfg_xyz[:, i]
        Bmatrix[:, ((5*i)-4):(5*i)] = [
            dNdZ 0.0 0.0 0.0 0.0
            0.0 dNdX 0.0 0.0 0.0
            0.0 0.0 dNdY 0.0 0.0
            0.0 dNdY dNdX 0.0 0.0
            dNdY 0.0 dNdZ 0.0 0.0
            dNdX dNdZ 0.0 0.0 0.0
            0.0 0.0 0.0 -dNdZ 0.0
            0.0 0.0 0.0 -dNdX 0.0
            0.0 0.0 0.0 -dNdY 0.0
            0.0 0.0 0.0 0.0 -dNdZ
            0.0 0.0 0.0 0.0 -dNdX
            0.0 0.0 0.0 0.0 -dNdY
        ]
    end
    return Bmatrix
end

function get_strain_disp_matrix_2D(
    shfg_xyz::Union{SMatrix,Matrix}
)::Matrix{Float64}
    num_shape_func::Int64 = size(shfg_xyz)[2]
    Bmatrix::Matrix{Float64} = zeros(Float64, 3, 2 * num_shape_func)
    for i in 1:num_shape_func
        dNdX, dNdY = shfg_xyz[:, i]
        Bmatrix[:, ((2*i)-1):(2*i)] = [
            dNdX 0.0
            0.0 dNdY
            dNdY dNdX
        ]
    end
    return Bmatrix
end

function get_thermal_conductivity_Bmat(
    shfg_xyz::Union{SMatrix,Matrix}
)::Matrix{Float64}
    num_shape_func::Int64 = size(shfg_xyz)[2]
    Bmatrix::Matrix{Float64} = zeros(Float64, 3, num_shape_func)
    for i in 1:num_shape_func
        dNdX, dNdY, dNdZ = shfg_xyz[:, i]
        Bmatrix[:, i:i] = [dNdZ, dNdX, dNdY ]
    end
    return Bmatrix
end

function get_B_matrix(
    shf_grad_xyz::Union{SMatrix,Matrix},
    analysis_type::DataType
)::Matrix{Float64}
    if analysis_type in (PlaneStress_2DFEA, PlaneStrain_2DFEA)
        return get_strain_disp_matrix_2D(shf_grad_xyz)
    elseif analysis_type in (Elastic_3DFEA, Thermo_Elastic_3DFEA, )
        return get_strain_disp_matrix_3D(shf_grad_xyz)
    elseif analysis_type === Piezo_Electric_3DFEA
        return get_piezo_Bmat(shf_grad_xyz)
    elseif analysis_type in (Thermal_Conductivity_3DFEA, )
        return get_thermal_conductivity_Bmat(shf_grad_xyz)
    elseif analysis_type in (Magneto_Electro_Elastic_3DFEA, )
        return get_MEE_Bmat(shf_grad_xyz)
    else
        @error "Encountered unknown analysis type :: $analysis_type, in B matrix evaluation."
    end
end

"""
    preapre_constraints_matrix_csc() -> SparseMatrixCSC

# Args
- **node_tags** :: Vector{T}
- **parent_child_node_tag_pairs** :: Vector{Tuple{T, T}}
- **dof_per_node** :: T

"""
function preapre_constraints_matrix_csc(
    analysis_type::DataType,
    node_tags::Vector{T},
    parent_child_node_tag_pairs::Vector{Tuple{T, T}};
    hdbc::Dict{T, Vector{Symbol}}=Dict{T, Vector{Symbol}}(),
)::SparseMatrixCSC{Float64, T} where {T}
    # ==============================
    # term: meaning
    # -------------
    # rnt: retaining node tag
    # pnt: parent node tag
    # cnt: child node tag
    # hdbc: homogenous Dirichlet boundary conditions
    # ==============================
    dof_quantities = get_nodal_dof(analysis_type)
    dof_per_node = length(dof_quantities)
    child_node_tags::Vector{T} = [i[2] for i in parent_child_node_tag_pairs]
    retaining_node_tags::Vector{T} = [i for i in node_tags if !(i in child_node_tags)]
    total_dof::T = dof_per_node*length(node_tags)
    retaining_dof::T = dof_per_node * length(retaining_node_tags)
    #
    row_indices::Vector{Int64} = T[]
    col_indices::Vector{Int64} = T[]
    for ant in node_tags
        append!(row_indices, (1+((ant-1) * dof_per_node)):(ant*dof_per_node))
        if !(ant in child_node_tags)
            a_rnt_index = findfirst(isequal(ant), retaining_node_tags)
        else
            a_cnt_index = findfirst(isequal(ant), child_node_tags)
            a_pnt = parent_child_node_tag_pairs[a_cnt_index][1]
            a_rnt_index = findfirst(isequal(a_pnt), retaining_node_tags)
        end
        append!(col_indices, (1+((a_rnt_index-1) * dof_per_node)):(a_rnt_index*dof_per_node))
    end
    if length(hdbc) > 0
        #
        # prepare columns where hdbc needs to be imposed.
        rm_column_indices::Vector{T} = T[]
        for (a_hdbc_nt, spec_dof) in hdbc
            a_hdbc_nt_index = findfirst(isequal(a_hdbc_nt), retaining_node_tags)
            if length(spec_dof)==1 && spec_dof[1]==:ALL
                dof_indices = collect(1:dof_per_node)
            else
                dof_indices = [kount for (kount, i) in enumerate(dof_quantities) if i in spec_dof]
            end
            append!(
                rm_column_indices,
                ((1+((a_hdbc_nt_index-1)*dof_per_node)):(a_hdbc_nt_index*dof_per_node))[dof_indices]
            )
        end
        # remove columns
        #
        row_indices_1::Vector{T} = T[]
        col_indices_1::Vector{T} = T[]
        for (ii, jj) in zip(row_indices, col_indices)
            if !(jj in rm_column_indices)
                # reduce jj by jj1 where jj1 is the number of elements in rm_column_indices which are < jj
                jj1 = length(findall(rm_column_indices .< jj))
                push!(row_indices_1, ii)
                push!(col_indices_1, jj - jj1)
            end
        end
        #
        retaining_dof -= length(rm_column_indices)
        #
        return sparse(row_indices_1, col_indices_1, 1.0, total_dof, retaining_dof)
    else
        return sparse(row_indices, col_indices, 1.0, total_dof, retaining_dof)
    end
end


text:
using MutualInformationImageRegistration
using Random, ImageFiltering, ComputationalResources, Test
using MutualInformationImageRegistration.FastHistograms
using JLD2

MAX_SHIFT = 11
padding = [-10, -10, 10, 10]

@testset "MutualInformationImageRegistration.jl" begin
    @testset "register without filtering" begin
        mi = MutualInformationContainer(
            create_fast_histogram(
                FastHistograms.FixedWidth(),
                FastHistograms.Arithmetic(),
                FastHistograms.NoParallelization(),
                [(0x00, 0xff, 8), (0x00, 0xff, 8)],
            ),
        )
        full_image = rand(UInt8, 500, 300)
        view(full_image, 300:330, 200:220) .= 0xff

        fixed = full_image[(300-10):(330+10), (200-10):(220+10)]

        buffer = Array{UInt8}(undef, (size(fixed) .+ (MAX_SHIFT * 2))...)

        for i = 1:1000
            expected_x = rand(-5:5)
            expected_y = rand(-5:5)

            shift, mm, mms = register!(
                mi,
                full_image,
                fixed,
                [300, 200, 330, 220] .+ padding .+
                [expected_x, expected_y, expected_x, expected_y],
                MAX_SHIFT,
                MAX_SHIFT,
                buffer,
            )

            # The shift we get out should be equal and opposite of the shift we applied
            @test shift == (-expected_x, -expected_y)
            if shift != (-expected_x, -expected_y)
                jldsave("register without filtering fail.jld2"; mi, full_image, fixed, buffer, expected_x, expected_y, shift, mm, mms)
                break
            end
        end
    end

    @testset "register with filtering" begin
        function prefilter!(img::Array{UInt8,2})
            buf = Float32.(img)
            buf = imfilter(
                CPU1(Algorithm.FIR()),
                buf,
                (centered(ones(2, 1) ./ 2), centered(ones(1, 2) ./ 2)),
                "reflect",
            )
            img .= round.(UInt8, buf)
            return nothing
        end

        mi = MutualInformationContainer(
            create_fast_histogram(
                FastHistograms.FixedWidth(),
                FastHistograms.Arithmetic(),
                FastHistograms.NoParallelization(),
                [(0x00, 0xff, 8), (0x00, 0xff, 8)],
            ),
        )
        full_image = rand(UInt8, 500, 300)
        view(full_image, 300:330, 200:220) .= 0xff

        fixed = full_image[(300-10):(330+10), (200-10):(220+10)]

        prefilter!(fixed)

        buffer = Array{UInt8}(undef, (size(fixed) .+ MAX_SHIFT * 2)...)

        for i = 1:1000
            expected_x = rand(-5:5)
            expected_y = rand(-5:5)

            shift, mm, mms = register!(
                mi,
                full_image,
                fixed,
                [300, 200, 330, 220] .+ padding .+
                [expected_x, expected_y, expected_x, expected_y],
                MAX_SHIFT,
                MAX_SHIFT,
                buffer;
                prefilter_frame_crop! = prefilter!,
            )

            # The shift we get out should be equal and opposite of the shift we applied
            @test shift == (-expected_x, -expected_y)
            if shift != (-expected_x, -expected_y)
                jldsave("register with filtering fail.jld2"; mi, full_image, fixed, buffer, expected_x, expected_y, shift, mm, mms)
                break
            end
        end
    end
end


text:
# Copyright 2014 Mitchell Kember. Subject to the MIT License.
# Google Code Jam Round 1A 2014
# Problem A: Charging Chaos

function nflips(l, aes, res)
	least = 9999999
	for i = 0:2^l-1
		xa = map(a -> a $ i, aes)
		if match(xa, res)
			no = nones(i)
			if no < least
				least = no
			end
		end
	end
	if least != 9999999
		return least
	end
end

function nones(n)
	c = 0
	while n != 0
		c += 1
		n &= (n - 1)
	end
	c
end

function match(l1, l2)
	for a in l1
		if !in(a, l2)
			return false
		end
	end
	true
end

function message(n)
	if n === nothing
		"NOT POSSIBLE"
	else
		n
	end
end

function parselist(l)
	xs = Array(Int, length(l))
	for i = 1:length(l)
		xs[i] = parseint(l[i], 2)
	end
	xs
end

function main()
	n = int(readline())
	for i = 1:n
		_, l = int(split(readline()))
		aes = parselist(split(readline()))
		res = parselist(split(readline()))
		a = nflips(l, aes, res)
		println("Case #$i: ", message(a))
	end
end

main()


text:
using IterativeClosestPoint
using BioStructures

ICP = IterativeClosestPoint

pdb_path = "C:\\Users\\Tom\\Downloads"
receptor1 = "Vmn1r78"
receptor2 = "Rhodopsin"

struc1 = read(joinpath(pdb_path,receptor1*".pdb"),PDB)
struc2 = read(joinpath(pdb_path,receptor2*".pdb"),PDB)

coords1 = hcat( [ struc1["A"].residues["$i"].atoms[" CA "].coords for i=1:length(struc1["A"].residues)]...)
coords2 = hcat( [ struc2["A"].residues["$i"].atoms[" CA "].coords for i=1:length(struc2["A"].residues)]...)

alignment = ICP.tm_align(coords1, coords2, )
tmscore = ICP.tm_score(coords1, coords2, alignment)


text:
struct tdvpIterable{TMPS <: AbstractMPS, TMPO <: AbstractMPO,T}
    initmps::TMPS
    mpo::TMPO
    Δt::T
end

struct tdvpState{TMPS <: AbstractMPS, TE <: AbstractTensor}
    mps::TMPS
    envs::Vector{TE}
end

function tdvpiterable(mps::TMPS, mpo::TMPO, Δt::T) where
        {TMPS <: AbstractMPS, TMPO <: AbstractMPO, T}
    cmps = normalize!(canonicalize(mps,1))
    return tdvpIterable{typeof(cmps), TMPO, T}(cmps , mpo, Δt)
end

function iterate(iter::tdvpIterable{TMPS}) where TMPS
    @unpack initmps, mpo = iter
    envs = initenvironments(initmps, mpo)
    TE = typeof(envs[1])
    state = tdvpState{TMPS,TE}(initmps, envs)
    return state, state
end

function iterate(iter::tdvpIterable, state::tdvpState)
    @unpack mps, envs = state
    @unpack mpo, Δt = iter
    tdvpsweep!(mps, mpo, envs, Δt)
    return state, state
end

leftenv(A::AbstractMPS{L,T,TE,TB}, O::MPO{L,T}, i) where
    {L,T,TE,TB}  = foldl(contractsites, zip(A[1:i], O[1:i], A[1:i]'), init=())

rightenv(A::AbstractMPS{L,T,TE,TB}, O::MPO{L,T}, i) where
    {L,T,TE,TB}  = foldr(contractsites, collect(zip(A[i:L], O[i:L], A[i:L]')), init=())

"""
    tdvp(mps::AbstractMPS, mpo::AbstractMPO, n, t, [;verbose=false, period=1, collectall=false])
returns the mps after applying the `mpo` for time `t` where `t` is split into `n` equal timesteps.
If `verbose=true`, information about the execution will be printed every `period` iteration.
If `collectall=true`, a tuple of times and states at those times will be returned.
"""
function tdvp(mps::AbstractMPS{L,T}, mpo::AbstractMPO{L,T}, n,t;
        verbose::Bool = false,
        period::Int = 1,
        collectall::Bool = false) where {L,T}
    Δt = t/n
    disp(state) = @printf("%5d \t| %.3e | %.3e\n",
                        state[2][1]-1,
                        state[1]/1e9,
                        (state[2][1]-1) * Δt)
    iter = tdvpiterable(mps, mpo, Δt)
    if collectall
        mpss = Vector{typeof(iter.initmps)}()
        iter = tee(iter, state -> push!(mpss, deepcopy(state.mps)))
    end
    iter = take(iter, n+1)
    iter = enumerate(iter)

    if verbose
        @printf("\tn \t| time (s)\t| t \n")
        iter = sample(iter, period)
        iter = stopwatch(iter)
        iter = tee(iter, disp)
        (_, (n, state)) = loop(iter)
    else
        (n, state) = loop(iter)
    end
    collectall && return ((1:n) .* Δt, normalize!.(mpss))
    return normalize!(state.mps)
end

function tdvpsweep!(mps::AbstractMPS{L}, mpo, envs, Δt) where L
    for i in 1:L-1
        onesitepropagate!(mps, mpo, i, envs, Δt/2)
        envsi = deepcopy(envs[i])
        updateenvironments!(envs, mps, mpo, i, true)
        zerositepropagate!(mps, envs[i], envsi, -Δt/2)
    end
    onesitepropagate!(mps, mpo, L, envs, Δt)
    envsi = deepcopy(envs[L-1])
    updateenvironments!(envs, mps, mpo, L)
    zerositepropagate!(mps, envsi, envs[L-1], -Δt/2)

    for i in L-1:-1:2
        onesitepropagate!(mps, mpo, i, envs, Δt/2)
        canonicalize!(mps,i-1)
        envsi = deepcopy(envs[i-1])
        updateenvironments!(envs, mps, mpo, i, false)
        zerositepropagate!(mps, envsi, envs[i-1], -Δt/2)
    end
    onesitepropagate!(mps,mpo,1, envs, Δt/2)
    canonicalize!(mps,1)
    return nothing
end

function onesitepropagate!(mps::CanonicalMPS{L}, mpo::AbstractMPO, i, envs, Δt) where L
    canonicalize!(mps, min(i, L-1))
    #absorb zero site
    centerlink(mps) ∈ (i,i-1) || error()
    if centerlink(mps) == i
        mps[i] = contractsite(mps[i], zerosite(mps))
    elseif centerlink(mps) == i-1
        mps[i] = contractsite(zerosite(mps), mps[i])
    end

    #propagate onesite
    heffmapfun = heffmapconst(mpo, envs, i)
    mps[i], info = exponentiate(
        heffmapfun,
        -1im * Δt,
        mps[i],
        ishermitian=true)

    #canonical form again
    if i == L
        U, S, mps[L] = onesitesvd2left(mps[L])
        mps.zerosite = contractsite(U, S)
    else
        mps[i], S, V = onesitesvd2right(mps[i])
        mps.zerosite = contractsite(S, V)
    end
    return nothing
end

function zerositepropagate!(mps::CanonicalMPS{L}, hl, hr, Δt) where {L}
    mps.zerosite, info = exponentiate(
        heffmap((hl,hr)),
        -1im * Δt,
        zerosite(mps),
        ishermitian=true)
    return nothing
end


text:
"""
    mveltsplot(setting; alpha = 0.05, showplot = true)

Generate MVE - LTS plot for visual detecting of regression outliers.

# Arguments
- `setting::RegressionSetting`: A regression setting object.
- `alpha::Float64`: Probability for quantiles of Chi-Squared statistic.
- `showplot::Bool`: Whether a plot is shown or only return statistics.

# Description
This is a method of combination of `lts` and `mve`. Regression residuals and robust distances obtained
by `mve` and `mve` are used to generate a plot. Despite this is a visual method, drawing a plot is not really
necessary. The algorithm divides the residuals-distances space into 4 parts, one for clean observations,
one for vertical outliers (y-space outliers), one for bad-leverage points (x-space outliers), and one for
good leverage points (observations far from the remaining of data in both x and y space).

# Output
- `["plot"]`: Generated plot object
- `["robust.distances"]`: Robust Mahalanobis distances
- `["scaled.residuals"]`: Scaled residuals of an `lts` estimate
- `["chi.squared"]`: Quantile of Chi-Squared distribution
- `["regular.points"]`: Array of indices of clean observations
- `["outlier.points"]`: Array of indices of y-space outliers (vertical outliers)
- `["leverage.points"]`: Array of indices of x-space outliers (bad leverage points)
- `["outlier.and.leverage.points"]`: Array of indices of xy-space outliers (good leverage points)

# References
Van Aelst, Stefan, and Peter Rousseeuw. "Minimum volume ellipsoid." Wiley
Interdisciplinary Reviews: Computational Statistics 1.1 (2009): 71-82.
"""
function mveltsplot(setting::RegressionSetting; alpha=0.05, showplot=true)
    ltsresult = lts(setting)
    mveresult = mve(setting.data)

    n, p = size(setting.data)
    indices = collect(1:n)

    chidist = Chisq(p)
    chicrit = sqrt(quantile(chidist, 1.0 - alpha))

    scaledresiduals = ltsresult["scaled.residuals"]
    robdistances = sqrt.(mveresult["squared.mahalanobis"])

    regularpoints = filter(i -> abs(scaledresiduals[i]) < 2.5 && robdistances[i] < chicrit, indices)
    outliers = filter(i -> scaledresiduals[i] > 2.5 && robdistances[i] < chicrit, indices)
    leverage = filter(i -> abs(scaledresiduals[i]) < 2.5 && robdistances[i] > chicrit, indices)
    outlierandleveragepoints = filter(i -> abs(scaledresiduals[i]) > 2.5 && robdistances[i] > chicrit, indices)

    scplot = nothing
    if showplot
        scplot = scatter(robdistances,
            scaledresiduals,
            legend=false,
            series_annotations=text.(1:n, :bottom),
            tickfont=font(10), guidefont=font(10), labelfont=font(10)
            )
        title!("MVE - LTS Plot")
        xlabel!("Robust distances")
        ylabel!("Scaled residuals")
        hline!([-2.5, 2.5])
        vline!([chicrit])
    end

    result = Dict()
    result["plot"] = scplot
    result["robust.distances"] = robdistances
    result["scaled.residuals"] = scaledresiduals
    result["chi.squared"] = chicrit^2
    result["regular.points"] = regularpoints
    result["outlier.points"] = outliers
    result["leverage.points"] = leverage
    result["outlier.and.leverage.points"] = outlierandleveragepoints
    return result
end

max_stars_repo_name: jbytecode/LinRegOutliers

max_issues_repo_name: jbytecode/LinRegOutliers

max_forks_repo_name: jbytecode/LinRegOutliers


text:
using BitBlit2D
using Test

@testset "BitBlit2D.jl" begin
    include("image-tests.jl")
end


text:
using LinearAlgebra: Adjoint, Transpose, UnitLowerTriangular, LowerTriangular, UnitUpperTriangular, UpperTriangular, AbstractTriangular

istransposed(x) = 'N'
istransposed(x::Adjoint) = 'C'
istransposed(x::Adjoint{<:Real}) = 'T'
istransposed(x::Transpose) = 'T'

untransposed(A) = A
untransposed(A::Adjoint) = adjoint(A)
untransposed(A::Transpose) = transpose(A)

_transpose(A) = A'
_transpose(A::Adjoint) = adjoint(A)
_transpose(A::Transpose) = transpose(A)

uplochar(::Union{LowerTriangular,UnitLowerTriangular}) = 'L'
uplochar(::Union{UpperTriangular,UnitUpperTriangular}) = 'U'
diagchar(::Union{LowerTriangular,UpperTriangular}) = 'N'
diagchar(::Union{UnitLowerTriangular,UnitUpperTriangular}) = 'U'
untranspose_flag(A::Adjoint, f::Bool = false) = untranspose_flag(parent(A), !f)
untranspose_flag(A::Transpose, f::Bool = false) = untranspose_flag(parent(A), !f)
untranspose_flag(A, f::Bool = false) = (A, f)

struct LU{T,M<:StridedMatrix{T},I}
  factors::M
  ipiv::Vector{I}
  info::I
end
function getrf_ipiv(A::StridedMatrix, ::Val{T}) where {T}
  M,N = size(A)
  Vector{T}(undef, min(M,N))
end

function _ipiv_rows!(A::LU, order::OrdinalRange, B::StridedVecOrMat)
  @inbounds for i ∈ order
    i ≠ A.ipiv[i] && LinearAlgebra._swap_rows!(B, i, A.ipiv[i])
  end
  B
end
_apply_ipiv_rows!(A::LU, B::StridedVecOrMat) = _ipiv_rows!(A, 1 : length(A.ipiv), B)
function _ipiv_cols!(A::LU, order::OrdinalRange, B::StridedVecOrMat)
  ipiv = A.ipiv
  @inbounds for i ∈ order
    i ≠ ipiv[i] && LinearAlgebra._swap_cols!(B, i, ipiv[i])
  end
  B
end
_apply_inverse_ipiv_cols!(A::LU, B::StridedVecOrMat) = _ipiv_cols!(A, length(A.ipiv) : -1 : 1, B)

for (name,BlasInt,suff) ∈ [
  ("mkl", :Int32, ""),
  ("openblas", :Int64, "_64_"),
  ("blis", :Int64, "_64_")
]
  uname = uppercase(name)
  lib = Symbol("lib", uname)
  fgemm = Symbol("gemm", name, '!')
  fgetrf = Symbol("getrf", name, '!')
  ftrsm = Symbol("trsm", name, '!')
  frdiv = Symbol("rdiv", name)
  fldiv = Symbol("ldiv", name)
  frdivbang = Symbol("rdiv", name, '!')
  fldivbang = Symbol("ldiv", name, '!')
  flu = Symbol("lu", name)
  flubang = Symbol("lu", name, '!')
  for (T,prefix) ∈ [(:Float32,'s'),(:Float64,'d')]
    fmgemm = QuoteNode(Symbol(prefix, "gemm", suff))
    fmgetrf = QuoteNode(Symbol(prefix, "getrf", suff))
    fmtrsm = QuoteNode(Symbol(prefix, "trsm", suff))
    @eval begin
      function $fgemm(
        C::AbstractMatrix{$T}, A::AbstractMatrix{$T}, B::AbstractMatrix{$T}, α = one($T), β = zero($T)
      )
        istransposed(C) === 'N' || ($fgemm(untransposed(C), _transpose(B), _transpose(A)); return C)
        transA = istransposed(A)
        transB = istransposed(B)
        pA     = untransposed(A);
        pB     = untransposed(B)
        M, N = size(C); K = size(B, 1)
        ldA = stride(pA, 2)
        ldB = stride(pB, 2)
        ldC = stride(C, 2)
        ccall(
          ($fmgemm, $lib), Cvoid,
          (Ref{UInt8}, Ref{UInt8}, Ref{$BlasInt}, Ref{$BlasInt}, Ref{$BlasInt}, Ref{$T}, Ref{$T},
           Ref{$BlasInt}, Ref{$T}, Ref{$BlasInt}, Ref{$T}, Ref{$T}, Ref{$BlasInt}),
          transA, transB, M, N, K, α, pA, ldA, pB, ldB, β, C, ldC
        )
        C
      end
    end
    name == "blis" && continue
    @eval begin
      function $fgetrf(
        A::StridedMatrix{$T}, ipiv::StridedVector{$BlasInt} = getrf_ipiv(A, Val($BlasInt))
      )
        M, N = size(A)
        info = Ref{$BlasInt}()
        ccall(
          ($fmgetrf,$lib), Cvoid, (Ref{$BlasInt},Ref{$BlasInt},Ptr{$T},Ref{$BlasInt},Ptr{$BlasInt},Ref{$BlasInt}),
          M, N, A, max(1,stride(A,2)), ipiv, info
        )
        LU(A, ipiv, info[])
      end
      function $ftrsm(
        B::StridedMatrix{$T}, α::$T, A::AbstractMatrix{$T}, side::Char
      )
        M, N = size(A)
        pA, transa = untranspose_flag(A)
        uplo = uplochar(pA)
        diag = diagchar(pA)
        ppA = parent(pA)
        ccall(
          ($fmtrsm, $lib), Cvoid, (Ref{UInt8}, Ref{UInt8}, Ref{UInt8}, Ref{UInt8}, Ref{$BlasInt}, Ref{$BlasInt}, Ref{$T}, Ptr{Float64}, Ref{$BlasInt}, Ptr{$T}, Ref{$BlasInt}),
          side, uplo, ifelse(transa, 'T', 'N'), diag, M % $BlasInt, N % $BlasInt, α, ppA, max(1,stride(ppA,2)), B, max(1,stride(B,2))
        )
        B
      end
    end
  end
  name == "blis" && continue
  @eval begin
    function $frdivbang(
      A::StridedMatrix{T}, B::AbstractTriangular{T}
    ) where {T <: Union{Float32,Float64}}
      $ftrsm(A, one(T), B, 'R')
    end
    function $fldivbang(
      A::StridedMatrix{T}, B::AbstractTriangular{T}
    ) where {T <: Union{Float32,Float64}}
      $ftrsm(A, one(T), B, 'L')
    end
    $flubang(A::StridedMatrix) = $fgetrf(A, getrf_ipiv(A, Val($BlasInt)))
    $flu(A::StridedMatrix) = $flubang(copy(A))
    function $frdivbang(A::AbstractMatrix, B::LU{<:Any,<:AbstractMatrix})
      $frdivbang($frdivbang(A, UpperTriangular(B.factors)), UnitLowerTriangular(B.factors))
      _apply_inverse_ipiv_cols!(B, A) # mutates `A`
    end
    function $fldivbang(A::LU{<:Any,<:AbstractMatrix}, B::AbstractMatrix)
      _apply_ipiv_rows!(A, B)
      $fldivbang($fldivbang(B, UnitLowerTriangular(A.factors)), UpperTriangular(A.factors))
    end

    $frdivbang(A::AbstractMatrix, B::AbstractMatrix) = $frdivbang(A, $flubang(B))
    $fldivbang(A::AbstractMatrix, B::AbstractMatrix) = $fldivbang($flubang(A), B)

    $frdiv(A::AbstractMatrix, B) = $frdivbang(copy(A), copy(B))
    $fldiv(A::AbstractMatrix, B) = $fldivbang(copy(A), copy(B))

  end
end
let BlasInt = :Int32
  for (T,prefix) ∈ [(:Float32,'s'),(:Float64,'d')]
    f = Symbol(prefix, "gemm_direct")
    @eval begin
      @inline function gemmmkl_direct!(C::AbstractMatrix{$T}, A::AbstractMatrix{$T}, B::AbstractMatrix{$T}, α = one($T), β = zero($T))
        istransposed(C) === 'N' || ($f(untransposed(C), _transpose(B), _transpose(A)); return C)
        transA = istransposed(A)
        transB = istransposed(B)
        pA     = untransposed(A);
        pB     = untransposed(B)
        M, N = size(C); K = size(B, 1)
        ldA = stride(pA, 2)
        ldB = stride(pB, 2)
        ldC = stride(C, 2)
        ccall(
          ($(QuoteNode(f)), libMKL), Cvoid,
          (Ref{UInt8}, Ref{UInt8}, Ref{$BlasInt}, Ref{$BlasInt}, Ref{$BlasInt}, Ref{$T}, Ref{$T},
           Ref{$BlasInt}, Ref{$T}, Ref{$BlasInt}, Ref{$T}, Ref{$T}, Ref{$BlasInt}, Ref{$BlasInt}),
          transA, transB, M, N, K, α, pA, ldA, pB, ldB, β, C, ldC, 0
        )
        C
      end
    end
  end
end


text:
"""
```julia
_posterior(
    prior::Distributions.DiscreteNonParametric,
    responses::Vector{<:AbstractResponse},
    items::Vector{<:AbstractResponse},
    )
```

#Description

It computes the posterior (not normalized) for a discrete non parametric distribution prior (`Distributions.DiscreteNonParametric`), vector of items `items` with vector of answers `responses`.
Items and responses must be sorted by idx (items) and item_idx (responses)
"""
function _posterior(
    prior::Distributions.DiscreteNonParametric,
    items::Vector{<:AbstractItem},
    responses::Vector{<:AbstractResponse},
    )
    return map( (x, w) ->
            mapreduce( (i, r) ->
            _likelihood(r, x, i)
            ,
            *,
            items,
            responses,
            )*w
        ,
        prior.support,
        prior.p
        )
end

function _update_posterior!(
    latent::Latent1D,
    items::Vector{<:AbstractItem}, #! must be sorted by idx
    responses::Vector{<:AbstractResponse}, #! must be sorted by item_idx
    )
    likelihood = _posterior(latent.prior, items, responses)
    normalizer = sum(likelihood)
    if normalizer > typemin(Float64)
        latent.posterior = Distributions.DiscreteNonParametric(latent.prior.support, likelihood ./ normalizer; check_args = false);
        latent.likelihood = normalizer
    else
        latent.posterior = Distributions.DiscreteNonParametric(latent.prior.support, likelihood; check_args = false);
    end
end

function update_posterior!(
    examinee::AbstractExaminee,
    items::Vector{<:AbstractItem},
    responses::Vector{<:AbstractResponse};
    already_sorted::Bool = false,
    kwargs...
    )
    if !already_sorted
        responses_filtered = filter( r -> r.examinee_idx == examinee.idx, responses)
        items_filtered = items[map(r -> r.item_idx, responses_filtered)]
    else
        responses_filtered = responses
        items_filtered = items
    end
    _update_posterior!(examinee.latent, items_filtered, responses_filtered)
    return nothing
end

function update_posterior!(
    examinees::Vector{<:AbstractExaminee},
    items::Vector{<:AbstractItem},
    responses::Vector{<:AbstractResponse};
    already_sorted::Bool = false,
    kwargs...
)
    map( e -> update_posterior!(e, items, responses; already_sorted = already_sorted, kwargs...), examinees)
    return nothing
end


text:
using SynapseClient

e = SynapseClient.Entity(name="Test object", description="I hope this works",
           annotations = Dict("foo"=>123, "nerds"=>["chris","jen","janey"], "annotations"=>"How confusing!"),
           properties  = Dict("annotations"=>"/repo/v1/entity/syn1234/annotations",
                              "md5"=>"cdef636522577fc8fb2de4d95875b27c",
                              "parentId"=>"syn1234"),
           concreteType="org.sagebionetworks.repo.model.Data")

syn = SynapseClient.Synapse(debug=true, skip_checks=true)
printentity(syn,e)

syn = SynapseClient.login()

file = get(syn, "syn1906479")
printentity(syn,file)
provenance = getprovenance(syn, "syn1906479")

project = get(syn, "syn1901847")
wiki = getwiki(syn, "syn1901847")
wiki2 = getwiki(syn, project)

@assert wiki == wiki2

team = getteam(syn,3341340)

evaluation = getevaluation(syn,2005090)

SynapseClient.check_for_updates()


text:
module GtkExtensions

# nice exports !
export show_iter,
    get_iter_at_position, text_view_window_to_buffer_coords, get_current_page_idx,
    set_current_page_idx, get_tab, set_position!, text_buffer_copy_clipboard, set_tab_label_text,
    MutableGtkTextIter, GtkTextIters,GtkEventBox, GtkCssProviderFromData!, GtkIconThemeAddResourcePath,
    GtkIconThemeGetDefault, index, style_css, text, PROPAGATE, INTERRUPT,
    end_iter, line_count,
    cursor_locations, gdk_window_get_origin, g_timeout_add, g_idle_add, delete_text, insert_text, insert,
    response, GdkKeySyms, offset, mutable, nonmutable, text_view_buffer_to_window_coords,
    line, show, scroll_to_iter, css_provider,
    push!, GtkCssProviderLeaf, GtkCssProvider, getbuffer, iter_nth_child, GtkIconThemeLoadIconForScale,
    expand_root, model, set_cursor_on_cell, expand, treepath, foreach, selection, selected, select_value

using Gtk

import ..Gtk: suffix, GtkCssProviderLeaf, GtkCssProvider
import Gtk: GtkTextIter, libgtk, iter_nth_child, selected
import Base: foreach, push!

const PROPAGATE = convert(Cint,false)
const INTERRUPT = convert(Cint,true)

import Gtk.GConstants: GdkModifierType
import Gtk.GdkKeySyms

include("MenuUtils.jl")

## TextIters

const MutableGtkTextIter = Gtk.GLib.MutableTypes.Mutable{GtkTextIter}
const GtkTextIters = Union{MutableGtkTextIter,GtkTextIter}
mutable(it::GtkTextIter) = Gtk.GLib.MutableTypes.mutable(it)

offset(it::GtkTextIters) = get_gtk_property(it,:offset,Integer)
line(it::GtkTextIters) = get_gtk_property(it,:line,Integer)+1#Gtk counts from zero
nonmutable(buffer::GtkTextBuffer,it::MutableGtkTextIter) = GtkTextIter(buffer,offset(it)+1)#this allows to convert to GtkTextBuffer without the -1 definition in Gtk.jl

getbuffer(it::GtkTextIter) = convert(GtkTextBuffer,
    ccall((:gtk_text_iter_get_buffer, libgtk),Ptr{GtkTextBuffer},(Ref{GtkTextIter},),it)
)
getbuffer(it::MutableGtkTextIter) = convert(GtkTextBuffer,
    ccall((:gtk_text_iter_get_buffer, libgtk),Ptr{GtkTextBuffer},(Ptr{GtkTextIter},),it)
)

function show_iter(it::MutableGtkTextIter,buffer::GtkTextBuffer,color::Int)
    Gtk.apply_tag(buffer, color > 0 ? "debug1" : "debug2",it, it+1)
end

function end_iter(buffer::Gtk.GtkTextBuffer)
    iter = Gtk.mutable(GtkTextIter)
    ccall((:gtk_text_buffer_get_end_iter,libgtk),Cvoid,(Ptr{Gtk.GObject},Ptr{GtkTextIter}),buffer,iter)
    return iter
end

line_count(buffer::GtkTextBuffer) = ccall((:gtk_text_buffer_get_line_count,libgtk),Cint,(Ptr{GObject},),buffer)

## TextView

get_iter_at_position(text_view::Gtk.GtkTextView,iter::MutableGtkTextIter,trailing,x::Int32,y::Int32) = ccall((:gtk_text_view_get_iter_at_position,libgtk),Cvoid,
	(Ptr{Gtk.GObject},Ptr{GtkTextIter},Ptr{Cint},Cint,Cint),text_view,iter,trailing,x,y)

function get_iter_at_position(text_view::Gtk.GtkTextView,x::Integer,y::Integer)
    buffer = get_gtk_property(text_view,:buffer,GtkTextBuffer)
    iter = mutable(GtkTextIter(buffer))
    get_iter_at_position(text_view::Gtk.GtkTextView,iter,C_NULL,Int32(x),Int32(y))
    return nonmutable(buffer,iter)
end

function text_view_window_to_buffer_coords(text_view::Gtk.GtkTextView,wintype::Integer,window_x::Integer,window_y::Integer)

	buffer_x = Gtk.mutable(Cint)
	buffer_y = Gtk.mutable(Cint)

	ccall((:gtk_text_view_window_to_buffer_coords,libgtk),Cvoid,
		(Ptr{Gtk.GObject},Cint,Cint,Cint,Ptr{Cint},Ptr{Cint}),text_view,Int32(wintype),window_x,window_y,buffer_x,buffer_y)

	return (buffer_x[],buffer_y[])
end
text_view_window_to_buffer_coords(text_view::Gtk.GtkTextView,window_x::Integer,window_y::Integer) = text_view_window_to_buffer_coords(text_view,2,window_x,window_y)

function text_view_buffer_to_window_coords(text_view::Gtk.GtkTextView,wintype::Integer,buffer_x::Integer,buffer_y::Integer)

	window_x = Gtk.mutable(Cint)
	window_y = Gtk.mutable(Cint)

	ccall((:gtk_text_view_buffer_to_window_coords,libgtk),Cvoid,
		(Ptr{Gtk.GObject},Cint,Cint,Cint,Ptr{Cint},Ptr{Cint}),text_view,Int32(wintype),buffer_x,buffer_y,window_x,window_y)

	return (window_x[],window_y[])
end
text_view_buffer_to_window_coords(text_view::Gtk.GtkTextView,buffer_x::Integer,buffer_y::Integer) = text_view_window_to_buffer_coords(text_view,0,buffer_x,buffer_y)

function cursor_locations(text_view::Gtk.GtkTextView)
    weak = Gtk.mutable(Gtk.GdkRectangle)
    strong = Gtk.mutable(Gtk.GdkRectangle)
    buffer = get_gtk_property(text_view,:buffer,GtkTextBuffer)
    iter = mutable( GtkTextIter(buffer, get_gtk_property(buffer,:cursor_position,Int)) )

    ccall((:gtk_text_view_get_cursor_locations,libgtk),Cvoid,(Ptr{Gtk.GObject},Ptr{GtkTextIter},Ptr{Gtk.GdkRectangle},Ptr{Gtk.GdkRectangle}),text_view,iter,strong,weak)
    return (iter,strong[],weak[])
end

scroll_to_iter(text_view::Gtk.GtkTextView,iter::GtkTextIter,within_margin::Number,use_align::Bool,xalign::Number,yalign::Number) = ccall((:gtk_text_view_scroll_to_iter,libgtk),Cint,
	(Ptr{Gtk.GObject},Ref{GtkTextIter},Cdouble,Cint,Cdouble,Cdouble),
    text_view,iter,within_margin,use_align,xalign,yalign)

scroll_to_iter(text_view::Gtk.GtkTextView,iter::MutableGtkTextIter,within_margin::Number,use_align::Bool,xalign::Number,yalign::Number) = ccall((:gtk_text_view_scroll_to_iter,libgtk),Cint,
	(Ptr{Gtk.GObject},Ptr{GtkTextIter},Cdouble,Cint,Cdouble,Cdouble),
    text_view,iter,within_margin,use_align,xalign,yalign)

scroll_to_iter(text_view::Gtk.GtkTextView,iter::GtkTextIters) = scroll_to_iter(text_view,iter,0.0,true,1.0,0.1)

# notebook
get_current_page_idx(notebook::Gtk.GtkNotebook) = ccall((:gtk_notebook_get_current_page,libgtk),Cint,
    (Ptr{Gtk.GObject},),notebook)+1
set_current_page_idx(notebook::Gtk.GtkNotebook,page_num::Int) = ccall((:gtk_notebook_set_current_page,libgtk),Cvoid,
    (Ptr{Gtk.GObject},Cint),notebook,page_num-1)

index(notebook::GtkNotebook) = get_current_page_idx(notebook)
index(notebook::GtkNotebook,i::Integer) = set_current_page_idx(notebook,i)
index(notebook::GtkNotebook, child::Gtk.GtkWidget) = pagenumber(notebook, child)+1

get_tab(notebook::Gtk.GtkNotebook,page_num::Int) = convert(Gtk.GtkWidget,ccall((:gtk_notebook_get_nth_page,libgtk),Ptr{Gtk.GObject},
	(Ptr{Gtk.GObject},Cint),notebook,page_num-1))

set_tab_label_text(notebook::Gtk.GtkNotebook,child,tab_text) = ccall((:gtk_notebook_set_tab_label_text,Gtk.libgtk),Cvoid,(Ptr{Gtk.GObject},
Ptr{Gtk.GObject},Ptr{UInt8}),notebook,child,tab_text)

popup_disble(notebook::Gtk.GtkNotebook) = ccall((:gtk_notebook_popup_disable ,Gtk.libgtk),
      Cvoid,
      (Ptr{Gtk.GObject},),
      notebook)
import Base.insert!
function insert!(w::Gtk.GtkNotebook, position::Integer, x::Union{Gtk.GtkWidget,Gtk.AbstractStringLike}, label::Union{Gtk.GtkWidget,Gtk.AbstractStringLike}, menu::Gtk.GtkWidget)
    ccall((:gtk_notebook_insert_page_menu,libgtk), Cint,
        (Ptr{GObject}, Ptr{Gtk.GObject}, Ptr{Gtk.GObject},Ptr{Gtk.GObject}, Cint),
        w, x, label, menu,position-1)+1
    w

end

## entry

function set_position!(editable::Gtk.Entry,position_)
    ccall((:gtk_editable_set_position,libgtk),Cvoid,(Ptr{Gtk.GObject},Cint),editable,position_)
end

#####  GtkClipboard #####

# Gtk.@gtktype GtkClipboard

# baremodule GdkAtoms
#     const NONE = 0x0000
#     const SELECTION_PRIMARY = 0x0001
#     const SELECTION_SECONDARY = 0x0002
#     const SELECTION_TYPE_ATOM = 0x0004
#     const SELECTION_TYPE_BITMAP = 0x0005
#     const SELECTION_TYPE_COLORMAP = 0x0007
#     const SELECTION_TYPE_DRAWABLE = 0x0011
#     const SELECTION_TYPE_INTEGER = 0x0013
#     const SELECTION_TYPE_PIXMAP = 0x0014
#     const SELECTION_TYPE_STRING = 0x001f
#     const SELECTION_TYPE_WINDOW = 0x0021
#     const SELECTION_CLIPBOARD = 0x0045
# end

# GtkClipboardLeaf(selection::UInt16) =  GtkClipboardLeaf(ccall((:gtk_clipboard_get,libgtk), Ptr{GObject},
#     (UInt16,), selection))
# GtkClipboardLeaf() = GtkClipboardLeaf(GdkAtoms.SELECTION_CLIPBOARD)
# clipboard_set_text(clip::GtkClipboard,text::AbstractString) = ccall((:gtk_clipboard_set_text,libgtk), Cvoid,
#     (Ptr{GObject}, Ptr{UInt8},Cint), clip, text, sizeof(text))
# clipboard_store(clip::GtkClipboard) = ccall((:gtk_clipboard_store,libgtk), Cvoid,
#     (Ptr{GObject},), clip)

# #note: this needs main_loops to run
# function clipboard_wait_for_text(clip::GtkClipboard)
#     ptr = ccall((:gtk_clipboard_wait_for_text,libgtk), Ptr{UInt8},
#         (Ptr{GObject},), clip)
#     return ptr == C_NULL ? "" : unsafe_string(ptr)
# end

# text_buffer_copy_clipboard(buffer::GtkTextBuffer,clip::GtkClipboard)  = ccall((:gtk_text_buffer_copy_clipboard, libgtk),Cvoid,
#     (Ptr{GObject},Ptr{GObject}),buffer,clip)

##
function GtkCssProviderFromData!(provider::GtkCssProvider;data=nothing,filename=nothing)
    source_count = (data!==nothing) + (filename!==nothing)
    @assert(source_count <= 1,
        "GtkCssProvider must have at most one data or filename argument")

    if data !== nothing
        Gtk.GError() do error_check

          ccall((:gtk_css_provider_load_from_data,libgtk), Bool,
            (Ptr{Gtk.GObject}, Ptr{UInt8}, Clong, Ptr{Ptr{Gtk.GError}}),
            provider, string(data), sizeof(data), error_check)
        end
    elseif filename !== nothing
        Gtk.GError() do error_check
          ccall((:gtk_css_provider_load_from_path,libgtk), Bool,
            (Ptr{Gtk.GObject}, Ptr{UInt8}, Ptr{Ptr{Gtk.GError}}),
            provider, string(filename), error_check)
        end
    end
    return provider
end

#TODO fix with Gtk.jl
GtkCssProviderLeaf() = GtkCssProviderLeaf(ccall((:gtk_css_provider_new,libgtk),Ptr{Gtk.GObject},()))

push!(context::GtkStyleContext, provider::GtkCssProvider, priority::Integer) =
  ccall((:gtk_style_context_add_provider, libgtk), Cvoid, (Ptr{GObject}, Ptr{GObject}, Cuint),
		 context, provider, priority)

function style_css(w::Gtk.GtkWidget,css::String)
    sc = Gtk.G_.style_context(w)
    provider = GtkCssProvider()
    push!(sc, GtkCssProviderFromData!(provider,data=css), 600)
end

function style_css(w::Gtk.GtkWidget,provider::GtkCssProvider)
    sc = Gtk.G_.style_context(w)
    push!(sc, provider, 600)
end

## Gdk

function gdk_window_get_origin(window)

	window_x = Gtk.mutable(Cint)
	window_y = Gtk.mutable(Cint)

	ccall((:gdk_window_get_origin,Gtk.libgdk),Cint,
		(Ptr{Gtk.GObject},Ptr{Cint},Ptr{Cint}),window,window_x,window_y)

	return (window_x[],window_y[])
end

gdk_keyval_name(val) = unsafe_string(
    ccall((:gdk_keyval_name,libgtk),Ptr{UInt8},(Cuint,),val),
true)

GtkIconThemeGetDefault() =  ccall((:gtk_icon_theme_get_default,Gtk.libgtk),Ptr{GObject},())

GtkIconThemeAddResourcePath(iconTheme,path::AbstractString) =  ccall(
                                                                   (:gtk_icon_theme_append_search_path,Gtk.libgtk),Cvoid,
                                                                   (Ptr{GObject}, Ptr{UInt8}),
                                                                   iconTheme,
                                                                   path);
function GtkIconThemeLoadIconForScale(iconTheme,icon_name::AbstractString, size::Integer, scale::Integer, flags::Integer)
    local pixbuf::Ptr{GObject}
    Gtk.GError() do error_check
        pixbuf = ccall((:gtk_icon_theme_load_icon_for_scale,Gtk.libgtk),
                   Ptr{GObject},
                   (Ptr{GObject},Ptr{UInt8},Cint,Cint,Cint,Ptr{Ptr{Gtk.GError}}),
                   iconTheme,string(icon_name),size,scale,flags,error_check)

        return pixbuf !== C_NULL
    end
    return convert(GdkPixbuf,pixbuf)
end

##GtkTreeStore
function insert(store::GtkTreeStore, it::GtkTreeIter, parent::GtkTreeIter, pos::Int)
    ccall((:gtk_tree_store_insert, Gtk.libgtk), Cvoid,  (Ptr{Gtk.GObject},Ptr{Gtk.GtkTreeIter},Ptr{Gtk.GtkTreeIter},Cint),
            store,it,parent,pos)
end

## update iter pointing to nth child n in 1:nchildren)
## return boolean
function iter_nth_child(treeModel::Gtk.GtkTreeModel, iter::Gtk.Mutable{Gtk.GtkTreeIter}, piter, n::Int)
  if (piter==nothing)
    ret = ccall((:gtk_tree_model_iter_nth_child, Gtk.libgtk), Cint,
        (Ptr{Gtk.GObject}, Ptr{GtkTreeIter}, Ptr{Gtk.GtkTreeIter}, Cint),
        treeModel, iter, C_NULL, n - 1) # 0-based
  else
    ret = ccall((:gtk_tree_model_iter_nth_child, Gtk.libgtk), Cint,
        (Ptr{Gtk.GObject}, Ptr{Gtk.GtkTreeIter}, Ptr{Gtk.GtkTreeIter}, Cint),
        treeModel, iter, Gtk.mutable(piter), n - 1) # 0-based
  end
    ret != 0
end

#GtkTreeView

function model(tree_view::Gtk.GtkTreeView)
    return convert(Gtk.GtkTreeStore,
        ccall((:gtk_tree_view_get_model, Gtk.libgtk),
        Ptr{Gtk.GObject},
        (Ptr{Gtk.GObject},),
        tree_view))
end

function selection(tree_view::Gtk.GtkTreeView)
    return convert(Gtk.GtkTreeSelection,
        ccall((:gtk_tree_view_get_selection,Gtk.libgtk),
        Ptr{Gtk.GObject},
        (Ptr{Gtk.GObject},),
        tree_view))
end

function set_cursor_on_cell(tree_view::Gtk.GtkTreeView, path::Gtk.GtkTreePath)
    return  ccall((:gtk_tree_view_set_cursor_on_cell , Gtk.libgtk),
                   Cvoid,
                  (Ptr{Gtk.GObject},Ptr{Gtk.GtkTreePath},Ptr{Gtk.GObject},Ptr{Gtk.GObject},Cint),
                  tree_view,path,C_NULL,C_NULL,false)
end

import Base.expand
function expand(tree_view::GtkTreeView,path::GtkTreePath)
    return  ccall((:gtk_tree_view_expand_to_path,libgtk),Cvoid,
                  (Ptr{Gtk.GObject},Ptr{Gtk.GtkTreePath}),
                  tree_view,path)
end

function treepath(path::AbstractString)
    ptr = ccall((:gtk_tree_path_new_from_string,libgtk),Ptr{GtkTreePath},
                  (Ptr{UInt8},),
                  string(path))
    if ptr != C_NULL
        return convert(GtkTreePath,ptr)
    else
        return GtkTreePath()
    end
end
expand_root(tree_view::GtkTreeView) = expand(tree_view,treepath("0"))

function selected(tree_view::GtkTreeView,list::GtkTreeStore)
    selmodel = selection(tree_view)
    if hasselection(selmodel)
        iter = selected(selmodel)
        return list[iter]
    end
    return nothing
end

#select the first entry that is equal to v
function select_value(tree_view::GtkTreeView,list::GtkTreeStore,v)
    selmodel = Gtk.G_.selection(tree_view)
    for i = 1:length(list)
        if list[i] == v
            partialsort!(selmodel, Gtk.iter_from_index(list, i))
            return
        end
    end
end

#GtkTreeModel
function foreach(model::Gtk.GtkTreeModel, f::Function, data)
  foreach_function = @cfunction($f, Cint, (Ptr{Gtk.GObject}, Ptr{Gtk.GtkTreePath}, Ptr{Gtk.GtkTreeIter}, Ptr{Cvoid}))
   ccall((:gtk_tree_model_foreach, Gtk.libgtk),      Cvoid,
                (Ptr{Gtk.GObject},Ptr{Cvoid}, Ptr{Cvoid}),
                model,foreach_function,pointer_from_objref(data))
end

global default_css_provider = convert(Ptr{Gtk.GObject}, 0)
function __init__()
    global default_css_provider = GtkCssProviderLeaf(
        ccall((:gtk_css_provider_get_default,libgtk),Ptr{Gtk.GObject},())
    )

    #import Base.show
    #show(io::IO, it::GtkTextIter) = println("GtkTextIter($(offset(it)))")

end

end#module

max_stars_repo_name: jonathanBieler/GtkExtensions

max_issues_repo_name: jonathanBieler/GtkExtensions

max_forks_repo_name: jonathanBieler/GtkExtensions


text:
"""
    GridData <: AbstractArray

Simulation data specific to a single grid.
"""
abstract type GridData{Y,X,R,T,N,I} <: AbstractArray{T,N} end

function (::Type{G})(d::GridData{Y,X,R,T,N}) where {G<:GridData,Y,X,R,T,N}
    args = init(d), mask(d), proc(d), opt(d), boundary(d), padval(d),
        source(d), dest(d), sourcestatus(d), deststatus(d)
    G{Y,X,R,T,N,map(typeof, args)...}(args...)
end
function ConstructionBase.constructorof(::Type{T}) where T<:GridData{Y,X,R} where {Y,X,R}
    T.name.wrapper{Y,X,R}
end

GridDataOrReps = Union{GridData, Vector{<:GridData}}

# Array interface
Base.size(d::GridData) = size(source(d))
Base.axes(d::GridData) = axes(source(d))
Base.eltype(d::GridData) = eltype(source(d))
Base.firstindex(d::GridData) = firstindex(source(d))
Base.lastindex(d::GridData) = lastindex(source(d))

# Getters
init(d::GridData) = d.init
mask(d::GridData) = d.mask
radius(d::GridData{<:Any,<:Any,R}) where R = R
proc(d::GridData) = d.proc
opt(d::GridData) = d.opt
boundary(d::GridData) = d.boundary
padval(d::GridData) = d.padval
source(d::GridData) = d.source
dest(d::GridData) = d.dest
sourcestatus(d::GridData) = d.sourcestatus
deststatus(d::GridData) = d.deststatus

gridsize(d::GridData) = size(init(d))
gridsize(A::AbstractArray) = size(A)
gridsize(nt::NamedTuple) = gridsize(first(nt))
gridsize(nt::NamedTuple{(),Tuple{}}) = 0, 0
gridaxes(d::GridData) = map(Base.OneTo, gridsize(d))
gridview(d::GridData) = view(parent(dest(d)), map(a -> a .+ radius(d), gridaxes(d))...)

"""
    ReadableGridData <: GridData

    ReadableGridData(grid::GridData)
    ReadableGridData{Y,X,R}(init::AbstractArray, mask, opt, boundary, padval)

Simulation data and storage passed to rules for each timestep.
"""
struct ReadableGridData{
    Y,X,R,T,N,I<:AbstractArray{T,N},M,P<:Processor,Op<:PerformanceOpt,Ov,PV,S,D,SSt,DSt
} <: GridData{Y,X,R,T,N,I}
    init::I
    mask::M
    proc::P
    opt::Op
    boundary::Ov
    padval::PV
    source::S
    dest::D
    sourcestatus::SSt
    deststatus::DSt
end
function ReadableGridData{Y,X,R}(
    init::I, mask::M, proc::P, opt::Op, boundary::Ov, padval::PV, source::S,
    dest::D, sourcestatus::SSt, deststatus::DSt
) where {Y,X,R,I<:AbstractArray{T,N},M,P,Op,Ov,PV,S,D,SSt,DSt} where {T,N}
    ReadableGridData{Y,X,R,T,N,I,M,P,Op,Ov,PV,S,D,SSt,DSt}(
        init, mask, proc, opt, boundary, padval, source, dest, sourcestatus, deststatus
    )
end
# Generate simulation data to match a ruleset and init array.
@inline function ReadableGridData{X,Y,R}(
    init::AbstractArray, mask, proc, opt, boundary, padval
) where {Y,X,R}
    # We add one extra row and column of status blocks so
    # we dont have to worry about special casing the last block
    if R > 0
        source = _addpadding(init, R, padval)
        dest = _addpadding(init, R, padval)
    else
        if opt isa SparseOpt
            opt = NoOpt()
        end
        source = deepcopy(init)
        dest = deepcopy(init)
    end
    sourcestatus, deststatus = _build_status(opt, source, R)

    grid = ReadableGridData{X,Y,R}(
        init, mask, proc, opt, boundary, padval, source, dest, sourcestatus, deststatus
    )
    _updatestatus!(grid)
    return grid
end

#=
Initialise the block status array.
This tracks whether anything has to be done in an area of the main array.
=#
_updatestatus!(grid::GridData) = _updatestatus!(opt(grid), grid)
function _updatestatus!(opt::SparseOpt, grid)
    blocksize = 2 * radius(grid)
    src = parent(source(grid))
    for i in CartesianIndices(src)
        # Mark the status block if there is a non-zero value
        if !can_skip(opt, src[i])
            bi = _indtoblock.(Tuple(i), blocksize)
            @inbounds sourcestatus(grid)[bi...] = true
            @inbounds deststatus(grid)[bi...] = true
        end
    end
    return nothing
end
_updatestatus!(opt, grid) = nothing

function _copygrid!(grid::GridData{<:Any,<:Any,R}, init) where R
    pad_axes = map(ax -> ax .+ R, axes(init))
    copyto!(parent(source(grid)), CartesianIndices(pad_axes), init, CartesianIndices(init))
    _updatestatus!(grid)
end

# Swap source and dest arrays. Allways returns regular SimData.
_swapsource(d::Tuple) = map(_swapsource, d)
function _swapsource(grid::GridData)
    src = grid.source
    dst = grid.dest
    @set! grid.dest = src
    @set! grid.source = dst
    srcstatus = grid.sourcestatus
    dststatus = grid.deststatus
    @set! grid.deststatus = srcstatus
    return @set grid.sourcestatus = dststatus
end

#=
Find the maximum radius required by all rules
Add padding around the original init array, offset into the negative
So that the first real cell is still 1, 1
=#
function _addpadding(init::AbstractArray{T,N}, r, padval) where {T,N}
    h, w = size(init)
    paddedsize = h + 4r, w + 2r
    paddedindices = -r + 1:h + 3r, -r + 1:w + r
    sourceparent = fill(convert(eltype(init), padval), paddedsize)
    source = OffsetArray(sourceparent, paddedindices...)
    # Copy the init array to the middle section of the source array
    for j in 1:size(init, 2), i in 1:size(init, 1)
        @inbounds source[i, j] = init[i, j]
    end
    return source
end

Base.parent(d::ReadableGridData) = parent(source(d))

@propagate_inbounds function Base.getindex(d::ReadableGridData, I...)
    getindex(source(d), I...)
end

function _build_status(opt::SparseOpt, source, r)
    hoodsize = 2r + 1
    blocksize = 2r
    nblocs = _indtoblock.(size(source), blocksize) .+ 1
    sourcestatus = zeros(Bool, nblocs)
    deststatus = zeros(Bool, nblocs)
    sourcestatus, deststatus
end
_build_status(opt::PerformanceOpt, init, r) = nothing, nothing

"""
    WritableGridData <: GridData

    WritableGridData(grid::GridData)

Passed to rules `<: `, and can be written to directly as
an array. This handles updates to SparseOpt() and writing to
the correct source/dest array.
"""
struct WritableGridData{
    Y,X,R,T,N,I<:AbstractArray{T,N},M,P<:Processor,Op<:PerformanceOpt,Ov,PV,S,D,SSt,DSt
} <: GridData{Y,X,R,T,N,I}
    init::I
    mask::M
    proc::P
    opt::Op
    boundary::Ov
    padval::PV
    source::S
    dest::D
    sourcestatus::SSt
    deststatus::DSt
end
function WritableGridData{Y,X,R}(
    init::I, mask::M, proc::P, opt::Op, boundary::Ov, padval::PV,
    source::S, dest::D, sourcestatus::SSt, deststatus::DSt
) where {Y,X,R,I<:AbstractArray{T,N},M,P,Op,Ov,PV,S,D,SSt,DSt} where {T,N}
    WritableGridData{Y,X,R,T,N,I,M,P,Op,Ov,PV,S,D,SSt,DSt}(
        init, mask, proc, opt, boundary, padval, source, dest, sourcestatus, deststatus
    )
end

Base.parent(d::WritableGridData) = parent(dest(d))
@propagate_inbounds function Base.getindex(d::WritableGridData, I...)
    getindex(source(d), I...)
end
@propagate_inbounds function Base.setindex!(d::WritableGridData, x, I...)
    _setdeststatus!(d, x, I...)
    dest(d)[I...] = x
end

_setdeststatus!(d::WritableGridData, x, I...) = _setdeststatus!(d::WritableGridData, opt(d), x, I...)
function _setdeststatus!(d::WritableGridData{Y,X,R}, opt::SparseOpt, x, I...) where {Y,X,R}
    blockindex = _indtoblock.(I .+ R, 2R)
    @inbounds deststatus(d)[blockindex...] = !(opt.f(x))
    return nothing
end
_setdeststatus!(d::WritableGridData, opt::PerformanceOpt, x, I...) = nothing


text:
using D3Trees
using JSON
using Test
using NBInclude
using Base64

children = [[2,3], [], [4], []]
text = ["one\n(second line)", "2", "III", "four"]
style = ["", "fill:red", "r:14", "opacity:0.7"]
link_style = ["", "stroke:blue", "", "stroke-width:10px"]
tooltip = ["pops", "up", "on", "hover"]
t1 = D3Tree(children,
           text=text,
           style=style,
           tooltip=tooltip,
           link_style=link_style)

@show json(t1)

stringmime(MIME("text/html"), t1)

t2 = D3Tree(children,
            text=text,
            style=style,
            tooltip=tooltip,
            link_style=link_style,
            init_expand=1000,
            init_duration=5000,
            svg_height=300
           )

# inchrome(t1)
# inbrowser(t1, "firefox")
# inchrome(t2)
# inbrowser(t2, "firefox")

@show D3Trees.children(D3TreeNode(t1, 1))
show(stdout, MIME("text/plain"), t1)

@nbinclude("../examples/hello.ipynb")

n = 1_000_000
println("generating $n children")
@time begin
    children = [[i+1] for i in 1:n-1]
    push!(children, [])
end

println("creating tree object")
@time t = D3Tree(children)

println("html string")
@time stringmime(MIME("text/html"), t)

println("AbstractTrees constructor")
@time t3 = D3Tree(t1)
@time t4 = D3Tree(t1, detect_repeat=false)
@test t3.text == t1.text
@test t3.tooltip == t1.tooltip
@test t3.style == t1.style
@test t3.link_style == t1.link_style

inbrowser(D3Tree([[2],[]]), `ls`)


text:
certifysubsetsum(subsetindices, set, sum_) = maximum(subsetindices) <= length(set) && sum([set[i] for i in subsetindices]) == sum_

function certifysubsetsum_long(subsetindices, set, sum_)
	s = 0
	for i in subsetindices
		if i > length(set)
			return false
		end
		s += set[i]
	end
	return s == sum_
end

function certifytsp(path, maxweight)

end


text:
export MtlArgumentEncoder, set_buffer!, set_buffers!, set_constant!

const MTLArgumentEncoder = Ptr{MtArgumentEncoder}

mutable struct MtlArgumentEncoder
    handle::MTLArgumentEncoder
end

function MtlArgumentEncoder(fun::MtlFunction, entry::Integer)
    handle = mtNewArgumentEncoderWithBufferIndexFromFunction(fun, entry-1)
    obj = MtlArgumentEncoder(handle)
    finalizer(unsafe_destroy!, obj)
    return obj
end

function unsafe_destroy!(enc::MtlArgumentEncoder)
    mtRelease(enc.handle)
end

Base.unsafe_convert(::Type{MTLArgumentEncoder}, enc::MtlArgumentEncoder) = enc.handle

Base.:(==)(a::MtlArgumentEncoder, b::MtlArgumentEncoder) = a.handle == b.handle
Base.hash(fun::MtlArgumentEncoder, h::UInt) = hash(mod.handle, h)

## properties

Base.propertynames(o::MtlArgumentEncoder) = (
    # identification
    #=:device, :label,=#
    # creation
    :encodedLength,
    # alignment
    :alignment,
)

function Base.getproperty(o::MtlArgumentEncoder, f::Symbol)
    if f === :encodedLength
        mtArgumentEncoderLength(o)
    elseif f === :alignment
        mtArgumentEncoderAlignment(o)
    else
        getfield(o, f)
    end
end

Base.sizeof(a::MtlArgumentEncoder) = Int(mtArgumentEncoderLength(a))

## operations

# NOTE: indices aren't 1-based here, because they map onto exact IDs in the metadata

function assign_argument_buffer!(enc::MtlArgumentEncoder, buf::MtlBuffer, offset::Integer=0)
    mtArgumentEncoderSetArgumentBufferWithOffset(enc, buf, offset)
end

function assign_argument_buffer!(enc::MtlArgumentEncoder, buf::MtlBuffer, offset::Integer, element::Integer)
    mtArgumentEncoderSetArgumentBufferWithOffsetForElement(enc, buf, offset, element)
end

set_buffer!(enc::MtlArgumentEncoder, buf::MtlBuffer, offset::Integer, index::Integer) =
    mtArgumentEncoderSetBufferOffsetAtIndex(enc, buf, offset, index)
set_buffers!(enc::MtlArgumentEncoder, bufs::Vector{<:MtlBuffer},
             offsets::Vector{Int}, indices::UnitRange{Int}) =
    mtArgumentEncoderSetBuffersOffsetsWithRange(enc, handle_array(bufs), offsets, indices)

function set_constant!(enc::MtlArgumentEncoder, val, index::Integer)
    dst = Base.bitcast(Ptr{typeof(val)}, mtArgumentEncoderConstantDataAtIndex(enc, index))
    Base.unsafe_store!(dst, val, 1)
    return
end


text:
"""
    adress(P::ECC.S256Point, compressed::Bool, testnet::Bool) -> String

Returns the Base58 Bitcoin address given an S256Point
Compressed is set to true if not provided.
Testnet is set to false by default.
"""
function address(P::T, compressed::Bool=true, testnet::Bool=false) where {T<:S256Point}
    s = point2sec(P, compressed)
    h160 = ripemd160(sha256(s))
    if testnet
        prefix = 0x6f
    else
        prefix = 0x00
    end
    result = pushfirst!(h160, prefix)
    return String(base58checkencode(result))
end

"""
    wif(pk::PrivateKey, compressed::Bool=true, testnet::Bool=false) -> String

Returns a PrivateKey in Wallet Import Format
Compressed is set to true if not provided.
Testnet is set to false by default.
"""
function wif(pk::PrivateKey, compressed::Bool=true, testnet::Bool=false)
    secret_bytes = int2bytes(pk.𝑒)
    if testnet
        prefix = 0xef
    else
        prefix = 0x80
    end
    result = pushfirst!(secret_bytes, prefix)
    if compressed
        return String(base58checkencode(push!(result, 0x01)))
    else
        return String(base58checkencode(result))
    end
end


text:
using PyPlot

x_lines = [-.5, 1.1, 0.6]
y_lines = [.6, -.3, 1.3]
plot(x_lines, y_lines)

x_square = [0,1,1,0,0]
y_square = [0,0,1,1,0]
plot(x_square, y_square)

axis("equal")
grid(true)
xlabel("x-coordinate")
ylabel("y-coordinate")
title("Sample plot of some objects");

f(x) = sin(3x) + 0.5sin(7x) - 0.1x
x = range(0, stop=10, length=200)      # Create x-vector, 200 points should make a good plot

plot(x, f.(x))
plot(x, f.(x .- 3), linewidth=3)
plot(x, (@. -0.7f(100/(10 + x^2))), color="m", linestyle="--")
legend(("f(x)", "f(x-3)", "-0.7f(x^2)"))
grid(true);

s = 0:0.001:1
r = @. .1 + s + 0.1sin(50π*s)
x = @. r*cos(10π*s)
y = @. r*sin(10π*s)
plot(x, y, "k")
axis("equal");


text:
# TIPS:
# 1
const DEFAULT = 42

global x = 42
function f(args)
	#function body
end
y = f(x::Int + 1)

# 9
function with_keyword(x; name::String = "Smith")
   # ...
end

# 16
type Person
	name::String
	height::Float64
	weight::Float64
end
# instead of:
# type Person
# 	name
# 	height
# 	weight
# end

# 18
# Use
file = open("file2.txt", "w")
a = "Julia"
b = "Rust"
write(file, a, " ", b)
# instead of:
write(file, "$a $b")

# TOOLS:
# @time

# linter:
Pkg.add("Lint")
using Lint
lintfile("performance.jl")
#
# typechecker:
Pkg.add("TypeCheck")
using TypeCheck
foo1(x::Int) = isprime(x) ? x: false
checkreturntypes(foo1) #> foo1(Int64)::Union(Bool,Int64)

max_stars_repo_path: Module 1/Chapter09/performance.jl

max_issues_repo_path: Module 1/Chapter09/performance.jl

max_forks_repo_path: Module 1/Chapter09/performance.jl


text:
include("../src/Calendar.jl")
using Calendar
using Base.Test

@assert isleapyear(2012)
@assert !isleapyear(2011)
@assert !isleapyear(2013)
@assert isleapyear(2020)
@assert !isleapyear(1900)
@assert !isleapyear(1901)
@assert isleapyear(1904)
@assert isleapyear(2000)

# daylight savings time
t = ymd_hms(2013, 3, 10, 1, 59, 59)
t2 = t + seconds(1)
@assert hour(t2) == 3
@assert minute(t2) == 0
@assert second(t2) == 0

@assert years(4) - hours(1) == CalendarDuration(4,0,0,-3600e3)
@assert hours(1) - years(4) == CalendarDuration(-4,0,0,3600e3)

t1 = ymd_hms(2013, 7, 30, 14, 22, 11)
t2 = ymd_hms(2013, 7, 30, 14, 22, 11.5)
d = t2 - t1
@assert d == seconds(0.5)
@assert d === FixedCalendarDuration(500.)
@assert repr(d) == "0.5 seconds"

@assert repr(years(1) + months(1)) == "1 year + 1 month"
@assert repr(years(2) + months(1)) == "2 years + 1 month"
@assert repr(years(1) + months(2)) == "1 year + 2 months"
@assert repr(years(2) + months(2)) == "2 years + 2 months"
@assert repr(days(3) + seconds(4)) == "3 days + 4 seconds"

f = "yyyy.MM.dd G 'at' HH:mm:ss zzz"
d = ymd_hms(2014,7,4,16,33,29)
s = "2014.07.04 AD at 16:33:29 EDT"
@test Calendar.parse(f, s) == d
@test format(f, d) == s

d1 = yDhms(2014,66,1,2,3)
d2 = ymd_hms(2014,3,7,1,2,3)
@test format(f, d1) == format(f, d2)


text:
# Nauty canonical form wrapper for julia.

"Julia wrapper for the Nauty C library."
module nauty

import LightGraphs

const LIB_FILE = "$(@__DIR__)" * "/minnautywrap"

# {{{ Types and structs

# {{{ Julia versions of two important structs from nauty.h

"""
Define a mutable copy of a struct at the same time as the immutable one and add constructors to both to allow easy conversion.

Required because mutable structs are not substitutable for C structs.

Author: Michael Eastwood
Source: https://discourse.julialang.org/t/passing-an-array-of-structures-through-ccall/5194/15
"""
macro define_mutable(immutable_struct)
   immutable_name = immutable_struct.args[2]
   mutable_name = Symbol(immutable_name, "_mutable")

   mutable_struct = copy(immutable_struct)
   mutable_struct.args[1] = true # set the mutability to true
   mutable_struct.args[2] = mutable_name

   constructors = quote
       function $mutable_name(x::$immutable_name)
           $mutable_name(ntuple(i->getfield(x, i), nfields($immutable_name))...)
       end
       function $immutable_name(x::$mutable_name)
           $immutable_name(ntuple(i->getfield(x, i), nfields($mutable_name))...)
       end
   end

   output = Expr(:block, immutable_struct, mutable_struct, constructors)
   esc(output)
end

const Nboolean = Cint

@define_mutable struct optionblk
    getcanon::Cint             #  make canong and canonlab?
    digraph::Nboolean          #  multiple edges or loops?
    writeautoms::Nboolean      #  write automorphisms?
    writemarkers::Nboolean     #  write stats on pts fixed, etc.?
    defaultptn::Nboolean       #  set lab,ptn,active for single cell?
    cartesian::Nboolean        #  use cartesian rep for writing automs?
    linelength::Cint           #  max chars/line (excl. '\n') for output
    outfile::Ptr{Void} #FILE *outfile;            #  file for output, if any
    userrefproc::Ptr{Void} # void (*userrefproc)       #  replacement for usual refine procedure
         #(graph*,int*,int*,int,int*,int*,set*,int*,int,int);
    userautomproc::Ptr{Void} # void (*userautomproc)     #  procedure called for each automorphism
         # (int,int*,int*,int,int,int);
    userlevelproc::Ptr{Void} # void (*userlevelproc)     #  procedure called for each level
         #(int*,int*,int,int*,statsblk*,int,int,int,int,int,int);
    usernodeproc::Ptr{Void} # void (*usernodeproc)      #  procedure called for each node
         #(graph*,int*,int*,int,int,int,int,int,int);
    usercanonproc::Ptr{Void} # Cint  (*usercanonproc)     #  procedure called for better labellings
         #(graph*,int*,graph*,int,int,int,int);
    invarproc::Ptr{Void} # void (*invarproc)         #  procedure to compute vertex-invariant
         #(graph*,int*,int*,int,int,int,int*,int,Nboolean,int,int);
    tc_level::Cint             #  max level for smart target cell choosing
    mininvarlevel::Cint        #  min level for invariant computation
    maxinvarlevel::Cint        #  max level for invariant computation
    invararg::Cint             #  value passed to (*invarproc)()
    dispatch::Ptr{Void} # dispatchvec *dispatch;    #  vector of object-specific routines
    schreier::Nboolean         #  use random schreier method
    extra_options::Ptr{Void} # void *extra_options;      #  arbitrary extra options
end

"""
# Constructor. Normally default options wouldn't be set by an inner
# constructor, but in this case the nauty manual requires that all
# optionblks be constructed by calling the macro.

Not an inner constructor any more because I don't want to override default constructor...
"""
const optionblk() = ccall((:defaultoptions_graph, LIB_FILE), optionblk, ())

struct statsblk
    grpsize1::Cdouble	#        /* size of group is */
    grpsize2::Cint	#           /*    grpsize1 * 10^grpsize2 */
    numorbits::Cint	#          /* number of orbits in group */
    numgenerators::Cint	#      /* number of generators found */
    errstatus::Cint	#          /* if non-zero : an error code */
    numnodes::Culong	#      /* total number of nodes */
    numbadleaves::Culong	#  /* number of leaves of no use */
    maxlevel::Cint	#                /* maximum depth of search */
    tctotal::Culong	#       /* total size of all target cells */
    canupdates::Culong	#    /* number of updates of best label */
    invapplics::Culong	#    /* number of applications of invarproc */
    invsuccesses::Culong	#  /* number of successful uses of invarproc() */
    invarsuclevel::Cint	#      /* least level where invarproc worked */
end

# }}}

"Julia automatically converts arrays to suitable Ptr{UInt64}s when passed in ccall."
const NautyGraph = Ptr{UInt64}

# }}}

# Interface:

"""
    canonical_form(g::LightGraphs.Graph)

Find the canonical graph, orbits, relabelling and orbits of `g`.

Returns:
- canonical_graph::NautyGraph
    - The canonical graph of the class of isomorphs that g is in
- labelling::Array{Cint}
    - How to relabel g such that it becomes canonical_graph
- orbits::Array{Cint}
    - See nauty documentation/automorphism theory
- partition::Array{Cint}
    - Probably garbage
"""
function canonical_form(g::GraphType) where GraphType <: LightGraphs.SimpleGraphs.AbstractSimpleGraph

    ingraph = lg_to_nauty(g)
    (num_vertices, num_setwords) = size(ingraph, 1, 2)

    # These don't need to be zero'd, I'm just doing it for debugging reasons.
    outgraph = similar(ingraph)
    labelling = Array{Cint}(num_vertices)
    partition = similar(labelling)
    orbits = similar(labelling)

    ccall((:canonical_form, LIB_FILE), Void,
          (NautyGraph, Cint, Cint, Ptr{Cint}, Ptr{Cint}, Ptr{Cint}, NautyGraph),
          ingraph, num_setwords, num_vertices, labelling, partition, orbits, outgraph)

    # Return everything nauty gives us.
    # I'm not sure that partition is meaningful...
    return outgraph, labelling, partition, orbits
end

# {{{ This doesn't work yet.

function densenauty(g, options::optionblk)
    # Return the canonical form in some way that is cheap to test equivalence of.

    ingraph = lg_to_nauty(g)
    (num_vertices, num_setwords) = size(ingraph, 1, 2)

    # These don't need to be zero'd, I'm just doing it for debugging reasons.
    outgraph = zeros(ingraph)
    labelling = zeros(Cint, num_vertices)
    partition = zeros(labelling)
    orbits = zeros(labelling)

    ccall((:densenauty_wrap, LIB_FILE), Void,
          (NautyGraph, Ptr{Cint}, Ptr{Cint}, Ptr{Cint},
           Ref{optionblk}, Cint, Cint, NautyGraph),
          ingraph, labelling, partition, orbits,
          Ref(options), num_setwords, num_vertices, outgraph)

    #= ccall((:densenauty_defaults_wrap, LIB_FILE), Void, =#
    #=       (NautyGraph, Ptr{Cint}, Ptr{Cint}, Ptr{Cint}, =#
    #=        Cint, Cint, NautyGraph), =#
    #=       ingraph, labelling, partition, orbits, =#
    #=       num_setwords, num_vertices, outgraph) =#

    # Return everything nauty gives us.
    # I'm not sure that partition is meaningful...
    return labelling, outgraph, partition, orbits
end
# }}}

# {{{ Helpers

"""
    lg_to_nauty(g::LightGraphs.SimpleGraphs.AbstractSimpleGraph)

Convert to nauty-compatible adjacency matrix (uint array).
"""
function lg_to_nauty(g::GraphType) where GraphType <: LightGraphs.SimpleGraphs.AbstractSimpleGraph
    # Nauty compatible adjacency matrix:
    #   An array of m*n WORDSIZE bitfields.
    #   Where n = num vertices, m = num_setwords = ((n-1) / WORDSIZE) + 1
    #   bitfields are boolean arrays packed from left. Only the first n bits of
    #   each m*WORDSIZE row are significant.

    # assume WORDSIZE = 64, can use nauty_check to confirm values are OK.
    WORDSIZE = 64
    num_vertices = LightGraphs.nv(g)
    num_setwords = div(num_vertices - 1, WORDSIZE) + 1

    # Initialise
    arr = BitArray(num_setwords * WORDSIZE, num_vertices)
    fill!(arr, false)

    # Columns and rows reversed because I care about the column/row layout of
    # arr.chunks, not arr.
    for (rowi, row) = enumerate(g.fadjlist)
        for value = row
            arr[end-value+1,rowi] = true
        end
    end

    # nauty_graph as a vector of UInt64s, just what Nauty wants.
    # For the purposes of ccall, an Array{T} can be reasonably safely treated as Ptr{T}
    return arr.chunks #, num_setwords, num_vertices
end

function label_to_adj(label)
    # change 64 to wordsize
    temp = BitArray(64,size(label,1))
    temp.chunks = label
    temparr = Array{Int64,2}(temp[end-size(label,1)+1:end,:])
    flipdim(temparr',2)
end

# }}}

# Experimental crap to remove...
function graph_receiver(g)
    ingraph = lg_to_nauty(g)
    (num_vertices, num_setwords) = size(ingraph, 1, 2)

    ccall((:graph_receiver, LIB_FILE), UInt64, (NautyGraph, Cint), ingraph, num_vertices * num_setwords)
end

#= # Fields =#

#= type Field{K} end =#

#= Base.convert{K}(::Type{Symbol}, ::Field{K}) = K =#
#= Base.convert(::Type{Field}, s::Symbol) = Field{s}() =#

#= macro f_str(s) =#
#=   :(Field{$(Expr(:quote, symbol(s)))}()) =#
#= end =#

#= typealias FieldPair{F<:Field, T} Pair{F, T} =#
#= const FieldPair{F<:Field, T} = Pair{F, T} =#

#= # Immutable `with` =#

#= for nargs = 1:5 =#
#=   args = [symbol("p$i") for i = 1:nargs] =#
#=   @eval with(x, $([:($p::FieldPair) for p = args]...), p::FieldPair) = =#
#=       with(with(x, $(args...)), p) =#
#= end =#

#= @generated function with{F, T}(x, p::Pair{Field{F}, T}) =#
#=   :($(x.name.primary)($([name == F ? :(p.second) : :(x.$name) =#
#=                          for name in fieldnames(x)]...))) =#
#= end =#

#using Base.Test
#
#@testset begin
#	"Convert the adjacency matrix of a directed graph into an undirected graph."
#	helper(x) = LightGraphs.Graph(x .| x')
#
#	# Two simple isomorphic graphs.
#	iso1a = helper(Array([0 1 1; 0 0 0; 0 0 0]))
#	iso1b = helper(Array([0 0 0; 1 0 1; 0 0 0]))
#
#	@test lg_to_nauty(iso1a) == Array{UInt64,1}([0x6000000000000000, 0x8000000000000000, 0x8000000000000000])
#	@test lg_to_nauty(iso1b) == Array{UInt64,1}([0x4000000000000000, 0xa000000000000000, 0x4000000000000000])
#
#	# Two simple isomorphic digraphs.
#	diso1a = LightGraphs.DiGraph(Array([0 1 1; 0 0 0; 0 0 0]))
#	diso1b = LightGraphs.DiGraph(Array([0 0 0; 1 0 1; 0 0 0]))
#
#	@test lg_to_nauty(diso1a) == Array{UInt64,1}([0x6000000000000000, 0x0000000000000000, 0x0000000000000000])
#	@test lg_to_nauty(diso1b) == Array{UInt64,1}([0x0000000000000000, 0xa000000000000000, 0x0000000000000000])
#
#	@test canonical_form(iso1a)[1] == canonical_form(iso1b)[1]
#
#	@test canonical_form(diso1a)[1] == canonical_form(diso1b)[1]
#end

end


text:
function basicDirectoryTests(comm::Comm{GID, PID, LID}) where {GID, PID, LID}
    n = 8
    nProc = numProc(comm)
    pid = myPid(comm)

    map = BlockMap(n*nProc, n, comm)

    dir = BasicDirectory{GID, PID, LID}(map)
    @test isa(dir, BasicDirectory{GID, PID, LID})

    dir = createDirectory(comm, map)
    @test isa(dir, BasicDirectory{GID, PID, LID})
    @test gidsAllUniquelyOwned(dir)
    @test (repeat(1:nProc, inner=n), repeat(1:n, outer=nProc)) == getDirectoryEntries(dir, map, Vector{GID}(1:n*nProc))
end


text:
@testset "fixedpoints" begin
    A = [0.5 0; 0 0.2]
    b = [1.0; 2.0]
    f(x) = A * x + b

    x0 = rand(2)

    Γx, x, sc = fixedpoint(f, x0)

    @test sc == :success
    @test norm(Γx-x) <= sqrt(eps(eltype(x)))

    As = SMatrix{2,2}(0.5,0,0,0.2)
    bs = SVector(1.0, 2.0)
    xs0 = SVector(0.0, 0.0)
    fs(xs) = As * xs + bs

    Γxs, xs, scs = fixedpoint(fs, xs0)
    @test scs == :success
end


text:
#
# Cache routines for Business Days precalculated days
#

"""
Holds caches for Holiday Calendars.

* Key = `HolidayCalendar` instance
* Value = instance of `HolidayCalendarCache`
"""
global _CACHE_DICT = Dict{HolidayCalendar, HolidayCalendarCache}()

function _getcachestate(hc::HolidayCalendar)
    return haskey(_CACHE_DICT, hc)
end

function _getholidaycalendarcache(hc::HolidayCalendar)
    return _CACHE_DICT[hc]
end

function checkbounds(hcc::HolidayCalendarCache, dt::Date)
    if !((hcc.dtmin <= dt) && (dt <= hcc.dtmax))
        error("Date out of cache bounds. Use initcache function with a wider time spread. Provided date: $(dt).")
    end
end

function _linenumber(hcc::HolidayCalendarCache, dt::Date)
    return Dates.days(dt) - Dates.days(hcc.dtmin) + 1
end

function isbday(hcc::HolidayCalendarCache, dt::Date)
    checkbounds(hcc, dt)
    return hcc.isbday_array[ _linenumber(hcc, dt) ]
end

function bdays(hcc::HolidayCalendarCache, dt0::Date, dt1::Date)
    # Computation is always based on next Business Days if given dates are not Business Days, inspired by Banking Account convention.
    dt0 = tobday(hcc.hc, dt0) # cache bounds are checked inside tobday -> isbday
    dt1 = tobday(hcc.hc, dt1) # cache bounds are checked inside tobday -> isbday

    return Day(convert(Int, hcc.bdayscounter_array[_linenumber(hcc, dt1)]) - convert(Int, hcc.bdayscounter_array[_linenumber(hcc, dt0)]))
end

# Be sure to use this function on a syncronized code (not multithreaded).
"""
    initcache(calendar, [d0], [d1])

Creates cache for a given Holiday Calendar. After calling this function, any call to `isbday`
function, or any function that uses `isbday`, will be optimized to use this cache.

You can pass `calendar` as an instance of `HolidayCalendar`, `Symbol` or `AbstractString`.
You can also pass `calendar` as an `AbstractArray` of those types.
"""
function initcache(hc::HolidayCalendar, d0::Date, d1::Date)
    isbday_array , bdayscounter_array = _createbdayscache(hc, d0, d1)
    _CACHE_DICT[hc] = HolidayCalendarCache(hc, isbday_array, bdayscounter_array, min(d0, d1), max(d0, d1))
end

initcache(calendar, d0::Date, d1::Date) = initcache(convert(HolidayCalendar, calendar), d0, d1)

# Defaults to d0 = 1st Jan 1980 , d1 = 20th dec 2150
function initcache(hc::HolidayCalendar)
    initcache(hc, Date(1980, 01, 01), Date(2150, 12, 20))
end

initcache(calendar) = initcache(convert(HolidayCalendar,calendar))
initcache{A<:AbstractArray}(calendars::A) = initcache(convert(Vector{HolidayCalendar}, calendars))

@vectorize_1arg HolidayCalendar initcache

# remove all elements from cache
function cleancache()
    for k in keys(_CACHE_DICT)
        delete!(_CACHE_DICT, k)
    end
end

"""
    cleancache([calendar])

Cleans cache for a given instance or list of `HolidayCalendar`, `Symbol` or `AbstractString`.
"""
function cleancache(hc::HolidayCalendar)
    if haskey(_CACHE_DICT, hc)
        delete!(_CACHE_DICT, hc)
    end
end

function cleancache(hc_vec::Vector{HolidayCalendar})
    for k in hc_vec
        if k in keys(_CACHE_DICT)
            delete!(_CACHE_DICT, k)
        end
    end
end

cleancache(calendar) = cleancache(convert(HolidayCalendar, calendar))
cleancache{A<:AbstractArray}(calendars::A) = cleancache(convert(Vector{HolidayCalendar}, calendars))

# Returns tuple
# tuple[1] = Array of Bool (isBday) , tuple[2] = Array of UInt32 (bdaycounter)
function _createbdayscache(hc::HolidayCalendar, d0::Date, d1::Date)

    const d0_ = min(d0, d1)
    const d1_ = max(d0, d1)

    const d0_rata = Dates.days(d0_)
    const d1_rata = Dates.days(d1_)

    # length of the cache arrays
    const len::Int = d1_rata - d0_rata + 1

    # This function uses UInt32 to store bdayscounter array
    # We need to check if we'll exceed typemax(UInt32)
    if len > convert(Int, typemax(UInt32))
        error("Maximum size allowed for bdays array is $(typemax(UInt32)). Required lenght of $(len).")
    end

    isbday_array = zeros(Bool, len)
    bdayscounter_array = zeros(UInt32, len)

    isbday_array[1] = isbday(hc, d0_)
    bdayscounter_array[1] = 0

    for i in 2:len
        isbday_array[i] = isbday(hc, d0_ + Dates.Day(i-1))
        bdayscounter_array[i] = bdayscounter_array[i-1] + isbday_array[i]
    end

    return isbday_array, bdayscounter_array
end


text:
include("mlf.jl")

args = [ (1.0, 1.0, 1.0),
         (0.5, 2.0, -0.2),
         (1.5, 1.5, 2.0),
         (2.5, 1.6, 3.0),
         (0.7, 1.1, 2.4) ]

@printf("%8s  %8s  %8s  %12s\n\n", "α", "β", "x", "E_αβ(x)")
for arg_list in args
    α, β, x = arg_list
    y = mlf(α, β, x, 10)
    @printf("%8.4f  %8.4f  %8.4f  %12.8f\n", α, β, x, y)
end


text:
struct LinearElasticity{T,dim,M}
    G::T
    K::T
    Eᵉ::SymmetricTensor{4,dim,T,M}
end

"""
    LinearElasticity{dim}(;G, K) where dim

Construct a `LinearElasticity` material based on the shear modulus `G` and the bulk modulus `K`.
"""
function LinearElasticity{dim}(;G::Float64, K::Float64) where dim
    λ = K - 2/3*G
    δ(i,j) = i==j ? 1.0 : 0.0
    f(i,j,k,l) = λ*δ(i,j)*δ(k,l) + G*(δ(i,k)*δ(j,l) + δ(i,l)*δ(j,k))
    E = SymmetricTensor{4, dim}(f)
    return LinearElasticity(G, K, E)
end

bulk_modulus(;E, ν) = E / 3.0(1.0-2.0ν)
shear_modulus(;E, ν) = E / 2.0(1.0+ν)


text:
@testset "'Design' ............................. " begin

  srand(1234)

  function simonsDesign(r1, n1, r, n)
    nvec = [[n1 for x1 in 0:r1]; [n for x1 in (r1 + 1):n1]]
    cvec = [[Inf for x1 in 0:r1]; [r for x1 in (r1 + 1):n1]]
    return Design(nvec, cvec)
  end

  # Simon's designs for beta = .2, alpha = .05, p1 = p0 +0.2
  # cf. Simon, R "Optimal Two-Stage Designs for Phase II Clinical Trials",
  # Controlled Clinical Trials 10:1-10 (1989). (p. 4)
  p0 = collect(linspace(.1, .7, 7))
  sd = [
    simonsDesign( 1, 10,  5, 29), # p0 = 0.1
    simonsDesign( 3, 13, 12, 43), # p0 = 0.2
    simonsDesign( 5, 15, 18, 46), # p0 = 0.3
    simonsDesign( 7, 16, 23, 46), # p0 = 0.4
    simonsDesign( 8, 15, 26, 43), # p0 = 0.5
    simonsDesign( 7, 11, 30, 43), # p0 = 0.6
    simonsDesign( 4,  6, 22, 27)  # p0 = 0.7
  ]

  @test typeof(sd[1]) == Design{typeof(1), typeof(Inf), NoParameters}

  @test DataFrames.DataFrame(sd[1]) == DataFrames.DataFrame(
        x1 = 0:interimsamplesize(sd[1]),
        n  = samplesize(sd[1]),
        c  = criticalvalue(sd[1])
    )

  @test interimsamplesize(sd[1]) == 10

  @test parameters(sd[1]) == NoParameters()

  @test samplesize(sd[1]) == [[10 for x1 in 0:1]; [29 for x1 in (1 + 1):10]]
  @test samplesize(sd[1], 1) == 10
  @test samplesize(sd[1], 2) == 29

  @test criticalvalue(sd[1]) == [[Inf for x1 in 0:1]; [5 for x1 in (1 + 1):10]]
  @test criticalvalue(sd[1], 1) == Inf
  @test criticalvalue(sd[1], 2) == 5

  supp_ = support(sd[1])
  for p in linspace(0, 1, 11)

    @test pdf(sd[1], 1, 2, p) == 0.0
    pdf_ = pdf.(sd[1], support(sd[1])[:,1], support(sd[1])[:,2], p)
    @test all(pdf_ .>= 0.0)
    @test isapprox(sum(pdf_), 1.0, atol = 0.00001)

  end

  @test power(sd[1], p0[1]) <= .05
  @test power(sd[1], p0[1] + .2) >= .8

  prior = p -> 1
  @test expectedpower(sd[1], prior, mcrv = p0[1] + .1) < power(sd[1], .8)
  @test expectedpower(sd[1], prior, mcrv = p0[1] + .1) > power(sd[1], .2)
  @test expectedpower(sd[1], 1, prior, mcrv = p0[1] + .1) < expectedpower(sd[1], 2, prior, mcrv = p0[1] + .1)

  @test isapprox(
    stoppingforfutility(sd[1], .2), pdf(sd[1], 0, 0, .2) + pdf(sd[1], 1, 0, .2), atol = .00001
  )

  @test !test(sd[1], 0, 0)
  @test test(sd[1], 5, 11)

  sim_ = simulate(sd[1], p0[1], 1000)
  @test mean(sim_[:rejectedH0]) < .05
  @test all(sim_[:n] .<= 29)
  @test all(sim_[:n] .>= 0)

  @test isapprox(quadgk(jeffreysprior(sd[1]), 0, 1)[1], 1.0, atol = 0.0001)

  save("test.jls", sd[1])
  @test true
  rm("test.jls")
  writecsv("test.csv", sd[1]; label = "test")
  @test true
  rm("test.csv")

end

max_stars_repo_name: kkmann/binaryTwoStageDesigns

max_issues_repo_name: kkmann/binaryTwoStageDesigns

max_forks_repo_name: kkmann/binaryTwoStageDesigns


text:
module LossFunctions

import Base.*
using Base.Cartesian
using Markdown
using RecipesBase
using InteractiveUtils: subtypes

import LearnBase.AggMode
import LearnBase.ObsDim
import LearnBase:
    AggregateMode,
    Loss, SupervisedLoss,
    DistanceLoss, MarginLoss,
    value, deriv, deriv2,
    isdistancebased, ismarginbased,
    isminimizable, isdifferentiable,
    istwicedifferentiable,
    isconvex, isstrictlyconvex,
    isstronglyconvex, isnemitski,
    isunivfishercons, isfishercons,
    islipschitzcont, islocallylipschitzcont,
    isclipable, isclasscalibrated, issymmetric

# supervised losses
include("supervised.jl")

# IO and plot recipes
include("printing.jl")
include("plotrecipes.jl")

# deprecations
@deprecate LogitProbLoss CrossEntropyLoss
@deprecate PinballLoss QuantileLoss
@deprecate OrdinalHingeLoss OrdinalMarginLoss{HingeLoss}
@deprecate ScaledDistanceLoss ScaledLoss
@deprecate ScaledMarginLoss ScaledLoss
@deprecate weightedloss(l, w) WeightedMarginLoss(l, w)
@deprecate scaled(l, λ) ScaledLoss(l, λ)
@deprecate value_deriv(l,y,ŷ) (value(l,y,ŷ), deriv(l,y,ŷ))
@deprecate value!(b,l,y,ŷ,args...)  (b .= value(l,y,ŷ,args...))
@deprecate deriv!(b,l,y,ŷ,args...)  (b .= deriv(l,y,ŷ,args...))
@deprecate deriv2!(b,l,y,ŷ,args...) (b .= deriv2(l,y,ŷ,args...))

export
    # loss API
    Loss,
    SupervisedLoss,
    MarginLoss,
    DistanceLoss,
    value, deriv, deriv2,
    isdistancebased, ismarginbased,
    isminimizable, isdifferentiable,
    istwicedifferentiable,
    isconvex, isstrictlyconvex,
    isstronglyconvex, isnemitski,
    isunivfishercons, isfishercons,
    islipschitzcont, islocallylipschitzcont,
    isclipable, isclasscalibrated, issymmetric,

    # relevant submodules
    AggMode, ObsDim,

    # margin-based losses
    ZeroOneLoss,
    LogitMarginLoss,
    PerceptronLoss,
    HingeLoss,
    L1HingeLoss,
    L2HingeLoss,
    SmoothedL1HingeLoss,
    ModifiedHuberLoss,
    L2MarginLoss,
    ExpLoss,
    SigmoidLoss,
    DWDMarginLoss,

    # distance-based losses
    LPDistLoss,
    L1DistLoss,
    L2DistLoss,
    PeriodicLoss,
    HuberLoss,
    EpsilonInsLoss,
    L1EpsilonInsLoss,
    L2EpsilonInsLoss,
    LogitDistLoss,
    QuantileLoss,
    LogCoshLoss,

    # other losses
    MisclassLoss,
    PoissonLoss,
    CrossEntropyLoss,

    # meta losses
    ScaledLoss,
    OrdinalMarginLoss,
    WeightedMarginLoss

end # module


text:
using ITensors,
      ITensorsGPU,
      LinearAlgebra, # For tr()
      CuArrays,
      Test

      # gpu tests!
@testset "cutrunctate" begin
    @test ITensorsGPU.truncate!(CuArrays.zeros(Float64, 10)) == (0., 0., CuArrays.zeros(Float64, 1))
    trunc = ITensorsGPU.truncate!(CuArray([1.0, 0.5, 0.1, 0.05]); absoluteCutoff=true, cutoff=0.2)
    @test trunc[1] ≈ 0.15
    @test trunc[2] ≈ 0.3
    @test trunc[3] == CuArray([1.0, 0.5])
    trunc = ITensorsGPU.truncate!(CuArray([0.5, 0.4, 0.1]); relativeCutoff=true, cutoff=0.2)
    @test trunc[1] ≈ 0.1
    @test trunc[2] ≈ 0.45
    @test trunc[3] == CuArray([0.5])
end # End truncate test


text:
module JML
export to_lexp
export runparser
export scoping_analysis
export global_scope
# export codegen
export to_julia
export global_scope, new_module_spec
export JMLCompilerError

include("Error.jl")
include("Parser.jl")
include("LExp.jl")
include("OpReduction.jl")
include("DExp.jl")
include("Lowering.jl")
include("JMLi.jl")
# include("Codegen.jl")

end # module

max_forks_repo_name: thautwarm/Rupy


text:
module MshReader

  include("MshFileReader.jl")

  export MshFileReader

end # module


text:
function add_to_log(tree::ResDirTree, text::String)
    lfile = logfile(tree);
    if isfile(lfile)
        st = open(lfile, "a");
        text = "\n$(text)\n";
        write(st, text);
        close(st);
    else
        text = "$(text)\n";
        write(lfile, text);
    end
end

max_stars_repo_name: josePereiro/MyResultsManager

max_issues_repo_name: josePereiro/MyResultsManager

max_forks_repo_name: josePereiro/MyResultsManager