text:
@doc raw"""
    usinc(θ::Real)

Unnormalized version of `sinc` function, i.e. $\operatorname{usinc}(θ) = \frac{\sin(θ)}{θ}$.
This is equivalent to `sinc(θ/π)`.
"""
@inline usinc(θ::Real) = θ == 0 ? one(θ) : isinf(θ) ? zero(θ) : sin(θ) / θ

@doc raw"""
    usinc_from_cos(x::Real)

Unnormalized version of `sinc` function, i.e. $\operatorname{usinc}(θ) = \frac{\sin(θ)}{θ}$,
computed from $x = cos(θ)$.
"""
@inline function usinc_from_cos(x::Real)
    return if x >= 1
        one(x)
    elseif x <= -1
        zero(x)
    else
        sqrt(1 - x^2) / acos(x)
    end
end

@doc raw"""
    nzsign(z[, absz])

Compute a modified `sign(z)` that is always nonzero, i.e. where
````math
\operatorname(nzsign)(z) = \begin{cases}
    1 & \text{if } z = 0\\
    \frac{z}{|z|} & \text{otherwise}
\end{cases}
````
"""
@inline function nzsign(z, absz=abs(z))
    psignz = z / absz
    return ifelse(iszero(absz), one(psignz), psignz)
end

allocate(p, s::Size{S}) where {S} = similar(p, S...)
allocate(p::StaticArray, s::Size{S}) where {S} = similar(p, maybesize(s))
allocate(p, ::Type{T}, s::Size{S}) where {S,T} = similar(p, T, S...)
allocate(p::StaticArray, ::Type{T}, s::Size{S}) where {S,T} = similar(p, T, maybesize(s))

"""
    eigen_safe(x)

Compute the eigendecomposition of `x`. If `x` is a `StaticMatrix`, it is
converted to a `Matrix` before the decomposition.
"""
@inline eigen_safe(x; kwargs...) = eigen(x; kwargs...)
@inline function eigen_safe(x::StaticMatrix; kwargs...)
    s = size(x)
    E = eigen!(Matrix(parent(x)); kwargs...)
    return Eigen(SizedVector{s[1]}(E.values), SizedMatrix{s...}(E.vectors))
end

"""
    log_safe(x)

Compute the matrix logarithm of `x`. If `x` is a `StaticMatrix`, it is
converted to a `Matrix` before computing the log.
"""
@inline log_safe(x) = log(x)
@inline function log_safe(x::StaticMatrix)
    s = Size(x)
    return SizedMatrix{s[1],s[2]}(log(Matrix(parent(x))))
end

# NOTE: workaround until https://github.com/JuliaLang/julia/pull/39973 or similar is merged
"""
    log_safe!(y, x)

Compute the matrix logarithm of `x`. If the eltype of `y` is real, then the imaginary part
of `x` is ignored, and a `DomainError` is raised if `real(x)` has no real logarithm.
"""
function log_safe!(Y, A)
    if eltype(Y) <: Real
        if ishermitian(A)
            eigenF = eigen(Symmetric(real(A)))
            i = findfirst(≤(0), eigenF.values)
            if i !== nothing
                throw(
                    DomainError(
                        eigenF.values[i],
                        "All eigenvalues must be positive to compute a real logarithm.",
                    ),
                )
            end
            mul!(Y, eigenF.vectors .* log.(eigenF.values'), eigenF.vectors')
        elseif istriu(A)
            i = findfirst(≤(0), @view(A[diagind(A)]))
            if i !== nothing
                throw(
                    DomainError(
                        A[i, i],
                        "All eigenvalues must be positive to compute a real logarithm.",
                    ),
                )
            end
            copyto!(Y, real(log(UpperTriangular(A))))
        else
            schurF = schur(convert(Matrix, real(A)))
            i = findfirst(x -> isreal(x) && real(x) ≤ 0, schurF.values)
            if i !== nothing
                throw(
                    DomainError(
                        schurF.values[i],
                        "All eigenvalues must be positive to compute a real logarithm.",
                    ),
                )
            end
            if istriu(schurF.T)
                mul!(Y, schurF.Z, real(log(UpperTriangular(schurF.T))) * schurF.Z')
            else
                schurS = schur(complex(schurF.T))
                Y .= real.(schurS.Z * log(UpperTriangular(schurS.T)) * schurS.Z')
                mul!(Y, schurF.Z * Y, schurF.Z')
            end
        end
    else
        copyto!(Y, log_safe(A))
    end
    return Y
end

"""
    mul!_safe(Y, A, B) -> Y

Call `mul!` safely, that is, `A` and/or `B` are permitted to alias with `Y`.
"""
mul!_safe(Y, A, B) = (Y === A || Y === B) ? copyto!(Y, A * B) : mul!(Y, A, B)

@doc raw"""
    realify(X::AbstractMatrix{T𝔽}, 𝔽::AbstractNumbers) -> Y::AbstractMatrix{<:Real}

Given a matrix $X ∈ 𝔽^{n × n}$, compute $Y ∈ ℝ^{m × m}$, where $m = n \operatorname{dim}_𝔽$,
and $\operatorname{dim}_𝔽$ is the [`real_dimension`](@ref) of the number field $𝔽$, using
the map $ϕ \colon X ↦ Y$, that preserves the matrix product, so that for all
$C,D ∈ 𝔽^{n × n}$,
````math
ϕ(C) ϕ(D) = ϕ(CD).
````
See [`realify!`](@ref) for an in-place version, and [`unrealify!`](@ref) to compute the
inverse of $ϕ$.
"""
function realify(X, 𝔽)
    n = LinearAlgebra.checksquare(X)
    nℝ = real_dimension(𝔽) * n
    Y = allocate(X, real(eltype(X)), nℝ, nℝ)
    return realify!(Y, X, 𝔽, n)
end
realify(X, ::typeof(ℝ)) = X

"""
    realify!(Y::AbstractMatrix{<:Real}, X::AbstractMatrix{T𝔽}, 𝔽::AbstractNumbers)

In-place version of [`realify`](@ref).
"""
realify!(Y, X, 𝔽)

@doc raw"""
    realify!(Y::AbstractMatrix{<:Real}, X::AbstractMatrix{<:Complex}, ::typeof(ℂ))

Given a complex matrix $X = A + iB ∈ ℂ^{n × n}$, compute its realified matrix
$Y ∈ ℝ^{2n × 2n}$, written
where
````math
Y = \begin{pmatrix}A & -B \\ B & A \end{pmatrix}.
````
"""
function realify!(Y, X, ::typeof(ℂ), n=LinearAlgebra.checksquare(X))
    for i in 1:n, j in 1:n
        Xr, Xi = reim(X[i, j])
        Y[i, j] = Y[n + i, n + j] = Xr
        Y[n + i, j] = Xi
        Y[i, n + j] = -Xi
    end
    return Y
end

@doc raw"""
    unrealify!(X::AbstractMatrix{T𝔽}, Y::AbstractMatrix{<:Real}, 𝔽::AbstractNumbers[, n])

Given a real matrix $Y ∈ ℝ^{m × m}$, where $m = n \operatorname{dim}_𝔽$, and
$\operatorname{dim}_𝔽$ is the [`real_dimension`](@ref) of the number field $𝔽$, compute
in-place its equivalent matrix $X ∈ 𝔽^{n × n}$. Note that this function does not check that
$Y$ has a valid structure to be un-realified.

See [`realify!`](@ref) for the inverse of this function.
"""
unrealify!(X, Y, 𝔽)

function unrealify!(X, Y, ::typeof(ℂ), n=LinearAlgebra.checksquare(X))
    for i in 1:n, j in 1:n
        X[i, j] = complex((Y[i, j] + Y[n + i, n + j]) / 2, (Y[n + i, j] - Y[i, n + j]) / 2)
    end
    return X
end
unrealify!(Y, X, ::typeof(ℝ), args...) = copyto!(Y, X)

@generated maybesize(s::Size{S}) where {S} = prod(S) > 100 ? S : :(s)

"""
    select_from_tuple(t::NTuple{N, Any}, positions::Val{P})

Selects elements of tuple `t` at positions specified by the second argument.
For example `select_from_tuple(("a", "b", "c"), Val((3, 1, 1)))` returns
`("c", "a", "a")`.
"""
@generated function select_from_tuple(t::NTuple{N,Any}, positions::Val{P}) where {N,P}
    for k in P
        (k < 0 || k > N) && error("positions must be between 1 and $N")
    end
    return Expr(:tuple, [Expr(:ref, :t, k) for k in P]...)
end

@doc raw"""
    vec2skew!(X, v, k)

create a skew symmetric matrix inplace in `X` of size $k\times k$ from a vector `v`,
for example for `v=[1,2,3]` and `k=3` this
yields
````julia
[  0  1  2;
  -1  0  3;
  -2 -3  0
]
````
"""
function vec2skew!(X, v)
    k = size(X)[1]
    size(X)[2] != k && error("X is of wrong size, expected ($k,$k) got $(size(X)).")
    n = div(k * (k - 1), 2)
    length(v) < n && error("The vector $(v) is too short, expected $(n) got $(length(v)).")
    m = 0
    X .= [i < j ? (m += 1; v[m]) : zero(eltype(v)) for i in 1:k, j in 1:k]
    X .= X - X'
    return X
end
function vec2skew(v, k)
    X = similar(v, k, k)
    vec2skew!(X, v)
    return X
end

@doc raw"""
    isnormal(x; kwargs...) -> Bool

Check if the matrix or number `x` is normal, that is, if it commutes with its adjoint:
````math
x x^\mathrm{H} = x^\mathrm{H} x.
````
By default, this is an equality check. Provide `kwargs` for `isapprox` to perform an
approximate check.
"""
function isnormal(x; kwargs...)
    (isdiag(x) || ishermitian(x)) && return true
    isempty(kwargs) && return x * x' == x' * x
    return isapprox(x * x', x' * x; kwargs...)
end
isnormal(::LinearAlgebra.RealHermSymComplexHerm; kwargs...) = true

"""
    ziptuples(a, b[, c[, d[, e]]])

Zips tuples `a`, `b`, and remaining in a fast, type-stable way. If they have different
lengths, the result is trimmed to the length of the shorter tuple.
"""
@generated function ziptuples(a::NTuple{N,Any}, b::NTuple{M,Any}) where {N,M}
    ex = Expr(:tuple)
    for i in 1:min(N, M)
        push!(ex.args, :((a[$i], b[$i])))
    end
    return ex
end
@generated function ziptuples(
    a::NTuple{N,Any},
    b::NTuple{M,Any},
    c::NTuple{L,Any},
) where {N,M,L}
    ex = Expr(:tuple)
    for i in 1:min(N, M, L)
        push!(ex.args, :((a[$i], b[$i], c[$i])))
    end
    return ex
end
@generated function ziptuples(
    a::NTuple{N,Any},
    b::NTuple{M,Any},
    c::NTuple{L,Any},
    d::NTuple{K,Any},
) where {N,M,L,K}
    ex = Expr(:tuple)
    for i in 1:min(N, M, L, K)
        push!(ex.args, :((a[$i], b[$i], c[$i], d[$i])))
    end
    return ex
end
@generated function ziptuples(
    a::NTuple{N,Any},
    b::NTuple{M,Any},
    c::NTuple{L,Any},
    d::NTuple{K,Any},
    e::NTuple{J,Any},
) where {N,M,L,K,J}
    ex = Expr(:tuple)
    for i in 1:min(N, M, L, K, J)
        push!(ex.args, :((a[$i], b[$i], c[$i], d[$i], e[$i])))
    end
    return ex
end


text:
module AssociativeArrays

export AssociativeArray

import Base: getindex, setindex!, size, similar
using Base.Test

immutable AssociativeArray{T,N,C<:Associative} <: AbstractArray{T,N}
    data::C
    dims::NTuple{N,Int}
end

AssociativeArray{K,V,N}(data::Associative{K,V}, dims::NTuple{N,Int}) = AssociativeArray{V, N, typeof(data)}(data, dims)
AssociativeArray(data::Associative, dims::Int) = AssociativeArray(data, (dims,))

# outer constructor infers the types from the data
# triangular dispatch constraints are always satisfied if you use this constructor
# these constructors use run time type functions so are unstable.
# function AssociativeArray(data, dims)
#     nt, t = eltype(data).types
#     a=length(nt.types)
#     n = max(1, a)
#     C = typeof(data)
#     return AssociativeArray{t, n, C}(data, dims)
# end
# #special case for the 1D case so you can index by plain integers rather than (Int,)
# function AssociativeArray(data, dims::Integer)
#     nt, t = eltype(data).types
#     return AssociativeArray{t, 1, typeof(data)}(data, (dims,))
# end

# this doesn't work because you can't pass a Dict you need to instantiate
# is with something in order to instantiate it.
# type SparseVec{T,N,C<:Associative}
#   data::C{T,N}
#   dims::N
# end

# this doesn't work because you can't nest parametric types.
# type SparseVec{T,N,C<:Associative{N,T}}
#     data::C
#     dims::N
# end

# satisfying the interface
size(sa::AssociativeArray) = sa.dims
@inline getindex(sa::AssociativeArray, indices...) = sa.data[indices]
@inline setindex!(sa::AssociativeArray, x, indices...) = setindex!(sa.data, x, indices)
similar{T}(sa::AssociativeArray, ::Type{T}, dims) = AssociativeArray(Dict{typeof(dims), T}(), dims)
#special 1d case
@inline getindex{T,C<:Associative}(sa::AssociativeArray{T,1,C}, index::Integer) = sa.data[index]
@inline setindex!{T,C<:Associative}(sa::AssociativeArray{T,1,C}, x, index::Integer) = setindex!(sa.data, x, index)

# TODO enable abbreviated printing.
# println(AssociativeArray(Dict{NTuple{2, Int}, Float64}(), (150,150)))

type SparseVector{T,C}
    data::C
    length::NTuple{1, Int}
end
end

ext: jl


text:
module Evolution

# Content

end

max_stars_repo_name: TheSquareRoot/Evolution

max_issues_repo_name: TheSquareRoot/Evolution

max_forks_repo_name: TheSquareRoot/Evolution


text:
"""
    AbstractSearch{DATA}

Branch and bound search interface in element of type DATA.

This interface provide an iterable that perform the search.

There is currently three types of search supported `BreadFirstAbstractSearch`,
`AbstractDepthFirstSearch` and `AbstractKeySearch`, each one processing the element of the
tree in a different order. When subtyping one of these, the following methods
must be implemented:
  - `root_element(::AbstractSearch)`: return the element with which the search is started
  - `process(::AbstractSearch, elem::DATA)`: return a symbol representing the action
        to perform with the element `elem` and an object of type `DATA` representing
        the state of the element after processing (may return `elem` unchanged).
  - `bisect(::AbstractSearch, elem::DATA)`: return two elements of type `DATA` build
        by bisecting `elem`

Subtyping `AbstractSearch` directly allows to have control over the order in which
the elements are process. To do this the following methods must be implemented:
  - `root_element(::AbstractSearch)`: return the first element to be processed. Use
        to build the initial tree.
  - `get_leaf_id!(::AbstractSearch, wt::BPTree)`: return the id of the next leaf that
        will be processed and remove it from the list of working leaves of `wt`.
  - `insert_leaf!(::AbstractSearch, wt::BPTree, leaf::BPLeaf)`: insert a leaf in the
        list of working leaves.

# Valid symbols returned by the process function
  - `:store`: the element is considered as final and is stored, it will not be
        further processed
  - `:bisect`: the element is bisected and each of the two resulting part will
        be processed
  - `:discard`: the element is discarded from the tree, allowing to free memory
"""
abstract type AbstractSearch{DATA} end

abstract type AbstractBreadthFirstSearch{DATA} <: AbstractSearch{DATA} end
abstract type AbstractDepthFirstSearch{DATA} <: AbstractSearch{DATA} end

# TODO should be smallest key first to match sort functions
"""
    AbstractKeySearch{DATA} <: AbstractSearch{DATA}

Interface to a branch and bound search that use a key function to decide which
element to process first. The search process first the element with the largest
key as computed by `keyfunc(ks::AbstractKeySearch, elem)`.

!!! warning
    Untested.
"""
abstract type AbstractKeySearch{DATA} <: AbstractSearch{DATA} end

"""
    root_element(search::AbstractSearch)

Return the initial element of the search. The `BPTree` will be build around it.

Can be define for custom searches that are direct subtype of `AbstractSearch`, default
behavior is to fetch the field `initial` of the search.
"""
root_element(search::AbstractSearch) = search.initial

"""
    get_leaf_id!(::AbstractSearch, wt::BPTree)

Return the id of the next leaf that will be processed and remove it from the
list of working leaves.

Must be define for custom searches that are direct subtype of `AbstractSearch`.
"""
get_leaf_id!(::AbstractBreadthFirstSearch, wt::BPTree) = popfirst!(wt.working_leaves)
get_leaf_id!(::AbstractDepthFirstSearch, wt::BPTree) = pop!(wt.working_leaves)
get_leaf_id!(::AbstractKeySearch, wt::BPTree) = popfirst!(wt.working_leaves)

"""
    insert_leaf!(::AbstractSearch, wt::BPTree, leaf::BPLeaf)

Insert the id of a new leaf that has been produced by bisecting an older leaf
into the list of working leaves.

Must be define for custom searches that are direct subtype of `AbstractSearch`.
"""
function insert_leaf!(::Union{AbstractBreadthFirstSearch{DATA}, AbstractDepthFirstSearch{DATA}},
                      wt::BPTree{DATA}, leaf::BPLeaf{DATA}) where {DATA}
    id = newid(wt)
    wt.leaves[id] = leaf
    push!(wt.working_leaves, id)
    return id
end

function insert_leaf!(::KS, wt::BPTree{DATA}, leaf::BPLeaf{DATA}) where {DATA, KS <: AbstractKeySearch{DATA}}
    id = newid(wt)
    wt.leaves[id] = leaf
    keys = keyfunc.(KS, wt.working_leaves)
    current = keyfunc(KS, leaf)

    # Keep the working_leaves sorted
    insert!(wt.working_leaves, searchsortedfirst(keys, current), id)
    return id
end

eltype(::Type{BPS}) where {DATA, BPS <: AbstractSearch{DATA}} = BPTree{DATA}
IteratorSize(::Type{BPS}) where {BPS <: AbstractSearch} = Base.SizeUnknown()

function iterate(search::AbstractSearch{DATA},
                 wt::BPTree=BPTree(root_element(search))) where {DATA}

    isempty(wt.working_leaves) && return nothing

    id = get_leaf_id!(search, wt)
    X = wt.leaves[id]
    action, newdata = process(search, data(X))
    if action == :store
        wt.leaves[id] = BPLeaf(newdata, X.parent, :final)
    elseif action == :bisect
        child1, child2 = bisect(search, newdata)
        leaf1 = BPLeaf(child1, id, :working)
        leaf2 = BPLeaf(child2, id, :working)
        id1 = insert_leaf!(search, wt, leaf1)
        id2 = insert_leaf!(search, wt, leaf2)
        wt.nodes[id] = BPNode(X, id1, id2)
        delete!(wt.leaves, id)
    elseif action == :discard
        discard_leaf!(wt, id)
    else
        error("Branch and bound: process function of the search object return " *
              "unknown action: $action for element $X. Valid actions are " *
              ":store, :bisect and :discard.")
    end
    return wt, wt
end


text:
let gene="ABL1"
    @test query_geneloc(gene) == ("chr9", "133589333", "133763062")
end

let gene = "TNFRSF14"
    @test query_gene("chr1","2490013") == "TNFRSF14"
end


text:
"""
    @test_modern

Use this to mark tests that should pass on modern (>=v1.4) versions of julia,
but should be marked as `@test_broken` on Julia 1.0.
For more complicated cases, or where it shouldn't be defined at all on julia versions,
you can use `if VERSION` directly.
"""
macro test_modern(expr...)
    if VERSION >= v"1.4"
        return :(@test $(expr...)) |> esc
    else
        return nothing
    end
end


text:
include("drawFlow.jl")

function main()

    # set up a block series to test the flow drawing method
    blocks= []

    push!(blocks,String("raw string code",))

    push!(blocks,if_block("y==4","code to run",true))

    push!(blocks,if_else_block("y==7","code (else if) to run",true))

    push!(blocks,else_block("(code else) to run"))

    drawFlow(blocks,"out.svg",1.0)
end

main()

max_stars_repo_name: aklapatch/CodeVisualized

max_issues_repo_name: aklapatch/CodeVisualized

max_forks_repo_name: aklapatch/CodeVisualized


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

using Markdown
using InteractiveUtils

# ╔═╡ e708e17d-617e-4a3f-bca1-164d22c3dffc
begin
    import Pkg
    Pkg.activate("..")
end

# ╔═╡ b0dbee5a-9881-48be-9fe8-56cc4f09f4dd
begin
    using Revise
    using Hydroelastics

    using MeshCat
    using GeometryBasics
    using CoordinateTransformations
end

# ╔═╡ 27b6e401-897c-459a-b6f6-575f3f89b0b5
begin
    vis = Visualizer()
    IJuliaCell(vis)
end

# ╔═╡ d4f9ca2a-1913-48a6-bed1-ab7633615752
begin
    setobject!(vis, HyperRectangle(Vec(0.0, 0, 0), Vec(1.0, 1, 1)))
    settransform!(vis, Translation(-0.5, -0.5, 0))
end

# ╔═╡ Cell order:
# ╠═e708e17d-617e-4a3f-bca1-164d22c3dffc
# ╠═b0dbee5a-9881-48be-9fe8-56cc4f09f4dd
# ╠═27b6e401-897c-459a-b6f6-575f3f89b0b5
# ╠═d4f9ca2a-1913-48a6-bed1-ab7633615752

max_stars_repo_name: ekzhang/hydroelastics

max_issues_repo_name: ekzhang/hydroelastics

max_forks_repo_name: ekzhang/hydroelastics


text:
export DoorKey

using Random

mutable struct DoorKey{W<:GridWorldBase} <: AbstractGridWorld
    world::W
    agent_pos::CartesianIndex{2}
    agent::Agent
end

function DoorKey(;n=8, agent_start_pos=CartesianIndex(2,2), rng=Random.GLOBAL_RNG)
    objects = (EMPTY, WALL, GOAL, Door(:yellow), Key(:yellow))
    world = GridWorldBase(objects, n, n)
    world[EMPTY, 2:n-1, 2:n-1] .= true
    world[WALL, [1,n], 1:n] .= true
    world[WALL, 1:n, [1,n]] .= true
    world[EMPTY, n-1, n-1] = false
    world[GOAL, n-1, n-1] = true

    split_idx = rand(3:n-2)
    world[WALL, 1:n, split_idx] .= true
    world[EMPTY, 1:n, split_idx] .= false
    door_idx = rand(2:n-1)
    world[Door(:yellow), door_idx, split_idx] = true
    world[WALL, door_idx, split_idx] = false

    key_pos = CartesianIndex(rand(2:n-2), rand(2:split_idx-1))
    if key_pos == CartesianIndex(2,2)
        key_pos = CartesianIndex(3,2)
    end
    world[EMPTY, key_pos] = false
    world[Key(:yellow), key_pos] = true

    DoorKey(world, agent_start_pos, Agent(;dir=RIGHT))
end

function (w::DoorKey)(::MoveForward)
    dir = get_dir(w.agent)
    dest = dir(w.agent_pos)

    if w.world[Key(:yellow), dest]
        if PICK_UP(w.agent, Key(:yellow))
            w.world[Key(:yellow), dest] = false
            w.world[EMPTY, dest] = true
        end
        w.agent_pos = dest
    elseif w.world[Door(:yellow), dest] && w.agent.inventory !== Key(:yellow)
        nothing
    elseif w.world[Door(:yellow), dest] && w.agent.inventory === Key(:yellow)
        w.agent_pos = dest
    elseif dest ∈ CartesianIndices((size(w.world, 2), size(w.world, 3))) && !w.world[WALL,dest]
        w.agent_pos = dest
    end

    w
end


text: A[end]

max_stars_repo_name: lidiapopescu/JuliaDT

max_issues_repo_name: lidiapopescu/JuliaDT

max_forks_repo_name: JuliaComputing/JuliaDT


text:
export GOT_Q, @k_sym
let h = unsafe_load(cglobal(:jl_exe_handle, Ptr{Void}))
    # Is Julia running embedded in q?
    global const GOT_Q = Libdl.dlsym_e(h, :b9) != C_NULL
end

if GOT_Q  # Get q C API from the current process
    macro k_sym(func)
        esc(func)
    end
else      # Load "c" DLL
    const SYS_CHAR = Dict(
        :Linux => 'l',
        :Darwin => 'm',
    )
    const SYS_ARCH = @sprintf("%c%d", SYS_CHAR[Sys.KERNEL], Sys.WORD_SIZE)
    const C_SO_PATH = joinpath(dirname(@__FILE__), SYS_ARCH, "c")
    const C_SO = Libdl.dlopen(C_SO_PATH,
                Libdl.RTLD_LAZY|Libdl.RTLD_DEEPBIND|Libdl.RTLD_GLOBAL)

    macro k_sym(func)
        :(($(esc(func)), C_SO_PATH))
    end
end


text:
using ConvOptDL
using ConvOptDL.Utils
using StatsBase
using LinearAlgebra: I
using ArgParse
using Serialization
using Pipe: @pipe
using Flux: update!, trainmode!
using Statistics: mean
using Logging: @info
using Dates: now
import Flux

function loss_log_softmax(logits, one_hot_vec)
    log_prob = log.(Flux.softmax(logits, dims = 1))
    loss = @pipe -(log_prob .* one_hot_vec) |> sum(_, dims = 1)
    mean(loss)
end

_format_batch(batch) = reshape(Float32.(batch), size(batch)[1:end-2]..., :)

function train!(loss, model, batch, opt)
    # https://fluxml.ai/Flux.jl/stable/training/training/#Custom-Training-loops-1
    ps = Flux.params(model)
    # smoothed onehot encoding
    # onehot is not differentiable, moving out of the do block
    query_onehot = @pipe ConvOptDL.Utils.onehot(batch.query_labels) |>
          reshape(_, batch.support_n_ways, :) |>
          _ .* (1 - 5f-2) .+ 5f-2 * (1 .- _) ./ (batch.support_n_ways - 1f0)
    support_onehot = ConvOptDL.Utils.onehot(batch.support_labels)

    local meta_loss
    gs = Flux.gradient(ps) do
        # (`feat_dim`, `n_ways`*`k_shots`, `num_tasks`)
        embed_query = reshape(
            model(_format_batch(batch.query_samples)),
            :,
            batch.n_query,
            size(batch),
        )
        embed_support = reshape(
            model(_format_batch(batch.support_samples)),
            :,
            batch.n_support,
            size(batch),
        )
        Q, p, G, h, A, b = crammer_svm(embed_support, support_onehot, batch)
        # solve QP
        X, λ, ν = solve_qp_batch(Q, p, G, h, A, b)
        α = @pipe X |>
              reshape(_, 1, batch.n_support, batch.support_n_ways, size(batch)) |>
              repeat(_, outer = (batch.n_query, 1, 1, 1)) |>
              permutedims(_, (3, 2, 1, 4))
        # logits
        logits = @pipe ConvOptDL.Utils.gram_matrix(embed_query, embed_support) |>
              reshape(_, 1, size(_)...) |>
              repeat(_, outer = (5, 1, 1, 1)) |>
              _ .* α |>
              sum(_, dims = 2) |>
              dropdims(_, dims = 2) |>
              reshape(_, batch.support_n_ways, :)
        meta_loss = loss(logits, query_onehot)
        return meta_loss
    end
    update!(opt, ps, gs)
    return meta_loss
end

function parse_opts()
    s = ArgParseSettings("Meta Convex Optimization Learning training")
    @add_arg_table! s begin
        "--batch-size"
        arg_type = Int64
        metavar = "INT"
        default = 8
        "--batches-per-episode"
        arg_type = Int64
        metavar = "INT"
        default = 200
        "--log-period"
        help = "print log for every given number of batches"
        arg_type = Int64
        metavar = "INT"
        default = 25
        "--num-episodes"
        arg_type = Int64
        metavar = "INT"
        default = 50
        "--learning-rate"
        arg_type = Float32
        default = 1f-1
        "-o"
        arg_type = String
        default = "model.jls"
        "data_file"
        arg_type = String
        default = "./test/test_data.jls"
    end
    args = parse_args(s)
    return args
end

if nameof(@__MODULE__) == :Main
    args = parse_opts()
    @info "Training Options: " args
    batch_size = args["batch-size"]
    batches_per_episode = args["batches-per-episode"]
    log_period = args["log-period"]
    num_episodes = args["num-episodes"]
    data_file = args["data_file"]
    out_model_file = args["o"]
    lr = args["learning-rate"]
    model_name, _ = splitext(out_model_file)
    dloader = FewShotDataLoader(data_file)
    opt = Flux.Optimise.Momentum(lr)
    model = resnet12()
    trainmode!(model)
    record = Dict()
    for episode = 1:num_episodes
        @info "episode $(episode) start" now()
        meta_losses = []
        total_time = 0
        for i = 1:batches_per_episode
            time_start = time()
            batch = sample(dloader, batch_size, support_n_ways = 5, support_k_shots = 2)
            meta_loss = train!(loss_log_softmax, model, batch, opt)
            time_end = time()
            push!(meta_losses, meta_loss)
            total_time += (time_end - time_start)
            if i % log_period == 0
                @info "[$(now())] last $(log_period) meta losses" meta_losses[end-log_period+1:end] mean(meta_losses)
                @info "[$(now())] average execution time" mean_exec_time = (total_time / i)
            end
        end
        @info "episode $(episode) end" now()
        record[episode] = meta_losses
        serialize("train_$(model_name)_meta_losses_$(episode).jls", meta_losses)
    end
    serialize("all_train_$(model_name)_meta_losses.jls", record)
    serialize(out_model_file, model)
end


text:
module Practice_Julia_Repositoria

greet() = print("Hello World!")
export greet

end # module


text:
println("hello world")
println("testing")
print("testing")

max_stars_repo_path: AlanZ/hello_world backup 1.jl

max_stars_repo_name: KorfLab/julie

max_issues_repo_path: AlanZ/hello_world backup 1.jl

max_issues_repo_name: KorfLab/julie

max_forks_repo_path: AlanZ/hello_world backup 1.jl

max_forks_repo_name: KorfLab/julie


text:
using FunctionalCollections
using Test

@testset "Persistent Lists" begin

    @testset "UnionFind" begin
        uf = PersistentDisjointSet{Int}()
        uf = union(uf,2,3)
        uf = union(uf,5,7)
        uf = union(uf,2,4)
        uf = union(uf,4,8)
        @assert uf[1] == 1
        @assert uf[1] != uf[7]
        @assert uf[3] == uf[8]
        @assert uf[5] == uf[7]
        @assert uf[3] != uf[7]
    end

end


text:
using Test, StableRNGs
using LighthouseFlux, Lighthouse, Flux, Random

using CairoMakie
CairoMakie.activate!(type="png")

mutable struct TestModel
    chain::Chain
end

(m::TestModel)(X) = m.chain(X)

Flux.@functor TestModel (chain,)

function LighthouseFlux.loss(m::TestModel, input_feature_batch, soft_label_batch)
    return Flux.mse(m(input_feature_batch), soft_label_batch)
end

@testset "learn!(::TestModel, ...)" begin
    mktempdir() do tmpdir
        rng = StableRNG(157)
        classes = ["class_$i" for i in 1:5]
        c, n = 5, 3
        Random.seed!(0)
        model = TestModel(Chain(Dense(4 * c, 2 * c),
                                Dense(2 * c, c), softmax))
        # assure that `testmode!` and `trainmode!` is being utilized correctly after training
        test_input = rand(StableRNG(42), Float32, 4 * c) # get test input
        y_pretrained = model(test_input) # test model output before being trained
        classifier = FluxClassifier(model, ADAM(0.1), classes)
        train_batches = [(rand(rng, 4 * c, n), rand(rng, 1, n)) for _ in 1:100]
        test_batches = [((rand(rng, 4 * c, n), rand(rng, 1, n)), (n * i - n + 1):(n * i))
                        for i in 1:10]
        possible_vote_labels = collect(0:length(classes))
        votes = [rand(rng, possible_vote_labels) for sample in 1:(n * 10), voter in 1:7]
        logger = Lighthouse.LearnLogger(joinpath(tmpdir, "logs"), "test_run")
        limit = 4
        let counted = 0
            upon_loss_decrease = Lighthouse.upon(logger,
                                                 "test_set_prediction/mean_loss_per_epoch";
                                                 condition=<, initial=Inf)
            callback = n -> begin
                upon_loss_decrease() do _
                    counted += n
                end
            end
            elected = majority.((rng,), eachrow(votes), (1:length(classes),))
            Lighthouse.learn!(classifier, logger, () -> train_batches, () -> test_batches,
                              votes, elected; epoch_limit=limit, post_epoch_callback=callback)
            # NOTE: the RNG chosen above just happens to allow this to work every time,
            # since the loss happens to actually "improve" on the random data each epoch
            @test counted == sum(1:limit)
        end
        for key in [
            "train/loss_per_batch"
            "train/forward_pass/time_in_seconds_per_batch"
            "train/forward_pass/gc_time_in_seconds_per_batch"
            "train/forward_pass/allocations_per_batch"
            "train/forward_pass/memory_in_mb_per_batch"
            "train/reverse_pass/time_in_seconds_per_batch"
            "train/reverse_pass/gc_time_in_seconds_per_batch"
            "train/reverse_pass/allocations_per_batch"
            "train/reverse_pass/memory_in_mb_per_batch"
            "train/update/time_in_seconds_per_batch"
            "train/update/gc_time_in_seconds_per_batch"
            "train/update/allocations_per_batch"
            "train/update/memory_in_mb_per_batch"
        ]
            @test length(logger.logged[key]) == length(train_batches) * limit
        end
        for key in [
            "test_set_prediction/loss_per_batch"
            "test_set_prediction/time_in_seconds_per_batch"
            "test_set_prediction/gc_time_in_seconds_per_batch"
            "test_set_prediction/allocations_per_batch"
            "test_set_prediction/memory_in_mb_per_batch"
        ]
            @test length(logger.logged[key]) == length(test_batches) * limit
        end
        for key in [
            "test_set_prediction/mean_loss_per_epoch"
            "test_set_evaluation/time_in_seconds_per_epoch"
            "test_set_evaluation/gc_time_in_seconds_per_epoch"
            "test_set_evaluation/allocations_per_epoch"
            "test_set_evaluation/memory_in_mb_per_epoch"
        ]
            @test length(logger.logged[key]) == limit
        end
        @test length(logger.logged["test_set_evaluation/metrics_per_epoch"]) == limit

        # test `testmode!` is correctly utiltized
        y_post_train = model(test_input)
        @test y_post_train != y_pretrained # check if model has trained
        @test y_post_train == model(test_input) # make sure model output is determinsitic

        # test onehot/onecold overloads
        classifier = FluxClassifier(model, ADAM(0.1), classes;
                                    onehot=(x -> fill(x, length(classes))), onecold=sum)
        @test Lighthouse.onehot(classifier, 3) == fill(3, length(classes))
        @test Lighthouse.onecold(classifier, [0.31, 0.43, 0.13]) == 0.87
    end
end

@testset "`evaluate_chain_in_debug_mode`" begin
    chain = Chain(Dense(4, 2), Dense(2, 1))
    input = Float32[
        1 2 3 4
        5 6 7 8
        8 7 6 5
        4 3 2 1
    ]
    @test @test_logs((:info, "Executing layer 1 / 2..."), (:info, (2, 4)),
                     (:info, "Executing layer 2 / 2..."), (:info, (1, 4)),
                     LighthouseFlux.evaluate_chain_in_debug_mode(chain, input)) ≈
          chain(input)
end


text:
module HeadlessChromium

using HttpServer, WebSockets
import HTTP, URIParser, Mustache, JSON

import Base: send, close

export Target, TimeoutError, send, close, stopChromium

global const DEBUG = false

### Set the const pointing to the Chromium executable

depsjl = joinpath(dirname(@__FILE__), "..", "deps", "deps.jl")

isfile(depsjl) || error("""HeadlessChromium not properly installed.
Please run
   Pkg.build(\"HeadlessChromium\")""")

include(depsjl)

### Includes + global var 'chromiumHandle' definition

global chromiumHandle = nothing

include("communication.jl")
include("send.jl")
include("misc.jl")

end # module


text:
### load all experiments, then combine them using the first experiment's 0 ppm, sampling rate, and spectral width.
function loadexperiments(paths::Vector{String};
    weights = ones(Float64, length(paths)),
    solvent_ppm_guess = 4.7,
    solvent_window_ppm = 0.1)

    M = length(paths)
    s_t_set = Vector{Vector{Complex{Float64}}}(undef, M)
    ν_0ppm_set = Vector{Float64}(undef, M)
    fs_set = Vector{Float64}(undef, M)
    SW_set = Vector{Float64}(undef, M)

    for n = 1:M

        s_t, S, hz2ppmfunc, ppm2hzfunc, ν_0ppm_set[n], fs_set[n], SW_set[n], α_0ppm, β_0ppm, λ_0ppm, Ω_0ppm,
        results_0ppm, dic, α_solvent, β_solvent, λ_solvent, Ω_solvent,
        results_solvent = loadspectrum(paths[n];
        solvent_ppm = solvent_ppm_guess,
        solvent_window_ppm = solvent_window_ppm)

        s_t_set[n] = weights[n] .* s_t
    end

    return s_t_set, ν_0ppm_set, fs_set, SW_set
end

function assemblemixture(paths::Vector{String};
    weights = ones(Float64, length(paths)),
    solvent_ppm = 4.7,
    solvent_window_ppm = 0.1)

    s_t_set, ν_0ppm_set, fs_set, SW_set = loadexperiments(paths;
    weights = weights,
    solvent_ppm_guess = solvent_ppm,
    solvent_window_ppm = solvent_window_ppm)

    # combine via average.
    s_t = sum(s_t_set) ./ length(s_t_set)

    return s_t, ν_0ppm_set[1], fs_set[1], SW_set[1],
        s_t_set, ν_0ppm_set, fs_set, SW_set # just have it.
end


text:
using Plots

#Useful functions

function hat(i, N, x)
    @assert 0 <= i < N
    h = 1 / (N-1) #Domain then goes from 0 to 1
    xim1 = h * (i-1)
    xi = h * i
    xip1 = h * (i+1)
    if x <= xim1
        return 0.
    elseif x <= xi
        return (x - xim1) / h
    elseif x <= xip1
        return -(x - xip1) / h
    else
        return 0.
    end
end

function sample(f, N)
    h = 1 / (N-1)
    return coeffs = Float64[f(h*i) for i in 0:N-1]
end

function evaluate(coeffs, x)
    N = length(coeffs)
    sum(coeffs[i+1] * hat(i, N, x) for i in 0:N-1)
end

function deriv(f)
    N = length(f)
    h = 1 / (N-1)
    deriv_coeffs = similar(f)
    for i in 1:N
        if i == 1
            deriv_coeffs[i] = (f[2] - f[1]) / h
        elseif i == N
            deriv_coeffs[i] = (f[N] - f[N-1]) / h
        else
            deriv_coeffs[i] = (f[i+1] - f[i-1]) / (2*h)
        end
    end
    return deriv_coeffs
end

function deriv2(f)
    N = length(f)
    h = 1 / (N-1)
    deriv2_coeffs = similar(f)
    for i in 1:N
        if i == 1
            deriv2_coeffs[i] = (f[3] - 2*f[2] + f[1]) / h^2
        elseif i == N
            deriv2_coeffs[i] = (f[N] - 2*f[N-1] + f[N-2]) / h^2
        else
            deriv2_coeffs[i] = (f[i+1] - 2*f[i] + f[i-1]) / h^2
        end
    end
    return deriv2_coeffs
end

function rk2step(f, y0, h)
    k0 = f(y0)
    y1 = y0 + (h/2)*k0
    k1 = f(y1)
    y2 = y0 + h*k1
    y2
end

function rk2integrate(f, y0, h, n)
    res = [y0]
    for i in 1:n
        y1 = rk2step(f, y0, h)
        push!(res, y1)
        y0 = y1
    end
    res
end

#Solving now the wave equation in a box.
#My state vector will be an array where the first index
#indicates position and the second index indicates time step

L=100
T=100

#Create the array

u = zeros(L,T)

#Boundary Conditions

for n=1:T
   u[1,n]=0.
   u[L,n]=0.
end

#Initial Conditions

f(x)=sin(2*pi*x)

#Discretization Parameter.
#It has a big effect on stability

s = 1

#Finite difference scheme

fi=sample(f,L)
u[:,1]=fi

for j=2:L-1
  u[j,2]=u[j,1]-1/2*s*(u[j+1,1]-2*u[j,1]+u[j-1,1])
end

for n=2:T-1
    for j=2:L-1
        u[j,n+1]=s*(u[j+1,n]+u[j-1,n]) + 2*(1.0-s)*u[j,n]-u[j,n-1]
    end
end

#Plot the solution

g(x,j)=evaluate(u[:,j],x)

u5(x)=g(x,5)
u10(x)=g(x,10)
u25(x)=g(x,25)
u27(x)=g(x,27)
u30(x)=g(x,30)
u50(x)=g(x,50)

labels = ["t=5" "t=10" "t=25" "t=27" "t=30" "t=50"]

plot([u5,u10,u25,u27,u30,u50],0,1,xlabel="x",ylabel="y",plot_title=`Hello`,label=labels)

#Energy of the solution (Conserved)- Using to check for stability

ud=similar(u)
udt=similar(u)

for n=1:T
   ud[:,n]=deriv(u[:,n])
end

for i=1:L
   udt[i,:]=deriv(u[i,:])
end

h(x,j)=evaluate(ud[:,j],x)
k(x,j)=evaluate(udt[:,j],x)

E=[]
ind=[]

for j=1:T

     p(x)=(h(x,j)^2+k(x,j)^2)

     Energy=rk2integrate(p,0.001,0.001,9000)[9000]

     push!(E,Energy)

     push!(ind,j)

end

plot(ind,E,ylabel="Energy",xlabel="Time step",label="Energy")

## Check for convergence

E=[]
ind=[]

for i=1:10

    T=10*i

    L=10*i

    u = zeros(L,T)

    for n=1:T
       u[1,n]=0.
       u[L,n]=0.
    end

    fi=sample(f,L)
    u[:,1]=fi

    for j=2:L-1
      u[j,2]=u[j,1]-1/2*s*(u[j+1,1]-2*u[j,1]+u[j-1,1])
    end

    for n=2:T-1
        for j=2:L-1
            u[j,n+1]=s*(u[j+1,n]+u[j-1,n]) + 2*(1.0-s)*u[j,n]-u[j,n-1]
        end
    end

    ud=similar(u)
    udt=similar(u)

    for n=1:T
       ud[:,n]=deriv(u[:,n])
    end

    for i=1:L
       udt[i,:]=deriv(u[i,:])
    end

    #Let's take a particular j, that is, a particular time step

    j=10

    g(x,j)=evaluate(u[:,j],x)
    h(x,j)=evaluate(ud[:,j],x)
    k(x,j)=evaluate(udt[:,j],x)

    p(x)=(h(x,j)^2+k(x,j)^2)

    Energy=rk2integrate(p,0.001,0.001,9000)[9000]

    push!(E,Energy)

    push!(ind,L)

end

plot(ind,ylabel="Energy",xlabel="Resolution",label="Energy")

max_stars_repo_name: JuliMaristany/WaveEquationBox

max_issues_repo_name: JuliMaristany/WaveEquationBox

max_forks_repo_name: JuliMaristany/WaveEquationBox


text:
using MAT
using Dates
using GLPK
using JuMP
using Statistics
using CartesianGeneticProgramming
using Cambrian
import Cambrian.mutate
using Plots
using Plots.PlotMeasures
using StatsPlots
include("../utils.jl")

function get_wave_t_in_range(param, min, max)
    indices_more = findall(>=(min), param)
    indices_less = findall(<(max), param)

    indices_to_keep = intersect(indices_more, indices_less)
    param = param[indices_to_keep]

    param, indices_to_keep
end

function fit_poly_inf_norm(x,y,n)
	m = Model(GLPK.Optimizer)

    @variable(m, alpha[1:n+1])
	@variable(m, t >= 0.)
	@constraint(m, [i=1:length(y)], sum(alpha[j]*x[i]^(j-1) for j in 1:n+1) - y[i] <= t);
	@constraint(m, [i=1:length(y)], -t <= sum(alpha[j]*x[i]^(j-1) for j in 1:n+1) - y[i]);
	@objective(m, Min, t)

	optimize!(m) # added new 2020-04-08

    return JuMP.value.(alpha)
end

f2(x, args) = (args[2] * x) .+ args[1]

const MATLAB_EPOCH = Dates.DateTime(-1,12,31)

date2num(d::Dates.DateTime) = Dates.value(d-MATLAB_EPOCH)/(1000*60*60*24)
num2date(n::Number) =  MATLAB_EPOCH + Dates.Millisecond(round(Int64, n*1000*60*60*24))

function plot_bathy(xs, profile; show_points=false, x_last_positive=NaN, x_first_negative=NaN, value_last_positive=NaN, value_first_negative=NaN, x_at_zero=NaN, title=NaN)
    p1 = plot(xs, profile, label="profile_t0")
    if show_points
        scatter!([x_last_positive, x_first_negative], [value_last_positive, value_first_negative], label="A, B")
        scatter!([x_at_zero], [0], label="Shoreline")
    end
    display(plot(p1, title=title))
end

function find_x_at_zero(profile, xs, i)
    t0_index_last_positive = findlast(profile.>0)
    if t0_index_last_positive == nothing t0_index_last_positive = 1 end

    t0_index_first_negative = findfirst(profile.<0)
    if t0_index_first_negative == nothing t0_index_first_negative = 1 end

    t0_x_last_positive = xs[t0_index_last_positive]
    t0_x_first_negative = xs[t0_index_first_negative]
    t0_value_last_positive = profile[t0_index_last_positive]
    t0_value_first_negative = profile[t0_index_first_negative]

    t0_coeffs = fit_poly_inf_norm([t0_value_last_positive, t0_value_first_negative], [t0_x_last_positive, t0_x_first_negative], 1)
    t0_x_at_zero = f2([0], t0_coeffs)[1]
    t0_x_at_zero
end

struct Day
    bathy0
    bathy1

    tp
    hs
    direction

    shoreline_change
    shoreline_change_class

    input_row
end

dir_path = get_dataset_path(;dataset_name="grand_popo")
kalman_bathy_vars = matread(string(dir_path, "Kalman_bathy.mat"))
bathy_xb = kalman_bathy_vars["xb"][1, :]# 1×213 Matrix{Float64}
# bathy_tn = kalman_bathy_vars["tn"]      # 1×1214 Matrix{Float64}
bathy_hbK = kalman_bathy_vars["hbK"]    # 213×1214 Matrix{Float64}

min_target_matdate = minimum(filter(!isnan, kalman_bathy_vars["tn"]))
max_target_matdate = maximum(filter(!isnan, kalman_bathy_vars["tn"]))

wave_vars = matread(string(dir_path, "Wave_GPP_ERA.mat"))
wave_t = wave_vars["t"]
wave_Tm = wave_vars["Tm"]
wave_hs = wave_vars["hs"]
wave_dir = wave_vars["dir"]

wave_t, indices = get_wave_t_in_range(wave_t, min_target_matdate, max_target_matdate)
wave_Tm = wave_Tm[indices]
wave_hs = wave_hs[indices]
wave_dir = wave_dir[indices]

bathy_dates = map(num2date, collect(min_target_matdate:max_target_matdate))  # collect() instead of using bathy_tn because bathy_tn includes NaN's, even though bathy_hbK has no NaN
wave_dates = map(num2date, wave_t)

days = []
n_offshore = 0
n_landward = 0
n_none = 0

all_tm = []
all_hs = []
all_dir = []
for i in 2:length(bathy_dates)
    _d0 = bathy_dates[i-1]
    _d1 = bathy_dates[i]

    _t0 = date2num(_d0)
    _t1 = date2num(_d1)

    _bathy_t0 = bathy_hbK[:, i-1]
    _bathy_t1 = bathy_hbK[:, i]

    _wave_t = deepcopy(wave_t)
    _wave_Tm = deepcopy(wave_Tm)
    _wave_dir = deepcopy(wave_dir)
    _wave_hs = deepcopy(wave_hs)

    _wave_t, _indices = get_wave_t_in_range(_wave_t, _t0, _t1)
    _wave_Tm = _wave_Tm[_indices]
    _wave_hs = _wave_hs[_indices]
    _wave_dir = _wave_dir[_indices]

    t0_x_at_zero = find_x_at_zero(_bathy_t0, bathy_xb, i)
    t1_x_at_zero = find_x_at_zero(_bathy_t1, bathy_xb, i)
    change = (t1_x_at_zero - t0_x_at_zero) / 100

    if change > 1 || change < -1
        println(i, "    ", change)
        exit()
    end

    if change > 0
        class = (1, 0, 0)  # (offshore, landward, none)
        global n_offshore
        n_offshore += 1
    elseif change < 0
        class = (0, 1, 0)
        global n_landward
        n_landward += 1
    else
        class = (0, 0, 1)
        global n_none
        n_none += 1
    end

    _wave_Tm = _wave_Tm ./ 100
    _wave_hs = _wave_hs ./ 10
    _wave_dir = _wave_dir ./ 100

    input_row = vcat(
        mean(_wave_Tm),
        minimum(_wave_Tm),
        maximum(_wave_Tm),
        std(_wave_Tm),
        mean(_wave_hs),
        minimum(_wave_hs),
        maximum(_wave_hs),
        std(_wave_hs),
        mean(_wave_dir),
        minimum(_wave_dir),
        maximum(_wave_dir),
        std(_wave_dir)
    )

    push!(all_tm, vcat(mean(_wave_Tm), minimum(_wave_Tm), maximum(_wave_Tm), std(_wave_Tm))...)
    push!(all_hs, vcat(mean(_wave_hs), minimum(_wave_hs), maximum(_wave_hs), std(_wave_hs))...)
    push!(all_dir, vcat(mean(_wave_dir), minimum(_wave_dir), maximum(_wave_dir), std(_wave_dir))...)

    day = Day(_bathy_t0, _bathy_t1, _wave_Tm, _wave_hs, _wave_dir, change, class, input_row)
    push!(days, day)
end

println("normalized(all_tm)     ", minimum(all_tm), "  ", maximum(all_tm))      # RAW: 0.007926686233555585     14.506896347331828
println("normalized(all_hs)     ", minimum(all_hs), "  ", maximum(all_hs))      # RAW: 0.0007081265658494296    2.6867510628524585
println("normalized(all_dir)     ", minimum(all_dir), "  ", maximum(all_dir))   # RAW: -13.212138751720701      17.384923209385946

println("n_offshore: ", n_offshore, ", n_landward: ", n_landward, ", n_none: ", n_none)  # n_offshore: 274, n_landward: 243, n_none: 696

function evaluate(ind::CGPInd, days::Vector{Any}, cfg::NamedTuple)
    errors = Vector{Float64}(undef, length(days))
    reset!(ind)
    for di in 1:length(days)
        day = days[di]

        new_change = process(ind, day.input_row)[1]
        error = (new_change*100 - day.shoreline_change*100)^2

        errors[di] = error
    end
    [mean(errors) * -1]
end

function make_configs(template)
    lst_n_cols = collect(20:10:50)
    lst_n_pop = collect(40:20:200)
    lst_f_recur = collect(0:0.1:1)
    println(length(lst_n_cols))
    println(length(lst_n_pop))
    println(length(lst_f_recur))
    configs = []
    for n_pop in lst_n_pop
        n_elite = 0.1 * n_pop
        for n_cols in lst_n_cols
            for f_recur in lst_f_recur
                config = deepcopy(template)
                two_arity = falses(length(config["functions"]))
                functions = Array{Function}(undef, length(config["functions"]))
                for i in eachindex(config["functions"])
                    fname = config["functions"][i]
                    if CGPFunctions.arity[fname] == 2
                        two_arity[i] = true
                    end
                    functions[i] = eval(Meta.parse(string("CGPFunctions.", fname)))
                end
                config["two_arity"] = two_arity
                config["functions"] = functions

                config["n_population"] = n_pop
                config["n_elite"] = n_elite
                config["columns"] = n_cols
                config["recur"] = f_recur
                config["id"] = string("auto_", "npop", n_pop, "ncols", n_cols, "frecur", f_recur)
                push!(configs, deepcopy(config))
            end
        end
    end
    configs
end

include("cfg_template.jl")
configs = make_configs(template)
println(length(configs))

cfg = configs[1]
mutate(i::CGPInd) = goldman_mutate(cfg, i)

for config in configs
    global cfg
    cfg = (; (Symbol(k)=>v for (k, v) in config)...)
    println(cfg)
    fit(i::CGPInd) = evaluate(i, days, cfg)
    e = CGPEvolution(cfg, fit)
    run!(e)
end

max_stars_repo_name: mahmoud-al-najar/cgp-morphcast

max_issues_repo_name: mahmoud-al-najar/cgp-morphcast

max_forks_repo_name: mahmoud-al-najar/cgp-morphcast


text:
module TestExamplesTutorialMissings
include("../examples/tutorial_missings.jl")
include("preamble.jl")

slow_xf_fullextrema = Zip(Count(), NotA(Missing)) |>
    Zip(xf_scanext(>), xf_scanext(<))

slow_xf_argextrema = slow_xf_fullextrema |>
    Map() do ((argmin, min), (argmax, max))
        (argmin, argmax)
    end

@testset "slow compilation" begin
    @test mapfoldl(slow_xf_fullextrema, right, [1.0, 3.0, -1.0, missing, 2.0]) ===
        ((3, -1.0), (2, 3.0))

    @test mapfoldl(slow_xf_argextrema, right, [1.0, 3.0, -1.0, missing, 2.0]) ===
        (3, 2)

    # Excluding `missing` from the input vector seems to slow down
    # compilation a lot.
    slow_test("slow compilation test with xf_fullextrema", 20) do
        @test mapfoldl(slow_xf_fullextrema, right, [1.0, 3.0, -1.0, 2.0]) ==
            ((3, -1.0), (2, 3.0))
    end
    slow_test("slow compilation test with xf_argextrema", 20) do
        @test mapfoldl(slow_xf_argextrema, right, [1.0, 3.0, -1.0, 2.0]) ==
            (3, 2)
    end
end

end  # module


text:
#########################################################
# Elementary action/properties of permutation actions

@inline @inbounds function Base.:^(n::Integer, p::Perm)
    if 1 <= n <= length(p.d)
        return oftype(n, p.d[n])
    end
    return n
end

@inline Base.:^(v::Tuple, p::Perm) = ntuple(i->v[i^p], length(v))
@inline function Base.:^(v::AbstractVector, p::Perm)
    res = similar(v);
    @inbounds for i in eachindex(v)
        res[i] = v[i^p]
    end
    return res
end

@inline fixes(p::GroupElem, pt, op=^) = op(pt, p) == pt
@inline fixes(p::Perm, v::AbstractVector, op=^) =
    all( i-> v[i] == v[op(i, p)], eachindex(v))

@inline Base.isone(p::Perm) = all(i->i^p == i, 1:degree(p))

fixedpoints(p::Perm, range=1:degree(p), op=^) = [i for i in range if fixes(p, i, op)]
nfixedpoints(p::Perm, range=1:degree(p), op=^) = count(i->fixes(p, i, op), range)

for (fname, findname) in [(:firstmoved, :findfirst), (:lastmoved, :findlast)]
    @eval begin
        function $fname(p::Generic.Perm{I}, op=^) where I
            k = $findname(i -> i != p.d[i], eachindex(p.d))
            k === nothing && return zero(I)
            return I(k)
        end
    end
end

#########################################################
# Misc functions that should go to AbstractAlgebra

import Base: one, conj
import AbstractAlgebra: degree

AbstractAlgebra.degree(p::Perm{I}) where I = I(length(p.d))

function Generic.emb(p::Generic.Perm{I}, n) where I
    return Generic.emb!(Perm(I(n)), p, 1:degree(p))
end

"""
    conj!(out::Perm, h::Perm, g::Perm)
Computes the conjugation action of `g` on `h` and stores the result in `out`.

The action is understood to be `h → g^-1*h*g`. `out` will be unaliased, if necessary.
"""
function Base.conj!(out::Perm, h::Perm, g::Perm)
    if out === h
        out = deepcopy(out)
    end
    @inbounds for i in 1:degree(g)
        out[g[i]] = g[h[i]]
    end
    return out
end

Base.conj(h::GroupElem, g::GroupElem) = conj!(h, h, g)
Base.:(^)(h::Perm, g::Perm) = conj(h,g)

AbstractAlgebra.degree(S::Generic.SymmetricGroup) = S.n
function AbstractAlgebra.gens(G::Generic.SymmetricGroup{I}) where I
    a = one(G)
    a.d[1], a.d[2] = 2, 1
    b = Perm(circshift(I(1):degree(G), -1))
    return [a, b]
end


text:
export reflect

"""
    reflect(P::LazySet{N}) where {N<:Real}

Concrete reflection of a convex set `P`, resulting in the reflected set `-P`.

### Note

This function requires that the list of constraints of the set `P` is
available, i.e. such that it can be written as
``P = \\{z ∈ ℝⁿ: ⋂ sᵢᵀz ≤ rᵢ, i = 1, ..., N\\}.``

This function can be used to implement the alternative definition of the
Minkowski Difference, which writes as
```math
A ⊖ B = \\{a − b | a ∈ A, b ∈ B\\} = A ⊕ (-B)
```
by calling `minkowski_sum(A, reflect(B))`.
"""
function reflect(P::LazySet)

    @assert applicable(constraints_list, P)  "this function " *
        "requires that the list of constraints is available, but it is not; " *
        "if the set is bounded, try overapproximating with an `HPolytope` first"

    F,g = tosimplehrep(P)
    T = isbounded(P) ? HPolytope : HPolyhedron
    return T(-F, g)
end


text:
struct Fab
    f::Function
    a::Float64
    b::Float64
end


text:
using DataFrames, Distributions, Random, LinearAlgebra, CSV, Plots; pyplot()
Random.seed!(0)

data = CSV.read("../data/L1L2data.csv")
xVals, yVals, n = data.X, data.Y, size(data)[1]

N = 10^6
alphaMin, alphaMax, betaMin, betaMax = 0, 5, 0, 5

alpha1, beta1, alpha2, beta2, bestL1Cost, bestL2Cost = 0.0,0.0,0.0,0.0,Inf,Inf
for _ in 1:N
    rAlpha,rBeta=rand(Uniform(alphaMin,alphaMax)),rand(Uniform(betaMin,betaMax))
    L1Cost = norm(rAlpha .+ rBeta*xVals - yVals,1)
    if L1Cost < bestL1Cost
        global alpha1 = rAlpha
        global beta1 = rBeta
        global bestL1Cost = L1Cost
    end
    L2Cost = norm(rAlpha .+ rBeta*xVals  - yVals)
    if L2Cost < bestL2Cost
        global alpha2 = rAlpha
        global beta2 = rBeta
        global bestL2Cost = L2Cost
    end
end

println("L1 line: $(round(alpha1,digits = 2)) + $(round(beta1,digits = 2))x")
println("L2 line: $(round(alpha2,digits = 2)) + $(round(beta2,digits = 2))x")

d = yVals - (alpha2 .+ beta2*xVals)
rectangle(x, y, d) = Shape(x .- [0,d,d,0,0], y .- [0,0,d,d,0])

p1 = scatter(xVals,yVals, c=:black, ms=5, label="")
p1 = plot!([0,10],[alpha1, alpha1 .+ beta1*10], c=:blue,label="L1 minimized")
for i in 1:n
    x,y = xVals[i],yVals[i]
    p1 = plot!([x, x], [y, alpha1 .+ beta1*x],color="black", label="")
end

p2 = scatter(xVals,yVals, c=:black, ms=5, label="")
p2 = plot!([0,10],[alpha2, alpha2 .+ beta2*10],c=:red,label="L2 minimized")
for i in 1:n
    x,y = xVals[i],yVals[i]
	p2 = plot!(rectangle(x,y,d[i]), fc=:gray, fa=0.5, label="")
end

p3 = scatter(xVals,yVals, c=:black, ms=5, label="")
p3 = plot!([0,10],[alpha1, alpha1 .+ beta1*10], c=:blue, label="L1 minimized")
p3 = plot!([0,10],[alpha2, alpha2 .+ beta2*10], c=:red, label="L2 minimized")

plot(p1, p2, p3, layout = (1,3),
	ratio=:equal, xlims=(0,10), ylims=(0,10),
	legend=:topleft, size=(1200, 400),
    	xlabel = "x", ylabel = "y")

max_stars_repo_name: Yoshinobu-Ishizaki/StatsWithJuliaBook

max_issues_repo_name: Yoshinobu-Ishizaki/StatsWithJuliaBook

max_forks_repo_name: Yoshinobu-Ishizaki/StatsWithJuliaBook


text:
using MAT
using Random
include("../src/HyperLocal.jl")
include("../include/FlowSeed.jl")
#
## Read in the matrix, if it isn't read in already
using SparseArrays
s = time()
@time M = matread("../data/processed/stackoverflow_answer_H.mat")
LabelMatrix = M["LabelMatrix"]
LabelNames = M["LabelNames"]
MainLabels = M["MainLabels"]
H = M["H"]
order = round.(Int64,vec(sum(H,dims=2)))
d = vec(sum(H,dims=1))
volA = sum(d)
m,n = size(H)
Ht = sparse(H')
toc = time()-s
println("Done loading things into memory in $toc seconds.")

##
Tconds = Vector{Float64}()
Tsizes = Vector{Int64}()
labels = Vector{Int64}()
for lab = 1:length(MainLabels)
    label = MainLabels[lab]
    T = findnz(LabelMatrix[:,label])[1]
    nT = length(T)
    condT, volT, cutT = tl_cond(H,T,d,1.0,volA,order)
    if condT < .2
        # println("$nT \t $condT \t"*LabelNames[label])
        push!(Tconds,condT)
        push!(labels,label)
        push!(Tsizes,nT)
    end
end

using Statistics
ntimes = 1
epsilon = 1.0
delta = 5000.0
outputmat = "Output_Stack/Set45_$(epsilon)_$(delta).mat"

mat = matread(outputmat)
hl_time = round.(mean(mat["hl_time"],dims = 2),digits = 2)
sl_time = round.(mean(mat["sl_time"],dims = 2),digits = 2)
zl_time = round.(mean(mat["zl_time"],dims = 2),digits = 2)
hl_size= round.(mean(mat["hl_size"],dims = 2),digits = 2)
r = round.(mean(mat["r_f1"],dims = 2),digits = 2)
hl = round.(mean(mat["hl_f1"],dims = 2),digits = 2)
b1 = round.(mean(mat["b1_f1"],dims = 2),digits = 2)
b2 = round.(mean(mat["b2_f1"],dims = 2),digits = 2)
nS = round.(mean(mat["newS"],dims = 2),digits = 2)
sl = round.(mean(mat["sl_f1"],dims = 2),digits = 2)
zl = round.(mean(mat["zl_f1"],dims = 2),digits = 2)

for i = 1:length(labels)
    lab = labels[i]
    println(" $(hl[i]) \t $(sl[i]) \t $(zl[i]) \t $(b1[i]) \t $(b2[i])\t $(r[i]) \t $(hl_time[i]) \t $(sl_time[i]) \t $(zl_time[i])\t $(Tconds[i])"*"\t"*LabelNames[lab])
end

max_stars_repo_name: nveldt/HypergraphFlowClustering

max_issues_repo_name: nveldt/HypergraphFlowClustering

max_forks_repo_name: nveldt/HypergraphFlowClustering


text:
export get_projector

function get_projector(constraint,comp_grid,special_operator_list::Array{String,1},A,TD_n::Tuple,TF::DataType)

  if constraint.set_type == "bounds"
    if constraint.app_mode[1] in ["matrix","tensor"]
      if constraint.TD_OP in special_operator_list
        P = x -> copyto!(x,A'*project_bounds!(A*x,constraint.min,constraint.max))
      else
        P = x -> project_bounds!(x,constraint.min,constraint.max)
      end
    else
      if constraint.TD_OP in special_operator_list
        P = x -> copyto!(x,A'*project_bounds!(reshape(A*x,TD_n),constraint.min,constraint.max,constraint.app_mode))
      else
        P = x -> project_bounds!(reshape(x,TD_n),constraint.min,constraint.max,constraint.app_mode)
      end
    end
  end

  if constraint.set_type == "prox_l1"
    if constraint.TD_OP in special_operator_list
      P = x -> copyto!(x,.5f0*A'*prox_l1!(A*x,constraint.max))
    else
      P = x -> prox_l1!(x,constraint.max)
    end
  end

  if constraint.set_type == "l1"
    if constraint.TD_OP in special_operator_list
      P = x -> copyto!(x,A'*project_l1_Duchi!(A*x,constraint.max))
    else
      P = x -> project_l1_Duchi!(x,constraint.max)
    end
  end

  if constraint.set_type == "l2"
    if constraint.TD_OP in special_operator_list
      P = x -> copyto!(x,A'*project_l2!(A*x,constraint.max))
    else
      P = x -> project_l2!(x,constraint.max)
    end
  end

  if constraint.set_type == "annulus"
    if constraint.TD_OP in special_operator_list
      P = x -> copyto!(x,A'*project_annulus!(A*x,constraint.min,constraint.max))
    else
      P =  x -> project_annulus!(x,constraint.min,constraint.max)
    end
  end

  if constraint.set_type == "subspace"
    if constraint.app_mode[1]  in ["matrix","tensor"]
      P = x -> project_subspace!(x,constraint.custom_TD_OP[1],constraint.custom_TD_OP[2])
    else
      P = x -> project_subspace!(reshape(x,comp_grid.n),constraint.custom_TD_OP[1],constraint.custom_TD_OP[2],constraint.app_mode)
    end
  end

  if constraint.set_type == "nuclear"
    if constraint.TD_OP in special_operator_list
      P = x -> copyto!(x,A'*project_nuclear!(reshape(A*x,TD_n),constraint.max,constraint.app_mode))
    else
      P = x -> project_nuclear!(reshape(x,TD_n),constraint.max,constraint.app_mode)
    end
  end

  if constraint.set_type == "rank"
    if constraint.TD_OP in special_operator_list
      P = x -> copyto!(x,A'*project_rank!(reshape(A*x,TD_n),constraint.max,constraint.app_mode[2]))
    else
      P = x -> project_rank!(reshape(x,TD_n),constraint.max,constraint.app_mode)
    end
  end

  if constraint.set_type == "histogram" #1D histogram only
    if constraint.TD_OP in special_operator_list
      P = x -> copyto!(x,A'*project_histogram_relaxed!(A*x,constraint.min,constraint.max)) #min and max should be vector valued
    else
      P = x -> project_histogram_relaxed!(x,constraint.min,constraint.max) #min and max should be vector valued
    end
  end

  if constraint.set_type == "cardinality"
    if constraint.app_mode[1]  in ["matrix","tensor"]
      if constraint.TD_OP in special_operator_list
        P = x -> copyto!(x,A'*project_cardinality!(A*x,convert(Integer,constraint.max)))
      else
        P = x -> project_cardinality!(x,convert(Integer,constraint.max))
      end
    else
      if constraint.TD_OP in special_operator_list
        P = x -> copyto!(x,A'*project_cardinality!(reshape(A*x,TD_n),convert(Integer,constraint.max)),constraint.app_mode)
      else
        P = x -> project_cardinality!(reshape(x,TD_n),convert(Integer,constraint.max),constraint.app_mode)
      end
    end
  end

return P
end

max_stars_repo_name: slimgroup/SetIntersectionProjection

max_issues_repo_name: slimgroup/SetIntersectionProjection

max_forks_repo_name: slimgroup/SetIntersectionProjection


text:
module Instruments

import Base.ref
export VanillaPayoff, ForwardPayoff, VanillaOption, AmericanOption

abstract Payoff
abstract Instrument

type VanillaPayoff <: Payoff
	typ
	strike
end

type ForwardPayoff <: Payoff
	typ
	strike
end

type StraddlePayoff <: Payoff
	typ
	strike
end

ref(p::VanillaPayoff, s::Real) = typ == :call ? max(s - p.strike, 0) : max(p.strike - s, 0)
ref(p::ForwardPayoff, s::Real) = typ == :long ? (s - p.strike) : (p.strike - s)
ref(p::StraddlePayoff, s::Real) = typ == :long ? max(s - p.strike, 0) + max(p.strike - s, 0) : -max(s - p.strike, 0) - max(p.strike - s, 0)

type VanillaOption
	payoff
	excercise
end

type AmericanOption
	payoff
	excercise
end

end


text:
# Load the required packages
begin
    using Base.Iterators
    using BNPRegressionGGA2021
    using CSV
    using DataFrames
    using Distributions
    using LinearAlgebra
    using Random
    using RCall
    const BNP = BNPRegressionGGA2021
end;

# Define the DGP
function simulate_sample(N0, N1)
    # Set the conditional density
    dy(xc) = Normal(0.2 * xc^3, 0.5)

    # Simulate a sample
    x0 = collect(LinRange(-3, 3, N0))
    y0 = [rand(dy(xc)) for xc in x0]
    X0 = rcopy(R"cbind(1, splines::bs($x0, df = 6, Boundary.knots = c(-3, 3)))")

    # Generate the grid points
    y_grid = LinRange(-3.0, 3.0, N1)
    x_grid = [-1, 1]
    y1, x1 = Iterators.product(y_grid, x_grid) |> x -> zip(x...) .|> collect
    X1 = rcopy(R"cbind(1, predict(splines::bs($x0, df = 6, Boundary.knots = c(-3, 3)), $x1))")
    return dy, y0, x0, X0, y1, x1, X1
end

# Simulate the data
begin
    Random.seed!(1);
    N0, N1 = 500, 50;
    dy, y0, x0, X0, y1, x1, X1 = simulate_sample(N0, N1);
end;

# Fit the model
begin
    mapping = [[1], collect(2:size(X0, 2))];
    m = BNP.DGSBPNormal(; y0, X0, y1, X1, mapping, update_γ = false);
    chainf, chainβ, chaing = BNP.sample!(m; mcmcsize = 4000);
end;

# Collect the estimated densities
df = DataFrame(y = y1, x = x1, fh = mean(chainf), f0 = pdf.(dy.(x1), y1))

# Save the results
CSV.write("data/walker-example.csv", df)

# Create a plot
R"""
# Munge the data
df <-
    readr::read_csv("data/walker-example.csv") |>
    dplyr::mutate(x = factor(x)) |>
    tidyr::pivot_longer(c(fh, f0))

# Plot the conditional posterior predictive densities
p <-
    df |>
    ggplot2::ggplot(ggplot2::aes(x = y, y = value, color = name)) +
    ggplot2::facet_grid(row = ggplot2::vars(x)) +
    ggplot2::geom_line()

# Save the plot in png format
"figures/walker-example.png" |>
    ggplot2::ggsave(plot = p, width = 5 * 8 / 9, height = 3)
"""


text:
"""
An approach that, given an instance, will merely convert that instance into a
`MathOptInterface.ModelLike`. For instance, this could encapsulate a
`AbstractOptimizer`, which ingests an instance and passes that instance to the
underlying solver through the `MathOptInterface` API.

    MOIBasedApproach{T}() where {T <: MOI.AbstractOptimizer}
    MOIBasedApproach{T}(factory::Function) where {T <: MOI.AbstractOptimizer}

The underlying `MathOptInterface` will either be inferred from `T`, or you can
optionally pass a `factory` function to the constructor. This `factory` must
respect the signature

    (instance::AbstractInstance, config::Dict{String,Any}) --> MathOptInterface.AbstractOptimizer.

Some pieces of the [functional interface](@ref functional_interface) are
implemented:

    function build_model(
        approach::MOIBasedApproach{T},
        inst::AbstractMOIBackedInstance,
        config::Dict{String,Any} = Dict{String,Any}(),
    ) where {T <: MathOptInterface.AbstractOptimizer}
"""
struct MOIBasedApproach{T<:MOI.AbstractOptimizer} <: AbstractApproach
    optimizer_factory::Union{Nothing,Function}
end
function MOIBasedApproach{T}() where {T<:MOI.AbstractOptimizer}
    return MOIBasedApproach{T}(nothing)
end

_get_factory(::Nothing, ::Type{T}) where {T<:MOI.AbstractOptimizer} = T
_get_factory(f::Function, ::Type{T}) where {T<:MOI.AbstractOptimizer} = f

function get_factory(
    approach::MOIBasedApproach{T},
) where {T<:MOI.AbstractOptimizer}
    return _get_factory(approach.optimizer_factory, T)
end

function build_model(
    approach::MOIBasedApproach{T},
    inst::AbstractMOIBackedInstance,
    config::Dict{String,Any} = Dict{String,Any}(),
) where {T<:MOI.AbstractOptimizer}
    factory = get_factory(approach)
    model = MOI.instantiate(factory)
    @assert model isa T
    MOI.copy_to(model, get_moi_model(inst), copy_names = false)
    return MOIModel(model)
end


text:
# Author: Pedro Henrique Pereira Braga
# Date: 2019-01-26
# Subject: Code to output HTML, PDF and MD from a Markdown document

# Add dependencies
using Pkg; Pkg.add.(["Plots", "DSP"])
using Pkg; Pkg.add.("Weave")
using Pkg; Pkg.add.(["Gadfly"])
using Pkg; Pkg.add.(["LaTeXStrings"])

#Use Weave
using Weave

# Weave the Document
## Output it to Markdown
weave(joinpath("", "WeavingDocumentsJl.jmd"),
      informat="markdown",
      out_path = :pwd, # This outputs the document to your current working directory
      doctype = "pandoc")

## to HTML
weave(joinpath("", "WeavingDocumentsJl.jmd"),
      informat="markdown",
      out_path = :pwd,
      doctype = "md2html")

## to PDF
weave(joinpath("", "WeavingDocumentsJl.jmd"),
      informat="markdown",
      out_path = joinpath("", "pdf"),
      doctype = "pandoc2pdf")

max_stars_repo_name: pedrohbraga/WeavingDocumentsInJulia

max_issues_repo_name: pedrohbraga/WeavingDocumentsInJulia

max_forks_repo_name: pedrohbraga/WeavingDocumentsInJulia


text:
include("../MNIST/MNIST.jl")

module MLPDemo

using MNIST, GenEvo

export runDemo

runDemo(ctx) = scope!(ctx) do this
    const batchSize = 1000

    println("Setting up training data.")

    data = loadSubset(ctx)

    println("MNIST Training sample count: $(numTrain(data))")
    println("MNIST Test sample count: $(numTest(data))")
    println("MNIST Feature size: $(featureSize(data))")

    validateSampleCount = numTrain(data) % batchSize
    if validateSampleCount < 10
        validateSampleCount += batchSize
    end
    const trainSampleCount = numTrain(data) - validateSampleCount

    println("Training sample count: $trainSampleCount")
    println("Validate sample count: $validateSampleCount")

    const trainImages = data.trainImages[seq(ctx, 0, trainSampleCount - 1), :]
    const trainLabels = data.trainLabels[seq(ctx, 0, trainSampleCount - 1), :]
    @assert dims(trainImages, 0) == trainSampleCount

    const validateImages = data.trainImages[seq(ctx, trainSampleCount, numTrain(data) - 1), :]
    const validateLabels = data.trainLabels[seq(ctx, trainSampleCount, numTrain(data) - 1), :]
    @assert dims(validateImages, 0) == validateSampleCount

    batcher = ForwardBatcher(ctx, trainImages, trainLabels, batchSize)

    println("Creating Network.")
    network = FFMLP(ctx, [featureSize(data), 100, 50, numClasses(data)])
    println(network)

    bestMSE = 2f0
    function fitness(arr)
        setWeights!(network, arr)
        batchInputs, batchTargets = current(batcher)
        output = predict!(network, batchInputs)
        mse = calculateMSE(ctx, output, batchTargets)
        bestMSE = min(bestMSE, mse)
        mse
    end

    println("Initializing algorithms.")

    const synth = Synthesizer(ctx)
    define!(synth, ScalarRule(count = weightCount(network), min = -0.7f0, max = 0.7f0, variationCount = 1))

    const algo = CrossEntropy(
        ctx,
        dnaSize!(synth),
        CalculateComparer(fitness),
        dna -> decodeAsContextArray(synth, dna),
        populationSize = 30,
        mutationChance = 0.03f0,
        mutationStrength = 0.1f0,
        selectionStdDev = 0.3f0,
        keepElitesRate = 0.05f0)

    next!(batcher)
    start!(algo)

    println("Starting optimization loop.")
    const maxIterations = 300
    t = time()
    sec = @elapsed for i in 1:maxIterations
        next!(batcher)
        lastTime = @elapsed step!(algo)
        now = time()
        if now - t >= 5.0
            println(i, ". ",lastTime," : " ,bestMSE)
            t = now
        end
        sync(ctx)
    end

    println("\nDONE.\n\nResults:\n")

    trainOutput = predict!(network, trainImages)

    println("Training set:")
    println("Accuracy on training data: $(accuracy(ctx, trainOutput, trainLabels))")

    validateOutput = predict!(network, validateImages)

    println("Validation set:")
    println("Accuracy on testing  data: $(accuracy(ctx, validateOutput, validateLabels))")

    println("Training time: $sec seconds")
end

function calculateMSE(ctx, output, target)
    scope!(ctx) do this
        dif = output .- target
        sq = dif .* dif
        s = sum(Float32, sq)
        s / elements(output)
    end
end

function accuracy(ctx, predicted, target)
    scope!(ctx) do this
        pMaxArray, pMaxIndex = imax(predicted, 1)
        tMaxArray, tMaxIndex = imax(target, 1)
        (100.0f0 * count(Float32, pMaxIndex .== tMaxIndex)) / elements(tMaxIndex)
    end
end

end

max_stars_repo_name: unbornchikken/julia-ml-proto

max_issues_repo_name: unbornchikken/julia-ml-proto

max_forks_repo_name: unbornchikken/julia-ml-proto


text:
"""
This module contains transcriptions of the tables from
[Stancil 1994](https://ui.adsabs.harvard.edu/abs/1994ApJ...430..360S/abstract), who calculated ff
and bf absorption coefficients for H₂⁺ and He₂⁺.

It also contains pre-constructed interpolation objects.
"""
module Stancil1994

using Interpolations: LinearInterpolation, Linear
Ts_He2plus      = [4200,6300,8400,12600,16800,25200,33600,50400]
K_He2plus_vals = 1e21 .* [0.3606, 2.6330, 7.1409, 21.097, 39.842, 87.960, 147.41, 293.57]
Ts_H2plus       = [3150,4200,5040,6300,8400,12600,16800,25200]
K_H2plus_vals  = 1e19 .* [0.9600,9.7683,29.997,89.599,265.32,845.01,1685.3,4289.5]

#Qs not used, copied for completeness
#Qs_He2plus = [1166.8,2591.9,4372.1,8078.8,11353.,16288.,19652.,23818.]
#Qs_H2plus  = [142.07,247.23,358.74,568.37,999.84,1968.2,2877.3,4313.9]

# they are provided in nm, we convert them to Å
# for some reason, the wavelenth table for H₂⁺ contains 295 nm, but that for He₂⁺ doesn't
λs_H2plus_ff = 10.0 .* [70, 80, 90, 100, 110, 120, 130, 140, 150, 160, 170, 180, 190, 200, 210, 220,
                        230, 240, 250, 260, 270, 280, 290, 295, 300, 350, 400, 450, 500, 600, 700,
                        800, 900, 1000, 2000, 3000, 4000, 5000, 11000, 15000, 20000]
λs_He2plus_ff = 10.0 .* [70, 80, 90, 100, 110, 120, 130, 140, 150, 160, 170, 180, 190, 200, 210,
                         220, 230, 240, 250, 260, 270, 280, 290, 300, 350, 400, 450, 500, 600, 700,
                         800, 900, 1000, 2000, 3000, 4000, 5000, 11000, 15000, 20000]
# for bf absorption, the wavelenths are the same for H₂⁺ and He₂⁺
λs_bf = 10.0 .* [50, 60, 70, 80, 90, 100, 110, 120, 130, 140, 150, 160, 170, 180, 190, 200, 210,
                 220, 230, 240, 250, 260, 270, 280, 290, 300, 350, 400, 450, 500, 600, 700, 800,
                 900, 1000, 2000, 3000, 4500, 5000, 10000, 15000, 20000]

# The "cross-sections" used in this function have units of cm^-5
σ_He2plus_ff_table = [7.08e-3 6.05e-3 5.45e-3 4.78e-3 4.41e-3 3.99e-3 3.75e-3 3.50e-3
                      0.0103  8.76e-3 7.84e-3 6.79e-3 6.20e-3 5.52e-3 5.15e-3 4.74e-3
                      0.0134  0.0113  0.0101  8.66e-3 7.87e-3 6.97e-3 6.48e-3 5.91e-3
                      0.0165  0.0139  0.0124  0.0106  9.65e-3 8.52e-3 7.89e-3 7.19e-3
                      0.0196  0.0165  0.0147  0.0126  0.0114  0.0101  9.30e-3 8.47e-3
                      0.0226  0.0190  0.0169  0.0145  0.0131  0.0116  0.0107  9.72e-3
                      0.0259  0.0217  0.0193  0.0166  0.0150  0.0132  0.0122  0.0111
                      0.0288  0.0242  0.0215  0.0185  0.0167  0.0147  0.0136  0.0124
                      0.0317  0.0267  0.0237  0.0204  0.0185  0.0163  0.0151  0.0137
                      0.0347  0.0292  0.0260  0.0223  0.0203  0.0179  0.0165  0.0151
                      0.0377  0.0317  0.0283  0.0243  0.0220  0.0195  0.0180  0.0164
                      0.0405  0.0341  0.0304  0.0262  0.0237  0.0210  0.0195  0.0178
                      0.0433  0.0365  0.0326  0.0280  0.0254  0.0225  0.0209  0.0191
                      0.0460  0.0388  0.0346  0.0298  0.0271  0.0240  0.0223  0.0204
                      0.0487  0.0411  0.0367  0.0316  0.0288  0.0255  0.0237  0.0217
                      0.0514  0.0434  0.0387  0.0335  0.0304  0.0270  0.0251  0.0230
                      0.0540  0.0457  0.0408  0.0352  0.0321  0.0285  0.0266  0.0244
                      0.0567  0.0479  0.0429  0.0371  0.0338  0.0301  0.0280  0.0257
                      0.0593  0.0502  0.0449  0.0389  0.0354  0.0316  0.0294  0.0271
                      0.0619  0.0524  0.0469  0.0407  0.0371  0.0331  0.0308  0.0284
                      0.0645  0.0546  0.0489  0.0424  0.0387  0.0346  0.0323  0.0297
                      0.0669  0.0567  0.0508  0.0441  0.0403  0.0360  0.0336  0.0310
                      0.0693  0.0587  0.0527  0.0457  0.0418  0.0374  0.0350  0.0323
                      0.0716  0.0607  0.0545  0.0474  0.0433  0.0388  0.0363  0.0335
                      0.0825  0.0702  0.0632  0.0551  0.0506  0.0455  0.0427  0.0396
                      0.0934  0.0798  0.0719  0.0630  0.0580  0.0524  0.0494  0.0460
                      0.1036  0.0887  0.0802  0.0706  0.0651  0.0591  0.0558  0.0522
                      0.1122  0.0964  0.0873  0.0771  0.0714  0.0650  0.0615  0.0577
                      0.1280  0.1106  0.1006  0.0894  0.0831  0.0762  0.0724  0.0683
                      0.1446  0.1255  0.1146  0.1024  0.0956  0.0882  0.0841  0.0797
                      0.1573  0.1371  0.1257  0.1129  0.1058  0.0980  0.0938  0.0892
                      0.1673  0.1465  0.1347  0.1215  0.1143  0.1063  0.1020  0.0973
                      0.1787  0.1570  0.1448  0.1312  0.1237  0.1155  0.1111  0.1064
                      0.2615  0.2361  0.2221  0.2067  0.1984  0.1894  0.1847  0.1796
                      0.3078  0.2825  0.2687  0.2537  0.2456  0.2371  0.2325  0.2278
                      0.3556  0.3301  0.3162  0.3014  0.2934  0.2850  0.2806  0.2760
                      0.4175  0.3909  0.3766  0.3612  0.3531  0.3445  0.3400  0.3353
                      0.5835  0.5607  0.5487  0.5360  0.5294  0.5225  0.5190  0.5153
                      0.6620  0.6406  0.6293  0.6176  0.6114  0.6051  0.6018  0.5985
                      0.7563  0.7359  0.7252  0.7140  0.7083  0.7023  0.6992  0.6961]

σ_He2plus_bf_table = [5.99e-6 1.15e-5 1.50e-5 1.86e-5 2.62e-5 2.18e-5 2.25e-5 2.31e-5
                      4.31e-2 8.49e-3 0.0112  0.0140  0.0153  0.0165  0.0171  0.0176
                      0.0523  0.0816  0.0973  0.1108  0.1166  0.1213  0.1232  0.1249
                      0.2503  0.2883  0.2956  0.2927  0.2876  0.2809  0.2770  0.2729
                      0.6587  0.6169  0.5716  0.5169  0.4869  0.4584  0.4446  0.4313
                      1.2197  1.0054  0.8757  0.7493  0.6879  0.6333  0.6081  0.5842
                      1.7949  1.3801  1.1630  0.9667  0.8758  0.7973  0.7617  0.7285
                      2.2739  1.6955  1.4083  1.1557  1.0409  0.9426  0.8984  0.8573
                      2.5948  1.9282  1.5999  1.3121  1.1816  1.0701  1.0200  0.9734
                      2.7699  2.0876  1.7438  1.4385  1.2986  1.1784  1.1241  1.0735
                      2.8293  2.1868  1.8478  1.5389  1.3950  1.2699  1.2129  1.1596
                      2.7978  2.2354  1.9185  1.6198  1.4775  1.3523  1.2948  1.2407
                      2.7012  2.2420  1.9582  1.6787  1.5418  1.4194  1.3625  1.3087
                      2.5652  2.2185  1.9747  1.7208  1.5924  1.4754  1.4203  1.3678
                      2.4108  2.1766  1.9764  1.7532  1.6360  1.5269  1.4750  1.4250
                      2.2515  2.1226  1.9663  1.7752  1.6702  1.5701  1.5217  1.4748
                      2.0956  2.0615  1.9468  1.7873  1.6944  1.6034  1.5586  1.5146
                      1.9478  1.9976  1.9221  1.7934  1.7127  1.6310  1.5899  1.5492
                      1.8099  1.9334  1.8946  1.7958  1.7275  1.6553  1.6182  1.5809
                      1.6822  1.8694  1.8646  1.7945  1.7382  1.6756  1.6425  1.6087
                      1.5643  1.8060  1.8322  1.7894  1.7449  1.6919  1.6628  1.6326
                      1.4557  1.7438  1.7985  1.7817  1.7488  1.7054  1.6805  1.6541
                      1.3561  1.6834  1.7641  1.7720  1.7506  1.7168  1.6962  1.6736
                      1.2648  1.6249  1.7291  1.7601  1.7497  1.7253  1.7088  1.6901
                      1.1815  1.5685  1.6935  1.7458  1.7457  1.7302  1.7177  1.7026
                      1.1056  1.5144  1.6579  1.7297  1.7391  1.7320  1.7234  1.7119
                      0.8192  1.2864  1.4950  1.6420  1.6905  1.7188  1.7266  1.7304
                      0.6390  1.1210  1.3675  1.5665  1.6448  1.7017  1.7233  1.7406
                      0.5156  0.9900  1.2562  1.4893  1.5890  1.6678  1.7003  1.7283
                      0.4267  0.8818  1.1552  1.4080  1.5221  1.6161  1.6566  1.6927
                      0.3145  0.7275  1.0010  1.2736  1.4051  1.5195  1.5711  1.6187
                      0.2506  0.6309  0.9016  1.1870  1.3314  1.4619  1.5226  1.5799
                      0.2079  0.5578  0.8202  1.1078  1.2583  1.3977  1.4638  1.5271
                      0.1772  0.4991  0.7498  1.0327  1.1842  1.3269  1.3955  1.4618
                      0.1551  0.4541  0.6941  0.9710  1.1219  1.2660  1.3360  1.4041
                      0.0716  0.2458  0.4036  0.6028  0.7194  0.8366  0.8958  0.9549
                      0.0454  0.1634  0.2743  0.4207  0.5041  0.5919  0.6368  0.6820
                      0.0344  0.1268  0.2151  0.3334  0.4015  0.4738  0.5110  0.5486
                      0.0290  0.1088  0.1861  0.2908  0.3516  0.4165  0.4501  0.4841
                      0.0134  0.0517  0.0900  0.1429  0.1742  0.2080  0.2257  0.2437
                      0.0104  0.0408  0.0711  0.1134  0.1384  0.1656  0.1798  0.1943
                      8.88e-3 0.0349  0.0611  0.0976  0.1194  0.1430  0.1554  0.1681]

σ_H2plus_bf_table = [7.34e-5 1.43e-4 2.04e-4 2.87e-4 3.89e-4 4.97e-4 5.49e-4 5.98e-4
                     0.0100  0.0150  0.0186  0.0230  0.0276  0.0319  0.0337  0.0353
                     0.1676  0.1965  0.2105  0.2215  0.2266  0.2246  0.2211  0.2163
                     0.8477  0.8199  0.7797  0.7183  0.6376  0.5477  0.5037  0.4622
                     2.1113  1.8166  1.6135  1.3823  1.1403  0.9157  0.8176  0.7313
                     3.4427  2.8069  2.4213  2.0136  1.6137  1.2616  1.1129  0.9845
                     4.3470  3.5155  3.0224  2.5070  2.0062  1.5685  1.3846  1.2262
                     4.6981  3.8841  3.3806  2.8402  2.3019  1.8203  1.6147  1.4358
                     4.6169  3.9763  3.5361  3.0358  2.5120  2.0239  1.8096  1.6202
                     4.2811  3.8840  3.5476  3.1272  2.6535  2.1858  1.9729  1.7809
                     3.8331  3.6850  3.4660  3.1434  2.7388  2.3082  2.1033  1.9138
                     3.3624  3.4344  3.3301  3.1098  2.7840  2.4020  2.2106  2.0287
                     2.9167  3.1670  3.1662  3.0438  2.7985  2.4707  2.2958  2.1244
                     2.5172  2.9031  2.9909  2.9573  2.7898  2.5175  2.3610  2.2019
                     2.1697  2.6538  2.8149  2.8600  2.7657  2.5495  2.4128  2.2682
                     1.8726  2.4243  2.6444  2.7573  2.7300  2.5682  2.4516  2.3222
                     1.6208  2.2160  2.4821  2.6522  2.6843  2.5725  2.4748  2.3597
                     1.4085  2.0289  2.3297  2.5473  2.6314  2.5653  2.4852  2.3839
                     1.2295  1.8617  2.1883  2.4453  2.5755  2.5517  2.4889  2.4012
                     1.0785  1.7129  2.0581  2.3479  2.5188  2.5347  2.4885  2.4144
                     0.9508  1.5804  1.9387  2.2553  2.4622  2.5145  2.4840  2.4230
                     0.8424  1.4623  1.8293  2.1680  2.4063  2.4920  2.4762  2.4277
                     0.7501  1.3568  1.7296  2.0865  2.3530  2.4701  2.4687  2.4327
                     0.6712  1.2626  1.6392  2.0118  2.3045  2.4525  2.4660  2.4433
                     0.6033  1.1783  1.5577  1.9443  2.2620  2.4415  2.4707  2.4626
                     0.5448  1.1029  1.4843  1.8834  2.2256  2.4369  2.4825  2.4901
                     0.3477  0.8233  1.2037  1.6442  2.0832  2.4343  2.5567  2.6400
                     0.2412  0.6469  1.0069  1.4533  1.9363  2.3659  2.5349  2.6658
                     0.1782  0.5283  0.8605  1.2951  1.7920  2.2630  2.4604  2.6222
                     0.1389  0.4487  0.7547  1.1780  1.6859  2.1963  2.4228  2.6179
                     0.0927  0.3396  0.6120  1.0109  1.5255  2.0874  2.3564  2.6031
                     0.0685  0.2749  0.5175  0.8894  1.3923  1.9700  2.2587  2.5323
                     0.0542  0.2335  0.4534  0.8033  1.2940  1.8800  2.1824  2.4760
                     0.0448  0.2041  0.4046  0.7334  1.2067  1.7876  2.0941  2.3965
                     0.0382  0.1813  0.3635  0.6701  1.1194  1.6816  1.9827  2.2828
                     0.0159  0.0901  0.1982  0.3951  0.7108  1.1448  1.3953  1.6590
                     0.0100  0.0596  0.1325  0.2699  0.4954  0.8132  1.0003  1.1999
                     6.88e-3 0.0425  0.0962  0.1994  0.3723  0.6216  0.7710  0.9325
                     6.41e-3 0.0400  0.0908  0.1889  0.3540  0.5932  0.7371  0.8932
                     3.56e-3 0.0229  0.0526  0.1110  0.2109  0.3582  0.4479  0.5463
                     2.50e-3 0.0161  0.0373  0.0790  0.1506  0.2567  0.3216  0.3929
                     1.90e-3 0.0123  0.0286  0.0607  0.1161  0.1982  0.2487  0.3042]

σ_H2plus_ff_table = [0.0174 0.0154 0.0142 0.0130 0.0116 0.0100 9.10e-3 8.08e-3
                     0.0280 0.0246 0.0227 0.0207 0.0184 0.0158 0.0143  0.0126
                     0.0394 0.0346 0.0319 0.0290 0.0257 0.0220 0.0199  0.0175
                     0.0514 0.0451 0.0416 0.0378 0.0336 0.0287 0.0259  0.0227
                     0.0640 0.0562 0.0519 0.0471 0.0418 0.0357 0.0322  0.0283
                     0.0770 0.0676 0.0624 0.0567 0.0504 0.0431 0.0389  0.0341
                     0.0903 0.0794 0.0733 0.0666 0.0592 0.0506 0.0457  0.0401
                     0.1040 0.0914 0.0843 0.0767 0.0682 0.0584 0.0527  0.0464
                     0.1177 0.1035 0.0956 0.0869 0.0773 0.0663 0.0599  0.0527
                     0.1317 0.1158 0.1070 0.0973 0.0866 0.0743 0.0672  0.0592
                     0.1456 0.1281 0.1184 0.1078 0.0960 0.0824 0.0746  0.0658
                     0.1597 0.1405 0.1299 0.1183 0.1054 0.0906 0.0821  0.0725
                     0.1737 0.1530 0.1414 0.1288 0.1149 0.0988 0.0896  0.0793
                     0.1877 0.1654 0.1529 0.1394 0.1243 0.1071 0.0972  0.0861
                     0.2017 0.1777 0.1644 0.1499 0.1338 0.1154 0.1048  0.0930
                     0.2156 0.1901 0.1759 0.1605 0.1433 0.1237 0.1125  0.0998
                     0.2294 0.2023 0.1873 0.1709 0.1527 0.1319 0.1201  0.1068
                     0.2431 0.2145 0.1987 0.1814 0.1622 0.1402 0.1277  0.1137
                     0.2568 0.2266 0.2099 0.1917 0.1716 0.1485 0.1354  0.1206
                     0.2703 0.2387 0.2211 0.2020 0.1809 0.1567 0.1430  0.1276
                     0.2836 0.2506 0.2322 0.2123 0.1902 0.1649 0.1506  0.1345
                     0.2969 0.2624 0.2433 0.2225 0.1994 0.1731 0.1582  0.1414
                     0.3100 0.2741 0.2542 0.2325 0.2086 0.1812 0.1657  0.1483
                     0.3165 0.2799 0.2596 0.2376 0.2131 0.1853 0.1695  0.1518
                     0.3230 0.2857 0.2650 0.2425 0.2177 0.1893 0.1732  0.1552
                     0.3858 0.3419 0.3176 0.2913 0.2621 0.2290 0.2103  0.1894
                     0.4451 0.3952 0.3677 0.3378 0.3048 0.2674 0.2463  0.2228
                     0.5011 0.4457 0.4152 0.3821 0.3456 0.3044 0.2812  0.2555
                     0.5539 0.4935 0.4603 0.4243 0.3848 0.3401 0.3150  0.2873
                     0.6511 0.5821 0.5442 0.5032 0.4583 0.4078 0.3795  0.3484
                     0.7388 0.6625 0.6207 0.5756 0.5262 0.4710 0.4402  0.4064
                     0.8186 0.7362 0.6911 0.6425 0.5895 0.5303 0.4975  0.4615
                     0.8918 0.8042 0.7563 0.7049 0.6488 0.5864 0.5518  0.5141
                     0.9596 0.8675 0.8173 0.7633 0.7047 0.6395 0.6036  0.5644
                     1.4600 1.3450 1.2830 1.2170 1.1460 1.0680 1.0260  0.9815
                     1.8050 1.6820 1.6170 1.5470 1.4740 1.3940 1.3510  1.3060
                     2.2000 2.0750 2.0090 1.9390 1.8650 1.7860 1.7450  1.7000
                     2.3130 2.1880 2.1220 2.0520 1.9790 1.9010 1.8590  1.8160
                     3.1800 3.0620 3.0000 2.9360 2.8690 2.7980 2.7610  2.7230
                     3.8110 3.7000 3.6420 3.5820 3.5200 3.4560 3.4220  3.3870
                     4.3230 4.2180 4.1640 4.1080 4.0500 3.9900 3.9580  3.9260]

# construct interpolations, converting cross-sections coeffs to cgs
σ_He2plus_ff = LinearInterpolation((λs_He2plus_ff, Ts_He2plus), σ_He2plus_ff_table .* 1e-39;
                                   extrapolation_bc=Linear())
σ_He2plus_bf = LinearInterpolation((λs_bf, Ts_He2plus), σ_He2plus_bf_table .* 1e-18;
                                   extrapolation_bc=Linear())
σ_H2plus_ff  = LinearInterpolation((λs_H2plus_ff, Ts_H2plus), σ_H2plus_ff_table .* 1e-39;
                                   extrapolation_bc=Linear())
σ_H2plus_bf  = LinearInterpolation((λs_bf, Ts_H2plus), σ_H2plus_bf_table .* 1e-18;
                                   extrapolation_bc=Linear())

# construct equilibrium constants.
# When Korg supports ionized molecules, we can use the Barklem and Collet values.
# For now, these are fine.
K_He2plus = LinearInterpolation(Ts_He2plus, K_He2plus_vals; extrapolation_bc=Linear())
K_H2plus  = LinearInterpolation(Ts_H2plus, K_H2plus_vals; extrapolation_bc=Linear())

end #module


text:
@graphs begin
	@node phase_2 = Start() -> p2_check
	@node p2_check = CheckStrategy("military") -> [n_true=p2_mili, n_false=p2_corr]

	# Military
	@node p2_mili = BinaryCondition("""
									The Shadow's victory points are less than the corruption points after the Free Peoples' have chosen whether to reveal.
									""") -> [n_true=p2_mili_change, n_false=p2_mili_end]
	@node p2_mili_change = SetStrategy("corruption") -> p2_mili_end

	@node p2_mili_end = End("End of Phase") -> []

	# Corruption
	@node p2_corr = BinaryCondition("""
									The corruption points are less than the Shadow's victory points after the Free Peoples' have chosen whether to reveal.
									""") -> [n_true=p2_corr_change, n_false=p2_corr_end]
	@node p2_corr_change = SetStrategy("military") -> p2_corr_end

	@node p2_corr_end = End("End of Phase") -> []

end

max_stars_repo_name: mvmorin/queller-bot

max_issues_repo_name: mvmorin/queller-bot

max_forks_repo_name: mvmorin/queller-bot


text:
using HeatEqCalcOld
sol = HeatEqCalcOld.calc_sol()
HeatEqCalcOld.save_sol(sol)

max_stars_repo_name: genkuroki/public

max_issues_repo_name: genkuroki/public

max_forks_repo_name: genkuroki/public


text:
zero_like(x) = x - x

type Beam
    """Complex beam parameter."""
    q

    """Wavelength."""
    λ
end

beam(waist_radius, offset, λ) = Beam(-offset + 1im * rayleigh_range(waist_radius, λ), λ)

rayleigh_range(b::Beam) = imag(b.q)
rayleigh_range(waist_radius, λ) = π * waist_radius^2 / λ

waist_offset(b::Beam) = -real(b.q)
waist_radius(b::Beam) = sqrt(rayleigh_range(b) * b.λ / π)

radius(b::Beam) = sqrt(-b.λ / (π * imag(1 / b.q)))
curvature(b::Beam) = 1 / real(1 / b.q)

abstract type Element end

struct Lens <: Element
    focal_length
end

struct FreeSpace <: Element
    distance
end

function expand_system(positions_and_elems, up_to_pos=nothing)
    last_pos = zero_like(up_to_pos == nothing ? positions_and_elems[1][1] : up_to_pos)
    system = Element[]
    for (p, e) in sort(positions_and_elems, by=a -> a[1])
        if up_to_pos != nothing && p > up_to_pos
            break
        end
        delta_pos = p - last_pos
        if delta_pos != zero_like(p)
            push!(system, FreeSpace(delta_pos))
        end
        last_pos = p

        push!(system, e)
    end
    if up_to_pos != nothing && last_pos < up_to_pos
        push!(system, FreeSpace(up_to_pos - last_pos))
    end
    system
end

function apply_abcd(beam::Beam, matrix)
    Beam((matrix[1, 1] * beam.q + matrix[1, 2]) /
        (matrix[2, 1] * beam.q + matrix[2, 2]), beam.λ)
end

function apply(beam::Beam, elem::Lens)
    apply_abcd(beam, [1 zero_like(elem.focal_length); (-1 / elem.focal_length) 1])
end

function apply(beam::Beam, elem::FreeSpace)
    apply_abcd(beam, [1 elem.distance; zero_like(1 / elem.distance) 1])
end

function propagate{E<:Element}(elements::Array{E, 1}, beam::Beam)
    reduce(apply, beam, elements)
end

export Beam, beam, Element, Lens, FreeSpace, expand_system, propagate


text:
using BayesLaplace
using Documenter

DocMeta.setdocmeta!(BayesLaplace, :DocTestSetup, :(using BayesLaplace); recursive=true)

makedocs(;
    modules=[BayesLaplace],
    authors="Patrick Altmeyer",
    repo="https://github.com/pat-alt/BayesLaplace.jl/blob/{commit}{path}#{line}",
    sitename="BayesLaplace.jl",
    format=Documenter.HTML(;
        prettyurls=get(ENV, "CI", "false") == "true",
        canonical="https://pat-alt.github.io/BayesLaplace.jl",
        assets=String[],
    ),
    pages=[
        "Home" => "index.md",
        "Tutorials" => [
            "Logistic Regression" => "tutorials/logit.md",
            "MLP" => "tutorials/mlp.md"
        ],
        "Reference" => "reference.md"
    ],
)

deploydocs(;
    repo="github.com/pat-alt/BayesLaplace.jl",
    devbranch="main",
)


text:
module MessagesTest

using ..Test

using HTTP.Messages
import HTTP.Messages.appendheader
import HTTP.URI
import HTTP.request
import HTTP: bytes

using HTTP: StatusError

using HTTP.MessageRequest: bodylength
using HTTP.MessageRequest: bodybytes
using HTTP.MessageRequest: unknown_length

using JSON

@testset "HTTP.Messages" begin

    @test bodylength(7) == unknown_length
    @test bodylength(UInt8[1,2,3]) == 3
    @test bodylength(view(UInt8[1,2,3], 1:2)) == 2
    @test bodylength("Hello") == 5
    @test bodylength(SubString("World!",1,5)) == 5
    @test bodylength(["Hello", " ", "World!"]) == 12
    @test bodylength(["Hello", " ", SubString("World!",1,5)]) == 11
    @test bodylength([SubString("Hello", 1,5), " ", SubString("World!",1,5)]) == 11
    @test bodylength([UInt8[1,2,3], UInt8[4,5,6]]) == 6
    @test bodylength([UInt8[1,2,3], view(UInt8[4,5,6],1:2)]) == 5
    @test bodylength([view(UInt8[1,2,3],1:2), view(UInt8[4,5,6],1:2)]) == 4
    @test bodylength(IOBuffer("foo")) == 3
    @test bodylength([IOBuffer("foo"), IOBuffer("bar")]) == 6

    @test bodybytes(7) == UInt8[]
    @test bodybytes(UInt8[1,2,3]) == UInt8[1,2,3]
    @test bodybytes(view(UInt8[1,2,3], 1:2)) == UInt8[1,2]
    @test bodybytes("Hello") == bytes("Hello")
    @test bodybytes(SubString("World!",1,5)) == bytes("World")
    @test bodybytes(["Hello", " ", "World!"]) == UInt8[]
    @test bodybytes([UInt8[1,2,3], UInt8[4,5,6]]) == UInt8[]

    req = Request("GET", "/foo", ["Foo" => "Bar"])
    res = Response(200, ["Content-Length" => "5"]; body="Hello", request=req)

    @test req.method == "GET"
    @test res.request.method == "GET"

    #display(req); println()
    #display(res); println()

    @test String(req) == "GET /foo HTTP/1.1\r\nFoo: Bar\r\n\r\n"
    @test String(res) == "HTTP/1.1 200 OK\r\nContent-Length: 5\r\n\r\nHello"

    @test header(req, "Foo") == "Bar"
    @test header(res, "Content-Length") == "5"
    setheader(req, "X" => "Y")
    @test header(req, "X") == "Y"

    ah(m::Message, h) = appendheader(m, SubString(h[1]) => SubString(h[2]))

    ah(req, "" => "Z")
    @test header(req, "X") == "YZ"

    ah(req, "X" => "more")
    @test header(req, "X") == "YZ, more"

    ah(req, "Set-Cookie" => "A")
    ah(req, "Set-Cookie" => "B")
    @test filter(x->first(x) == "Set-Cookie", req.headers) ==
        ["Set-Cookie" => "A", "Set-Cookie" => "B"]

    @test Messages.httpversion(req) == "HTTP/1.1"
    @test Messages.httpversion(res) == "HTTP/1.1"

    raw = String(req)
    #@show raw
    req = parse(Request,raw)
    #display(req); println()
    @test String(req) == raw

    req = parse(Request, raw * "xxx")
    @test String(req) == raw

    raw = String(res)
    #@show raw
    res = parse(Response,raw)
    #display(res); println()
    @test String(res) == raw

    res = parse(Response,raw * "xxx")
    @test String(res) == raw

    for sch in ["http", "https"]
        for m in ["GET", "HEAD", "OPTIONS"]
            @test request(m, "$sch://httpbin.org/ip").status == 200
        end
        try
            request("POST", "$sch://httpbin.org/ip")
            @test false
        catch e
            @test isa(e, StatusError)
            @test e.status == 405
        end
    end

#=
    @sync begin
        io = Base.BufferStream()
        @async begin
            for i = 1:100
                sleep(0.1)
                write(io, "Hello!")
            end
            close(io)
        end
        yield()
        r = request("POST", "http://httpbin.org/post", [], io)
        @test r.status == 200
    end
=#

    for sch in ["http", "https"]
        for m in ["POST", "PUT", "DELETE", "PATCH"]

            uri = "$sch://httpbin.org/$(lowercase(m))"
            r = request(m, uri)
            @test r.status == 200
            body = r.body

            io = Base.BufferStream()
            r = request(m, uri, response_stream=io)
            @test r.status == 200
            @test read(io) == body
        end
    end

    for sch in ["http", "https"]
        for m in ["POST", "PUT", "DELETE", "PATCH"]

            uri = "$sch://httpbin.org/$(lowercase(m))"
            io = Base.BufferStream()
            r = request(m, uri, response_stream=io)
            @test r.status == 200
        end

        r = request("POST", "$sch://httpbin.org/post",
                   ["Expect" => "100-continue"], "Hello")
        @test r.status == 200
        r = JSON.parse(String(r.body))
        @test r["data"] == "Hello"
    end

    mktempdir() do d
        cd(d) do

            n = 50
            io = open("result_file", "w")
            r = request("GET", "http://httpbin.org/stream/$n",
                        response_stream=io)
            @show filesize("result_file")
            i = 0
            for l in readlines("result_file")
                x = JSON.parse(l)
                @test i == x["id"]
                i += 1
            end
            @test i == n
        end
    end

end

end # module MessagesTest