text:
#=
--- Day 7: Handy Haversacks ---

You land at the regional airport in time for your next flight. In fact, it looks like you'll even have time to grab some food: all flights are currently delayed due to issues in luggage processing.

Due to recent aviation regulations, many rules (your puzzle input) are being enforced about bags and their contents; bags must be color-coded and must contain specific quantities of other color-coded bags. Apparently, nobody responsible for these regulations considered how long they would take to enforce!

For example, consider the following rules:

light red bags contain 1 bright white bag, 2 muted yellow bags.
dark orange bags contain 3 bright white bags, 4 muted yellow bags.
bright white bags contain 1 shiny gold bag.
muted yellow bags contain 2 shiny gold bags, 9 faded blue bags.
shiny gold bags contain 1 dark olive bag, 2 vibrant plum bags.
dark olive bags contain 3 faded blue bags, 4 dotted black bags.
vibrant plum bags contain 5 faded blue bags, 6 dotted black bags.
faded blue bags contain no other bags.
dotted black bags contain no other bags.

These rules specify the required contents for 9 bag types. In this example, every faded blue bag is empty, every vibrant plum bag contains 11 bags (5 faded blue and 6 dotted black), and so on.

You have a shiny gold bag. If you wanted to carry it in at least one other bag, how many different bag colors would be valid for the outermost bag? (In other words: how many colors can, eventually, contain at least one shiny gold bag?)

In the above rules, the following options would be available to you:

    A bright white bag, which can hold your shiny gold bag directly.
    A muted yellow bag, which can hold your shiny gold bag directly, plus some other bags.
    A dark orange bag, which can hold bright white and muted yellow bags, either of which could then hold your shiny gold bag.
    A light red bag, which can hold bright white and muted yellow bags, either of which could then hold your shiny gold bag.

So, in this example, the number of bag colors that can eventually contain at least one shiny gold bag is 4.

How many bag colors can eventually contain at least one shiny gold bag? (The list of rules is quite long; make sure you get all of it.)
=#

"""
    Array[
            bag_color => [(Int,bag_color), (Int,bag_color), ...],
            bag_color => [(Int,bag_color), (Int,bag_color), ...],
            bag_color => [(Int,bag_color), (Int,bag_color), ...]
        ]
"""

function read_input()
    rules = []
    open("./input", "r") do f
        while ! eof(f)
            s::String = readline(f)
            bag_color_rule = split(s," bags contain ")[1]
            contained_bags = split(s," bags contain ")[2]

            list_contained = split(contained_bags,",")
            dict_value = []
            for contained in list_contained

                if contained[1] == ' '
                    contained = contained[2:end]
                end

                # gambi because muita mão
                qtd =  match(r"[\d]*",contained).match

                bag_quantity = tryparse(Int,qtd)
                if bag_quantity == nothing
                    break end

                bag_color = match(r"[^\d]*.?(?= bag)",contained).match[2:end]
                push!(dict_value,(bag_quantity,bag_color))
            end

            push!(rules, Dict(bag_color_rule => dict_value))
        end
    end
    return rules
end

bag_dict_list = read_input()
parent_bags = Set()

function find_parent_bags(bag_color)
    for rule in bag_dict_list
        for (bag, tuples) in rule   # rule has only one element, but i used a for loop in order to extract its key,value pair
            for t in tuples
                if bag_color == t[2]
                    push!(parent_bags, bag)
                    find_parent_bags(bag)
                    break
                end
            end
        end
    end
end

find_parent_bags("shiny gold")
println("There are ", length(parent_bags), " bags that can eventually contain shiny gold bags.")

function check_if_rule_repeat()
    colors_set=Set()
    colors_list=[]
    for rule in bag_dict_list
        for (bag, tuple) in rule
            push!(colors_set, bag)
            push!(colors_list, bag)
        end
    end
    if length(colors_set) != length(colors_list)
        return true
    else
        return false
    end
end

println("Do rules reapeat? ", check_if_rule_repeat() ? "yes" : "no")
println("Number of rules: ", length(bag_dict_list))

#=
--- Part Two ---

It's getting pretty expensive to fly these days - not because of ticket prices, but because of the ridiculous number of bags you need to buy!

Consider again your shiny gold bag and the rules from the above example:

    faded blue bags contain 0 other bags.
    dotted black bags contain 0 other bags.
    vibrant plum bags contain 11 other bags: 5 faded blue bags and 6 dotted black bags.
    dark olive bags contain 7 other bags: 3 faded blue bags and 4 dotted black bags.

So, a single shiny gold bag must contain 1 dark olive bag (and the 7 bags within it) plus 2 vibrant plum bags (and the 11 bags within each of those): 1 + 1*7 + 2 + 2*11 = 32 bags!

Of course, the actual rules have a small chance of going several levels deeper than this example; be sure to count all of the bags, even if the nesting becomes topologically impractical!

Here's another example:

shiny gold bags contain 2 dark red bags.
dark red bags contain 2 dark orange bags.
dark orange bags contain 2 dark yellow bags.
dark yellow bags contain 2 dark green bags.
dark green bags contain 2 dark blue bags.
dark blue bags contain 2 dark violet bags.
dark violet bags contain no other bags.

In this example, a single shiny gold bag must contain 126 other bags.

How many individual bags are required inside your single shiny gold bag?
=#

# because no rule repeats, we can create a dict where each key is the rule
bag_dict = Dict()
function gen_bag_dict()
    for rule in bag_dict_list
        for (key,value) in rule
            bag_dict[key] = value
        end
    end
end
gen_bag_dict()

function count_bags_inside_bag(bag)

    if !haskey(bag_dict, bag)
        return 0
    end

    count = 0
    for tuple in bag_dict[bag]
        count += tuple[1] + tuple[1]*count_bags_inside_bag(tuple[2])
    end

    return count
end

count = count_bags_inside_bag("shiny gold")
println("Shiny gold bags can have ", count, " bags inside it.")


text:
"""
    append_gs!(model::JuMP.Model, x, ng::Int)

Append inequality constraints (auto-generated)

# Arguments
    `model::JuMP.Model`: optimization model
    `x::JuMP.variables`: optimization model variables
    `ng::Int`: number of inequality constraints
"""
function append_gs!(model::JuMP.Model, x, ng::Int)
    for ig = 1:ng
        if ig == 1
            @NLconstraint(model, -Inf <= g1(x...) <= 0)
        elseif ig == 2
            @NLconstraint(model, -Inf <= g2(x...) <= 0)
        elseif ig == 3
            @NLconstraint(model, -Inf <= g3(x...) <= 0)
        elseif ig == 4
            @NLconstraint(model, -Inf <= g4(x...) <= 0)
        elseif ig == 5
            @NLconstraint(model, -Inf <= g5(x...) <= 0)
        elseif ig == 6
            @NLconstraint(model, -Inf <= g6(x...) <= 0)
        elseif ig == 7
            @NLconstraint(model, -Inf <= g7(x...) <= 0)
        elseif ig == 8
            @NLconstraint(model, -Inf <= g8(x...) <= 0)
        elseif ig == 9
            @NLconstraint(model, -Inf <= g9(x...) <= 0)
        elseif ig == 10
            @NLconstraint(model, -Inf <= g10(x...) <= 0)
        elseif ig == 11
            @NLconstraint(model, -Inf <= g11(x...) <= 0)
        elseif ig == 12
            @NLconstraint(model, -Inf <= g12(x...) <= 0)
        elseif ig == 13
            @NLconstraint(model, -Inf <= g13(x...) <= 0)
        elseif ig == 14
            @NLconstraint(model, -Inf <= g14(x...) <= 0)
        elseif ig == 15
            @NLconstraint(model, -Inf <= g15(x...) <= 0)
        elseif ig == 16
            @NLconstraint(model, -Inf <= g16(x...) <= 0)
        elseif ig == 17
            @NLconstraint(model, -Inf <= g17(x...) <= 0)
        elseif ig == 18
            @NLconstraint(model, -Inf <= g18(x...) <= 0)
        elseif ig == 19
            @NLconstraint(model, -Inf <= g19(x...) <= 0)
        elseif ig == 20
            @NLconstraint(model, -Inf <= g20(x...) <= 0)
        elseif ig == 21
            @NLconstraint(model, -Inf <= g21(x...) <= 0)
        elseif ig == 22
            @NLconstraint(model, -Inf <= g22(x...) <= 0)
        elseif ig == 23
            @NLconstraint(model, -Inf <= g23(x...) <= 0)
        elseif ig == 24
            @NLconstraint(model, -Inf <= g24(x...) <= 0)
        elseif ig == 25
            @NLconstraint(model, -Inf <= g25(x...) <= 0)
        elseif ig == 26
            @NLconstraint(model, -Inf <= g26(x...) <= 0)
        elseif ig == 27
            @NLconstraint(model, -Inf <= g27(x...) <= 0)
        elseif ig == 28
            @NLconstraint(model, -Inf <= g28(x...) <= 0)
        elseif ig == 29
            @NLconstraint(model, -Inf <= g29(x...) <= 0)
        elseif ig == 30
            @NLconstraint(model, -Inf <= g30(x...) <= 0)
        elseif ig == 31
            @NLconstraint(model, -Inf <= g31(x...) <= 0)
        elseif ig == 32
            @NLconstraint(model, -Inf <= g32(x...) <= 0)
        elseif ig == 33
            @NLconstraint(model, -Inf <= g33(x...) <= 0)
        elseif ig == 34
            @NLconstraint(model, -Inf <= g34(x...) <= 0)
        elseif ig == 35
            @NLconstraint(model, -Inf <= g35(x...) <= 0)
        elseif ig == 36
            @NLconstraint(model, -Inf <= g36(x...) <= 0)
        elseif ig == 37
            @NLconstraint(model, -Inf <= g37(x...) <= 0)
        elseif ig == 38
            @NLconstraint(model, -Inf <= g38(x...) <= 0)
        elseif ig == 39
            @NLconstraint(model, -Inf <= g39(x...) <= 0)
        elseif ig == 40
            @NLconstraint(model, -Inf <= g40(x...) <= 0)
        elseif ig == 41
            @NLconstraint(model, -Inf <= g41(x...) <= 0)
        elseif ig == 42
            @NLconstraint(model, -Inf <= g42(x...) <= 0)
        elseif ig == 43
            @NLconstraint(model, -Inf <= g43(x...) <= 0)
        elseif ig == 44
            @NLconstraint(model, -Inf <= g44(x...) <= 0)
        elseif ig == 45
            @NLconstraint(model, -Inf <= g45(x...) <= 0)
        elseif ig == 46
            @NLconstraint(model, -Inf <= g46(x...) <= 0)
        elseif ig == 47
            @NLconstraint(model, -Inf <= g47(x...) <= 0)
        elseif ig == 48
            @NLconstraint(model, -Inf <= g48(x...) <= 0)
        elseif ig == 49
            @NLconstraint(model, -Inf <= g49(x...) <= 0)
        elseif ig == 50
            @NLconstraint(model, -Inf <= g50(x...) <= 0)
        elseif ig == 51
            @NLconstraint(model, -Inf <= g51(x...) <= 0)
        elseif ig == 52
            @NLconstraint(model, -Inf <= g52(x...) <= 0)
        elseif ig == 53
            @NLconstraint(model, -Inf <= g53(x...) <= 0)
        elseif ig == 54
            @NLconstraint(model, -Inf <= g54(x...) <= 0)
        elseif ig == 55
            @NLconstraint(model, -Inf <= g55(x...) <= 0)
        elseif ig == 56
            @NLconstraint(model, -Inf <= g56(x...) <= 0)
        elseif ig == 57
            @NLconstraint(model, -Inf <= g57(x...) <= 0)
        elseif ig == 58
            @NLconstraint(model, -Inf <= g58(x...) <= 0)
        elseif ig == 59
            @NLconstraint(model, -Inf <= g59(x...) <= 0)
        elseif ig == 60
            @NLconstraint(model, -Inf <= g60(x...) <= 0)
        elseif ig == 61
            @NLconstraint(model, -Inf <= g61(x...) <= 0)
        elseif ig == 62
            @NLconstraint(model, -Inf <= g62(x...) <= 0)
        elseif ig == 63
            @NLconstraint(model, -Inf <= g63(x...) <= 0)
        elseif ig == 64
            @NLconstraint(model, -Inf <= g64(x...) <= 0)
        elseif ig == 65
            @NLconstraint(model, -Inf <= g65(x...) <= 0)
        elseif ig == 66
            @NLconstraint(model, -Inf <= g66(x...) <= 0)
        elseif ig == 67
            @NLconstraint(model, -Inf <= g67(x...) <= 0)
        elseif ig == 68
            @NLconstraint(model, -Inf <= g68(x...) <= 0)
        elseif ig == 69
            @NLconstraint(model, -Inf <= g69(x...) <= 0)
        elseif ig == 70
            @NLconstraint(model, -Inf <= g70(x...) <= 0)
        elseif ig == 71
            @NLconstraint(model, -Inf <= g71(x...) <= 0)
        elseif ig == 72
            @NLconstraint(model, -Inf <= g72(x...) <= 0)
        elseif ig == 73
            @NLconstraint(model, -Inf <= g73(x...) <= 0)
        elseif ig == 74
            @NLconstraint(model, -Inf <= g74(x...) <= 0)
        elseif ig == 75
            @NLconstraint(model, -Inf <= g75(x...) <= 0)
        elseif ig == 76
            @NLconstraint(model, -Inf <= g76(x...) <= 0)
        elseif ig == 77
            @NLconstraint(model, -Inf <= g77(x...) <= 0)
        elseif ig == 78
            @NLconstraint(model, -Inf <= g78(x...) <= 0)
        elseif ig == 79
            @NLconstraint(model, -Inf <= g79(x...) <= 0)
        elseif ig == 80
            @NLconstraint(model, -Inf <= g80(x...) <= 0)
        elseif ig == 81
            @NLconstraint(model, -Inf <= g81(x...) <= 0)
        elseif ig == 82
            @NLconstraint(model, -Inf <= g82(x...) <= 0)
        elseif ig == 83
            @NLconstraint(model, -Inf <= g83(x...) <= 0)
        elseif ig == 84
            @NLconstraint(model, -Inf <= g84(x...) <= 0)
        elseif ig == 85
            @NLconstraint(model, -Inf <= g85(x...) <= 0)
        elseif ig == 86
            @NLconstraint(model, -Inf <= g86(x...) <= 0)
        elseif ig == 87
            @NLconstraint(model, -Inf <= g87(x...) <= 0)
        elseif ig == 88
            @NLconstraint(model, -Inf <= g88(x...) <= 0)
        elseif ig == 89
            @NLconstraint(model, -Inf <= g89(x...) <= 0)
        elseif ig == 90
            @NLconstraint(model, -Inf <= g90(x...) <= 0)
        elseif ig == 91
            @NLconstraint(model, -Inf <= g91(x...) <= 0)
        elseif ig == 92
            @NLconstraint(model, -Inf <= g92(x...) <= 0)
        elseif ig == 93
            @NLconstraint(model, -Inf <= g93(x...) <= 0)
        elseif ig == 94
            @NLconstraint(model, -Inf <= g94(x...) <= 0)
        elseif ig == 95
            @NLconstraint(model, -Inf <= g95(x...) <= 0)
        elseif ig == 96
            @NLconstraint(model, -Inf <= g96(x...) <= 0)
        elseif ig == 97
            @NLconstraint(model, -Inf <= g97(x...) <= 0)
        elseif ig == 98
            @NLconstraint(model, -Inf <= g98(x...) <= 0)
        elseif ig == 99
            @NLconstraint(model, -Inf <= g99(x...) <= 0)
        elseif ig == 100
            @NLconstraint(model, -Inf <= g100(x...) <= 0)
        end
    end
end

ext: jl

max_stars_repo_name: Yuricst/NLPSaUT

max_issues_repo_name: Yuricst/NLPSaUT

max_forks_repo_name: Yuricst/NLPSaUT


text:
wrap_api_call(spec::SpecFunc, args; with_func_ptr = false) = wrap_return(
    :($(spec.name)($((with_func_ptr ? [args; first(func_ptrs(spec))] : args)...))),
    spec.return_type,
    nice_julian_type(spec.return_type),
)

function wrap_enumeration_api_call(spec::SpecFunc, exs::Expr...; free = [])
    if must_repeat_while_incomplete(spec)
        if !isempty(free)
            free_block = quote
                if _return_code == VK_INCOMPLETE
                    $(map(x -> :(Libc.free($(x.name))), free)...)
                end
            end
            [:(@repeat_while_incomplete $(Expr(:block, exs..., free_block)))]
        else
            [:(@repeat_while_incomplete $(Expr(:block, exs...)))]
        end
    else
        exs
    end
end

function wrap(spec::SpecFunc; with_func_ptr = false)
    p = init_wrapper_func(spec)

    count_ptr_index = findfirst(x -> (is_length(x) || is_size(x)) && x.requirement == POINTER_REQUIRED, children(spec))
    queried_params = getindex(children(spec), findall(is_implicit_return, children(spec)))
    if !isnothing(count_ptr_index)
        count_ptr = children(spec)[count_ptr_index]
        queried_params =
            getindex(children(spec), findall(x -> x.len == count_ptr.name && !x.is_constant, children(spec)))

        first_call_args = map(@λ(begin
            &count_ptr => count_ptr.name
            GuardBy(in(queried_params)) => :C_NULL
            x => vk_call(x)
        end), children(spec))

        i = 0
        second_call_args = map(@λ(begin
            :C_NULL && Do(i += 1) => queried_params[i].name
            x => x
        end), first_call_args)

        p[:body] = concat_exs(
            initialize_ptr(count_ptr),
            wrap_enumeration_api_call(
                spec,
                wrap_api_call(spec, first_call_args; with_func_ptr),
                (is_length(count_ptr) ? initialize_array : initialize_ptr).(queried_params, count_ptr)...,
                wrap_api_call(spec, second_call_args; with_func_ptr),
                ;
                free = filter(is_data, queried_params),
            )...,
            wrap_implicit_return(spec, queried_params; with_func_ptr),
        )

        args = filter(!in(vcat(queried_params, count_ptr)), children(spec))

        ret_type = @match length(queried_params) begin
            if any(is_data_with_retrievable_size, queried_params)
            end => Expr(:curly, :Tuple, nice_julian_type.([unique(len.(queried_params)); queried_params])...)
            1 => nice_julian_type(first(queried_params))
            _ => Expr(
                :curly,
                :Tuple,
                (nice_julian_type(param) for param ∈ queried_params)...,
            )
        end
    elseif !isempty(queried_params)
        call_args = map(@λ(begin
            x && GuardBy(in(queried_params)) => x.name
            x => vk_call(x)
        end), children(spec))

        p[:body] = concat_exs(
            map(initialize_ptr, queried_params)...,
            wrap_api_call(spec, call_args; with_func_ptr),
            wrap_implicit_return(spec, queried_params; with_func_ptr),
        )

        args = filter(!in(filter(x -> x.requirement ≠ POINTER_REQUIRED, queried_params)), children(spec))

        ret_type = @match length(queried_params) begin
            1 => nice_julian_type(first(queried_params))
            _ => Expr(:curly, :Tuple, (nice_julian_type(param) for param ∈ queried_params)...)
        end
    else
        p[:short] = true
        p[:body] = :($(wrap_api_call(spec, map(vk_call, children(spec)); with_func_ptr)))

        args = children(spec)

        ret_type = nice_julian_type(spec.return_type)
    end

    add_func_args!(p, spec, args; with_func_ptr)
    p[:return_type] = wrap_return_type(spec, ret_type)
    p
end

"""
Extend functions that create (or allocate) one or several handles,
by exposing the parameters of the associated CreateInfo structures.
`spec` must have one or several CreateInfo arguments.
"""
function extend_handle_constructor(spec::CreateFunc; with_func_ptr = false)
    p_func = wrap(spec.func)
    @assert !isnothing(spec.create_info_param) "Cannot extend handle constructor with no create info parameter."
    p_info = add_constructor(spec.create_info_struct)

    args = [p_func[:args]; p_info[:args]]
    kwargs = [p_func[:kwargs]; p_info[:kwargs]]

    info_expr = reconstruct_call(p_info; is_decl = false)
    info_index = findfirst(==(p_info[:name]), innermost_type.(type.(p_func[:args])))
    deleteat!(args, info_index)

    func_call_args = Vector{Any}(name.(p_func[:args]))
    func_call_args[info_index] = info_expr

    if with_func_ptr
        append!(args, func_ptr_args(spec.func))
        append!(func_call_args, name.(func_ptrs(spec.func)))
    end

    body = reconstruct_call(Dict(:name => p_func[:name], :args => func_call_args, :kwargs => name.(p_func[:kwargs])))

    Dict(
        :category => :function,
        :name => p_func[:name],
        :args => args,
        :kwargs => kwargs,
        :short => true,
        :body => body,
    )
end


text:
module FastAIDatasets

using DataDeps
using FilePathsBase
using FilePathsBase: filename
import FileIO
using FileTrees
using MLDataPattern
using MLDataPattern: splitobs
import LearnBase
using Colors
using FixedPointNumbers

include("fastaidatasets.jl")

function __init__()
    initdatadeps()
end

include("containers.jl")
include("transformations.jl")
include("load.jl")

export
    # reexports from MLDataPattern
    splitobs,
    getobs,
    nobs,

    # container transformations
    mapobs,
    filterobs,
    groupobs,
    joinobs,

    # primitive containers
    FileDataset,
    TableDataset,

    # utilities
    isimagefile,
    loadfile,
    filename,

    # datasets
    DATASETS,
    loadtaskdata,
    datasetpath

end  # module


text:
using Printf

test = ["I", "III", "IX", "IVI", "IIM",
        "CMMDXL", "icv", "cDxLiV", "MCMLD", "ccccccd",
        "iiiiiv", "MMXV", "MCMLXXXIV", "ivxmm", "SPQR"]
for rnum in test
    @printf("%15s → %s\n", rnum, try parseroman(rnum) catch "not valid" end)
end

max_stars_repo_name: ethansaxenian/RosettaDecode

max_issues_repo_name: ethansaxenian/RosettaDecode

max_forks_repo_name: ethansaxenian/RosettaDecode


text:
module TestIoFuncs

using Base.Test
using DataArrays
using DataFrames
using Dates

using TimeData

println("Running IO function tests")

dirPath = Pkg.dir("TimeData")
filePath = joinpath(dirPath, "data", "logRet.csv")
td = readTimedata(filePath)

end


text:
using MeasureBase
using Distributions

@testset "distributions" begin
    d = Rice(ν=2.0, σ=0.5)
    logdensity(d, 1.0)

    samples = [rand(d) for _ in 1:100_000]
    ms = mean(samples)
    vs = var(samples)
    @test isapprox(ms, mean(d); atol=1e-2)
    @test isapprox(vs, var(d); atol=1e-2)

    dhsnr = Rice(ν=5.0, σ=1e-3)
    μ = sqrt(dhsnr.ν^2 - dhsnr.σ^2)
    @test isapprox(logdensity(dhsnr, 4.9), logpdf(Normal(μ, dhsnr.σ), 4.9); atol=1e-1)

    dvm1 = CPVonMises(μ=0.0, σ=0.1)
    dvm2 = CPVonMises(μ=0.0, κ=1/0.1^2)

    # Take off the normalization piece
    ds = VonMises(0.0, 1/0.1^2)
    norm = -log(ds.I0κx) - log(2π)
    @test logpdf(ds, 0.5) - norm ≈ logdensity(dvm2, 0.5)

    @test logdensity(dvm1, 0.5) ≈ logdensity(dvm2, 0.5)
end


text:
"""
    half_saturation(max, half, x)

Half satration curve.
"""
half_saturation(max, half, x) = max/(oneunit(half/x) + half/x)

"""
    reserve_drain!(o::AbstractOrgan, column, drain)

Generalised reserve drain for any flux column (ie :gro),
specified by its `Symbol` or `Int` index.
"""
@inline reserve_drain!(o::AbstractOrgan, col, drain) = reserve_drain!(has_reserves(o), o, col, drain)
@inline reserve_drain!(::HasCNE, o, col, drain) = begin
    ΘE, J = θE(o), flux(o)
    J_CN = -drain * (1 - θ) # fraction on drain from C and N reserves
    @inbounds J[:C,col] = J_CN / y_E_C(o)
    @inbounds J[:N,col] = J_CN / y_E_N(o)
    @inbounds J[:E,col] = -drain * ΘE
    nothing
end
@inline reserve_drain!(::HasCN, o, col, drain) = begin
    J = flux(o)
    @inbounds J[:C,col] = -drain / y_E_C(o)
    @inbounds J[:N,col] = -drain / y_E_N(o)
    nothing
end

"""
    sum_flux!(du, organs::Tuple)

Sum flux matrix and write to `du`.
"""
sum_flux!(du, organs::Tuple) = begin
    offset_apply(sum_flux!, du, organs, 0)
    du
end
sum_flux!(du, organ::Organ, offset::Int) = begin
    J = flux(organ)
    z = zero(J[1, 1])
    for i in 1:size(J, 1)
        s = z
        for j in 1:size(J, 2)
            s += J[i, j]
        end
        du[i + offset] = s
    end
    offset + size(J, 1)
end


text:
function write_cstring(io, str::String)
    write(io, string(str), 0x00)
end

function write_elem(io::IO, elem::Element)
    ret = 0
    ret += write(io, bson_code(elem.value), elem.key, 0x00)
    ret += write_value(io, elem.value)
    ret
end

function write_elist(io::IO, elist::Vector{Element})
    ret = 0
    head = position(io)
    ret += write(io, Int32(0))
    [write_elem(io, x) for x in elist]
    ret += write(io, 0x00)
    tail = position(io)
    seek(io, head)
    ret += write(io, Int32(tail - head))
    seek(io, tail)
    ret
end

write_value(io::IO, x::Missing) = 0
write_value(io::IO, x::Nothing) = 0
write_value(io::IO, x::Maxkey) = 0
write_value(io::IO, x::Minkey) = 0
write_value(io::IO, x::Bool) = write(io, x)
write_value(io::IO, x::Int32) = write(io, x)
write_value(io::IO, x::Int64) = write(io, x)
write_value(io::IO, x::Float64) = write(io, x)
write_value(io::IO, x::DateTime) = write(io, value(x) - UNIXEPOCH)
write_value(io::IO, x::BSONSymbol) = write_value(io, x.value)
write_value(io::IO, x::Timestamp) = write(io, x.value)
write_value(io::IO, x::ObjectId) = write(io, x.value)
write_value(io::IO, x::Decimal128) = write(io, x.value)
write_value(io::IO, x::DBPointer) = write_value(io, x.ref) + write_value(io, x.id)
write_value(io::IO, x::Code) = write_value(io, x.code)
write_value(io::IO, x::BSONRegex) = write_cstring(io, x.pattern) + write_cstring(io, x.options)
write_value(io::IO, x::Document) = write_elist(io, x.elist)
write_value(io::IO, x::BSONArray) = write_elist(io, x.elist)

function write_value(io::IO, x::String)
    @assert isvalid(x)
    write(io, Int32(sizeof(x) + 1), x, 0x00)
end

function write_value(io::IO, x::Binary)
    len = Int32(length(x.bytes))
    ret = 0
    if x.subtype == 0x02
        ret += write(io, Int32(len + sizeof(len)), x.subtype, len)
    else
        ret += write(io, len, x.subtype)
    end
    ret += write(io, x.bytes)
    ret
end

function write_value(io::IO, x::CodeWithScope)
    ret = 0
    head = position(io)
    ret += write(io, Int32(0))
    ret += write_value(io, x.code)
    ret += write_value(io, x.scope)
    tail = position(io)
    seek(io, head)
    ret += write(io, Int32(tail - head))
    seek(io, tail)
    ret
end

# high level functions

function write_bson(io::IO, doc; codec::Union{AbstractBSONCodec,Symbol} = :bson)
    codec = codec isa Symbol ? _codec[codec] : codec
    write_value(io, encode(codec, doc))
end

function write_bson(doc; kwargs...)
    io = IOBuffer(read = false, write = true)
    write_bson(io, doc; kwargs...)
    take!(io)
end

write_bson(filename::AbstractString, doc; kwargs...) = open(x -> write_bson(x, doc; kwargs...), filename, "w")


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

using Markdown
using InteractiveUtils

# ╔═╡ 249f32b0-573c-11ec-234c-57a7c9b9fc62
3+4

# ╔═╡ Cell order:
# ╠═249f32b0-573c-11ec-234c-57a7c9b9fc62


text:
Base.promote_rule(::Type{DataValueArray{T,N}}, ::Type{Array{T,N}}) where {T,N} = DataValueArray{T,N}
Base.promote_rule(::Type{Array{T,N}}, ::Type{DataValueArray{T,N}}) where {T,N} = DataValueArray{T,N}
Base.promote_rule(::Type{DataValueArray{T,N}}, ::Type{Array{S,N}}) where {T,S,N} = DataValueArray{Base.promote_type(T,S),N}
Base.promote_rule(::Type{Array{S,N}}, ::Type{DataValueArray{T,N}}) where {T,S,N} = DataValueArray{Base.promote_type(T,S),N}


text:
using Documenter, AtomDescribe

makedocs(;
    modules=[AtomDescribe],
    format=Documenter.HTML(),
    pages=[
        "Home" => "index.md",
    ],
    repo="https://github.com/yaozhenghangma/AtomDescribe.jl/blob/{commit}{path}#L{line}",
    sitename="AtomDescribe.jl",
    authors="Yaozhenghang Ma",
    assets=String[],
)

deploydocs(;
    repo="github.com/yaozhenghangma/AtomDescribe.jl",
)


text: precompile(_package_hook, (Base.PkgId,))


text:
##### MEAN-VARIANCE OPTIMIZATION
#This is not supposed to produce great results, since mean and covariance
#is supposed to be very exact before good weights are found
using JuMP
using Gurobi
using StatsBase
using DataFrames
using CSV
include("SupportFunction.jl")
include("DataLoad.jl")
println("Leeeeroooy Jenkins")
#Esben's path
path = "$(homedir())/Documents/GitHub/Thesis/Data/IndexDataDiff"

#path = "/Users/SkipperAfRosenborg/Google Drive/DTU/10. Semester/Thesis/GitHubCode/Thesis/Data/IndexData/"
#path = "/zhome/9f/d/88706/SpecialeCode/Thesis/Data/IndexData/"

trainingSize = 240
possibilities = 5
industries = ["NoDur", "Durbl", "Manuf", "Enrgy", "HiTec", "Telcm", "Shops", "Hlth", "Utils", "Other"]
industriesTotal = length(industries)

modelMatrix = zeros(industriesTotal, possibilities)
noDurModel = [1 0 1 1 1]
testModel = [0 1 0 0 0]
modelMatrix[1, :] = noDurModel
for i=2:industriesTotal
    modelMatrix[i, :] = noDurModel
end

XArrays = Array{Array{Float64, 2}}(industriesTotal)
YArrays = Array{Array{Float64, 2}}(industriesTotal)

XArrays, YArrays = generateXandYs(industries, modelMatrix)

nGammas = 5
standY = YArrays[1]
SSTO = sum((standY[i]-mean(standY[:]))^2 for i=1:length(standY))
lambdaValues = log.(logspace(100, SSTO/2, nGammas))
nRows = size(standY)[1]
amountOfModels = nGammas^4

modelConfig = zeros(amountOfModels, 4)
counter = 1
for l1 = 1:nGammas
    for l2 = 1:nGammas
        for l3 = 1:nGammas
            for l4 = 1:nGammas
                modelConfig[counter,:] = [lambdaValues[l1] lambdaValues[l2] lambdaValues[l3] lambdaValues[l4]]
                counter += 1
            end
        end
    end
end
modelConfig

#Initialization of parameters
w1N = repeat([0.1], outer = 11) #1/N weights
gamma = 10 #risk aversion
validationPeriod = 5
PMatrix = zeros(nRows-trainingSize, amountOfModels)
return1NMatrix = zeros(nRows-trainingSize, amountOfModels)
returnCEOMatrix = zeros(nRows-trainingSize, amountOfModels)
returnPerfectMatrix = zeros(nRows-trainingSize)

bestModelAmount = 5
bestModelConfigs = zeros(bestModelAmount, 4)
bestModelIndexes = zeros(bestModelAmount)

weightsPerfect = zeros(nRows-trainingSize, 10)
weightsCEO     = zeros(nRows-trainingSize, 10)

#Establishing perfect results in order to avoid doing same mean-variance calculation over and over
for t=1:(nRows-trainingSize-2)
    trainingXArrays, trainingYArrays, validationXRows, validationY, OOSXArrays, OOSYArrays, OOSRow, OOSY = createDataSplits(XArrays, YArrays, t, trainingSize)
    valY = zeros(10)
    for i = 1:10
        valY[i] = validationY[i][1]
    end
    weightsPerfect[t, :], returnPerfectMatrix[t] = findPerfectResults(trainingXArrays, valY, valY, gamma)
end

println("Starting CEO Validation loop")
for t=1:50# (nRows-trainingSize-2)
    println("Time $t/50")
    trainingXArrays, trainingYArrays, validationXRows, validationY, OOSXArrays, OOSYArrays, OOSRow, OOSY = createDataSplits(XArrays, YArrays, t, trainingSize)

    if t <= validationPeriod
        for m = 1:amountOfModels
            betaArray, U = @time(runCEO(trainingXArrays, trainingYArrays, modelConfig[m, :], gamma))
            expectedReturns = generateExpectedReturns(betaArray, trainingXArrays, trainingYArrays, validationXRows)

            #Need to send OOSRow to mean-variance optimization to get "perfect information" since validationY is the values in OOSRow[1:10]
            valY = zeros(10)
            for i = 1:10
                valY[i] = validationY[i][1]
            end
            #return1N, returnCEO, wStar = performMVOptimization(expectedReturns, U, gamma, OOSRow[1][1:10], valY)
            return1N, returnCEO, wStar = performMVOptimization(expectedReturns, U, gamma, valY, valY)
            weightsCEO[t, 1:10]     = wStar
            return1NMatrix[t, m]      = return1N
            returnCEOMatrix[t, m]     = returnCEO
            returnPerfect = returnPerfectMatrix[t]
            println("1N returns is $return1N, returnPerfect is $returnPerfect and returnCEO is $returnCEO")
            PMatrix[t, m] = calculatePvalue(return1N, returnPerfect, returnCEO)
            #trackReturn(returnCEOTotal, returnCEO)
        end
    elseif t == (validationPeriod+1)
        modelMeans = mean(PMatrix, 1)
        for i = 1:bestModelAmount
            bestModelIndex = indmax(modelMeans)
            bestModelConfig[i] = modelConfig[bestModelIndex, :]
            bestModelIndexes[i] = bestModelIndex
            #place 0 at index place now and take new max as "next best model"
            modelMeans[bestModelIndex] = 0
        end

        for m = 1:bestModelAmount
            betaArray, U = @time(runCEO(trainingXArrays, trainingYArrays, bestModelConfigs[m, :]))
            expectedReturns = generateExpectedReturns(betaArray, trainingXArrays, trainingYArrays, validationXRows)

            #Need to send OOSRow to mean-variance optimization to get "perfect information" since validationY is the values in OOSRow[1:10]
            valY = zeros(10)
            for i = 1:10
                valY[i] = validationY[i][1]
            end
            return1N, returnPerfect, returnCEO, wStar = performMVOptimization(expectedReturns, U, gamma, valY, valY)
            weightsCEO[t, 1:10]     = wStar
            return1NMatrix[t, bestModelIndexes[m]]      = return1N
            returnCEOMatrix[t, bestModelIndexes[m]]     = returnCEO
            returnPerfect = returnPerfectMatrix[t]
            println("1N returns is $return1N, returnPerfect is $returnPerfect and returnCEO is $returnCEO")
            PMatrix[t, bestModelIndexes[m]] = calculatePvalue(return1N, returnPerfect, returnCEO)
        end
    else
        for m = 1:bestModelAmount
            betaArray, U = @time(runCEO(trainingXArrays, trainingYArrays, bestModelConfigs[m, :]))
            expectedReturns = generateExpectedReturns(betaArray, trainingXArrays, trainingYArrays, validationXRows)

            #Need to send OOSRow to mean-variance optimization to get "perfect information" since validationY is the values in OOSRow[1:10]
            valY = zeros(10)
            for i = 1:10
                valY[i] = validationY[i][1]
            end
            return1N, returnCEO, wStar = performMVOptimization(expectedReturns, U, gamma, valY, valY)
            weightsCEO[t, 1:10]     = wStar
            return1NMatrix[t, bestModelIndexes[m]]      = return1N
            returnCEOMatrix[t, bestModelIndexes[m]]     = returnCEO
            returnPerfect = returnPerfectMatrix[t]
            println("1N returns is $return1N, returnPerfect is $returnPerfect and returnCEO is $returnCEO")
            PMatrix[t, bestModelIndexes[m]] = calculatePvalue(return1N, returnPerfect, returnCEO)
        end

    end

end
Array(return1NMatrix[1:200, 1])
returnPerfectMatrix[1:200, 1]
returnCEOMatrix[1:200, 1]
combinedPortfolios = hcat(returnPerfectMatrix[1:200, 1], return1NMatrix[1:200, 1], returnCEOMatrix[1:200, 1], PMatrix[1:200, 1])
writedlm("returnPvalueOutcome1to200.csv", combinedPortfolios, ",")

function createDataSplits(XArrays, YArrays, t, trainingSize)
    trainingXArrays = Array{Array{Float64, 2}}(industriesTotal)
    trainingYArrays = Array{Array{Float64, 2}}(industriesTotal)

    validationXRows = Array{Array{Float64, 2}}(industriesTotal)
    validationY     = Array{Array{Float64, 2}}(industriesTotal)

    OOSXArrays      = Array{Array{Float64, 2}}(industriesTotal)
    OOSYArrays      = Array{Array{Float64, 2}}(industriesTotal)

    OOSRow          = Array{Array{Float64, 2}}(industriesTotal)
    OOSY            = Array{Array{Float64, 2}}(industriesTotal)
    for i = 1:industriesTotal
        trainingXArrays[i] = XArrays[i][t:(t-1+trainingSize), :]
        trainingYArrays[i] = YArrays[i][t:(t-1+trainingSize), :]

        validationXRows[i] = XArrays[i][(t+trainingSize):(t+trainingSize), :]
        validationY[i]     = YArrays[i][(t+trainingSize):(t+trainingSize), :]

        OOSXArrays[i]      = XArrays[i][(t+1):(t+trainingSize), :]
        OOSYArrays[i]      = YArrays[i][(t+1):(t+trainingSize), :]

        OOSRow[i]          = XArrays[i][(t+trainingSize+1):(t+trainingSize+1), :]
        OOSY[i]            = YArrays[i][(t+trainingSize+1):(t+trainingSize+1), :]
    end
    return trainingXArrays, trainingYArrays, validationXRows, validationY, OOSXArrays, OOSYArrays, OOSRow, OOSY
end

function calculatePvalue(return1N, returnPerfect, returnCEO)
    PValue = 1 - (returnCEO-returnPerfect)/(return1N-returnPerfect)
    return PValue
end

function findPerfectResults(trainingXArrays, Xrow, Yvalues, gamma)
    indexes = 10
    Sigma =  cov(trainingXArrays[1][:,1:10])

    #A=U^(T)U where U is upper triangular with real positive diagonal entries
    F = lufact(Sigma)
    U = F[:U]  #Cholesky factorization of Sigma

    M = JuMP.Model(solver = GurobiSolver(OutputFlag = 0))
    @variables M begin
            w[1:indexes]
            u[1:indexes]
            z
            y
    end

    @objective(M,Min, gamma*y - Yvalues'*w)
    @constraint(M, 0 .<= w)
    @constraint(M, sum(w[i] for i=1:indexes) == 1)
    @constraint(M, norm([2*U'*w;y-1]) <= y+1)
    solve(M)
    wPerfect = getvalue(w)
    forecastRow = (exp10(Xrow')-1)*100
    periodPerfectReturn = forecastRow*wPerfect

    return wPerfect, periodPerfectReturn
end

function performMVOptimization(expectedReturns, U, gamma, Xrow, Yvalues)
    indexes = 10

    M = JuMP.Model(solver = GurobiSolver(OutputFlag = 0))
    @variables M begin
            w[1:indexes]
            u[1:indexes]
            z
            y
    end

    @objective(M,Min, gamma*y - expectedReturns'*w)
    @constraint(M, 0 .<= w)
    @constraint(M, sum(w[i] for i=1:indexes) == 1)
    @constraint(M, norm([2*U'*w;y-1]) <= y+1)
    solve(M)
    wStar = getvalue(w)

    forecastRow = (exp(Xrow)-1)*100

    periodReturn = forecastRow'*wStar
    period1NReturn = forecastRow'*w1N

    return period1NReturn, periodReturn, wStar
end

function generateExpectedReturns(betaArray, trainingXArrays, trainingYArrays, validationXRows)
    industriesTotal = 10

    expectedReturns = zeros(10)
    bRows = 240
    bCols = floor(Int,zeros(industriesTotal))
    for i = 1:industriesTotal
        bCols[i] = size(trainingXArrays[i])[2]
    end

    for i = 1:10
        errorArray = (trainingXArrays[i]*betaArray[i,1:bCols[i]]-trainingYArrays[i])*(1/bRows)
        errorSum = sum(errorArray)
        expectedReturns[i] = validationXRows[i][:]'*betaArray[i,1:bCols[i]]
    end
    return expectedReturns
end

function runCEO(trainingXArrays, trainingYArrays, modelConfigRow, gamma)
    industriesTotal = 10
    periodMean = zeros(10)
    for i = 1:10
        periodMean[i] = mean(trainingXArrays[1][:,i])
    end
    l1, l2, l3, l4 = modelConfigRow
    env = Gurobi.Env()

    maxPredictors = 1412
    bCols = floor(Int,zeros(industriesTotal))
    Sigma =  cov(trainingXArrays[1][:,1:10])

    #A=U^(T)U where U is upper triangular with real positive diagonal entries
    F = lufact(Sigma)
    U = F[:U]  #Cholesky factorization of Sigma

    for i = 1:industriesTotal
        bCols[i] = size(trainingXArrays[i])[2]
    end
    gamma = gamma
    bRows = size(trainingXArrays[1])[1]
    M = JuMP.Model(solver = GurobiSolver(env, OutputFlag = 0, Threads=(nprocs()-1)))

    @variables M begin
            b[1:industriesTotal, 1:maxPredictors] #Beta in LASSO; industry i and predictors 1:1304
            z[1:industriesTotal, 1:maxPredictors] #auxilliary for 1norm
            q[1:industriesTotal] #auxilliary for norm
            f[1:industriesTotal, 1:bRows] #expected return for industry i at time t
            w[1:industriesTotal, 1:bRows] #weights for industry i at time t
            y[1:bRows] #auxilliary variables for MV optimization
    end

    #@objective(M,Min,sum(0.5*q[i] + l2*ones(maxPredictors)'*z[i,:] for i=1:industriesTotal)+ l3*(sum(gamma*y[t]-f[:,t]'*w[:,t] for t=1:bRows)) - l4*(sum(w[i,t]*trainingYArrays[i][t] for i=1:10, t=1:bRows)))

    @objective(M,Min,sum(l1*q[i] + l2*ones(maxPredictors)'*z[i,:] for i=1:industriesTotal)+ l3*(sum(gamma*y[t]-sum(f[i,t] for i=1:10) for t=1:bRows)) - l4*(sum(w[i,t]*trainingYArrays[i][t] for i=1:10, t=1:bRows)))

    for i = 1:industriesTotal
        for t = 1:bRows
            @constraint(M, (periodMean[i]+1000)*w[i,t] >= f[i,t])
        end
    end
    for i = 1:industriesTotal
        @constraint(M, norm( [1-q[i];2*(trainingXArrays[i]*b[i,1:bCols[i]]-trainingYArrays[i])] ) <= 1+q[i]) #second order cone constraint SOC
    end
    @constraint(M,  b .<= z)
    @constraint(M, -z .<= b)

    for i = 1:10
        errorArray = (trainingXArrays[i]*b[i,1:bCols[i]]-trainingYArrays[i])*(1/bRows)
        errorSum = sum(errorArray)
        for t = 1:bRows
            prediction = trainingXArrays[i][t, :]'*b[i,1:bCols[i]]
            @constraint(M, prediction + errorSum >= f[i, t]) ##== originally
        end
    end

    @constraint(M, 0 .<= w)
    for t = 1:bRows
        @constraint(M, sum(w[i,t] for i=1:industriesTotal) == 1)
        @constraint(M, norm([2*U'*w[:,t];y[t]-1]) <= y[t]+1)
    end

    @time(solve(M))
    betaArray = getvalue(b)
    #Insert shrinking here
    return betaArray, U
end

function generateXandYs(industries, modelMatrix)
    industriesTotal = length(industries)
    for i = 1:industriesTotal
        industry = industries[i]
        mainData = loadIndexDataLOGReturn(industry, path)

        #path = "/zhome/9f/d/88706/SpecialeCode/"
        #path = "/Users/SkipperAfRosenborg/Google Drive/DTU/10. Semester/Thesis/GitHubCode/"
        dateAndRecession = Array(mainData[:,end-1:end])
        mainDataArr = Array(mainData[:,1:end-2])

        colNames = names(mainData)

        nRows = size(mainDataArr)[1]
        nCols = size(mainDataArr)[2]

        #fileName = path*"/Results/IndexData/LassoTests/240-1/2401_Shrink_"
        #fileName = path*"/Results/Test/Parallel/"
        #fileName = path*"Results/IndexData/LassoTest/"*industry*"/"*string(trainingSize)*"-1/"*string(trainingSize)*"_"

        #Reset HPC path
        #path = "/zhome/9f/d/88706/SpecialeCode/Thesis/ML"
        #cd(path)
        fileName = path*"Results/IndexData/CEOTest/"*industry*"/"*string(trainingSize)*"-1/"*string(trainingSize)*"_"

        mainYarr = mainDataArr[:, nCols:nCols]

        mainXarr = mainDataArr[:,1:nCols-1]

        VIX = modelMatrix[i, 1]
        raw = modelMatrix[i, 2]
        expTrans = modelMatrix[i, 3]
        timeTrans = modelMatrix[i, 4]
        TA = modelMatrix[i,5]

        if raw == 1
            fileName = fileName*"Raw"
            mainXarr = mainXarr[:,1:10]
        else
            if VIX == 1
                fileName = fileName*"VIX_Macro"
                mainXarr = mainXarr[:,1:end-1]
            else
                fileName = fileName*"Macro"
            end
        end

        # Transform with time elements
        if timeTrans == 1
            fileName = fileName*"Time"
            mainXarr = expandWithTime3612(mainXarr)
        end

        # Transform #
        if expTrans == 1
            fileName = fileName*"Exp"
            mainXarr = expandWithTransformations(mainXarr)
        end

        if TA == 1
            fileName = fileName*"TA"
            mainXarr = expandWithMAandMomentum(mainXarr, mainYarr, (nCols-1))
        end

        # Standardize
        standX = zScoreByColumn(mainXarr)
        standY = mainYarr

        XArrays[i] = standX
        YArrays[i] = standY
    end
    return XArrays, YArrays
end

"""
Keeps track of returns for both the proposed method and 1/N.
Needs the overall array as well as the two next returns to push
"""

function trackPortfolioReturns(ReturnOverall, individualReturn, newPeriodReturn, Return1NOverall, Individual1NReturn, period1NReturn)
    currentLength = length(ReturnOverall)
    ReturnOverall = push!(ReturnOverall, ReturnOverall[currentLength]*newPeriodReturn)
    individualReturn = push!(individualReturn, newPeriodReturn)

    current1NLength = length(Return1NOverall)
    Return1NOverall = push!(Return1NOverall, Return1NOverall[current1NLength]*period1NReturn)
    Individual1NReturn = push!(Individual1NReturn, period1NReturn)
end

function trackPortfolioWeights(weights, periodWeights, iteration)
    weights[iteration, 1:10] = periodWeights
end

max_stars_repo_name: skipperafrosenborg/Thesis

max_issues_repo_name: skipperafrosenborg/Thesis

max_forks_repo_name: skipperafrosenborg/Thesis


text:
function make_buttons(cur_observer::Observable)
  cur_buttons = dom"div"(
    dom"div"(
      dom"a"(
        dom"i"(
          attributes=Dict(:class => "fab fa-github-alt cs-icon js-github-button")
        ),
        attributes=Dict(
          :href => "https://github.com/djsegal/Tetris.jl",
          :target => "_blank"
        )
      )
    ),
    dom"div"(
      dom"a"(
        dom"i"(
          attributes=Dict(:class => "fas fa-volume-down cs-icon js-music-button")
        ),
        events=Dict("click" => @js function (cur_event)
          if $cur_observer[] == "music"
            $cur_observer[] = "music-music";
          else
            $cur_observer[] = "music";
          end
        end)
      )
    ),
    dom"div"(
      dom"a"(
        dom"i"(
          attributes=Dict(:class => "fas fa-pause cs-icon js-pause-button")
        ),
        events=Dict("click" => @js function (cur_event)
          if $cur_observer[] == "pause"
            $cur_observer[] = "pause-pause";
          else
            $cur_observer[] = "pause";
          end
        end)
      )
    ),
    dom"div"(
      dom"a"(
        dom"i"(
          attributes=Dict(:class => "fas fa-times cs-icon")
        ),
        events=Dict("click" => @js function (cur_event)
          $cur_observer[] = "suspend";

          is_sure = window.confirm("Are you sure you want to restart?");

          if is_sure
            $cur_observer[] = "restart";
          end
        end)
      )
    ),
    attributes=Dict(:class => "cs-icon-list")
  )

  cur_buttons
end


text:
"""
Latex code

For notation: set up a `Dict{String, String}`. No need for code.
"""
module LatexLH

using DocStringExtensions

export write_figure_slide, figure_slide

include("table.jl");
include("parameter_table.jl")
include("beamer.jl")
# include("latex/symbol_table.jl")

end


text:
include("euler/euler.jl")
using .Grid: max_path_sum

function compute(path::String)::Int
    triangle = map(x -> parse.(Int, split(x, " ")), split(read(path, String), "\n"))
    return max_path_sum(triangle)
end

max_stars_repo_name: Dynortice/Project-Euler

max_issues_repo_name: Dynortice/Project-Euler

max_forks_repo_name: Dynortice/Project-Euler


text:
# # Logistic regression
using DataFrames
using Plots
using RDatasets
using Convex
using SCS

# This is an example logistic regression using `RDatasets`'s iris data.
# Our goal is to gredict whether the iris species is versicolor
# using the sepal length and width and petal length and width.
iris = dataset("datasets", "iris");
iris[1:10,:]

# We'll define `Y` as the outcome variable: +1 for versicolor, -1 otherwise.
Y = [species == "versicolor" ? 1.0 : -1.0 for species in iris.Species]

# We'll create our data matrix with one column for each feature
# (first column corresponds to offset).
X = hcat(ones(size(iris, 1)), iris.SepalLength, iris.SepalWidth, iris.PetalLength, iris.PetalWidth);

# Now to solve the logistic regression problem.
n, p = size(X)
beta = Variable(p)
problem = minimize(logisticloss(-Y.*(X*beta)))
solve!(problem, SCS.Optimizer(verbose=false))

# Let's see how well the model fits.
using Plots
logistic(x::Real) = inv(exp(-x) + one(x))
perm = sortperm(vec(X*beta.value))
plot(1:n, (Y[perm] .+ 1)/2, st=:scatter)
plot!(1:n, logistic.(X*beta.value)[perm])


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

import Test

Test.@testset "Package ValueShapes" begin
    include("test_value_shape.jl")
    include("test_value_accessor.jl")
    include("test_scalar_shape.jl")
    include("test_array_shape.jl")
    include("test_const_value_shape.jl")
    include("test_named_tuple_shape.jl")
    include("test_functions.jl")
    include("test_distributions.jl")
    include("test_const_value_dist.jl")
    include("test_named_tuple_dist.jl")
    include("test_reshaped_dist.jl")
end # testset


text:
@defcomp sealeveldamages begin
    regions = Index()

    SLRDAMAGES = Variable(index=[time, regions])

    slrmultiplier = Parameter(index=[regions])
    slrelasticity = Parameter(index=[regions])
    slrdamlinear = Parameter(index=[regions])
    slrdamquadratic = Parameter(index=[regions])
    TOTALSLR = Parameter(index=[time])
    YGROSS = Parameter(index=[time, regions])

    function run_timestep(p, v, d, t)
        #Define function for SLRDAMAGES
        for r in d.regions
            if is_first(t)
                v.SLRDAMAGES[t,r] = 0.
            else
                v.SLRDAMAGES[t,r] = 100. * p.slrmultiplier[r] * (p.TOTALSLR[t-1] * p.slrdamlinear[r] + p.TOTALSLR[t-1]^2 * p.slrdamquadratic[r]) * (p.YGROSS[t-1,r] / p.YGROSS[TimestepIndex(1),r])^(1/p.slrelasticity[r])
            end
        end
    end
end


text:
println("\nRunLargeTest: Demonstrates arrays of components and timing.")

using Modia
using Modia.Electric

@static if VERSION < v"0.7.0-DEV.2005"
  using Base.Test
else
  using Test
end

@testset "RunLargeTests" begin

@model LPfilter begin
  R=Resistor(R=100.0)
  C=Capacitor(C=2.5E-3, v = Float(start=0.0))
  V=ConstantVoltage(V=10.0)
  ground=Ground()
@equations begin
  connect(V.n, ground.p)
  connect(V.p, R.p)
  connect(R.n, C.p)
  connect(C.n, V.n)
  end
end

nFilters=100
@model ManyFilters begin
  F = [LPfilter() for i in 1:nFilters]
end

@time simulate(ManyFilters, 2, logTiming=true, storeEliminated=false, removeSingularities=false)

end # testset
nothing


text:
#=
@doc """
Facilities to help map GUI input / expression into a function that aggregates columns.
We would store that function and use it whenever we need to aggregate a column.

Three ways to instantiate the DataFrameAggr
* Symbol. e.g. `:mean`. It's lowered into `mean(:_)`. See below for the meaning of `:_`
* Expr. e.g. `:( mean(:_,:wcol) )`, `:( quantile( :_, 0.75 ) )`
   symbols (quoted) are interpreted as the named column.
   The underscore `:_` column is special -- it's interpreted on-the-fly as the column that
   needs to be aggregated. So multiple columns can use the same expression `mean(_,:wcol)` and
   they will all use the same weight (from `wcol`), but produces target specific mean.
* Function: (most likely lambda) e.g. `df -> quantile( df[:col], 0.5 )`. It is expected
   to generate either a value (typically a scalar consistent to the column type), or
   a 1-row, 1-col dataframe.
""" ->
=#

DataFrameAggrCache = Dict{Any,Function}()

defaultAggr( ::Type{} ) = :uniqvalue
defaultAggr{T<:Real}( ::Type{T} ) = :sum
defaultAggr{T}( ::Type{Array{T,1}} ) = :unionall

function liftAggrSpecToFunc( c::Symbol, dfa::UTF8String )
    if haskey( DataFrameAggrCache, (c, dfa ) )
        return DataFrameAggrCache[ (c, dfa ) ]
    end
    ret = liftAggrSpecToFunc( c, parse( dfa ) )
    DataFrameAggrCache[ (c, dfa) ] = ret
end
liftAggrSpecToFunc( c::Symbol, dfa::ASCIIString ) = liftAggrSpecToFunc( c, utf8( dfa ) )

function liftAggrSpecToFunc( c::Symbol, dfa::Union{ Function, Symbol, Expr } )
    if typeof( dfa ) == Function
        return dfa
    end
    if haskey( DataFrameAggrCache, (c, dfa) )
        return DataFrameAggrCache[ (c, dfa ) ]
    end
    # "mean" or "Module.aggrfunc"
    if typeof( dfa ) == Symbol || typeof( dfa ) == Expr && dfa.head == :(.) && all( x->typeof(x)==Symbol, dfa.args )
        funnameouter = gensym( "DFAggr" )
        code = :(
            function $funnameouter( _df_::AbstractDataFrame )
            end )
        push!( code.args[2].args, Expr( :call, dfa, Expr( :call, :getindex, :_df_, QuoteNode(c) ) ) )
        ret = eval( code )
    else # expr
        if !Base.Meta.isexpr( dfa, :call )
            error( string( dfa ) * " does not look like an aggregator function")
        end
        # disallow mutating function
        fname = dfa.args[1]
        if Base.Meta.isexpr( dfa, :curly )
            error( "DataFrameAggr: curly not supported")
        elseif !( typeof( fname ) <: Symbol ) && !Base.Meta.isexpr( fname, :(.) )
            error( "DataFrameAggr: only simple function name please")
        end
        if contains( string(fname), "!" )
            error( string(fname) * " seems to have side effects" )
        end

        # replace _ with _df_[$c], and then leverage @with

        # before we do that, note that
        # (A) in DataFramesMeta, macro converts :x to Expr( :quote, :x ) but
        #     :( :x ) or parse( ":x" ) it is actually QuoteNode( :x ),
        #     so we need to do a little conversion in order to leverage that package
        cdfa = deepcopy( dfa )
        convertExpression!( cdfa, c )

        membernames = Dict{Union{Symbol,Expr}, Symbol}()
        cdfa = DataFramesMeta.replace_syms(cdfa, membernames)
        funargs = map(x -> :( getindex( _df_, $(x)) ), collect(keys(membernames)))
        funnameouter = gensym("DFAggr")
        funname = gensym()
        code = quote
            function $funnameouter( _df_::AbstractDataFrame )
                function $funname($(collect(values(membernames))...))
                    $cdfa
                end
                $funname($(funargs...))
            end
        end
        ret = eval( code )
    end
    DataFrameAggrCache[ (c, dfa) ] = ret
end

# used by aggregation lifting and calcpivot lifting
function convertExpression!( ex::Expr, column_ctx::Symbol = Symbol("") )
    for i in 1:length( ex.args )
        a = ex.args[i]
        if typeof( a ) == QuoteNode
            if a.value == :_ && column_ctx != Symbol("")
                ex.args[i] = Expr( :quote, column_ctx )
            else
                ex.args[i] = Expr( :quote, a.value )
            end
        elseif typeof( a ) == Expr
            convertExpression!( a )
        end
    end
end

function unionall( x::AbstractDataArray )
    l = dropna( x )
    t = eltype( eltype( x ) )
    s = Set{t}()
    for el in l
        push!( s, el... )
    end
    collect( s )
end

function unionall( x::Array )
    t = eltype( eltype( x ) )
    s = Set{t}()
    for el in x
        push!( s, el... )
    end
    collect( s )
end

function uniqvalue( x::AbstractDataArray; skipna::Bool=true )
    lvls = DataArrays.levels(x)
    if skipna
        l = dropna( lvls )
        if length(l) == 1
            return l[1]
        end
        return NA
    end
    if length(lvls) == 1
        return lvls[1]
    end
    return NA
end

function uniqvalue{T<:AbstractString}( x::Union{ Array{T}, DataArray{T}, PooledDataArray{T} }; skipna::Bool=true, skipempty::Bool=true )
    lvls = DataArrays.levels(x)
    if skipna
        l = dropna( lvls )
        if skipempty
            emptyidx = findfirst( l, "" )
            if length( l ) == 1 && emptyidx == 0
                return l[1]
            elseif length( l ) == 2 && emptyidx != 0
                if emptyidx == 1
                    return l[2]
                else
                    return l[1]
                end
            end
        elseif length( l ) == 1
            return l[1]
        end
        return NA
    end
    if skipempty
        emptyidx = findfirst( lvls, "" )
        if length( lvls ) == 1 && emptyidx == 0
            return lvls[1] # could be NA
        elseif length( lvls ) == 2 && emptyidx != 0
            if emptyidx == 1
                return lvls[2]
            else
                return lvls[1]
            end
        end
    elseif length( lvls ) == 1
        return lvls[1]
    end
    return NA
end

immutable CalcPivot
    spec::Expr
    by::Array{Symbol,1}
    CalcPivot( x::UTF8String, by::Array{Symbol,1}=Symbol[] ) = CalcPivot( parse(x), by )
    CalcPivot( x::ASCIIString, by::Array{Symbol,1}=Symbol[] ) = CalcPivot( parse(utf8(x)), by )
    CalcPivot( x::UTF8String, by::Symbol ) = CalcPivot( parse(x), Symbol[ by ] )
    CalcPivot( x::ASCIIString, by::Symbol ) = CalcPivot( parse(utf8(x)), Symbol[ by ] )
    function CalcPivot( x::Expr, by::Symbol )
        CalcPivot( x, Symbol[ by ] )
    end
    function CalcPivot( ex::Expr, by::Array{Symbol,1}=Symbol[] )
        if !Base.Meta.isexpr( ex, :call )
            error( string( ex ) * " does not look like an aggregator function")
        end
        # disallow mutating function
        fname = ex.args[1]
        if Base.Meta.isexpr( ex, :curly )
            error( "CalcPivot: curly not supported")
        elseif !( typeof( fname ) <: Symbol ) && !Base.Meta.isexpr( fname, :(.) )
            error( "CalcPivot: only simple function name please")
        end
        if contains( string(fname), "!" )
            error( string(fname) * " seems to have side effects" )
        end

        if fname == :topnames # ensure we have the name in "by"
            # expect the first argument is a QuoteNode, or Expr( :quote, symbol )
            if typeof( ex.args[2] ) == QuoteNode
                name_col = ex.args[2].value
            elseif Base.Meta.isexpr( ex.args[2], :quote )
                name_col = ex.args[2].args[1]
            else
                throw( "topnames: 1st argument expects a symbol (name column)")
            end
            if !in( name_col, by )
                new( ex, Symbol[ by..., name_col ] )
            else
                new( ex, by )
            end
        else
            new( ex, by )
        end
    end
end

CalcPivotFuncCache = Dict{ Any, Function }()
CalcPivotAggrDepCache = Dict{ Any, Array{Symbol,1} }()

# lifting a CalcPivot is a more complicated beast compared to aggregation.
# It is basically the fancier "apply" component of the # split-apply-combine strategy,
# with the "apply" also doing its nested split-apply-combine steps if necessary.

function liftCalcPivotToFunc( ex::Expr, by::Array{Symbol,1} )
    if haskey( CalcPivotFuncCache, (ex, by ) )
        return CalcPivotFuncCache[ (ex, by ) ]
    end

    cex = deepcopy( ex )
    convertExpression!( cex )
    funnameouter = gensym("calcpvt")
    funname = gensym()

    membernames = Dict{Union{Symbol,Expr}, Symbol}()
    cex = DataFramesMeta.replace_syms(cex, membernames)
    # keys are the columns. values are the unique gensyms

    if !isempty( by ) # micro split-apply-combine
        funargs = map(x -> :( getindex( byf, $(Meta.quot(x))) ), collect(keys(membernames)))
        aggregates = setdiff( collect( keys( membernames ) ), by )
        CalcPivotAggrDepCache[ (ex, by ) ] = aggregates
        # basically we want to do
        # byf = by(_df_, by, lambdadf -> DataFrame(b = aggrfuncs[:b](lambdadf), c = aggrfuncs[:c](lambdadf), ... )
        aggr_args = Any[]
        lambdasym = gensym( "lambdadf" )
        for a in aggregates
            push!( aggr_args, Expr( :kw, a, Expr( :call, Expr( :ref, :aggrfuncs, QuoteNode( a ) ), lambdasym ) ) )
        end

        bycolsexpr=Expr( :vcat, map( _->QuoteNode(_), by )... ) # [ :a, :b, :c ... ]
        aggrcode = Expr( :->, lambdasym, Expr( :call, DataFrame, aggr_args... ) )
        code = :(
            function $funnameouter( _df_::AbstractDataFrame, c::Symbol; kwargs... ) # we need kwargs here for aggregate specs
                aggrfuncs = Dict{Symbol,Function}()
                for (aggrc,spec) in kwargs
                    # this would throw if spec is not compliant, as usual
                    aggrfuncs[aggrc] = TermWin.liftAggrSpecToFunc(aggrc, spec )
                end
                # if kwargs doesn't have the required columns, it'll throw in this line
                byf = DataFrames.by( _df_, $bycolsexpr, $aggrcode )
                function $funname($(collect(values(membernames))...))
                    $cex
                end
                vals = $funname($(funargs...))
                ret = byf[ $bycolsexpr ]
                ret[ c ] = vals
                ret # This has all the "by" columns and the result, named by c
                # combine is done by the caller
            end
        )
    else # much simplier, we are doing line-by-line apply.
        # Common use case: a simple bucketing, or a line by line ranking.
        # a "pure" reading on empty "by" array may mean we aggregate everything in the
        # dataframe into one row, and then attach this row's value (a constant) into every row.
        # However, I have not known any legitimate use case for this.
        CalcPivotAggrDepCache[ (ex, by )] = Symbol[]
        funargs = map(x -> :( getindex( _df_, $(x)) ), collect(keys(membernames)))
        code = :(
            function $funnameouter( _df_::AbstractDataFrame ) # we need kwargs here for aggregate specs
                function $funname($(collect(values(membernames))...))
                    $cex
                end
                $funname($(funargs...))
                # the creation of column is done by the caller
            end
        )
    end
    ret = eval( code )
    CalcPivotFuncCache[ (ex, by ) ] = ret
end

function cut_categories{S<:Real, T<:Real}( ::Type{S}, breaks::Vector{T};
    boundedness = :unbounded,
    leftequal=true, # t1 <= x < t2 or t1 < x <= t2?
    absolute=false, # t1 <= |x| < t2?
    rank=true, # add a rank to the string output for easier sorting?
    ranksep = ". ", # "1. t1 <= x < t2"?
    label = "", # if not compact, what label do we use for x?
    compact=(label==""), # <t1, [t1,t2), t2+. Further shortened for integer intervals with length=1
    reverse=false, # reverse the rank from the largest first?
    # the following format the boundary numbers
    # see Formatting.jl
    prefix="", suffix="", scale=1, precision=-1,
    commas=false,stripzeros=(precision==-1),parens=false,
    mixedfraction=false,autoscale=:none,conversion=""
    )
    n = length(breaks)
    breakstrs = UTF8String[]
    function formatter(_)
        prefix * format( _*scale,
            precision=precision,
            commas=commas,
            stripzeros=stripzeros,
            parens=parens,
            mixedfraction=mixedfraction,
            autoscale=autoscale,
            conversion=conversion ) * suffix
    end
    for b in breaks
        push!( breakstrs, formatter( b ) )
    end
    if boundedness == :unbounded
        ncategories = n + 1
    elseif boundedness == :bounded
        ncategories = n-1
    else
        ncategories = n
    end
    pool = Array(UTF8String, ncategories )
    if rank
        rankwidth = length(string(ncategories))
    end
    if !rank
        rankprefixfunc = _->""
    elseif reverse
        rankprefixfunc = j -> format( n+2-j, width=rankwidth ) * ranksep
    else
        rankprefixfunc = j -> format( j, width=rankwidth ) * ranksep
    end
    if compact
        if S <: Integer && T <: Integer && scale == 1
            # we use 1...5, 6, 7...10, 11+etc.
            if leftequal
                breakminus1strs = UTF8String[]
                for b in breaks
                    push!( breakminus1strs, formatter( b-1 ) )
                end
                poolindexshift = -1
                if boundedness in [ :unbounded, :boundedabove ]
                    pool[1] = rankprefixfunc(1) * "≤" * breakminus1strs[1]
                    poolindexshift = 0
                end
                for i in 2:n
                    if breaks[i-1] == breaks[i]-1
                        pool[i+poolindexshift] = rankprefixfunc(i+poolindexshift)*breakstrs[i-1]
                    else
                        pool[i+poolindexshift] = rankprefixfunc(i+poolindexshift)*breakstrs[i-1]*"…"*breakminus1strs[i]
                    end
                end
                if boundedness in [ :unbounded, :boundedbelow ]
                    pool[n+1+poolindexshift] = rankprefixfunc(n+1+poolindexshift)*breakstrs[n]*"+"
                end
            else
                breakplus1strs = UTF8String[]
                for b in breaks
                    push!( breakminus1strs, formatter( b+1 ) )
                end
                poolindexshift = -1
                if boundedness in [ :unbounded, :boundedabove ]
                    pool[1] = rankprefixfunc(1) *"≤ "* breakstrs[1]
                    poolindexshift = 0
                end
                for i in 2:n
                    if breaks[i-1]+1 == breaks[i]
                        pool[i+poolindexshift] = rankprefixfunc(i+poolindexshift)*breakstrs[i]
                    else
                        pool[i+poolindexshift] = rankprefixfunc(i+poolindexshift)*breakplus1strs[i-1]*"…"*breakstrs[i]
                    end
                end
                if boundedness in [ :unbounded, :boundedbelow ]
                    pool[n+1+poolindexshift] = rankprefixfunc(n+1+poolindexshift)*breakplus1strs[n]*"+"
                end
            end
        else # by the way, we don't show absolute in compact
            if leftequal
                brackL = "["
                brackR = ")"
                compareL = utf8( "<" )
                compareR = "≥"
            else
                brackL = "("
                brackR = "]"
                compareL = "≤"
                compareR = utf8(">")
            end
            poolindexshift = -1
            if boundedness in [ :unbounded, :boundedabove ]
                pool[1] = rankprefixfunc(1) * compareL * breakstrs[1]
                poolindexshift = 0
            end
            for i in 2:n
                if i == 2 && boundedness in [ :boundedbelow, :bounded ]
                    pool[i+poolindexshift] = rankprefixfunc(i+poolindexshift) * "[" * breakstrs[i-1]* "," *breakstrs[i] * brackR
                elseif i == n && boundedness in [ :boundedabove, :bounded ]
                    pool[i+poolindexshift] = rankprefixfunc(i+poolindexshift) * brackL * breakstrs[i-1]* "," *breakstrs[i] * "]"
                else
                    pool[i+poolindexshift] = rankprefixfunc(i+poolindexshift) * brackL * breakstrs[i-1]* "," *breakstrs[i] * brackR
                end
            end
            if boundedness in [ :unbounded, :boundedbelow ]
                pool[n+1+poolindexshift] = rankprefixfunc(n+1+poolindexshift) * compareR * breakstrs[n]
            end
        end
    else
        if absolute
            label2 = "|"*label*"|"
        else
            label2 = label
        end
        if leftequal
            compareL = " ≤ "
            compareR = utf8(" < ")
        else
            compareL = utf8(" < ")
            compareR = " ≤ "
        end
        poolindexshift = -1
        if boundedness in [ :unbounded, :boundedabove ]
            pool[1] = rankprefixfunc( 1 ) * label2 * compareR * breakstrs[1]
            poolindexshift = 0
        end
        for i in 2:n
            if i == 2 && boundedness in [ :boundedbelow, :bounded ]
                pool[i+poolindexshift] = rankprefixfunc(i+poolindexshift) * breakstrs[i-1] * " ≤ " * label2 * compareR * breakstrs[i]
            elseif i == n && boundedness in [ :boundedabove, :bounded ]
                pool[i+poolindexshift] = rankprefixfunc(i+poolindexshift) * breakstrs[i-1] * compareL * label2 * " ≤ " * breakstrs[i]
            else
                pool[i+poolindexshift] = rankprefixfunc(i+poolindexshift) * breakstrs[i-1] * compareL * label2 * compareR * breakstrs[i]
            end
        end
        if boundedness in [ :unbounded, :boundedbelow ]
            pool[n+1+poolindexshift] = rankprefixfunc(n+1+poolindexshift) * breakstrs[n] * compareL * label2
        end
    end
    return pool
end

# boundedness:
#    unbounded    gives n+1 categories for n breaks.
#    boundedbelow gives n   categories for n breaks. Values below min will be NA
#    boundedabove gives n   categories for n breaks. Values above max will be NA
#    bounded      gives n-1 categories for n breaks. Values below min or above max will be NA
function discretize{S<:Real, T<:Real}(x::AbstractArray{S,1}, breaks::Vector{T};
    boundedness = :unbounded,
    bucketstrs = UTF8String[], # if provided, all of below will be ignored. length must be length(breaks)+1
    leftequal=true, # t1 <= x < t2 or t1 < x <= t2?
    absolute=false, # t1 <= |x| < t2?
    rank=true, # add a rank to the string output for easier sorting?
    ranksep = ". ", # "1. t1 <= x < t2"?
    label = "", # if not compact, what label do we use for x?
    compact=(label==""), # <t1, [t1,t2), t2+. Further shortened for integer intervals with length=1
    reverse=false, # reverse the rank from the largest first?
    # the following format the boundary numbers
    # see Formatting.jl
    prefix="", suffix="", scale=1, precision=-1,
    commas=false,stripzeros=(precision==-1),parens=false,
    mixedfraction=false,autoscale=:none,conversion=""
    )
    if !issorted(breaks)
        sort!(breaks)
    end
    refs = fill(zero(DataArrays.DEFAULT_POOLED_REF_TYPE), length(x))
    n = length(breaks)
    if absolute
        x2 = abs(x)
    else
        x2 = x
    end

    if boundedness == :unbounded
        below_min_mult = 1
        above_max_mult = 1
        ref_shift = 1
        ncategories = length( breaks ) + 1
    elseif boundedness == :boundedbelow
        below_min_mult = 0
        above_max_mult = 1
        ref_shift = 0
        ncategories = length( breaks )
    elseif boundedness == :boundedabove
        below_min_mult = 1
        above_max_mult = 0
        ref_shift = 1
        ncategories = length( breaks )
    elseif boundedness == :bounded
        below_min_mult = 0
        above_max_mult = 0
        ref_shift = 0
        ncategories = length( breaks ) - 1
    end

    if ncategories < 1
        error( "Too few categories. Change boundedness or add breaks")
    end

    if leftequal
        for i in 1:length(x)
            if isna( x, i )
                refs[i] = 0
            elseif x2[i] <  breaks[1]
                refs[i] = below_min_mult
            elseif x2[i] > breaks[end]
                refs[i] = (n+ref_shift) * above_max_mult
            elseif x2[i] == breaks[end]
                if boundedness in [ :bounded, :boundedabove ]
                    refs[i] = ncategories
                else
                    refs[i] = n+ref_shift
                end
            else
                refs[i] = searchsortedlast(breaks, x2[i]) + ref_shift
            end
        end
    else
        for i in 1:length(x)
            if isna( x, i )
                refs[i] = 0
            elseif x2[i] < breaks[1]
                refs[i] = below_min_mult
            elseif x2[i] > breaks[end]
                refs[i] = (n+ref_shift) * above_max_mult
            else
                refs[i] = searchsortedfirst(breaks, x2[i])
            end
        end
    end

    if length( bucketstrs ) != 0
        if length( bucketstrs ) != ncategories
            error( "bucketstrs expected to have size " * string( ncategories ) *
                 ". Got " * string( length( bucketstrs ) ) )
        end
        if maximum( refs ) > ncategories
            maxref = maximum( refs )
            s = "ncategories < max refs \n maxref="  * string( maxref ) * "\n ncategories=" * string( ncategories )
            s *= "\n buckets = " * string( bucketstrs )
            idx = findfirst( refs, maxref )
            s *= "\n Example x = " * string(x2[idx] )
            s *= "\n breaks" * string( breaks )
            error( s )
        end
        return DataArrays.PooledDataArray(DataArrays.RefArray(refs), bucketstrs )
    end
    pool = cut_categories( S, breaks,
        boundedness = boundedness,
        leftequal   = leftequal,
        absolute    = absolute,
        rank        = rank,
        ranksep     = ranksep,
        label       = label,
        compact     = compact,
        reverse     = reverse,
        prefix=prefix, suffix=suffix, scale=scale, precision=precision,
        commas=commas,stripzeros=stripzeros,parens=parens,
        mixedfraction=mixedfraction,autoscale=autoscale,conversion=conversion
        )

    DataArrays.PooledDataArray(DataArrays.RefArray(refs), pool)
end

# quantile-based auto-breaks
# weighted quantile is not implemented
# use scale=100.0, suffix="%", to express the quantiles in percentages
function discretize{S<:Real}(x::AbstractArray{S,1}; quantiles = Float64[], ngroups::Int = 4, kwargs ... )
    if length( quantiles ) != 0
        if any( _ -> _ < 0.0 || _ > 1.0 , quantiles )
            error( "illegal quantile numbers outside [0,1]")
        end
        if !issorted(quantiles)
            sort!(quantiles)
        end
        if quantiles[1] != 0.0
            insert!( quantiles, 1, 0.0 )
        end

        if quantiles[end] != 1.0
            push!( quantiles, 1.0 )
        end
        bucketstrs = cut_categories( Float64, quantiles; boundedness = :bounded, kwargs... )
        discretize( x, quantile( x, quantiles ); bucketstrs = bucketstrs, boundedness = :bounded, kwargs... )
    else
        qs = [0:ngroups]/ngroups
        bucketstrs = cut_categories( Float64, qs; boundedness = :bounded, kwargs... )
        discretize( x, quantile( x, qs); bucketstrs = bucketstrs, boundedness = :bounded, kwargs... )
    end
end

# names are expected to be unique
# n is the maximum rank number to report. Actual outcome may depend on existence of a tie, and dense option
function topnames{S<:AbstractString,T<:Real}( name::AbstractArray{S,1}, measure::AbstractArray{T,1}, n::Int;
    absolute=false,
    ranksep = ". ",
    dense = true, # if there is a tie in the 2nd place, do we do "1,2,2,4", or "1,2,2,3"
    tol = 0,  # if absolute, what is the smallest contribution that we would consider
    others = "Others",
    parens = false # put parentheses around names with negative measure?
    )

    if absolute
        df = DataFrame( name = name, measure = measure, absmeasure = abs(measure) )
        if tol > 0 # filter out too small names
            dfsorted = sort( df[ df[ :absmeasure ] >= tol ], cols = [:absmeasure, :measure ], rev=[true,true] )
        else
            dfsorted = sort( df, cols = [ :absmeasure, :measure ], rev = [ true, true ] )
        end
    else
        df = DataFrame( name = name, measure = measure )
        dfsorted = sort( df, cols = [ :measure ], rev = [ true ] )
    end

    rankcount = 1
    rankwidth = length( string( n ) )
    nr = nrow( dfsorted )

    if !absolute
        pool = UTF8String[]
        refs = fill(zero(DataArrays.DEFAULT_POOLED_REF_TYPE), nr )
        lastval  = zero( T )
        lastrank = 0
        for r in 1:nr
            if isna( dfsorted[:measure], r )
                continue
            else
                val = dfsorted[ r, :measure ]
                if lastrank != 0 && lastval == val # tie
                    push!( pool, format( lastrank, width=rankwidth ) * ranksep * dfsorted[ r, :name ] )
                    refs[r] = length( pool )
                    if !dense
                        rankcount += 1
                    end
                elseif rankcount > n
                    break
                else
                    push!( pool, format( rankcount, width=rankwidth ) * ranksep * dfsorted[ r, :name ] )
                    lastrank = rankcount
                    lastval = val
                    refs[r] = length( pool )
                    rankcount += 1
                end
            end
        end
        dfsorted[ :rankstr ] = DataArrays.PooledDataArray(DataArrays.RefArray(refs), pool)
        rdf = dfsorted[ [:name, :rankstr ] ]
        jdf = join( df, rdf, on = :name, kind = :left )
    else
        rankedflag = fill( zero( Bool ), nr )
        lastval  = zero( T )
        lastrank = 0
        for r in 1:nr
            if isna( dfsorted[:measure], r )
                continue
            else
                val = dfsorted[ r, :measure ]
                if lastrank != 0 && lastval == val # tie
                    rankedflag[r] = true
                    if !dense
                        rankcount += 1
                    end
                elseif rankcount > n
                    break
                else
                    rankedflag[r] = true
                    lastrank = rankcount
                    lastval = val
                    rankcount += 1
                end
            end
        end
        dfsorted[ :rankedflag ] = rankedflag
        dfsorted2 = sort( dfsorted, cols = [ :measure ], rev = [ true ] )
        rankstr = DataArray(UTF8String,nr)
        rankcount = 1
        lastval  = zero( T )
        lastrank = 0
        for r in 1:nr
            if isna( dfsorted2[ :measure ], r )
                continue
            elseif dfsorted2[ r, :rankedflag ]
                val = dfsorted2[ r, :measure ]
                if lastrank != 0 && lastval == val # tie
                    if parens && val < 0
                        rankstr[r] = format( lastrank, width=rankwidth ) * ranksep * "("*dfsorted2[r,:name]*")"
                    else
                        rankstr[r] = format( lastrank, width=rankwidth ) * ranksep * dfsorted2[r,:name]
                    end
                    if !dense
                        rankcount += 1
                    end
                elseif rankcount > n
                    break
                else
                    if parens && val < 0
                        rankstr[r] = format( rankcount, width=rankwidth ) * ranksep * "("*dfsorted2[r,:name]*")"
                    else
                        rankstr[r] = format( rankcount, width=rankwidth ) * ranksep * dfsorted2[r,:name]
                    end
                    lastrank = rankcount
                    lastval = val
                    rankcount += 1
                end
            end
        end
        dfsorted2[ :rankstr ] = rankstr
        jdf = join( df, dfsorted2[ [:name, :rankstr] ], on = :name, kind=:left )
    end

    # replace NA with "others"
    ret = DataArrays.PooledDataArray( jdf[ :rankstr ] )
    push!( ret.pool, others )
    poollen = length( ret.pool )
    for i = 1:length( ret.refs )
        if ret.refs[i] == 0
            ret.refs[i] = poollen
        end
    end
    ret
end

import DataFrames.describe
export describe

function describe{T}( io, dv::Array{T,1} )
    describe( io, DataArray( dv ) )
end


text:
using Metatheory
using Metatheory.Rewriters
using Test
# using SymbolicUtils

array_theory = @theory x y f g M N begin
    #map(f,x)[n:m] = map(f,x[n:m]) # but does NOT commute with filter
    map(f,fill(x,N))            == fill(apply(f,x), N) # hmm
    # cumsum(fill(x,N))           == collect(x:x:(N*x))
    fill(x,N)[y]                --> x
    length(fill(x,N))           --> N
    reverse(reverse(x))         --> x
    sum(fill(x,N))              --> x * N
    map(f,reverse(x))           == reverse(map(f, x))
    filter(f,reverse(x))        == reverse(filter(f,x))
    reverse(fill(x,N))          == fill(x,N)
    filter(f, fill(x,N))        == (if apply(f, x); fill(x,N) else fill(x,0) end)
    filter(f, filter(g, x))     == filter(fand(f,g), x) # using functional &&
    cat(fill(x,N),fill(x,M))    == fill(x,N + M)
    cat(map(f,x), map(f,y))     == map(f, cat(x,y))
    map(f, cat(x,y))            == cat(map(f,x), map(f,y))
    map(f,map(g,x))             == map(f ∘ g, x)
    reverse( cat(x,y) )         == cat(reverse(y), reverse(x))
    map(f,x)[y]                 == apply(f,x[y])
    apply(f ⋅ g, x)             == apply(f, apply(g, x))

    map(f, reduce(g, x))        == mapreduce(f, g, x)
    map(f, foldl(g, x))         == mapfoldl(f, g, x)
    map(f, foldr(g, x))         == mapfoldr(f, g, x)
end

normalize_theory = @theory x y f g begin
    fand(f, g)  --> ((x) -> f(x) && g(x))
    apply(f, x) => Expr(:call, f, x)
end

fold_theory = @theory x y begin
    x::Number * y::Number => x*y
    x::Number + y::Number => x+y
    x::Number / y::Number => x/y
    # etc...
end

function expr_replace(env)
    function _expr_replace()
    end
    return _expr_replace
end

# function inline_lambda(ex::Expr)
#     exprhead(ex) != :call && return ex
#     fun = operation(ex)
#     !(fun isa Expr) && return ex
#     exprhead(fun) != :(->) && return ex

#     f_arg = arguments(fun)[1]
#     !(f_arg isa Symbol) && return ex
#     f_body = arguments(fun)[2]

#     args = arguments(ex)
#     length(args) != 1 && return ex

#     args[1]
#     body = Metatheory.Syntax.rmlines(substitute(f_body, Dict(f_arg => args[1]); fold=false))
#     if exprhead(body) == :block && length(arguments(body)) == 1
#         return arguments(body)[1]
#     else
#         return body
#     end
# end

# inline_lambda(x) = x

params = SaturationParams()

function stream_optimize(ex)
    g = EGraph(ex)
    rep = saturate!(g, array_theory, params)
    @info rep
    ex = extract!(g, astsize;) # TODO cost fun with asymptotic complexity
    ex = Fixpoint(Postwalk(Chain(vcat(normalize_theory, fold_theory))))(ex)
    return ex
end

build_fun(ex) = eval(:(()->$ex))

@testset "Stream Fusion" begin
    ex = :( map(x -> 7 * x, fill(3,4)))
    opt = stream_optimize(ex)
    @test opt == :(fill(((x-> 7x))(3), 4))

    ex = :( map(x -> 7 * x, fill(3,4) )[1])
    opt = stream_optimize(ex)
    @test opt == :((x->7x)(3))
end


text:
function test_binarysearch_all()
    print(rpad("Algorithm.BinarySearch...", 50, ' '))

    @test [binarysearch(x -> x, i, 1, 100) for i in 1:100] == [(i,i) for i in 1:100]
    @test [binarysearch(x -> 2x, i, 0, 10) for i in 2:10] == [(1,1); (1,2); (2,2); (2,3); (3,3); (3,4); (4,4); (4,5); (5,5)]

    @test_throws DomainError binarysearch(x -> x, 0, 1, 10)
    @test_throws DomainError binarysearch(x -> x, 11, 1, 10)

    println("PASS")
end


text:
module Diagnostics

export
    Average, TimeSeries, FieldMaximum,
    CFL, AdvectiveCFL, DiffusiveCFL, NaNChecker,
    run_diagnostic

using Oceananigans
using Oceananigans.Operators

using Oceananigans: AbstractDiagnostic

include("nan_checker.jl")
include("average.jl")
include("time_series.jl")
include("field_maximum.jl")
include("cfl.jl")

end


text:
using Printf
import Flux: σ
using ModelingToolkit
using GalacticOptim
using Optim
using DiffEqFlux
using NeuralPDE
using Quadrature, Cubature, Cuba
using Plots

@parameters t,x,y
@variables c(..)
@derivatives Dt'~t
@derivatives Dxx''~x
@derivatives Dyy''~y
@derivatives Dx'~x
@derivatives Dy'~y

# Parameters

u = 1.0
v = 1.0
v_vector = (u, v)
R = 0
D = 0 # 0.1 # diffusion
t_max = 2.0
x_min = -1.0
x_max = 1.0
y_min = -1.0
y_max = 1.0

# exp(-(x^2+y^2)/0.1)
# div(v_vector * c) = dx(uc) + dy(vc)
# Equations, initial and boundary conditions
eqs = [ Dt(c(t,x,y)) ~ D * (Dxx(c(t,x,y)) + Dyy(c(t,x,y))) - (u*Dx(c(t,x,y)) + v*Dy(c(t,x,y))) + R]

bcs = [
        c(0, x, y) ~ exp(-(x^2+y^2)/0.1), #cos(π*x) * cos(π*y) + 1.0,
        c(t, x_min, y) ~ c(t, x_max, y),
        c(t, x, y_min) ~ c(t, x, y_max)
]

# Space and time domains
domains = [t ∈ IntervalDomain(0.0,t_max),
        x ∈ IntervalDomain(x_min,x_max),
        y ∈ IntervalDomain(y_min,y_max)
]

# Discretization
nx = 32
# ny = 32
dx = (x_max-x_min) / (nx - 1)
# dy = (y_max-y_min) / (ny -1)
dt = 0.01

# Neural network
dim = length(domains)
output = length(eqs)
hidden = 32

chain = FastChain( FastDense(dim, hidden, tanh),
                    FastDense(hidden, hidden, tanh),
                    FastDense(hidden, 1))

# strategy = GridTraining(dx=[dt,dx,dy])
strategy = StochasticTraining(900)

discretization = PhysicsInformedNN(chain, strategy=strategy)

pde_system = PDESystem(eqs, bcs, domains, [t,x,y], [c])
prob = discretize(pde_system,discretization)

cb = function (p,l)
    println("Current loss is: $l")
    return false
end

res = GalacticOptim.solve(prob,ADAM();cb=cb,maxiters=100)

# Plots

phi = discretization.phi

initθ = discretization.initθ

acum =  [0;accumulate(+, length.(initθ))]
sep = [acum[i]+1 : acum[i+1] for i in 1:length(acum)-1]
minimizers = [res.minimizer[s] for s in sep]
ts,xs,ys = [domain.domain.lower:dx:domain.domain.upper for domain in domains]

# Plot initial condition
c_predict = reshape([ phi([0, x, y], res.minimizer)[1] for x in xs for y in ys], length(xs), length(ys))
plot(xs, ys, c_predict)

# Animate
anim = @animate for (i, t) in enumerate(0:dt:t_max)
    @info "Animating frame $i..."
    c_predict = reshape([phi([t, x, y], res.minimizer)[1] for x in xs for y in ys], length(xs), length(ys))
    title = @sprintf("Advection-diffusion t = %.3f", t)
    heatmap(xs, ys, c_predict, label="", title=title , xlims=(-1, 1), ylims=(-1, 1), color=:thermal, clims=(0, 1))
end

gif(anim, "advection_diffusion_2d_pinn.gif", fps=15)

# Plot correct solution
using JLD2
file = jldopen("advection_diffusion/simulation/cosine_advection_diffusion.jld2")
iterations = parse.(Int, keys(file["timeseries/t"]))

anim = @animate for (i, iter) in enumerate(iterations)
    @info "Animating frame $i..."
    Hx = file["grid/Hx"]
    x = file["grid/xC"][1+Hx:end-Hx]
    t = file["timeseries/t/$iter"]
    c = file["timeseries/c/$iter"][:]

    title = @sprintf("Advection-diffusion t = %.3f", t)
    p = plot(x, c .+ 1, linewidth=2, title=title, label="Oceananigans",
             xlabel="x", ylabel="Tracer", xlims=(-1, 1), ylims=(0, 2))

    c_predict = reshape([phi([t, x], res.minimizer)[1] for x in xs], length(xs))

    plot!(p, x, c_predict, linewidth=2, label="Neural PDE")
end

gif(anim, "tanh_hidden32_it100_comparison.gif", fps=15)

max_stars_repo_name: drozzy/NeuralOceananigans

max_issues_repo_name: drozzy/NeuralOceananigans

max_forks_repo_name: drozzy/NeuralOceananigans


text:
using OceanWaveSpectralFitting
using Test, FiniteDifferences
import OceanWaveSpectralFitting.WhittleLikelihoodInference: grad_sdf, hess_sdf, grad_acv, hess_acv

## define useful functions
# Function for approximating the gradient using finite differences.
approx_gradient(func, x) = reinterpret(ComplexF64,jacobian(central_fdm(5, 1), func, x)[1])
approx_gradient_uni(func, x) = complex.(permutedims(jacobian(central_fdm(5, 1), func, x)[1]))
# Function for approximating the hessian using finite differences.
function approx_hessian(grad, x)
    hess = zeros(ComplexF64, size(grad(x),1), length(x)*(length(x)+1)÷2)
    count = 1
    @views for i ∈ 1:length(x)
        hess[:,count:count+length(x)-i] = approx_gradient(y->grad(y)[:,i], x)[:,i:end]
        count += length(x)-i+1
    end
    return hess
end
function approx_hessian_uni(grad, x)
    full_hess = approx_gradient(grad, x)
    hess = ones(ComplexF64, size(full_hess,1)*(size(full_hess,1)+1)÷2)
    for i ∈ 1:size(full_hess,1), j ∈ 1:i
        hess[WhittleLikelihoodInference.indexLT(i,j,size(full_hess,1))] = full_hess[i,j]
    end
    return hess
end

@testset "OceanWaveSpectralFitting.jl" begin
    include("models/univariate/JONSWAP_test.jl")
    include("models/univariate/generaljonswap_test.jl")
    include("models/univariate/gaussian_test.jl")
    include("models/multivariate/JS_BWG_HNE_test.jl")
    include("models/multivariate/JS_BWG_HNE_DL_test.jl")
    include("models/multivariate/JS_WG_HNE_DL_test.jl")
    include("models/multivariate/GJS_BWG_HNE_DL_test.jl")
    include("models/multivariate/Gaussian_WG_HNE_DL_test.jl")
end


text:
using Base.Filesystem

include("../src/Triangulation.jl")

# create_nodes
nodes = Triangulation.create_nodes(10, [0., 1.], [0., 1.])

# in_circumcircle
Triangulation.in_circumcircle([1, 2, 3], 4, nodes)

# is_clockwise
Triangulation.is_clockwise([0., 0.], [1.,0.], [0., 1.])

# create_supertriangle
Triangulation.create_supertriangle(nodes)

# require_flip
AB = Triangulation.Edge([1, 2], [])
BC = Triangulation.Edge([2, 3], [])
CD = Triangulation.Edge([3, 4], [])
BD = Triangulation.Edge([2, 4], [])
AD = Triangulation.Edge([1, 4], [])
edges = [AB, BC, CD, BD, AD]
triangles = [Triangulation.Triangle([1, 4, 5]), Triangulation.Triangle([2, 3, 4])]
width = 1.
height = width/tan(pi/3)
nodes = [[0., -height], [width, 0.],  [0., height], [-width, 0.]]
Triangulation.require_flip(1, 2, triangles, edges, nodes) == true

# compute_angle
nodes_2 = [[0., 0.], [1., 0.], [0., tan(pi/3)]]
edge1 = Triangulation.Edge([1, 2], [])
edge2 = Triangulation.Edge([1, 3], [])
Triangulation.compute_angle(edge1, edge2, nodes_2)

# get_tri_nodes
edges = [
    Triangulation.Edge([1, 2], []),
    Triangulation.Edge([2, 3], []),
    Triangulation.Edge([3, 1], [])
]
triangles = [Triangulation.Triangle([1, 2, 3])]
Triangulation.get_tri_node_ids(1, triangles, edges)

# inside_triangle
Triangulation.inside_triangle([1, 2, 3], 4, nodes)

# get_triangulation_nodes
Triangulation.get_triangulation_nodes(edges)

# flip_triangulation!
AB = Triangulation.Edge([1, 2], [1, 0])
BC = Triangulation.Edge([2, 3], [2, 0])
CD = Triangulation.Edge([3, 4], [2, 0])
BD = Triangulation.Edge([2, 4], [1, 2])
AD = Triangulation.Edge([1, 4], [1, 0])
edges = Dict{Int, Triangulation.Edge}(1=>AB, 2=>BC, 3=>CD, 4=>BD, 5=>AD)
triangles = Dict{Int, Triangulation.Triangle}(
    1 => Triangulation.Triangle([1, 4, 5]),
    2 => Triangulation.Triangle([2, 3, 4])
)
width = 1.
height = width/tan(pi/3)
nodes = [[0., -height], [width, 0.],  [0., height], [-width, 0.]]
Triangulation.flip_triangulation!(triangles, edges, nodes, true, "./_tmp_.txt", init_flip_step=0)

# plot_step
Triangulation.plot_step("./_tmp_flip_1.txt", false)

# triangulate
nodes = Triangulation.create_nodes(5, [0., 1.], [0., 1.])
Triangulation.triangulate(nodes, false, "./_tmp_.txt")

# parse_args
Triangulation.parse_args("-i 5 --long=3 --long2 3")

# print_help
Triangulation.print_help()

# get_arg_value
Triangulation.get_arg_value(["c", "d"], Dict("c" => "ok", "d" => "ok"))

# get_main_parameters
Triangulation.get_main_parameters("-c -o myOutput -h -p -s")
Triangulation.get_main_parameters("--create-nodes --output=myOutput --help --plot --save-steps")
Filesystem.touch("./_tmp_test.txt")
Triangulation.get_main_parameters("-f ./_tmp_test.txt")
Triangulation.get_main_parameters("-l ./_tmp_test*.txt")

# main
Triangulation.main(["--help"])
Triangulation.main(["-l ./results*.txt"])

# julia_main
Triangulation.julia_main()

# remove tmp files
Filesystem.rm("./_tmp_test.txt")
Filesystem.rm("./_tmp_flip_1.txt")
Filesystem.rm("./_tmp_flip_1.png")


text:
"""
    HullSearch(tolerance::Float64)

HullSearch, SGSV Combine HullReach, performs over-approximated reachability analysis to compute the over-approximated output reachable set for a network.

# Problem requirement
1. Network: any depth, any activation that is monotone
2. Input: `HPolytope`
3. Output: `HPolytope`

# Return
`BasicResult`

# Method
Search and Reachability

# Property
Sound but not complete.
"""
@with_kw struct HullTree
    tolerance::Float64 = 1.0
end

# This is the main function
function solve(solver::HullTree, problem::Problem)
    input = problem.input
    stack = Vector{Hyperrectangle}(undef, 0)
    push!(stack, input)
    #count = 1
    while !isempty(stack)
        interval = popfirst!(stack)
        reach = forward_network(solver, problem.network, interval)
        if issubset(reach, problem.output)
            continue
        else
            if get_largest_width(interval) > solver.tolerance
                sections = bisect(interval)
                for i in 1:2
                    if isborder(sections[i], problem.input)
                        push!(stack, sections[i])
                        #count += 1
                    end
                end
            else
                return BasicResult(:violated)
            end
        end
    end
    #print("\n$(count)\n")
    return BasicResult(:holds)
end

#to determine whether x has intersection with any border of y
function isborder(x::Hyperrectangle, y::Hyperrectangle)
    x_lower, x_upper = low(x), high(x)
    y_lower, y_upper = low(y), high(y)
    for i in 1:lastindex(x_lower)
        if x_lower[i] == y_lower[i] || x_upper[i] == y_upper[i]
            return true
        end
    end
    return false
end

function forward_layer(solver::HullTree, L::Layer, input::Hyperrectangle)
    (W, b, act) = (L.weights, L.bias, L.activation)
    center = zeros(size(W, 1))
    gamma  = zeros(size(W, 1))
    for j in 1:size(W, 1)
        node = Node(W[j,:], b[j], act)
        center[j], gamma[j] = forward_node(solver, node, input)
    end
    return Hyperrectangle(center, gamma)
end

function forward_node(solver::HullTree, node::Node, input::Hyperrectangle)
    output    = node.w' * input.center + node.b
    deviation = sum(abs.(node.w) .* input.radius)
    βmax = node.act(output + deviation)
    βmin = node.act(output - deviation)
    return ((βmax + βmin)/2, (βmax - βmin)/2)
end


text:
function readinput(filename)
    return map(x -> split(x, " = "), readlines(filename))
end

function masks(mask)
    and, or, x = UInt(0), UInt(0), UInt(0)
    for c in mask
        and, or, x = [and, or, x] .<<= 1
        if c == '0'
            and += 1
        elseif c == '1'
            or += 1
        else
            x += 1
        end
    end
    return (~and, or, x)
end

function allcombinations(addr, floating)
    if floating == 0
        return addr
    else
        attempt = UInt(1)
        while true
            if attempt == floating & attempt
                return cat(allcombinations(addr, floating &= ~attempt), allcombinations(xor(addr, attempt), floating &= ~attempt); dims=1)
            end
            attempt <<= 1
        end
    end
end

function memoryallocation(instructions)
    P1, P2 = Dict(), Dict()
    and = or = x = 0
    for instruction in instructions
        cmd, val = instruction[1], instruction[2]
        if cmd == "mask"
            and, or, x = masks(val)
        else
            addr, val = parse(UInt, match(r"^mem\[(\d+)\]$", cmd)[1]), parse(UInt, val)
            P1[addr] = (val | or) & and
            for loc in allcombinations(addr | or, x)
                P2[loc] = val
            end
        end
    end
    return P1, P2
end

instructions = readinput("d14.in")
P1, P2 = memoryallocation(instructions)

println("P1: ", sum(values(P1)))
println("P2: ", sum(values(P2)))

max_stars_repo_name: eduellery/adventofcode

max_issues_repo_name: eduellery/adventofcode

max_forks_repo_name: eduellery/adventofcode


text:
module CustomAlphabetTest

using Hashids
using Test

const TESTDATA = [
    [0],
    [1, 2, 3],
    [10000000000, 0, 0, 0, 999999999999999],
    [827364827638428762424, 0x410498019283098083],
    [340282366920938463463374607431768211456]
]

function test_with_custom_alphabet(alphabet::String)
    config = Hashids.configure(alphabet=alphabet)
    for numbers in TESTDATA
        @test decode(config, encode(config, numbers)) == numbers
    end
end

@testset "Custom Alphabet" begin

@testset "Worst Alphabet" begin
    test_with_custom_alphabet("cCsSfFhHuUiItT01")
end  # of "Worst Alphabet"

@testset "Containing Spaces" begin
    test_with_custom_alphabet("cCsSfFhH uUiItT01")
end  # of "Containing Spaces"

@testset "a Half be Separators" begin
    test_with_custom_alphabet("abdegjklCFHISTUc")
end  # of "a Half be Separators"

@testset "2 Separators" begin
    test_with_custom_alphabet("abdegjklmnopqrSF")
end  # of "2 Separators"

@testset "No Separators" begin
    test_with_custom_alphabet("abdegjklmnopqrvwxyzABDEGJKLMNOPQRVWXYZ1234567890")
end  # of "No Separators"

@testset "Long Alphabet" begin
    test_with_custom_alphabet(raw"abcdefghijklmnopqrstuvwxyzABCDEFGHIJKLMNOPQRSTUVWXYZ1234567890`~!@#$%^&*()-_=+\\|'\";:/?.>,<{[}]")
end  # of "Long Alphabet"

@testset "Weird Alphabet" begin
    test_with_custom_alphabet(raw"`~!@#$%^&*()-_=+\\|'\";:/?.>,<{[}]")
end  # of "Weird Alphabet"

@testset "Unicode Chars" begin
    test_with_custom_alphabet("😀😁😂🤣😃😄😅😆😉😊😋😎😍😘🥰😗😙😚")
end  # of "Unicode Chars"

@testset "Complex Unicode Chars" begin
    test_with_custom_alphabet("🤺👩🏿🦳🛁👩🏻🦷🤦♂🐁☝🏼✍🏾👉🏽🇸🇰❤🍭")
end  # of "Complex Unicode Chars"

end  # of "Custom Alphabet"

end  # module