text:
"""
    nonparametricbootstrap([rng::AbstractRNG,] nsamp::Integer, m::LinearMixedModel;
                           use_threads=false, blup_method=ranef, β=coef(morig))

Perform `nsamp` nonparametric bootstrap replication fits of `m`, returning a `MixedModelBootstrap`.

The default random number generator is `Random.GLOBAL_RNG`.

`GeneralizedLinearMixedModel` is currently unsupported.

# Named Arguments
`use_threads` determines whether or not to use thread-based parallelism.

!!! note
    Note that `use_threads=true` may not offer a performance boost and may even
    decrease peformance if multithreaded linear algebra (BLAS) routines are available.
    In this case, threads at the level of the linear algebra may already occupy all
    processors/processor cores. There are plans to provide better support in coordinating
    Julia- and BLAS-level threads in the future.

!!! warning
    The PRNG shared between threads is locked using [`Threads.SpinLock`](@ref), which
    should not be used recursively. Do not wrap `nonparametricbootstrap` in an outer `SpinLock`.

`hide_progress` can be used to disable the progress bar. Note that the progress
bar is automatically disabled for non-interactive (i.e. logging) contexts.

`blup_method` provides options for how/which group-level effects are passed for resampling.
The default `ranef` uses the shrunken conditional modes / BLUPs. Unshrunken estimates from
ordinary least squares (OLS) can be used with `olsranef`. There is no shrinkage of the
group-level estimates with this approach, which means singular estimates can be avoided.
However, if the design matrix for the random effects is rank deficient (e.g., through the use
of `MixedModels.fulldummy` or missing cells in the data), then this method will fail.
See [`olsranef`](@ref) and `MixedModels.ranef` for more information.

# Method

The method implemented here is based on the approach given in Section 3.2 of:
Carpenter, J.R., Goldstein, H. and Rasbash, J. (2003).
A novel bootstrap procedure for assessing the relationship between class size and achievement.
Journal of the Royal Statistical Society: Series C (Applied Statistics), 52: 431-443.
https://doi.org/10.1111/1467-9876.00415
"""
function nonparametricbootstrap(
    rng::AbstractRNG,
    n::Integer,
    morig::LinearMixedModel{T};
    use_threads::Bool=false,
    hide_progress=false,
    blup_method=ranef,
    β=coef(morig),
) where {T}
    # XXX should we allow specifying betas and blups?
    #     if so, should we use residuals computed based on those or the observed ones?
    βsc, θsc = similar(morig.β), similar(morig.θ)
    p, k = length(βsc), length(θsc)
    m = deepcopy(morig)

    β_names = (Symbol.(fixefnames(morig))..., )
    rank = length(β_names)

    blups = blup_method(morig)
    resids = residuals(morig)#, blups)
    reterms = morig.reterms
    scalings = inflation_factor(morig, blups, resids)
    # we need arrays of these for in-place operations to work across threads
    m_threads = [m]
    βsc_threads = [βsc]
    θsc_threads = [θsc]

    if use_threads
        Threads.resize_nthreads!(m_threads)
        Threads.resize_nthreads!(βsc_threads)
        Threads.resize_nthreads!(θsc_threads)
    end
    # we use locks to guarantee thread-safety, but there might be better ways to do this for some RNGs
    # see https://docs.julialang.org/en/v1.3/manual/parallel-computing/#Side-effects-and-mutable-function-arguments-1
    # see https://docs.julialang.org/en/v1/stdlib/Future/index.html
    rnglock = Threads.SpinLock()
    samp = replicate(n; use_threads=use_threads, hide_progress=hide_progress) do
        tidx = use_threads ? Threads.threadid() : 1
        model = m_threads[tidx]
        local βsc = βsc_threads[tidx]
        local θsc = θsc_threads[tidx]
        lock(rnglock)
        model = resample!(rng, model; β=β, blups=blups, resids=resids, scalings=scalings)
        unlock(rnglock)
        refit!(model)
        (
         objective = model.objective,
         σ = model.σ,
         β = NamedTuple{β_names}(fixef!(βsc, model)),
         se = SVector{p,T}(stderror!(βsc, model)),
         θ = SVector{k,T}(getθ!(θsc, model)),
        )
    end
    MixedModelBootstrap(
        samp,
        deepcopy(morig.λ),
        getfield.(morig.reterms, :inds),
        copy(morig.optsum.lowerbd),
        NamedTuple{Symbol.(fnames(morig))}(map(t -> (t.cnames...,), morig.reterms)),
    )
end

function nonparametricbootstrap(nsamp::Integer, m::LinearMixedModel, args...; kwargs...)
    return nonparametricbootstrap(Random.GLOBAL_RNG, nsamp, m, args...; kwargs...)
end

function nonparametricbootstrap(rng::AbstractRNG, n::Integer,
                                morig::GeneralizedLinearMixedModel;
                                use_threads::Bool=false) where {T}

    throw(ArgumentError("GLMM support is not yet implemented"))
end

resample!(model::LinearMixedModel, args...; kwargs...) =
    resample!(Random.GLOBAL_RNG, model, args...; kwargs...)

"""
    resample!([rng::AbstractRNG,] model::LinearMixedModel;
              β=coef(model), blups=ranef(model), resids=residuals(model,blups)
              scalings=inflation_factor(model))

Simulate and install a new response using resampling at the observational and group level.

At both levels, resampling is done with replacement. At the observational level,
this is resampling of the residuals, i.e. comparable to the step in the classical
nonparametric bootstrap for OLS regression. At the group level, samples are done
jointly for all  terms ("random intercept and random slopes", conditional modes)
associated with a particular level of a blocking factor. For example, the
predicted slope and intercept for a single participant or item are kept together.
This clumping in the resampling procedure is necessary to preserve the original
correlation structure of the slopes and intercepts.

In addition to this resampling step, there is also an inflation step. Due to the
shrinkage associated with the random effects in a mixed model, the variance of the
conditional modes / BLUPs / random intercepts and slopes is less than the variance
estimated by the model and displayed in the model summary or via `MixedModels.VarCorr`.
This shrinkage also impacts the observational level residuals. To compensate for this,
the resampled residuals and groups are scale-inflated so that their standard deviation
matches that of the estimates in the original model. The default inflation factor is
computed using [`inflation_factor`](@ref) on the model.

See also [`nonparametricbootstrap`](@ref) and `MixedModels.simulate!`.

!!! warning
    This method has serious limitations for singular models because resampling from
    a distribution with many zeros (e.g. the random effects components with zero variance)
    will often generate new data with even less variance.

# Reference
The method implemented here is based on the approach given in Section 3.2 of:
Carpenter, J.R., Goldstein, H. and Rasbash, J. (2003).
A novel bootstrap procedure for assessing the relationship between class size and achievement.
Journal of the Royal Statistical Society: Series C (Applied Statistics), 52: 431-443.
https://doi.org/10.1111/1467-9876.00415
"""
function resample!(rng::AbstractRNG, model::LinearMixedModel{T};
                   β=coef(model),
                   blups=ranef(model),
                   resids=residuals(model,blups),
                   scalings=inflation_factor(model)) where {T}

    reterms = model.reterms
    y = response(model) # we are now modifying the model

    # sampling with replacement
    sample!(rng, resids, y; replace=true)

    # inflate these to be on the same scale as the empirical variation instead of the MLE
    y .*= last(scalings)

    for (inflation, re, trm) in zip(scalings[2:end], blups, reterms)
        npreds, ngrps = size(re)
        # sampling with replacement
        samp = sample(rng, 1:ngrps, ngrps; replace=true)
         # allocate now to avoid allocating later
         # while taking advantage of LowerTriangular lmul!
        newre = re[:, samp]

        # this just multiplies the Z matrices by the BLUPs
        # and add that to y
        mul!(y, trm, lmul!(inflation, newre), one(T), one(T))
        # XXX inflation is resampling invariant -- should we move it out?
    end

    mul!(y, model.X, β, one(T), one(T))

    # mark model as unfitted
    model.optsum.feval = -1

    return model
end


text:
module NodeMap

using DataFrames
using CSV

using LightGraphs
import LightGraphs.loadgraph
import GraphIO.EdgeList.EdgeListFormat

function loadedgelist(io::IO, gname::String)
    srcs = Vector{String}()
    dsts = Vector{String}()
    while !eof(io)
        line = strip(chomp(readline(io)))
        if !startswith(line, "#") && (line != "")
            r = r"([\w-]+)[\s,]+([\w-]+)"
            src_s, dst_s = match(r, line).captures
            push!(srcs, src_s)
            push!(dsts, dst_s)
        end
    end
    vxset = unique(vcat(srcs, dsts))
    vxdict = Dict{String,Int}()
    for (v, k) in enumerate(vxset)
        vxdict[k] = v
    end

    n_v = length(vxset)
    g = LightGraphs.DiGraph(n_v)
    for (u, v) in zip(srcs, dsts)
        add_edge!(g, vxdict[u], vxdict[v])
    end
    return g, vxdict
end

loadgraph(io::IO, gname::String, ::EdgeListFormat) = loadedgelist(io, gname)

function create_node_map(graph_file)
    G, node_map = loadgraph(graph_file, EdgeListFormat())
    df = DataFrame(original_id = [keys(node_map)...], julia_id = [values(node_map)...])
    save_file = splitext(graph_file)[1] * "-nodemap.csv"
    CSV.write(save_file, df)
end

end

ext: jl

max_stars_repo_name: dataspider/momo

max_issues_repo_name: dataspider/momo

max_forks_repo_name: dataspider/momo


text:
function mirror(A,d::Int=3)
    w = div(d-1,2) #width of window, d is an ODD integer
    mA = [A[(w+1):-1:2];A;A[lastindex(A)-1:-1:(length(A)-w)]]
end

function find_extrema(A::Array,d::Int=3)
    mA = mirror(A)
    count = 0
    maxs = Int[]
    mins = Int[]
    exts = Int[]
    for i in 1:lastindex(mA)-d+1
        win = view(mA,i:i+d-1)
        if (win[2] > win[1] && win[2] > win[3])
                count+=1
                push!(maxs,i)
                push!(exts,i)
        elseif (win[2] < win[1] && win[2] < win[3])
                count+=1
                push!(mins,i)
                push!(exts,i)
        end
    end
    maxs, mins,exts
end

function find_extrema_count(A::Array,d::Int=3)
    dd = diff(sign.(diff(A)))
    return length(dd[dd .!= 0]) + 2
end

function find_extrema_minmax(A::Array,d::Int=3)
    dd = diff(sign.(diff(A)))
    return findall(dd .< 0), findall(dd .> 0)
end

function stream_minmax(env1, env2, A, d::Int)
    a = mirror(A, d)
    N = length(a)
    upper = UInt32[] #buffer for indices
    lower = UInt32[]
    push!(upper,1)
    push!(lower,1)
    w = d
    for i in 2:lastindex(a)
        if i >= w
            env1[i-w+1] = a[upper[1]]
            env2[i-w+1] = a[lower[1]]
        end

        if a[i] > a[i-1]
            pop!(upper)
            # remove maxima from buffer that are less than our new one
            while isempty(upper) != true
                if a[i] <= a[upper[end]]
                    break
                end
                pop!(upper)
            end
        else
            pop!(lower)
            # remove minima from buffer that are greater than our new one
            while isempty(lower) != true
                if a[i] >= a[lower[end]]
                    break
                end
                pop!(lower)
            end
        end
        push!(upper,i)
        push!(lower,i)
        if i == w + upper[1]
            popfirst!(upper)
        elseif i == w + lower[1]
            popfirst!(lower)
        end
        env1[end-w] = a[upper[1]]
        env2[end-w] = a[lower[1]]
    end
    return nothing
end

function moving_average(A, d::Int=3)
    A = mirror(A,d)
    T = typeof(one(eltype(A))/1)
    ret = Vector{T}(undef, length(A) - d + 1)
    id = 1 / d
    s = sum(view(A, 1:d))
    ret[1] = s * id
    @inbounds for n = 1:length(ret)-1
        s += A[n+d] - A[n]
        ret[n+1] = s * id
    end
    return ret
end


text:
push!(LOAD_PATH, "../../DeepCorrection")
push!(LOAD_PATH, "../../DeepCorrection/src")
push!(LOAD_PATH, "../../VICAS/src/")
push!(LOAD_PATH, "../TrajectoryViz")

include("../TrajectoryViz/trajViz.jl")

headingSet = deg2rad.([30 * i - 360 for i in 1:23])

function sensitivity(;
	numInt::Int64=1,
	numSamples::Int64=10000,
	CASType::Symbol=:correctedSector,
	policyFile::AbstractString="",
	printOption::Bool=false,
	correction_weight::Float64=0.5)

	problem_file = "../../DeepCorrection/logs/" * policyFile * "/final_problem.jld"
	weights_file = "../../DeepCorrection/logs/" * policyFile * "/final_weights.jld"
	policy = restore(problem_file=problem_file, weights_file=weights_file)
	policy.env.correction_weight = correction_weight

	x0 = 0.0
	y0 = 0.0
	x_dest = 10000.0
	y_dest = 0.0
	ownship = initOwnship(x0, y0, x_dest, y_dest, CASType, policy)

	sensitiveRate = 0
	for t in 1 : numSamples
		intCoords = Vector{Float64}()
		for i in 1 : numInt
			intruder_x = rand(-1000 : 50 : 1000)
			y = sqrt(SENSING_RANGE^2 - intruder_x^2)
			intruder_y = rand(trunc(- y) : 50 : trunc(y))
			intruder_heading = rand(headingSet)
			intCoords = vcat(intCoords, [intruder_x, intruder_y, intruder_heading])
		end
		getState!(ownship, numInt, intCoords)
		getAdvisory!(ownship)
		if ownship.advisory != :COC
			sensitiveRate += 1
		end
		if printOption == true
			println(intCoords, ", ", ownship.advisory)
		end
	end
	return sensitiveRate / numSamples
end

max_stars_repo_name: parachutel/DeepCorrectionCAS


text:
using Test, LinearMaps, LinearAlgebra, SparseArrays

@testset "transpose/adjoint" begin
    CS = @inferred LinearMap{ComplexF64}(cumsum, x -> reverse(cumsum(reverse(x))), 10; ismutating=false)
    for transform in (adjoint, transpose)
        @test transform(CS) != CS
        @test CS != transform(CS)
        @test transform(transform(CS)) == CS
        @test LinearMaps.MulStyle(transform(CS)) === LinearMaps.MulStyle(CS)
    end
    @test occursin("10×10 LinearMaps.TransposeMap{$(eltype(CS))}", sprint((t, s) -> show(t, "text/plain", s), transpose(CS)))
    @test occursin("10×10 LinearMaps.AdjointMap{$(eltype(CS))}", sprint((t, s) -> show(t, "text/plain", s), adjoint(CS)))
    @test !(transpose(CS) == adjoint(CS))
    @test !(adjoint(CS) == transpose(CS))
    M = Matrix(CS)
    @test M == LowerTriangular(ones(10,10)) == mul!(copy(M), CS, 1, true, false)
    x = rand(ComplexF64, 10); w = rand(ComplexF64, 10)
    X = rand(ComplexF64, 10, 3); W = rand(ComplexF64, 10, 3)
    α, β = rand(ComplexF64, 2)
    for transform in (adjoint, transpose)
        @test convert(AbstractMatrix, transform(CS)) == transform(M)
        @test sparse(transform(CS)) == transform(M)
        @test transform(CS) * x ≈ transform(M)*x
        @test mul!(copy(w), transform(CS), x) ≈ transform(M)*x
        @test mul!(copy(W), transform(CS), X) ≈ transform(M)*X
        @test mul!(copy(w), transform(CS), x, α, β) ≈ transform(M)*x*α + w*β
        @test mul!(copy(W), transform(CS), X, α, β) ≈ transform(M)*X*α + W*β
    end
    for transform1 in (adjoint, transpose), transform2 in (identity, adjoint, transpose)
        @test transform1(transform1(CS)) === CS
        @test transform2(transform1(transform2(CS))) === transform1(CS)
        @test transform2(transform1(CS)) * x ≈ transform2(transform1(M))*x
        @test mul!(copy(w), transform2(transform1(CS)), x) ≈ transform2(transform1(M))*x
        @test mul!(copy(W), transform2(transform1(CS)), X) ≈ transform2(transform1(M))*X
        @test mul!(copy(M), transform2(transform1(CS)), 2) ≈ transform2(transform1(M))*2
        @test mul!(copy(w), transform2(transform1(CS)), x, α, β) ≈ transform2(transform1(M))*x*α + w*β
        @test mul!(copy(W), transform2(transform1(CS)), X, α, β) ≈ transform2(transform1(M))*X*α + W*β
        @test mul!(copy(M), transform2(transform1(CS)), 2, α, β) ≈ transform2(transform1(M))*(2*α) + M*β
    end

    id = @inferred LinearMap(identity, identity, 10; issymmetric=true, ishermitian=true, isposdef=true)
    for transform in (adjoint, transpose)
        @test transform(id) == id
        @test id == transform(id)
        for prop in (issymmetric, ishermitian, isposdef)
            @test prop(transform(id))
        end
    end
end


text:
using ArgParse
using ConvHMM
using BSON
using LinearAlgebra
using MarkovModels
using Random
using YAML

function getargs()
    s = ArgParseSettings()
    @add_arg_table! s begin
        "--prior-strength", "-p"
            help = "strength of the prior, i.e. pseudo-counts"
            arg_type = Int
            default = 1
        "--filter-order", "-f"
            help = "order of the filter"
            arg_type = Int
            default = 5
        "--n-utts", "-n"
            help = "use `n` random utterances to compute the statistics"
            arg_type = Int
            default = 100
        "--seed", "-s"
            help = "seed of the random generator"
            arg_type = Int
            default = 1234
         "--verbose", "-v"
            help = "enable debug messages"
            action = :store_true
        "uttids"
            help = "list of training utterances"
            arg_type = String
            required = true
        "feadir"
            help = "features directory"
            arg_type = String
            required = true
        "conf-hmm"
            help = "HMM configuration file"
            arg_type = String
            required = true
        "units"
            help = "list of \"unit\" for which to create a pdf"
            arg_type = String
            required = true
        "outdir"
            help = "output directory where will be stored the HMMs and the emissions"
            arg_type = String
            required = true
    end
    args = parse_args(s)
end

const args = getargs()
const uttids = readlines(args["uttids"])
const feadir = args["feadir"]
const confhmm = args["conf-hmm"]
const units = readlines(args["units"])
const outdir = args["outdir"]

run(`mkdir -p $outdir`)

if args["verbose"] ENV["JULIA_DEBUG"] = Main end

########################################################################
# Create the HMMs

function getstate(fsm, sid, pdfindex)
    if sid == "start"
        return initstate(fsm)
    elseif sid == "end"
        return finalstate(fsm)
    else
        if sid in keys(fsm.states)
            return fsm.states[sid]
        else
            return addstate!(fsm, pdfindex = pdfindex)
        end
    end
end

@info "creating HMMs..."

conf = open(confhmm, "r") do f YAML.load(f) end

const hmms = Dict{String, FSM}()
pdfcount = 0
for unit in units
    label, type = split(unit)
    unitconf = conf[type]

    @debug "creating HMM for unit $unit with topology type $type"

    fsm = FSM()
    for arc in unitconf["topology"]
        src = getstate(fsm, arc["start_id"], pdfcount + 1)
        if isemitting(src) && src.pdfindex > pdfcount && ! unitconf["sharedpdf"]
            global pdfcount += 1
        end

        dest = getstate(fsm, arc["end_id"], pdfcount + 1)
        if isemitting(dest) && dest.pdfindex > pdfcount && ! unitconf["sharedpdf"]
            global pdfcount += 1
        end

        link!(fsm, src, dest, log(arc["trans_prob"]))
    end
    if unitconf["sharedpdf"] global pdfcount += 1 end

    hmms[label] = fsm |> weightnormalize!
end
@debug "number of created pdf: $pdfcount"

bson(joinpath(outdir, "hmms.bson"), Dict(:hmms => hmms))

########################################################################
# Create the emissions

function getstats(rng, uttids, n, D)
    x̂ = zeros(D)
    x̂² = zeros(D)
    N = 0.
    for i in 1:n
        path = joinpath(feadir, uttids[rand(rng, 1:length(uttids))] * ".bson")
        data = BSON.load(path)[:data]
        N += size(data, 2)
        x̂ .+= dropdims(sum(data, dims = 2), dims = 2)
        x̂² .+= dropdims(sum(data.^2, dims = 2), dims = 2)
    end
    μ̂ = x̂ ./ N
    μ̂, (x̂² ./ N) - μ̂.^2
end

const S = pdfcount
const K = args["filter-order"]
const strength = args["prior-strength"]

# Check the dimension of the features
data = BSON.load(joinpath(feadir, uttids[1] * ".bson"))[:data]
const D = size(data, 1)

seed = args["seed"]
@debug "seeding the random generator with $seed"
const rng = Random.seed!(seed)

@info "estimate the statistics from the training data..."
μ̂, σ̂² = getstats(rng, uttids, args["n-utts"], D)

@info "creating the emissions"
μ₀s = [begin
    μ₀ = zeros(K+1)
    μ₀[end] = μ̂[d]
    μ₀
end for d in 1:D]
Σ₀s = [begin
    Σ₀ = Matrix{eltype(μ₀s[1])}(I, K+1, K+1) ./ strength
    Σ₀[end, end] = σ̂²[d]
    Σ₀
end for d in 1:D]
a₀ = strength
b₀ = strength
emissions = [[ARNormal1D(K, μ₀s[d], Σ₀s[d], a₀, b₀) for i in 1:S] for d in 1:D];

save(joinpath(outdir, "emissions_0.bson"), emissions)

max_stars_repo_name: BUTSpeechFIT/ConvHMM

max_issues_repo_name: BUTSpeechFIT/ConvHMM

max_forks_repo_name: BUTSpeechFIT/ConvHMM


text:
module Fread

using RCall, Feather

export fread

const PKG_NOT_INSTALLED_ERR_MSG = "{data.table} and/or {arrow} package(s) not installed\n run 'Fread.install_pkgs()' or install {data.table} and {arrow} manually in R using 'install.packages(c('data.table', 'arrow'))'"

"""
    Fread.install_pkgs()

Install the R packages needed namely data.table and arrow
"""
function install_pkgs(repos = "https://cran.rstudio.com")
    run(`R -e "install.packages(c('data.table', 'arrow'),repos='$repos')"`)
    pkgs_installed()
end

"""
    Fread.pkgs_installed()

Checks if all required R packages are installed
"""
function pkgs_installed()
    R"""
        ml = memory.limit()
        # TODO better
        memory.limit(max(ml, ml*2, 2e14))
        datatable_installed <- require(data.table)
        arrow_installed <- require(arrow)
    """
    @rget datatable_installed
    @rget arrow_installed
    (datatable_installed = datatable_installed, arrow_installed = arrow_installed)
end

"""
    fread(path, feather_out_path = path*".feather"; params...)

Reads a CSV from `path` and create a DataFrame. An intermediate feather file
is created at `feather_out_path`

* `path` input CSV path
* `feather_out_path` intermediate feather file path
* `params` parameters passed to data.table
"""
fread(path,  feather_out_path = path*".feather"; params...) = begin
    csv_to_feather(path, feather_out_path; params...)

    Feather.read(feather_out_path)
end

"""
    csv_to_feather(path, outpath; params...)

Converts a CSV file to feather format
"""
csv_to_feather(path, outpath; params...) = begin
    if !all(pkgs_installed())
        throw(ErrorException(PKG_NOT_INSTALLED_ERR_MSG))
    end

    @rput params
    R"""
        ml = memory.limit()
        # TODO better
        memory.limit(max(ml, ml*2, 2e14))
        if(!require(data.table)) {
            stop("data.table not installed")
        }
        if(!require(arrow)) {
            stop("arrow R package not installed")
        }
        #res = data.table::fread($path)
        res = do.call(data.table::fread, c(list($path), params))
        arrow::write_feather(res, $outpath)
    """
end

"""
    csv_to_parquet(path, outpath; params...)

Converts a CSV file to parquet format
"""
csv_to_parquet(path, outpath; params...) = begin
    if !all(pkgs_installed())
        throw(ErrorException(PKG_NOT_INSTALLED_ERR_MSG))
    end

    @rput params
    R"""
        ml = memory.limit()
        # TODO better
        memory.limit(max(ml, ml*2, 2e14))
        if(!require(data.table)) {
            stop("data.table not installed")
        }
        if(!require(arrow)) {
            stop("arrow R package not installed")
        }
        #res = data.table::fread($path)
        res = do.call(data.table::fread, c(list($path), params))
        arrow::write_parquet(res, $outpath)
    """
end

end # module


text:
@testset "Exceptions" begin

function oob_kernel(X)
    X[0] = 1
    nothing
end

RA = ROCArray(ones(Float32, 4))
try
    wait(@roc oob_kernel(RA))
catch err
    @test err isa AMDGPU.KernelException
    if err isa AMDGPU.KernelException
        @test err.exstr !== nothing
        @test occursin("BoundsError", err.exstr)
    end
end

end


text:
using CGMethod1D
using Ferrite

grid = generate_grid(Line, (10,))

dim = 1
order = 3
ip = Lagrange{dim, RefLine, order}()
qr = QuadratureRule{dim, RefLine}(order)
# cellvalues = CellScalarValues(qr, ip)

# dh = DofHandler(grid)
# push!(dh, :p, 1, ip)
# # push!(dh, :p, 1)
# close!(dh)

# K = create_sparsity_pattern(dh);

# using UnicodePlots
# fill!(K.nzval, 1.0)
# spy(K; height = 15)

# ch = ConstraintHandler(dh);

# ∂Ω = union(getfaceset.((grid, ), ["left", "right", "top", "bottom"])...);

# dbc = Dirichlet(:u, ∂Ω, (x, t) -> 0)
# add!(ch, dbc);

# close!(ch)
# update!(ch, 0.0);

# function doassemble(cellvalues::CellScalarValues{dim}, K::SparseMatrixCSC, dh::DofHandler) where {dim}

#     n_basefuncs = getnbasefunctions(cellvalues)
#     Ke = zeros(n_basefuncs, n_basefuncs)
#     fe = zeros(n_basefuncs)

#     f = zeros(ndofs(dh))
#     assembler = start_assemble(K, f)

#     @inbounds for cell in CellIterator(dh)

#         fill!(Ke, 0)
#         fill!(fe, 0)

#         reinit!(cellvalues, cell)

#         for q_point in 1:getnquadpoints(cellvalues)
#             dΩ = getdetJdV(cellvalues, q_point)

#             for i in 1:n_basefuncs
#                 v  = shape_value(cellvalues, q_point, i)
#                 ∇v = shape_gradient(cellvalues, q_point, i)
#                 fe[i] += v * dΩ
#                 for j in 1:n_basefuncs
#                     ∇u = shape_gradient(cellvalues, q_point, j)
#                     Ke[i, j] += (∇v ⋅ ∇u) * dΩ
#                 end
#             end
#         end

#         assemble!(assembler, celldofs(cell), fe, Ke)
#     end
#     return K, f
# end

# K, f = doassemble(cellvalues, K, dh);

# apply!(K, f, ch)
# u = K \ f;


text:
using DifferentiableNAS
using Flux
using CUDA
using Distributions: Bernoulli

@testset "DARTS MixedOp" begin
    steps = 4
    k = floor(Int, steps^2/2+3*steps/2)
    num_ops = length(PRIMITIVES)
    mo = MixedOp(1,"1-2",4,1)  |> gpu
    @test length(Flux.params(mo).order |> cpu) > 0
    data = rand(Float32,8,8,4,2)  |> gpu
    α = rand(Float32, num_ops)
    @test size(data) == size(mo(data, α))
    g = gradient((x,α) -> sum(mo(x, α)), data, α)
    @test size(data) == size(g[1])
    @test size(α) == size(g[2])
end

@testset "DARTS Cell" begin
    steps = 4
    k = floor(Int, steps^2/2+3*steps/2)
    num_ops = length(PRIMITIVES)
    data = rand(Float32,8,8,4,2) |> gpu
    cell = Cell(4, 4, 1, false, false, 4, 4, 1) |> gpu
    @test length(Flux.params(cell).order) > 0
    as = [2e-3*(rand(num_ops).-0.5) |> f32   for _ in 1:k]
    @test size(data) == size(cell(data, data, as))
    grad = gradient((x1, x2, αs) -> sum(cell(x1, x2, αs)), data, data, as)
    @test size(data) == size(grad[1])
    @test size(as) == size(grad[3])
end

@testset "DARTS Model" begin
    steps = 4
    k = floor(Int, steps^2/2+3*steps/2)
    num_ops = length(PRIMITIVES)
    data = rand(Float32,8,8,4,2) |> gpu
    masked_αs = [2e-3*(rand(num_ops).-0.5).*rand(Bernoulli(),num_ops) |> f32 for _ in 1:k]
    m = DARTSModel(num_cells=4, channels=4, track_acts = true)
    @test length(Flux.params(m).order) > 1
    @test length(Flux.params(m.cells).order) > 1
    test_image = rand(Float32, 32, 32, 3, 1) |> gpu
    grad = gradient(x->sum(m(x)), test_image)
    @test size(test_image) ==  size(grad[1])
    loss(m, x) = sum(m(x))
    gws = gradient(Flux.params(m.cells)) do
        sum(m(test_image))
    end
    @test typeof(gws[Flux.params(m.cells)[1]]) != Nothing
    a = all_αs(m)
    gαs = gradient(a) do
        sum(m(test_image))
    end
    @test typeof(gαs[a[1]]) != Nothing
    @test length(m.activations.currentacts) > 0
    key = collect(keys(m.activations.currentacts))[1]
    sample1 = m.activations.currentacts[key]
    out = m(test_image, [masked_αs, masked_αs])
    sample2 = m.activations.currentacts[key]
    @test sample1 != sample2
end
@testset "DARTS Model BN" begin
    steps = 4
    k = floor(Int, steps^2/2+3*steps/2)
    num_ops = length(PRIMITIVES)
    data = rand(Float32,8,8,4,2) |> gpu
    m = DARTSModelBN(num_cells=4, channels=4)
    @test length(Flux.params(m).order) > 1
    @test length(Flux.params(m.cells).order) > 1
    test_image = rand(Float32, 32, 32, 3, 1) |> gpu
    grad = gradient(x->sum(m(x)), test_image)
    @test size(test_image) ==  size(grad[1])
    gws = gradient(Flux.params(m.cells)) do
        sum(m(test_image))
    end
    @test typeof(gws[Flux.params(m.cells)[1]]) != Nothing
    a = all_αs(m)
    gαs = gradient(a) do
        sum(m(test_image))
    end
    @test typeof(gαs[a[1]]) != Nothing
end

@testset "DARTS Eval Cell" begin
    steps = 4
    k = floor(Int, steps^2/2+3*steps/2)
    num_ops = length(PRIMITIVES)
    data = rand(Float32,8,8,4,2) |> gpu
    αs = [2e-3*(rand(num_ops).-0.5) |> f32  for _ in 1:k]
    cell = EvalCell(4, 4, 1, false, false, 4, 4, αs, PRIMITIVES)
    @test length(Flux.params(cell).order) > 0
    @test size(data) == size(cell(data, data, 0.4f0))
    grad = gradient((x1, x2) -> sum(cell(x1, x2, 0.4f0)), data, data)
    @test size(data) == size(grad[1])
end

@testset "DARTS Eval Model" begin
    steps = 4
    k = floor(Int, steps^2/2+3*steps/2)
    num_ops = length(PRIMITIVES)
    normal = [2e-3*(rand(num_ops).-0.5) |> f32 for _ in 1:k]
    reduce = [2e-3*(rand(num_ops).-0.5) |> f32 for _ in 1:k]
    m = DARTSEvalModel(normal, reduce)
    @test length(Flux.params(m).order) > 1
    @test length(Flux.params(m.cells).order) > 1
    test_image = rand(Float32, 32, 32, 3, 1) |> gpu
    grad = gradient(x->sum(m(x)), test_image)
    @test size(test_image) ==  size(grad[1])
    loss(m, x) = sum(m(x, 0.4f0))
    gws = gradient(Flux.params(m.cells)) do
        sum(m(test_image))
    end
    @test typeof(gws[Flux.params(m.cells)[1]]) != Nothing
end

@testset "DARTS Eval Aux Model" begin
    steps = 4
    k = floor(Int, steps^2/2+3*steps/2)
    num_ops = length(PRIMITIVES)
    normal = [2e-3*(rand(num_ops).-0.5) |> f32  for _ in 1:k]
    reduce = [2e-3*(rand(num_ops).-0.5) |> f32  for _ in 1:k]
    m = DARTSEvalAuxModel(normal, reduce)
    @test length(Flux.params(m).order) > 1
    @test length(Flux.params(m.cells).order) > 1
    test_image = rand(Float32, 32, 32, 3, 1) |> gpu
    out, out_aux = m(test_image, true)
    @test size(out,1) == size(out_aux,1)
    grad = gradient(x->sum(m(x, true)[1]), test_image)
    @test size(test_image) ==  size(grad[1])
    loss(m, x) = sum(m(x, true, 0.4f0))
    gws = gradient(Flux.params(m.cells)) do
        sum(sum(m(test_image, true)))
    end
    @test typeof(gws[Flux.params(m.cells)[1]]) != Nothing
end

max_stars_repo_name: huangyxi/DifferentiableNAS

max_issues_repo_name: huangyxi/DifferentiableNAS

max_forks_repo_name: huangyxi/DifferentiableNAS


text:
function Base.print(io::IO, f::ReducerFunctionAndFoldable)
    show(IOContext(io, :limit => false), f)
end

function Base.show(io::IO, ::MIME"text/plain", f::ReducerFunctionAndFoldable)
    n = length(methods(f))
    m = n > 1 ? "methods" : "method"
    print(io, f)
    print(io, " (generic function with $n $m)")
end


text:
const DISABLE_FM_PARALLEL = true
FM_M = Matrix{Float64}
FM_V = Vector{Float64}
FM_QR = Matrix{Float64}

function maketable1(panel::DataFrame; trialrunonly::Bool=trialrunonly,
    runpanela::Bool=true,
    runpanelb::Bool=true,
    runpanelc::Bool=true,
    runpaneld::Bool=true,
    )::Nothing
  #load the data we will work with

  sort!(panel, [:permno, :date])
  apanel::DataFrame = panel[panel.annflag,:]

  ################common data
  cqports::PortfolioConstructions = leverageconstructions()
  T1_FNAMES::Vector{Symbol} = CQUART_FNAMES
  T1_FVALS::Vector{Symbol} = CQUART_FVALS
  T1_FGROUPS::Vector{Symbol} = CQUART_FGROUPS
  T1_FNAMES_ALT::Vector{Symbol} = CQUART_FNAMES_ALT
  T1_FGROUPS_ALT::Vector{Symbol} = CQUART_FGROUPS_ALT

  ################Table A
  T1_PANEL_A_ROWKEYS::Vector{Symbol} = [:mean, :std, :stdxt]

  T1_PANEL_A_ROWLABEL::Dict = Dict(
    :mean=>"Mean",
    :std=>"SD",
    :stdxt=>"SD(X,T)")

  T1_PANEL_A_IDX::Dict = Dict(
    :mean=>1,
    :std=>2,
    :stdxt=>3)
  T1_PANEL_A_NAME::String = "t1panela-$(REPLICATION_TYPE[])-$(OUT_SUFFIX[])"

  runpanela && tablexpanelA(apanel,
    fvals = T1_FVALS,
    rowkeys=T1_PANEL_A_ROWKEYS,
    rowlabel=T1_PANEL_A_ROWLABEL,
    idx=T1_PANEL_A_IDX,
    tablename=T1_PANEL_A_NAME)

  #####################TABLE B
  T1_PANEL_B_ROWKEYS::Vector{Symbol} = [
    :leverageQmeans,
    :quartilenum,
    :LLsd,
    :sd,
    :Δsd,
    :LLret12m,
    :ret12m,
    :Δret12m]

  T1_PANEL_B_ROWLABEL::Dict = Dict(
    :leverageQmeans=>"Leverage Q Means",
    :quartilenum=>"Quartile \\#",
    :LLsd=>"\\\\ SD Net Lagged\$^2\$",
    :sd=>"SD Net Lead",
    :Δsd=>"SD Net 2Y Delta",
    :LLret12m=>"\\\\ Compound Raw Lagged\$^2\$",
    :ret12m=>"Compound Raw Lead",
    :Δret12m=>"Compound Raw 2Y Delta")

  T1_PANEL_B_IDX::Dict = Dict(
    T1_PANEL_B_ROWKEYS[i]=>i for i ∈ 1:(length(T1_PANEL_B_ROWKEYS)))

  #the below is to help with future abstraction if necessary
  T1_PANEL_B_VALCOLIDX::Dict= Dict(f=>f for f ∈ T1_FVALS)

  #these are pairs of lagged and unlagged fields
  T1_PANEL_B_VOLFIELDS::NTuple{2,Symbol} = (:LLvol252d, :vol252d)
  T1_PANEL_B_RETFIELDS::NTuple{2,Symbol} = (:LLret12m, :ret12m)

  T1_PANEL_B_NAME::String = "t1panelb-$(REPLICATION_TYPE[])-$(OUT_SUFFIX[])"

  runpanelb && tablexpanelB(apanel,
    fvals = T1_FVALS,
    rowkeys=T1_PANEL_B_ROWKEYS,
    rowlabel=T1_PANEL_B_ROWLABEL,
    idx=T1_PANEL_B_IDX,
    valcolidx = T1_PANEL_B_VALCOLIDX,
    tablename=T1_PANEL_B_NAME,
    volfields=T1_PANEL_B_VOLFIELDS,
    retfields=T1_PANEL_B_RETFIELDS)

  ###################PANEL C############

  T1_PANEL_C_ROWKEYS::Vector{Symbol} = [:constant, :plusxmkt, :plussmbhml, :plusrmwcma, :plusumd]

  T1_PANEL_C_ROWLABEL::Dict = Dict(
    :constant=>"Constant",
    :plusxmkt=>"+XMKT",
    :plussmbhml=>"+SMB+HML",
    :plusrmwcma=>"+RMW+CMA",
    :plusumd=>"+UMD")

  T1_PANEL_C_CONTROLS::Dict = Dict(
    :constant=>CONTROL0,
    :plusxmkt=>[:mkt],
    :plussmbhml=>[:mkt, :P_ff3m_smb,	:P_ff3m_hml],
    :plusrmwcma=>[:mkt, :P_ff5_smb,	:P_ff5_hml,	:P_ff5_rmw,	:P_ff5_cma],
    :plusumd=>[:mkt, :P_ff5_smb,	:P_ff5_hml,	:P_ff5_rmw,	:P_ff5_cma, :P_ff3m_umd])

  T1_PANEL_C_FOCALS::Dict = Dict(
    :constant=>T1_FNAMES,
    :plusxmkt=>T1_FNAMES,
    :plussmbhml=>T1_FNAMES,
    :plusrmwcma=>T1_FNAMES,
    :plusumd=>T1_FNAMES)

  T1_PANEL_C_NAME::String = "t1panelc-$(REPLICATION_TYPE[])-$(OUT_SUFFIX[])"

  T1_PANEL_C_MIN_ROWS::Int = 24

  runpanelc && tablexpanelC(panel,
    fvals=T1_FVALS,
    fgroups=T1_FGROUPS,
    Fnames = T1_FNAMES,
    rowkeys=T1_PANEL_C_ROWKEYS,
    rowlabel=T1_PANEL_C_ROWLABEL,
    controls=T1_PANEL_C_CONTROLS,
    focals=T1_PANEL_C_FOCALS,
    tablename=T1_PANEL_C_NAME,
    minrowspergroup = T1_PANEL_C_MIN_ROWS
  )

  #########################PANEL D#########################
  T1_PANEL_D_CONTROLS::Dict = Dict(
    :quartiles=>CONTROL0,
    :quartiles4ctrl=>T1_CONTROL4,
    :quartiles10ctrl=>T1_CONTROL10,
    :values=>CONTROL0,
    :values4ctrl=>T1_CONTROL4,
    :values10ctrl=>T1_CONTROL10,
    :quartiles_ac=>CONTROL0,
    :quartiles_lret=>T1_CONTROL10)

  T1_PANEL_D_RET::Dict = Dict(
    :quartiles=>:ret,
    :quartiles4ctrl=>:ret,
    :quartiles10ctrl=>:ret,
    :values=>:ret,
    :values4ctrl=>:ret,
    :values10ctrl=>:ret,
    :quartiles_ac=>:ret,
    :quartiles_lret=>:lret)

  T1_PANEL_D_ROWKEYS::Vector{Symbol} = [:quartiles, :quartiles4ctrl, :quartiles10ctrl, :values,
    :values4ctrl, :values10ctrl, :quartiles_ac, :quartiles_lret]

  T1_PANEL_D_ROWLABEL::Dict = Dict(
    :quartiles=>"Leverage Quartiles",
    :quartiles4ctrl=>"... +4 vars",
    :quartiles10ctrl=>"... +10 vars",
    :values=>"Leverage Values",
    :values4ctrl=>"... +4 vars",
    :values10ctrl=>"... +10 vars",
    :quartiles_ac=>"\\midrule Asset-Controlled Quartiles",
    :quartiles_lret=>"Log Returns +10 vars")

  T1_PANEL_D_FOCALS::Dict = Dict(
    :quartiles=>T1_FGROUPS,
    :quartiles4ctrl=>T1_FGROUPS,
    :quartiles10ctrl=>T1_FGROUPS,
    :values=>T1_FVALS,
    :values4ctrl=>T1_FVALS,
    :values10ctrl=>T1_FVALS,
    :quartiles_ac=>T1_FGROUPS_ALT,
    :quartiles_lret=>T1_FGROUPS)

  T1_PANEL_D_NAME::String = "t1paneld-$(REPLICATION_TYPE[])-$(OUT_SUFFIX[])"
  T1_PANEL_D_MIN_ROWS::Int = 24

  runpaneld && tablexpanelD(panel,
      fvals=T1_FVALS,
      fgroups=T1_FGROUPS,
      ret = T1_PANEL_D_RET,
      rowkeys=T1_PANEL_D_ROWKEYS,
      rowlabel=T1_PANEL_D_ROWLABEL,
      controls=T1_PANEL_D_CONTROLS,
      focals=T1_PANEL_D_FOCALS,
      tablename=T1_PANEL_D_NAME,
      minrowspergroup = T1_PANEL_D_MIN_ROWS
    )

  return nothing
end

function tablexpanelA(panel::DataFrame;
    fvals::Vector{Symbol} = error("PANELA: fvals is required"),
    rowkeys::Vector{Symbol} = error("PANELA: rowkeys is required"),
    rowlabel::Dict = error("PANELA: rowlabel is required"),
    idx::Dict = error("PANELA: idx is required"),
    tablename::String = error("PANELA: table name is required"),
    tablepath::String = TABLE_PATH,
    decimals::Int = DEFAULT_DECIMALS)
  #local cqports::CQuartileConstructions = leverageconstructions()

  #NOTE: These are of the form mknegcash ~ α + β*mkt ...

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

  #make the table
  colnames::Vector{Vector{String}} = [(string).(fvals)]
  descrownames::Vector{String} = (s->rowlabel[s]).(rowkeys)

  desccontent::Vector{Vector{String}} = (i->Vector{String}(undef, length(fvals))).(1:length(rowkeys))
  for c ∈ 1:length(fvals)
    local sspanels::GroupedDataFrame
    local ngroups::Int

    f::Symbol = fvals[c]
    allfields::Vector{Symbol} = [:date, :permno, f]
    spanel::SubDataFrame = view(panel, completecases(panel[!,allfields]), allfields)

    cmean::Float64 = mean(spanel[!, f])
    cstd::Float64 = std(spanel[!, f])

    cstdx::Float64 = 0.0
    ngroups = 0
    sspanels = groupby(spanel, :date)
    for sspanel ∈ sspanels
      cstdxi::Float64 = std(sspanel[!,f])
      if isfinite(cstdxi)
        cstdx += cstdxi
        ngroups+=1
      end
    end
    cstdx /= ngroups

    cstdt::Float64 = 0.0
    ngroups = 0
    sspanels = groupby(spanel, :permno)
    for sspanel ∈ groupby(spanel, :permno)
      cstdti::Float64 = std(sspanel[!,f])
      if isfinite(cstdti)
        cstdt += cstdti
        ngroups+=1
      end
    end
    cstdt /= ngroups

    desccontent[idx[:mean]][c] = brnum(cmean, decimals)
    desccontent[idx[:std]][c] = "($(num2str(cstd, decimals)))"
    desccontent[idx[:stdxt]][c] = "($(num2str(cstdx, decimals)), $(num2str(cstdt, decimals)))"
  end

  tt::String = textable(colnames=colnames, descrownames=descrownames, desccontent=desccontent)

  writetextable(tt, path=tablepath, outname="$tablename.tex")
end

function tablexpanelB(panel::DataFrame;
    fvals::Vector{Symbol} = error("PANELB: fvals is required"),
    #fgroups::Vector{Symbol} = error("PANELB: fgroups is required"),
    rowkeys::Vector{Symbol} = error("PANELB: rowkeys is required"),
    rowlabel::Dict = error("PANELB: rowlabel is required"),
    idx::Dict = error("PANELB: idx is required"),
    valcolidx::Dict = error("PANELB: valcolidx is required"),
    tablename::String = error("PANELB: table name is required"),
    volfields::NTuple{2, Symbol}=error("PANELB: volfields is required"),
    retfields::NTuple{2, Symbol}=error("PANELB: retfields is required"),
    tablepath::String = TABLE_PATH,
    decimals::Int = 1)
  local cqports::CQuartileConstructions = leverageconstructions()

  #(s->winsorize!(panel[!,s], prop = .999, sorted=false)).([retfields...; volfields...;])
  #println(describe(panel))

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

  #will be using this a lot
  @inline skipmean(v::AbstractVector{<:MFloat64})::Float64 = mean(skipmissing(v))
  @inline n2s(x::T where T<:MFloat64, d::Real=decimals) = num2str(x, d, scalefactor=100.)
  @inline br(x::T where T<:MFloat64, d::Real=decimals) = brnum(x, d, scalefactor=100.)

  #=println("Test neg cash: (")
  sdf::SubDataFrame = view(panel, (g->(!ismissing(g)) && g==1).(panel[!, cqports[1].Fgroup]),:)
  print("$(skipmean(sdf[!, T1_PANEL_B_RETFIELDS[1]])), ")
  sdf = view(panel, (g->(!ismissing(g)) && g==4).(panel[!, cqports[1].Fgroup]),:)
  print("$(skipmean(sdf[!, T1_PANEL_B_RETFIELDS[1]])))")=#

  #make the table
  colnames::Vector{Vector{String}} = [[(s->string(s)).(fvals)...;]]
  widthcolnames::Vector{Vector{Int}} = [[(s->2).(fvals)...;]]
  descrownames::Vector{String} = (s->rowlabel[s]).(rowkeys)
  alignmentstring::String = string(" l", join(["r@{ }r" for i ∈ 1:length(fvals)]))

  #pre-allocate for content
  numcols::Int = length(fvals)*2
  desccontent::Vector{Vector{String}} = (
    i->Vector{String}(undef, numcols)).(1:length(rowkeys))

  volscale::Float64 = REPLICATION_TYPE[] == :daily ? 252^0.5 : 12^0.5

  for i ∈ 1:length(fvals)
    c::Int = 1+(i-1)*2

    local sspanelq1::SubDataFrame
    local sspanelq4::SubDataFrame
    local allfields::Vector{Symbol}

    #acquire the fields we will use
    f::Symbol = fvals[i]
    name::Symbol = Symbol("P_", f)
    Fgroup::Symbol = cqports[name].Fgroup
    #println("f: $f, name: $name, Fgroup: $Fgroup")
    basicfields::Vector{Symbol} = [:date; :permno; Fgroup]
    local spanel::SubDataFrame = view(panel, completecases(panel[!,basicfields]),:)

    ###Value Rows
    Fval::Symbol = valcolidx[f]
    #allfields = [basicfields; Fval]
    sspanelq1 = view(spanel, spanel[!,Fgroup] .== 1, :)
    sspanelq4 = view(spanel, spanel[!,Fgroup] .== 4, :)

    q1val::Float64 = mean(skipmissing(sspanelq1[!, Fval]))
    q4val::Float64 = mean(skipmissing(sspanelq4[!, Fval]))
    desccontent[idx[:leverageQmeans]][c] = (
      "$(brnum(q1val, 2))")
    desccontent[idx[:leverageQmeans]][c+1] = ("$(brnum(q4val, 2))")

    desccontent[idx[:quartilenum]][c] = ("\\text{Q1}")
    desccontent[idx[:quartilenum]][c+1] = ("\\text{Q4}")

    ###SD Rows
    #allfields = [basicfields; T1_PANEL_B_VOLFIELDS]
    sspanelq1 = view(spanel, spanel[!,Fgroup] .== 1, :)
    sspanelq4 = view(spanel, spanel[!,Fgroup] .== 4, :)

    q1LLvol::Float64 = skipmean(sspanelq1[!, volfields[1]]) * volscale
    q4LLvol::Float64 = skipmean(sspanelq4[!, volfields[1]]) * volscale
    desccontent[idx[:LLsd]][c] = ("$(n2s(q1LLvol))")
    desccontent[idx[:LLsd]][c+1] = ("$(n2s(q4LLvol))")

    q1vol::Float64 = skipmean(sspanelq1[!, volfields[2]]) * volscale
    q4vol::Float64 = skipmean(sspanelq4[!, volfields[2]]) * volscale
    desccontent[idx[:sd]][c] = ("$(n2s(q1vol))")
    desccontent[idx[:sd]][c+1] = ("$(n2s(q4vol))")

    Δq1vol::Float64 = skipmean(sspanelq1[!, volfields[2]] .-
      sspanelq1[!, volfields[1]]) * volscale
    Δq4vol::Float64 = skipmean(sspanelq4[!, volfields[2]] .-
      sspanelq4[!, volfields[1]]) * volscale
    desccontent[idx[:Δsd]][c] = ("$(br(Δq1vol))")
    desccontent[idx[:Δsd]][c+1] = ("$(br(Δq4vol))")

    ###Raw Ret Rows
    #allfields = [basicfields; T1_PANEL_B_RETFIELDS]
    sspanelq1 = view(spanel, spanel[!,Fgroup] .== 1, :)
    sspanelq4 = view(spanel, spanel[!,Fgroup] .== 4, :)

    q1LLret::Float64 = skipmean(sspanelq1[!, retfields[1]])
    q4LLret::Float64 = skipmean(sspanelq4[!, retfields[1]])
    desccontent[idx[:LLret12m]][c] = ("$(n2s(q1LLret))")
    desccontent[idx[:LLret12m]][c+1] = ("$(n2s(q4LLret))")

    q1ret::Float64 = skipmean(sspanelq1[!, retfields[2]])
    q4ret::Float64 = skipmean(sspanelq4[!, retfields[2]])
    desccontent[idx[:ret12m]][c] = ("$(n2s(q1ret))")
    desccontent[idx[:ret12m]][c+1] = ("$(n2s(q4ret))")

    Δq1ret::Float64 = skipmean(sspanelq1[!, retfields[2]] .-
      sspanelq1[!, retfields[1]])
    Δq4ret::Float64 = skipmean(sspanelq4[!, retfields[2]] .-
      sspanelq4[!, retfields[1]])
    desccontent[idx[:Δret12m]][c] = ("$(br(Δq1ret))")
    desccontent[idx[:Δret12m]][c+1] = ("$(br(Δq4ret))")

  end

  tt::String = textable(colnames=colnames,
    descrownames=descrownames,
    desccontent=desccontent,
    widthcolnames=widthcolnames,
    alignmentstring=alignmentstring)

  writetextable(tt, path=tablepath, outname="$tablename.tex")
end

function tablexpanelC(panel::DataFrame;
    fvals::Vector{Symbol} = error("PANELC: fvals is required"),
    Fnames::Vector{Symbol} = error("PANELC: Fnames is required"),
    rowkeys::Vector{Symbol} = error("PANELC: rowkeys is required"),
    rowlabel::Dict = error("PANELC: rowlabel is required"),
    focals::Dict = error("PANELC: focals is required"),
    controls::Dict = error("PANELC: controls is required"),
    tablename::String = error("PANELC: table name is required"),
    fgroups::Vector{Symbol} = error("PANELC: fgroups is required"),
    tablepath::String = TABLE_PATH,
    decimals::Int = DEFAULT_DECIMALS,
    minrowspergroup::Int = error("minrowspergroup is required"))

  #local cqports::CQuartileConstructions = leverageconstructions()

  #println(Fnames)

  local parallel::Bool = PARALLEL[] && (!DISABLE_FM_PARALLEL)

  #################################
  fss::Vector{FactorSpecification} = Vector{FactorSpecification}()
  sizehint!(fss, length(rowkeys) * length(fgroups))

  #need to re-shape for this exercise
  idcols::Vector{Symbol} = unique([:date; values(controls)...;])
  measurecols::Vector{Symbol} = Fnames
  allcols::Vector{Symbol} = unique([:date; values(controls)...; values(focals)...;])
  minipanel::DataFrame = unique(panel[!,allcols])
  minipanel = melt(minipanel, idcols, measurecols, variable_name = :portname, value_name = :ret)
  sort!(minipanel, [:portname, :date])

  #println("allcols: $allcols")
  #println("minipanelun: ", unique(minipanel[!, :portname]))
  #println(describe(minipanel))
  #println("Size minipanel: $(size(minipanel))")

  #make the regression specifications
  for rname ∈ rowkeys
    sallcols::Vector{Symbol} = [:ret; :date; :portname; controls[rname]]
    sminipanel = view(minipanel, completecases(minipanel[!,sallcols]), sallcols)
    #println("sminipanelun: ", describe(sminipanel))
    fs = FactorSpecification(sminipanel, :ret, :date, :portname,
      Ffactors = nothing,
      Fcharacteristics=nothing,
      Fcontrols = controls[rname],
      name=Symbol("row$(rname)"),
      skipcrosssectiondfcheck=true,
      minrowspergroup=minrowspergroup)
    #lock(spinlock)
    push!(fss, fs) #only doing time series regressions
    #unlock(spinlock)
  end

  #run the regressions
  #WARNING: When BLAS paralllelization issues are solved replace with parallel=PARALLEL[]
  frs::Vector{FactorRegression} = Vector{FactorRegression}()

  #generate and push all the regressions
  #each call to manyFF pulls a row of regressions
  for i ∈ 1:length(fss)
    sfrs::Vector{FactorRegression} = manyFF(fss[i], FM_M, FM_V,
      qrtype=FM_QR, parallel=parallel)

    push!(frs, sfrs...)
  end

  #println((fr->fr.name).(frs))

  #make the table
  colnames::Vector{Vector{String}} = [
    (string).(fvals), (i->"Alpha (T-stat)").(1:length(fvals))]
  descrownames::Vector{String} = (s->rowlabel[s]).(rowkeys)

  desccontent::Vector{Vector{String}} = factordesccontent(
    fr::FactorRegression->coeftstat(fr, i=1), #fss=fss,
    frs=frs, rowkeys = rowkeys, colindex=focals)

  tt::String = textable(colnames=colnames, descrownames=descrownames, desccontent=desccontent)

  writetextable(tt, path=tablepath, outname="$tablename.tex")
end

function tablexpanelD(panel::DataFrame;
    fvals::Vector{Symbol} = error("fvals is required"),
    ret::Dict = error("ret is required"),
    rowkeys::Vector{Symbol} = error("rowkeys is required"),
    rowlabel::Dict = error("rowlabel is required"),
    focals::Dict = error("focals is required"),
    controls::Dict = error("controls is required"),
    tablename::String = error("table name is required"),
    fgroups::Vector{Symbol} = error("fgroups is required"),
    tablepath::String = TABLE_PATH,
    minrowspergroup::Int = error("minrowspergroup is required"))

  local parallel::Bool = PARALLEL[] && (!DISABLE_FM_PARALLEL)

  #################################
  fss::Vector{FactorSpecification} = Vector{FactorSpecification}()
  sizehint!(fss, length(rowkeys) * length(fgroups))

  #make the regression specifications
  for rname ∈ rowkeys
    for focal ∈ focals[rname]
      fs::FactorSpecification = FactorSpecification(panel, ret[rname], :date, :permno,
        Ffactors=nothing,
        Fcharacteristics = [focal],
        Fcontrols = controls[rname],
        minrowspergroup=minrowspergroup,
        skiptimeseriesdfcheck=true,
        name=Symbol("row$(rname)_col$(focal)"))

      #lock(spinlock)
      push!(fss, fs) #only doing crossectional regressions
      #unlock(spinlock)
    end
  end

  #Run the regressions
  #WARNING: When BLAS paralllelization issues are solved replace with parallel=PARALLEL[]

  frs::Vector{FactorRegression} = Vector{FactorRegression}(undef, length(fss))
  for i ∈ 1:length(fss)
    frs[i] = characteristicFM(fss[i], FM_M, FM_V,
      parallel=parallel, qrtype=FM_QR)
  end

  #make the table
  colnames::Vector{Vector{String}} = [
    (string).(fvals), (i->"Coef (T-stat)").(1:length(fvals))]
  descrownames::Vector{String} = (s->rowlabel[s]).(rowkeys)
  desccontent::Vector{Vector{String}} = factordesccontent(#fss=fss,
    frs=frs,  rowkeys = rowkeys, colindex=focals)

  tt::String = textable(colnames=colnames,
    descrownames=descrownames, desccontent=desccontent)

  writetextable(tt, path=tablepath, outname="$tablename.tex")
end

function coeftstat(fr::FactorRegression; decimals::Int=DEFAULT_DECIMALS, i::Int=2,
  )::String
  λᵢ::Float64 = fr.λᵢ[i] * 100
  tstat::Float64 = fr.λᵢ[i]/fr.σᵢ[i]

  outstring = "$(num2str(λᵢ, decimals))\\;($(num2str(tstat, decimals)))"

  return brnum(λᵢ, numstring=outstring)
end

#this is a kernel function to generate the content matrix under a number of different scenarios
function factordesccontent(f::Function = coeftstat;
    #fss::Vector{FactorSpecification} = error("fss is required for factordesccontent"),
    frs::Vector{FactorRegression} = error("frs is required for factordesccontent"),
    rowkeys::Vector{Symbol} = error("Rowkeys is required"),
    colindex::Dict = error("colindex is required for factordesccontent"))

    Nrows::Int = length(rowkeys)
    Ncols::Int = length(collect(values(colindex))[1])

    #pre-allocate
    desccontent::Vector{Vector{String}} = (i->Vector{String}(undef, Ncols)).(1:Nrows)

    #convenience dictionaries
    frsidx = Dict(fr.name=>fr for fr ∈ frs)

    rowcolnameidx = Dict(
      (rname::Symbol, focal::Symbol)=>Symbol("row$(rname)_col$(focal)")
      for rname::Symbol ∈ rowkeys for focal::Symbol ∈ colindex[rname])

    #write the content matrix
    for (r,rname) ∈ enumerate(rowkeys)
      for (c, focal) ∈ enumerate(colindex[rname])
        fr::FactorRegression = frsidx[rowcolnameidx[(rname, focal)]]
        desccontent[r][c] = f(fr)
      end
    end

    return desccontent
  end

#creates the blue/red numbers used in much fo the tables
@inline function brnum(f::Real,
    decimals::Int = DEFAULT_DECIMALS;
    Ints::Bool = false,
    scalefactor=1.,
    numstring::MString = num2str(f, decimals, Ints=Ints, scalefactor=scalefactor))

  local outstring::String

  if f≥0
    outstring = "\\textcolor{blue}{$numstring}"
  else
    outstring = "\\textcolor{red}{$numstring}"
  end

  return outstring
end

brnum(::Missing, args...; keyargs...)::MString = missing

max_stars_repo_name: clintonTE/CCA

max_issues_repo_name: clintonTE/CCA

max_forks_repo_name: clintonTE/CCA


text:
module ConfusionMatrices

using StatsBase

export ConfusionMatrix, CombinedConfusionMatrix

include("ConfusionMatrix.jl")
include("CombinedConfusionMatrix.jl")

end # module


text:
module LogNumbers

import Base
import Random

export LogNumber, Log,
    floattype, logsumexp

include("types.jl")

# IO, literals
Base.show(io::IO, x::Union{LogFloat32,LogFloat64}) = print(io, "log", '\"', exp(logvalue(x)), '\"')
Base.show(io::IO, x::AbstractLogNumber{F}) where {F} = print(io, "Log(", F, ", ", exp(logvalue(x)), ")")

include("literal_macro.jl")

# Comparison, floating point stuff
Base.:(==)(x::AbstractLogNumber, y::AbstractLogNumber) = logvalue(x) == logvalue(y)
Base.isequal(x::AbstractLogNumber, y::AbstractLogNumber) = isequal(logvalue(x), logvalue(y))
Base.hash(x::AbstractLogNumber, h) = hash(logvalue(x), hash(typeof(x), h))

Base.isapprox(x::AbstractLogNumber, y::AbstractLogNumber; args...) =
    isapprox(logvalue(x), logvalue(y); args...)
Base.eps(::Type{<:AbstractLogNumber{F}}) where F = Base.eps(F) # is this the right thing?

Base.:<(x::AbstractLogNumber, y::AbstractLogNumber) = logvalue(x) < logvalue(y)
Base.:<=(x::AbstractLogNumber, y::AbstractLogNumber) = logvalue(x) <= logvalue(y)
Base.isless(x::AbstractLogNumber, y::AbstractLogNumber) = isless(logvalue(x), logvalue(y))

Base.isinf(x::AbstractLogNumber) = isinf(logvalue(x)) && logvalue(x) > 0
Base.isnan(x::AbstractLogNumber) = isnan(logvalue(x))

# Arithmetic etc.
# See https://en.wikipedia.org/wiki/Log_probability for the formulae and precautions

function Base.:+(x::AbstractLogNumber, y::AbstractLogNumber)
    y, x = minmax(x, y)
    iszero(x) && return x
    isinf(y) && return y
    LogNumber(logvalue(x) + log1p(exp(logvalue(y) - logvalue(x))))
end

function Base.:-(x::AbstractLogNumber, y::AbstractLogNumber)
    m = max(x, y)               # preserver order to automatically throw DomainError
    iszero(m) && return m
    LogNumber(logvalue(x) + log1p(-exp(logvalue(y) - logvalue(x))))
end

Base.:*(x::AbstractLogNumber, y::AbstractLogNumber) = LogNumber(logvalue(x) + logvalue(y))
Base.:/(x::AbstractLogNumber, y::AbstractLogNumber) = LogNumber(logvalue(x) - logvalue(y))
Base.log(x::AbstractLogNumber) = LogNumber(log(logvalue(x)))

# Summing values in log space
# See the following pages for explanation:
# http://www.nowozin.net/sebastian/blog/streaming-log-sum-exp-computation.html
# https://www.xarg.org/2016/06/the-log-sum-exp-trick-in-machine-learning/

logsumexp(xs) = logsumexp(xs, eltype(xs))

function logsumexp(xs, ::Type{<:AbstractLogNumber{F}}) where {F}
    α, r = mapfoldr(logvalue, expsum_update, xs, init = (-infty(F), zero(F)))
    LogNumber(log(r) + α)
end

function logsumexp(xs, ::Type{F}) where {F}
    α, r = mapfoldr(identity, expsum_update, xs, init = (-infty(F), zero(F)))
    log(r) + α
end

function expsum_update(x, acc)
    α, r = acc

    if x <= α
        return (α, r + exp(x - α))
    else
        return (x, r * exp(α - x) + one(r))
    end
end

# Random value generation
Random.rand(rng::Random.AbstractRNG, ::Random.SamplerType{L}) where {L<:PrimitiveLogNumber} =
    reinterpret(L, rand(rng, floattype(L)))
Random.rand(rng::Random.AbstractRNG, ::Random.SamplerType{WrappedLogNumber{F}}) where {F<:AbstractFloat} =
    WrappedLogNumber{F}(rand(rng, F))

end


text:
module TestingInteractive

# TODO: This is a workaround for CI with GTK since precompilation fails with headless systems
# Remove after https://github.com/JuliaGraphics/Gtk.jl/issues/346 is resolved
if isinteractive()
  @info "is interactive"
  using Gtk
else
  @info "is not interactive"
end

ask() = ask_dialog("Hallo")

end # module

TestingInteractive.ask()


text:
module JungleHelperRollingRock

using ..Ahorn, Maple

@mapdef Entity "JungleHelper/RollingRock" RollingRock(x::Integer, y::Integer, cracked::Bool=false,
    spriteXmlName::String="", debrisSpriteDirectory::String="JungleHelper/RollingRock", flag::String="", rollingSpeed::Number=100.0, fallingSpeed::Number=200.0, instantFalling::Bool=false)

const placements = Ahorn.PlacementDict(
    "Rolling Rock (Jungle Helper)" => Ahorn.EntityPlacement(
        RollingRock
    )
)

function Ahorn.selection(entity::RollingRock)
    x, y = Ahorn.position(entity)
    spriteName = get(entity.data, "cracked", false) ? "boulder_cracked" : "boulder"
    sprite = "JungleHelper/RollingRock/$(spriteName)"
    return Ahorn.getSpriteRectangle(sprite, x, y)
end

function Ahorn.render(ctx::Ahorn.Cairo.CairoContext, entity::RollingRock, room::Maple.Room)
    spriteName = get(entity.data, "cracked", false) ? "boulder_cracked" : "boulder"
    sprite = "JungleHelper/RollingRock/$(spriteName)"
    Ahorn.drawSprite(ctx, sprite, 0, 0)
end

end

max_stars_repo_name: Jackal-Celeste/JungleHelper

max_issues_repo_name: Jackal-Celeste/JungleHelper

max_forks_repo_name: Jackal-Celeste/JungleHelper


text:
@testset "Test SplineSparse" begin

    @testset "full" begin
        s = SplineSparse(1, 2, 3, 1:4, 1:2)
        want = [1 2 0; 0 3 4]
        @test full(s) == want

        s = SplineSparse(2, 2, 6, 1:8, [1, 5, 2, 4])
        want = [1 2 0 0 3 4
                0 5 6 7 8 0]

        @test full(s) == want
    end

    @testset "row_kron" begin
        s1 = SplineSparse(1,2, 3, rand(1:10, 12), rand(1:2, 6))
        s2 = SplineSparse(1,2, 4, rand(1:10, 12), rand(1:3, 6))
        want = row_kron(full(s1), full(s2))
        s12 = row_kron(s1, s2)

        @test full(s12) == want
        @test full(row_kron(s12, s1)) == row_kron(full(s12), full(s1))
        @test full(row_kron(s12, s12)) == row_kron(full(s12), full(s12))

        Base.zero(::Type{String}) = ""

        s1 = SplineSparse(1,2, 2, ["a", "b", "c", "d"], [1, 1])
        s2 = SplineSparse(1,3, 3, map(string, 1:6), [1, 1])

        want = ["a1" "a2" "a3" "b1" "b2" "b3"
                "c4" "c5" "c6" "d4" "d5" "d6"]
        have = row_kron(s1, s2)

        @test full(have) == want

    end

    @testset "getindex" begin
        s = SplineSparse(1, 2, 10, rand(12), rand(1:9, 6))
        full_s = full(s)

        for r in 1:size(s, 1)
            for c in 1:size(s, 1)
                @test s[r, c] == full_s[r, c]
            end
        end
    end

    @testset "*" begin
        s1 =  SplineSparse(1, 2, 3, 1:4, [1, 2])
        s2 =  SplineSparse(1, 3, 4, 1:6, [2, 1])
        s12 = SplineSparse(2, 3, 12, [1, 2, 3, 2, 4, 6, 12, 15, 18, 16, 20, 24],
                           [2, 6, 5, 9])

        x1 = rand(3)
        x2 = rand(4)
        x12 = rand(12)

        @test s1*x1 == full(s1)*x1
        @test s2*x2 == full(s2)*x2
        @test s12*x12 == full(s12)*x12
    end

end


text:
#!/usr/bin/env julia

# RainFARM
# Stochastic downscaling following
# D'Onofrio et al. 2014, J of Hydrometeorology 15 , 830-843 and
# Rebora et. al 2006, JHM 7, 724
# Includes orographic corrections

# Implementation in Julia language

# Copyright (c) 2016, Jost von Hardenberg - ISAC-CNR, Italy

using RainFARM
using ArgParse
using Printf

function parse_commandline()
    s = ArgParseSettings()

    @add_arg_table! s begin
        "--slope", "-s"
            help = "spatial spectral slope"
            arg_type = Float64
            default = 0.0
        "--nens", "-e"
            help = "number of ensemble members"
            arg_type = Int
            default = 1
        "--nf", "-n"
            help = "Subdivisions for downscaling"
            arg_type = Int
            default = 2
        "--region", "-r", "-R"
            help = "Indices for region to cutout imin/imax/jmin/jmax"
            arg_type = AbstractString
            default = "0/0/0/0"
        "--weights", "-w", "--weight"
            help = "Weights file"
            arg_type = AbstractString
            default = ""
        "--outfile", "-o", "--out"
            help = "Output filename radix"
            arg_type = AbstractString
            default = "rainfarm"
        "infile"
            help = "The input file to downscale"
            arg_type = AbstractString
            required = true
        "--varname", "-v"
            help = "Input variable name"
            arg_type = AbstractString
            default = ""
        "--global", "-g"
            action = :store_true
            help = "conserve precipitation over full domain"
        "--conv", "-c"
            action = :store_true
            help = "conserve precipitation using convolution"
        "--kmin", "-k"
            help = "Minimum wavenumber"
            arg_type = Int
            default = 1
    end

    s.description = "RainFARM downscaling: creates NENS realizations, downscaling INFILE, increasing spatial resolution by a factor NF. The slope is computed automatically unless specified. \ua0 Weights can be created with rfweights.jl"

    return parse_args(s)
end

args = parse_commandline()
nf = args["nf"]
filenc = args["infile"]
weightsnc = args["weights"]
nens = args["nens"]
varname = args["varname"]
fnbase = args["outfile"]
sx = args["slope"]
fglob = args["global"]
fsmooth = args["conv"]
region = args["region"]
kmin = args["kmin"]

imin = parse(Int, split(region, "/")[1])
imax = parse(Int, split(region, "/")[2])
jmin = parse(Int, split(region, "/")[3])
jmax = parse(Int, split(region, "/")[4])

println("Downscaling ", filenc)

(pr, lon_mat, lat_mat, varname)=read_netcdf2d(filenc, varname)

# Crete the fine grid
(lon_f, lat_f)=lon_lat_fine(lon_mat, lat_mat, nf)

ns=size(lon_f)
println("Output size: ", ns)

if imin==0
    imin = 1
    jmin = 1
    if length(ns)==2
        imax = ns[1]; jmax = ns[2]
    else
        imax = ns[1]; (jmax,) = size(lat_f)
    end
else
    println("Cutout region: ", imin, "/", imax, "/", jmin, "/", jmax)
end

if length(ns)==2
    lon_fr = lon_f[imin:imax,jmin:jmax]
    lat_fr = lat_f[imin:imax,jmin:jmax]
else
    lon_fr = lon_f[imin:imax]
    lat_fr = lat_f[jmin:jmax]
end

if sx==0.
    # Compute fft3d and slope
    (fxp, ftp) = fft3d(pr)
    sx = fitslope(fxp, kmin=kmin)
    println("Computed spatial spectral slope: ", sx)
else
    println("Fixed spatial spectral slope: ", sx)
end

if fglob
    println("Conserving only global precipitation")
end

if weightsnc!=""
    println("Using weights file ", weightsnc)
    (ww, lon_mat2, lat_mat2) = read_netcdf2d(weightsnc, "")
else
    ww = 1.
end
# Downscaling
for iens = 1:nens
    @printf("Realization %d\n", iens)
    @time rd=rainfarm(pr, sx, nf, ww, fglob=fglob, fsmooth=fsmooth, verbose=true)
    fname = @sprintf("%s_%04d.nc", fnbase, iens)
    write_netcdf2d(fname, rd[imin:imax,jmin:jmax,:], lon_fr, lat_fr, varname, filenc)
end


text:
module ModuleRotation
using MultipleViewGeometry.ModuleOperators, MultipleViewGeometry.ModuleMathAliases
using LinearAlgebra
export rotx, roty, rotz, rotxyz, rodrigues2matrix
include("rotations.jl")
end


text:
import FFNetwork
using Test

@testset "FFNetwork.jl" begin

    println();
    printstyled("Training with MNIST,\n This could take a few minutes and migth require internet!!!"; color = :green);
    println();

    FFNetwork.loadMNIST();
    traindata, testdata = FFNetwork.get_MNIST_train_and_test_data();

    println();
    println("Train set length: 50000");
    println("Mini_batch_size 20");
    println("Eta 3.0");
    println("Test set length: 10000");
    println();

    println("Running 10 epochs ");
    println("Negative is Good, Positive is Bad");
    println();

    println("Using ReLU function as act_function ");
    println();
    FFNetwork.act_function(z) = FFNetwork.ReLU(z);
    FFNetwork.act_function_prime(z) = FFNetwork.ReLU_prime(z);
    netReLU = FFNetwork.Network(28*28,16,10);
    FFNetwork.train!(netReLU, FFNetwork.data.(traindata) , FFNetwork.data.(testdata) , 10 , 20, 3.0; verbose = true);

    println("Using Sigmoid function as act_function ");
    println();
    FFNetwork.act_function(z) = FFNetwork.sigmoid(z);
    FFNetwork.act_function_prime(z) = FFNetwork.sigmoid_prime(z);
    netSigmoid = FFNetwork.Network(28*28,16,10);
    FFNetwork.train!(netSigmoid, FFNetwork.data.(traindata) , FFNetwork.data.(testdata) , 10, 20, 3.0; verbose = true);

    sigmoidacc = FFNetwork.evaluate!(netSigmoid, FFNetwork.get_MNIST_img_data());
    ReLUacc = FFNetwork.evaluate!(netReLU, FFNetwork.get_MNIST_img_data());
    println();

    println("Evaluation over whole MNIST data");
    println("Sigmoid accuracy $sigmoidacc");
    println("ReLU accuracy $ReLUacc");
    println();

    @test 1==1;

end


text:
"""
   openxlsxtemplate(filepath::AbstractString) :: XLSXFile

Open an Excel file as template for editing and saving to another file with `XLSX.writexlsx`.
"""
@inline openxlsxtemplate(filepath::AbstractString) :: XLSXFile = open_or_read_xlsx(filepath, true, true, true)

"""
    writexlsx(output_filepath, xlsx_file; [rewrite])
"""
function writexlsx(output_filepath::AbstractString, xf::XLSXFile; rewrite::Bool=false)
    @assert is_writable(xf) "XLSXFile instance is not writable."
    @assert !isopen(xf) "Can't save an open XLSXFile."
    @assert all(values(xf.files)) "Some internal files were not loaded into memory. Did you use `XLSX.openxlsxtemplate` to open this file?"
    if !rewrite
        @assert !isfile(output_filepath) "Output file $output_filepath already exists."
    end

    update_worksheets_xml!(xf)

    xlsx = ZipFile.Writer(output_filepath)

    # write XML files
    for f in keys(xf.files)
        if f == "xl/sharedStrings.xml"
            # sst will be generated below
            continue
        end

        io = ZipFile.addfile(xlsx, f)
        EzXML.print(io, xf.data[f])
    end

    # write binary files
    for f in keys(xf.binary_data)
        io = ZipFile.addfile(xlsx, f)
        ZipFile.write(io, xf.binary_data[f])
    end

    if !isempty(get_sst(xf))
        io = ZipFile.addfile(xlsx, "xl/sharedStrings.xml")
        print(io, generate_sst_xml_string(get_sst(xf)))
    end

    close(xlsx)

    nothing
end

get_worksheet_internal_file(ws::Worksheet) = "xl/" * get_relationship_target_by_id(get_workbook(ws), ws.relationship_id)
get_worksheet_xml_document(ws::Worksheet) = get_xlsxfile(ws).data[ get_worksheet_internal_file(ws) ]

function set_worksheet_xml_document!(ws::Worksheet, xdoc::EzXML.Document)
    xf = get_xlsxfile(ws)
    filename = get_worksheet_internal_file(ws)
    @assert haskey(xf.data, filename) "Internal file not found for $(ws.name)."
    xf.data[filename] = xdoc
end

function generate_sst_xml_string(sst::SharedStrings) :: String
    @assert sst.is_loaded "Can't generate XML string from a Shared String Table that is not loaded."
    buff = IOBuffer()

    # TODO: <sst count="89"
    print(buff, """
<?xml version="1.0" encoding="UTF-8" standalone="yes"?><sst uniqueCount="$(length(sst))" xmlns="http://schemas.openxmlformats.org/spreadsheetml/2006/main">
""")

    for s in sst.formatted_strings
        print(buff, s)
    end

    print(buff, "</sst>")
    return String(take!(buff))
end

function update_worksheets_xml!(xl::XLSXFile)
    buff = IOBuffer()

    wb = get_workbook(xl)
    for i in 1:sheetcount(wb)
        sheet = getsheet(wb, i)
        doc = get_worksheet_xml_document(sheet)
        xroot = EzXML.root(doc)

        # check namespace and root node name
        @assert get_default_namespace(xroot) == SPREADSHEET_NAMESPACE_XPATH_ARG[1][2] "Unsupported Spreadsheet XML namespace $(get_default_namespace(xroot))."
        @assert EzXML.nodename(xroot) == "worksheet" "Malformed Excel file. Expected root node named `worksheet` in worksheet XML file."

        # forces a document copy to avoid crash: munmap_chunk(): invalid pointer
        EzXML.print(buff, doc)
        doc_copy = EzXML.parsexml(String(take!(buff)))

        # deletes all elements under sheetData
        child_nodes = EzXML.find(EzXML.root(doc_copy), "/xpath:worksheet/xpath:sheetData/xpath:row", SPREADSHEET_NAMESPACE_XPATH_ARG)
        for c in child_nodes
            EzXML.unlink!(c)
        end
        c = nothing
        child_nodes = nothing

        # updates sheetData
        sheetData_node = EzXML.findfirst(EzXML.root(doc_copy), "/xpath:worksheet/xpath:sheetData", SPREADSHEET_NAMESPACE_XPATH_ARG)
        spans_str = string(column_number(get_dimension(sheet).start), ":", column_number(get_dimension(sheet).stop))

        # iterates over WorksheetCache cells and write the XML
        for r in eachrow(sheet)
            ordered_column_indexes = sort(collect(keys(r.rowcells)))

            row_node = EzXML.addelement!(sheetData_node, "row")
            row_node["r"] = string(row_number(r))
            row_node["spans"] = spans_str

            # add cells to row
            for c in ordered_column_indexes
                cell = getcell(r, c)
                c_element = EzXML.addelement!(row_node, "c")

                c_element["r"] = cell.ref.name

                if cell.datatype != ""
                    c_element["t"] = cell.datatype
                end

                if cell.style != ""
                    c_element["s"] = cell.style
                end

                if cell.formula != ""
                    f_node = EzXML.addelement!(c_element, "f")
                    EzXML.setnodecontent!(f_node, cell.formula)
                end

                if cell.value != ""
                    v_node = EzXML.addelement!(c_element, "v")
                    EzXML.setnodecontent!(v_node, cell.value)
                end
            end
        end

        # updates worksheet dimension
        dimension_node = EzXML.findfirst(EzXML.root(doc_copy), "/xpath:worksheet/xpath:dimension", SPREADSHEET_NAMESPACE_XPATH_ARG)
        dimension_node["ref"] = string(get_dimension(sheet))

        set_worksheet_xml_document!(sheet, doc_copy)
    end

    nothing
end

function setdata!(ws::Worksheet, cell::Cell)
    @assert is_writable(get_xlsxfile(ws)) "XLSXFile instance is not writable."
    @assert !isnull(ws.cache) "Can't write data to a Worksheet with empty cache."
    cache = get(ws.cache)

    r = row_number(cell)
    c = column_number(cell)

    if !haskey(cache.cells, r)
        push!(cache.rows_in_cache, r)
        sort!(cache.rows_in_cache)
        cache.cells[r] = Dict{Int, Cell}()

        for i in 1:length(cache.rows_in_cache)
            if cache.rows_in_cache[i] == r
                cache.row_index[r] = i
            end
        end
    end
    cache.cells[r][c] = cell

    # update worksheet dimension
    ws_dimension = get_dimension(ws)

    top = row_number(ws_dimension.start)
    left = column_number(ws_dimension.start)

    bottom = row_number(ws_dimension.stop)
    right = column_number(ws_dimension.stop)

    if r < top || c < left
        top = min(r, top)
        left = min(c, left)
        set_dimension!(ws, CellRange(top, left, bottom, right))
    elseif r > bottom || c > right
        bottom = max(r, bottom)
        right = max(c, right)
        set_dimension!(ws, CellRange(top, left, bottom, right))
    end

    nothing
end

function setdata!(ws::Worksheet, ref::CellRef, v::Union{Int, Float64})
    cell = Cell(ref, "", "", string(v), "")
    setdata!(ws, cell)
end

function setdata!(ws::Worksheet, ref::CellRef, str::AbstractString)
    sst_ind = add_shared_string!(get_workbook(ws), str)
    cell = Cell(ref, "s", "", string(sst_ind), "")
    setdata!(ws, cell)
end

function setdata!(ws::Worksheet, ref::CellRef, m::Missings.Missing)
    # no-op
    nothing
end

function setdata!(ws::Worksheet, ref::CellRef, b::Bool)
    value_str = b ? "1" : "0"
    cell = Cell(ref, "b", "", value_str, "")
    setdata!(ws, cell)
end

const DEFAULT_DATE_numFmtId = 14
const DEFAULT_TIME_numFmtId = 20
const DEFAULT_DATETIME_numFmtId = 22

function setdata!(ws::Worksheet, ref::CellRef, date::Date)
    wb = get_workbook(ws)
    date_style_index = styles_get_cellXf_with_numFmtId(wb, DEFAULT_DATE_numFmtId)
    if date_style_index == -1
        # adds default style <xf applyNumberFormat="1" borderId="0" fillId="0" fontId="0" numFmtId="14" xfId="0"/>
        date_style_index = styles_add_cell_xf(wb, Dict("applyNumberFormat"=>"1", "borderId"=>"0", "fillId"=>"0", "fontId"=>"0", "numFmtId"=>"14", "xfId"=>"0"))
    end

    value_str = string(date_to_excel_value(date, isdate1904(get_xlsxfile(ws))))
    cell = Cell(ref, "", string(date_style_index), value_str, "")
    setdata!(ws, cell)
end

function setdata!(ws::Worksheet, ref::CellRef, t::Dates.Time)
    wb = get_workbook(ws)
    time_style_index = styles_get_cellXf_with_numFmtId(wb, DEFAULT_TIME_numFmtId)
    if time_style_index == -1
        # adds default style <xf applyNumberFormat="1" borderId="0" fillId="0" fontId="0" numFmtId="20" xfId="0"/>
        time_style_index = styles_add_cell_xf(wb, Dict("applyNumberFormat"=>"1", "borderId"=>"0", "fillId"=>"0", "fontId"=>"0", "numFmtId"=>"20", "xfId"=>"0"))
    end

    value_str = string(time_to_excel_value(t))
    cell = Cell(ref, "", string(time_style_index), value_str, "")
    setdata!(ws, cell)
end

function setdata!(ws::Worksheet, ref::CellRef, t::Dates.DateTime)
    wb = get_workbook(ws)
    datetime_style_index = styles_get_cellXf_with_numFmtId(wb, DEFAULT_DATETIME_numFmtId)
    if datetime_style_index == -1
        # adds default style <xf applyNumberFormat="1" borderId="0" fillId="0" fontId="0" numFmtId="22" xfId="0"/>
        datetime_style_index = styles_add_cell_xf(wb, Dict("applyNumberFormat"=>"1", "borderId"=>"0", "fillId"=>"0", "fontId"=>"0", "numFmtId"=>"22", "xfId"=>"0"))
    end

    value_str = string(datetime_to_excel_value(t, isdate1904(get_xlsxfile(ws))))
    cell = Cell(ref, "", string(datetime_style_index), value_str, "")
    setdata!(ws, cell)
end

setdata!(ws::Worksheet, ref_str::AbstractString, value::CellValue) = setdata!(ws, CellRef(ref_str), value)
setdata!(ws::Worksheet, ref_str::AbstractString, value) = error("Unsupported datatype $(typeof(value)) for writing data to Excel file. Supported data types are $(CellValue).")

Base.setindex!(ws::Worksheet, v, ref) = setdata!(ws, ref, v)

const DEFAULT_EXCEL_TEMPLATE = joinpath(dirname(@__FILE__), "..", "data", "blank.xlsx")

function open_default_template() :: XLSXFile
    @assert isfile(DEFAULT_EXCEL_TEMPLATE) "Couldn't find template file $DEFAULT_EXCEL_TEMPLATE."
    return openxlsxtemplate(DEFAULT_EXCEL_TEMPLATE)
end

function writetable!(sheet::Worksheet, data, columnnames; anchor_cell::CellRef=CellRef("A1"))

    # read dimensions
    col_count = length(data)
    @assert col_count == length(columnnames) "Column count mismatch between `data` ($col_count columns) and `columnnames` ($(length(columnnames)) columns)."
    @assert col_count > 0 "Can't write table with no columns."
    row_count = length(data[1])
    if col_count > 1
        for c in 2:col_count
            @assert length(data[c]) == row_count "Row count mismatch between column 1 ($row_count rows) and column $c ($(length(data[c])) rows)."
        end
    end

    anchor_row = row_number(anchor_cell)
    anchor_col = column_number(anchor_cell)

    # write table header
    for c in 1:col_count
        target_cell_ref = CellRef(anchor_row, c + anchor_col - 1)
        sheet[target_cell_ref] = string(columnnames[c])
    end

    # write table data
    for r in 1:row_count, c in 1:col_count
        target_cell_ref = CellRef(r + anchor_row, c + anchor_col - 1)
        sheet[target_cell_ref] = data[c][r]
    end
end

function rename!(ws::Worksheet, name::AbstractString)
    xf = get_xlsxfile(ws)
    @assert is_writable(xf) "XLSXFile instance is not writable."

    # updates XML
    xroot = xmlroot(xf, "xl/workbook.xml")
    for node in EzXML.eachelement(xroot)
        if EzXML.nodename(node) == "sheets"

            for sheet_node in EzXML.eachelement(node)
                if sheet_node["name"] == ws.name
                    # assign new name
                    sheet_node["name"] = name
                    break
                end
            end

            break
        end
    end

    # updates the new name in the worksheet instance
    ws.name = name
    nothing
end

const FILEPATH_SHEET_TEMPLATE = joinpath(dirname(@__FILE__), "..", "data", "sheet_template.xml")

addsheet!(xl::XLSXFile, name::AbstractString="") :: Worksheet = addsheet!(get_workbook(xl), name)

"""
    addsheet!(workbook, [name]) :: Worksheet

Create a new worksheet with named `name`.
If `name` is not provided, a unique name is created.

"""
function addsheet!(wb::Workbook, name::AbstractString="") :: Worksheet

    xf = get_xlsxfile(wb)
    @assert is_writable(xf) "XLSXFile instance is not writable."

    @assert isfile(FILEPATH_SHEET_TEMPLATE) "Couldn't find template file $FILEPATH_SHEET_TEMPLATE."

    if name == ""
        # name was not provided. Will find a unique name.
        i = 1
        current_sheet_names = sheetnames(wb)
        while true
            name = "Sheet$i"
            if !in(name, current_sheet_names)
                # found a unique name
                break
            end
            i += 1
        end
    end

    # generate sheetId
    current_sheet_ids = [ ws.sheetId for ws in wb.sheets ]
    sheetId = max(current_sheet_ids...) + 1

    xdoc = EzXML.readxml(FILEPATH_SHEET_TEMPLATE)

    # generate a unique name for the XML
    local xml_filename::String
    i = 1
    while true
        xml_filename = "xl/worksheets/sheet$i.xml"
        if !in(xml_filename, keys(xf.files))
            break
        end
        i += 1
    end

    # adds doc do XLSXFile
    xf.files[xml_filename] = true # is read
    xf.data[xml_filename] = xdoc

    # adds workbook-level relationship
    # <Relationship Id="rId1" Target="worksheets/sheet1.xml" Type="http://schemas.openxmlformats.org/officeDocument/2006/relationships/worksheet"/>
    rId = add_relationship!(wb, xml_filename[4:end], "http://schemas.openxmlformats.org/officeDocument/2006/relationships/worksheet")

    # creates Worksheet instance
    ws = Worksheet(xf, sheetId, rId, name, CellRange("A1:A1"))

    # creates a mock WorksheetCache
    itr = SheetRowStreamIterator(ws)
    zip_io, reader = open_internal_file_stream(xf, "[Content_Types].xml") # could be any file
    state = SheetRowStreamIteratorState(zip_io, reader, true, true, 0)
    close(state)
    ws.cache = WorksheetCache(CellCache(), Vector{Int}(), Dict{Int, Int}(), itr, state)

    push!(wb.sheets, ws)

    # updates workbook xml
    xroot = xmlroot(xf, "xl/workbook.xml")
    for node in EzXML.eachelement(xroot)
        if EzXML.nodename(node) == "sheets"

            #<sheet name="Sheet1" r:id="rId1" sheetId="1"/>
            sheet_element = EzXML.addelement!(node, "sheet")
            sheet_element["name"] = name
            sheet_element["r:id"] = rId
            sheet_element["sheetId"] = string(sheetId)

            break
        end
    end

    return ws
end


text:
import Base.isempty

#=
isempty : check if given Assoc is empty.

Note: Associ can be considered empty even if there are mapping for potential or past values.
=#
function isempty(A::Assoc)
    return isempty(A.A)
end

########################################################
# D4M: Dynamic Distributed Dimensional Data Model
# Architect: Dr. Jeremy Kepner (kepner@ll.mit.edu)
# Software Engineer: Alexander Chen (alexc89@mit.edu)
########################################################


text:
#=
@author : Spencer Lyon

@date: 07/09/2014
=#

using PyPlot

const r = 0.05
const beta = 1.0 / (1.0 + r)
const T = 20  # Time horizon
const S = 5   # Impulse date
const sigma1 = 0.15
const sigma2 = 0.15

function time_path(permanent=false)
    w1 = zeros(T+1)
    w2 = zeros(T+1)
    b = zeros(T+1)
    c = zeros(T+1)

    if permanent === false
        w2[S+2] = 1.0
    else
        w1[S+2] = 1.0
    end

    for t=2:T
        b[t+1] = b[t] - sigma2 * w2[t]
        c[t+1] = c[t] + sigma1 * w1[t+1] + (1 - beta) * sigma2 * w2[t+1]
    end

    return b, c
end

function main()
    fix, axes = subplots(2, 1)
    plt.subplots_adjust(hspace=0.5)
    p_args = {:lw=> 2, :alpha => 0.7}

    L = 0.175

    for ax in axes
        ax[:grid](alpha=0.5)
        ax[:set_xlabel]("Time")
        ax[:set_ylim](-L, L)
        ax[:plot]((S, S), (-L, L), "k-", lw=0.5)
    end

    ax = axes[1]
    b, c = time_path(false)
    ax[:set_title]("impulse-response, transitory income shock")
    ax[:plot](0:T, c, "g-", label="consumption"; p_args...)
    ax[:plot](0:T, b, "b-", label="debt"; p_args...)
    ax[:legend](loc="upper right")

    ax = axes[2]
    b, c = time_path(true)
    ax[:set_title]("impulse-response, permanent income shock")
    ax[:plot](0:T, c, "g-", label="consumption"; p_args...)
    ax[:plot](0:T, b, "b-", label="debt"; p_args...)
    ax[:legend](loc="lower right")

    return Void
end


text:
using SlurmAutoDetectWorkers
using Test
using Distributed
using ExistingProcessManagers

@testset "SlurmAutoDetectWorkers.jl" begin
    @test SlurmAutoDetectWorkers._foo(3) == 6
end


text:
module SurrealNumbers

using DocStringExtensions

export Surreal, left, right, SurrealForm, SurrealSet
export issurreal, isordinal, generation
export ∅

include("util.jl")
include("base.jl")
include("sets/SurrealSet.jl")
include("sets/subtypes.jl")
include("show.jl")
include("SurrealForm.jl")
include("builtins.jl")

end


text:
"""
Demonstration of solution of SYK equation by Matsubara frequency/imaginary time
alternation method using the discrete Lehmann representation

 The SYK equation is the nonlinear Dyson equation corresponding to self-energy

     Σ(τ) = J^2 * G(τ) * G(τ) * G(β-τ).

 We solve the Dyson equation self-consistently by a weighted
 fixed point iteration, with weight w assigned to the new iterate
 and weight 1-w assigned to the previous iterate. The self-energy
 is evaluated in the imaginary time domain, and each linear Dyson
 equation, corresponding to fixed self-energy, is solved in the
 Matsubara frequency domain, where it is diagonal.

 To solve the equation with a desired chemical potential μ,
 we pick a number nmu>1, and solve a sequence intermediate
 problems to obtain a good initial guess. First we solve the
 equation with chemical potential μ_0 = 0, then use this
 solution as an initial guess for the fixed point iteration with
 μ_1 = μ/nmu, then use this the solution as an initial guess
 for the fixed point iteration with μ = 2*μ/nmu, and so on,
 until we reach μ_{nmu} = μ.
"""

using Lehmann

diff(a, b) = maximum(abs.(a - b)) # return the maximum deviation between a and b

conformal_tau(τ, β) = -π^(1 / 4) / sqrt(2β) * 1 / sqrt(sin(π * τ / β))

function syk_sigma_dlr(d, G_x, J = 1.0)

    tau_k = d.τ # DLR imaginary time nodes
    tau_k_rev = d.β .- tau_k # Reversed imaginary time nodes

    G_k = dlr2tau(d, G_x) # G at tau_k
    G_k_rev = dlr2tau(d, G_x, tau_k_rev) # G at beta - tau_k

    # for i in 1:length(G_k)
    #     println("$(d.τ[i])    $(real(G_k[i]))     $(imag(G_k[i]))")
    # end

    # for i = 1:length(tau_k)
    #     G_k[i] = conformal_tau(tau_k[i], d.β)
    #     G_k_rev[i] = conformal_tau(tau_k_rev[i], d.β)
    # end

    Sigma_k = J .^ 2 .* G_k .^ 2 .* G_k_rev # SYK self-energy in imaginary time
    Sigma_x = tau2dlr(d, Sigma_k) # DLR coeffs of self-energy

    # println("sigma diff: ", diff(Sigma_k, dlr2tau(d, Sigma_x)))
    # for i in 1:length(Sigma_x)
    #     println("$(d.τ[i])    $(real(Sigma_x[i]))     $(imag(Sigma_x[i]))")
    # end
    # exit(0)
    return Sigma_x
end

function solve_syk_with_fixpoint_iter(d, mu, tol = d.rtol, mix = 0.3, maxiter = 100)

    Sigma_q = zeros(length(d)) # Initial guess
    G_q = zeros(ComplexF64, length(d))
    for iter in 1:maxiter

        G_q .= -1 ./ (d.ωn * 1im .- mu .+ Sigma_q) # Solve Dyson
        G_x = matfreq2dlr(d, G_q) # Get DLR coeffs
        Sigma_x_new = syk_sigma_dlr(d, G_x)
        Sigma_q_new = dlr2matfreq(d, Sigma_x_new)

        # println(diff(Sigma_q_new, Sigma_q))
        if maximum(abs.(Sigma_q_new .- Sigma_q)) < tol
            break
        end
        Sigma_q = mix * Sigma_q_new + (1 - mix) * Sigma_q # Linear mixing
    end
    return G_q
end

d = DLRGrid(Euv = 5.0, β = 1000.0, isFermi = true, rtol = 1e-14) # Initialize DLR object
G_q = solve_syk_with_fixpoint_iter(d, 0.0)

G_q = matfreq2tau(d, G_q)
println(G_q)
for i in 1:length(G_q)
    println("$(d.τ[i])    $(real(G_q[i]))     $(imag(G_q[i]))")
end


text:
############### print_human_readable(X::AbstractMultiScaleArray
function add_separator!(first_element,toprint,level)
    if first_element
        toprint[level][end] = "+"*toprint[level][end]
    end
    return false
end

function take_first_char(s::Union{String,SubString{String}},N::Int64)     s[1:min(length(s),N)]     end

function toprint_AbstractMultiScaleArray!(toprint,X::AbstractMultiScaleArray,fields,levelmax,n_char_per_name,level = 1)
    if length(toprint) < level
        push!(toprint,Vector{String}())
    end
    first_element = true
    if level<levelmax
        for x in X.nodes
            if fieldname(typeof(x),1) == :nodes
                push!(toprint[level],  "|"*take_first_char(split(string(typeof(x)),"{")[1],n_char_per_name)  )
                first_element = add_separator!(first_element,toprint,level)
                toprint_AbstractMultiScaleArray!(toprint,x,fields,levelmax,n_char_per_name,level+1)
            else
                if fields == nothing
                    push!(toprint[level],  take_first_char(split(string(typeof(x)),"{")[1],n_char_per_name)  )
                    first_element = add_separator!(first_element,toprint,level)
                else
                    if !(typeof(x) <: AbstractMultiScaleArrayLeaf)
                        push!(toprint[level],  take_first_char(split(string(typeof(x)),"{")[1],n_char_per_name)  )
                        first_element = add_separator!(first_element,toprint,level)
                    else
                        push!(toprint[level],  join([take_first_char(string(field),n_char_per_name)*": "*string(getfield(x,Symbol(field))) for field in collect(fields)],", ")  )
                        first_element = add_separator!(first_element,toprint,level)
                    end
                end
            end
        end
    end
end

function print_human_readable(X::AbstractMultiScaleArray ; n_char_per_name = 6, fields = nothing, levelmax=Inf, n_item_max_per_levels = Inf) # fields = nothing OR [field1,field2,...]
#     if typeof(X) <: AbstractMultiScaleArrayLeaf
#         println(X)
#     else
        toprint = Vector{Vector{String}}([[]]) # one vector per row, one string per element
        toprint_AbstractMultiScaleArray!(toprint,X,fields,levelmax,n_char_per_name)
        toprint = map(x -> join(x,"; "),toprint)
    #     toprint[1] = join(split(toprint[1]," "),"; ")
        Done = false
        while !Done
            Done = true # temporary
            for level1 in length(toprint):-1:2
                toprint0 = split(toprint[level1-1],"|")
                toprint1 = split(toprint[level1  ],"+")
                toprint0_ , toprint1_ = "", ""
                for (x0, x1) in zip(toprint0,toprint1)
                    if length(x0) > length(x1)
                        x1=join([x1,repeat(" ",length(x0)-length(x1))],"")
                        Done = false
                    elseif length(x0) < length(x1)
                        x0=join([x0,repeat(" ",length(x1)-length(x0))],"")
                        Done = false
                    end
                    toprint0_ = join([toprint0_, x0],"|")
                    toprint1_ = join([toprint1_, x1],"+")
                end
                toprint0_ = toprint0_[2:end]
                toprint1_ = toprint1_[2:end]
                toprint[level1-1] = toprint0_
                toprint[level1  ] = toprint1_
            end
        end
        # handle "ula; +         |Popu" cases
        #              ^^^^^^^^^^^
        for level in length(toprint):-1:1
            while occursin("+ ", toprint[level])
                Move = findfirst("+ ", toprint[level])[1]
                To   = findnext(" |", toprint[level], Move[end])[1]
                toprint[level] = join([toprint[level][1:(Move-1)],  toprint[level][(Move+1):To] , toprint[level][Move]  , toprint[level][(To+1):end]])
                for level_ in 1:(level-1)
                    toprint[level_] = join([toprint[level_][1:(Move-1)],  toprint[level_][(To+1):(To+To-Move)],  toprint[level_][Move:To],  toprint[level_][(To+To-Move+1):end]])
                end
            end
        end
        println(join(toprint,"\n"))
#     end
end


text:
using LazyViews
using Base.Test

# write your own tests here
@test 1 == 1

const fs = [
  x -> x,
  x -> 2*x,
  x -> sin(x),
  x -> x^3 + 3*x - 2,
  x -> sqrt(x)
];

const _NFPTESTS = 5000;
pdone = 10;

println("Floating point, $_NFPTESTS tests...");
for n=1:_NFPTESTS
  d   =   rand(1:5);
  ns  =   rand(1:25, d);
  a   =   rand(ns...);
  for f in fs
    mv    =   MapView{Float64, d}(f, a);
    lv    =   LazyView{Float64, d}(f, a);
    lvlv  =   LazyView{Float64, d}(f, lv);
    @test(mv == map(f, a));
    @test(lv == map(f, a));
    @test(map(f, a) == mv);
    @test(map(f, a) == lv);
    @test(mv == lv);
    @test(lvlv == map(f, map(f, a)));
    for (x, y) in zip(mv, map(f, a))
      @test(x == y);
    end
  end
  if (n % (_NFPTESTS / 10)) == 0
    print("$pdone%... ");
    pdone += 10;
  end
end
print_with_color(:green, "\nTEST PASSED.\n");

const cs = rand(-1e4:1e-9:1e4, 100);
const _SMTESTS = 100;
pdone = 10;

println("Scaler multiple, $_SMTESTS tests...");
for n=1:_SMTESTS
  d   =   rand(1:5);
  ns  =   rand(1:25, d);
  a   =   rand(ns...);
  for c in cs
    cv  =   CView{Float64, d}(c, a);
    lv  =   LazyView{Float64, d}(x -> c*x, a);
    @test(cv == map(x -> c*x, a));
    @test(cv == c * a);
    @test(c * a == cv);
    @test(cv == lv);
    @test((c * (c * a)) == (c * lv));
    @test((c * (c * a)) == (lv * c));
    for (x, y) in zip(cv, map(x -> c*x, a))
      @test(x == y);
    end
  end
  if (n % (_SMTESTS / 10)) == 0
    print("$pdone%... ");
    pdone += 10;
  end
end
print_with_color(:green, "\nTEST PASSED.\n");


text:
module MIMEFileExtensionsTests

using Documenter
using MIMEFileExtensions
using Test

macro s_str(x::AbstractString)
    QuoteNode(Symbol(x))
end

function test_list()
    dict = Dict(row.mime => row.fileexts for row in MIMEFileExtensions.list())
    @test [s"txt", s"log"] ⊆ dict[s"text/plain"]
    @test [s"jpeg", s"jpg"] ⊆ dict[s"image/jpeg"]
end

function test_fileexts_from_mime()
    @test [s"txt", s"log"] ⊆ fileexts_from_mime(s"text/plain")
    @test [s"jpeg", s"jpg"] ⊆ fileexts_from_mime("image/jpeg")
    @test [s"png"] ⊆ fileexts_from_mime(MIME"image/png"())
end

function test_mimes_from_fileext()
    @test [s"text/plain"] ⊆ MIMEFileExtensions.mimes_from_fileext(s"txt")
    @test [s"image/jpeg"] ⊆ MIMEFileExtensions.mimes_from_fileext("jpg")
end

function test_doctest()
    doctest(MIMEFileExtensions; manual = false)
end

end  # module MIMEFileExtensionsTests


text:
module RecSys

using ParallelSparseMatMul
using MAT
using Blobs

include("chunks/matrix.jl")
using .MatrixBlobs

if isless(Base.VERSION, v"0.5.0-")
    using SparseVectors
else
    nonzeroinds = Base.SparseArrays.nonzeroinds
    nonzeros = Base.SparseArrays.nonzeros
end

import Base: zero
import Blobs: save, load

export FileSpec, DlmFile, MatFile, SparseMat, SparseBlobs, DenseBlobs, read_input
export ALSWR, train, recommend, rmse, zero
export ParShmem, ParBlob
export save, load, clear, localize!

typealias RatingMatrix SparseMatrixCSC{Float64,Int64}
typealias SharedRatingMatrix ParallelSparseMatMul.SharedSparseMatrixCSC{Float64,Int64}
typealias InputRatings Union{RatingMatrix,SharedRatingMatrix,SparseMatBlobs}
typealias InputIdMap Union{Vector{Int64}, SharedVector{Int64}}
typealias ModelFactor Union{Matrix{Float64}, SharedArray{Float64,2}, DenseMatBlobs{Float64}}

abstract FileSpec
abstract Inputs
abstract Model

abstract Parallelism
type ParShmem <: Parallelism end
type ParBlob <: Parallelism end

if (Base.VERSION >= v"0.5.0-")
using Base.Threads
type ParThread <: Parallelism end
export ParThread
else
threadid() = 1
macro threads(x)
end
end

function tstr()
    t = time()
    string(Libc.strftime("%Y-%m-%dT%H:%M:%S",t), Libc.strftime("%z",t)[1:end-2], ":", Libc.strftime("%z",t)[end-1:end])
end

# enable logging only during debugging
#using Logging
##const logger = Logging.configure(filename="recsys.log", level=DEBUG)
#const logger = Logging.configure(level=DEBUG)
#macro logmsg(s)
#    quote
#        debug("[", myid(), "-", threadid(), "] ", $(esc(s)))
#    end
#end
macro logmsg(s)
end
#macro logmsg(s)
#    quote
#        info(tstr(), " [", myid(), "-", threadid(), "] ", $(esc(s)))
#    end
#end

include("inputs/inputs.jl")
include("models/models.jl")

include("als-wr.jl")
include("utils.jl")

end


text:
module Day01

using AdventOfCode2019

function day01(input::String = readInput(joinpath(@__DIR__, "..", "data", "day01.txt")))
    masses = parse.(Int, split(input))
    fuel = masses .÷ 3 .- 2

    solution2 = 0
    for elem in fuel
        while elem > 0
            solution2 += elem
            elem = elem ÷ 3 - 2
        end
    end

    return [fuel |> sum, solution2]
end

end # module


text:
module QuantizedNoise

export randqn

import Random: GLOBAL_RNG

"""
    randqn([rng,] T::Type{<:Integer}, dims::Integer...; σ=1.0, µ=0.0)
    randqn([rng,] T::Type{<:Integer}, dims::Integer...; std=1.0, mean=0.0)
    randqn([rng,] ::Type{<:Complex{T}}, dims::Integer...; σ=1.0, µ=0.0) where T<:Integer
    randqn([rng,] ::Type{<:Complex{T}}, dims::Integer...; std=1.0, mean=0.0) where T<:Integer

Generate samples from the normal distribution `𝒩(μ,σ²)` quantized and clamped
to the Integer subtype `T`.  Note that the `σ` and `μ` are properties of the
distribution from which the samples are drawn.  They are not necessarily
properties of the returned samples.  For example, quantization will alter the
standard deviation of the output values, so `std(randqn(Int8, 10^6, std=x))` may
be quite different from `x`, especially as `x` approaches (and exceeds)
`typemax(T)`.  Keyword arguments `std` and `mean` may be used instead of `σ` and
`μ`, resp.  If both forms are used, the single character keyword argument takes
precedence.

For complex types, the standard deviation `σ` or `std` may be given as a `Real`
or `Complex`.  If `σ` is given as a `Real`, it specifies the standard deviation
of the complex values (i.e. each real/imaginary component will have a standard
deviarion of `σ/sqrt(2)`). If given as a `Complex`, each real/imaginary
component of `σ` will be the standard deviation of the corresponding components
of the output.  Note that `σ=1.0` is equivalent to `σ=(1+1im)/sqrt(2)`.
Likewise, if `µ` is a `Real`, the complex mean will be `µ+µ*im`.
"""
function randqn(rng, T::Type{<:Integer}, dims::Integer...;
                    std=1.0, σ=std, mean=0.0, µ=mean)
    clamp.(round.(Integer, σ.*Base.randn(rng, Float32, dims...).+µ), T)
end

function randqn(T::Type{<:Integer}, dims::Integer...;
                    std=1.0, σ=std, mean=0.0, μ=mean)
    randqn(GLOBAL_RNG, T, dims...; σ=σ, μ=μ)
end

function randqn(rng, ::Type{<:Complex{T}}, dims::Integer...;
                    std=1.0, σ=std, mean=1.0, µ=mean) where T<:Integer
    σr, σi = (typeof(σ) <: Real) ? (σ, σ)./sqrt(2) : (real(σ), imag(σ))
    µr, µi = (typeof(µ) <: Real) ? (µ, µ) : (real(µ), imag(µ))

    complex.(randqn(rng, T, dims...; σ=σr, μ=μr),
             randqn(rng, T, dims...; σ=σi, μ=μi))
end

function randqn(T::Type{<:Complex{<:Integer}}, dims::Integer...;
                    std=1.0, σ=std, mean=0.0, μ=mean)
    randqn(GLOBAL_RNG, T, dims...; σ=σ, μ=μ)
end

end


text:
# © 2016 Massachusetts Institute of Technology.  See LICENSE file for details.
immutable EvaluatedInstance{T} <: Instance
  storage::Dict{ANode,T}
end

evaluate(n::ANode...) = evaluate(collect(ANode, n))
evaluate(T::Type, n::ANode...) = evaluate(T, collect(ANode, n))
evaluate(pool::StoragePool, n::ANode...) = evaluate(pool, collect(ANode, n))

evaluate{A<:ANode}(n::Vector{A}) = evaluate(Array{Float64}, n)
evaluate{A<:ANode}(T::Type, n::Vector{A}) = evaluate(NullPool(T), n)

function evaluate{A<:ANode}(pool::StoragePool, n::Vector{A})
  sorted = toposort(n)

  T = atype(pool)
  dict = Dict{ANode,T}()

  refcnt = Dict{ANode,Vector{Int}}()
  for a in n
    if sym(a) == :noop continue end
    get!(refcnt, output_node(a), Int[0])[1] += 1
  end

  for a in sorted
    for b in a.input
      get!(refcnt, output_node(b), Int[0])[1] += 1
    end
  end

  for i in sorted
    @assert refcnt[output_node(i)][1] > 0

    if sym(i) == :input
      error("Input nodes are not allowed here (try passing an array to the ANode constructor).")
    elseif sym(i) == :const || sym(i) == :load
      if isa(arg(i), T)
        dict[i] = arg(i)
      else
        get!(refcnt, i, Int[0])[1] += 1
        t = dict[i] = allocate(pool, i.size...)
        copy!(t, arg(i))
      end
    else
      if mutates(i)
        dict[i] = dict[i.input[1]]
        do_forward!(i, map(j->dict[j], i.input)...)
      else
        dict[i] = allocate(pool, i.size...)
        do_forward!(i, dict[i], map(j->dict[j], i.input)...)
      end
    end

    for a in i.input
      aa = output_node(a)
      @assert refcnt[aa][1] > 0

      refcnt[aa][1] -= 1
      if refcnt[aa][1] == 0 && arg(aa) != dict[aa]
        deallocate(pool, dict[aa])
        delete!(dict, aa)
      end
    end
  end

  for (n,r) in refcnt
    if r[1] == 0
      continue
    elseif r[1] == 1
      r[1] -= 1
      if dict[n] !== arg(n)
        deallocate(pool, dict[n])
      end
    else
      error("invalid refcnt")
    end
  end

  return EvaluatedInstance(dict)
end