text:
#
# Date created: 2022-03-06
# Author: aradclif
#
#
############################################################################################

#### Clean working draft
# function vargmax(A; dims=nothing)
#     dims === nothing && return _vargmax(A)
#     isone(ndims(A)) && return [_vargmax(A)]
#     @assert length(dims) == 1
#     axes_arg = axes(A)
#     axes_out = Base.setindex(axes_arg, length_one_axis(axes_arg[dims]), dims)
#     out = similar(A, CartesianIndex{ndims(A)}, axes_out)
#     Base.Cartesian.@nif 5 d -> (d ≤ ndims(A) && dims == d) d -> begin
#         Rpre = CartesianIndices(ntuple(i -> axes_arg[i], d - 1))
#         Rpost = CartesianIndices(ntuple(i -> axes_arg[i + d], ndims(A) - d))
#     end d -> begin
#         Rpre = CartesianIndices(axes_arg[1:dims-1])
#         Rpost = CartesianIndices(axes_arg[dims+1:end])
#     end
#     _vargmax_dims!(out, Rpre, 1:size(A, dims), Rpost, A)
# end

# function _vargmax_dims!(out, Rpre, is, Rpost, A)
#     s = typemin(eltype(A))
#     for Ipost ∈ Rpost, Ipre ∈ Rpre
#         indmax = 0
#         maxval = s
#         @turbo for i ∈ is
#             newmax = A[Ipre, i, Ipost] > maxval
#             maxval = newmax ? A[Ipre, i, Ipost] : maxval
#             indmax = newmax ? i : indmax
#         end
#         out[Ipre, 1, Ipost] = CartesianIndex(Ipre, indmax, Ipost)
#     end
#     return out
# end

# function _vargmax_dims_post!(out, is, Rpost, A)
#     s = typemin(eltype(A))
#     for j ∈ eachindex(Rpost) #Ipost ∈ Rpost
#         indm = 0
#         mval = s
#         Ipost = Rpost[j]
#         for i ∈ is
#             newm = >(A[i, Ipost], mval)
#             mval = newm ? A[i, Ipost] : mval
#             indm = newm ? i : indm
#         end
#         out[1, Ipost] = CartesianIndex(indm, Ipost)
#     end
#     out
# end

# function _vargmax_dims_ipost!(out, Rpre, is, Ipost, A)
#     s = typemin(eltype(A))
#     for Ipre ∈ Rpre
#         indm = 0
#         mval = s
#         @turbo for i ∈ is
#             newm = >(A[Ipre, i, Ipost], mval) #A[Ipre, i, Ipost] > maxval
#             mval = newm ? A[Ipre, i, Ipost] : mval
#             indm = newm ? i : indm
#         end
#         out[Ipre, 1, Ipost] = CartesianIndex(Ipre, indm, Ipost)
#     end
#     return out
# end

# function _vargmax_prepost!(out, Ipre, is, Ipost, A)
#     indm = 0
#     mval = typemin(eltype(A))
#     @turbo for i ∈ is
#         newm = >(A[Ipre, i, Ipost], mval) #A[Ipre, i, Ipost] > maxval
#         mval = newm ? A[Ipre, i, Ipost] : mval
#         indm = newm ? i : indm
#     end
#     out[Ipre, 1, Ipost] = CartesianIndex(Ipre, indm, Ipost)
# end

# function _vargmax(A::AbstractArray{T, N}) where {T, N}
#     indmax = 0
#     maxval = typemin(T)
#     @turbo for i ∈ eachindex(A)
#         newmax = A[i] > maxval
#         maxval = newmax ? A[i] : maxval
#         indmax = newmax ? i : indmax
#     end
#     CartesianIndices(A)[indmax]
# end

# # Not ideal, but makes return type stable
# function _vargmax(A::AbstractVector{T}) where {T}
#     indmax = 0
#     maxval = typemin(T)
#     @turbo for i ∈ eachindex(A)
#         newmax = >(A[i], maxval) #A[i] $op maxval
#         maxval = newmax ? A[i] : maxval
#         indmax = newmax ? i : indmax
#     end
#     indmax
# end

# # needs customization analogous to vmap, vmapreduce
# function vargmax(f, A::AbstractArray{T, N}) where {T, N}
#     indmax = 0
#     maxval = typemin(T)
#     @turbo for i ∈ eachindex(A)
#         φ = f(A[i])
#         newmax = >(ϕ, maxval)
#         maxval = newmax ? φ : maxval
#         indmax = newmax ? i : indmax
#     end
#     maxval, CartesianIndices(A)[indmax]
# end
# function vargmax(f, A::AbstractVector{T}) where {T}
#     indmax = 0
#     maxval = typemin(T)
#     @turbo for i ∈ eachindex(A)
#         φ = f(A[i])
#         newmax = >(ϕ, maxval)
#         maxval = newmax ? φ : maxval
#         indmax = newmax ? i : indmax
#     end
#     maxval, indmax
# end

# # # In fact, this is more than the Base findmax provides. It simply returns:
# # (A[argmax(A, dims)], argmax(A, dims))
# function vfindmax(A; dims=nothing)
#     dims === nothing && return _vfindmax(A)
#     isone(ndims(A)) && (begin vout, iout = _vfindmax(A) end; return [vout], [iout])
#     @assert length(dims) == 1
#     axes_arg = axes(A)
#     axes_out = Base.setindex(axes_arg, length_one_axis(axes_arg[dims]), dims)
#     Iout = similar(A, CartesianIndex{ndims(A)}, axes_out)
#     Vout = similar(A, axes_out)
#     Base.Cartesian.@nif 5 d -> (d ≤ ndims(A) && dims == d) d -> begin
#         Rpre = CartesianIndices(ntuple(i -> axes_arg[i], d - 1))
#         Rpost = CartesianIndices(ntuple(i -> axes_arg[i + d], ndims(A) - d))
#     end d -> begin
#         Rpre = CartesianIndices(axes_arg[1:dims-1])
#         Rpost = CartesianIndices(axes_arg[dims+1:end])
#     end
#     _vfindmax_dims!(Vout, Iout, Rpre, 1:size(A, dims), Rpost, A)
# end

# function _vfindmax_dims!(Vout, Iout, Rpre, is, Rpost, A)
#     s = typemin(eltype(A))
#     for Ipost ∈ Rpost, Ipre ∈ Rpre
#         indmax = 0
#         maxval = s
#         @turbo for i ∈ is
#             newmax = >(A[Ipre, i, Ipost], maxval)
#             maxval = newmax ? A[Ipre, i, Ipost] : maxval
#             indmax = newmax ? i : indmax
#         end
#         Vout[Ipre, 1, Ipost] = maxval
#         Iout[Ipre, 1, Ipost] = CartesianIndex(Ipre, indmax, Ipost)
#     end
#     return Vout, Iout
# end

# function _vfindmax(A::AbstractArray{T, N}) where {T, N}
#     indmax = 0
#     maxval = typemin(T)
#     @turbo for i ∈ eachindex(A)
#         newmax = >(A[i], maxval)
#         maxval = newmax ? A[i] : maxval
#         indmax = newmax ? i : indmax
#     end
#     maxval, CartesianIndices(A)[indmax]
# end
# function _vfindmax(A::AbstractVector{T}) where {T}
#     indmax = 0
#     maxval = typemin(T)
#     @turbo for i ∈ eachindex(A)
#         newmax = >(A[i], maxval)
#         maxval = newmax ? A[i] : maxval
#         indmax = newmax ? i : indmax
#     end
#     maxval, indmax
# end

# # Simpler, but not as efficient
# function vfindmax(A; dims=nothing)
#     dims === nothing && return _vfindmax(A)
#     @assert length(dims) == 1
#     Imax = vargmax(A, dims)
#     A[Imax], Imax
# end
# function _vfindmax(A::AbstractArray)
#     Imax = vargmax(A)
#     A[Imax], Imax
# end

#### vargmax, vargmin

for (name, f1, fdims) ∈ zip((:vargmax, :vargmin), (:_vargmax, :_vargmin), (:_vargmax_dims!, :_vargmin_dims!))
    @eval function $name(A; dims=nothing)
        dims === nothing && return $f1(A)
        isone(ndims(A)) && return [$f1(A)]
        @assert length(dims) == 1
        axes_arg = axes(A)
        axes_out = Base.setindex(axes_arg, length_one_axis(axes_arg[dims]), dims)
        out = similar(A, CartesianIndex{ndims(A)}, axes_out)
        Base.Cartesian.@nif 5 d -> (d ≤ ndims(A) && dims == d) d -> begin
            Rpre = CartesianIndices(ntuple(i -> axes_arg[i], d - 1))
            Rpost = CartesianIndices(ntuple(i -> axes_arg[i + d], ndims(A) - d))
        end d -> begin
            Rpre = CartesianIndices(axes_arg[1:dims-1])
            Rpost = CartesianIndices(axes_arg[dims+1:end])
        end
        $fdims(out, Rpre, 1:size(A, dims), Rpost, A)
    end
    @eval $name(A, ::Colon) = $name(A)
    # @eval $name(A; dims::Colon) = $name(A)
    # A neat idea, but poor performance due to much memory allocation
    # @eval function $name(A::AbstractArray, region::NTuple{N, Int}) where {N}
    #     if N == 1
    #         return $name(A, dims=region[1])
    #     elseif N == 2
    #         Iout1 = $name(A, dims=region[1])
    #         Iout2 = $name(A[Iout1], dims=region[2])
    #         return Iout1[Iout2]
    #     else
    #         I = $name(A, dims=region[1])
    #         for d = 2:N
    #             Ĩ = $name(A[I], dims=region[d])
    #             I = I[Ĩ]
    #         end
    #         return I
    #     end
    # end
end
for (op, name, init) ∈ zip((:>, :<), (:_vargmax_dims!, :_vargmin_dims!), (:typemin, :typemax))
    @eval function $name(out, Rpre, is, Rpost, A)
        s = $init(eltype(A))
        @inbounds for Ipost ∈ Rpost, Ipre ∈ Rpre
            indm = 0
            mval = s
            @turbo for i ∈ is
                newm = $op(A[Ipre, i, Ipost], mval)
                mval = newm ? A[Ipre, i, Ipost] : mval
                indm = newm ? i : indm
            end
            out[Ipre, 1, Ipost] = CartesianIndex(Ipre, indm, Ipost)
        end
        return out
    end
end
for (op, name, init) ∈ zip((:>, :<), (:_vargmax, :_vargmin), (:typemin, :typemax))
    @eval function $name(A::AbstractArray{T, N}) where {T, N}
        indm = 0
        mval = $init(T)
        @turbo for i ∈ eachindex(A)
            newm = $op(A[i], mval)
            mval = newm ? A[i] : mval
            indm = newm ? i : indm
        end
        # just the linear indm is nice, but one actually wants:
        CartesianIndices(A)[indm]
    end
    # Not ideal, but makes return type stable
    @eval function $name(A::AbstractVector{T}) where {T}
        indm = 0
        mval = $init(T)
        @turbo for i ∈ eachindex(A)
            newm = $op(A[i], mval) #A[i] $op maxval
            mval = newm ? A[i] : mval
            indm = newm ? i : indm
        end
        indm
    end
end

#### vtargmax, vtargmin : threaded version
for (name, f1, fdims) ∈ zip((:vtargmax, :vtargmin), (:_vtargmax, :_vtargmin), (:_vtargmax_dims!, :_vtargmin_dims!))
    @eval function $name(A; dims=nothing)
        dims === nothing && return $f1(A)
        isone(ndims(A)) && return [$f1(A)]
        @assert length(dims) == 1
        axes_arg = axes(A)
        axes_out = Base.setindex(axes_arg, length_one_axis(axes_arg[dims]), dims)
        out = similar(A, CartesianIndex{ndims(A)}, axes_out)
        Base.Cartesian.@nif 5 d -> (d ≤ ndims(A) && dims == d) d -> begin
            Rpre = CartesianIndices(ntuple(i -> axes_arg[i], d - 1))
            Rpost = CartesianIndices(ntuple(i -> axes_arg[i + d], ndims(A) - d))
        end d -> begin
            Rpre = CartesianIndices(axes_arg[1:dims-1])
            Rpost = CartesianIndices(axes_arg[dims+1:end])
        end
        $fdims(out, Rpre, 1:size(A, dims), Rpost, A)
    end
    @eval $name(A, ::Colon) = $name(A)
end
for (op, name, init) ∈ zip((:>, :<), (:_vargmax_prepost!, :_vargmin_prepost!), (:typemin, :typemax))
    @eval function $name(out, Ipre, is, Ipost, A)
        indm = 0
        mval = $init(eltype(A))
        @turbo for i ∈ is
            newm = $op(A[Ipre, i, Ipost], mval)
            mval = newm ? A[Ipre, i, Ipost] : mval
            indm = newm ? i : indm
        end
        out[Ipre, 1, Ipost] = CartesianIndex(Ipre, indm, Ipost)
        return out
    end
end
for (f, name) ∈ zip((:_vargmax_prepost!, :_vargmin_prepost!), (:_vtargmax_dims!, :_vtargmin_dims!))
    @eval function $name(out, Rpre, is, Rpost, A)
        if length(Rpre) > length(Rpost)
            @batch for Ipre ∈ Rpre, Ipost ∈ Rpost
                $f(out, Ipre, is, Ipost, A)
            end
        else
            @batch for Ipost ∈ Rpost, Ipre ∈ Rpre
                $f(out, Ipre, is, Ipost, A)
            end
        end
        return out
    end
end
for (op, name, init) ∈ zip((:>, :<), (:_vtargmax, :_vtargmin), (:typemin, :typemax))
    @eval function $name(A::AbstractArray{T, N}) where {T, N}
        indm = 0
        mval = $init(T)
        @tturbo for i ∈ eachindex(A)
            newm = $op(A[i], mval)
            mval = newm ? A[i] : mval
            indm = newm ? i : indm
        end
        CartesianIndices(A)[indm]
    end
    @eval function $name(A::AbstractVector{T}) where {T}
        indm = 0
        mval = $init(T)
        @tturbo for i ∈ eachindex(A)
            newm = $op(A[i], mval)
            mval = newm ? A[i] : mval
            indm = newm ? i : indm
        end
        indm
    end
end


text:
using DelaySSAdocs
using Test

@testset "DelaySSAdocs.jl" begin
    # Write your tests here.
end

ext: jl


text:
function RinexHeader(arg0::jdouble, arg1::SatelliteSystem, arg10::JString, arg11::JString, arg12::Vector3D, arg13::jdouble, arg14::Vector2D, arg15::Vector3D, arg16::JString, arg17::Vector3D, arg18::Vector3D, arg19::jdouble, arg2::JString, arg20::Vector3D, arg21::Vector3D, arg22::JString, arg23::jdouble, arg24::AbsoluteDate, arg25::AbsoluteDate, arg26::jint, arg27::List, arg28::List, arg29::List, arg3::JString, arg30::jint, arg31::jint, arg32::jint, arg33::jint, arg4::JString, arg5::JString, arg6::JString, arg7::JString, arg8::JString, arg9::JString)
    return RinexHeader((jdouble, SatelliteSystem, JString, JString, Vector3D, jdouble, Vector2D, Vector3D, JString, Vector3D, Vector3D, jdouble, JString, Vector3D, Vector3D, JString, jdouble, AbsoluteDate, AbsoluteDate, jint, List, List, List, JString, jint, jint, jint, jint, JString, JString, JString, JString, JString, JString), arg0, arg1, arg10, arg11, arg12, arg13, arg14, arg15, arg16, arg17, arg18, arg19, arg2, arg20, arg21, arg22, arg23, arg24, arg25, arg26, arg27, arg28, arg29, arg3, arg30, arg31, arg32, arg33, arg4, arg5, arg6, arg7, arg8, arg9)
end

function RinexHeader(arg0::jdouble, arg1::SatelliteSystem, arg10::JString, arg11::JString, arg12::Vector3D, arg13::jdouble, arg14::Vector2D, arg15::Vector3D, arg16::Vector3D, arg17::Vector3D, arg18::jdouble, arg19::AbsoluteDate, arg2::JString, arg20::AbsoluteDate, arg21::jint, arg22::jint, arg3::JString, arg4::JString, arg5::JString, arg6::JString, arg7::JString, arg8::JString, arg9::JString)
    return RinexHeader((jdouble, SatelliteSystem, JString, JString, Vector3D, jdouble, Vector2D, Vector3D, Vector3D, Vector3D, jdouble, AbsoluteDate, JString, AbsoluteDate, jint, jint, JString, JString, JString, JString, JString, JString, JString), arg0, arg1, arg10, arg11, arg12, arg13, arg14, arg15, arg16, arg17, arg18, arg19, arg2, arg20, arg21, arg22, arg3, arg4, arg5, arg6, arg7, arg8, arg9)
end

function get_agency_name(obj::RinexHeader)
    return jcall(obj, "getAgencyName", JString, ())
end

function get_antenna_azimuth(obj::RinexHeader)
    return jcall(obj, "getAntennaAzimuth", jdouble, ())
end

function get_antenna_b_sight(obj::RinexHeader)
    return jcall(obj, "getAntennaBSight", Vector3D, ())
end

function get_antenna_height(obj::RinexHeader)
    return jcall(obj, "getAntennaHeight", jdouble, ())
end

function get_antenna_number(obj::RinexHeader)
    return jcall(obj, "getAntennaNumber", JString, ())
end

function get_antenna_phase_center(obj::RinexHeader)
    return jcall(obj, "getAntennaPhaseCenter", Vector3D, ())
end

function get_antenna_reference_point(obj::RinexHeader)
    return jcall(obj, "getAntennaReferencePoint", Vector3D, ())
end

function get_antenna_type(obj::RinexHeader)
    return jcall(obj, "getAntennaType", JString, ())
end

function get_antenna_zero_direction(obj::RinexHeader)
    return jcall(obj, "getAntennaZeroDirection", Vector3D, ())
end

function get_approx_pos(obj::RinexHeader)
    return jcall(obj, "getApproxPos", Vector3D, ())
end

function get_center_mass(obj::RinexHeader)
    return jcall(obj, "getCenterMass", Vector3D, ())
end

function get_clk_offset(obj::RinexHeader)
    return jcall(obj, "getClkOffset", jint, ())
end

function get_eccentricities(obj::RinexHeader)
    return jcall(obj, "getEccentricities", Vector2D, ())
end

function get_interval(obj::RinexHeader)
    return jcall(obj, "getInterval", jdouble, ())
end

function get_leap_seconds(obj::RinexHeader)
    return jcall(obj, "getLeapSeconds", jint, ())
end

function get_leap_seconds_day_num(obj::RinexHeader)
    return jcall(obj, "getLeapSecondsDayNum", jint, ())
end

function get_leap_seconds_future(obj::RinexHeader)
    return jcall(obj, "getLeapSecondsFuture", jint, ())
end

function get_leap_seconds_week_num(obj::RinexHeader)
    return jcall(obj, "getLeapSecondsWeekNum", jint, ())
end

function get_list_applied_dcbs(obj::RinexHeader)
    return jcall(obj, "getListAppliedDCBS", List, ())
end

function get_list_applied_pcvs(obj::RinexHeader)
    return jcall(obj, "getListAppliedPCVS", List, ())
end

function get_marker_name(obj::RinexHeader)
    return jcall(obj, "getMarkerName", JString, ())
end

function get_marker_number(obj::RinexHeader)
    return jcall(obj, "getMarkerNumber", JString, ())
end

function get_marker_type(obj::RinexHeader)
    return jcall(obj, "getMarkerType", JString, ())
end

function get_observation_code(obj::RinexHeader)
    return jcall(obj, "getObservationCode", JString, ())
end

function get_observer_name(obj::RinexHeader)
    return jcall(obj, "getObserverName", JString, ())
end

function get_phase_shift_corrections(obj::RinexHeader)
    return jcall(obj, "getPhaseShiftCorrections", List, ())
end

function get_receiver_number(obj::RinexHeader)
    return jcall(obj, "getReceiverNumber", JString, ())
end

function get_receiver_type(obj::RinexHeader)
    return jcall(obj, "getReceiverType", JString, ())
end

function get_receiver_version(obj::RinexHeader)
    return jcall(obj, "getReceiverVersion", JString, ())
end

function get_rinex_version(obj::RinexHeader)
    return jcall(obj, "getRinexVersion", jdouble, ())
end

function get_satellite_system(obj::RinexHeader)
    return jcall(obj, "getSatelliteSystem", SatelliteSystem, ())
end

function get_signal_strength_unit(obj::RinexHeader)
    return jcall(obj, "getSignalStrengthUnit", JString, ())
end

function get_t_first_obs(obj::RinexHeader)
    return jcall(obj, "getTFirstObs", AbsoluteDate, ())
end

function get_t_last_obs(obj::RinexHeader)
    return jcall(obj, "getTLastObs", AbsoluteDate, ())
end


text:
module UnrollingAverages

export unroll

using Combinatorics
using LinearAlgebra
using Statistics

include("utilities.jl")
include("NaturalSeriesUnrolling.jl")

end


text:
# Autogenerated wrapper script for Xorg_xorgproto_jll for any
## Global variables
PATH = ""
LIBPATH = ""
LIBPATH_env = "LD_LIBRARY_PATH"
LIBPATH_default = ""

"""
Open all libraries
"""
function __init__()
    global artifact_dir = abspath(artifact"Xorg_xorgproto")

    # Initialize PATH and LIBPATH environment variable listings
    global PATH_list, LIBPATH_list
    # Filter out duplicate and empty entries in our PATH and LIBPATH entries
    filter!(!isempty, unique!(PATH_list))
    filter!(!isempty, unique!(LIBPATH_list))
    global PATH = join(PATH_list, ':')
    global LIBPATH = join(vcat(LIBPATH_list, [joinpath(Sys.BINDIR, Base.LIBDIR, "julia"), joinpath(Sys.BINDIR, Base.LIBDIR)]), ':')

end  # __init__()


text:
module Mixed

using ExtractMacro
using ..Common
using ..Interface

export GraphMixed

import ..Interface: energy, delta_energy, neighbors, update_cache!

struct GraphMixed{ET} <: SimpleGraph{ET}
    N::Int
    graphs::Vector{AbstractGraph}
end

"""
    GraphMixed(graphs::AbstractGraph...)

A `SimpleGraph` mixing 2 or more base graphs into a single model.
The energy function of a mixed graph is just the sum of the energy functions
of the individual graphs in `graphs`.

***Usage Example***:

    X = RRRMC.GraphMixed(RRRMC.GraphRRG(10, 3),
                         RRRMC.GraphSK(10))

    X = RRRMC.GraphMixed(RRRMC.GraphSK(31),
                         RRRMC.GraphPercStep(31, 30),
                         RRRMC.GraphSAT(31, 3, 4.2))
"""
function GraphMixed(graphs::AbstractGraph...)
    length(graphs) >= 2 || throw(ArgumentError("at least 2 base graphs needed"))
    N = getN(graphs[1])
    all(getN(x) == N for x in graphs) || throw(ArgumentError("same N for all graphs required"))
    ET = typeof(sum(zero(energy_type(x)) for x in graphs))
    gs = AbstractGraph[x for x in graphs]
    return GraphMixed{ET}(N, gs)
end

function energy(X::GraphMixed{ET}, C::Config) where ET
    sum(energy(x, C) for x in X.graphs)::ET
end

function delta_energy(X::GraphMixed{ET}, C::Config, move::Int) where ET
    sum(delta_energy(x, C, move) for x in X.graphs)::ET
end

function neighbors(X::GraphMixed, i::Int)
    union(neighbors.(X.graphs, i)...)
end

function update_cache!(X::GraphMixed, C::Config, move::Int)
    foreach(x->update_cache!(x, C, move), X.graphs)
end

end #module


text:
##
## Fermi-Ulam model
##
## Tomáš Kalvoda, FIT CTU, 2015
##

##
## Related articles
##
## E. Fermi, On the Origin of the Cosmic Radiation, Phys. Rev. 75 (8), 1949, 1169-1174.
##
## S. Ulam, On the Statistical Properties of Dynamical Systems, Proc. 4th Berkeley Symp. Math. Stat. Probab. 3, 1961, 315-320.
##

using FixedPointNumbers, Colors, Images

using ProgressMeter

# parameters of the model

const period = 1.0
const L = 1.3

wall_prime(t) = 0.2 * ( abs((t % period) - period/2) - period/4 )

# computation of trajectory

function fermi(t, v, n)
    output = Array(Float64, (n,2))
    for j in 1:n
        t = (t + (2*L/v)) % period
        v = abs(-v + 2*wall_prime(t))
        output[j,:] = [t v]
    end
    return output
end

# draw the trajectory

function draw_trajectory(data, image, tox, toy, color)
    height, width = size(image)
    for j in 1:size(data)[1]
        x,y = (tox(data[j,1]), toy(data[j,2]))
        if 1 <= x <= width && 1 <= y <= height
            image[y,x,:] = color
        end
    end
end

# draw the phase space

function phase_space_image(;width=800, height=1600, tmin=0., tmax=period, vmin=1., vmax=100., deltav=1.1, nmax=10000)
    img = Array(Float64, (height,width,3))
    fill!(img, 1.)
    tox(t) = convert(Int64, floor(width*(t-tmin)/(tmax-tmin)))
    toy(v) = convert(Int64, floor(height*(v-vmin)/(vmax-vmin)))
    number_of_trajectories = convert(Int64, ceil((vmax-vmin)/deltav))
    getcolor(j) = convert(RGB, HSV(360.*j/number_of_trajectories, 1, 1))
    v = vmin + deltav
    j = 1

    info("Computing trajectories...")
    p = Progress(number_of_trajectories, 1)

    while v <= vmax
        trajectory = fermi(tmin+rand()*(tmax-tmin), v, nmax)
        if isodd(j)
            clr = getcolor(j)
        else
            clr = getcolor(number_of_trajectories-j)
        end
        draw_trajectory(trajectory, img, tox, toy, [clr.r, clr.g, clr.b])
        v += deltav
        j += 1
        next!(p)
    end

    return flipud(img)
end

max_stars_repo_name: kalvotom/julia-examples

max_issues_repo_name: kalvotom/julia-examples

max_forks_repo_name: kalvotom/julia-examples


text:
# This file was formerly a part of Julia. License is MIT: https://julialang.org/license

include("UnixDomainCM.jl")

addprocs(UnixDomainCM(parse(Int,ARGS[1])))
resp = pmap(x -> myid() *2, [1:nworkers()])
println(resp)

max_stars_repo_name: vtjnash/Examples

max_issues_repo_name: JuliaArchive/Examples

max_forks_repo_name: JuliaArchive/Examples


text:
function _precompile_()
    ccall(:jl_generating_output, Cint, ()) == 1 || return nothing
    precompile(Makie.backend_display, (GLBackend, Scene))

    activate!()
    precompile(refreshwindowcb, (GLFW.Window, Screen))
    p = plot(rand(10))
    # Mimic `display(p)` without actually creating a display
    function insertplotstype(scene)
        for elem in scene.plots
            inserttype(scene, elem)
        end
        foreach(s-> insertplotstype(s), scene.children)
    end
    function inserttype(scene, @nospecialize(x))
        if isa(x, Combined)
            if isempty(x.plots)
                precompile(insert!, (Screen, typeof(scene), typeof(x)))
            else
                foreach(x.plots) do x
                    inserttype(scene, x)
                end
            end
        else
            precompile(insert!, (Screen, typeof(scene), typeof(x)))
        end
    end
    scene = p.figure.scene
    insertplotstype(scene)
end


text:
module DABenchmarks

using Distributed
using LShapedSolvers
using HydroModels
using BenchmarkTools
using RecipesBase
using LaTeXStrings
using Colors
using Statistics
using Colors
using Printf
using Logging

export
    DABenchmark,
    save_benchmark,
    load_benchmark,
    warmup_benchmarks!,
    run_benchmarks!,
    save_results!,
    load_results!,
    load_results,
    Sequential,
    Scaling

mutable struct DABenchmark{BType}
    names::Vector{String}
    benchmarks::BenchmarkGroup
    medians::Matrix{Float64}
    results::BenchmarkGroup

    function (::Type{DABenchmark})(names::Vector{String}, benchmarks::BenchmarkGroup, btype::Symbol)
        return new{btype}(names, benchmarks, zeros(length(names),1), BenchmarkGroup())
    end

    function (::Type{DABenchmark})(names::Vector{String}, benchmarks::BenchmarkGroup, results::BenchmarkGroup, btype::Symbol)
        return new{btype}(names, benchmarks, zeros(length(names),1), results)
    end
end

function Sequential(solvernames::Vector{String}, solvers::Vector, nscenarios::Integer, nsamples::Integer; timeout::Int = 100)
    length(solvernames) == length(solvers) || error("Inconsistent number of solvernames/solvers")
    benchmarks = BenchmarkGroup(solvernames)
    for (solvername, solver) in zip(solvernames, solvers)
        solve_time = @elapsed begin
            da = load_model(nscenarios)
            plan!(da, optimsolver = solver)
        end
        max_time = min((nsamples+1)*solve_time, timeout)
        benchmarks[solvername] = @benchmarkable plan!(da,optimsolver=$solver) seconds=max_time samples=nsamples setup=(da = load_model($nscenarios))
    end
    dab = DABenchmark(solvernames, benchmarks, :seq)
    return dab
end

function Scaling(solvernames::Vector{String}, solvers::Vector, nscenarios::Integer, nsamples::Integer; timeout::Int = 1000)
    length(solvernames) == length(solvers) || error("Inconsistent number of solvernames/solvers")
    benchmarks = BenchmarkGroup()
    b = benchmarks[string(nworkers())] = BenchmarkGroup(solvernames)
    for (solvername, solver) in zip(solvernames, solvers)
        solve_time = @elapsed begin
            da = load_model(nscenarios)
            plan!(da, optimsolver = solver)
        end
        max_time = min((nsamples+1)*solve_time, timeout)
        b[solvername] = @benchmarkable plan!(da,optimsolver=$solver) seconds=max_time samples=nsamples setup=(da = load_model($nscenarios))
    end
    dab = DABenchmark(solvernames, benchmarks, :scaling)
    return dab
end

function load_model(nscenarios::Integer)
    dayahead_data = HydroModels.NordPoolDayAheadData("data/plantdata.csv", "data/spotprices.csv", 1, 35.0, 60.0)
    sampler = DayAheadSampler(dayahead_data)
    dayahead_model = DayAheadModel(dayahead_data, sampler, nscenarios, [:Skellefteälven])
    return dayahead_model
end

function collect_estimates!(da_benchmark::DABenchmark{:seq})
    medians = time(median(da_benchmark.results))
    da_benchmark.medians = zeros(length(da_benchmark.names))
    for (j,solvername) in enumerate(da_benchmark.names)
        da_benchmark.medians[j] = medians[solvername] / 1e9
    end
end

function collect_estimates!(da_benchmark::DABenchmark{:scaling})
    medians = time(median(da_benchmark.results))
    nworkers = sort(parse.(Int, keys(da_benchmark.results)))
    da_benchmark.medians = zeros(length(keys(da_benchmark.results)), length(da_benchmark.names))
    for (i,w) in enumerate(nworkers)
        for (j,solvername) in enumerate(da_benchmark.names)
            da_benchmark.medians[i,j] = medians[string(w)][solvername] / 1e9
        end
    end
end

function collect_estimates!(da_benchmark::DABenchmark{:async})
    medians = time(median(da_benchmark.results))
    κs = sort(parse.(Float64, keys(da_benchmark.results[first(da_benchmark.names)])))
    da_benchmark.medians = zeros(length(da_benchmark.names), length(κs))
    for (i,solvername) in enumerate(da_benchmark.names)
        for (j,κ) in enumerate(κs)
            da_benchmark.medians[i,j] = medians[solvername][string(κ)] / 1e9
        end
    end
end

function collect_estimates!(da_benchmark::DABenchmark{:bundle})
    medians = time(median(da_benchmark.results))
    bundles = sort(parse.(Int, keys(da_benchmark.results[first(da_benchmark.names)])))
    da_benchmark.medians = zeros(length(da_benchmark.names), length(bundles))
    for (i,solvername) in enumerate(da_benchmark.names)
        for (j,b) in enumerate(bundles)
            da_benchmark.medians[i,j] = medians[solvername][string(b)] / 1e9
        end
    end
end

function save_benchmark(da_benchmark::DABenchmark, filename::String)
    BenchmarkTools.save(filename, da_benchmark.benchmarks)
end

function load_benchmark(filename::String)
    benchmarks = BenchmarkTools.load(filename)
    solvernames::Vector{String} = collect(keys(first(benchmarks)))
    da_benchmark = DABenchmark(solvernames,benchmarks)
    return da_benchmark
end

function warmup_benchmarks!(da_benchmark::DABenchmark)
    warmup(da_benchmark.benchmarks)
    nothing
end

function run_benchmarks!(da_benchmark::DABenchmark)
    da_benchmark.results = run(da_benchmark.benchmarks, verbose=true)
    collect_estimates!(da_benchmark)
    return da_benchmark.results
end

function save_results!(da_benchmark::DABenchmark, filename::String)
    BenchmarkTools.save(filename, da_benchmark.results)
    nothing
end

function load_results!(da_benchmark::DABenchmark, filename::String)
    da_benchmark.results = BenchmarkTools.load(filename)[1]
    collect_estimates!(da_benchmark)
    return da_benchmark.results
end

function load_results(filename::String, ::Val{:async}) where btype
    results = BenchmarkTools.load(filename)[1]
    solvernames::Vector{String} = collect(keys(results))
    da_benchmark = DABenchmark(solvernames,BenchmarkGroup(), :async)
    load_results!(da_benchmark, filename)
    return da_benchmark
end

function load_results(filename::String, ::Val{:bundle}) where btype
    results = BenchmarkTools.load(filename)[1]
    solvernames::Vector{String} = collect(keys(results))
    da_benchmark = DABenchmark(solvernames,BenchmarkGroup(), :bundle)
    load_results!(da_benchmark, filename)
    return da_benchmark
end

function load_results(filename::String, ::Val{btype}) where btype
    results = BenchmarkTools.load(filename)[1]
    solvernames::Vector{String} = collect(keys(first(values(results))))
    da_benchmark = DABenchmark(solvernames,BenchmarkGroup(),btype)
    load_results!(da_benchmark, filename)
    return da_benchmark
end

KTH_colors = [RGB(25/255,84/255,166/255),
              RGB(157/255,16/255,45/255),
              RGB(98/255,146/255,46/255),
              RGB(36/255,160/255,216/255),
              RGB(228/255,54/255,62/255),
              RGB(176/255,201/255,43/255),
              RGB(216/255,84/255,151/255),
              RGB(250/255,185/255,25/255),
              RGB(101/255,101/255,108/255),
              RGB(189/255,188/255,188/255),
              RGB(227/255,229/255,227/255)]

@recipe function f(benchmark::DABenchmark{:scaling})
    nworkers = sort(parse.(Int, keys(benchmark.results)))
    increment = std(benchmark.medians)
    tmin = 0.9*minimum(benchmark.medians)
    tmax = 1.05*maximum(benchmark.medians)

    linewidth --> 2
    tickfontsize := 10
    tickfontfamily := "sans-serif"
    guidefontsize := 10
    guidefontfamily := "sans-serif"
    legend := :topright
    xticks := sort(nworkers)
    xlabel := "Number of Cores P"
    ylabel := "Computation Time T [s]"
    ylims --> (tmin,tmax)
    yticks --> tmin:increment:tmax
    yformatter := (d) -> @sprintf("%.1f",d)
    color_palette := KTH_colors

    for (i,solver) in enumerate(benchmark.names)
        @series begin
            label --> solver
            nworkers, benchmark.medians[:,i]
        end
    end
end

@recipe function f(benchmark::DABenchmark{:async})
    κs = sort(parse.(Float64, keys(benchmark.results[first(benchmark.names)])))
    increment = std(benchmark.medians)
    tmin = 0.9*minimum(benchmark.medians)
    tmax = 1.05*maximum(benchmark.medians)

    linewidth --> 2
    tickfontsize := 10
    tickfontfamily := "sans-serif"
    guidefontsize := 10
    guidefontfamily := "sans-serif"
    legend := :topright
    xticks := sort(κs)
    xlabel := L"\kappa"
    ylabel := "Computation Time T [s]"
    ylims --> (tmin,tmax)
    yticks --> tmin:increment:tmax
    yformatter := (d) -> @sprintf("%.1f",d)
    color_palette := KTH_colors

    for (i,solver) in enumerate(benchmark.names)
        @series begin
            label --> solver
            κs, benchmark.medians[i,:]
        end
    end
end

@recipe function f(benchmark::DABenchmark{:bundle})
    bundles = sort(parse.(Int, keys(benchmark.results[first(benchmark.names)])))
    increment = std(benchmark.medians)
    tmin = 0.9*minimum(benchmark.medians)
    tmax = 1.05*maximum(benchmark.medians)

    linewidth --> 2
    tickfontsize := 10
    tickfontfamily := "sans-serif"
    guidefontsize := 10
    guidefontfamily := "sans-serif"
    legend := :topright
    xticks := sort(bundles)
    xlabel := "# Cuts in Bundle B"
    ylabel := "Computation Time T [s]"
    ylims --> (tmin,tmax)
    yticks --> tmin:increment:tmax
    yformatter := (d) -> @sprintf("%.1f",d)
    color_palette := KTH_colors

    for (i,solver) in enumerate(benchmark.names)
        @series begin
            label --> solver
            bundles, benchmark.medians[i,:]
        end
    end
end

end#module

max_stars_repo_name: martinbiel/mpc-reproducibility

max_issues_repo_name: martinbiel/mpc-reproducibility

max_forks_repo_name: martinbiel/mpc-reproducibility


text:
"""
Global choice of form of sylvester matrix.
"""
global const sylves = sylves2
global const extract_sylves = extract_sylves2

"""
  gcdroot(z::Vector{Number} [, tol=1e-10]) -> PolyZeros, bkerror

  Calculate the multiplicity structure and initial root
  approximation of a real or complex polynomial p,

        p(x) = a_1 x^n + a_2 x^n-1 + ... + a_n x + a_n+1,

  given by an (n+1)-row vector p = (a_1, ..., a_n+1)

  METHOD : details can be found in a paper of the author  Zhonggang Zeng

  INPUT :  p = polynomial coefficients.
           tol = zero remainder threshold;
                 default = 10^{-10}, set internally
                 if no value is specified.

  OUTPUT : z = roots of polynomial p
           l = corresponding multiplicities
           bkerr = backward error

  CALL :   It is best to call pzero with only one argument as
               `z = gcdroot(p)`.

"""
function gcdroot(p::AbstractVector{T}, tol::S = 1e-10; gamma=S(100), delta = S(100), thresh=tol*100, drop=S(5e-5)) where {T<:Number,S<:AbstractFloat}
    # gamma     residual growth factor, default = 100
    # delta     threshold growth factor, default = 100
    # thresh    zero singular threshold, default = 100*tol
    # drop      try Gauß-Newton if sigma(k) < drop * sigma(k-1)

    E = one(T); RE = real(E); Z = zero(T); RZ = real(Z)

    a = findfirst(!iszero, p)
    a !== nothing || throw(ArgumentError("gcdroot of zero polynomial not allowed"))
    # make the polynomial monic
    k = length(p)
    n = k - a
    p = f = view(p, a:k) / p[a] # working copy of the monic polynomial
    g = deriv(f) / n            #  g(x) = f'(x) / degree(f)
    nf = norm(f, Inf)           # the largest coefficient

    mx = n; wtol = tol; s0 = RZ; s = RZ; wthrh = thresh

    z = complex(T)[]        # all different roots
    l = Int64[]             # the multiplicities of those roots as known so far
    nroots = 0              # invariant: length(z) == length(l) == nroots

    k = n;                  # the degree of working polynomial
    while k >= 1
        m = 1 # ensure this variable is used as loop index and value preserved after break
        if k == 1     # the polynomial is linear, GCD = 1
            h = E; u = f; v = E
        else
            for mm = 1:k
                m = mm
                A = sylves(f, g, m)
                s0 = s  # save previous sigma

                s, x = zminsv(A, tol)
                if m > 1 && s > s0 + norm(f, Inf) * eps()
                    m -= 1
                    break
                end
                # @printf("s. value %g,%g,%g\n", m, k, s)
                if x[1] == 0
                    x[1] = eps()
                end
                # I don't understand why s is compared to nf.
                if s < wthrh * nf || m == mx || s < drop * s0
                    # h = gcd(f, g); u = f / h; v = g / h
                    h, u, v, res0, res, sm = gcd_refinement(x, m, f, g, A)
                    if m == length(u) - 1 && (res < wtol || m == mx)
                        wtol = max(wtol, res * gamma)	     # increase tolerance by factor
                        wthrh = max(wthrh, (s / nf) * delta) # increase threshold by factor
                        break # for mm
                    end
                end
            end # for mm
        end
        # println("after m-loop: m = $m mx=$mx k=$k n=$n s=$s thr=$(wthrh*nf) drop=$(drop*s0) ")
        if k == n			# the root values of u contain all roots of f
            z = roots(u)	# u has only simple roots
            # assert length(z) == m
            length(z) == m || error("number of roots of u $(length(z)) is not $m")
            nroots = m
            l = ones(Int64, nroots)
            if m == 1
                # only one root value with multiplicity n
                l[1] = k
                k = 0
            end
            # verify roots of u are roots of f: horner_analysis(z, f)
        else
            t = roots(u)	# u has only simple roots - they should be in z
            jj = 0
            ix = Vector(1:nroots)
            for j = 1:m
                tj = t[j]
                tz, jm = findmin(abs.(z[ix] .- tj))	# find root closest to tj
                jj = ix[jm]
                deleteat!(ix, jm)
                l[jj] += 1                     # increase multiplicity of that root
                # println("k = $k m = $m: z[$jj] = $(z[jj]) tz = $tz")
                tz <= 0.1 * (abs(z[jj])+0.01) || error("inacceptable root, k = $k tz = $tz")
            end
            if m == 1
                l[jj] += k - 1
                k = 0
            end
        end
        k = max(k - m, 0)
        if k <= 0
            return finish(z, l, p)
        end
        f = h               # f = gcd(f, f') -- all positive multiplicities decreased by 1
        g = deriv(f) / k
        nf = norm(f, Inf)
        mx = m
    end # k-while loop
end

function finish(z, l, p)
    f = polymult(z, l)
    w = [ max(abs(pj), 1) for pj in p]
    # println("f = $(size(f)) p = $(size(p)) w = $(size(w))")
    bkerr = maximum(abs, (f - p) ./ w)

    # println("finish: bkerr = $bkerr")
    # display([f p (f-p) ./ w])

    PolyZeros(z, l), bkerr
end

function gcd_refinement(x, m, f, g, A)
    u0, v0 = extract_sylves(x)
    # println("after extract_sylves: x = $x\nu0 = $u0\nv0 = $v0")
    # println("refinement u0: $(norm(convolute(u0,g) - convolute(v0,f))/norm(f))")
    B = cauchymt(u0, length(f) - length(u0) + 1)
    g0 = scalsq(B, f)	# scaled least squares solver g0 = div(f, u0)
    g0 = g0 / g0[1]		# first approximation of gcd(p, p')
    # printresidue("g0", g0, u0, v0, f, g)
    h, u, v, res0, res = gcdgn(f, g, g0, u0, v0)	# refinement of g, u, v by G-N
    # println("refinement u:  $(norm(convolute(u,g) - convolute(v,f))/norm(f))")
    # printresidue("h ", h, u, v, f, g)
    sm = norm(A * [u; v]) / norm([u;v])
    h, u, v, res0, res, sm
end

function printresidue(text, g0, u0, v0, f, g)
    uuu = convolute(g0, u0) - f	# should be zero
    vvv = convolute(g0, v0) - g	# should be zero

    # println("$text p(x)  $f")
    # println("$text p'(x)  $g")
    println("$text u(x)  $u0")
    println("$text v(x)  $v0")
    # println("$text h(x)  $g0")
    # println("$text Δp(x)  $uuu")
    # println("$text Δp'(x) $vvv")
end

function analysis_sv(k, m, A, x, s)

    #println("svcheck: $k-$m norm(Ax) = $(norm(A*x)) s = $s")

end

function scalerows!(A::Array{<:Number,2})
    s = map(i->2.0^-exponent(norm(A[i,:])), 1:size(A,1))
    scale!(s, A)
    s
end

scale!(s, A) = lmul!(Diagonal(s), A)


text:
# Test the default case
rank_data = rankinequalities(m, active_constr = true)
@test isa(rank_data, Vector{Dict})

# Test varying the max_ranks to work as expected
rank_data = rankinequalities(m, max_ranks = 1, active_constr = true)
@test isa(rank_data, Vector{Dict}) && length(rank_data) == 1
@test_throws ErrorException rankinequalities(m, max_ranks = 0, active_constr = true)


text:
# SPDX-License-Identifier: MIT

"""
    struct CompiledFSM{T<:Semifield}
        π        # vector of initial probabilities
        ω        # vector of final probabilities
        A        # matrix of transition probabilities
        Aᵀ       # transpose of `A`
        pdfmap   # mapping state -> pdfindex
    end

Compiled FSM: matrix/vector format of an FSM used by inference
algorithms. Note that in this form, every state is an emitting state.
"""
struct CompiledFSM{T<:Semifield}
    π::AbstractVector{T}
    ω::AbstractVector{T}
    A::AbstractMatrix{T}
    Aᵀ::AbstractMatrix{T}
    pdfmap::AbstractVector{PdfIndex}
end

fsmtype(::CompiledFSM{T}) where T<:Semifield = T

"""
    compile(fsm; allocator = spzeros)

Compile `fsm` into a inference-friendly format [`CompiledFSM`](@ref).
`allocator` is a function analogous to `zeros` which creates an
N-dimensional array and fills it with zero elements.
"""
function compile(fsm::FSM{T}; allocator = spzeros) where T
    allstates = collect(states(fsm))
    S = length(allstates)

    # Initial states' probabilities.
    π = allocator(T, S)
    for s in filter(isinit, allstates) π[s.id] = s.initweight end

    # Final states' probabilities.
    ω = allocator(T, S)
    for s in filter(isfinal, allstates) ω[s.id] = s.finalweight end

    # Transition matrix.
    A = allocator(T, S, S)
    Aᵀ = allocator(T, S, S)
    for src in allstates
        for link in links(fsm, src)
            A[src.id, link.dest.id] = link.weight
            Aᵀ[link.dest.id, src.id] = link.weight
        end
    end

    pdfmap = [s.pdfindex for s in sort(allstates, by = p -> p.id)]

    CompiledFSM{T}(π, ω, A, Aᵀ, pdfmap)
end

"""
    gpu(cfsm)

Move the compiled fsm `cfsm` to GPU.
"""
function gpu(cfsm::CompiledFSM{T}) where T
    if ! issparse(cfsm.π)
        return CompiledFSM{T}(
            CuArray(cfsm.π),
            CuArray(cfsm.ω),
            CuArray(cfsm.A),
            CuArray(cfsm.Aᵀ),
            CuArray(cfsm.pdfmap)
        )
    end

    A = CuSparseMatrixCSC(cfsm.A)
    Aᵀ = CuSparseMatrixCSC(cfsm.Aᵀ)
    return CompiledFSM{T}(
        CuSparseVector(cfsm.π),
        CuSparseVector(cfsm.ω),
        CuSparseMatrixCSR(Aᵀ.colPtr, Aᵀ.rowVal, Aᵀ.nzVal, A.dims),
        CuSparseMatrixCSR(A.colPtr, A.rowVal, A.nzVal, A.dims),
        CuArray(cfsm.pdfmap)
    )
end

function Base.convert(::Type{CompiledFSM{NT}}, cfsm) where NT <: Semifield
    CompiledFSM(
        copyto!(similar(cfsm.π, NT), cfsm.π),
        copyto!(similar(cfsm.ω, NT), cfsm.ω),
        copyto!(similar(cfsm.A, NT), cfsm.A),
        copyto!(similar(cfsm.Aᵀ, NT), cfsm.Aᵀ),
        copyto!(similar(cfsm.pdfmap, Int), cfsm.pdfmap)
    )
end
Base.convert(::Type{T}, cfsm::T) where T <: CompiledFSM = cfsm

max_stars_repo_name: MartinKocour/MarkovModels

max_issues_repo_name: MartinKocour/MarkovModels

max_forks_repo_name: MartinKocour/MarkovModels


text:
using BallroomSkatingSystem
using Documenter

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

makedocs(;
    modules = [BallroomSkatingSystem],
    authors = "Sören Dobberschütz and contributors",
    repo = "https://github.com/sdobber/BallroomSkatingSystem.jl/blob/{commit}{path}#{line}",
    sitename = "BallroomSkatingSystem.jl",
    format = Documenter.HTML(;
        prettyurls = get(ENV, "CI", "false") == "true",
        canonical = "https://sdobber.github.io/BallroomSkatingSystem.jl",
        assets = String[]
    ),
    pages = [
        "Home" => "index.md",
        "API Reference" => "api.md",
        "Examples" => "examples.md",
    ]
)

deploydocs(;
    repo = "github.com/sdobber/BallroomSkatingSystem.jl"
)


text:
include("jo_invertible_matrix.jl")
include("jo_lagrange_interp_1d.jl")


text:
function solution()
    values = readlines("smalltest.txt")
    values = parse.(Int, values)
    windowindices = [1, 2, 3]
    prevsum = -Inf
    count = 0
    while windowindices[3] <= length(values)
        currsum = values[windowindices[1]] + values[windowindices[2]] + values[windowindices[3]]
        if prevsum > 0 && currsum > prevsum
            count += 1
        end
        windowindices .+= 1
        prevsum = currsum
    end
    println(count)
end

solution()


text:
@testset "Test `Indices` constructors" begin
    @test Miller([-1, 0, 1]) == Miller(-1, 0, 1)
    @test ReciprocalMiller(-1, 0, 1) == ReciprocalMiller([-1, 0, 1])
    @test MillerBravais([-2, 1, 1, 3]) == MillerBravais([-2, 1, 1, 3])
    @test ReciprocalMillerBravais([-2, 1, 1, 3]) == ReciprocalMillerBravais([-2, 1, 1, 3])
end

@testset "Test `@m_str` & `@mb_str`" begin
    @test m"[-1, 0, 1]" == Miller([-1, 0, 1])
    @test m"<2, -1, -1, 3>" == MillerBravais([2, -1, -1, 3])
    @test m"[  2,-1, -1,3  ]" == MillerBravais([2, -1, -1, 3])
    @test m"(-1, 0, 1  )" == ReciprocalMiller([-1, 0, 1])
    @test m"{  -1,0, 1}" == ReciprocalMiller([-1, 0, 1])
    @test m"(1, 0, -1, 0)" == ReciprocalMillerBravais(1, 0, -1, 0)
    @test m"{  1,0, -1, 0}" == ReciprocalMillerBravais(1, 0, -1, 0)
end

@testset "Test conversion between real `Miller` and `MillerBravais`" begin
    millerreal = [[1, 1, 1], [-1, 0, 1], [0, -1, 1], [-1, -1, 1], [1, 0, 1], [0, 1, 1]]
    millerbravaisreal = [
        [1, 1, -2, 3],
        [-2, 1, 1, 3],
        [1, -2, 1, 3],
        [-1, -1, 2, 3],
        [2, -1, -1, 3],
        [-1, 2, -1, 3],
    ]
    for (x, y) in zip(millerreal, millerbravaisreal)
        m = Miller(x)
        mb = MillerBravais(y)
        @test convert(MillerBravais, m) == mb
        @test convert(Miller, mb) == m
        @test convert(Miller, m) == m
        @test convert(MillerBravais, mb) == mb
        @test_throws MethodError convert(ReciprocalMiller, m)
        @test_throws MethodError convert(ReciprocalMillerBravais, mb)
        @test_throws MethodError convert(ReciprocalMillerBravais, m)
        @test_throws MethodError convert(ReciprocalMiller, mb)
    end
end

@testset "Test conversion between reciprocal `Miller` and `MillerBravais`" begin
    millerreciprocal =
        [[1, 0, 0], [0, 1, 0], [1, -1, 0], [-1, 0, 0], [0, -1, 0], [-1, 1, 0]]
    millerbravaisreciprocal = [
        [1, 0, -1, 0],
        [0, 1, -1, 0],
        [1, -1, 0, 0],
        [-1, 0, 1, 0],
        [0, -1, 1, 0],
        [-1, 1, 0, 0],
    ]
    for (x, y) in zip(millerreciprocal, millerbravaisreciprocal)
        m = ReciprocalMiller(x)
        mb = ReciprocalMillerBravais(y)
        @test convert(ReciprocalMillerBravais, m) == mb
        @test convert(ReciprocalMiller, mb) == m
        @test convert(ReciprocalMiller, m) == m
        @test convert(ReciprocalMillerBravais, mb) == mb
        @test_throws MethodError convert(Miller, m)
        @test_throws MethodError convert(MillerBravais, mb)
        @test_throws MethodError convert(MillerBravais, m)
        @test_throws MethodError convert(Miller, mb)
    end
end

# From https://ssd.phys.strath.ac.uk/wp-content/uploads/Crystallographic_maths.pdf
@testset "Test conversion between reciprocal `Miller` and `MillerBravais`" begin
    millerreciprocal =
        [[1, 1, 0], [1, -2, 0], [-2, 1, 0], [-1, -1, 0], [-1, 2, 0], [2, -1, 0]]
    millerbravaisreciprocal = [
        [1, 1, -2, 0],
        [1, -2, 1, 0],
        [-2, 1, 1, 0],
        [-1, -1, 2, 0],
        [-1, 2, -1, 0],
        [2, -1, -1, 0],
    ]
    for (x, y) in zip(millerreciprocal, millerbravaisreciprocal)
        m = ReciprocalMiller(x)
        mb = ReciprocalMillerBravais(y)
        @test convert(ReciprocalMillerBravais, m) == mb
        @test convert(ReciprocalMiller, mb) == m
        @test convert(ReciprocalMiller, m) == m
        @test convert(ReciprocalMillerBravais, mb) == mb
        @test_throws MethodError convert(Miller, m)
        @test_throws MethodError convert(MillerBravais, mb)
        @test_throws MethodError convert(MillerBravais, m)
        @test_throws MethodError convert(Miller, mb)
    end
end
@test convert(MillerBravais, Miller([1, 0, 0])) == MillerBravais([2, -1, -1, 0])

@testset "Test `family`" begin
    @testset "Test 1" begin
        m = Miller(-1, 0, 1)
        @test Set(family(m)) == Set(
            Miller.([[1, 1, 1], [-1, 0, 1], [0, -1, 1], [-1, -1, 1], [1, 0, 1], [0, 1, 1]]),
        )
        mb = MillerBravais(-1, -1, 2, 3)
        @test Set(family(mb)) == Set(
            MillerBravais.([
                [1, 1, -2, 3],
                [-2, 1, 1, 3],
                [1, -2, 1, 3],
                [-1, -1, 2, 3],
                [2, -1, -1, 3],
                [-1, 2, -1, 3],
            ]),
        )
    end
    @testset "Test 2" begin
        m = ReciprocalMiller(0, 1, 0)
        @test Set(family(m)) == Set(
            ReciprocalMiller.([
                [1, 0, 0],
                [0, 1, 0],
                [1, -1, 0],
                [-1, 0, 0],
                [0, -1, 0],
                [-1, 1, 0],
            ]),
        )
        mb = ReciprocalMillerBravais(0, -1, 1, 0)
        @test Set(family(mb)) == Set(
            ReciprocalMillerBravais.([
                [1, 0, -1, 0],
                [0, 1, -1, 0],
                [1, -1, 0, 0],
                [-1, 0, 1, 0],
                [0, -1, 1, 0],
                [-1, 1, 0, 0],
            ]),
        )
    end
    @testset "Test 3" begin
        m = ReciprocalMiller(1, 1, 0)
        @test Set(family(m)) == Set(
            ReciprocalMiller.([
                [1, 1, 0],
                [1, -2, 0],
                [-2, 1, 0],
                [-1, -1, 0],
                [-1, 2, 0],
                [2, -1, 0],
            ]),
        )
        mb = ReciprocalMillerBravais(-2, 1, 1, 0)
        @test Set(family(mb)) == Set(
            ReciprocalMillerBravais.([
                [1, 1, -2, 0],
                [1, -2, 1, 0],
                [-2, 1, 1, 0],
                [-1, -1, 2, 0],
                [-1, 2, -1, 0],
                [2, -1, -1, 0],
            ]),
        )
    end
end


text:
module TestIterators

import DataFrames: DataFrame
import Volcanito: columns_tuple, tuple_iterator
import Test: @testset, @test

@testset "Iterators" begin
    df = DataFrame(a = [1, 2, 3], b = [2, 3, 4])

    @test columns_tuple(df, (:a, :b)) == (df[:a], df[:b])

    tuples = collect(tuple_iterator(df, (:a, :b)))
    @test tuples == [(1, 2), (2, 3), (3, 4)]
end

end


text:
using PyCall
using Conda

mutable struct OptimizerCmaEsDeap
    opt::Any
    lambda_::Any
    dim::Any
    chiN::Any
    mu::Any
    weights::Any
    mueff::Any
    cc::Any
    cs::Any
    ps::Any
    pc::Any
    centroid::Any
    update_count::Any
    ccov1::Any
    ccovmu::Any
    C::Any
    sigma::Any
    damps::Any
    diagD::Any
    B::Any
    BD::Any
    genomes::Any

    function OptimizerCmaEsDeap(individual_size::Int, configuration::OrderedDict)
        scriptdir = @__DIR__
        pushfirst!(PyVector(pyimport("sys")."path"), scriptdir)
        optimizer_deap = pyimport("cma_es_deap")
        opt = optimizer_deap.OptimizerCmaEsDeap(individual_size, configuration)

        lambda_ = opt.strategy.lambda_
        dim = opt.strategy.dim
        chiN = opt.strategy.chiN
        mu = opt.strategy.mu
        weights = opt.strategy.weights
        mueff = opt.strategy.mueff
        cc = opt.strategy.cc
        cs = opt.strategy.cs
        ps = opt.strategy.ps
        pc = opt.strategy.pc
        centroid = opt.strategy.centroid
        update_count = opt.strategy.update_count
        ccov1 = opt.strategy.ccov1
        ccovmu = opt.strategy.ccovmu
        C = opt.strategy.C
        sigma = opt.strategy.sigma
        damps = opt.strategy.damps
        diagD = opt.strategy.diagD
        B = opt.strategy.B
        BD = opt.strategy.BD
        eigenvectors = opt.strategy.eigenvectors
        indx = opt.strategy.indx
        genomes = zeros(opt.strategy.lambda_, individual_size)

        optimizer = new(opt, lambda_, dim, chiN, mu, weights, mueff, cc, cs, ps, pc, centroid, update_count, ccov1, ccovmu, C, sigma, damps, diagD, B, BD, genomes)

        return optimizer, eigenvectors, indx .+ 1
    end
end

function ask(optimizer)

    genomes_list, population_size, individual_size, strategy = optimizer.opt.ask()

    # The genomes need to be reshaped into a MxN matrix.
    for i = 1:population_size
        for j = 1:individual_size
            optimizer.genomes[i, j] = (genomes_list[i])[j]
        end
    end

    return optimizer.genomes, strategy.randoms
end

function tell(optimizer, rewards)

    strategy = optimizer.opt.tell(rewards)

    optimizer.centroid = strategy.centroid
    optimizer.ps = strategy.ps
    optimizer.update_count = strategy.update_count
    optimizer.pc = strategy.pc
    optimizer.C = strategy.C
    optimizer.sigma = strategy.sigma
    optimizer.diagD = strategy.diagD
    optimizer.B = strategy.B
    optimizer.BD = strategy.BD

    return strategy.eigenvectors, strategy.indx .+ 1

end


text:
using GLMdiagnostics
using Test
using Distributions
using Random
using DataFrames
using StatsBase
using GLM
using Gadfly
using Statistics
using Compose

###########################################
# - Create simulated data to test functions
###########################################

# Set reproducible seed
Random.seed!(123)

# Simulate data from a Normal distribution
d = Normal(1, 2)
y = rand(d, 100)

# Assign grouping variables (5 groups of 20)
x = repeat([1, 2, 3, 4, 5],
            inner = 20,
            outer = 1)

df = hcat(y, x)
df = convert(DataFrame, df)
colnames = ["y","x"]
rename!(df, Symbol.(colnames))

##########################################
# - Run linear model (one-way ANOVA)
##########################################

# Run model
modAov = lm(@formula(y ~ 1 + x), df)

# Make qqplot for ANOVA model
anova_qqplot = plot_qq(modAov)

max_stars_repo_name: guysutton/GLMdiagnostics

max_issues_repo_name: guysutton/GLMdiagnostics

max_forks_repo_name: guysutton/GLMdiagnostics


text:
for frange in (mrange,srange)
    @testset "reverse-$frange" begin
        @test reverse(reverse(frange(1, 10))) == 1:10
        @test reverse(reverse(frange(typemin(Int), typemax(Int)))) == typemin(Int):typemax(Int)
        @test reverse(reverse(frange(typemin(Int), step=2, stop=typemax(Int)))) == typemin(Int):2:typemax(Int)
    end
end

@testset "reverse!" begin
    r = StepMRange(1, 1, 10)
    @test reverse!(r) == StepMRange(10, -1, 1)
    @test reverse!(r) == StepMRange(1, 1, 10)

    r = LinMRange(1, 10, 10)
    @test reverse!(r) == LinMRange(10, 1, 10)
    @test reverse!(r) == LinMRange(1, 10, 10)

    r = StepMRangeLen(1, 1, 10)
    @test reverse!(r) == StepMRangeLen(10, -1, 10)
    @test reverse!(r) == StepMRangeLen(1, 1, 10)
end


text:
using MetadataArrays, Test

@testset "MetadataVector" begin
    v = [1, 3, 5, 4]
    sv = MetadataVector(v, "Numbers")
    @test length(sv) == 4
    @test collect(sv) == parent(sv) == copy(sv)
    @test metadata(sv) == "Numbers"
    s_view1 = view(sv, 1:3)
    s_view2 = view(s_view1, 1:2)
    @test metadata(s_view1) == metadata(s_view2) == "Numbers"
    sv[3] = 50
    @test all(sv .== [1, 3, 50, 4])
    @test eltype(sv) == Int
    @test metadata(similar(sv)) == metadata(sv)
end

@testset "MetadataArray" begin
    v = [1 2; 3 4; 5 4]
    sv = MetadataArray(v, Dict(1 => "Not at all", 10 => "A lot"))
    @test size(sv) == size(v)
    @test sv[6] == parent(sv)[6]
    @test collect(sv) == parent(sv) == v == copy(sv)
    @test metadata(sv) == Dict(1 => "Not at all", 10 => "A lot")
    sv[3, 2] = 50
    @test all(sv .== [1 2; 3 4; 5 50])
    @test parent(sv[1:2, 1]) == parent(sv)[1:2, 1]
    @test metadata(sv[1:2, 1]) == metadata(sv)
    @test metadata(rand(2)) === nothing
    @test eltype(sv) == Int

    s = MetadataArray(rand(4), "test")
    sr = reshape(s, 2, 2)
    @test parent(sr) == reshape(parent(s), 2, 2)
    @test metadata(sr) == metadata(s)
    @test axes(sr) == (Base.OneTo(2), Base.OneTo(2))
    s2 = similar(s, (2, 3))
    @test size(s2) == (2, 3)
    @test metadata(s2) == metadata(s)

    x1 = rand(4, 4)
    m1 = MetadataArray(x1, "something")
    @test IndexStyle(m1) == IndexLinear()
    x2 = view(x1, 1:2, 1:2)
    m2 = MetadataArray(x2, "something")
    @test IndexStyle(m2) == IndexCartesian()
end


text:
function formulation(
    i::UniformMatroidInstance{Float64, Maximise},
    f::DefaultLinearFormulation
)
    dim = dimension(i)

    model = Model(f.solver)
    x = @variable(model, [1:dim], Bin)
    @objective(model, Max, sum(x[j] * i.rewards[j] for j in 1:dim))
    @constraint(model, sum(x) <= i.m)

    set_silent(model)

    return model, x
end

function solve(
    i::UniformMatroidInstance{Float64, Maximise},
    f::DefaultLinearFormulation
)
    m, x = formulation(i, f)
    optimize!(m)
    return UniformMatroidSolution(i, findall(JuMP.value.(x) .>= 0.5))
end


text:
module Cache

["Caching of docstrings and package metadata."]

export getraw, getparsed, getmeta, clear!, objects

using Compat; import Compat.String

import ..Utilities
import ..Collector: Collector, PackageData, ModuleData
import ..Formats

include("types.jl")
include("utilities.jl")
include("storage.jl")
include("interface.jl")

end


text:
# Adapted from:
# SOSDEMO1 --- Sum of Squares Test
# Section 3.1 of SOSTOOLS User's Manual

facts("SOSDEMO1") do
for solver in sdp_solvers
context("With solver $(typeof(solver))") do
  @polyvar x y

  m = Model(solver = solver)

  # Example 2.4 of
  # P. Parrilo and A. Jadbabaie
  # "Approximation of the joint spectral radius using sum of squares."
  # Linear Algebra and its Applications, Elsevier, 2008, 428, 2385-2402

  p = 2*x^4 + 2*x^3*y - x^2*y^2 + 5*y^4

  soscon = @polyconstraint m p >= 0

  status = solve(m)

  @fact status --> :Optimal

  q = getslack(soscon)
  Q = getmat(q)
  @fact issymmetric(Q) --> true
  @fact isapprox(Q[1, 1], 2, rtol=1e-5) --> true
  @fact isapprox(Q[1, 2], 1, rtol=1e-5) --> true
  @fact isapprox(Q[3, 3], 5, rtol=1e-5) --> true
  @fact abs(Q[2, 3]) < 1e-5 --> true
  @fact isapprox(Q[2, 2] + 2Q[1, 3], -1, rtol=1e-5) --> true
  sosdec = SOSDecomposition(q)
  @fact isapprox(sum(sosdec.ps.^2), p; rtol=1e-4) --> true

  M = Model(solver = solver)

  p = 4*x^4*y^6 + x^2 - x*y^2 + y^2

  soscon = @polyconstraint M p >= 0

  status = solve(M)

  @fact status --> :Optimal

  # p should be MatPolynomial([1, 0, -1/2, 0, -1, 1, 0, -2/3, 0, 4/3, 0, 0, 2, 0, 4], [y, x, x*y, x*y^2, x^2*y^3])
end; end; end

max_stars_repo_name: polislizarralde/jsos

max_issues_repo_name: polislizarralde/jsos

max_forks_repo_name: polislizarralde/jsos


text:
using Query
using DataFrames

df = DataFrame(name=["John", "Sally", "Kirk"], age=[23., 42., 59.], children=[3,5,2])

q = @from i in df begin
    @select i.name=>i.children
    @collect Dict
end

println(q)


text:
import Term: Panel, TextBox
import Term.Layout: PlaceHolder

@testset "\e[34mPanel - no content" begin
    for style in ("default", "red", "on_blue")
        testpanel(Panel(fit = true, style = style), 3, 2)
        testpanel(Panel(), 88, 2)
        testpanel(Panel(width = 12, height = 4, style = style), 12, 4)
    end
end

@testset "\e[34mPanel - fit basic shape" begin
    testpanel(Panel("MYTEXT"; width = 60, height = 10), 60, 10)
    testpanel(Panel("MYTEXT"; width = 60), 60, 3)
    testpanel(Panel("MYTEXT"), 12, 3)

    testpanel(Panel("MYTEXT"; width = 60, height = 10, padding = (4, 4, 2, 2)), 60, 10)
    testpanel(Panel("MYTEXT"; width = 60, padding = (4, 4, 2, 2)), 60, 7)
    testpanel(Panel("MYTEXT"; padding = (4, 4, 2, 2)), 16, 7)
end

@testset "\e[34mPanel - fit overflow" begin
    testpanel(
        Panel(
            Panel("MYTEXT"^50; width = 60, height = 10) /
            Panel("MYTEXT"^50; width = 60, height = 10),
        ),
        console_width() - 1,
        nothing,
    )
end

@testset "\e[34mPANEL - fit - measure" begin
    for style in ("default", "red", "on_blue")
        # ----------------------------- text only content ---------------------------- #
        _kw = (fit = true, style = style)
        testpanel(Panel("t"; _kw...), 7, 3)
        testpanel(Panel("test"; _kw...), 10, 3)
        testpanel(Panel("1234\n123456789012"; _kw...), 18, 4)
        testpanel(Panel("나랏말싸미 듕귁에 달아"; _kw...), 28, 3)
        testpanel(Panel("나랏말싸미 듕귁에 달아\n1234567890123456789012"; _kw...), 28, 4)
        testpanel(Panel("."^500; _kw...), displaysize(stdout)[2] - 1, nothing)

        # ------------------------------- nested panels ------------------------------ #
        testpanel(Panel(Panel("test"; _kw...); _kw...), 16, 5)

        testpanel(
            Panel(Panel(Panel("."; _kw...); _kw...); fit = true, style = style),
            19,
            7,
        )
    end
end

@testset "\e[34mPANEL - nofit - measure" begin
    for style in ("default", "red", "on_blue")
        testpanel(Panel("t"; style = style, fit = false), 88, 3)
        testpanel(Panel("test"; style = style, fit = false), 88, 3)
        testpanel(Panel("1234\n123456789012"; style = style, fit = false), 88, 4)
        testpanel(Panel("나랏말싸미 듕귁에 달아"; style = style, fit = false), 88, 3)
        testpanel(
            Panel("나랏말싸미 듕귁에 달아\n1234567890123456789012"; style = style, fit = false),
            88,
            4,
        )
        for justify in (:left, :center, :right)
            # ----------------------------- text only content ---------------------------- #
            _kw = (justify = justify, style = style)
            _nofit = (; fit = false, _kw...)
            testpanel(Panel("."^1500; _nofit...), 88, 21)

            # ------------------------------- nested panels ------------------------------ #
            testpanel(Panel(Panel("test"); fit = true, _kw...), 16, 5)

            testpanel(Panel(Panel(Panel("."); _nofit...); _nofit...), 88, 7)

            testpanel(Panel(Panel("."^250); _nofit...), 88, 6)

            testpanel(Panel(Panel("test"; _kw...); fit = false), 88, 5)

            testpanel(Panel(Panel(Panel("."; _kw...); _kw...); fit = false), 88, 7)

            testpanel(Panel(Panel("."^250; _kw...); fit = false), 88, 6)

            testpanel(Panel(Panel("t1"; _kw...), Panel("t2"; _kw...); fit = false), 88, 8)

            _kw = (fit = false, width = 30, height = 8)
            testpanel(Panel(Panel("test"; width = 22); _kw...), 30, 8)

            testpanel(Panel(Panel("test"; width = 42); _kw...), 30, 8)

            testpanel(
                Panel(
                    Panel("test"; width = 42, height = 12);
                    fit = false,
                    width = 30,
                    height = 8,
                ),
                30,
                8,
            )
        end
    end
end

@testset "\e[34mPANEL - FIT - measure" begin
    testpanel(Panel("t"; fit = true), 7, 3)
    testpanel(Panel("test"; fit = true), 10, 3)
    testpanel(Panel("1234\n123456789012"; fit = true), 18, 4)
    testpanel(Panel("나랏말싸미 듕귁에 달아"; fit = true), 28, 3)
    testpanel(Panel("나랏말싸미 듕귁에 달아\n1234567890123456789012"; fit = true), 28, 4)
    for justify in (:left, :center, :right)
        # ----------------------------- text only content ---------------------------- #
        _kw = (fit = true, justify = justify)
        testpanel(Panel("."^1500; _kw...), console_width() - 1, nothing)

        # ------------------------------- nested panels ------------------------------ #
        testpanel(Panel(Panel("test"); _kw...), 16, 5)

        testpanel(Panel(Panel(Panel("."); _kw...); _kw...), 19, 7)

        testpanel(Panel(Panel("."^250); _kw...), console_width() - 1, nothing)

        testpanel(Panel(Panel("test"; justify = justify); fit = true), 16, 5)

        testpanel(
            Panel(Panel(Panel("."; justify = justify); justify = justify); fit = true),
            19,
            7,
        )

        testpanel(
            Panel(Panel("."^250; justify = justify); fit = true),
            console_width() - 1,
            nothing,
        )

        testpanel(
            Panel(
                Panel("t1"; justify = justify),
                Panel("t2"; justify = justify);
                fit = true,
            ),
            14,
            8,
        )

        _kw = (fit = true, width = 30, height = 8)
        testpanel(Panel(Panel("test"; width = 22); _kw...), 30, 8)

        testpanel(Panel(Panel("test"; width = 42); _kw...), 48, 8)

        testpanel(Panel(Panel("test"; width = 42, height = 12); _kw...), 48, 14)
    end
end

@testset "PANEL - centered title style" begin
    @test string(
        Panel(title = "test", title_justify = :left, title_style = "italic red"),
    ) ==
          "\e[22m╭──── \e[3m\e[31mtest\e[23m\e[39m\e[22m\e[22m ────────────────────────────────────────────────────────────────────────────╮\e[22m\e[0m\e[22m\n\e[22m╰──────────────────────────────────────────────────────────────────────────────────────╯\e[22m\e[0m"
    @test string(
        Panel(title = "test", title_justify = :center, title_style = "italic red"),
    ) ==
          "\e[22m╭────────────────────────────────────────\e[22m \e[3m\e[31mtest\e[23m\e[39m\e[22m\e[22m ────────────────────────────────────────╮\e[22m\e[0m\n\e[22m╰──────────────────────────────────────────────────────────────────────────────────────╯\e[22m\e[0m"
    @test string(
        Panel(title = "test", title_justify = :right, title_style = "italic red"),
    ) ==
          "\e[22m╭───────────────────────────────────────────────────────────────────────────── \e[3m\e[31mtest\e[23m\e[39m\e[22m\e[22m ───╮\e[22m\e[0m\e[22m\n\e[22m╰──────────────────────────────────────────────────────────────────────────────────────╯\e[22m\e[0m"

    @test string(Panel(title = "test", title_justify = :left)) ==
          "\e[22m╭──── test\e[22m ────────────────────────────────────────────────────────────────────────────╮\e[22m\e[0m\e[22m\n\e[22m╰──────────────────────────────────────────────────────────────────────────────────────╯\e[22m\e[0m"
    @test string(Panel(title = "test", title_justify = :center)) ==
          "\e[22m╭────────────────────────────────────────\e[22m test\e[22m ────────────────────────────────────────╮\e[22m\e[0m\n\e[22m╰──────────────────────────────────────────────────────────────────────────────────────╯\e[22m\e[0m"
    @test string(Panel(title = "test", title_justify = :right)) ==
          "\e[22m╭───────────────────────────────────────────────────────────────────────────── test\e[22m ───╮\e[22m\e[0m\e[22m\n\e[22m╰──────────────────────────────────────────────────────────────────────────────────────╯\e[22m\e[0m"

    p = Panel(
        title = "test",
        title_justify = :left,
        subtitle = "aaaaa",
        subtitle_justify = :right,
        width = 22,
    )
    @test string(p) ==
          "\e[22m╭──── test\e[22m ──────────╮\e[22m\e[0m\e[22m\n\e[22m╰────────── aaaaa\e[22m ───╯\e[22m\e[0m\e[22m\e[0m"
end

@testset "PANEL - compare to string" begin
    pts = """
Lorem{red} ipsum dolor s{/red}it amet, consectetur adipiscing elit,
ed do eiusmod tempor {bold blue}incididu{underline}nt ut labore et dolore magna aliqua. Ut enim ad minim
veniam, quis nos{/underline}trud exercitation ullamco laboris nisi ut aliquip ex
ea commodo consequat. Duis aute {/bold blue}{red on_black}irure dolor in reprehenderit
in voluptate velit esse cillum dolore eu fugiat nulla
pariatur. Excepteur sint occaecat{/red on_black} cupidatat non proident,
sunt in {green}culpa qui officia{/green} deserunt mollit anim
id est laborum."""

    pts = replace(pts, "\n" => "")

    pts2 = replace(
        """
Lorem ipsum {bold}dolor sit{/bold} amet, consectetur adipiscing elit,
ed do e{red}iusmod tempor incididunt{/red} ut {bold}labore et {underline}dolore{/underline} magna aliqua.{/bold} Ut enim ad minim
veniam, quis{green} nostrud exercitation {on_black}ullamco laboris nisi ut aliquip ex {/on_black}
ea commodo consequat.{blue} Duis aute irure dolor in{/blue} reprehenderit
in voluptate velit{/green} esse {italic}cillum dolore{/italic}{red} eu{/red}{italic green} fugiat {/italic green}nulla
pariatur. Excepteur{red} sint{/red}{blue} occaecat cupidatat {/blue}non proident,
sunt in culpa qui {italic}officia{/italic} deserunt mollit anim
id est laborum.""",
        "\n" => "",
    )

    @test string(Panel(pts2; fit = false, width = 44)) ==
          "\e[22m╭──────────────────────────────────────────╮\e[22m\n\e[0m\e[22m│\e[22m\e[0m  Lorem ipsum \e[1mdolor sit\e[22m amet,             \e[0m\e[22m│\e[22m\e[0m\n\e[0m\e[22m│\e[22m\e[0m  consectetur adipiscing elit,ed do       \e[0m\e[22m│\e[22m\e[0m\n\e[0m\e[22m│\e[22m\e[0m  e\e[31miusmod tempor incididunt\e[39m\e[22m ut \e[1mlabore     \e[0m\e[22m│\e[22m\e[0m\n\e[0m\e[22m│\e[22m\e[0m  et \e[4mdolore\e[24m\e[1m magna aliqua.\e[22m Ut enim ad      \e[0m\e[22m│\e[22m\e[0m\n\e[0m\e[22m│\e[22m\e[0m  minimveniam, quis\e[32m nostrud               \e[0m\e[22m│\e[22m\e[0m\n\e[0m\e[22m│\e[22m\e[0m  exercitation \e[40mullamco laboris nisi ut    \e[0m\e[22m│\e[22m\e[0m\n\e[0m\e[22m│\e[22m\e[0m  aliquip ex \e[49m\e[32mea commodo consequat.\e[34m Duis   \e[0m\e[22m│\e[22m\e[0m\n\e[0m\e[22m│\e[22m\e[0m  aute irure dolor in\e[39m\e[32m reprehenderit in    \e[0m\e[22m│\e[22m\e[0m\n\e[0m\e[22m│\e[22m\e[0m  voluptate velit\e[39m\e[22m esse \e[3mcillum dolore\e[23m\e[22m\e[31m eu\e[39m\e[22m\e[3m\e[32m   \e[0m\e[22m│\e[22m\e[0m\n\e[0m\e[22m│\e[22m\e[0m  fugiat \e[23m\e[22m\e[39m\e[3mnulla pariatur. Excepteur\e[31m sint\e[39m\e[3m\e[34m   \e[0m\e[22m│\e[22m\e[0m\n\e[0m\e[22m│\e[22m\e[0m  occaecat cupidatat \e[39m\e[3mnon proident, sunt   \e[0m\e[22m│\e[22m\e[0m\n\e[0m\e[22m│\e[22m\e[0m  in culpa qui \e[3mofficia\e[23m\e[3m deserunt mollit    \e[0m\e[22m│\e[22m\e[0m\n\e[0m\e[22m│\e[22m\e[0m  anim id est laborum.                    \e[0m\e[22m│\e[22m\e[0m\n\e[22m╰──────────────────────────────────────────╯\e[22m\e[0m"

    @test string(Panel(pts2; fit = false, width = 51)) ==
          "\e[22m╭─────────────────────────────────────────────────╮\e[22m\n\e[0m\e[22m│\e[22m\e[0m  Lorem ipsum \e[1mdolor sit\e[22m amet, consectetur        \e[0m\e[22m│\e[22m\e[0m\n\e[0m\e[22m│\e[22m\e[0m  adipiscing elit,ed do e\e[31miusmod tempor           \e[0m\e[22m│\e[22m\e[0m\n\e[0m\e[22m│\e[22m\e[0m  incididunt\e[39m\e[22m ut \e[1mlabore et \e[4mdolore\e[24m\e[1m magna aliqua.\e[22m   \e[0m\e[22m│\e[22m\e[0m\n\e[0m\e[22m│\e[22m\e[0m  Ut enim ad minimveniam, quis\e[32m nostrud           \e[0m\e[22m│\e[22m\e[0m\n\e[0m\e[22m│\e[22m\e[0m  exercitation \e[40mullamco laboris nisi ut aliquip   \e[0m\e[22m│\e[22m\e[0m\n\e[0m\e[22m│\e[22m\e[0m  ex \e[49m\e[32mea commodo consequat.\e[34m Duis aute irure       \e[0m\e[22m│\e[22m\e[0m\n\e[0m\e[22m│\e[22m\e[0m  dolor in\e[39m\e[32m reprehenderit in voluptate velit\e[39m\e[22m      \e[0m\e[22m│\e[22m\e[0m\n\e[0m\e[22m│\e[22m\e[0m  esse \e[3mcillum dolore\e[23m\e[22m\e[31m eu\e[39m\e[22m\e[3m\e[32m fugiat \e[23m\e[22m\e[39m\e[3mnulla pariatur.   \e[0m\e[22m│\e[22m\e[0m\n\e[0m\e[22m│\e[22m\e[0m  Excepteur\e[31m sint\e[39m\e[3m\e[34m occaecat cupidatat \e[39m\e[3mnon          \e[0m\e[22m│\e[22m\e[0m\n\e[0m\e[22m│\e[22m\e[0m  proident, sunt in culpa qui \e[3mofficia\e[23m\e[3m deserunt   \e[0m\e[22m│\e[22m\e[0m\n\e[0m\e[22m│\e[22m\e[0m  mollit anim id est laborum.                    \e[0m\e[22m│\e[22m\e[0m\n\e[22m╰─────────────────────────────────────────────────╯\e[22m\e[0m"

    @test string(Panel(pts2; fit = false, width = 67)) ==
          "\e[22m╭─────────────────────────────────────────────────────────────────╮\e[22m\n\e[0m\e[22m│\e[22m\e[0m  Lorem ipsum \e[1mdolor sit\e[22m amet, consectetur adipiscing elit,ed     \e[0m\e[22m│\e[22m\e[0m\n\e[0m\e[22m│\e[22m\e[0m  do e\e[31miusmod tempor incididunt\e[39m\e[22m ut \e[1mlabore et \e[4mdolore\e[24m\e[1m magna         \e[0m\e[22m│\e[22m\e[0m\n\e[0m\e[22m│\e[22m\e[0m  aliqua.\e[22m Ut enim ad minimveniam, quis\e[32m nostrud exercitation      \e[0m\e[22m│\e[22m\e[0m\n\e[0m\e[22m│\e[22m\e[0m  \e[40mullamco laboris nisi ut aliquip ex \e[49m\e[32mea commodo consequat.\e[34m       \e[0m\e[22m│\e[22m\e[0m\n\e[0m\e[22m│\e[22m\e[0m  Duis aute irure dolor in\e[39m\e[32m reprehenderit in voluptate velit\e[39m\e[22m      \e[0m\e[22m│\e[22m\e[0m\n\e[0m\e[22m│\e[22m\e[0m  esse \e[3mcillum dolore\e[23m\e[22m\e[31m eu\e[39m\e[22m\e[3m\e[32m fugiat \e[23m\e[22m\e[39m\e[3mnulla pariatur. Excepteur\e[31m sint\e[39m\e[3m\e[34m    \e[0m\e[22m│\e[22m\e[0m\n\e[0m\e[22m│\e[22m\e[0m  occaecat cupidatat \e[39m\e[3mnon proident, sunt in culpa qui \e[3mofficia\e[23m\e[3m     \e[0m\e[22m│\e[22m\e[0m\n\e[0m\e[22m│\e[22m\e[0m  deserunt mollit anim id est laborum.                           \e[0m\e[22m│\e[22m\e[0m\n\e[22m╰─────────────────────────────────────────────────────────────────╯\e[22m\e[0m"

    circle = """
          oooo
       oooooooooo
      oooooooooooo
      oooooooooooo
       oooooooooo
          oooo    """
    p = Panel(circle; fit = true, padding = (0, 0, 0, 0))
    @test string(p) ==
          "\e[22m╭────────────╮\e[22m\n\e[0m\e[22m│\e[22m\e[0m    oooo    \e[0m\e[22m│\e[22m\e[0m\n\e[0m\e[22m│\e[22m\e[0m oooooooooo \e[0m\e[22m│\e[22m\e[0m\n\e[0m\e[22m│\e[22m\e[0moooooooooooo\e[0m\e[22m│\e[22m\e[0m\n\e[0m\e[22m│\e[22m\e[0moooooooooooo\e[0m\e[22m│\e[22m\e[0m\n\e[0m\e[22m│\e[22m\e[0m oooooooooo \e[0m\e[22m│\e[22m\e[0m\n\e[0m\e[22m│\e[22m\e[0m    oooo    \e[0m\e[22m│\e[22m\e[0m\n\e[22m╰────────────╯\e[22m\e[0m"

    p = Panel(
        "test"^25;
        width = 49,
        height = 22,
        box = :HEAVY,
        justify = :center,
        title = "description",
        title_justify = :center,
        title_style = "bold bright_blue",
        padding = (2, 4, 4, 2),
        fit = false,
    )
    @test size(p.measure) == (49, 22)

    p = Panel(
        "test"^25;
        width = 49,
        height = 22,
        box = :SIMPLE,
        justify = :center,
        title = "description",
        title_justify = :center,
        title_style = "bold bright_blue",
        padding = (2, 4, 4, 2),
        fit = false,
    )
    @test size(p.measure) == (49, 22)

    p = Panel(
        PlaceHolder(60, 30);
        width = 49,
        height = 22,
        box = :HEAVY,
        justify = :center,
        title = "description",
        title_justify = :center,
        title_style = "bold bright_blue",
        padding = (2, 0, 2, 0),
        fit = false,
    )
    @test size(p.measure) == (49, 22)
end

@testset "\e[34mPanel + renderables" begin
    testpanel(Panel(RenderableText("x" .^ 5)), 11, 3)

    testpanel(Panel(RenderableText("x" .^ 500); fit = false), 88, 9)

    testpanel(Panel(RenderableText("x" .^ 5); fit = true), 11, 3)

    testpanel(
        Panel(RenderableText("x" .^ 500); fit = true),
        displaysize(stdout)[2] - 1,
        nothing,
    )
end

@testset "\e[34mPANEL - titles" begin
    style = "red"
    for fit in (true, false)
        for justify in (:left, :center, :right)
            testpanel(
                Panel(
                    "."^50;
                    title = "test",
                    title_style = style,
                    title_justify = justify,
                    subtitle = "subtest",
                    subtitle_style = style,
                    subtitle_justify = justify,
                    fit = fit,
                ),
                fit ? nothing : 88,
                3,
            )
        end
    end
end

@testset "\e[34mPanel - padding" begin
    p = Panel("."^24; padding = (4, 4, 2, 2), fit = false)
    testpanel(p, 88, 7)

    p = Panel("."^24; padding = (4, 4, 2, 2), fit = true)
    testpanel(p, 34, 7)
end


text:
# ==========================================================================================================
# Backend funcion

formula(model::MixedModel) = model.formula

lhs_no_intercept(lhs::String) = Set([lhs])
lhs_no_intercept(lhs) = Set(filter!(!=("(Intercept)"), lhs))

"""
    calcdof(model::LinearMixedModel)

Calculate degree of freedom of factors and residuals for linear mixed effect models
DOF of residuals are estimated by between-within method: dofᵢ = nobsᵢ - dofᵢ₋₁ - nfixᵢ

Reference algorithm: [GLMM FAQ](https://bbolker.github.io/mixedmodels-misc/glmmFAQ.html#why-doesnt-lme4-display-denominator-degrees-of-freedomp-values-what-other-options-do-i-have) for details.
"""
function calcdof(model::LinearMixedModel)
    randoms = collect(formula(model).rhs) # ranef formula terms
    fixs = popfirst!(randoms) # fixef formula terms
    nfac = length(fixs.terms) # number of factors
    isempty(randoms) && (return repeat([nobs(model) - nfac], nfac))
    reterms = reverse(model.reterms) # Vector of ReMat
    fixef = model.feterm.cnames # fixef name
    nfix = length(fixef) # number of fix effects

    # Determine affected fix effects for each random effects
    affectfixef = Dict{String, Set{String}}()
    for pair in randoms
        lhs = coefnames(pair.lhs)
        lhs = lhs_no_intercept(lhs)
        rhs = coefnames(pair.rhs, Val(:anova))
        haskey(affectfixef, rhs) ? union!(affectfixef[rhs], lhs) : push!(affectfixef, rhs => lhs)
    end

    # Determines if fix effects vary at each random effects group
    within = zeros(Bool, nfix, size(reterms, 1) + 2)
    for (reid, remat) in enumerate(reterms)
        levels = unique(remat.refs)
        within[:, reid + 1] .= map(1:nfix) do fixid
            all(levels) do level
                vals = view(model.Xymat.wtxy, findall(==(level), remat.refs), fixid)
                val = first(vals)
                for i in vals
                    i != val && (return false)
                end
                true
            end
        end
    end
    within[:, 1] = repeat([false], nfix) # population level
    within[1, 1] = true
    within[:, end] = repeat([true], nfix) # observation level

    # Turn affected fix effects into true
    ranef = ["", coefnames.(getproperty.(reterms, :trm), Val(:anova))..., ""]
    for (key, value) in affectfixef
        id = findfirst(in(value), fixef)
        isnothing(id) || (within[id, findfirst(==(key), ranef)] = true)
    end
    # Find first true for each fix effect, aka level
    level = mapslices(within, dims = 2) do fixef
        findfirst(fixef)
    end
    # Number of fix effects per level
    nfixperlv = zeros(Int, length(ranef))
    for i in level
        nfixperlv[i] += 1
    end

    # dfᵢ = nobsᵢ - dfᵢ₋₁ - nfixᵢ
    n0 = 1
    n = nobs(model)
    nobsperlv = (n0, length.(getproperty.(reterms, :levels))..., n) # number of random effects observation per level
    ndiff = _diff(nobsperlv)
    dfperlv = ndiff .- nfixperlv[2:end]
    dfperlv = (last(dfperlv), dfperlv...)

    # Assign df to each factors based on the level
    assign = asgn(fixs)
    offset = 0
    intercept = first(assign) == 1
    intercept || (offset = 1)
    dfr = ntuple(last(assign) - offset) do i
        dfperlv[level[findlast(==(i + offset), assign)]] # use the last one to avoid wrong assign of the first arguments
    end
    df = dof(assign)
    intercept || popfirst!(df)
    tuple(df...), dfr
end

nestedmodels(::Type{<: LinearMixedModel}, f::FormulaTerm, tbl;
                null::Bool = true,
                wts = [],
                contrasts = Dict{Symbol, Any}(),
                progress::Bool = true,
                REML::Bool = false) =
    nestedmodels(fit(LinearMixedModel, f, tbl; wts, contrasts, progress, REML); null)

function nestedmodels(model::LinearMixedModel; null::Bool = true, kwargs...)
    f = formula(model)
    range = null ? (0:length(f.rhs[1].terms) - 1) : (1:length(f.rhs[1].terms) - 1)
    assign = asgn(first(f.rhs))
    REML = model.optsum.REML
    pivot = model.feterm.piv
    X = copy(model.X)
    y = copy(model.y)
    T = promote_type(Float64, eltype(y))
    size(X, 2) > 0 && (X = pivot == collect(1:size(X, 2)) ? X : X[:, pivot])
    reterms = deepcopy(model.reterms)
    sqrtwts = copy(model.sqrtwts)
    models = map(range) do id
        subf = subformula(f.lhs, f.rhs, id)
        cnames = coefnames(first(subf.rhs))
        # modify X by assign
        select = findall(x->x <= id, assign)
        subX = X[:, select]
        feterms = MixedModels.FeTerm{T}[]
        push!(feterms, MixedModels.FeTerm(subX, isa(cnames, String) ? [cnames] : collect(cnames)))
        feterm = only(feterms)
        Xy = MixedModels.FeMat(feterm, vec(y))
        reweight!(Xy, sqrtwts)
        A, L = createAL(reterms, Xy)
        lbd = foldl(vcat, lowerbd(c) for c in reterms)
        θ = foldl(vcat, getθ(c) for c in reterms)
        optsum = OptSummary(θ, lbd, :LN_BOBYQA, ftol_rel = T(1.0e-12), ftol_abs = T(1.0e-8))
        fill!(optsum.xtol_abs, 1.0e-10)
        lmm = LinearMixedModel(
            subf,
            reterms,
            Xy,
            feterm,
            sqrtwts,
            model.parmap,
            (n = length(y), p = feterm.rank, nretrms = length(reterms)),
            A,
            L,
            optsum,
            )
        fit!(lmm; REML, progress = true)
        lmm
    end
    (models..., model)
end

# For nestedmodels
isnullable(::LinearMixedModel) = true

# Specialized dof_residual
dof_residual(aov::AnovaResult{<: MixedModel, FTest}) = aov.tests.resdof


text:
# Autogenerated wrapper script for gdk_pixbuf_jll for aarch64-apple-darwin
export gdk_pixbuf_loaders_dir, gdk_pixbuf_query_loaders, libgdkpixbuf

using Glib_jll
using JpegTurbo_jll
using libpng_jll
using Libtiff_jll
JLLWrappers.@generate_wrapper_header("gdk_pixbuf")
JLLWrappers.@declare_file_product(gdk_pixbuf_loaders_dir)
JLLWrappers.@declare_library_product(libgdkpixbuf, "@rpath/libgdk_pixbuf-2.0.0.dylib")
JLLWrappers.@declare_executable_product(gdk_pixbuf_query_loaders)
function __init__()
    JLLWrappers.@generate_init_header(Glib_jll, JpegTurbo_jll, libpng_jll, Libtiff_jll)
    JLLWrappers.@init_file_product(
        gdk_pixbuf_loaders_dir,
        "lib/gdk-pixbuf-2.0/2.10.0/loaders",
    )

    JLLWrappers.@init_library_product(
        libgdkpixbuf,
        "lib/libgdk_pixbuf-2.0.0.dylib",
        RTLD_LAZY | RTLD_DEEPBIND,
    )

    JLLWrappers.@init_executable_product(
        gdk_pixbuf_query_loaders,
        "bin/gdk-pixbuf-query-loaders",
    )

    JLLWrappers.@generate_init_footer()
end  # __init__()


text:
"""
    select!(node::Node, args..., <keyword arguments>)
    select(node::Node, args..., <keyword arguments>)

Delete all attributes not selected in `args...`, and optionally apply transformations on the
fly on the selected variables. This function works similarly to [`transform!`](@ref) except
it keeps only the selected variables, while [`transform!`](@ref) add new variables.

See the documentation of [`transform!`](@ref) for more details on the format of `args` and on
how to use the arguments.

This function adds one more form to `args...`: a variable name to just select a variable.

# Examples

```julia
file = joinpath(dirname(dirname(pathof(MultiScaleTreeGraph))),"test","files","simple_plant.mtg")
mtg = read_mtg(file)

select!(mtg, :Length => (x -> x / 10) => :length_m, :Width, ignore_nothing = true)
```
"""
function select!(
    mtg::Node,
    args...;
    scale = nothing,
    symbol = nothing,
    link = nothing,
    filter_fun = nothing,
    ignore_nothing = false
)

    check_filters(mtg, scale = scale, symbol = symbol, link = link)

    keep_var = []

    for nc in args

        if nc isa Symbol || nc isa String
            # `new name` form, i.e. just selecting the variable
            push!(keep_var, nc)
        elseif last(nc) isa Symbol
            # `function => new name` form, i.e. `node -> sum(descendants(node, :var)) => :newvar`.
            # or `Name => new name` form, i.e. :x => :y.
            push!(keep_var, last(nc))
        elseif last(nc) isa Pair && (last(last(nc)) isa Symbol || last(last(nc)) isa String)
            # `Name => function => new name` form, i.e. :x => sqrt => :x_sq
            push!(keep_var, last(last(nc)))
        elseif nc isa Base.Callable
            # If we provide just a function, it is just a transformation to apply on the node
            # directly. No need to handle this here (it is handled by transform!)
            continue
        else
            # `Name => function` form, i.e. :x => sqrt
            # ?NOTE: Here the function takes one or more attributes as input
            col_idx, fun = nc
            newname, col_idx = col_name_from_call(col_idx, fun)
            push!(keep_var, newname)
        end
    end

    # Transform the MTG if any transformation is needed:
    transform!(
        mtg,
        args...;
        scale = scale,
        symbol = symbol,
        link = link,
        filter_fun = filter_fun,
        ignore_nothing = ignore_nothing
    )

    # Remove all un-selected attributes from the MTG:
    traverse!(
        mtg,
        node -> [pop!(node, attr) for attr in setdiff(keys(node.attributes), keep_var)]
    )

end

function select(
    mtg::Node,
    args...;
    scale = nothing,
    symbol = nothing,
    link = nothing,
    filter_fun = nothing,
    ignore_nothing = false
)
    new_mtg = deepcopy(mtg)

    select!(
        new_mtg,
        args...;
        scale = scale,
        symbol = symbol,
        link = link,
        filter_fun = filter_fun,
        ignore_nothing = ignore_nothing
    )

    return new_mtg
end


text:
module SpringCollab2020CustomSnow

using ..Ahorn, Maple

@mapdef Effect "SpringCollab2020/CustomSnow" CustomSnow(only::String="*", exclude::String="", colors::String="FFFFFF,FFFFFF")

placements = CustomSnow

function Ahorn.canFgBg(effect::CustomSnow)
    return true, true
end

end


text:
module Patches

using AWS
using Dates
using HTTP
using JSON
using GitHub
using Mocking
using OrderedCollections: LittleDict

version = v"1.1.0"
status = 200
headers = Pair[
    "x-amz-id-2" => "x-amz-id-2",
    "x-amz-request-id" => "x-amz-request-id",
    "Date" => "Tue, 16 Jun 2020 21:29:18 GMT",
    "x-amz-bucket-region" => "us-east-1",
    "Content-Type" => "application/xml",
    "Transfer-Encoding" => "chunked",
    "Server" => "AmazonS3"
]

body = """
    <?xml version=\"1.0\" encoding=\"UTF-8\"?>\n
    <ListBucketResult xmlns=\"http://s3.amazonaws.com/doc/2006-03-01/\">
        <Name>sample-bucket</Name>
        <Prefix></Prefix>
        <Marker></Marker>
        <MaxKeys>1000</MaxKeys>
        <IsTruncated>false</IsTruncated>
        <Contents>
            <Key>test.txt</Key>
            <LastModified>2020-06-16T21:37:34.000Z</LastModified>
            <ETag>"d41d8cd98f00b204e9800998ecf8427e"</ETag>
            <Size>0</Size>
            <Owner>
                <ID>id</ID>
                <DisplayName>matt.brzezinski</DisplayName>
            </Owner>
            <StorageClass>STANDARD</StorageClass>
        </Contents>
    </ListBucketResult>
    """

function _response(; version::VersionNumber=version, status::Int64=status, headers::Array=headers, body::String=body)
    response = HTTP.Messages.Response()

    response.version = version
    response.status = status
    response.headers = headers
    response.body = Vector{UInt8}(body)

    return response
end

function _aws_http_request_patch(response::HTTP.Messages.Response=_response())
    return @patch AWS._http_request(::AWS.AbstractBackend, request::Request) = response
end

_cred_file_patch = @patch function dot_aws_credentials_file()
    return ""
end

_config_file_patch = @patch function dot_aws_config_file()
    return ""
end

_assume_role_patch = function (
    op;
    access_key="access_key",
    secret_key="secret_key",
    session_token="token",
    role_arn="arn:aws:sts:::assumed-role/role-name",
)
    @patch function AWSServices.sts(op, params; aws_config)
        return Dict(
            "$(op)Result" => Dict(
                "Credentials" => Dict(
                    "AccessKeyId" => access_key,
                    "SecretAccessKey" => secret_key,
                    "SessionToken" => session_token,
                    "Expiration" => string(now(UTC)),
                ),
                "AssumedRoleUser" => Dict(
                    "Arn" => "$(role_arn)/$(params["RoleSessionName"])",
                ),
            ),
        )
    end
end

_github_tree_patch = @patch function tree(repo, tree_obj; kwargs...)
    if tree_obj == "master"
        return Tree("test-sha", HTTP.URI(), [Dict("path"=>"apis", "sha"=>"apis-sha")], false)
    else
        return Tree("test-sha", HTTP.URI(), [Dict("path"=>"test-2020-01-01.normal.json")], false)
    end
end

_instance_metadata_timeout_patch = @patch function HTTP.request(method::String, url; kwargs...)
    throw(HTTP.ConnectionPool.ConnectTimeout("169.254.169.254", "80"))
end

# This patch causes `HTTP.request` to return all of its keyword arguments
# except `require_ssl_verification` and `response_stream`. This is used to
# test which other options are being passed to `HTTP.Request` inside of
# `_http_request`.
_http_options_patch = @patch function HTTP.request(args...; kwargs...)
    options = Dict(kwargs)
    delete!(options, :require_ssl_verification)
    delete!(options, :response_stream)
    return options
end

end


text:
using StatsBase: sample
using EvoTrees: sigmoid, logit
using MLJBase

##################################################
### Regression - small data
##################################################
features = rand(10_000) .* 5 .- 2
X = reshape(features, (size(features)[1], 1))
Y = sin.(features) .* 0.5 .+ 0.5
Y = logit(Y) + randn(size(Y))
Y = sigmoid(Y)
y = Y
X = MLJBase.table(X)

# @load EvoTreeRegressor
# linear regression
tree_model = EvoTreeRegressor(max_depth=5, η=0.05, nrounds=10)
# logistic regression
tree_model = EvoTreeRegressor(loss=:logistic, max_depth=5, η=0.05, nrounds=10)
# quantile regression
# tree_model = EvoTreeRegressor(loss=:quantile, α=0.75, max_depth=5, η=0.05, nrounds=10)

mach = machine(tree_model, X, y)
train, test = partition(eachindex(y), 0.7, shuffle=true); # 70:30 split
fit!(mach, rows=train, verbosity=1)

mach.model.nrounds += 10
fit!(mach, rows=train, verbosity=1)

# predict on train data
pred_train = predict(mach, selectrows(X,train))
mean(abs.(pred_train - selectrows(Y,train)))

# predict on test data
pred_test = predict(mach, selectrows(X,test))
mean(abs.(pred_test - selectrows(Y,test)))

@test MLJBase.iteration_parameter(EvoTreeRegressor) == :nrounds

##################################################
### classif - categorical target
##################################################
X, y = @load_crabs

tree_model = EvoTreeClassifier(max_depth=4, η=0.05, λ=0.0, γ=0.0, nrounds=10)

# @load EvoTreeRegressor
mach = machine(tree_model, X, y)
train, test = partition(eachindex(y), 0.7, shuffle=true); # 70:30 split
fit!(mach, rows=train, verbosity=1)

mach.model.nrounds += 50
fit!(mach, rows=train, verbosity=1)

pred_train = predict(mach, selectrows(X,train))
pred_train_mode = predict_mode(mach, selectrows(X,train))
cross_entropy(pred_train, selectrows(y, train)) |> mean
sum(pred_train_mode .== y[train]) / length(y[train])

pred_test = predict(mach, selectrows(X,test))
pred_test_mode = predict_mode(mach, selectrows(X,test))
cross_entropy(pred_test, selectrows(y, test)) |> mean
sum(pred_test_mode .== y[test]) / length(y[test])
pred_test_mode = predict_mode(mach, selectrows(X,test))

##################################################
### count
##################################################
features = rand(10_000, 10)
# features = rand(100, 10)
X = features
Y = rand(UInt8, size(X, 1))
𝑖 = collect(1:size(X,1))

# train-eval split
𝑖_sample = sample(𝑖, size(𝑖, 1), replace = false)
train_size = 0.8
𝑖_train = 𝑖_sample[1:floor(Int, train_size * size(𝑖, 1))]
𝑖_eval = 𝑖_sample[floor(Int, train_size * size(𝑖, 1)) + 1:end]

X_train, X_eval = X[𝑖_train, :], X[𝑖_eval, :]
Y_train, Y_eval = Y[𝑖_train], Y[𝑖_eval]

# @load EvoTreeRegressor
tree_model = EvoTreeCount(
    loss=:poisson, metric=:poisson,
    nrounds=10,
    λ = 0.0, γ=0.0, η=0.1,
    max_depth = 6, min_weight = 1.0,
    rowsample=0.5, colsample=0.5, nbins=32)

X = MLJBase.table(X)
X = MLJBase.matrix(X)

# typeof(X)
mach = machine(tree_model, X, Y)
train, test = partition(eachindex(Y), 0.8, shuffle=true); # 70:30 split
fit!(mach, rows=train, verbosity=1, force=true)

mach.model.nrounds += 10
fit!(mach, rows=train, verbosity=1)

pred = predict(mach, selectrows(X,train))
pred_mean = predict_mean(mach, selectrows(X,train))
pred_mode = predict_mode(mach, selectrows(X,train))
# pred_mode = predict_median(mach, selectrows(X,train))

##################################################
### Gaussian - Larger data
##################################################
features = rand(10_000, 10)
X = features
Y = rand(size(X, 1))
𝑖 = collect(1:size(X,1))

# train-eval split
𝑖_sample = sample(𝑖, size(𝑖, 1), replace = false)
train_size = 0.8
𝑖_train = 𝑖_sample[1:floor(Int, train_size * size(𝑖, 1))]
𝑖_eval = 𝑖_sample[floor(Int, train_size * size(𝑖, 1)) + 1:end]

X_train, X_eval = X[𝑖_train, :], X[𝑖_eval, :]
Y_train, Y_eval = Y[𝑖_train], Y[𝑖_eval]

# @load EvoTreeRegressor
tree_model = EvoTreeGaussian(
    loss=:gaussian, metric=:gaussian,
    nrounds=10,
    λ = 0.0, γ=0.0, η=0.1,
    max_depth = 6, min_weight = 1.0,
    rowsample=0.5, colsample=0.5, nbins=32)

X = MLJBase.table(X)

# typeof(X)
mach = machine(tree_model, X, Y)
train, test = partition(eachindex(Y), 0.8, shuffle=true); # 70:30 split
fit!(mach, rows=train, verbosity=1, force=true)

mach.model.nrounds += 10
fit!(mach, rows=train, verbosity=1)

pred = predict(mach, selectrows(X,train))
pred_mean = predict_mean(mach, selectrows(X,train))
pred_mode = predict_mode(mach, selectrows(X,train))
# pred_mode = predict_median(mach, selectrows(X,train))
mean(abs.(pred_mean - selectrows(Y,train)))

q_20 = quantile.(pred, 0.20)
q_20 = quantile.(pred, 0.80)

report(mach)

############################
# Added in response to #92 #
############################

# tests that `update` handles data correctly in the case of a cold
# restatrt:

X = MLJBase.table(rand(5,2))
y = rand(5)
model = EvoTreeRegressor()
data = MLJBase.reformat(model, X, y);
f, c, r = MLJBase.fit(model, 2, data...);
model.λ = 0.1
MLJBase.update(model, 2, f, c, data...);


text:
import GR

import Random
srand(seed) = Random.seed!(seed)

function main()
    srand(37)
    y = randn(20, 500)

    GR.setviewport(0.1, 0.95, 0.1, 0.95)
    GR.setcharheight(0.020)
    GR.settextcolorind(82)
    GR.setfillcolorind(90)
    GR.setfillintstyle(1)

    for x in 1:5000
        GR.clearws()
        GR.setwindow(x, x+500, -200, 200)
        GR.fillrect(x, x+500, -200, 200)
        GR.setlinecolorind(0);  GR.grid(50, 50, 0, -200, 2, 2)
        GR.setlinecolorind(82); GR.axes(50, 50, x, -200, 2, 2, -0.005)
        y = hcat(y, randn(20))
        for i in 1:20
            GR.setlinecolorind(980 + i)
            s = cumsum(reshape(y[i,:], x+500))
            GR.polyline([x:x+500;], s[x:x+500])
        end
        GR.updatews()
    end
end

main()


text:
@model function simulate_data(::Type{T} = Float64;
    DOdata = missing,
    ODE,
    times,
    nDays,
    GPP_daily_mu,
    GPP_daily_sigma,
    ER_daily_mu,
    ER_daily_sigma,
    k600_meanlog,
    k600_sdlog,
    K02_conv,
    err_obs_iid_sigma_scale,
    err_proc_iid_sigma_scale
    ) where {T}

    # Priors on daily GPP and ER, k600, observation error, and process error
    GPP_daily ~ MvNormal(FillArrays.Fill(GPP_daily_mu, nDays), GPP_daily_sigma)
    ER_daily ~ MvNormal(FillArrays.Fill(ER_daily_mu, nDays), ER_daily_sigma)
    k600 ~ LogNormal(k600_meanlog, k600_sdlog)
    kGas = K02_conv*k600
    err_obs_iid_sigma ~ truncated(Cauchy(0, err_obs_iid_sigma_scale), 0, Inf)
    err_proc_iid_sigma ~ truncated(Cauchy(0, err_proc_iid_sigma_scale), 0, Inf)

    p = Dict(
        "PAR" => ODE.p["PAR"],
        "OSat" => ODE.p["OSat"],
        "day" => ODE.p["day"],
        "GPP_daily"=> GPP_daily,
        "ER_daily" => ER_daily,
        "kGas" => kGas,
        "Q" => ODE.p["Q"],
        "V" => ODE.p["V"])

    if DOdata === missing
        DOdata = Vector{T}(undef, length(times))
        DOdata[1] = ODE.u0[1]
    end

    DO_true = Vector{T}(undef, (length(times)))
    DO_true[1] ~ Normal(DOdata[1], err_obs_iid_sigma)

    for i = 2:length(times)
        problem = remake(ODE, p = p, u0 = [DO_true[i-1]], tspan = (times[i-1], times[i]))
        predicted = solve(problem, alg = Trapezoid())
        DO_true[i] ~ Normal(predicted(times[i])[1], err_proc_iid_sigma)
        DOdata[i] ~ Normal(DO_true[i], err_obs_iid_sigma)
    end

    return DOdata
end

max_stars_repo_name: elizabethjmohr/mesocosm-metabolism

max_issues_repo_name: elizabethjmohr/mesocosm-metabolism

max_forks_repo_name: elizabethjmohr/mesocosm-metabolism


text:
# This file contains functions to compute posterior inference results

using StatsBase

# posterior quantile intervals

doc"""
    calcquantileint(data::Matrix, lower::Real=2.5,upper::Real=97.5)

Computes posterior quantile intervals, 2.5th and 97.5th quantiles as default.
"""
function calcquantileint(data::Matrix, lower::Real=2.5,upper::Real=97.5)

    # find dim. for data
    dim = minimum(size(data))

    # set nbr of intervals
    intervals = zeros(dim, 2)

    # transform data to column-major order if necessary
    if size(data)[1] > size(data)[2]
        data = data'
    end

    # calc intervals over all dimensions
    for i = 1:dim
        intervals[i,:] = quantile(data[i,:], [lower/100 upper/100])
    end

    # return intervals
    return intervals

end

doc"""
    calcquantileint(data::Matrix, lower::Real=2.5,upper::Real=97.5)

Computes posterior quantile intervals, 2.5th and 97.5th quantiles as default.
"""
function calcquantileint(data::Vector, lower::Real=2.5,upper::Real=97.5)

    # set nbr of intervals
    intervals = zeros(1, 2)

    # calc intervals over all dimensions
    intervals[1,:] = quantile(data[1,:], [lower/100 upper/100])

    # return intervals
    return intervals

end

# RMSE

doc"""
    RMSE(theta_true::Vector, theta_est::Matrix)

Computes the root mean square error for the parameter estimations.
"""
function RMSE(theta_true::Vector, Theta_est::Matrix)
  rmse = copy(theta_true)
  for i = 1:length(rmse)
    rmse[i] = sqrt(mean((theta_true[i]-Theta_est[i,:]).^2))
  end
  return rmse
end

max_stars_repo_name: SamuelWiqvist/abcmcmcrickermodel

max_issues_repo_name: SamuelWiqvist/abcmcmcrickermodel

max_forks_repo_name: SamuelWiqvist/abcmcmcrickermodel


text:
const AbstractArraySet = Union{CartesianProductArray,
                               ConvexHullArray,
                               IntersectionArray,
                               MinkowskiSumArray,
                               UnionSetArray
                              }

function Base.getindex(X::AbstractArraySet, i)
    return getindex(array(X), i)
end

function Base.length(X::AbstractArraySet)
    return length(array(X))
end


text:
using GeometricFlux
using GraphSignals
using Graphs
using SparseArrays
using Plots
using GraphPlot
using Clustering
using Cairo, Compose

clusters = [1, 1, 1, 1, 1, 1, 1, 1, 2, 2, 1, 1, 1, 1, 2, 2, 1, 1, 2, 1, 2, 1, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2]

int2col_str(x::Int) = x==1 ? "lightblue" : "red"

g = smallgraph(:karate)
fg = FeaturedGraph(g)
vectors = node2vec(fg; walks_per_node=10, len=80, p=1.0, q=1.0)
R = kmeans(vectors, 2)

learned_clusters = copy(assignments(R))
# ensure that the cluster containing node 1 is cluster 1
if assignments(R)[1] != 1
    learned_clusters = [i == 1 ? 2 : 1 for i in assignments(R)]
end

output_plot_name = "karateclub.pdf"
draw(
    PDF(output_plot_name, 16cm, 16cm),
    gplot(g,
        nodelabel=map(string, 1:34),
        nodefillc=[int2col_str(learned_clusters[i]) for i in 1:34],
        nodestrokec=["white" for _ in 1:34]
    )
)

incorrect = sum(learned_clusters .!= clusters)
println(incorrect, " incorrect cluster labelings")
println("Drawn graph to ", output_plot_name)