text:
"""
    audit_library_name(dag::ExprDAG, options::Dict)
Check whether the library name is valid.
"""
function audit_library_name(dag::ExprDAG, options::Dict)
    general_options = get(options, "general", Dict())
    if haskey(general_options, "library_name")
        try
            Meta.parse(general_options["library_name"])
        catch err
            general_options["library_name"] = "libxxx"
            @warn "failed to parse the `library_name` as a valid Julia expression, default library name: \":libxxx\" is being used instead."
        end
    else
        general_options["library_name"] = "libxxx"
        @warn "default libname: \":libxxx\" is being used, did you forget to set `library_name` in the toml file?"
    end
    if haskey(options, "general")
        merge!(options["general"], general_options)
    else
        options["general"] = general_options
    end
    return nothing
end

"""
    sanity_check(dag::ExprDAG, options::Dict)
Check whether the C code uses bad code-styles listed below:
* typedef an identifier that has the same name as a tag-type "A" but the underlying type
is actually something else. For example,
```c
struct foo {
  int x;
};

struct bar {
  int y;
};

typedef struct foo bar;
typedef struct bar foo;
```
* use the same name for structs and functions. For example,
```c
struct post {
    char m;
};

void post(post m);
```
* two identifiers with same name are actually different things.

It is intended to treat these as "bug"s of the upstream C code.
"""
function sanity_check(dag::ExprDAG, options::Dict)
    for (tk, tv) in dag.tags, (idk, idv) in dag.ids
        idk == tk || continue

        tn = dag.nodes[tv]
        idn = dag.nodes[idv]
        @assert !is_tagtype(idn)

        # in case the audit pass is executed after the `Skip` marking
        is_hardskip(idn) && continue

        if is_typedef(idn)
            ty = getCanonicalType(getTypedefDeclUnderlyingType(idn.cursor))
            if !is_same(getTypeDeclaration(ty), dag.nodes[tv].cursor)
                ifile, iline, icol = get_file_line_column(idn.cursor)
                tfile, tline, tcol = get_file_line_column(tn.cursor)
                error("sanity check failed. [REASON]: typedef an identifier $(idn.cursor) at $ifile:$iline:$icol that has the same name as a tag-type $(tn.cursor) at $tfile:$tline:$tcol. The C code use the same name for different things, which is not a good code-style, please fix it in the upstream.")
            end
        elseif is_function(idn)
            ifile, iline, icol = get_file_line_column(idn.cursor)
            tfile, tline, tcol = get_file_line_column(tn.cursor)
            error("sanity check failed. [REASON]: use the same name $tk for struct $(tn.cursor) at $tfile:$tline:$tcol and function $(idn.cursor) at $ifile:$iline:$icol. The C code use the same name for structs and functions, which is not a good code-style, please fix it in the upstream.")
        else
            error("sanity check failed. please file an issue to Clang.jl.")
        end
    end

    # check whether all duplicated identifiers are actually the same thing (is this even legal in C?)
    for node in dag.nodes
        is_dup_identifier(node) || continue
        for (idk, idv) in dag.ids
            idn = dag.nodes[idv]
            ty1 = getCanonicalType(getTypedefDeclUnderlyingType(node.cursor))
            ty2 = getCanonicalType(getTypedefDeclUnderlyingType(idn.cursor))
            if is_same(getTypeDeclaration(ty1), getTypeDeclaration(ty2))
                file1, line1, col1 = get_file_line_column(node.cursor)
                file2, line2, col2 = get_file_line_column(idn.cursor)
                error("sanity check failed. [REASON]: the identifier $(node.cursor) at $file1:$line1:$col1 and the identifier $(idn.cursor) at $file2:$line2:$col2 have the same name but are actually different things. This is not a good code-style, please fix it in the upstream.")
            end
        end
    end

    return nothing
end

"""
    report_default_tag_types(dag::ExprDAG, options::Dict)
Report those tag-types labeled `StructDefault`/`UnionDefault`/`EnumDefault`.
"""
function report_default_tag_types(dag::ExprDAG, options::Dict)
    default_types = filter(dag.nodes) do node
        node.type isa StructDefault ||
        node.type isa UnionDefault ||
        node.type isa EnumDefault
    end
    if !isempty(default_types) && options["Audit_log"]
        @info "found default tag-types: " default_types
    end
    return nothing
end


text:
abstract type AbstractFacebook <: Name end
struct FacebookV1 <: AbstractFacebook end
struct FacebookV2 <: AbstractFacebook end
struct FacebookV3 <: AbstractFacebook end
struct FacebookV4 <: AbstractFacebook end
struct FacebookV5 <: AbstractFacebook end

# mandatory methods
format(::Type{<:AbstractFacebook}) = TabularData
task(::Type{<:AbstractFacebook}) = Regression
source(::Type{<:AbstractFacebook}) = "https://archive.ics.uci.edu/ml/datasets/Facebook+Comment+Volume+Dataset"

function downloadlink(::Type{<:AbstractFacebook})
    return "https://archive.ics.uci.edu/ml/machine-learning-databases/00363/Dataset.zip"
end

function preprocess(N::Type{<:AbstractFacebook})
    i = parse(Int, string(nameof(N))[end])
    istarget = 54
    iscategorical = vcat(35, 38, 39, 40:53)
    header = false
    files = [
        "Features_Variant_$(i).csv",
        "Features_TestSet.csv",
    ]

    return function (path)
        unzip(path; rename = basename, skip = f -> !in(basename(f.name), files))
        files = joinpath.(dirname(path), files)
        preprocess_csv(N, files[1], :train; istarget, iscategorical, header)
        preprocess_csv(N, files[2], :test; istarget, iscategorical, header)
        return
    end
end

nattributes(::Type{<:AbstractFacebook}) = (54, )
ninstances(::Type{FacebookV1}) = (40949, 0, 10044)
ninstances(::Type{FacebookV2}) = (81312, 0, 10044)
ninstances(::Type{FacebookV3}) = (121098, 0, 10044)
ninstances(::Type{FacebookV4}) = (160424, 0, 10044)
ninstances(::Type{FacebookV5}) = (199030, 0, 10044)

# optional methods
function checksum(::Type{<:AbstractFacebook})
    return "f073dc6c33ee8704cd4eb05d982974debce255028a08babf6d4b0785d9ee31f8"
end

name(::Type{<:AbstractFacebook}) = "Facebook Comment Volume Dataset Data Set"
author(::Type{<:AbstractFacebook}) = ["Kamaljot Singh"]

function citation(::Type{<:AbstractFacebook})
    return """
    @inproceedings{Sing1503:Comment,
        author = 'Kamaljot Singh and Ranjeet Kaur Sandhu and Dinesh Kumar',
        title = 'Comment Volume Prediction Using Neural Networks and Decision Trees',
        booktitle = 'IEEE UKSim-AMSS 17th International Conference on Computer Modelling and
        Simulation, UKSim2015 (UKSim2015)',
        address = 'Cambridge, United Kingdom',
        days = 25,
        month = mar,
        year = 2015,
    }
    """
end

ext: jl

max_stars_repo_name: VaclavMacha/DatasetProvider

max_issues_repo_name: VaclavMacha/DatasetProvider

max_forks_repo_name: VaclavMacha/DatasetProvider


text:
using JLD
using BenchmarkTools
using Dates
using Base.Sys
using CellListsBenchmarks
using Plots

gettime(x) = getfield(x, :time)
gettime_nogc(x) = getfield(x, :time) - getfield(x, :gctime)

files = [
    "n100-d2-r0.01.jld",
    "n100-d3-r0.01.jld",
    "n200-d4-r0.01.jld",
    "n300-d5-r0.01.jld"
]
dir_brute_force = "output/brute_force_2021-04-09T08:28:23+0300/"
dir_cell_lists = "output/cell_lists_2021-04-09T08:28:29+0300/"

for file in files
    lb = load(joinpath(dir_brute_force, file))["benchmarks"]
    lc = load(joinpath(dir_cell_lists, file))["benchmarks"]

    d = getfield(lb[1], :d)
    r = getfield(lc[1], :r)
    ns = getfield.(lb, :n)
    tb = getfield.(lb, :trials)
    tc = getfield.(lc, :trials)

    mb = [median.(t) for t in tb]
    mc = [median.(t) for t in tc]

    meanb = [mean(gettime.(v)) for v in mb]
    meanc = [mean(gettime.(v)) for v in mc]

    p1 = plot(legend=:topleft)
    plot!(p1, ns, meanb, markershape=:circle, markersize=2, linestyle=:dash, label="brute force")
    plot!(p1, ns, meanc, markershape=:circle, markersize=2, linestyle=:dash, label="cell lists")

    p2 = plot(legend=:topleft)
    plot!(p2, ns, meanb./meanc, markershape=:circle, markersize=2, linestyle=:dash, label="brute force/cell lists")

    directory = joinpath("figures", "cell_lists_vs_brute_force")
    if !isdir(directory)
        mkpath(directory)
    end

    savefig(p1, joinpath(directory, "d$d-r$r-n$(ns[end]).svg"))
    savefig(p2, joinpath(directory, "d$d-r$r-n$(ns[end])-ratio.svg"))
end

max_stars_repo_name: jaantollander/cell-lists-benchmarks

max_issues_repo_name: jaantollander/cell-lists-benchmarks

max_forks_repo_name: jaantollander/cell-lists-benchmarks


text:
"""
"""
struct FCycle <: AbstractCycle end

"""
"""
struct FkCycle <: AbstractCycle
    k::Int
end

""
function _core_cycle!(problem, level, ::FCycle)
    process_fine!(problem, level)
    _descend!(problem, level, FCycle())
    process_coarse!(problem, level)

    process_uncoarse!(problem, level)
    _descend!(problem, level, VCycle())
    process_coarse!(problem, level)
    return nothing
end

""
@inline function _repeatcycle!(problem, level, cycle_t::FkCycle)
    for i in 1:(cycle_t.k)
        process_uncoarse!(problem, level)
        _descend!(problem, level, VCycle())
        process_coarse!(problem, level)
    end
    return nothing
end

""
function _core_cycle!(problem, level, cycle_t::FkCycle)
    process_fine!(problem, level)
    _descend!(problem, level, cycle_t)
    process_coarse!(problem, level)

    _repeatcycle!(problem, level, cycle_t)
    return nothing
end


text:
# Parts of this file are based on:
# 1. https://github.com/JuliaPlots/Plots.jl/blob/master/src/backends/web.jl
# 2. https://github.com/JuliaGraphics/Luxor.jl/blob/master/src/Luxor.jl
# 3. https://github.com/tpapp/DefaultApplication.jl/blob/master/src/DefaultApplication.jl

import FileIO

# """
# """
# function open_browser_window(filename::Nothing, a::AbstractDict = ENV)
#     @warn("no filename to open")
#     return filename
# end

# """
# """
# function open_browser_window(filename::AbstractString, a::AbstractDict = ENV)
#     extension = filename_extension(filename)
#     is_svg_file = extension == ".svg"
#     is_png_file = extension == ".png"
#     is_svg_or_png_file = is_svg_file || is_png_file
#     is_ijulia = isdefined(Main, :IJulia) && Main.IJulia.inited
#     if is_ijulia && is_svg_or_png_file
#         Main.IJulia.clear_output(true)
#         if is_svg_file
#             open(filename) do f
#                 display(
#                     "image/svg+xml",
#                     readstring(f),
#                     )
#             end
#         elseif is_png_file
#             # We use Base.invokelatest to avoid world age errors
#             Base.invokelatest(
#                 display,
#                 "image/png",
#                 FileIO.load(filename),
#                 )
#         end
#     elseif (is_make_examples(a)) ||
#             (is_make_docs(a)) ||
#             (is_runtests(a) && !open_plots_during_tests(a))
#         @debug(
#             string(
#                 "Skipping opening file: ",
#                 filename,
#                 )
#             )
#     else
#         @debug(string("Opening file ",filename,))
#         if Sys.isapple()
#             try
#                 run(`open $(filename)`)
#             catch e
#                 @warn(string("ignoring error: "), e)
#             end
#         elseif Sys.islinux() || Sys.isbsd()
#             try
#                 run(`xdg-open $(filename)`)
#             catch e
#                 @warn(string("ignoring error: "), e)
#             end
#         elseif Sys.iswindows()
#             try
#                 run(`$(ENV["COMSPEC"]) /c start "" "$(filename)"`)
#             catch e
#                 @warn(string("ignoring error: "), e)
#             end
#         else
#             @warn(
#                 string(
#                     "unknown operating system; could not open file ",
#                     filename,
#                     )
#                 )
#         end
#     end
#     return filename
# end


text:
@testset "Attrition tests" begin

piseed = floor( Int, pi * 1_000_000 )
eseed = floor( Int, exp( 1 ) * 1_000_000 )
s2seed = floor( Int, sqrt( 2 ) * 1_000_000 )

mpSim = ManpowerSimulation( "sim" )
setSimulationLength!( mpSim, 300 )

attribute = Attribute( "attribute A" )
addInitialAttributeValue!( attribute, "value A", 1 )
addSimulationAttribute!( mpSim, attribute )

attrition = Attrition( "default" )
setAttritionRate!( attrition, 0.1 )
setAttritionPeriod!( attrition, 12 )
addSimulationAttrition!( mpSim, attrition )

node = BaseNode( "node A" )
addNodeRequirement!( node, "attribute A", "value A" )
addSimulationBaseNode!( mpSim, node )

recruitment = Recruitment( "EW" )
setRecruitmentSchedule!( recruitment, 600 )
setRecruitmentTarget!( recruitment, "node A" )
setRecruitmentFixed!( recruitment, 1000 )
addSimulationRecruitment!( mpSim, recruitment )

@testset "Basic attrition test" begin
    @test verifySimulation!( mpSim )
    run( mpSim, saveConfig=false )
    report = nodeFluxReport( mpSim, 12, :out, "active" )["active"]
    @test !all( report[:, "active => external"] .== 0 ) &&
        !all( report[:, Symbol( "attrition: node A => external" )] .== 0 )
end  # @testset "Basic attrition test"

@testset "Piecewise attrition test" begin
    attrition = Attrition( "default" )
    setAttritionPeriod!( attrition, 12 )
    setAttritionCurve!( attrition, Dict( 0.0 => 0.1, 120.0 => 0.0 ) )
    addSimulationAttrition!( mpSim, attrition )
    @test verifySimulation!( mpSim )
    run( mpSim, saveConfig=false )
    report = nodeFluxReport( mpSim, 12, :out, "active" )["active"]
    @test all( report[12:end, "active => external"] .== 0 )
end  # @testset "Piecewise attrition test"

@testset "Segregated attrition test" begin
    attribute = Attribute( "attribute A" )
    addInitialAttributeValue!( attribute, "value A", 1 )
    addInitialAttributeValue!( attribute, "value B", 1 )
    addSimulationAttribute!( mpSim, attribute )

    attrition = Attrition( "default" )
    setAttritionRate!( attrition, 0 )
    setAttritionPeriod!( attrition, 12 )
    addSimulationAttrition!( mpSim, attrition )

    attrition = Attrition( "attrition B" )
    setAttritionRate!( attrition, 0.1 )
    setAttritionPeriod!( attrition, 12 )
    addSimulationAttrition!( mpSim, attrition )

    node = BaseNode( "node B" )
    addNodeRequirement!( node, "attribute A", "value B" )
    setNodeAttritionScheme!( node, "attrition B" )
    addSimulationBaseNode!( mpSim, node )

    recruitment = Recruitment( "EW" )
    setRecruitmentSchedule!( recruitment, 600 )
    setRecruitmentTarget!( recruitment, "node B" )
    setRecruitmentFixed!( recruitment, 1000 )
    addSimulationRecruitment!( mpSim, recruitment )

    @test verifySimulation!( mpSim )
    run( mpSim, saveConfig=false )
    report = nodeFluxReport( mpSim, 12, :out, "node A", "node B" )
    @test all( report["node A"][:, Symbol( "node A => other" )] .== 0 )
end  # @testset "Segregated attrition test"

mpSim = ManpowerSimulation( "sim" )
setSimulationLength!( mpSim, 300 )

clearSimulationAttributes!( mpSim )
attribute = Attribute( "attribute A" )
addInitialAttributeValue!( attribute, "value A", 1 )
addSimulationAttribute!( mpSim, attribute )

clearSimulationBaseNodes!( mpSim )
node = BaseNode( "node A" )
addNodeRequirement!( node, "attribute A", "value A" )
addSimulationBaseNode!( mpSim, node )

clearSimulationRecruitment!( mpSim )
recruitment = Recruitment( "EW" )
setRecruitmentSchedule!( recruitment, 600 )
setRecruitmentTarget!( recruitment, "node A" )
setRecruitmentFixed!( recruitment, 1000 )
addSimulationRecruitment!( mpSim, recruitment )

@testset "Distribution test I" begin
    attrition = Attrition( "default" )
    setAttritionRate!( attrition, 0.1 )
    setAttritionPeriod!( attrition, 12 )
    addSimulationAttrition!( mpSim, attrition )

    @test verifySimulation!( mpSim )
    testDistribution = truncated( Exponential( -12 / log(0.9) ), 0, 300.0 )

    run( mpSim, saveConfig=false, seed=piseed )
    attritionTimes = Vector{Float64}( DataFrame( DBInterface.execute(
        mpSim.simDB,
        "SELECT * FROM Transitions_Sim WHERE targetNode IS NULL" ) )[:,
        :timeIndex] )
    report = nodeFluxReport( mpSim, 12, :out, "active" )["active"]
    pval = pvalue( ExactOneSampleKSTest( attritionTimes, testDistribution ) )
    @test pval > 0.05
    @test all( report[2:end, "active => external"] .==
        [116, 81, 81, 77, 53, 58, 45, 57, 52, 41, 37, 31, 28, 23, 23, 15, 19,
        13, 18, 9, 18, 14, 9, 11, 4] )

    run( mpSim, saveConfig=false, seed=piseed )
    attritionTimes = Vector{Float64}( DataFrame( DBInterface.execute(
        mpSim.simDB,
        "SELECT * FROM Transitions_Sim WHERE targetNode IS NULL" ) )[:,
        :timeIndex] )
    report = nodeFluxReport( mpSim, 12, :out, "active" )["active"]
    @test pval ==
        pvalue( ExactOneSampleKSTest( attritionTimes, testDistribution ) )
    @test all( report[2:end, "active => external"] .==
        [116, 81, 81, 77, 53, 58, 45, 57, 52, 41, 37, 31, 28, 23, 23, 15, 19,
        13, 18, 9, 18, 14, 9, 11, 4] )

    run( mpSim, saveConfig=false, seed=piseed, sysEnt=true )
    attritionTimes = Vector{Float64}( DataFrame( DBInterface.execute(
        mpSim.simDB,
        "SELECT * FROM Transitions_Sim WHERE targetNode IS NULL" ) )[:,
        :timeIndex] )
    report = nodeFluxReport( mpSim, 12, :out, "active" )["active"]
    @test pval !=
        pvalue( ExactOneSampleKSTest( attritionTimes, testDistribution ) )
    @test report[2:end, "active => external"] !=
        [116, 81, 81, 77, 53, 58, 45, 57, 52, 41, 37, 31, 28, 23, 23, 15, 19,
        13, 18, 9, 18, 14, 9, 11, 4]

    run( mpSim, saveConfig=false, seed=eseed )
    attritionTimes = Vector{Float64}( DataFrame( DBInterface.execute(
        mpSim.simDB,
        "SELECT * FROM Transitions_Sim WHERE targetNode IS NULL" ) )[:,
        :timeIndex] )
    report = nodeFluxReport( mpSim, 12, :out, "active" )["active"]
    @test pvalue( ExactOneSampleKSTest( attritionTimes, testDistribution ) ) >
        0.05
    @test all( report[2:end, "active => external"] .==
        [112, 85, 76, 83, 51, 60, 60, 52, 40, 39, 34, 36, 21, 26, 22, 20, 20,
        14, 16, 15, 15, 8, 7, 14, 7] )

    run( mpSim, saveConfig=false, seed=s2seed )
    attritionTimes = Vector{Float64}( DataFrame( DBInterface.execute(
        mpSim.simDB,
        "SELECT * FROM Transitions_Sim WHERE targetNode IS NULL" ) )[:,
        :timeIndex] )
    report = nodeFluxReport( mpSim, 12, :out, "active" )["active"]
    @test pvalue( ExactOneSampleKSTest( attritionTimes, testDistribution ) ) >
        0.05
    @test all( report[2:end, "active => external"] .==
        [118, 91, 70, 63, 69, 55, 56, 50, 35, 44, 37, 36, 27, 21, 16, 24, 18,
        14, 7, 14, 8, 10, 10, 12, 7] )
end  # @testset "Distribution test I"

@testset "Distribution test II" begin
    attrition = Attrition( "default" )
    setAttritionPeriod!( attrition, 12 )
    setAttritionCurve!( attrition, Dict( 0.0 => 0.1, 120.0 => 0.0, 180.0 => 0.25, 240.0 => 0.0 ) )
    addSimulationAttrition!( mpSim, attrition )

    @test verifySimulation!( mpSim )
    testDistribution1 = truncated( Exponential( -12 / log(0.9) ), 0, 120.0 )
    p1 = 1 - 0.9^10
    testDistribution2 = truncated( Exponential( -12 / log(0.75) ), 180.0,
        240.0 )
    p2 = 0.9^10 * ( 1 - 0.75^5 )
    p = [p1, p2]
    p ./= sum( p )
    testDistribution = MixtureModel( [testDistribution1, testDistribution2],
        p )

    run( mpSim, saveConfig=false, seed=piseed )
    attritionTimes = Vector{Float64}( DataFrame( DBInterface.execute(
        mpSim.simDB,
        "SELECT * FROM Transitions_Sim WHERE targetNode IS NULL" ) )[:,
        :timeIndex] )
    report = nodeFluxReport( mpSim, 12, :out, "active" )["active"]
    pval = pvalue( ExactOneSampleKSTest( attritionTimes, testDistribution ) )
    @test pval > 0.05
    @test all( report[2:end, "active => external"] .==
        [116, 81, 81, 77, 53, 58, 45, 57, 52, 41, 0, 0, 0, 0, 0, 87, 61, 46,
        40, 29, 0, 0, 0, 0, 0] )

    run( mpSim, saveConfig=false, seed=piseed )
    attritionTimes = Vector{Float64}( DataFrame( DBInterface.execute(
        mpSim.simDB,
        "SELECT * FROM Transitions_Sim WHERE targetNode IS NULL" ) )[:,
        :timeIndex] )
    report = nodeFluxReport( mpSim, 12, :out, "active" )["active"]
    @test pval ==
        pvalue( ExactOneSampleKSTest( attritionTimes, testDistribution ) )
    @test all( report[2:end, "active => external"] .==
        [116, 81, 81, 77, 53, 58, 45, 57, 52, 41, 0, 0, 0, 0, 0, 87, 61, 46,
        40, 29, 0, 0, 0, 0, 0] )

    run( mpSim, saveConfig=false, seed=piseed, sysEnt=true )
    attritionTimes = Vector{Float64}( DataFrame( DBInterface.execute(
        mpSim.simDB,
        "SELECT * FROM Transitions_Sim WHERE targetNode IS NULL" ) )[:,
        :timeIndex] )
    report = nodeFluxReport( mpSim, 12, :out, "active" )["active"]
    @test pval !=
        pvalue( ExactOneSampleKSTest( attritionTimes, testDistribution ) )
    @test report[2:end, "active => external"] !=
        [116, 81, 81, 77, 53, 58, 45, 57, 52, 41, 0, 0, 0, 0, 0, 87, 61, 46,
        40, 29, 0, 0, 0, 0, 0]

    run( mpSim, saveConfig=false, seed=eseed )
    attritionTimes = Vector{Float64}( DataFrame( DBInterface.execute(
        mpSim.simDB,
        "SELECT * FROM Transitions_Sim WHERE targetNode IS NULL" ) )[:,
        :timeIndex] )
    report = nodeFluxReport( mpSim, 12, :out, "active" )["active"]
    @test pvalue( ExactOneSampleKSTest( attritionTimes, testDistribution ) ) >
        0.05
    @test all( report[2:end, "active => external"] .==
        [112, 85, 76, 83, 51, 60, 60, 52, 40, 39, 0, 0, 0, 0, 0, 90, 58, 47,
        42, 24, 0, 0, 0, 0, 0] )

    run( mpSim, saveConfig=false, seed=s2seed )
    attritionTimes = Vector{Float64}( DataFrame( DBInterface.execute(
        mpSim.simDB,
        "SELECT * FROM Transitions_Sim WHERE targetNode IS NULL" ) )[:,
        :timeIndex] )
    report = nodeFluxReport( mpSim, 12, :out, "active" )["active"]
    @test pvalue( ExactOneSampleKSTest( attritionTimes, testDistribution ) ) >
        0.05
    @test all( report[2:end, "active => external"] .==
        [118, 91, 70, 63, 69, 55, 56, 50, 35, 44, 0, 0, 0, 0, 0, 96, 54, 46,
        24, 30, 0, 0, 0, 0, 0] )
end  # @testset "Distribution test II"

end  # @testset "Attrition tests"


text:
# Automatically generated using Clang.jl wrap_c, version 0.0.0

export
    unix,
    linux,
    AVFILTER_FLAG_DYNAMIC_INPUTS,
    AVFILTER_FLAG_DYNAMIC_OUTPUTS,
    AVFILTER_FLAG_SLICE_THREADS,
    AVFILTER_FLAG_SUPPORT_TIMELINE_GENERIC,
    AVFILTER_FLAG_SUPPORT_TIMELINE_INTERNAL,
    AVFILTER_FLAG_SUPPORT_TIMELINE,
    AVFILTER_THREAD_SLICE,
    AVFILTER_CMD_FLAG_ONE,
    AVFILTER_CMD_FLAG_FAST,
    AVFilterPad,
    AVFilter_,
    AVFilterGraphInternal,
    AVFilterGraph,
    AVFilterInternal,
    AVFilterContext,
    AVFilterFormats,
    ANONYMOUS_1,
    AVFILTER_AUTO_CONVERT_ALL,
    AVFILTER_AUTO_CONVERT_NONE,
    AVFilterInOut,
    AV_BUFFERSINK_FLAG_PEEK,
    AV_BUFFERSINK_FLAG_NO_REQUEST,
    AVBufferSinkParams,
    AVABufferSinkParams,
    ANONYMOUS_2,
    AV_BUFFERSRC_FLAG_NO_CHECK_FORMAT,
    AV_BUFFERSRC_FLAG_PUSH,
    AV_BUFFERSRC_FLAG_KEEP_REF,
    AVBufferSrcParameters,
    LIBAVFILTER_VERSION_MAJOR,
    LIBAVFILTER_VERSION_MINOR,
    LIBAVFILTER_VERSION_MICRO,
    LIBAVFILTER_BUILD

const unix = 1
const linux = 1
const AVFILTER_FLAG_DYNAMIC_INPUTS = 1 << 0
const AVFILTER_FLAG_DYNAMIC_OUTPUTS = 1 << 1
const AVFILTER_FLAG_SLICE_THREADS = 1 << 2
const AVFILTER_FLAG_SUPPORT_TIMELINE_GENERIC = 1 << 16
const AVFILTER_FLAG_SUPPORT_TIMELINE_INTERNAL = 1 << 17
const AVFILTER_FLAG_SUPPORT_TIMELINE = AVFILTER_FLAG_SUPPORT_TIMELINE_GENERIC | AVFILTER_FLAG_SUPPORT_TIMELINE_INTERNAL
const AVFILTER_THREAD_SLICE = 1 << 0
const AVFILTER_CMD_FLAG_ONE = 1
const AVFILTER_CMD_FLAG_FAST = 2

struct AVFilterPad
end

struct AVFilter_
    name::Cstring
    description::Cstring
    inputs::Ptr{AVFilterPad}
    outputs::Ptr{AVFilterPad}
    priv_class::Ptr{AVClass}
    flags::Cint
    preinit::Ptr{Cvoid}
    init::Ptr{Cvoid}
    init_dict::Ptr{Cvoid}
    uninit::Ptr{Cvoid}
    query_formats::Ptr{Cvoid}
    priv_size::Cint
    flags_internal::Cint
    next::Ptr{AVFilter}
    process_command::Ptr{Cvoid}
    init_opaque::Ptr{Cvoid}
    activate::Ptr{Cvoid}
end

struct AVFilterGraphInternal
end

struct AVFilterGraph
    av_class::Ptr{AVClass}
    filters::Ptr{Ptr{AVFilterContext}}
    nb_filters::UInt32
    scale_sws_opts::Cstring
    resample_lavr_opts::Cstring
    thread_type::Cint
    nb_threads::Cint
    internal::Ptr{AVFilterGraphInternal}
    opaque::Ptr{Cvoid}
    execute::Ptr{avfilter_execute_func}
    aresample_swr_opts::Cstring
    sink_links::Ptr{Ptr{AVFilterLink}}
    sink_links_count::Cint
    disable_auto_convert::UInt32
end

struct AVFilterInternal
end

struct AVFilterContext
    av_class::Ptr{AVClass}
    filter::Ptr{AVFilter_}
    name::Cstring
    input_pads::Ptr{AVFilterPad}
    inputs::Ptr{Ptr{AVFilterLink}}
    nb_inputs::UInt32
    output_pads::Ptr{AVFilterPad}
    outputs::Ptr{Ptr{AVFilterLink}}
    nb_outputs::UInt32
    priv::Ptr{Cvoid}
    graph::Ptr{AVFilterGraph}
    thread_type::Cint
    internal::Ptr{AVFilterInternal}
    command_queue::Ptr{AVFilterCommand}
    enable_str::Cstring
    enable::Ptr{Cvoid}
    var_values::Ptr{Cdouble}
    is_disabled::Cint
    hw_device_ctx::Ptr{AVBufferRef}
    nb_threads::Cint
    ready::UInt32
    extra_hw_frames::Cint
end

struct AVFilterFormats
end

# begin enum ANONYMOUS_1
const ANONYMOUS_1 = Cint
const AVFILTER_AUTO_CONVERT_ALL = 0 |> Int32
const AVFILTER_AUTO_CONVERT_NONE = -1 |> Int32
# end enum ANONYMOUS_1

struct AVFilterInOut
    name::Cstring
    filter_ctx::Ptr{AVFilterContext}
    pad_idx::Cint
    next::Ptr{AVFilterInOut}
end

const AV_BUFFERSINK_FLAG_PEEK = 1
const AV_BUFFERSINK_FLAG_NO_REQUEST = 2

struct AVBufferSinkParams
    pixel_fmts::Ptr{AVPixelFormat}
end

struct AVABufferSinkParams
    sample_fmts::Ptr{AVSampleFormat}
    channel_layouts::Ptr{Int64}
    channel_counts::Ptr{Cint}
    all_channel_counts::Cint
    sample_rates::Ptr{Cint}
end

# begin enum ANONYMOUS_2
const ANONYMOUS_2 = UInt32
const AV_BUFFERSRC_FLAG_NO_CHECK_FORMAT = 1 |> UInt32
const AV_BUFFERSRC_FLAG_PUSH = 4 |> UInt32
const AV_BUFFERSRC_FLAG_KEEP_REF = 8 |> UInt32
# end enum ANONYMOUS_2

struct AVBufferSrcParameters
    format::Cint
    time_base::AVRational
    width::Cint
    height::Cint
    sample_aspect_ratio::AVRational
    frame_rate::AVRational
    hw_frames_ctx::Ptr{AVBufferRef}
    sample_rate::Cint
    channel_layout::UInt64
end

const LIBAVFILTER_VERSION_MAJOR = 7
const LIBAVFILTER_VERSION_MINOR = 40
const LIBAVFILTER_VERSION_MICRO = 101

# Skipping MacroDefinition: LIBAVFILTER_VERSION_INT AV_VERSION_INT ( LIBAVFILTER_VERSION_MAJOR , LIBAVFILTER_VERSION_MINOR , LIBAVFILTER_VERSION_MICRO )
# Skipping MacroDefinition: LIBAVFILTER_VERSION AV_VERSION ( LIBAVFILTER_VERSION_MAJOR , LIBAVFILTER_VERSION_MINOR , LIBAVFILTER_VERSION_MICRO )

const LIBAVFILTER_BUILD = LIBAVFILTER_VERSION_INT

# Skipping MacroDefinition: LIBAVFILTER_IDENT "Lavfi" AV_STRINGIFY ( LIBAVFILTER_VERSION )
# Skipping MacroDefinition: FF_API_OLD_FILTER_OPTS_ERROR ( LIBAVFILTER_VERSION_MAJOR < 8 )
# Skipping MacroDefinition: FF_API_LAVR_OPTS ( LIBAVFILTER_VERSION_MAJOR < 8 )
# Skipping MacroDefinition: FF_API_FILTER_GET_SET ( LIBAVFILTER_VERSION_MAJOR < 8 )
# Skipping MacroDefinition: FF_API_NEXT ( LIBAVFILTER_VERSION_MAJOR < 8 )


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

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

module CyclotomicBinaryForms
using Nemo, Counts, NumberTheory

export ModuleCyclotomicBinaryForms
export isA206864, F206864, I206864, L206864
export isA206942, F206942, I206942, L206942
export isA293654, F293654, I293654, L293654
export isA296095, F296095, I296095, L296095
export V299214, L299214
export isA299498, F299498, I299498, L299498
export isA299733, L299733
export isA299928, F299928, I299928, L299928
export isA299929, F299929, I299929, L299929
export isA299930, F299930, I299930, L299930
export isA325143, F325143, I325143, L325143
export isA325145, F325145, I325145, L325145

"""
E. Fouvry, C. Levesque, M. Waldschmidt,
[Representation of integers by cyclotomic binary forms](https://arxiv.org/pdf/1712.09019.pdf), arXiv:1712.09019 [math.NT], 2017.

* isA206864, F206864, I206864, L206864, isA206942, F206942, I206942, L206942, isA293654, F293654, I293654, L293654, isA296095, F296095, I296095, L296095, V299214, L299214, isA299498, F299498, I299498, L299498, isA299733, L299733, isA299928, F299928, I299928, L299928, isA299929, F299929, I299929, L299929, isA299930, F299930, I299930, L299930, isA325143, F325143, I325143, L325143, isA325145, F325145, I325145, L325145
"""
const ModuleCyclotomicBinaryForms = ""

#A206942########################################################################
"""
Is ``n`` a numbers of the form ``Phi_k(m)`` with ``k > 2`` and ``|m| > 1``
where ``Phi_k(m)`` denotes the ``k``-th cyclotomic polynomial evaluated at ``m``.
"""
function isA206942(n)
    if n < 3 return false end
    R, x = PolynomialRing(ZZ, "x")
    K = floor(Int, 5.383*log(n)^1.161) # Bounds from
    M = floor(Int, 2*sqrt(n/3))        # Fouvry & Levesque & Waldschmidt
    for k in 3:K
		c = cyclotomic(k, x)
        for m in 2:M
            n == subst(c, m) && return true
        end
    end
    return false
end

"""
Filter the integers which are A206942 and <= ``n``.
"""
F206942(n) = filter(isA206942, 1:n)

"""
Iterate over the first ``n`` integers which are A206942.
"""
I206942(n) = takeFirst(isA206942, n)

"""
Return the list of the first ``n`` integers which are A206942.
"""
L206942(n) = collect(I206942(n))

#A206864########################################################################

"""
Is ``n`` a prime of the form ``Phi_k(m)`` with ``k > 2`` and ``|m| > 1``
where ``Phi_k(m)`` denotes the ``k``-th cyclotomic polynomial evaluated at ``m``.
"""
isA206864(n) = isPrime(n) && isA206942(n)

"""
Filter the integers which are A206864 and <= ``n``.
"""
F206864(n) = (j for j in 1:n if isA206864(j))

"""
Iterate over the first ``n`` integers which are A206864.
"""
I206864(n) = takeFirst(isA206864, n)

"""
Return the list of the first ``n`` integers which are A206864.
"""
L206864(n) = collect(I206864(n))

#A299498########################################################################

"""
Is ``n`` primitively represented by a cyclotomic binary forms?
"""
function isA299498(n)
    isPrimeTo(n, k) = gcd(ZZ(n), ZZ(k)) == ZZ(1)
    R, x = PolynomialRing(ZZ, "x")
    K = floor(Int, 5.383*log(n)^1.161) # Bounds from
    M = floor(Int, 2*sqrt(n/3))  # Fouvry & Levesque & Waldschmidt
    N = QQ(n)

    for k in 3:K
        e = Int(eulerphi(ZZ(k)))
        c = cyclotomic(k, x)
        for m in 1:M, j in m+1:M if isPrimeTo(m, j)
            N == m^e*subst(c, QQ(j,m)) && return true
    end end end
    return false
end

"""
Filter the integers which are A299498 and <= ``n``.
"""
F299498(n) = (j for j in 1:n if isA299498(j))

"""
Iterate over the first ``n`` integers which are A299498.
"""
I299498(n) = takeFirst(isA299498, n)

"""
Return the list of the first ``n`` integers which are A299498.
"""
L299498(n) = collect(I299498(n))

"""
Count the number of cyclotomic binary forms which primitively represent ``n``.
"""
function countA299498(n)
    if n < 3 return 0 end
    isPrimeTo(n, k) = gcd(ZZ(n), ZZ(k)) == ZZ(1)
    R, x = PolynomialRing(ZZ, "x")
    K = floor(Int, 5.383*log(n)^1.161) # Bounds from
    M = floor(Int, 2*sqrt(n/3))  # Fouvry & Levesque & Waldschmidt
    N = QQ(n); count = 0

    for k in 3:K
        e = Int(eulerphi(ZZ(k)))
        c = cyclotomic(k, x)
        for m in 1:M, j in m+1:M if isPrimeTo(m, j)
            N == m^e*subst(c, QQ(j,m)) && (count += 1)
    end end end
    count
end
# [countA299498(n) for n in 1:100] |> println

#A299928########################################################################

"""
Is ``n`` represented by a cyclotomic binary form f(x, y) where x and y are prime numbers and 0 < y < x?
"""
function isA299928(n)
    R, z = PolynomialRing(ZZ, "z")
    K = floor(Int, 5.383*log(n)^1.161) # Bounds from
    M = floor(Int, 2*sqrt(n/3))  # Fouvry & Levesque & Waldschmidt
    N = QQ(n)
    P(u) = (p for p in u:M if isprime(ZZ(p)))
    for k in 3:K
        e = Int(eulerphi(ZZ(k)))
        c = cyclotomic(k, z)
        for y in P(2), x in P(y+1)
            N == y^e*subst(c, QQ(x, y)) && return true
        end
    end
    return false
end

"""
Filter the integers which are A299928 and <= ``n``.
"""
F299928(n) = (j for j in 1:n if isA299928(j))

"""
Iterate over the first ``n`` integers which are A299928.
"""
I299928(n) = takeFirst(isA299928, n)

"""
Return the list of the first ``n`` integers which are A299928.
"""
L299928(n) = collect(I299928(n))

#A299929########################################################################

"""
Is ``n`` a prime represented by a cyclotomic binary form f(x, y) where x and y are prime numbers and 0 < y < x?
"""
isA299929(n) = isPrime(n) && isA299928(n)

"""
Filter the integers which are A299929 and <= ``n``.
"""
F299929(n) = (j for j in 1:n if isA299929(j))

"""
Iterate over the first ``n`` integers which are A299929.
"""
I299929(n) = takeFirst(isA299929, n)

"""
Return the list of the first ``n`` integers which are A299929.
"""
L299929(n) = collect(I299929(n))

#A299930########################################################################

"""
Is ``n`` a prime represented by a cyclotomic binary form f(x, y) with x and y odd prime numbers and x > y.
"""
function isA299930(n)
    !isprime(ZZ(n)) && return false
    R, z = PolynomialRing(ZZ, "z")
    K = floor(Int, 5.383*log(n)^1.161)
    M = floor(Int, 2*sqrt(n/3))
    N = QQ(n)
    P(u) = (p for p in u:M if isprime(ZZ(p)))
    for k in 3:K
        e = Int(eulerphi(ZZ(k)))
        c = cyclotomic(k, z)
        for y in P(3), x in P(y+2)
            N == y^e*subst(c, QQ(x, y)) && return true
    end end
    return false
end

"""
Filter the integers which are A299930 and <= ``n``.
"""
F299930(n) = (j for j in 1:n if isA299930(j))

"""
Iterate over the first ``n`` integers which are A299930.
"""
I299930(n) = takeFirst(isA299930, n)

"""
Return the list of the first ``n`` integers which are A299930.
"""
L299930(n) = collect(I299930(n))

#A296095########################################################################

"""
Is ``n`` represented by a cyclotomic binary form?
"""
function isA296095(n)
    n < 3 && return false
    R, z = PolynomialRing(ZZ, "z")
    N = QQ(n)
    # Bounds from Fouvry & Levesque & Waldschmidt
    logn = log(n)^1.161
    K = floor(Int, 5.383*logn)
    M = floor(Int, 2*(n/3)^(1/2))
    k = 3
    while true
        c = cyclotomic(k, z)
        e = Int(eulerphi(ZZ(k)))
        if k == 7
            K = ceil(Int, 4.864*logn)
            M = ceil(Int, 2*(n/11)^(1/4))
        end
        for y in 2:M, x in 1:y
            N == y^e*subst(c, QQ(x, y)) && return true
        end
        k += 1
        k > K && break
    end
    return false
end

"""
Filter the integers which are A296095 and <= ``n``.
"""
F296095(n) = (j for j in 1:n if isA296095(j))

"""
Iterate over the first ``n`` integers which are A296095.
"""
I296095(n) = takeFirst(isA296095, n)

"""
Return the list of the first ``n`` integers which are A296095.
"""
L296095(n) = collect(I296095(n))

#A293654########################################################################

"""
Is ``n`` unrepresentable by a cyclotomic binary forms?
"""
function isA293654(n)
    if n < 3 return true end
    R, x = PolynomialRing(ZZ, "x")
    K = floor(Int, 5.383*log(n)^1.161) # Bounds from
    M = floor(Int, 2*sqrt(n/3)) # Fouvry & Levesque & Waldschmidt
    N = QQ(n); count = 0

    for k in 3:K
        e = Int(eulerphi(ZZ(k)))
        c = cyclotomic(k, x)
        for m in 1:M, j in 0:M if max(j, m) > 1
            N == m^e*subst(c, QQ(j,m)) && return false
    end end end
    return true
end

"""
Filter the integers which are A293654 and <= ``n``.
"""
F293654(n) = (j for j in 1:n if isA293654(j))

"""
Iterate over the first ``n`` integers which are A293654.
"""
I293654(n) = takeFirst(isA293654, n)

"""
Return the list of the first ``n`` integers which are A293654.
"""
L293654(n) = collect(I293654(n))

#A325143########################################################################

"""
Is ``n`` a prime represented by a cyclotomic binary form?
"""
function isA325143(n)
    (n < 3 || !isprime(ZZ(n))) && return false
    R, x = PolynomialRing(ZZ, "x")
    K = floor(Int, 5.383*log(n)^1.161) # Bounds from
    M = floor(Int, 2*sqrt(n/3)) # Fouvry & Levesque & Waldschmidt
    N = QQ(n)

    for k in 3:K
        e = Int(eulerphi(ZZ(k)))
        c = cyclotomic(k, x)
        for m in 1:M, j in 0:M if max(j, m) > 1
            N == m^e*subst(c, QQ(j,m)) && return true
    end end end
    return false
end

"""
Filter the integers which are A325143 and <= ``n``.
"""
F325143(n) = (j for j in 1:n if isA325143(j))

"""
Iterate over the first ``n`` integers which are A325143.
"""
I325143(n) = takeFirst(isA325143, n)

"""
Return the list of the first ``n`` integers which are A325143.
"""
L325143(n) = collect(I325143(n))

#A325145########################################################################

"""
Is ``n`` a prime unrepresentable by a cyclotomic binary form?
"""
isA325145(n) = isprime(ZZ(n)) && ! isA325143(n)

"""
Filter the integers which are A325145 and <= ``n``.
"""
F325145(n) = (j for j in 1:n if isA325145(j))

"""
Iterate over the first ``n`` integers which are A325145.
"""
I325145(n) = takeFirst(isA325145, n)

"""
Return the list of the first ``n`` integers which are A325145.
"""
L325145(n) = collect(I325145(n))

#A299214########################################################################

"""
Return the number of representations of ``n`` by cyclotomic binary forms.
"""
function V299214(n)
    if n < 3 return 0 end
    R, x = PolynomialRing(ZZ, "x")
    K = floor(Int, 5.383*log(n)^1.161) # Bounds from
    M = floor(Int, 2*sqrt(n/3))  # Fouvry & Levesque & Waldschmidt
    N = QQ(n); count = 0

    for k in 3:K
        e = Int(eulerphi(ZZ(k)))
        c = cyclotomic(k, x)
        for m in 1:M, j in 0:M if max(j, m) > 1
            N == m^e*subst(c, QQ(j,m)) && (count += 1)
    end end end
    4*count
end

"""
Return the initial list of V299214 of length len .
"""
L299214(len) = [V299214(i) for i in 1:len]

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

"""
Is n a prime represented in more than one way by cyclotomic binary forms f(x,y) with x and y prime numbers and y < x?
"""
function isA299733(n)
    if n < 3 || !isprime(ZZ(n)) return false end
    R, x = PolynomialRing(ZZ, "x")
    K = floor(Int, 5.383*log(n)^1.161) # Bounds from
    M = floor(Int, 2*sqrt(n/3)) # Fouvry & Levesque & Waldschmidt
    N = QQ(n); multi = 0

    for k in 3:K
        e = Int(eulerphi(ZZ(k)))
        c = cyclotomic(k, x)
        for m in 2:M if isprime(ZZ(m))
            for j in m:M if isprime(ZZ(j))
                if N == m^e*subst(c, QQ(j,m)) multi += 1
    end end end end end end
    multi > 1
end

"""
Return the list of the first ``n`` integers which are A299733.
"""
L299733(n) = [i for i in 1:n if isA299733(i)]

# L299733(1000) |> println

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

using Test, SeqTests

function test()
    @testset "CyclotomicForms" begin
        if is_oeis_installed()

            L = [L206942, L206864, L299498, L299928, L299929, L299930,
                 L299214, L296095, L293654]
            SeqTest(L, 'L')
        end
    end
end

function demo()
    println("\nIntegers which are A206942 and <= 32:")
    for n in 1:32
        isA206942(n) && println(n)
    end

    println("\nIterate over the first 32 integers which are A206942:")
    for i in I206942(32)
        print("$i, ")
    end

    println("\n\nFilter the integers which are A206942 and <= 60:")
    F206942(60) |> println

    println("\nReturn the list of the first 32 integers which are A206942:")
    L206942(32) |> println
end

"""
L299214(100) :: 1.822155 seconds (7.22 M allocations: 185.560 MiB, 36.76% gc time)
"""
function perf()
    @time L299214(100)
end

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

main()

end # module


text:
# Autogenerated wrapper script for Libgcrypt_jll for x86_64-w64-mingw32
export libgcrypt

using Libgpg_error_jll
JLLWrappers.@generate_wrapper_header("Libgcrypt")
JLLWrappers.@declare_library_product(libgcrypt, "libgcrypt-20.dll")
function __init__()
    JLLWrappers.@generate_init_header(Libgpg_error_jll)
    JLLWrappers.@init_library_product(
        libgcrypt,
        "bin\\libgcrypt-20.dll",
        RTLD_LAZY | RTLD_DEEPBIND,
    )

    JLLWrappers.@generate_init_footer()
end  # __init__()


text:
function lglnodes(N,a,b)

# %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
# %
# % lglnodes.m
# %
# % Computes the Legendre-Gauss-Lobatto nodes, weights and the LGL Vandermonde
# % matrix. The LGL nodes are the zeros of (1-x^2)*P'_N(x). Useful for numerical
# % integration and spectral methods.
# %
# % Reference on LGL nodes and weights:
# %   C. Canuto, M. Y. Hussaini, A. Quarteroni, T. A. Tang, "Spectral Methods
# %   in Fluid Dynamics," Section 2.3. Springer-Verlag 1987
# %
# % Written by Greg von Winckel - 04/17/2004
# % Contact: gregvw@chtm.unm.edu
# %
# %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
#
# % Truncation + 1
	N1=N+1;

	#% Use the Chebyshev-Gauss-Lobatto nodes as the first guess
	x=cos(pi*[0:N]/N);

	# The Legendre Vandermonde Matrix
	P=zeros(N1,N1);

	# 	% Compute P_(N) using the recursion relation
	# 	% Compute its first and second derivatives and
	# 	% update x using the Newton-Raphson method.

	xold=2;

	while maximum(abs(x-xold))>eps()

	    xold=x;

	    P[:,1]=1;    P[:,2]=x;

	    for k=2:N
		P[:,k+1]=( (2*k-1)*x.*P[:,k]-(k-1)*P[:,k-1] )/k;
	    end

	    x=xold-( x.*P[:,N1]-P[:,N] )./( N1*P[:,N1] );

	end

	w=2./(N*N1*P[:,N1].^2)/2.*(b-a);
	x = (x+1)*(b-a)/2.+a
	return x,w,P
end

max_stars_repo_name: flavioluiz/PortHamiltonian

max_issues_repo_name: flavioluiz/PortHamiltonian

max_forks_repo_name: flavioluiz/PortHamiltonian


text:
using DimensionalData, Test
using DimensionalData.LookupArrays, DimensionalData.Dimensions

using DimensionalData.LookupArrays: _set

a = [1 2; 3 4]
dimz = (X(143.0:2.0:145.0; lookup=Sampled(order=ForwardOrdered()), metadata=Metadata(Dict(:meta => "X"))),
        Y(-38.0:2.0:-36.0; lookup=Sampled(order=ForwardOrdered()), metadata=Metadata(Dict(:meta => "Y"))))
da = DimArray(a, dimz; name=:test)

a2 = [1 2 3 4
      3 4 5 6
      4 5 6 7]
dimz2 = (Dim{:row}(10.0:10.0:30.0), Dim{:column}(-2:1.0:1.0))
da2 = DimArray(a2, dimz2; name=:test2)

s = DimStack(da2, DimArray(2a2, dimz2; name=:test3))

@testset "Array fields" begin
    @test parent(set(da2, fill(9, 3, 4))) == fill(9, 3, 4)
    # A differently sized array can't be set
    @test_throws ArgumentError parent(set(da2, [9 9; 9 9])) == [9 9; 9 9]
end

@testset "DimStack fields" begin
    @test_throws ArgumentError set(s, (x=a2, y=3a2))
    @test_throws ArgumentError set(s, (test2=a2, test3=hcat(a2, a2)))
    s2 = set(s, (test2=a2, test3=3a2))
    @test keys(s2) == (:test2, :test3)
    @test values(s2) == (a2, 3a2)
end

@testset "DimStack Dimension" begin
    @test typeof(dims(set(s, row=X, column=Z))) <: Tuple{<:X,<:Z}
    @test typeof(dims(set(s, row=X(), column=Z()))) <: Tuple{<:X,<:Z}
    @test typeof(dims(set(s, row=:row2, column=:column2))) <: Tuple{<:Dim{:row2},<:Dim{:column2}}
    @test typeof(dims(set(s, :column => Ti(), :row => Z))) <: Tuple{<:Z,<:Ti}
    @test typeof(dims(set(s, Dim{:row}(Y()), Dim{:column}(X())))) <: Tuple{<:Y,<:X}
    @test typeof(dims(set(s, (Dim{:row}(Y), Dim{:column}(X))))) <: Tuple{<:Y,<:X}
    @test index(set(s, Dim{:row}([:x, :y, :z])), :row) == [:x, :y, :z]
end

@testset "DimArray dim Dimension" begin
    @test typeof(dims(set(da, X=:a, Y=:b))) <: Tuple{<:Dim{:a},<:Dim{:b}}
    @test typeof(dims(set(da2, Dim{:row}(Y()), Dim{:column}(X())))) <: Tuple{<:Y,<:X}
    @test typeof(dims(set(da, X => Ti(), Y => Z()))) <: Tuple{<:Ti,<:Z}
    @test typeof(dims(set(da2, :column => Ti(), :row => Z()))) <: Tuple{<:Z,<:Ti}
    @test typeof(dims(set(da2, (Dim{:row}(Y), Dim{:column}(X))))) <: Tuple{<:Y,<:X}
    @test typeof(dims(set(da2, row=X, column=Z))) <: Tuple{<:X,<:Z}
    @test typeof(dims(set(da2, row=X(), column=Z()))) <: Tuple{<:X,<:Z}
    @test typeof(dims(set(da2, row=:row2, column=:column2))) <: Tuple{<:Dim{:row2},<:Dim{:column2}}
    @test index(set(da2, Dim{:row}([:x, :y, :z])), :row) == [:x, :y, :z]
end

@testset "Dimension index" begin
    @test index(set(da2, :column => [:a, :b, :c, :d], :row => 4:6)) ==
        (4:6, [:a, :b, :c, :d])
    @test index(set(s, :column => 10:5:20, :row => 4:6)) == (4:6, 10:5:20)
    @test step.(span(set(da2, :column => 10:5:20, :row => 4:6))) == (1, 5)
end

@testset "dim lookup" begin
    @test lookup(set(da2, :column => NoLookup(), :row => Sampled(sampling=Intervals(Center())))) ==
        (Sampled(10.0:10.0:30.0, ForwardOrdered(), Regular(10.0), Intervals(Center()), NoMetadata()), NoLookup(Base.OneTo(4)))
    @test lookup(set(da2, column=NoLookup())) ==
        (Sampled(10.0:10.0:30.0, ForwardOrdered(), Regular(10.0), Points(), NoMetadata()), NoLookup(Base.OneTo(4)))
    @test lookup(set(da2, :column => NoLookup(), :row => Sampled())) ==
        (Sampled(10.0:10.0:30.0, ForwardOrdered(), Regular(10.0), Points(), NoMetadata()), NoLookup(Base.OneTo(4)))
    cat_da = set(da, X=NoLookup(), Y=Categorical())
    @test index(cat_da) ==
        (NoLookup(Base.OneTo(2)), Categorical(-38.0:2.0:-36.0, Unordered(), NoMetadata()))
    cat_da_m = set(dims(cat_da, Y), X(DimensionalData.AutoIndex(); metadata=Dict()))
    @test metadata(cat_da_m) == Dict()

    @testset "span" begin
        # TODO: should this error? the span step doesn't match the index step
        @test span(set(da2, row=Irregular(10, 12), column=Regular(9.9))) ==
            (Irregular(10, 12), Regular(9.9))
        @test _set(Sampled(), AutoSpan()) == Sampled()
        @test _set(Sampled(), Irregular()) == Sampled(; span= Irregular())
    end

    interval_da = set(da, X=Intervals(), Y=Intervals())
    @testset "locus" begin
        @test_throws ArgumentError set(da2, (End(), Center()))
        @test locus(set(interval_da, X(End()), Y(Center()))) == (End(), Center())
        @test locus(set(interval_da, X=>End(), Y=>Center())) == (End(), Center())
        @test locus(set(interval_da, X=End, Y=Center)) == (End(), Center())
        @test locus(set(da, Y=Center())) == (Center(), Center())
        @test _set(Intervals(Center()), Start()) == Intervals(Start())
        @test _set(Intervals(Center()), AutoLocus()) == Intervals(Center())
        @test _set(Points(), Center()) == Points()
        @test_throws ArgumentError _set(Points(), Start())
        @test_throws ArgumentError _set(Points(), End())
    end

    @testset "sampling" begin
        @test sampling(interval_da) == (Intervals(), Intervals())
        @test sampling(set(da, (X(Intervals(End())), Y(Intervals(Start()))))) ==
            (Intervals(End()), Intervals(Start()))
        @test _set(Sampled(), AutoSampling()) == Sampled()
        @test _set(Sampled(), Intervals()) == Sampled(; sampling=Intervals())
        @test _set(Points(), AutoSampling()) == Points()
        @test _set(AutoSampling(), Intervals()) == Intervals()
        @test _set(AutoSampling(), AutoSampling()) == AutoSampling()
    end

    @testset "order" begin
        uda = set(da, Y(Unordered()))
        @test order(uda) == (ForwardOrdered(), Unordered())
        @test order(set(uda, X=ReverseOrdered())) == (ReverseOrdered(), Unordered())
    end

end

@testset "metadata" begin
    @test metadata(set(Sampled(), Metadata(Dict(:a=>1, :b=>2)))).val == Dict(:a=>1, :b=>2)
    dax = set(da, X => Metadata(Dict(:a=>1, :b=>2)))
    @test metadata(dims(dax), X).val == Dict(:a=>1, :b=>2)
    @test metadata(dims(dax), Y).val == Dict(:meta => "Y")
    dax = set(da, X => Metadata(Dict(:a=>1, :b=>2)))
    @test metadata(dims(dax, X)).val == Dict(:a=>1, :b=>2)
    dax = set(da2, row=Metadata(Dict(:a=>1, :b=>2)))
    @test metadata(dims(dax, :row)).val == Dict(:a=>1, :b=>2)
    dax = set(da2, column=Metadata(Dict(:a=>1, :b=>2)))
    @test metadata(dims(dax, :column)).val == Dict(:a=>1, :b=>2)
end

@testset "all dim fields" begin
    md = Metadata(Dict(:a=>1, :b=>2))
    dax = set(da, X(20:-10:10; metadata=md))
    x = dims(dax, X)
    @test index(x) === 20:-10:10
    @test order(x) === ReverseOrdered()
    @test span(x) === Regular(-10)
    @test lookup(x) == Sampled(20:-10:10, ReverseOrdered(), Regular(-10), Points(), md)
    @test metadata(x).val == Dict(:a=>1, :b=>2)
end

@testset "errors with set" begin
    @test_throws ArgumentError set(da, Sampled())
    @test_throws ArgumentError set(s, Categorical())
    @test_throws ArgumentError set(da, Irregular())
    @test_throws ArgumentError set(da, X=7)
    @test_throws ArgumentError _set(dims(da, X), X(7))
    @test_throws ArgumentError set(da, notadimname=Sampled())
end

@testset "_set nothing" begin
    @test _set(nothing, nothing) == nothing
    @test _set(1, nothing) == 1
    @test _set(nothing, 2) == 2
end


text:
export IWrapper, Partial, unwrap, value

"""
    value(x)

Get the `value` from a wrapper instance.
"""
value(x) = x
NiLangCore.chfield(x::T, ::typeof(value), y::T) where T = y

"""
    IWrapper{T} <: Real

IWrapper{T} is a wrapper of for data of type T.
It will forward `>, <, >=, <=, ≈` operations.
"""
abstract type IWrapper{T} <: Real end
NiLangCore.chfield(x, ::Type{T}, v) where {T<:IWrapper} = (~T)(v)
Base.eps(::Type{<:IWrapper{T}}) where T = Base.eps(T)

"""
    unwrap(x)

Unwrap a wrapper instance (recursively) to get the content value.
"""
unwrap(x::IWrapper) = unwrap(value(x))
unwrap(x) = x

for op in [:>, :<, :>=, :<=, :isless, :(==)]
    @eval Base.$op(a::IWrapper, b::IWrapper) = $op(unwrap(a), unwrap(b))
    @eval Base.$op(a::IWrapper, b::Real) = $op(unwrap(a), b)
    @eval Base.$op(a::IWrapper, b::AbstractFloat) = $op(unwrap(a), b)
    @eval Base.$op(a::Real, b::IWrapper) = $op(a, unwrap(b))
    @eval Base.$op(a::AbstractFloat, b::IWrapper) = $op(a, unwrap(b))
end

"""
Partial{FIELD, T, T2} <: IWrapper{T2}

Take a field `FIELD` without dropping information.
This operation can be undone by calling `~Partial{FIELD}`.
"""
struct Partial{FIELD, T, T2} <: IWrapper{T2}
    x::T
    function Partial{FIELD,T,T2}(x::T) where {T,T2,FIELD}
        new{FIELD,T,T2}(x)
    end
    function Partial{FIELD,T,T2}(x::T) where {T<:Complex,T2,FIELD}
        new{FIELD,T,T2}(x)
    end
end
Partial{FIELD}(x::T) where {T,FIELD} = Partial{FIELD,T,typeof(getfield(x,FIELD))}(x)
Partial{FIELD}(x::T) where {T<:Complex,FIELD} = Partial{FIELD,T,typeof(getfield(x,FIELD))}(x)

@generated function (_::Type{Inv{Partial{FIELD}}})(x::Partial{FIELD}) where {FIELD}
    :(x.x)
end

function NiLangCore.chfield(hd::Partial{FIELD}, ::typeof(value), val) where FIELD
    chfield(hd, Val(:x), chfield(hd.x, Val(FIELD), val))
end

@generated function value(hv::Partial{FIELD}) where FIELD
    :(hv.x.$FIELD)
end

function Base.zero(x::T) where T<:Partial
    zero(T)
end

function Base.zero(x::Type{<:Partial{FIELD,T}}) where {FIELD, T}
    Partial{FIELD}(Base.zero(T))
end
Base.show(io::IO, gv::Partial{FIELD}) where FIELD = print(io, "$(gv.x).$FIELD")
Base.show(io::IO, ::MIME"plain/text", gv::Partial) = Base.show(io, gv)


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

mutable struct ProviderResourceType <: SwaggerModel
    resourceType::Any # spec type: Union{ Nothing, String } # spec name: resourceType
    locations::Any # spec type: Union{ Nothing, Vector{String} } # spec name: locations
    aliases::Any # spec type: Union{ Nothing, Vector{Alias} } # spec name: aliases
    apiVersions::Any # spec type: Union{ Nothing, Vector{String} } # spec name: apiVersions
    capabilities::Any # spec type: Union{ Nothing, String } # spec name: capabilities
    properties::Any # spec type: Union{ Nothing, Dict{String, String} } # spec name: properties

    function ProviderResourceType(;resourceType=nothing, locations=nothing, aliases=nothing, apiVersions=nothing, capabilities=nothing, properties=nothing)
        o = new()
        validate_property(ProviderResourceType, Symbol("resourceType"), resourceType)
        setfield!(o, Symbol("resourceType"), resourceType)
        validate_property(ProviderResourceType, Symbol("locations"), locations)
        setfield!(o, Symbol("locations"), locations)
        validate_property(ProviderResourceType, Symbol("aliases"), aliases)
        setfield!(o, Symbol("aliases"), aliases)
        validate_property(ProviderResourceType, Symbol("apiVersions"), apiVersions)
        setfield!(o, Symbol("apiVersions"), apiVersions)
        validate_property(ProviderResourceType, Symbol("capabilities"), capabilities)
        setfield!(o, Symbol("capabilities"), capabilities)
        validate_property(ProviderResourceType, Symbol("properties"), properties)
        setfield!(o, Symbol("properties"), properties)
        o
    end
end # type ProviderResourceType

const _property_map_ProviderResourceType = Dict{Symbol,Symbol}(Symbol("resourceType")=>Symbol("resourceType"), Symbol("locations")=>Symbol("locations"), Symbol("aliases")=>Symbol("aliases"), Symbol("apiVersions")=>Symbol("apiVersions"), Symbol("capabilities")=>Symbol("capabilities"), Symbol("properties")=>Symbol("properties"))
const _property_types_ProviderResourceType = Dict{Symbol,String}(Symbol("resourceType")=>"String", Symbol("locations")=>"Vector{String}", Symbol("aliases")=>"Vector{Alias}", Symbol("apiVersions")=>"Vector{String}", Symbol("capabilities")=>"String", Symbol("properties")=>"Dict{String, String}")
Base.propertynames(::Type{ ProviderResourceType }) = collect(keys(_property_map_ProviderResourceType))
Swagger.property_type(::Type{ ProviderResourceType }, name::Symbol) = Union{Nothing,eval(Base.Meta.parse(_property_types_ProviderResourceType[name]))}
Swagger.field_name(::Type{ ProviderResourceType }, property_name::Symbol) =  _property_map_ProviderResourceType[property_name]

function check_required(o::ProviderResourceType)
    true
end

function validate_property(::Type{ ProviderResourceType }, name::Symbol, val)
end


text:
socket = connect("localhost",256)
write(socket, "hello socket world")
close(socket)

max_stars_repo_name: ethansaxenian/RosettaDecode

max_issues_repo_name: ethansaxenian/RosettaDecode

max_forks_repo_name: ethansaxenian/RosettaDecode


text:
"""
    JoeCopula{d,T}

Fields:
  - θ::Real - parameter

Constructor

    JoeCopula(d, θ)

The [Joe](https://en.wikipedia.org/wiki/Copula_(probability_theory)#Most_important_Archimedean_copulas) copula in dimension ``d`` is parameterized by ``\\theta \\in [1,\\infty)``. It is an Archimedean copula with generator :

```math
\\phi(t) = 1 - \\left(1 - e^{-t}\\right)^{\\frac{1}{\\theta}}
```
"""
struct JoeCopula{d,T} <: ArchimedeanCopula{d}
    θ::T
end
JoeCopula(d,θ) = θ >= 1 ? JoeCopula{d,typeof(θ)}(θ) : @error "Theta must be greater than one."
ϕ(  C::JoeCopula,          t) = 1-(1-exp(-t))^(1/C.θ)
ϕ⁻¹(C::JoeCopula,          t) = -log(1-(1-t)^C.θ)
τ(C::JoeCopula) = 1 - 4sum(1/(k*(2+k*C.θ)*(C.θ*(k-1)+2)) for k in 1:1000) # 446 in R copula.
radial_dist(C::JoeCopula) = Sibuya(1/C.θ)


text:
"""
    IsingDistribution

Represents the Ising distribution ``P(x) = 1/Z exp(x' J x)``, where ``J`` is a lower
triangular matrix.

"""
mutable struct IsingDistribution <: AbstractBinaryVectorDistribution
    # W # these are the n(n+1)/2 "true parameters."
    J::LowerTriangular # This is the lower triangular matrix used in calculations
    # I::Vector{Int} # moving to metadata
    metadata::Dict{Any,Any}
    cache::Dict{Any,Any}

    function IsingDistribution(J; kwargs...)
        if size(J,1) != size(J,2)
            error("IsingDistribution: J must be square")
        end
        Jsym = (J + J') / 2
        n = size(J,1)
        new(LowerTriangular(Jsym), Dict(kwargs), Dict(:pdf=>spzeros(2^n), :F=>spzeros(2^n)))
    end
    function IsingDistribution(W::Vector, kwargs...)
        IsingDistribution(v2lt(W), kwargs...)
    end
    # function IsingDistribution(J::Matrix, theta::Vector; kwargs...)
    #     if size(J,1) != size(J,2)
    #         error("Ising Distribution: J must be square")
    #     end
    #     if size(J,1) != length(theta)
    #         error("Ising Distribution: J and theta must have compatible size")
    #     end
    #     Jsym = (J + J') / 2
    #     Jsym = Jsym - Diagonal(Jsym)
    #     new(Jsym, theta, Dict(kwargs), Dict(:pdf=>spzeros(2^length(theta)), :energy=>spzeros(2^length(theta))))
    # end
    # function IsingDistribution(Jtilde::Matrix; kwargs...)
    #     theta = diag(Jtilde)
    #     Jnodiag = Jtilde - Diagonal(Jtilde)
    #     IsingDistribution(Jnodiag, theta; kwargs...)
    # end
    # function
end

n_bits(P::IsingDistribution) = size(P.J,1)

function show(io::IO, P::IsingDistribution)
    println(io, "Ising Distribution")
    println(io, "N_neurons: $(n_bits(P))")
    # println(io, "Indices:   $(P.I)")
    show_metadata(io, P)
end

# ==(P1::IsingDistribution, P2::IsingDistribution) = (P1.J == P2.J) && (P1.theta == P2.theta)
==(P::IsingDistribution, Q::IsingDistribution) = P.J == Q.J

################################################################################
#### Miscellaneous computations/internal functions
################################################################################
"""
    initialize_cache!(ID)

Checks if `ID.cache[:pdf]` and `ID.cache[:F]` exist; if not, initializes those values to empty sparse vectors of length `2^n_bits(ID)`
"""
function initialize_cache!(ID)
    if !haskey(ID.cache, :pdf)
        ID.cache[:pdf] = spzeros(2^n_bits(ID))
    end
    if !haskey(ID.cache, :F)
        ID.cache[:F] = spzeros(2^n_bits(ID))
    end
end

"""
    F_Ising(ID, x)

The "pseudo-energy" of state `x`, given by `x' J x`. Stores value in cache.
"""
function F_Ising(ID::IsingDistribution, x::AbstractVector{Bool})
    idx = 1 + _binary_to_int(x)
    if ID.cache[:F][idx] == 0.0
        ID.cache[:F][idx] = dot(x, ID.J * x)
    end
    return ID.cache[:F][idx]
end
# """
#     _get_energies(ID)
#
# Returns an array of the Ising energy of each state.
# """
# _get_energies(ID::IsingDistribution) = [_E_Ising(ID, digits(Bool,x,2,n_bits(ID))) for x in 0:(2^n_bits(ID) - 1)]
get_F(ID::IsingDistribution) = [F_Ising(ID, digits(Bool,x,2,n_bits(ID))) for x in 0:(2^n_bits(ID) - 1)]
function get_Z(ID::IsingDistribution)
    if !haskey(ID.cache, :Z)
        if n_bits(ID) > ISING_METHOD_THRESHOLD
            warn("Computing Z for Ising PDF with $(n_bits(ID)) neurons. This might take a while, and numerical accuracy is not guaranteed. (Warning only displays once)", once=true, key=METH_THRESH_WARN_KEY)
        end
        # ID.cache[:Z] = sum_kbn([exp(-_E_Ising(ID, digits(Bool, k, 2, n_bits(ID)))) for k = 0:(2^n_bits(ID) - 1)])
        ID.cache[:Z] = sum_kbn(exp.(get_F(ID)))
    end
    return ID.cache[:Z]
end

################################################################################
#### Implementation of methods for distributions
################################################################################
function entropy(P::IsingDistribution)
    if !haskey(P.metadata, :entropy)
        energies = get_F(P)
        # P.metadata[:entropy] = log(get_Z(P)) + sum_kbn([exp(-_E_Ising(P,digits(Bool,2,x,n_bits(P)))) for x in 0:(2^n_bis(P) - 1)])
        P.metadata[:entropy] = log(get_Z(P)) - sum_kbn(exp.(energies) .* energies) / get_Z(P)
    end
    return P.metadata[:entropy]
end
function entropy2(P::IsingDistribution)
    if !haskey(P.metadata, :entropy2)
        energies = get_F(P)
        P.metadata[:entropy2] = log2(get_Z(P)) - sum_kbn(exp.(energies) .* energies) / (get_Z(P) * log2(e))
    end
    return P.metadata[:entropy2]
end

function pdf(ID::IsingDistribution, x)
    if length(x) != n_bits(ID)
        error("IsingDistribution pdf: out of domain error")
    else
        idx = 1 + _binary_to_int(x)
        if ID.cache[:pdf][idx] == 0.0
            ID.cache[:pdf][idx] = exp(F_Ising(ID,x)) / get_Z(ID)
        end
        return ID.cache[:pdf][idx]
    end
end

# function random(ID::IsingDistribution, n_samples=1, force_gibbs=false)
#     if n_bits(ID) <= ISING_METHOD_THRESHOLD && !force_gibbs
#         return _random_exact(ID, n_samples)
#     else
#         return _random_gibbs(ID, n_samples)
#     end
# end
function random(ID::IsingDistribution, n_samples=1)
    burnin = _get_burnin(ID)
    ind_steps = _get_ind_steps(ID)
    N_steps = burnin + (n_samples - 1) * ind_steps

    d_i = rand(1:n_bits(ID), N_steps) # sequence of bits to update
    u_r = rand(N_steps) # sequence of uniform random samples from [0,1]
    X_out = falses(n_bits(ID), n_samples)
    x = rand(Bool, n_bits(ID)) # the current state

    sample = 1
    next_sample_step = burnin

    for step = 1:N_steps
        # compute the conditional probability of bit d_i[step] given the other bits
        # E_act = dot(x, -0.5*ID.J[d_i[step],:]) + ID.theta[d_i[step]]
        E_act = F_Ising(ID, x)
        P_act = 1 / (1 + exp(E_act))
        x[d_i[step]] = P_act < u_r[step]
        if step == next_sample_step
            X_out[:,sample] = x
            sample += 1
            next_sample_step += ind_steps
        end
    end
    return X_out
end
_get_burnin(ID::IsingDistribution) = get!(ID.metadata, :burnin, 100*n_bits(ID))
_get_ind_steps(ID::IsingDistribution) = get!(ID.metadata, :ind_steps, 10*n_bits(ID))

"""
    lt2v(V, n)

Return the lower triangular matrix `L` such that `L[:]`, skipping the entries above the
diagonal, is `V`.
"""
function v2lt(V, n=Int((-1 + sqrt(1 + 8*length(V)))/2))
    v = start(V)
    L = zeros(n,n)
    for i = 1:n
        for j = i:n
            L[j,i], v = next(V, v)
        end
    end
    return LowerTriangular(L)
end

lt2v(L) = vcat([[L[j,i] for j = i:size(L,1)] for i = 1:size(L,2)]...)

max_stars_repo_name: sekunder/BinaryVectorProbability

max_issues_repo_name: sekunder/BinaryVectorProbability

max_forks_repo_name: sekunder/BinaryVectorProbability