text:
"""
GLOBAL_PAGE_VARS

Dictionary of variables assumed to be set for the entire website. Entries have the format
KEY => PAIR where KEY is a string (e.g.: "author") and PAIR is a pair where the first element is
the default value for the variable and the second is a tuple of accepted possible (super)types
for that value. (e.g.: "THE AUTHOR" => (String, Nothing))

DEVNOTE: marked as constant for perf reasons but can be modified since Dict.
"""
const GLOBAL_PAGE_VARS = PageVars()

"""
$(SIGNATURES)

Convenience function to allocate default values of the global site variables. This is called once,
when JuDoc is started.
"""
@inline function def_GLOBAL_PAGE_VARS!()::Nothing
    empty!(GLOBAL_PAGE_VARS)
    GLOBAL_PAGE_VARS["author"]      = Pair("THE AUTHOR",   (String, Nothing))
    GLOBAL_PAGE_VARS["date_format"] = Pair("U dd, yyyy",   (String,))
    GLOBAL_PAGE_VARS["prepath"]     = Pair("",             (String,))
    # these must be defined for the RSS file to be generated
    GLOBAL_PAGE_VARS["website_title"] = Pair("",    (String,))
    GLOBAL_PAGE_VARS["website_descr"] = Pair("",    (String,))
    GLOBAL_PAGE_VARS["website_url"]   = Pair("",    (String,))
    # if set to false, nothing rss will be considered
    GLOBAL_PAGE_VARS["generate_rss"]  = Pair(true,  (Bool,))
    return nothing
end

"""
CODE_SCOPE

Page-related struct to keep track of the code blocks that have been evaluated.
"""
mutable struct CodeScope
    rpaths::Vector{SubString}
    codes::Vector{SubString}
end
CodeScope() = CodeScope(String[], String[])

"""Convenience function to add a code block to the code scope."""
function push!(cs::CodeScope, rpath::SubString, code::SubString)::Nothing
    push!(cs.rpaths, rpath)
    push!(cs.codes, code)
    return nothing
end

"""Convenience function to (re)start a code scope."""
function reset!(cs::CodeScope)::Nothing
    cs.rpaths = []
    cs.codes  = []
    return nothing
end

"""Convenience function to purge code scope from head"""
function purgefrom!(cs::CodeScope, head::Int)
    cs.rpaths = cs.rpaths[1:head-1]
    cs.codes  = cs.codes[1:head-1]
    return nothing
end

"""
LOCAL_PAGE_VARS

Dictionary of variables copied and then set for each page (through definitions). Entries have the
same format as for `GLOBAL_PAGE_VARS`.

DEVNOTE: marked as constant for perf reasons but can be modified since Dict.
"""
const LOCAL_PAGE_VARS = PageVars()

"""
$(SIGNATURES)

Convenience function to allocate default values of page variables. This is called every time a page
is processed.
"""
@inline function def_LOCAL_PAGE_VARS!()::Nothing
    # NOTE `jd_code` is the only page var we KEEP (stays alive)
    code_scope = get(LOCAL_PAGE_VARS, "jd_code_scope") do
        Pair(CodeScope(),  (CodeScope,)) # the "Any" is lazy but easy
    end
    empty!(LOCAL_PAGE_VARS)

    # Local page vars defaults
    LOCAL_PAGE_VARS["title"]      = Pair(nothing, (String, Nothing))
    LOCAL_PAGE_VARS["hasmath"]    = Pair(true,    (Bool,))
    LOCAL_PAGE_VARS["hascode"]    = Pair(false,   (Bool,))
    LOCAL_PAGE_VARS["date"]       = Pair(Date(1), (String, Date, Nothing))
    LOCAL_PAGE_VARS["lang"]       = Pair("julia", (String,)) # default lang for indented code
    LOCAL_PAGE_VARS["reflinks"]   = Pair(true,    (Bool,))   # whether there are reflinks or not
    # Table of contents controls
    LOCAL_PAGE_VARS["mintoclevel"] = Pair(1,      (Int,)) # set to 2 to ignore h1
    LOCAL_PAGE_VARS["maxtoclevel"] = Pair(100,    (Int,))
    # CODE EVALUATION
    #
    LOCAL_PAGE_VARS["reeval"]        = Pair(false,  (Bool,)) # whether to always re-evals all on pg
    LOCAL_PAGE_VARS["freezecode"]    = Pair(false,  (Bool,)) # no-reevaluation of the code
    LOCAL_PAGE_VARS["showall"]       = Pair(false,  (Bool,)) # like a notebook on each cell
    # NOTE: when using literate, `literate_only` will assume that it's the only source of
    # code, so if it doesn't see change there, it will freeze the code to avoid an eval, this will
    # cause problems if there's more code on the page than from just the call to \literate
    # in such cases set literate_only to false.
    LOCAL_PAGE_VARS["literate_only"] = Pair(true,       (Bool,))
    # the jd_* should not be assigned externally
    LOCAL_PAGE_VARS["jd_code_scope"] = code_scope
    LOCAL_PAGE_VARS["jd_code_head"]  = Pair(Ref(0),     (Ref{Int},))
    LOCAL_PAGE_VARS["jd_code_eval"]  = Pair(Ref(false), (Ref{Bool},)) # toggle reeval
    LOCAL_PAGE_VARS["jd_code"]       = Pair("",         (String,))    # just the script

    # RSS 2.0 item specs:
    # only title, link and description must be defined
    #
    #     title       -- rss_title // fallback to title
    # (*) link        -- [automatically generated]
    #     description -- rss // rss_description NOTE: if undefined, no item generated
    #     author      -- rss_author // fallback to author
    #     category    -- rss_category
    #     comments    -- rss_comments
    #     enclosure   -- rss_enclosure
    # (*) guid        -- [automatically generated from link]
    #     pubDate     -- rss_pubdate // fallback date // fallback jd_ctime
    # (*) source      -- [unsupported assumes for now there's only one channel]
    #
    LOCAL_PAGE_VARS["rss"]             = Pair("", (String,))
    LOCAL_PAGE_VARS["rss_description"] = Pair("", (String,))
    LOCAL_PAGE_VARS["rss_title"]       = Pair("",      (String,))
    LOCAL_PAGE_VARS["rss_author"]      = Pair("",      (String,))
    LOCAL_PAGE_VARS["rss_category"]    = Pair("",      (String,))
    LOCAL_PAGE_VARS["rss_comments"]    = Pair("",      (String,))
    LOCAL_PAGE_VARS["rss_enclosure"]   = Pair("",      (String,))
    LOCAL_PAGE_VARS["rss_pubdate"]     = Pair(Date(1), (Date,))

    # page vars used by judoc, should not be accessed or defined
    LOCAL_PAGE_VARS["jd_ctime"]  = Pair(Date(1), (Date,))   # time of creation
    LOCAL_PAGE_VARS["jd_mtime"]  = Pair(Date(1), (Date,))   # time of last modification
    LOCAL_PAGE_VARS["jd_rpath"]  = Pair("",      (String,)) # local path to file src/[...]/blah.md

    # If there are GLOBAL vars that are defined, they take precedence
    local_keys = keys(LOCAL_PAGE_VARS)
    for k in keys(GLOBAL_PAGE_VARS)
        k in local_keys || continue
        LOCAL_PAGE_VARS[k] = GLOBAL_PAGE_VARS[k]
    end
    return nothing
end

"""
PAGE_HEADERS

Keep track of seen headers. The key is the refstring, the value contains the title,
the occurence number for the first appearance of that title and the level (1, ..., 6).
"""
const PAGE_HEADERS = LittleDict{AS,Tuple{AS,Int,Int}}()

"""
$(SIGNATURES)

Empties `PAGE_HEADERS`.
"""
@inline function def_PAGE_HEADERS!()::Nothing
    empty!(PAGE_HEADERS)
    return nothing
end

"""
PAGE_FNREFS

Keep track of name of seen footnotes; the order is kept as it's a list.
"""
const PAGE_FNREFS = String[]

"""
$(SIGNATURES)

Empties `PAGE_FNREFS`.
"""
@inline function def_PAGE_FNREFS!()::Nothing
    empty!(PAGE_FNREFS)
    return nothing
end

"""
PAGE_LINK_DEFS

Keep track of link def candidates
"""
const PAGE_LINK_DEFS = LittleDict{String,String}()

"""
$(SIGNATURES)

Empties `PAGE_LINK_DEFS`.
"""
@inline function def_PAGE_LINK_DEFS!()::Nothing
    empty!(PAGE_LINK_DEFS)
    return nothing
end

"""
GLOBAL_LXDEFS

List of latex definitions accessible to all pages. This is filled when the config file is read
(via manager/file_utils/process_config).
"""
const GLOBAL_LXDEFS = LittleDict{String, LxDef}()

"""
EMPTY_SS

Convenience constant for an empty substring, used in LXDEFS.
"""
const EMPTY_SS = SubString("")

"""
$(SIGNATURES)

Convenience function to allocate default values of global latex commands accessible throughout
the site. See [`resolve_lxcom`](@ref).
"""
@inline function def_GLOBAL_LXDEFS!()::Nothing
    empty!(GLOBAL_LXDEFS)
    # hyperreferences
    GLOBAL_LXDEFS["\\eqref"]    = LxDef("\\eqref",    1, EMPTY_SS)
    GLOBAL_LXDEFS["\\cite"]     = LxDef("\\cite",     1, EMPTY_SS)
    GLOBAL_LXDEFS["\\citet"]    = LxDef("\\citet",    1, EMPTY_SS)
    GLOBAL_LXDEFS["\\citep"]    = LxDef("\\citep",    1, EMPTY_SS)
    GLOBAL_LXDEFS["\\label"]    = LxDef("\\label",    1, EMPTY_SS)
    GLOBAL_LXDEFS["\\biblabel"] = LxDef("\\biblabel", 2, EMPTY_SS)
    GLOBAL_LXDEFS["\\toc"]      = LxDef("\\toc",      0, EMPTY_SS)
    GLOBAL_LXDEFS["\\tableofcontents"] = LxDef("\\tableofcontents", 0, EMPTY_SS)
    # inclusion
    GLOBAL_LXDEFS["\\input"]      = LxDef("\\input",      2, EMPTY_SS)
    GLOBAL_LXDEFS["\\output"]     = LxDef("\\output",     1, EMPTY_SS)
    GLOBAL_LXDEFS["\\codeoutput"] = LxDef("\\codeoutput", 1, subs("@@code_output \\output{#1}@@"))
    GLOBAL_LXDEFS["\\textoutput"] = LxDef("\\textoutput", 1, EMPTY_SS)
    GLOBAL_LXDEFS["\\textinput"]  = LxDef("\\textinput",  1, EMPTY_SS)
    GLOBAL_LXDEFS["\\show"]       = LxDef("\\show",       1, EMPTY_SS)
    GLOBAL_LXDEFS["\\figalt"]     = LxDef("\\figalt",     2, EMPTY_SS)
    GLOBAL_LXDEFS["\\fig"]        = LxDef("\\fig",        1, subs("\\figalt{}{#1}"))
    GLOBAL_LXDEFS["\\file"]       = LxDef("\\file",       2, subs("[#1]()"))
    GLOBAL_LXDEFS["\\tableinput"] = LxDef("\\tableinput", 2, EMPTY_SS)
    GLOBAL_LXDEFS["\\literate"]   = LxDef("\\literate",   1, EMPTY_SS)
    # text formatting
    GLOBAL_LXDEFS["\\underline"] = LxDef("\\underline", 1,
                            subs("~~~!#1~~~"))
    GLOBAL_LXDEFS["\\textcss"]   = LxDef("\\underline", 2,
                            subs("~~~!#2~~~"))
    return nothing
end

#= ==========================================
Convenience functions related to the jd_vars
============================================= =#

"""
$(SIGNATURES)

Convenience function taking a `DateTime` object and returning the corresponding formatted string
with the format contained in `GLOBAL_PAGE_VARS["date_format"]`.
"""
jd_date(d::DateTime)::AS = Dates.format(d, GLOBAL_PAGE_VARS["date_format"].first)

"""
$(SIGNATURES)

Checks if a data type `t` is a subtype of a tuple of accepted types `tt`.
"""
check_type(t::DataType, tt::NTuple{N,DataType} where N)::Bool = any(<:(t, tᵢ) for tᵢ ∈ tt)

"""
$(SIGNATURES)

Take a var dictionary `dict` and update the corresponding pair. This should only be used internally
as it does not check the validity of `val`. See [`write_page`](@ref) where it is used to store a
file's creation and last modification time.
"""
set_var!(d::PageVars, k::K, v) where K = (d[k] = Pair(v, d[k].second); nothing)

#= =================================================
set_vars, the key function to assign site variables
==================================================== =#

"""
$(SIGNATURES)

Given a set of definitions `assignments`, update the variables dictionary `jd_vars`. Keys in
`assignments` that do not match keys in `jd_vars` are ignored (a warning message is displayed).
The entries in `assignments` are of the form `KEY => STR` where `KEY` is a string key (e.g.:
"hasmath") and `STR` is an assignment to evaluate (e.g.: "=false").
"""
function set_vars!(jd_vars::PageVars, assignments::Vector{Pair{String,String}})::PageVars
    # if there's no assignment, cut it short
    isempty(assignments) && return jd_vars
    # process each assignment in turn
    for (key, assign) ∈ assignments
        # at this point there may still be a comment in the assignment
        # string e.g. if it came from @def title = "blah"
        # so let's strip


text:
using GLMakie, NetCDF, ColorSchemes
GLMakie.activate!()
using YAXArrays.Datasets: open_mfdataset

cs = collect(keys(colorschemes))
cmapIdx = Node(1)
cmapIdx
cmap = @lift(cs[$cmapIdx]) # @lift(Reverse(cs[$cmapIdx]))
ncolors = 101
colors = @lift(to_colormap($cmap, ncolors))
# some transparencies in the colors
g(x) = 0.5 + 0.5*tanh((x+3)/3)
alphas = g.(LinRange(-10,10, ncolors))
cmap_alpha = @lift(RGBAf.($colors, alphas))

dataset = open_mfdataset(joinpath(@__DIR__, "ETOPO1_halfdegree.nc"))
lon = dataset.lon.values
lat = dataset.lat.values
data = dataset.ETOPO1avg.data[:,:]
mindata = minimum(data)
front = data[:, end]
back = data[:, 1]
left = data[1, :]
right = data[end, :]

lower = [Point3f(i,-90, mindata) for i in lon]
upper = [Point3f(i,-90, front[idx]) for (idx,i) in enumerate(lon)]

lowerB = [Point3f(i,90, mindata) for i in lon]
upperB = [Point3f(i,90, back[idx]) for (idx,i) in enumerate(lon)]

lowerL = [Point3f(-180, i, mindata) for i in lat]
upperL = [Point3f(-180,i, left[idx]) for (idx,i) in enumerate(lat)]

lowerR = [Point3f(180, i, mindata) for i in lat]
upperR = [Point3f(180,i, right[idx]) for (idx,i) in enumerate(lat)]

fig = Figure()
ax = LScene(fig[1, 1], scenekw=(show_axis=false,))
mini, maxi = extrema(data)
data_norm = ((data .- mini) ./ (maxi-mini)) .* 20
crange = (((-6500,5500) .- mini) ./ (maxi-mini) .* 20)
pltobj = surface!(ax, lon, lat, data_norm, colorrange=crange, colormap=[:black, :brown, :green], backlight = 1.0f0)
color = :black
display(fig)
# surface!(ax, lon, lat, fill(0, size(data)); colorrange = (0,1),
    # lowclip = color, shading = false, transparency = true)
# objb = band!(ax, lower, upper; color = color)
# band!(ax, lowerL, upperL; color = color)
# band!(ax, lowerR, upperR; color = color)
# band!(ax, lowerB, upperB; color = color)
using RadeonProRender, GeometryBasics, Colors, Makie
using ReferenceTests, Colors
using RPRMakie
RPR = RadeonProRender

isdefined(Main, :context) && RPR.release(context)
context = RPR.Context()
RPRMakie.to_rpr_scene(context, ax.scene)
fb_size = (1500, 1500)
frame_buffer = RPR.FrameBuffer(context, RGBA, fb_size)
frame_bufferSolved = RPR.FrameBuffer(context, RGBA, fb_size)
set!(context, RPR.RPR_AOV_COLOR, frame_buffer)
set_standard_tonemapping!(context)
begin
    clear!(frame_buffer)
    RPR.rprContextSetParameterByKey1u(context, RPR.RPR_CONTEXT_ITERATIONS, 2)
    RPR.render(context)
    RPR.rprContextResolveFrameBuffer(context, frame_buffer, frame_bufferSolved, false)
    RPR.save(frame_bufferSolved, "test.png")
end

ext: jl


text:
using Flux
using Flux: onehot
using Flux.Data.Sentiment
using Flux.Batches: Tree, leaves

traintrees = Sentiment.train()

# Get the raw labels and phrases as separate trees.
labels  = map.(x -> x[1], traintrees)
phrases = map.(x -> x[2], traintrees)

# All tokens in the training set.
tokens = vcat(map(leaves, phrases)...)

# Count how many times each token appears.
freqs = Dict{String,Int}()
for t in tokens
  freqs[t] = get(freqs, t, 0) + 1
end

# Replace singleton tokens with an "unknown" marker.
# This roughly cuts our "alphabet" of tokens in half.
phrases = map.(t -> get(freqs, t, 0) == 1 ? "UNK" : t, phrases)

# Our alphabet of tokens.
alphabet = unique(vcat(map(leaves, phrases)...))

# One-hot-encode our training data with respect to the alphabet.
phrases_e = map.(t -> t == nothing ? t : onehot(t, alphabet), phrases)
labels_e  = map.(t -> onehot(t, 0:4), labels)

train = map.(tuple, phrases_e, labels_e)

max_stars_repo_name: domluna/model-zoo

max_issues_repo_name: domluna/model-zoo

max_forks_repo_name: domluna/model-zoo


text:
# fixit: also test @rput net2

@testset "rexport and sexp tests" begin

useape = false; # these extra tests would require the installation of ape by the CI platform,
# but ape itself requires a number of R packages.
# would make the tests taken longer. overkill.
# use it locally to run more checks when desired.

# on a tree with some edge lengths missing
s = "(A,(B:1.0,(C:1.0,D:1.0):1.0):1.0);";
tree1 = PhyloPlots.rexport(readTopology(s)); # for testing rexport
net2 = readTopology(s); # for testing sexp on clean network
if useape
## check for correct unrooted topology and for identical lengths
## of (unrooted) internal edges, but requires t
    R"library(ape)"
    R"tree2 = read.tree(text = $s)";
    @test convert(Bool, R"dist.topo($tree1, tree2, method='score') == 0")
    @test convert(Bool, R"dist.topo($net2, tree2, method='score') == 0")
end
R"""
tree1r = list(Nnode=3L, edge=matrix(c(5,5,6,6,7,7, 1,6,2,7,3,4), 6,2),
             tip.label=c("A","B","C","D"), edge.length=c(NA,1,1,1,1,1))
storage.mode(tree1r$edge) <- "integer"
class(tree1r) = "phylo"
tree2r = tree1r[c(2,1,4,3)]
class(tree2r) = "phylo"
""";
@test convert(Bool, R"isTRUE(all.equal($tree1,tree1r))")
@test convert(Bool, R"isTRUE(all.equal($net2, tree2r))")
@rput net2
@test convert(Bool, R"isTRUE(all.equal(net2, tree2r))")

# network, h=1, missing gamma values
s = "(((A:4.0,(B:1.0)#H1:1.1):0.5,(C:0.6,#H1:1.0):1.0):3.0,D:5.0);";
phy1 = PhyloPlots.rexport(readTopology(s));
net2 = readTopology(s);
if useape # needs ape version > 4.1 with read.evonet (not in 4.1)
    R"phyr = read.evonet(text = $s)"
    @test convert(Bool, R"dist.topo($phy1, phyr, method='score') == 0")
    @test convert(Bool, R"dist.topo($net2, phyr, method='score') == 0")
end
R"""
phy1r = list(Nnode=5L, edge=matrix(c(5,5,6,6,7,8,8,9, 6,4,8,7,3,1,9,2), 8,2),
             tip.label=c("A","B","C","D"), edge.length=c(3,5,.5,1,.6,4,1.1,1),
             reticulation=matrix(c(7L,9L),1,2), reticulation.length=1)
storage.mode(phy1r$edge) <- "integer"
class(phy1r) = c("evonet","phylo")
phy2r = phy1r[c(2,6,1,4,5,3)]
class(phy2r) = c("evonet","phylo")
""";
@test convert(Bool, R"isTRUE(all.equal($phy1, phy1r))")
@test convert(Bool, R"isTRUE(all.equal($net2, phy2r))")
@rput net2
@test convert(Bool, R"isTRUE(all.equal( net2, phy2r))")

# network, h=1, mainTree=true; minor hybrid edge length missing
s = "(((A:4.0,(B)#H1:1.1::0.9):0.5,(C:0.6,#H1):1.0):3.0,D:5.0);";
stree = "(((A:4.0,B):0.5,C:1.6):3.0,D:5.0);"; # main tree
phy1 = PhyloPlots.rexport(readTopology(s); mainTree=true);
net2 = readTopology(s);
tree2 = majorTree(readTopology(s));
if useape
    R"tree2r = read.tree(text = $stree)";
    @test convert(Bool, R"dist.topo($phy1, tree2r, method='score') == 0")
    @test convert(Bool, R"dist.topo($tree2, tree2r, method='score') == 0")
end
R"""
tree1r = list(Nnode=3L, edge=matrix(c(5,5,6,6,7,7, 6,4,7,3,1,2), 6,2),
             tip.label=c("A","B","C","D"), edge.length=c(3,5,.5,1.6,4,NA))
storage.mode(tree1r$edge) <- "integer"
class(tree1r) = "phylo"
phy2r = list(edge=matrix(c(5,5,6,6,7,8,8,9, 6,4,8,7,3,1,9,2), 8,2),
             Nnode=5L, edge.length=c(3,5,.5,1,.6,4,1.1,NA),
             reticulation=matrix(c(7L,9L),1,2),
             reticulation.gamma=0.1,
             tip.label=c("A","B","C","D"))
storage.mode(phy2r$edge) <- "integer"
class(phy2r) = c("evonet","phylo")
""";
@test convert(Bool, R"isTRUE(all.equal($phy1, tree1r))")
@test convert(Bool, R"isTRUE(all.equal($net2, phy2r))")
@rput net2
@test convert(Bool, R"isTRUE(all.equal( net2, phy2r))")

# on a network, h=1, with useEdgeLength=false
s = "(((A:1.0,(B:.5)#H1:.2::0.9):.8,(C:1.5,#H1:.01::0.1):2):.6,D:5.0);";
phy1 = PhyloPlots.rexport(readTopology(s), useEdgeLength=false);
if useape # needs ape version > 4.1 with read.evonet (not in 4.1), also dist.topo for networks
    R"phyr = read.evonet(text = $s)"
    @test convert(Bool, R"dist.topo($phy1, phyr) == 0")
end
R"""
phy1r = list(Nnode=5L, edge=matrix(c(5,5,6,6,7,8,8,9, 6,4,8,7,3,1,9,2), 8,2),
             tip.label=c("A","B","C","D"), reticulation=matrix(c(7L,9L),1,2),
             reticulation.gamma=0.1)
storage.mode(phy1r$edge) <- "integer"
class(phy1r) = c("evonet","phylo")
""";
@rput phy1
@test convert(Bool, R"isTRUE(all.equal(phy1, phy1r))")

# on a network with h=2 hybridizations
s = "(((Ag,(#H1:7.159::0.056,((Ak,(E:0.08,#H2:0.0::0.004):0.023):0.078,(M:0.0)#H2:::0.996):2.49):2.214):0.026,
      (((((Az:0.002,Ag2:0.023):2.11,As:2.027):1.697)#H1:0.0::0.944,Ap):0.187,Ar):0.723):5.943,(P,20):1.863,165);";
phy1 = PhyloPlots.rexport(readTopology(s));
net2 = readTopology(s);
R"""
phy1r = list(Nnode=14L, edge=matrix(c(13,13,13,14,14,15,15,16,16,17,17,18,18,19,20,20,21,21,22,23,24,25,25,26,26,
             15,14,12,10,11,18,16,17,9,24,8,1,19,20,21,23,2,22,3,4,25,26,7,5,6), 25,2),
             tip.label=c("Ag","Ak","E","M","Az","Ag2","As","Ap","Ar","P","20","165"),
             edge.length=c(5.943,1.863,NA,NA,NA,.026,.723,.187,NA,0,NA,NA,2.214,2.490,.078,NA,NA,
                           .023,.08,0,1.697,2.11,2.027,.002,.023),
             reticulation=matrix(c(19,22, 24,23),2,2),
             reticulation.gamma=c(0.056,0.004),
             reticulation.length=c(7.159,0.0))
storage.mode(phy1r$edge) <- "integer"
storage.mode(phy1r$reticulation) <- "integer"
class(phy1r) = c("evonet","phylo")
phy2r = phy1r[c(2,7,1,4,5,6,3)]
class(phy2r) = c("evonet","phylo")
""";
@test convert(Bool, R"identical($phy1, phy1r)")
@test convert(Bool, R"identical($net2, phy2r)")

end

max_stars_repo_name: cecileane/PhyloPlots

max_issues_repo_name: cecileane/PhyloPlots

max_forks_repo_name: cecileane/PhyloPlots


text:
using Test
using Unitful

#EmptyDiffusivity
dif = EmptyDiffusivity()
@test_throws ErrorException dif(T=400.0)

#ConstantDiffusivity
D = 4.0
dif = ConstantDiffusivity(D)
@test dif(T=400.0,mu=1e3) == D

#StokesDiffusivity
r = 0.47e-9
dif = StokesDiffusivity(r)
D = dif(T=500.0,mu=.891e-3)
@test kB*500.0/(6*pi*.891e-3*D) ≈ r rtol=1e-3


text:
using Test
if !@isdefined(interpret)
    include("../src/boolfuck.jl")
end

# brainfucktoboolfuck = let commands = Dict((
#     '+' => ">[>]+<[+<]>>>>>>>>>[+]<<<<<<<<<",
#     '-' => ">>>>>>>>>+<<<<<<<<+[>+]<[<]>>>>>>>>>[+]<<<<<<<<<",
#     '<' => "<<<<<<<<<",
#     '>' => ">>>>>>>>>",
#     ',' => ">,>,>,>,>,>,>,>,<<<<<<<<",
#     '.' => ">;>;>;>;>;>;>;>;<<<<<<<<",
#     '[' => ">>>>>>>>>+<<<<<<<<+[>+]<[<]>>>>>>>>>[+<<<<<<<<[>]+<[+<]",
#     ']' => ">>>>>>>>>+<<<<<<<<+[>+]<[<]>>>>>>>>>]<[+<]"
# ))
# (brainfuck) -> join(get(commands, command, "") for command in brainfuck)
# end

@testset "Hello, world!" begin
    helloworld = """;;;+;+;;+;+;
        +;+;+;+;;+;;+;
        ;;+;;+;+;;+;
        ;;+;;+;+;;+;
        +;;;;+;+;;+;
        ;;+;;+;+;+;;
        ;;;;;+;+;;
        +;;;+;+;;;+;
        +;;;;+;+;;+;
        ;+;+;;+;;;+;
        ;;+;;+;+;;+;
        ;;+;+;;+;;+;
        +;+;;;;+;+;;
        ;+;+;+;"""

    @test interpret(helloworld) == "Hello, world!\n"
    @test compile(helloworld)() == "Hello, world!\n"
end

@testset "Multiply" begin
    multiply = ">,>,>,>,>,>,>,>,>>,>,>,>,>,>,>,>,<<<<<<<<+<<<<<<<<+[>+]<[<]>>>>>>>>>[+<<<<<<<<[>]+<[+<]>>>>>>>>>>>>>>>>>>+<<<<<<<<+[>+]<[<]>>>>>>>>>[+<<<<<<<<[>]+<[+<]>>>>>>>>>+<<<<<<<<+[>+]<[<]>>>>>>>>>[+]>[>]+<[+<]>>>>>>>>>[+]>[>]+<[+<]>>>>>>>>>[+]<<<<<<<<<<<<<<<<<<+<<<<<<<<+[>+]<[<]>>>>>>>>>]<[+<]>>>>>>>>>>>>>>>>>>>>>>>>>>>+<<<<<<<<+[>+]<[<]>>>>>>>>>[+<<<<<<<<[>]+<[+<]>>>>>>>>>+<<<<<<<<+[>+]<[<]>>>>>>>>>[+]<<<<<<<<<<<<<<<<<<<<<<<<<<[>]+<[+<]>>>>>>>>>[+]>>>>>>>>>>>>>>>>>>+<<<<<<<<+[>+]<[<]>>>>>>>>>]<[+<]<<<<<<<<<<<<<<<<<<+<<<<<<<<+[>+]<[<]>>>>>>>>>[+]+<<<<<<<<+[>+]<[<]>>>>>>>>>]<[+<]>>>>>>>>>>>>>>>>>>>;>;>;>;>;>;>;>;<<<<<<<<"
    m = compile(multiply)
    @test m(join(Char[4, 5])) == join(Char[20])
end

@testset "Fibonacci Numbers" begin
    fibsbelow = let fibcache = [1, 1]
        function fib(n::Integer)
            while length(fibcache) < n
                push!(fibcache, fibcache[end-1] + fibcache[end])
            end
            view(fibcache, 1:n)
        end
    end
    julia_f(n) = join(fibsbelow(n), ", ")

    fibonacci = brainfucktoboolfuck(",>+>>>>++++++++++++++++++++++++++++++++++++++++++++>++++++++++++++++++++++++++++++++<<<<<<[>[>>>>>>+>+<<<<<<<-]>>>>>>>[<<<<<<<+>>>>>>>-]<[>++++++++++[-<-[>>+>+<<<-]>>>[<<<+>>>-]+<[>[-]<[-]]>[<<[>>>+<<<-]>>[-]]<<]>>>[>>+>+<<<-]>>>[<<<+>>>-]+<[>[-]<[-]]>[<<+>>[-]]<<<<<<<]>>>>>[++++++++++++++++++++++++++++++++++++++++++++++++.[-]]++++++++++<[->-<]>++++++++++++++++++++++++++++++++++++++++++++++++.[-]<<<<<<<<<<<<[>>>+>+<<<<-]>>>>[<<<<+>>>>-]<-[>>.>.<<<[-]]<<[>>+>+<<<-]>>>[<<<+>>>-]<<[<+>-]>[<+>-]<<<-]")
    f = compile(fibonacci)

    for n in rand(1:9, 3)
        @test julia_f(n) == f(String(Char[n]))
    end
end

max_stars_repo_name: stellartux/Esolangs

max_issues_repo_name: stellartux/Esolangs

max_forks_repo_name: stellartux/Esolangs


text:
@static if VERSION >= v"1.3"

    """
        (f::AbstractGP)(x::AbstractVector)

    Construct a `FiniteGP` representing the projection of `f` at `x`.
    """
    (f::AbstractGP)(x...) = FiniteGP(f, x...)
else

    """
        (f::Union{GP, CompositeGP})(x::AbstractVector)

    Construct a `FiniteGP` representing the projection of `f` at `x`.
    """
    (f::GP)(x...) = FiniteGP(f, x...)
    (f::CompositeGP)(x...) = FiniteGP(f, x...)
end


text:
using BinDeps
using Compat

@BinDeps.setup

nettle = library_dependency("nettle", aliases = ["libnettle","libnettle-4-6"])

@windows_only begin
  using WinRPM
  provides(WinRPM.RPM, "libnettle-4-6", nettle, os = :Windows )
end

@osx_only begin
  using Homebrew
  provides( Homebrew.HB, "nettle", nettle, os = :Darwin )
end

provides( AptGet, "libnettle4", nettle )
provides( Yum, "nettle", nettle )

julia_usrdir = normpath(JULIA_HOME*"/../") # This is a stopgap, we need a better built-in solution to get the included libraries
libdirs = AbstractString["$(julia_usrdir)/lib"]
includedirs = AbstractString["$(julia_usrdir)/include"]
env = @compat Dict("HOGWEED_LIBS" => "-L$(libdirs[1]) -L$(BinDeps.libdir(nettle)) -lhogweed -lgmp",
       "NETTLE_LIBS" => "-L$(libdirs[1]) -L$(BinDeps.libdir(nettle)) -lnettle -lgmp",
       "LD_LIBRARY_PATH" => join([libdirs[1];BinDeps.libdir(nettle)],":"))

provides( Sources,
          URI("http://www.lysator.liu.se/~nisse/archive/nettle-2.7.1.tar.gz"),
          SHA="bc71ebd43435537d767799e414fce88e521b7278d48c860651216e1fc6555b40",
          nettle )
provides( BuildProcess,
          Autotools(lib_dirs = libdirs,
                    include_dirs = includedirs,
                    env = env,
                    configure_options = ["--disable-openssl", "--libdir=$(BinDeps.libdir(nettle))"]),
          nettle )

@compat @BinDeps.install Dict(:nettle => :nettle)


text:
"""
Calculate the exact survival signature of `system` in regards to the given `types` and structure function `φ`.

This evaluates all possible combinations of network states and is time consuming for large networks. An optional
`preprocessor` can be used to exclude parts of the survival signature beforhand.
"""
function survivalsignature(
    system::Any,
    types::Dict{Int64,Array{Int64,1}},
    φ::Function,
    preprocessor = nothing,
)
    Φ, _ = prepare_survival_signature(types)

    preprocessor !== nothing && (preprocessor(Φ, system))

    indices = findall(x -> isinf(x), Φ)
    entries = @showprogress pmap(indices) do index
        return exactentry(index, system, types, φ)
    end

    for (index, entry) in zip(indices, entries)
        Φ[index] = entry
    end

    return Φ
end

function survivalsignature(
    system::Any,
    types::Dict{Int64,Array{Int64,1}},
    φ::Function,
    samples::Int64,
    limit::Float64 = 0.001,
    preprocessor = nothing,
)
    Φ, components_per_type = prepare_survival_signature(types)
    preprocessor !== nothing && (preprocessor(Φ, system))

    cov = zeros(size(Φ))

    indices = findall(x -> isinf(x), Φ)
    entries = @showprogress pmap(indices) do index
        n_c = numberofcombinations(components_per_type, index)
        if n_c <= samples
            return exactentry(index, system, types, φ), 0
        else
            return approximateentry(index, system, types, φ, samples, limit)
        end
    end

    for (index, entry) in zip(indices, entries)
        s, c = entry
        Φ[index] = s
        cov[index] = c
    end

    return Φ, cov
end

function exactentry(index::CartesianIndex, system, types, φ)
    combinations = []
    for (type, components) ∈ types
        push!(combinations, subsets(components, index[type] - 1))
    end

    results = map(
        x -> φ(system, Iterators.flatten(x) |> collect),
        Iterators.product(combinations...),
    )

    return count(results) / length(results)
end

function approximateentry(
    index::CartesianIndex,
    system,
    types::Dict{Int,Array{Int,1}},
    φ::Function,
    samples::Int,
    limit::Float64,
)
    working = 0
    Φ = 0
    cov = Inf

    for n = 1:samples
        # Generate a random network state for the current index
        x = []
        for (type, components) ∈ types
            shuffle!(components)
            append!(x, components[1:index[type]-1])
        end

        # Evaluate the structure function and update the current approximation of Φ
        φ(system, x) && (working += 1; Φ = working / n)
        Φ > 0 && (cov = √((Φ - Φ^2) / n) / Φ)

        # Break if the cov is below the defined threshold
        Φ != 1 && cov < limit && break

    end

    return Φ, cov
end

function numberofcombinations(components_per_type::Array{Int}, index::CartesianIndex)
    combinations = 1
    for (c, i) in zip(components_per_type, Tuple(index))
        combinations *= BigCombinatorics.Binomial(c - 1, i - 1)
    end
    return combinations
end

function prepare_survival_signature(types)
    components_per_type = ones(Int, length(types))
    for (type, components) in types
        components_per_type[type] += length(components)
    end

    Φ = ones(components_per_type...) * Inf

    return Φ, components_per_type
end


text:
using Assignment09
using Test

@testset "Assignment09" begin

@testset "Question 1" begin
    @test typeof(question1) <: AbstractString
end

@testset "Question 2" begin
    @test typeof(question2) <: AbstractString
end

@testset "Question 3" begin
    @test typeof(question3) <: AbstractString
end

@testset "Question 4" begin
    @test typeof(question4) <: AbstractString
end

@testset "Question 5" begin
    @test typeof(question5) <: AbstractString
end

end # Assignment09


text:
module test_Plot5

import ModiaMath

# Desired:
#   using Test
#
# In order that Test needs not to be defined in the user environment, it is included via ModiaMath:
using ModiaMath.Test
t = range(0.0, stop=10.0, length=100)

series = Dict{AbstractString,Any}()
series["time"] = t
series["phi"]  = sin.(t)
series["phi2"] = 0.5 * sin.(t)
series["w"]    = cos.(t)
series["w2"]   = 0.6 * cos.(t)
series["r"]    = hcat(0.5 * sin.(t), 0.2 * sin.(t), 0.3 * cos.(t))

time = ModiaMath.RealScalar(:time, unit="s")
phi  = ModiaMath.RealScalar(:phi,  unit="rad")
phi2 = ModiaMath.RealScalar(:phi2, unit="rad")
w    = ModiaMath.RealScalar(:w,    unit="rad/s")
w2   = ModiaMath.RealScalar(:w2,   unit="rad/s")
r    = ModiaMath.RealSVector3(:r,  unit="m")
var  = Dict{Symbol,Any}(:time => time, :phi => phi, :phi2 => phi2, :w => w, :w2 => w2, :r => r)

result = ModiaMath.ResultWithVariables(series, var, "")

ModiaMath.plot(result, [ (:phi, :r)      (:phi, :phi2, :w);
                         (:w, :w2, :phi2) (:phi, :w)      ], heading="Matrix of plots")

# Print result variables
println("\n... result variables = ", result)
end


text:
module QuerySynthesisStrategies

using Memento
using Reexport
using MLKernels

import JuMP
import MLLabelUtils
import LIBSVM
import SVDD

import ..QueryStrategies:
    QueryStrategy,
    HybridPQs,
    MissingLabelTypeException,
    QueryStrategy,

    get_query_objects,
    initialize_qs,
    qs_score

include("query_synthesis_base.jl")
include("query_synthesis_utils.jl")

include("Optimization/Optimization.jl")
@reexport using .Optimization

include("TestQss.jl")
include("RandomQss.jl")
include("RandomOutlierQss.jl")
include("DecisionBoundaryQss.jl")
include("NaiveExplorativeMarginQss.jl")
include("ExplorativeMarginQss.jl")

export
    QuerySynthesisStrategy,
    DataBasedQss, ModelBasedQss, HybridQss,
    # BaselineQss:
    TestQss, RandomQss,
    # DataBasedQss:
    # ModelBasedQss:
    DecisionBoundaryQss, RandomOutlierQss,
    # HybridQss:
    ExplorativeMarginQss, NaiveExplorativeMarginQss,

    get_query_object,
    estimate_margin_epsilon, estimate_limit_epsilon,
    check_epsilon, check_limits,
    extrema_arrays,
    rand_in_hypercube,
    estimate_boundary_shift_epsilon, data_limit_epsilon, data_boundaries
end


text:
function build_biophysical_dictionary(path_to_biophysical_constants_file::String)::Dict{String,Any}

    # check the path - is it legit?
    check_file_existence(path_to_biophysical_constants_file)

    # TODO: we need to check if the biophysical_constants_file is a TOML file
    # TODO: is this a TOML file?

    # Load the model dictionary -
    model_dictionary = TOML.parsefile(path_to_biophysical_constants_file)

    # initialize -
    biophysical_dictionary = Dict{String,Any}()

    # convert -
    list_of_constant_dictionaries = model_dictionary["biophysical_constants"]
    for (key_string, local_dictionary) in list_of_constant_dictionaries

        # get the value -
        value = parse(Float64, local_dictionary["value"])

        # cache -
        biophysical_dictionary[key_string] = value
    end

    # return -
    return biophysical_dictionary
end

function build_ec_data_dictionary(path_to_ec_data_file::String)::Dict{String,Any}

    # check the path - is it legit?
    check_file_existence(path_to_ec_data_file)

    # TODO: we need to check if the path_to_ec_data_file is a TOML file
    # TODO: is this a TOML file?

    # load the ec model dictionary -
    ec_data_dictionary = TOML.parsefile(path_to_ec_data_file)

    # return -
    return ec_data_dictionary
end

function build_control_constants_dictionary(path_to_ec_data_file::String)::Dict{String,Any}

    # check the path - is it legit?
    check_file_existence(path_to_ec_data_file)

    # TODO: we need to check if the path_to_ec_data_file is a TOML file
    # TODO: is this a TOML file?

    # load the ec model dictionary -
    control_parameter_data_dictionary = TOML.parsefile(path_to_ec_data_file)

    # return -
    return control_parameter_data_dictionary
end

function build_static_fba_problem_object(data_dictionary::Dict{String,Any}; is_minimum_flag::Bool=true)::VLStaticFBAProblem

    # takes stuff from the dd and bundles it up as a static problem object -

    # TODO: check = do we have the required fields in the dd?
    # TODO: error?

    # create a problem, and populate its fields -
    problem = VLStaticFBAProblem()
    problem.stoichiometric_matrix = data_dictionary["stoichiometric_matrix"]
    problem.flux_bounds_array = data_dictionary["flux_bounds_array"]
    problem.species_bounds_array = data_dictionary["species_bounds_array"]
    problem.obj_coeff_vector = data_dictionary["objective_coefficient_array"]
    problem.is_minimum_flag = is_minimum_flag
    return problem
end

function build_dynamic_fba_problem_object(data_dictionary::Dict{String,Any}; is_minimum_flag::Bool=true)::VLDynamicFBAProblem

    # takes stuff from the dd and bundles it up as a static problem object -

    # TODO: check = do we have the required fields in the dd?
    # TODO: error?

    # create a problem, and populate its fields -
    # for now, return blank problem -
    return VLDynamicFBAProblem()
end

max_stars_repo_name: varnerlab/CFMG

max_issues_repo_name: varnerlab/CFMG

max_forks_repo_name: varnerlab/CFMG


text:
using Test
import sbp
using sbp.StaggeredRK4: staggered_rk4

y = zeros(1)

@testset "Integrate cubic monomial" begin

        e = ones(1)

        function f(y::AbstractArray, t::Float64)
                return e * t ^ 3
        end

        function g(y::AbstractArray, t::Float64)
                return e * t ^ 3
        end

        nt = 10
        dt = 0.1
        u = zeros(1)
        v = zeros(1)
        U = zeros(nt+1)
        V = zeros(nt+1)

        ans_u = zeros(nt+1)
        ans_v = zeros(nt+1)
        for k in 1:(nt+1)
                ans_u[k] = 0.25 * ((k - 1) * dt) ^ 4
                ans_v[k] = 0.25 * ((k - 0.5) * dt) ^ 4
        end

        v[1] = ans_v[1]

        tk=0.0
        V[1] = v[1]
        for k in 1:nt
                tk = (k - 1) * dt
                @time u, v = staggered_rk4(f, g, u, v, tk, dt)
                U[k+1] = u[1]
                V[k+1] = v[1]
        end
        @test isapprox(ans_u, U)
        @test isapprox(ans_v, V)
end

max_stars_repo_name: ooreilly/sbpjl

max_issues_repo_name: ooreilly/sbpjl

max_forks_repo_name: ooreilly/sbpjl


text:
using OpenTelemetry
using Test

@testset "OpenTelemetry.jl" begin
    code = []
    is_code = false
    for line in readlines(joinpath(@__DIR__, "..", "README.md"))
        if line == "```julia"
            is_code = true
        elseif line == "```"
            is_code = false
            code_str = """
            begin
            $(join(code, "\n"))
            end
            """
            println("Executing the following code:")
            println(code_str)
            eval(Meta.parse(code_str))
            empty!(code)
        elseif is_code
            push!(code, line)
        end
    end
end


text:
using DiscreteFiltering
using OrdinaryDiffEq: ODEProblem, solve, QNDF
using LinearAlgebra
using Plots
using SparseArrays

## Filter
h₀ = 0.03
h(x) = h₀ * (1 - 1 / 2 * cos(x))
dh(x) = 0.0 # h₀ / 2 * sin(x)
α(x) = 1 / 3 * dh(x) * h(x)
filter = TopHatFilter(h)

## Domain
a = 0.0
b = 2π
domain = PeriodicIntervalDomain(a, b)
# domain = ClosedIntervalDomain(a, b)

## Discretization
N = 500
M = N
x = discretize(domain, M)
ξ = discretize(domain, N)
Δx = (b - a) / N

Δx^2 / maximum(abs.(α.(x)))

## Time
T = 1.0
t = T

## Plot filter
plot(x, h)
ylims!((0, ylims()[2]))

## Get matrices
C = advection_matrix(domain, N)
D = diffusion_matrix(domain, N)
W = filter_matrix(filter, domain, M, N)
R = reconstruction_matrix(filter, domain, M, N)
A = spdiagm(α.(ξ))

## Exact solutions
u₀(ξ) = sin(ξ) + 0.6cos(5ξ) + 0.04sin(20(ξ - 1))

## Discrete initial conditions
uₕ = u₀.(ξ, 0.0)
plot(ξ, uₕ, label = "Discretized")

## Extension to a non-linear case: Burgers equation
ν = 0.03
duₕ(uₕ, p, t) = -uₕ .* (C * uₕ) + ν * D * uₕ
prob = ODEProblem(duₕ, uₕ, (0, T))
sol = solve(prob, QNDF(), abstol = 1e-6, reltol = 1e-4)

function duₕ_allbar(uₕ_allbar, p, t)
    uₕ = R * uₕ_allbar
    -W * (uₕ .* (C * uₕ)) + ν * W * D * uₕ
end
prob_allbar = ODEProblem(duₕ_allbar, W * uₕ, (0, T))
sol_allbar = solve(prob_allbar, QNDF(), abstol = 1e-6, reltol = 1e-4)

plot(ξ, uₕ, label = "Initial")
plot!(x, W * uₕ, label = "Initial discretized-then-filtered")
#plot!(ξ, sol(t), label = "Discretized")
plot!(x, W * sol(t), label = "W * Discretized")
plot!(x, sol_allbar(t), label = "Discretized-then-filtered")
ylims!(minimum(uₕ), maximum(uₕ))


text:
using Revise
using NNFEM
using PyPlot
using LinearAlgebra
using ADCME
using ADCMEKit
using DelimitedFiles
using ProgressMeter
using PyCall
using MAT
using JLD2

NT = 50
Δt = 2.0/NT

FILE = splitdir(pathof(NNFEM))[1]*"/../deps/Data/"
node, elem = meshread(FILE*"twoholes.msh")
node *= 100
xmax = maximum(node[:,1])
ymax = maximum(node[:,2])

elements = []
prop = Dict("name"=> "PlaneStressIncompressibleRivlinSaunders", "rho"=> 1.522,  "C1"=>0.162, "C2"=>5.9e-3)

for i = 1:size(elem,1)
    nodes = node[elem[i,:], :]
    elnodes = elem[i,:]
    ngp = 3
    push!(elements, FiniteStrainContinuum(nodes, elnodes, prop, ngp))
end

Edge_Traction_Data = Array{Int64}[]
for i = 1:length(elements)
    elem = elements[i]
    for k = 1:4
        if elem.coords[k,2]<1e-5 && elem.coords[k+1>4 ? 1 : k+1,2]<1e-5
            push!(Edge_Traction_Data, Int64[i, k, 1])
        end
    end
end
Edge_Traction_Data = hcat(Edge_Traction_Data...)'|>Array

EBC = zeros(Int64, size(node, 1), 2)
FBC = zeros(Int64, size(node, 1), 2)
g = zeros(size(node, 1), 2)
f = zeros(size(node, 1), 2)
for i = 1:size(node, 1)
    if node[i,1]>xmax-1e-5 || node[i,2]>ymax-1e-5
        EBC[i,:] .= -1
    end
end

domain = Domain(node, elements, 2, EBC, g, FBC, f, Edge_Traction_Data)

# Set initial condition: at reset
Dstate = zeros(domain.neqs)
state = zeros(domain.neqs)
velo = zeros(domain.neqs)
acce = zeros(domain.neqs)
EBC_func = nothing
FBC_func = nothing
Body_func = nothing
# Construct Edge_func
function Edge_func_hyperelasticity(x, y, t, idx)
  return [zeros(length(x)) 0.01*ones(length(x)) * sin(2π*t)]
end

globaldata = GlobalData(state, Dstate, velo, acce, domain.neqs, EBC_func, FBC_func,Body_func, Edge_func_hyperelasticity)

assembleMassMatrix!(globaldata, domain)
SolverInitial!(Δt, globaldata, domain)

updateStates!(domain, globaldata)
ω = EigenMode(Δt, globaldata, domain)
@show "stable time step is ", 0.8 * 2/ω, " current time step is ", Δt
@showprogress for i = 1:NT
    global globaldata, domain = GeneralizedAlphaSolverStep(globaldata, domain, Δt)
end

p = visualize_von_mises_stress(domain)
saveanim(p, raw"C:\Users\kaila\Desktop\ADCMEImages\NNFEM\visualize_von_mises_stress.gif")

visualize_von_mises_stress(domain, 10)
savefig(raw"C:\Users\kaila\Desktop\ADCMEImages\NNFEM\visualize_von_mises_stress_50.png")

p = visualize_displacement(domain, scale_factor = 20.0)
saveanim(p, raw"C:\Users\kaila\Desktop\ADCMEImages\NNFEM\visualize_displacement.gif")

visualize_mesh(domain.nodes, getElems(domain))
visualize_mesh(domain)
savefig(raw"C:\Users\kaila\Desktop\ADCMEImages\NNFEM\visualize_mesh.png")

visualize_boundary(domain)
savefig(raw"C:\Users\kaila\Desktop\ADCMEImages\NNFEM\visualize_boundary.png")

d_ = hcat(domain.history["state"]...)'|>Array
visualize_total_deformation_on_scoped_body(d_, domain; scale_factor=20)
saveanim(p, raw"C:\Users\kaila\Desktop\ADCMEImages\NNFEM\visualize_total_deformation_on_scoped_body.gif")

visualize_von_mises_stress_on_scoped_body(d_, domain; scale_factor=20)
saveanim(p, raw"C:\Users\kaila\Desktop\ADCMEImages\NNFEM\visualize_von_mises_stress_on_scoped_body.gif")


text:
using DrWatson
@quickactivate :FractalDimension # uses DynamicalSystems, PyPlot
include(srcdir("style.jl"))
using3D()

function mackieglass(xnew, x, p, n)
    τ = length(x) # total delay, in units of δt
    β, γ, n, δt = p
    # x[1] is latest x
    xnew[1] = δt*(β*x[L]/(1 + x[L]^n) - γ*x[1]) + x[1]
    xnew[2:τ] .= view(x, 1:τ-1)
end

L = 300
p0 = [2, 1.0, 9.65, 0.01]
u0 = 0.1rand(L)

ds = DiscreteDynamicalSystem(mackieglass, u0, p0)
tr = trajectory(ds, 200000; Ttr = 10000, dt = 10)
x, y = tr[:, 1], tr[:, L]
figure()
plot(x,y, marker = "o", ls = "None", alpha = 0.5)
#
# t = estimate_delay(x, "mi_min")
# E = delay_afnn(x, t, 1:10)
# plot(1:10, E)

max_stars_repo_name: Datseris/FractalDimension

max_issues_repo_name: Datseris/FractalDimension

max_forks_repo_name: Datseris/FractalDimension


text:
#=
Core processing file for this project containing general context agnostic tools:
* data conversion
* aggregation
* convenience tools for DateTime, e.g. ensuring yearly averages
=#

# TODO: Add examples to all docstrings

shortname(attrib) = attrib["CF_name"]*" [$(attrib["units"])]"
longname(attrib) = attrib["long_name"]*" [$(attrib["units"])]"

export shortname, longname

export totalaverage
totalaverage(x) = spacemean(timemean(x))

using Distances, StaticArrays
export findnearest, SVector, split_at_longitude

"""
    split_at_longitude(A, λ) → A1, A2
Split `A` at given longitude, taking into account the periodic nature of the longitude.
It is expected that input `A` covers all longitudes from 0 to 360.
"""
function split_at_longitude(A, λ)
    λi = findnearest(λ, dims(A, Lon))
    # here 179 explicitly assumes longitudinal spacing of 1 degree.
    # TODO: change it to use internal spacing, and do 180 - δλ
    i1 = mod1.(λi:λi+179, size(A, Lon))
    i2 = setdiff(1:size(A, Lon), i1)
    return A[Lon(i1)], A[Lon(i2)]
end

"""
    findnearest(val, A)
Return the index of `A` which has value nearest to `val`.
"""
function findnearest(val, A)
    i = 1
    d = evaluate(Euclidean(), val, A[i])
    @inbounds for j in 1:length(A)
        dd = evaluate(Euclidean(), val, A[j])
        if dd < d
            i = j
            d = dd
        end
    end
    return i
end

export timeseries

_correct(::Coord, r) = r
_correct(::Coord, r::Real) = SVector(0, r)
_correct(::Lat, r) = r[2]
_correct(::Lat, r::Real) = r

"""
    timeseries(F, g)
Return the timeseries of field `F` at point `g` (either latitude or (lon, lat)
vector). Works with both global `F` as well as zonally-averaged.
"""
function timeseries(F, r)
    if hasdim(F, Lon)
        p = _correct(Coord(), r)
        i = findnearest(p[1], dims(F, 1))
        j = findnearest(p[2], dims(F, 2))
        given = SVector(dims(F, 1)[i], dims(F, 2)[j])
        a = Array(F[i, j, :])
    else
        p = _correct(dims(F, 1), r) #assumes spatial dim is first dim
        i = findnearest(p, dims(F, 1))
        given = dims(F, 1)[i]
        a = Array(F[i, :])
    end
    return a, given
end

using Dates: TimeType

"""
    monthly_insolation(t::TimeType, args...)
Average `insolation(τ, ...)` for `τ` in the month of given time.
"""
function monthly_insolation(t::TimeType, args...)
    d = monthspan(t)
    mean(insolation(τ, args...) for τ in d)
end

export linreg
"""
    linreg(x, y) -> s, o, fit, x
Return best linear fit so that y ≈ s*x + o.
If `x` is a date vector, `1:length(x)` is used instead.
"""
linreg(x::AbstractArray{<:Dates.AbstractDateTime}, y) = linreg(1:length(x), y)
function linreg(x, y)
    X = ones(length(x), 2)
    X[:, 1] .= x
    c = X\y
    s, o = c
    f = s .* x .+ o
    return s, o, f, x
end

max_stars_repo_name: Datseris/EarthAlbedoSymmetry

max_issues_repo_name: Datseris/EarthAlbedoSymmetry

max_forks_repo_name: Datseris/EarthAlbedoSymmetry


text:
# -------------------------------------------------------
#                        Generic
# -------------------------------------------------------
function slice_matrix(mc::AbstractDQMC, slice::Int, power::Float64=1.)
  a = mc.a

  res = Matrix{geltype(mc)}(I, mc.p.flv*mc.l.sites, mc.p.flv*mc.l.sites)

  @mytimeit a.to "slice_matrix" begin
  if power > 0
    multiply_B_left!(mc, slice, res)
  else
    multiply_B_inv_left!(mc, slice, res)
  end
  end
  return res
end

# -------------------------------------------------------
#  				No checkerboard
# -------------------------------------------------------
# Beff(slice) = exp(−1/2∆τT)exp(−1/2∆τT)exp(−∆τV(slice))
function slice_matrix!(mc::DQMC_CBFalse, slice::Int, power::Float64=1., Bl::AbstractMatrix=mc.s.Bl)
  eT = mc.l.hopping_matrix_exp
  eTinv = mc.l.hopping_matrix_exp_inv
  eV = mc.s.eV
  tmp = mc.s.tmp

  interaction_matrix_exp!(mc, slice, power, eV)

  if power > 0
    mul!(tmp, eT, eV)
    mul!(Bl, eT, tmp)
    # return eT * eT * eV
  else
    mul!(tmp, eTinv, eTinv)
    mul!(Bl, eV, tmp)
    # return eV * eTinv * eTinv
  end
  nothing # slice_matrix now in Bl
end

function multiply_B_left!(mc::DQMC_CBFalse, slice::Int, M::AbstractMatrix)
  @mytimeit mc.a.to "multiply_B (combined)" begin
  @mytimeit mc.a.to "multiply_B_left!" begin
  slice_matrix!(mc, slice, 1., mc.s.Bl)
  mul!(mc.s.tmp, mc.s.Bl, M)
  M .= mc.s.tmp
  end #timeit
  end #timeit
  nothing
end
function multiply_B_right!(mc::DQMC_CBFalse, slice::Int, M::AbstractMatrix)
  @mytimeit mc.a.to "multiply_B (combined)" begin
  @mytimeit mc.a.to "multiply_B_right!" begin
  slice_matrix!(mc, slice, 1., mc.s.Bl)
  mul!(mc.s.tmp, M, mc.s.Bl)
  M .= mc.s.tmp
	end #timeit
  end #timeit
  nothing
end
function multiply_B_inv_right!(mc::DQMC_CBFalse, slice::Int, M::AbstractMatrix)
  @mytimeit mc.a.to "multiply_B (combined)" begin
  @mytimeit mc.a.to "multiply_B_inv_right!" begin
  slice_matrix!(mc, slice, -1., mc.s.Bl)
  mul!(mc.s.tmp, M, mc.s.Bl)
  M .= mc.s.tmp
	end #timeit
  end #timeit
  nothing
end
function multiply_B_inv_left!(mc::DQMC_CBFalse, slice::Int, M::AbstractMatrix)
  @mytimeit mc.a.to "multiply_B (combined)" begin
  @mytimeit mc.a.to "multiply_B_inv_left!" begin
  slice_matrix!(mc, slice, -1., mc.s.Bl)
  mul!(mc.s.tmp, mc.s.Bl, M)
  M .= mc.s.tmp
	end #timeit
  end #timeit
  nothing
end
function multiply_daggered_B_left!(mc::DQMC_CBFalse, slice::Int, M::AbstractMatrix)
  @mytimeit mc.a.to "multiply_B (combined)" begin
  @mytimeit mc.a.to "multiply_daggered_B_left!" begin
  slice_matrix!(mc, slice, 1., mc.s.Bl)
  mul!(mc.s.tmp, adjoint(mc.s.Bl), M) # ctranspose
  M .= mc.s.tmp
  end #timeit
  end #timeit
  nothing
end

# -------------------------------------------------------
#  				Assaad checkerboard
# -------------------------------------------------------
function multiply_B_left!(mc::AbstractDQMC{CBAssaad}, slice::Int, M::AbstractMatrix{T}) where T<:Number
  l = mc.l
  s = mc.s
  a = mc.a
  @mytimeit a.to "multiply_B (combined)" begin
  @mytimeit a.to "multiply_B_left!" begin

  @mytimeit a.to "Bleft_construct_eV" begin
  interaction_matrix_exp!(mc,slice,1.,s.eV)
  end
  @mytimeit a.to "Bleft_mult_eV" begin
  mul!(s.tmp, s.eV, M)
  M .= s.tmp
  end
  mul!(s.tmp, l.chkr_mu, M)
  M .= s.tmp

  @mytimeit a.to "Bleft_mult_eT" begin
  mul!(s.tmp, l.chkr_hop_half[2], M)
  M .= s.tmp
  mul!(s.tmp, l.chkr_hop[1], M)
  M .= s.tmp
  mul!(s.tmp, l.chkr_hop_half[2], M)
  M .= s.tmp
  end
  end
  end #timeit
  nothing
end

function multiply_B_right!(mc::AbstractDQMC{CBAssaad}, slice::Int, M::AbstractMatrix{T}) where T<:Number
  l = mc.l
  s = mc.s
  a = mc.a
  @mytimeit a.to "multiply_B (combined)" begin
  @mytimeit a.to "multiply_B_right!" begin

  interaction_matrix_exp!(mc,slice,1.,s.eV)
  mul!(s.tmp, M, l.chkr_hop_half[2])
  M .= s.tmp
  mul!(s.tmp, M, l.chkr_hop[1])
  M .= s.tmp
  mul!(s.tmp, M, l.chkr_hop_half[2])
  M .= s.tmp

  mul!(s.tmp, M, l.chkr_mu)
  M .= s.tmp
  mul!(s.tmp, M, s.eV)
  M .= s.tmp
  end
  end #timeit
  nothing
end

function multiply_B_inv_left!(mc::AbstractDQMC{CBAssaad}, slice::Int, M::AbstractMatrix{T}) where T<:Number
  l = mc.l
  s = mc.s
  a = mc.a
  @mytimeit a.to "multiply_B (combined)" begin
  @mytimeit a.to "multiply_B_inv_left!" begin

  interaction_matrix_exp!(mc, slice, -1., s.eV)
  mul!(s.tmp, l.chkr_hop_half_inv[2], M)
  M .= s.tmp
  mul!(s.tmp, l.chkr_hop_inv[1], M)
  M .= s.tmp
  mul!(s.tmp, l.chkr_hop_half_inv[2], M)
  M .= s.tmp

  mul!(s.tmp, l.chkr_mu_inv, M)
  M .= s.tmp
  mul!(s.tmp, s.eV, M)
  M .= s.tmp
  end
  end #timeit
  nothing
end

function multiply_B_inv_right!(mc::AbstractDQMC{CBAssaad}, slice::Int, M::AbstractMatrix{T}) where T<:Number
  l = mc.l
  s = mc.s
  a = mc.a
  @mytimeit a.to "multiply_B (combined)" begin
  @mytimeit a.to "multiply_B_inv_right!" begin

  interaction_matrix_exp!(mc, slice, -1., s.eV)
  mul!(s.tmp, M, s.eV)
  M .= s.tmp
  mul!(s.tmp, M, l.chkr_mu_inv)
  M .= s.tmp

  mul!(s.tmp, M, l.chkr_hop_half_inv[2])
  M .= s.tmp
  mul!(s.tmp, M, l.chkr_hop_inv[1])
  M .= s.tmp
  mul!(s.tmp, M, l.chkr_hop_half_inv[2])
  M .= s.tmp
  end
  end #timeit
  nothing
end

function multiply_daggered_B_left!(mc::AbstractDQMC{CBAssaad}, slice::Int, M::AbstractMatrix{T}) where T<:Number
  l = mc.l
  s = mc.s
  a = mc.a
  @mytimeit a.to "multiply_B (combined)" begin
  @mytimeit a.to "multiply_daggered_B_left!" begin

  interaction_matrix_exp!(mc, slice, 1., s.eV)
  mul!(s.tmp, l.chkr_hop_half_dagger[2], M)
  M .= s.tmp
  mul!(s.tmp, l.chkr_hop_dagger[1], M)
  M .= s.tmp
  mul!(s.tmp, l.chkr_hop_half_dagger[2], M)
  M .= s.tmp

  # s.eV == adjoint(s.eV) and l.chkr_mu == adjoint(s.chkr_mu)
  mul!(s.tmp, l.chkr_mu, M)
  M .= s.tmp
  mul!(s.tmp, s.eV, M)
  M .= s.tmp
  end
  end #timeit
  nothing
end

# -------------------------------------------------------
#  				Generic checkerboard
# -------------------------------------------------------
function multiply_B_left!(mc::AbstractDQMC{CBGeneric}, slice::Int, M::AbstractMatrix{T}) where T<:Number
  l = mc.l
  p = mc.p
  s = mc.s
  a = mc.a
  @timeit a.to "multiply_B (combined)" begin
  @timeit a.to "multiply_B_left!" begin

  # @timeit a.to "Bleft_construct_eV" begin
  interaction_matrix_exp!(mc,slice,1.,s.eV)
  # end
  # @timeit a.to "Bleft_mult_eV" begin
  mul!(s.tmp, s.eV, M)
  M .= s.tmp
  # end
  mul!(s.tmp, l.chkr_mu, M)
  M .= s.tmp

  # @timeit a.to "Bleft_mult_eT" begin
  @inbounds @views begin
    for i in reverse(1:l.n_folded)
      mul!(s.tmp, l.chkr_hop_half_folded_rev[i], M)
      M .= s.tmp
    end
    mul!(s.tmp, l.chkr_hop[1], M)
    M .= s.tmp
    for i in 1:l.n_folded
      mul!(s.tmp, l.chkr_hop_half_folded[i], M)
      M .= s.tmp
    end
  end
  # end
  end
  end #timeit
  nothing
end

function multiply_B_right!(mc::AbstractDQMC{CBGeneric}, slice::Int, M::AbstractMatrix{T}) where T<:Number
  l = mc.l
  p = mc.p
  s = mc.s
  a = mc.a
  @mytimeit a.to "multiply_B (combined)" begin
  @mytimeit a.to "multiply_B_right!" begin

  @inbounds @views begin
    for i in reverse(1:l.n_folded)
      mul!(s.tmp, M, l.chkr_hop_half_folded[i])
      M .= s.tmp
    end
    mul!(s.tmp, M, l.chkr_hop[1])
    M .= s.tmp
    for i in 1:l.n_folded
      mul!(s.tmp, M, l.chkr_hop_half_folded_rev[i])
      M .= s.tmp
    end
  end

  interaction_matrix_exp!(mc,slice,1.,s.eV)
  mul!(s.tmp, M, l.chkr_mu)
  M .= s.tmp
  mul!(s.tmp, M, s.eV)
  M .= s.tmp
  end
  end #timeit
  nothing
end

function multiply_B_inv_left!(mc::AbstractDQMC{CBGeneric}, slice::Int, M::AbstractMatrix{T}) where T<:Number
  l = mc.l
  p = mc.p
  s = mc.s
  a = mc.a
  @mytimeit a.to "multiply_B (combined)" begin
  @mytimeit a.to "multiply_B_inv_left!" begin

  @inbounds @views begin
    for i in reverse(1:l.n_folded)
      mul!(s.tmp, l.chkr_hop_half_inv_folded_rev[i], M)
      M .= s.tmp
    end
    mul!(s.tmp, l.chkr_hop_inv[1], M)
    M .= s.tmp
    for i in 1:l.n_folded
      mul!(s.tmp, l.chkr_hop_half_inv_folded[i], M)
      M .= s.tmp
    end
  end

  interaction_matrix_exp!(mc, slice, -1., s.eV)
  mul!(s.tmp, l.chkr_mu_inv, M)
  M .= s.tmp
  mul!(s.tmp, s.eV, M)
  M .= s.tmp
  end
  end #timeit
  nothing
end

function multiply_B_inv_right!(mc::AbstractDQMC{CBGeneric}, slice::Int, M::AbstractMatrix{T}) where T<:Number
  l = mc.l
  p = mc.p
  s = mc.s
  a = mc.a
  @mytimeit a.to "multiply_B (combined)" begin
  @mytimeit a.to "multiply_B_inv_right!" begin

  interaction_matrix_exp!(mc, slice, -1., s.eV)
  mul!(s.tmp, M, s.eV)
  M .= s.tmp
  mul!(s.tmp, M, l.chkr_mu_inv)
  M .= s.tmp

  @inbounds @views begin
    for i in reverse(1:l.n_folded)
      mul!(s.tmp, M, l.chkr_hop_half_inv_folded[i])
      M .= s.tmp
    end
    mul!(s.tmp, M, l.chkr_hop_inv[1])
    M .= s.tmp
    for i in 1:l.n_folded
      mul!(s.tmp, M, l.chkr_hop_half_inv_folded_rev[i])
      M .= s.tmp
    end
  end
  end
  end #timeit
  nothing
end

function multiply_daggered_B_left!(mc::AbstractDQMC{CBGeneric}, slice::Int, M::AbstractMatrix{T}) where T<:Number
  l = mc.l
  p = mc.p
  s = mc.s
  a = mc.a
  @mytimeit a.to "multiply_B (combined)" begin
  @mytimeit a.to "multiply_daggered_B_left!" begin

  @inbounds @views begin
    for i in reverse(1:l.n_folded)
      mul!(s.tmp, l.chkr_hop_half_dagger_folded_rev[i], M)
      M .= s.tmp
    end
    mul!(s.tmp, l.chkr_hop_dagger[1], M)
    M .= s.tmp
    for i in 1:l.n_folded
      mul!(s.tmp, l.chkr_hop_half_dagger_folded[i], M)
      M .= s.tmp
    end
  end

  interaction_matrix_exp!(mc, slice, 1., s.eV)
  # s.eV == adjoint(s.eV) and l.chkr_mu == adjoint(s.chkr_mu)
  mul!(s.tmp, l.chkr_mu, M)
  M .= s.tmp
  mul!(s.tmp, s.eV, M)
  M .= s.tmp
  end
  end #timeit
  nothing
end

max_stars_repo_name: crstnbr/dqmc

max_issues_repo_name: crstnbr/dqmc

max_forks_repo_name: crstnbr/dqmc


text:
#!/usr/bin/env julia

using Merly
using JSON

server = Merly.app()

@route POST "/" begin
    "Gotcha! Here's your data back: $body"
end

server.start("localhost", 80)

max_stars_repo_name: espectalll/gtk-rest-examples

max_issues_repo_name: espectalll/gtk-rest-examples

max_forks_repo_name: espectalll/gtk-rest-examples


text:
### conversions.jl
#
# Copyright (C) 2016, 2017 Mosè Giordano.
#
# Maintainer: Mosè Giordano <mose AT gnu DOT org>
# Keywords: uncertainty, error propagation, physics
#
# This file is a part of Measurements.jl.
#
# License is MIT "Expat".
#
### Commentary:
#
# This file defines conversion and promotion rules involving Measurement
# objects.
#
### Code:

Base.convert(::Type{Measurement{T}}, a::Irrational) where {T<:AbstractFloat} =
    measurement(T(a))
Base.convert(::Type{Measurement{T}}, a::Rational{<:Integer}) where {T<:AbstractFloat} =
    measurement(T(a))
Base.convert(::Type{Measurement{T}}, a::Real) where {T<:AbstractFloat} = measurement(T(a))

Base.convert(::Type{Measurement{T}}, a::Measurement{T}) where {T<:AbstractFloat} = a
function Base.convert(::Type{Measurement{T}},
                      a::Measurement{<:AbstractFloat}) where {T<:AbstractFloat}
    newder = empty_der2(zero(T))
    for tag in keys(a.der)
        newder = Derivatives(newder, (T(tag[1]), T(tag[2]), tag[3])=>T(a.der[tag]))
    end
    Measurement(T(a.val), T(a.err), a.tag, newder)
end

Base.convert(::Type{Measurement}, a::Measurement) = a
Base.convert(::Type{Measurement}, a::Rational{<:Integer}) = measurement(a)
Base.convert(::Type{Measurement}, a::Real) = measurement(a)
Base.convert(::Type{Signed}, a::Measurement) = convert(Signed, a.val)

Base.promote_rule(::Type{Measurement{T}}, ::Type{S}) where {T<:AbstractFloat, S<:Real} =
    Measurement{promote_type(T, S)}
Base.promote_rule(::Type{Measurement{T}},
                  ::Type{Measurement{S}}) where {T<:AbstractFloat, S<:AbstractFloat} =
    Measurement{promote_type(T, S)}


text:
# Check whether the test is being run as standalone
if ~@isdefined Test
    using Test
end

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

@testset "hierarchy matrix" begin
    # Makes sure the hierarchy matrix generator is working as expected for a select set of orbital systems.
    # Fully nested Keplerians
    ϵ_nested::Matrix{Float64} = [-1  1  0  0;
                                 -1 -1  1  0;
                                 -1 -1 -1  1;
                                 -1 -1 -1 -1;]
    @test hierarchy([4,1,1,1]) == ϵ_nested

    # Double Binary
    ϵ_double::Matrix{Float64} = [-1  1  0  0;
                                  0  0 -1  1;
                                 -1 -1  1  1;
                                 -1 -1 -1 -1;]
    @test hierarchy([4,2,1]) == ϵ_double

    # Double binary + top level body
    ϵ_asym::Matrix{Float64} = [-1  1  0  0  0;
                                0  0 -1  1  0;
                               -1 -1  1  1  0;
                               -1 -1 -1 -1  1;
                               -1 -1 -1 -1 -1;]
    @test hierarchy([5,2,1,1]) == ϵ_asym

    # Test of specific ambiguous initialization case of a 'symmetric fully-nested' (top) and the double binary + tertiary binary (bottom)
    ϵ_sym_full::Matrix{Float64} = [-1  1  0  0  0  0;
                                    0  0  0 -1  1  0;
                                   -1 -1  1  0  0  0;
                                    0  0  0 -1 -1  1;
                                   -1 -1 -1  1  1  1;
                                   -1 -1 -1 -1 -1 -1;]

    ϵ_asym_multi::Matrix{Float64}=[-1  1  0  0  0  0;
                                    0  0 -1  1  0  0;
                                   -1 -1  1  1  0  0;
                                    0  0  0  0 -1  1;
                                   -1 -1 -1 -1  1  1;
                                   -1 -1 -1 -1 -1 -1;]
    h_sym_full,h_asym_multi = hierarchy.([[6,2,2,1,2],[6,2,2,1]])

    @test ϵ_sym_full == h_sym_full
    @test ϵ_asym_multi == h_asym_multi
end


text:
using KernelFunctions
using MLDataUtils
using Zygote
using Flux
using Distributions, LinearAlgebra
using Plots

N = 100 #Number of samples
μ = randn(2, 2) # Random Centers
xgrid = range(-3, 3; length=100) # Create a grid
Xgrid = hcat(collect.(Iterators.product(xgrid, xgrid))...)' #Combine into a 2D grid
y = rand([-1, 1], N) # Select randomly between the two classes
X = zeros(N, 2)
X[y .== 1, :] = rand(MvNormal(μ[:, 1], I), count(y .== 1))' #Attribute samples from class 1
X[y .== -1, :] = rand(MvNormal(μ[:, 2], I), count(y .== -1))' # Attribute samples from class 2

k = SqExponentialKernel(2.0) # Create kernel function
function f(x, k, λ)
    return kernelmatrix(k, x, X; obsdim=1) *
           inv(kernelmatrix(k, X; obsdim=1) + exp(λ[1]) * I) *
           y
end # Optimal prediction f
svmloss(y, ŷ) = ŷ -> sum(maximum.(0.0, 1 - y * ŷ)) - λ * norm(ŷ)(f(X, k, λ)) # Total svm loss with regularisation
pred = f(Xgrid, k, λ) #Compute prediction on a grid
contourf(xgrid, xgrid, pred)
scatter!(eachcol(X)...; color=y, lab="data")


text:
"""
    https://github.com/berkeleydeeprlcourse/homework/blob/60b3ef08c2dca3961efb63b32683bb588571f226/sp17_hw/hw3/dqn_utils.py#L49
"""

linear_interpolation(l, r, alpha) = l + alpha * (r - l)

struct PiecewiseSchedule
    endpoints
    interpolation
end

PiecewiseSchedule(endpoints) = PiecewiseSchedule(endpoints, linear_interpolation)

function value(schedule::PiecewiseSchedule, t)
    for ((l_t, l),(r_t,r)) in zip(schedule.endpoints[1:end-1],schedule.endpoints[2:end])
        if l_t <= t && t < r_t
            alpha = (t - l_t) / (r_t - l_t)
            return schedule.interpolation(l, r, alpha)
        end
    end
    return schedule.endpoints[end][2]
end


text:
"""
Definition of continuation problem.

Numerical continuation algorithms find curves in ``\\mathbb R^{N}``
implicitly defined by

```math
H(u) = 0
```

where ``H: \\mathbb R^{N} \\to \\mathbb R^{N-1}``.

A continuation problem type (a subtype of
`AbstractContinuationProblem`) defines problems for such algorithms to
solve by providing how to compute:

* ``H(u)``: [`residual`](@ref), [`residual!`](@ref)
* its derivative ``\\partial H /  \\partial u``: [`residual_jacobian!`](@ref)
* an initial guess ``u_0``: [`get_u0`](@ref)
* and computation cache: [`get_prob_cache`](@ref).
"""
abstract type AbstractContinuationProblem end

"""
    AbstractProblemCache{P <: AbstractContinuationProblem}

Cache for computing ``H`` and its Jacobian.
"""
abstract type AbstractProblemCache{P <: AbstractContinuationProblem} end

"""
    residual(u, cache) ↦ H

Compute ``H(u) \\in \\mathbb R^{N - 1}`` (aka in-place computation).
The definition of ``H`` is specified by `cache`.

The name `residual` of the function is came from the problem we are to
solve: i.e., find the set of ``u`` such that ``H(u) = 0``.  Thus, the
vector returned by `residual(u, cache)` is considered to be a
residual.

# Arguments
- `u::AbstractVector` (size: `(N,)`)
- `cache::AbstractProblemCache`
"""
function residual end  # TODO: remove it (make it internal)

"""
    residual!(H, u, cache) ↦ H

Compute ``H(u) \\in \\mathbb R^{N - 1}`` and store the result in `H`
(aka out-of-place computation).
The definition of ``H`` is specified by `cache`.

See also: [`residual`](@ref)

# Arguments
- `H::AbstractVector` (size: `(N - 1,)`)
- `u::AbstractVector` (size: `(N,)`)
- `cache::AbstractProblemCache`
"""
function residual! end

"""
    residual_jacobian!(H, J, u, cache) ↦ (H, J)

Compute ``H(u)`` and its Jacobian ``\\partial H /  \\partial u``.

# Arguments
- `H::AbstractVector` (size: `(N - 1,)`) ``= H(u)``
- `J::AbstractMatrix` (size: `(N - 1, N)`) ``= \\partial H /  \\partial u``
- `u::AbstractVector` (size: `(N,)`)
- `cache::AbstractProblemCache`
"""
function residual_jacobian! end

"""
    isindomain(u, cache) :: Bool

# Arguments
- `u::AbstractVector` (size: `(N,)`)
- `cache::AbstractProblemCache`
"""
function isindomain end

"""
    get_prob_cache(prob::AbstractContinuationProblem) :: AbstractProblemCache
"""
function get_prob_cache end
# TODO: Rename AbstractProblemCache to ProblemCache and use it instead
# of get_prob_cache.

"""
    get_u0(prob::AbstractContinuationProblem) ↦ u0
"""
function get_u0 end

abstract type AbstractContinuationCache{PC <: AbstractProblemCache} end


text:
function splitInt(i::Int)
    a = []
    while 0 < i
        pushfirst!(a, i % 10)
        i = i ÷ 10
    end
    return a
end

lower = parse(Int, ARGS[1])
upper = parse(Int, ARGS[2])

# Rule 1: must contain at least one consecutive bit string
containsConsecutive(a) = (a[1] == a[2]) || (a[2] == a[3]) ||
    (a[3] == a[4]) || (a[4] == a[5]) || (a[5] == a[6])

# Rule 2: digits in the string must be in non-decreasing order
increasing(a) = a[1] <= a[2] <= a[3] <= a[4] <= a[5] <= a[6]

# Rule 3 (part 2): only pairs count, not triples or longer matches.
function containsPair(i::Int)
    if i == 0 return false end
    current = i % 10
    matches = 0
    while (i ÷ (10^matches)) % 10 == current
        matches = matches + 1
    end
    return matches == 2 || containsPair(i ÷ (10^matches))
end

function countCombinations(lower::Int, upper::Int)
    combinations = Set()
    for x in lower:upper
        a = splitInt(x)
        if (containsConsecutive(a) && increasing(a))
            push!(combinations, x)
        end
    end
    return combinations
end

combinations = countCombinations(lower, upper)
println("Part 1: $(length(combinations))") # 925

grouped = filter(containsPair, combinations)
println("Part 2: $(length(grouped))")
# This one turned out to be more challenging for me than it should have.
# I had misinterpreted the specification for part 2. I thought that
# the groups could not be length 3, but that two groups of 2 next to each
# other would be OK somehow. As in, "567777" would be a satisfying assignment.
# I was wrong.
# I'm kind of OK with my algorithm to match pairs. It does not easily
# generalize to arbitrary strings.


text:
const ScalarPotential{T, N, S, AT} = Union{ElectricPotential{T, N, S, AT}, WeightingPotential{T, N, S, AT}, PointTypes{T, N, S, AT}, EffectiveChargeDensity{T, N, S, AT}, DielectricDistribution{T, N, S, AT}}

ScalarPotential(::ElectricPotential, data, grid) = ElectricPotential(data, grid)
ScalarPotential(::WeightingPotential, data, grid) = WeightingPotential(data, grid)
ScalarPotential(::PointTypes, data, grid) = PointTypes(data, grid)
ScalarPotential(::EffectiveChargeDensity, data, grid) = EffectiveChargeDensity(data, grid)

get_axes_type(::Type{Tuple{AT1, AT2, AT3}}) where {AT1, AT2, AT3} = (AT1, AT2, AT3)

function getindex(ep::P, grid::Grid{T, N, S}) where {T, N, S, P <: ScalarPotential{T, N, S}}
    gridsize::Tuple = size(grid)
    data::Array{T, N} = zeros(T, gridsize)
    ep_itp::Interpolations.Extrapolation{T, N} = interpolated_scalarfield(ep)
    PT = (S == Cylindrical ? CylindricalPoint : CartesianPoint)
    for i1 in eachindex(grid[1])
        for i2 in eachindex(grid[2])
            for i3 in eachindex(grid[3])
                data[i1, i2, i3] = get_interpolation(ep_itp, PT{T}(grid[i1, i2, i3]), S)
            end
        end
    end
    return P(data, grid)
end

function get_2π_potential(sp::ScalarPotential{T, 3, Cylindrical}, axφ::DiscreteAxis{AT, :periodic, :periodic}, int::Interval{:closed, :open, AT}) where {T, AT}
    @assert int.right != 0 "Right boundary of φ interval is not allowed to be 0"
    l::Int = length( axφ )
    Δφ::AT = width(int)
    new_int::Interval{:closed, :open, AT} = Interval{:closed, :open, AT}(0, 2π)
    n::Int = Int(round(T(2π) / Δφ, sigdigits = 6))
    new_ticks::Vector{AT} = Vector{AT}(undef, l * n)
    new_pot = Array{eltype(sp.data), 3}(undef, size(sp, 1), l * n, size(sp, 3))
    for idx_n in 1:n
        for il in 1:l
            idx::Int = il + (idx_n - 1) * l
            new_ticks[idx] = (idx_n - 1) * Δφ + axφ.ticks[il]
            new_pot[:, idx, :] = sp[:, il, :]
        end
    end
    new_axφ = DiscreteAxis{AT, :periodic, :periodic}( new_int, new_ticks )
    new_axes = (sp.grid[1], new_axφ, sp.grid[3])
    new_grid = Grid{AT, 3, Cylindrical, typeof(new_axes)}( new_axes )
    ScalarPotential(sp, new_pot, new_grid)
end

function get_2π_potential(sp::ScalarPotential{T, 3, Cylindrical}, axφ::DiscreteAxis{AT, :reflecting, :reflecting}, int::Interval{:closed, :closed, AT}) where {T, AT}
    @assert int.right != 0 "Right boundary of φ interval is not allowed to be 0"
    l::Int = 2 * length( axφ ) - 2
    Δφ::AT = width(int)
    new_int::Interval{:closed, :open, AT} = Interval{:closed, :open, AT}(int.left, int.right + Δφ)
    new_ticks::Vector{AT} = Vector{AT}(undef, l)
    new_ticks[1:length( axφ )] = collect(axφ.ticks)
    new_pot = Array{eltype(sp.data), 3}(undef, size(sp, 1), l, size(sp, 3))
    new_pot[:, 1:length( axφ ), :] = sp.data[:, :, :]
    for i in 1:(length(axφ) - 2)
        new_ticks[length(axφ) + i] = 2 * int.right - axφ.ticks[length(axφ) - i]
        new_pot[:, length(axφ) + i, :] = sp[:, length(axφ) - i, :]
    end
    new_axφ = DiscreteAxis{AT, :periodic, :periodic}( new_int, new_ticks )
    new_axes = (sp.grid[1], new_axφ, sp.grid[3])
    new_grid = Grid{AT, 3, Cylindrical, typeof(new_axes)}( new_axes )
    ScalarPotential(sp, new_pot, new_grid)
end

function extend_2D_to_3D_by_n_points(sp::ScalarPotential{T, 3, Cylindrical}, axφ::DiscreteAxis{AT, :reflecting, :reflecting}, int::Interval{:closed, :closed, AT},
        n_points_in_φ::Int ) where {T, AT}
    new_int::Interval{:closed, :open, AT} = Interval{:closed, :open, AT}(0, 2π)
    new_ticks::Vector{AT} = Vector{AT}(undef, n_points_in_φ)
    new_pot = Array{eltype(sp.data), 3}(undef, size(sp, 1), n_points_in_φ, size(sp, 3))
    Δφ::AT = 2π / n_points_in_φ
    for i in 1:n_points_in_φ
        new_ticks[i] = (i - 1) * Δφ
        new_pot[:, i, :] = sp[:, 1, :]
    end
    new_axφ = DiscreteAxis{AT, :periodic, :periodic}( new_int, new_ticks )
    new_axes = (sp.grid[1], new_axφ, sp.grid[3])
    new_grid = Grid{AT, 3, Cylindrical, typeof(new_axes)}( new_axes )
    ScalarPotential(sp, new_pot, new_grid)
end

function get_2π_potential(sp::ScalarPotential{T, 3, Cylindrical}; n_points_in_φ::Union{Missing, Int} = missing) where {T}
    axφ::DiscreteAxis{T} = sp.grid[2]
    int::Interval = axφ.interval
    if int.right == 0 && length(axφ) == 1 # 2D
        if ismissing(n_points_in_φ)
            error(ArgumentError, ": First Argument is a 2D potential (only 1 point in φ). User has to set the keyword `n_points_in_φ::Int` (e.g. `n_points_in_φ = 18`) in order to get a 3D potential.")
        else
            return extend_2D_to_3D_by_n_points(sp, axφ, int, n_points_in_φ)
        end
    else
        if isclosedset(axφ.interval)
            sp_tmp = get_2π_potential(sp, axφ, int)
            return get_2π_potential(sp_tmp, sp_tmp.grid[2], sp_tmp.grid[2].interval)
        else
            return get_2π_potential(sp, axφ, int)
        end
    end
end


text:
abstract type DynamicSuperRes <: SuperRes end

dim(model :: DynamicSuperRes) = dim(model.static)

function psi(model :: DynamicSuperRes, theta :: Vector{Float64})
        # This function computes the direct problem for a single point theta
        d = dim(model)
        return vcat([psi(model.static, theta[1:d] + t * theta[d+1:end]) for t in model.times]...)
end

function dpsi(model :: DynamicSuperRes, theta :: Vector{Float64})
        # This function computes the gradient of psi
        d = dim(model)
        dx = vcat([dpsi(model.static, theta[1:d] + t * theta[d+1:end]) for t in model.times]...)
        dv = vcat([t * dpsi(model.static, theta[1:d] + t * theta[d+1:end]) for t in model.times]...)
        return [dx dv]
end


text:
module GtkInteract

## TODO:
## * more layout options. (Did formlayout). Do `flex`
## * [DONE] menubars, toolbars
## * GtkMenuButton and toolbar (one isn't there, other isn't working...)
## * [DONE] use formlayout in MainWindow

using Gtk
using Reactive
using DataStructures
using Requires
using Colors
import Cairo

## selectively import pieces of `Interact`
import Interact
import Interact: Button,  button
import Interact: ToggleButton, togglebutton
import Interact: Slider, slider
import Interact: Options, dropdown, radiobuttons, togglebuttons
import Interact: Checkbox, checkbox
import Interact: Textbox, textbox
import Interact: Widget, InputWidget
import Interact: widget, signal
import Reactive: foreach, sampleon, value
import Gtk: destroy

## exports (most widgets of `Interact` and the modified `@manipulate` macro)
export slider, button, checkbox, togglebutton, dropdown, radiobuttons, selectlist, textbox, textarea, togglebuttons
export @manipulate
#export button
export buttongroup, cairographic, cairoimagesurface, label, progress
export icon, tooltip, separator
export mainwindow
export foreach, value
export destroy
export eatn

## Add a non-exclusive set of buttons
## Code is basically the Options code of Interact
type VectorOptions{view,T} <: InputWidget{T} # XXX This is a poor name, but it isn't exported XXX
    signal
    label::AbstractString
    values
    options::OrderedDict{AbstractString, T}
end

function VectorOptions{T}(view::Symbol, options::OrderedDict{AbstractString, T};
                          label = "",
                          value=T[],
                          signal=Signal(value))
    VectorOptions{view, T}(signal, label, value, options)
end

function VectorOptions{T}(view::Symbol, options::AbstractArray{T};
                    kwargs...)
    opts = OrderedDict{AbstractString, T}()
    map(v -> opts[string(v)] = v, options)
    VectorOptions(view, opts; kwargs...)
end

function VectorOptions{K, V}(view::Symbol, options::Associative{K, V};
                             kwargs...)
    opts = OrderedDict{AbstractString, V}()
    map(v->opts[string(v[1])] = v[2], options)
    VectorOptions(view, opts; kwargs...)
end

"""

A `buttongroup` is like `togglebuttons` only one can select 0, 1, or more of the items. (Non-exclusive selection.)

"""
buttongroup(opts; kwargs...) = VectorOptions(:ButtonGroup, opts; kwargs...)

"""

The selectlist constructor is just `Interact.select`

"""
selectlist(args...;kwargs...) = Interact.select(args...; kwargs...)

## Output widgets

"""

Output widgets are different from input widgets, in that they do not have signals propogate when they
are changed. Rather, they are for display.

Values are `push!`ed onto the widget to update the display.
"""
abstract OutputWidget <: Widget

Interact.signal(w::OutputWidget) = w.signal
"""
CairoGraphic.

for a plot window

Replace plot via `push!(cg, winston_plot_object)`
"""
type CairoGraphic <: OutputWidget
    width::Int
    height::Int
    signal
    value
    obj
end

function cairographic(;width::Int=480, height::Int=400)
    widget = CairoGraphic(width, height,nothing, nothing, nothing)
    widget.signal = Signal(widget)
    widget
end

# """
# CairoImageSurface.

# For displaying an in-memory image

# Replace the image via push!(cs, img::AbstractMatrix)
# """

type CairoImageSurface{T<:Union{UInt32,RGB24,ARGB32}} <: OutputWidget
    width::Int
    height::Int
    surf::Cairo.CairoSurface{T}
    signal
    value
    obj
end

function cairoimagesurface{T<:Union{UInt32,RGB24,ARGB32}}(buf::Matrix{T}; width::Int=size(buf,1), height::Int=size(buf,2))
    surf = Cairo.CairoImageSurface(buf)
    widget = CairoImageSurface(width, height, surf, nothing, nothing, nothing)
    widget.signal = Signal(widget)
    widget
end

function Base.push!(cs::CairoImageSurface, data::AbstractMatrix)
    copy!(cs.surf.data, data)
    if cs.obj != nothing
        Gtk.draw(cs.obj)
    end
    cs
end

type ImmerseFigure <: OutputWidget
    width::Int
    height::Int
    signal
    value
    obj
    toolbar
    cnv
end

"""

Add area for an `Immerse` graphic

"""
function immersefigure(;width::Int=480, height::Int=400)
    widget = ImmerseFigure(width, height, nothing, nothing, nothing,nothing, nothing)
    widget.signal = Signal(widget)
    widget
end

"""
Textarea for output

Replace text via `push!(obj, value)`
"""
type Textarea{T <: AbstractString} <: OutputWidget
    width::Int
    height::Int
    signal
    value::T
    buffer
    obj
end

function textarea(;width::Int=480, height::Int=400, value::AbstractString="")
    widget = Textarea(width, height, Signal(Any),  value, nothing, nothing)
    widget.signal = Signal(widget)
    widget
end

textarea(value; kwargs...) = textarea(value=value, kwargs...)

type Label <: OutputWidget
    signal
    value::AbstractString
    obj
end
"""
label.

Like text area, but is clearly not editable and allows for PANGO markup.
Replace text via `push!(obj, value)`
"""
function label(;value="")
    widget = Label(Signal(Any), string(value), nothing)
    widget.signal = Signal(widget)
    widget
end
label(lab; kwargs...) = label(value=lab, kwargs...)

type Image <: OutputWidget
    signal
    value::AbstractString
    obj
end

"""

Show an image from a file

"""
function image(;value="")
    widget = Image(Signal(Any), string(value), nothing)
    widget.signal = Signal(widget)
    widget
end
image(fname; kwargs...) = image(value=fname, kwargs...)

type Progress <: OutputWidget
    signal
    value
    range::Range
    obj
end

## Progress creates a progress bar
##
## `push!` values onto it where `value` is within `[first(range), last(range)]`
## function progress(args...)
##     widget = Progress(args...)
##     widget.signal = Signal(widget)
##     widget
## end
function progress(;label="", value=0, range=0:100)
    widget = Progress(nothing, value, range, nothing)
    widget.signal = Signal(widget)
    widget
end

## Decorative widgets
abstract DecorativeWidget <: Widget

## An icon
type Icon <: DecorativeWidget
    stock_id
    tile
    obj
end

"""

Add a static icon to a display or
add an icon to a GUI element, such as a button.

```
icon("window-close", nothing) # creates an image to be added to a GUI
icon("window-close", btn) # adds icon to button
button("button") |> icon("window-close")  # Curried form
```

Idea is `icon(id, button("text"))` will decorate button

A list of [icons](https://docs.google.com/spreadsheets/d/1HavJQRPpMuq-N0GoN1wJR-9KEGXpKy3-NEPpZZkUGJY/pub?output=html).
"""
icon(id::AbstractString, tile) = Icon(id, tile, nothing)
icon(id::AbstractString) = tile -> Icon(id, tile, nothing)

## tooltip
immutable Tooltip <: DecorativeWidget
    text
    tile
end
"""

Add a tooltip to a widget

```
tooltip("message", widget)
widget |> tooltip("message")
```
"""
tooltip(msg::AbstractString, tile::Widget) = Tooltip(msg, tile)
tooltip(msg::AbstractString) = tile -> tooltip(msg, tile)

## Separator
immutable Separator <: DecorativeWidget
    orient
end
"""

A simple separator between text. [Currently implemented in an old school way...]

"""
separator(orient::Symbol=:horizontal) = Separator(orient)
##################################################
### Container(s) and Layout
##
## We have some exports here to trim down. These kinda follow the layout of Escher
## but Gtk is not as rich as HTML
##
## Our main containers are: window, vbox, hbox, tabs
##
## our attributes work on children and are pad, align, grow, width, height

export size,
       width, #minwidth, maxwidth,
       height, #minheight, maxheight,
       vskip, hskip,
       grow, shrink, #flex,
       align, halign, valign,
       padding,
       vbox, hbox,
       tabs,
       formlayout,
       toolbar,
       menu,
       window,
       messagebox, confirmbox, inputbox,
       openfile, savefile, selectdir

## We follow Escher layout here
immutable Length{unit}
    value::Float64
end

abstract Layout <: Widget
abstract LayoutAttribute <: Layout

immutable Size <: LayoutAttribute
    tile
    w_prefix
    h_prefix
    w_px::Length
    h_px::Length
end

Base.size(w::Int, h::Int, lyt) = Size(lyt, "", "", Length{:px}(w), Length{:px}(h))
Base.size(w::Int, h::Int) = lyt -> size(w, h, lyt)

"Width of widget"
function width(prefix::AbstractString, px::Int, lyt::Layout)
    Size(lyt, prefix, "", Length{:px}(px), Length{:px}(-1))
end
width(prefix::AbstractString, px::Int) = lyt -> width(prefix, px, lyt)
width(px::Int, w) = width("", px, w)
width(px::Int) = w -> width(px, w)
#minwidth(w, x...) = width("min", w, x...)
#maxwidth(w, x...) = width("max", w, x...)

"Height of widget"
function height(prefix::AbstractString, px::Int, lyt::Layout)
    Size(lyt, "", prefix, Length{:px}(-1), Length{:px}(px))
end
height(prefix::AbstractString, px::Int) = lyt -> height(prefix, px, lyt)
height(px::Int, w) = height("", px, w)
height(px::Int) = w -> height(px, w)
#minheight(w, x...) = height("min", w, x...)
#maxheight(w, x...) = height("max", w, x...)

## pus int space between tiles
vskip(y) = size(0, y, empty)
vskip() = size(0,0, empty)              # for toolbars
hskip(y) = size(y, 0, empty)

## grow, shrink, flex
immutable Grow <: LayoutAttribute
    tile
    factor
    direction
end

## factor in [0,1]
"""

Have widget expand if space is available

* `grow(widget)` expand in both horizontal and vertical directions
* `grow(:horizontal, widget)` (also `:vertical`). Expand in horizontal direction, no directive for vertical
* `grow(factor, widget)` use factor in `[0,1]` to determine growth. Currently Only 0 and 1 are supported.

"""
grow(factor::Real, tile) = Grow(tile, factor, [:horizontal,:vertical])
grow(factor::Real) = tile -> grow(factor, tile)
grow(tile) = grow(1.0, tile)
grow(direction::Vector{Symbol}, factor::Real, tile) = Grow(tile, factor, direction)
grow(direction::Symbol, factor::Real, tile) = grow([direction;], factor, tile)
grow(direction::Symbol, tile) = grow([direction;], 1, tile)

immutable Shrink <: LayoutAttribute
    tile
    factor
    direction
end
## factor in [0,1]
"""

Have widget not expand when space is available

* `shrink(widget)`
* `shrink(:horizontal, widget)` (also `:vertical`) to not expand in a given direction
* `shrink(factor, widget)` use factor to determine shrinking.

"""
shrink(factor::Real, tile) = Shrink(tile, factor, [:horizontal,:vertical])
shrink(factor::Real) = tile -> shrink(factor, tile)
shrink(tile) = shrink(1.0, tile)
shrink(direction::Vector{Symbol}, factor::Real, tile) = Shrink(tile, factor, direction)
shrink(direction::Symbol, factor::Real, tile) = shrink([direction;], factor, tile)
shrink(direction::Symbol, tile) = shrink([direction;], 1, tile)

## immutable Flex <: LayoutAttribute
##     tile
##     factor
## end
## ## factor in [0,1]
## XXX implement
## Ignore the width (in a hbox) or height (in a vbox) and stretch or shrink a tile to fill / distribute remaining space
## flex(factor::Real, tile) = Flex(tile, factor)
## flex(factor::Real) = tile -> flex(factor, tile)
## flex(tile) = flex(1.0, tile)

## alignment
Alignments =Dict(:fill      => Gtk.GConstants.GtkAlign.FILL,
                 :axisstart => Gtk.GConstants.GtkAlign.START,
                 :start     => Gtk.GConstants.GtkAlign.START,
                 :axisend   => Gtk.GConstants.GtkAlign.END,
                 :end       => Gtk.GConstants.GtkAlign.END,
                 :center    => Gtk.GConstants.GtkAlign.CENTER,
                 :baseline  => Gtk.GConstants.GtkAlign.BASELINE)

immutable Align <: LayoutAttribute
    tile
    halign
    valign
end

"""

Align child if more space is needed than requested.

* `align(h,v,lyt)` where `h` and `v` are in `[:fill, :axisstart, :start, :axisend, :end, :center,:baseline]`

"""
align(h::Symbol, v::Symbol, lyt) = Align(lyt, Alignments[h], Alignments[v])
align(h::Symbol, v::Symbol) = lyt -> align(h, v, lyt)

"""
Align horizonatally. See `align`.
"""
halign(a::Symbol, lyt) = align(a, :fill, lyt)
halign(a::Symbol) = lyt -> halign(a, lyt)

"""
Align vertically. See `align`.
"""
valign(a::Symbol, lyt) = align(:fill, a, lyt)
valign(a::Symbol) = lyt -> valign(a, lyt)

## padding
immutable Pad <: LayoutAttribute
    tile
    sides
    len::Length
end
## what to put for sides?
padcontent(len, tile, sides=[:left, :right, :top, :bottom]) = Pad(tile, sides, Length{:px}(len))

"""

Add padding to a widget.

* `padding(n, widget)` add n pixels of space around each side
* `padding([:left, :right, :top, :bottom], n, widget)` specify which sides with a vector.

"""
padding(len::Int, widget) = padcontent(len, widget)
padding(len::Int) = widget -> padding(len, widget)
padding(sides::AbstractVector, len, tile) =
    padcontent(sides, len, Container(tile))

##################################################
## Boxes
type FlowContainer <: Layout
    obj
    direction
    children
end

"""

Vertical box for packing in children

```
vbox(child1, child2, ...)
```

Use attributes before packing to get spacing:

```
vbox(grow(child1), padding(10, child2))
```

"""
vbox(children...) = FlowContainer(nothing, "vertical", Any[children...])
vbox() = FlowContainer(nothing, "vertical", Any[])
"""

Horizontal box container. See `vbox`.

"""
hbox(children...) = FlowContainer(nothing, "horizontal", Any[children...])
hbox() = FlowContainer(nothing, "horizontal", Any[])

"""

An empty container for spacing purposes

"""
const empty = vbox()

immutable FormLayout <: Layout
    children
end

"""

Form layout container. There are two columns. The children are
arranged in the second. Their `label` property runs down the first.

Example:
```
sl = slider(1:10, label="slider")
rb = radiobuttons(["one", "two", "three"], label="radio buttons")
btn = button("click me")
fl = formlayout(sl, rb, separator(), btn)
fl |> padding(5) |> window(title="formlayout example")
```
"""
formlayout(children...) = FormLayout(Any[children...;])
formlayout() = FormLayout(Any[])
## can push children on if convenient
## fl = formlayout
## push!(fl, button("one"))
Base.push!(lyt::Layout, child) = (push!(lyt.children, child); nothing)
Base.append!(lyt::Layout, children) = (append!(lyt.children, children); nothing)

## Tabs...
immutable Tabs <: Layout
    children
    labels
    initial::Int
end

"""

Use a notebook to organize pages

The tab labels are specified at construction time using "pairs:"

*  `tabs("label"=>tile, "label1"=>tile, ...; selected=1)`

"""
function tabs(tiles...; selected::Int=1)
    labels = Any[label for (label, child) in tiles]
    children = Any[child for (label, child) in tiles]
    Tabs(children, labels, selected)
end
tabs(;selected::Int=1) = Tabs(Any[], Any[], selected)
function Base.push!(lyt::Tabs, child)
    push!(lyt.labels, child[1])
    push!(lyt.children, child[2])
    nothing
end
function Base.append!(lyt::Tabs, children::Pair...)
    append!(lyt.labels, [child[1] for child in children])
    append!(lyt.children, [child[2] for child in children])
    nothing
end

##################################################
## Toolbar
immutable Toolbar <: Layout
    children
end

"""

Toolbar container. Holds buttons, togglebuttons, separators

"""
toolbar(children...) = Toolbar(Any[children...;])
toolbar() = Toolbar(Any[])

## Menubar
immutable Menu <: Layout
    label
    children
end

"""

A menu bar. Children are items (buttons, separators, or other menu items)

Example:
```
btn1 = button("one")
btn2 = button("two")

mb =menu(menu(btn1, separator(), btn2,
menu(btn1, btn2, label="submenu"),
label="File"),
menu(btn1, btn2, label="Edit"))

do_something(args...) = println("hi")
map(do_something, btn1)            # give some action
w = window(mb, grow(label("space")))
```

Note: toggle buttons can be supported when Gtk does.

"""
menu(children...; label="") = Menu(label, Any[children...;])
menu(; label="") = Menu(label, Any[])

## XXX This isn't working, but should be
immutable MenuButton <: Layout
    label
    children
end
"""

A menu button is used like a combobox, but we can put in menu items. It can also be in a toolbar.

"""
menubutton(children...; label::AbstractString="") = MenuButton(label, Any[children...;])

type Window <: Layout
    width::Int
    height::Int
    title
    children
    obj
end

    """
A parentless container, like MainWindow, but less fuss.

Child widgets are packed into a `vbox`.

* `window(child1, child2, ...; title="some title")`

If no children are passed in, `window(kwargs...)` creates a window
object onto which children may be pushed or appended. The window will render when
its `display` method is called.

```
w = window(title="Title");
push!(w, button("button"))
append!(w, [slider(1:10), slider(1:20)])
w
## to destroy
destroy(w)
```

"""
window(children...; width::Int=-1, height::Int=-1, title::AbstractString="") = Window(width, height, title, Any[children...;], nothing)
Base.push!(widget::Window, child) = (push!(widget.children, child); nothing)
Base.append!(widget::Window, children) = (append!(widget.children, children); nothing)

"""

Close a GUI window.

"""
function Gtk.destroy(widget::Window)
    if widget.obj != nothing
        destroy(widget.obj)
        widget.obj = nothing
    end
end
##
## MainWindow
##
## A top-level window
type MainWindow <: Layout
    width::Int
    height::Int
    title::AbstractString
    children
    window                              # main @GtkWindow. Used by destroy()
    out                                 # output widget
    refs::Vector{Signal}                # signals to prevent from garbage collection
end

"""

Mainwindow for manipulate or easy layout of widgets with label.

The main window uses a "form" layout. That is, a two-column
format. The first column to hold labels, the second the controls. The
label comes from the label property of the control specified through
the `label=` keyword argument.

Example

```
sl = slider(1:10, label="slider")
rb = radiobuttons(["one", "two", "three"], label="radio")
w = mainwindow(title="a main window");
append!(w, [sl, rb])
w                                       # when displayed, creates window.  Call destroy(w) to close.
```

"""
function mainwindow(children...;width::Int=300, height::Int=200, title::AbstractString="")
    widget = MainWindow(width, height, title, Any[children...;],
                        nothing, nothing, Signal[]) # window, outputwidget
    widget
end

Base.push!(widget::MainWindow, w) = (push!(widget.children, w); nothing)
Base.append!(widget::MainWindow, ws) = (append!(widget.children, ws); nothing)
function Gtk.destroy(widget::MainWindow)
    if widget.window != nothing
        destroy(widget.window)
        widget.window = nothing
    end
end
Base.display(widget::MainWindow) = Gtk.showall(gtk_widget(widget))

function closerefs!{S<:Signal}(refs::Vector{S})
    for r in refs
        close(r)
    end
    empty!(refs)
end

##################################################
    ## dialogs
"""

Dialogs are modal so that they block until there is a returned value

"""
abstract Dialog <: Widget

immutable MessageBox <: Dialog
    msg
    style
end

"""

A simple modal message with a dismiss button.

* `msg`: a string
* `style`: a symbol. One of `:info`, `:warn`, `:error`

```
messagebox("my message", :info)
```
"""
messagebox(msg::AbstractString, style::Symbol=:ask) = MessageBox(msg, style) |> gtk_widget

immutable ConfirmBox <: Dialog
    msg
end

"""

A simple modal message with a yes/no pair of buttons. Returns boolean.

* `msg`: a string

```
confirmbox("my message")
```
"""
confirmbox(msg::AbstractString) = ConfirmBox(msg)  |> gtk_widget

immutable InputBox <: Dialog
    msg
    default
end

"""

Modal dialog to gather user input. Returns `""` for a cancelled query.

Example:
```
inputbox("What is your name", "John Doe")
```
"""
inputbox(msg::AbstractString, default::AbstractString="")  = InputBox(msg, default)  |> gtk_widget

immutable OpenFile <: Dialog
    title
end

"""

File open dialog

Example
```
fname = openfile("Select a file")
```
"""
openfile(title::AbstractString="Select a file") = OpenFile(title)  |> gtk_widget

immutable SaveFile <: Dialog
    title
end

"""

Modal dialog to select a filename for saving a file

Example
```
fname = savefile("Select a file")
length(fname) > 0 && "save the file to fname ..."
```
"""
savefile(title::AbstractString="Save file") = SaveFile(title)  |> gtk_widget

immutable SelectDir <: Dialog
    title
end

"""

Modal dialog to select a directory.

Example
```
d = selectdir("Select a directory...")
```
"""
selectdir(title::AbstractString="Save file") = SelectDir(title)  |> gtk_widget

##################################################
## Typography
## If these are useful, they could easily be expanded
immutable Bold <: Layout label end
"""
Make bold text in label
"""
bold(label::AbstractString) = Bold(label)

immutable Emph <: Layout label end
"""
Make emphasized test in a label
"""
emph(label::AbstractString) = Emph(label)

immutable Code <: Layout label end
"""
Use typewriter font in text for a label
"""
code(label::AbstractString) = Code(label)

##################################################
## Manipulate

## Modifications for @manipulate

## We add these output widgets to `widget`
widget_dict = Dict{Symbol, Function}()
widget_dict[:plot]=cairographic
widget_dict[:text] =textarea
widget_dict[:label]=label
widget_dict[:progress]=progress

function widget(x::Symbol, args...)
    widget_dict[x]()
end

## Main changes come from needing to pass through a parent container in order to "display" objects
## in "display_widgets" we just use `push!` though we could define `display` methods but passing in the
## parent container makes that awkward.
function display_widgets(win, widgetvars)
    map(v -> Expr(:call, esc(:push!), win, esc(v)), widgetvars)
end

## In `@manipulate` the macro builds up an expression, a. The `display` method is used to access the runtime
## value. Here we add a special type to couple the expression with a parent container, the main window.

type ManipulateWidget
    a
    w
end

macro manipulate(expr)
    if expr.head != :for
        error("@manipulate syntax is @manipulate for ",
              " [<variable>=<domain>,]... <expression> end")
    end
    block = expr.args[2]
    if expr.args[1].head == :block
        bindings = expr.args[1].args
    else
        bindings = [expr.args[1]]
    end
    syms = Interact.symbols(bindings)

    ## Modifications
    w = mainwindow(title="@manipulate")
    a = Expr(:let, Expr(:block,
                        display_widgets(w, syms)...,
                        esc(Interact.map_block(block, syms))),
             map(Interact.make_widget, bindings)...)
             b = Expr(:call, :ManipulateWidget, a, w)
    b
end

## connnect up Reactive with GtkInteract
Base.push!(w::Interact.InputWidget, value) = push!(w.signal, value)
Base.push!(w::OutputWidget, value::Interact.Signal) = push!(w, Reactive.value(value))

"""

The `Reactive.map(function, signals...)` function takes the signals
passed to it and passes their values to the function each time a
signal is changed through the GUI. (Signals are not emitted when a
value is `push!`ed to them)

This is a convenience extending the values to input widgets, as
otherwise, their signal must be pulled out before passing to `map`
using the `Interact.signal` function.

```
sl = slider(1:10)
sl1 = slider(1:10)
append!(mainwindow(), [sl, sl1])
map(println, sl, sl1)
```

"""

Base.map(f::Function, ws::Interact.InputWidget...) = map(f, map(Interact.signal, ws)...)

"""
Eat the first n function calls, defaulting to n=1
"""
function eatn(f, n::Int=1)
    ctr = 0
    (args... ) -> begin
        if ctr >= n
            f(args...)
        end
        ctr = ctr + 1
        nothing
    end
end

"""

We suppress the initial call when mapping to a button. This allows
something like `destroy` to be used when a button is clicked. Set
`initial_call` to `true` to force the initial call.

"""
function Base.map(f::Function, b::Interact.Button, initial_call::Bool=false)
    ctr = 0
    map(b.signal) do args...
        ctr >= 1 && f()
        ctr = ctr + 1
    end
    initial_call && f()
end

"""

The `sampleon` function of `Reactive.jl` passes on the values of the signals (as a
tuple) when the button is pressed. This can be used when the effects
of changing a GUI element can take a long time to propogate. In that
case, a GUI can be come unresponsive. This allows the effects to be
called only when the button is pressed.

```
sl = slider(1:10)
sl1 = slider(1:10)
btn = button("click me")
w = mainwindow()
append!(w, [sl, sl1, btn])

map(println, sampleon(btn, sl, sl1))
```

This function is a convenience for the following composition:

`Reactive.sampleon(btn.signal, map(tuple, map(Interact.signal, ws...)))`
"""
function Reactive.sampleon(btn::Button, ws::Interact.InputWidget...)
    Reactive.sampleon(btn.signal, map(tuple, ws...))
end

"""

Find current value of a widget. Short hand for `Reactive.value(w.signal)`.

"""
value(w::Widget) = Reactive.value(w.signal)

## load Gtk specific things
include("Gtk/gtkwidget.jl")

end # module


text:
arg_dict = Dict{String, Any}(
"physics" => "Euler",
"dimension" => 2,
"use_DG" => true,
"operator_type" => "SBPOmega",
"run_type" => 1,
"real_time" => true,
"order" => 1,
"Flux_name" => "IRFlux",
"IC_name" => "ICUnsteadyVortex",
"volume_integral_type" => 2,
"face_integral_type" => 2,
"addFaceIntegrals" => true,
"addVolumeIntegrals" => true,
"Volume_flux_name" => "IRFlux",
"FaceElementIntegral_name" => "ECFaceIntegral",
"vortex_x0" => 5.0,
"numBC" => 0,
"delta_t" => 0.0005,
"t_max" => 0.50,
#"smb_name" => "unsteadyvortexp2_.smb",
#"dmg_name" => "unsteadyvortexp2_.dmg",
"smb_name" => "SRCMESHES/smallp4_.smb",
"dmg_name" => ".null",
"output_freq" => 50,
"write_vis" => false,
"do_postproc" => true,
"exact_soln_func" => "ICUnsteadyVortex",
"solve" => true,
"write_entropy" => true,
"write_entropy_freq" => 1,
)


text:
import Base: eof, close
export BinningMap, close, value, position, eof, advance!

type BinningMap
    reader::SequenceReader
    binSize::Int64
    position::Int64
    value::Float64
    skipDup::Bool
end

function BinningMap(reader::SequenceReader, binSize; skipDup=true)
    fm = BinningMap(reader, binSize, 0, 0.0, skipDup)
    advance!(fm)
    fm
end
close(fm::BinningMap) = close(fm.reader)
value(fm::BinningMap) = fm.value
position(fm::BinningMap) = fm.position
eof(fm::BinningMap) = fm.position <= 0

function advance!(fm::BinningMap)
    fm.position = floor((fm.reader.position-1)/fm.binSize) + 1
    binEnd = fm.position*fm.binSize

    # Fill in the bin
    fm.value = 0.0
    lastPos = -1
    while fm.reader.position != -1 && fm.reader.position <= binEnd
        if !fm.skipDup || fm.reader.position != lastPos
            fm.value += 1
        end
        lastPos = fm.reader.position
        advance!(fm.reader)
    end
end