text:
using OutOfSampleBootstraps
using Base.Test

n = 100
R = 30
P = n - R
k = 4

y = randn(n)
x = [ones(n) randn(n, k-1)]

estimates = Array(Float64, k, P)
errors = Array(Float64, P)

## Verify that recursive_ols produces the right estimates.
b1 = recursive_ols(y, x, R)
@test_approx_eq(b1[:,1], x[1:R,:] \ y[1:R])
@test_approx_eq(b1[:,end], x[1:end-1,:] \ y[1:end-1])
## results after reordering too
shuffle = rand(1:n, n)
b1_shuffled = recursive_ols(y, x, R, shuffle)
@test_approx_eq(b1_shuffled[:,1], x[shuffle[1:R],:] \ y[shuffle[1:R]])
@test_approx_eq(b1_shuffled[:,end], x[shuffle[1:end-1],:] \ y[shuffle[1:end-1]])

## Make sure that different functions to get the recursive OLS
## estimates return the same values.
recursive_ols!(estimates, y, x)
@test_approx_eq(estimates, b1)
recursive_ols!(estimates, y, x, shuffle)
@test_approx_eq(estimates, b1_shuffled)
recursive_ols!(estimates, errors, y, x)
@test_approx_eq(estimates, b1)
recursive_ols!(estimates, errors, y, x, shuffle)
@test_approx_eq(estimates, b1_shuffled)

## Make sure basic bootstraps run
b = Array(Float64, 199) ## Initialize storage. Each call will
                        ## overwrite the previous results
recursive_bootstrap!(b, :naive, Block(:moving, 15), y, x, R)
recursive_bootstrap!(b, :cs07_ols, Block(:moving, 15), y, x, R)
recursive_bootstrap!(b, :mine_ols, Block(:moving, 15), y, x, R)


text:
function trinary_to_decimal(str)

end

ext: jl


text:
"""
# Game

Encapsulates the state of a Quoridor game.
"""
type Game <: AbstractGraph{Loc, Edge{Loc}}
    walls::Vector{Wall}
    p1::PlayerState
    p2::PlayerState
    active::Player
end

initstate(::Type{Val{Player1}}) = PlayerState(Loc(5,1))
initstate(::Type{Val{Player2}}) = PlayerState(Loc(5,9))

Game() = Game(Wall[],
              initstate(Val{Player1}),
              initstate(Val{Player2}),
              rand(Player))

"return a copy of `game`, nothing shared with the original"
copy(game::Game) = Game(copy(walls), p1, p2, turn)

active(game::Game) = game.active
inactive(game::Game) = other(active(game))

"return the current list of walls"
walls(game::Game) = game.walls

getindex(game::Game, player::Player) = ifelse(player == Player1,
                                              game.p1, game.p2)

occupied(g::Game, l::Loc) = loc(g[Player1]) == l || loc(g[Player2]) == l

# applying and un-applying moves to games

function apply!(game::Game, pm::PlayerMove)
    # TODO implement
end

function apply!(game::Game, wp::WallPlacement)
    # TODO implement
end

# note that these are not time-invariant: they assume the move being
# taken back is the most recent move made to bring the game to it's
# current state. could avoid this by keeping a stack of made moves on
# the game

function takeback!(game::Game, pm::PlayerMove)
    # TODO implement
end

function takeback!(game::Game, wp::WallPlacement)
    # TODO implement
end

apply(game::Game, move::Move) = apply!(copy(Game), move)
takeback(game::Game, move::Move) = takeback!(copy(Game), move)

"Test if two `Loc`s are neighbors in the context of a `Game`"
function neighbors(game::Game, l1::Loc, l2::Loc)
    !any_separates(walls(game), l1, l2)
end

# generate the set of all legal moves for a given game
function legalmoves(game::Game)
    ret = Move[] # TODO this is going to be slow, should eventually
                 # replace with a struct that distinguishes between
                 # wall placements and player moves and return that,
                 # or work up an encoding of both types of moves in a
                 # single bitstype and return a Vector of that
    append!(ret, playermoves(game))
    append!(ret, wallplacements(game))
    return ret
end

playermoves(game::Game) = out_edges(loc(game[active(game)]), game)

function wallplacements(game::Game)
    # note that some of these placements will be illegal because they
    # block a player's path to the endzone, let's plan to filter these
    # out later because i think doing so now would be prohibitively
    # expensive? and i dont feel like writing the graph generation
    # code yet

    # TODO this is slow due to allocating the new list of placements,
    # should come back and rewrite it as an iterator but i don't feel
    # like it right now
    placements = WallPlacement[]
    for possibility in ALL_WALLS
        placement = WallPlacement(possibility)
        if !any_blocks(walls(game), placement)
            push!(placements, placement)
        end
    end
    return placements
end


text:
function Base.show(io::IO, op_model::DynamicSystem)
    println(io, "System()")
end

function Base.show(io::IO, op_model::DynamicSimulation)
    println(io, "Simulation()")
end


text:
#=
This file is auto-generated. Do not edit.
=#
"""
    mutable struct TwoPartCost <: OperationalCost
        variable::VariableCost
        fixed::Float64
        forecasts::InfrastructureSystems.Forecasts
        internal::InfrastructureSystemsInternal
    end

Data Structure Operational Cost Data in two parts: fixed and variable cost.

# Arguments
- `variable::VariableCost`: variable cost
- `fixed::Float64`: fixed cost, validation range: (0, nothing), action if invalid: warn
- `forecasts::InfrastructureSystems.Forecasts`: internal forecast storage
- `internal::InfrastructureSystemsInternal`: power system internal reference, do not modify
"""
mutable struct TwoPartCost <: OperationalCost
    "variable cost"
    variable::VariableCost
    "fixed cost"
    fixed::Float64
    "internal forecast storage"
    forecasts::InfrastructureSystems.Forecasts
    "power system internal reference, do not modify"
    internal::InfrastructureSystemsInternal
end

function TwoPartCost(variable, fixed, forecasts=InfrastructureSystems.Forecasts(), )
    TwoPartCost(variable, fixed, forecasts, InfrastructureSystemsInternal(), )
end

function TwoPartCost(; variable, fixed, forecasts=InfrastructureSystems.Forecasts(), )
    TwoPartCost(variable, fixed, forecasts, )
end

# Constructor for demo purposes; non-functional.
function TwoPartCost(::Nothing)
    TwoPartCost(;
        variable=VariableCost((0.0, 0.0)),
        fixed=0.0,
        forecasts=InfrastructureSystems.Forecasts(),
    )
end

"""Get TwoPartCost variable."""
get_variable(value::TwoPartCost) = value.variable
"""Get TwoPartCost fixed."""
get_fixed(value::TwoPartCost) = value.fixed

InfrastructureSystems.get_forecasts(value::TwoPartCost) = value.forecasts
"""Get TwoPartCost internal."""
get_internal(value::TwoPartCost) = value.internal

"""Set TwoPartCost variable."""
set_variable!(value::TwoPartCost, val::VariableCost) = value.variable = val
"""Set TwoPartCost fixed."""
set_fixed!(value::TwoPartCost, val::Float64) = value.fixed = val

InfrastructureSystems.set_forecasts!(value::TwoPartCost, val::InfrastructureSystems.Forecasts) = value.forecasts = val
"""Set TwoPartCost internal."""
set_internal!(value::TwoPartCost, val::InfrastructureSystemsInternal) = value.internal = val


text:
struct SE3PDController{T, W, G<:SE3PDGains}
    base::BodyID
    body::BodyID
    trajectory::Base.RefValue{T}
    weight::W
    gains::Base.RefValue{G}
end

function SE3PDController(base::BodyID, body::BodyID, trajectory, weight, gains::SE3PDGains)
    SE3PDController(base, body, Ref(trajectory), weight, Ref(gains))
end

function (controller::SE3PDController)(t, state::MechanismState)
    Href, Tref, Tdref = controller.trajectory[](t, Val(2))
    H = relative_transform(state, Href.from, Href.to)
    T = transform(relative_twist(state, controller.body, controller.base), inv(transform_to_root(state, Tref.frame)))
    return Tdref + pd(controller.gains[], H, Href, T, Tref)
end


text:
module InternetCafe

using ..Ahorn, Maple

const placements = Ahorn.PlacementDict(
    "Internet Café" => Ahorn.EntityPlacement(
        Maple.InternetCafe
    )
)

spriteBack = "objects/wavedashtutorial/building_back"
spriteLeft = "objects/wavedashtutorial/building_front_left"
spriteRight = "objects/wavedashtutorial/building_front_right"

# We don't need the back texture, it fits inside
function Ahorn.selection(entity::Maple.InternetCafe)
    x, y = Ahorn.position(entity)

    return Ahorn.coverRectangles([
        Ahorn.getSpriteRectangle(spriteLeft, x, y, jx=0.5, jy=1.0),
        Ahorn.getSpriteRectangle(spriteRight, x, y, jx=0.5, jy=1.0),
    ])
end

function Ahorn.render(ctx::Ahorn.Cairo.CairoContext, entity::Maple.InternetCafe, room::Maple.Room)
    Ahorn.drawSprite(ctx, spriteBack, 0, 0, jx=0.5, jy=1.0)
    Ahorn.drawSprite(ctx, spriteLeft, 0, 0, jx=0.5, jy=1.0)
    Ahorn.drawSprite(ctx, spriteRight, 0, 0, jx=0.5, jy=1.0)
end

end


text:
__precompile__()

module CPLEX

    if VERSION >= v"0.7.0-DEV.3382"
        using Libdl
    end

    if is_apple()
        Libdl.dlopen("libstdc++",Libdl.RTLD_GLOBAL)
    end

    if isfile(joinpath(dirname(@__FILE__),"..","deps","deps.jl"))
        include("../deps/deps.jl")
    else
        error("CPLEX not properly installed. Please run Pkg.build(\"CPLEX\")")
    end

    ### imports
    import Base: convert, unsafe_convert, show, copy

    # Standard LP interface
    importall MathProgBase.SolverInterface

    # exported functions
    # export is_valid,
    #        set_logfile,
    #        get_error_msg,
    #        read_model,
    #        write_model,
    #        get_sense,
    #        set_sense!,
    #        get_obj,
    #        set_obj!,
    #        set_warm_start!,
    #        free_problem,
    #        close_CPLEX,
    #        get_param_type,
    #        set_param!,
    #        get_param,
    #        add_vars!,
    #        add_var!,
    #        get_varLB,
    #        set_varLB!,
    #        get_varUB,
    #        set_varUB!,
    #        set_vartype!,
    #        get_vartype,
    #        num_var,
    #        add_constrs!,
    #        add_constrs_t!,
    #        add_rangeconstrs!,
    #        add_rangeconstrs_t!,
    #        num_constr,
    #        get_constr_senses,
    #        set_constr_senses!,
    #        get_rhs,
    #        set_rhs!,
    #        get_constrLB,
    #        get_constrUB,
    #        set_constrLB!,
    #        set_constrUB!,
    #        get_nnz,
    #        get_constr_matrix,
    #        set_sos!,
    #        add_qpterms!,
    #        add_diag_qpterms!,
    #        add_qconstr!,
    #        optimize!,
    #        get_objval,
    #        get_solution,
    #        get_reduced_costs,
    #        get_constr_duals,
    #        get_constr_solution,
    #        get_infeasibility_ray,
    #        get_unbounded_ray,
    #        get_status,
    #        get_status_code,
    #        setcallbackcut,
    #        cbcut,
    #        cblazy,
    #        cbget_mipnode_rel,
    #        cbget_mipsol_sol,
    #        cplex_model

    using Compat

    include("cpx_common.jl")
    include("cpx_env.jl")
    v = version()
    if startswith(v,"12.6")
        include("full_defines_126.jl")
        include("cpx_params_126.jl")
    elseif startswith(v,"12.7.1")
        include("full_defines_1271.jl")
        include("cpx_params_1271.jl")
    elseif startswith(v,"12.7")
        include("full_defines_127.jl")
        include("cpx_params_127.jl")
    elseif startswith(v,"12.8")
        include("full_defines_1280.jl")
        include("cpx_params_1280.jl")
    else
        error("Unsupported CPLEX version $v. Only 12.6, 12.7 and 12.8 are currently supported.")
    end
    include("cpx_model.jl")
    include("cpx_params.jl")
    include("cpx_vars.jl")
    include("cpx_constrs.jl")
    include("cpx_quad.jl")
    include("cpx_solve.jl")
    include("cpx_callbacks.jl")
    include("cpx_highlevel.jl")

    include("CplexSolverInterface.jl")
    include("MOIWrapper.jl")

    # These are undocumented JuMP extensions for CPLEX which
    # will need to be hosted in a separate package for Julia 0.6 and later.
    if isdir(Pkg.dir("JuMP")) && VERSION < v"0.6-"
        try
            eval(current_module(), Expr(:import,:JuMP))
            include("JuMPfunctions.jl")
        end
    end
end


text:
if VERSION >= v"1.2.0-DEV.573"
    include("snoopi.jl")
end

using SnoopCompile
using JLD
using SparseArrays
using Test

# issue #26
logfile = joinpath(tempdir(), "anon.log")
@snoopc logfile begin
    map(x->x^2, [1,2,3])
end
data = SnoopCompile.read(logfile)
pc = SnoopCompile.parcel(reverse!(data[2]))
@test length(pc[:Base]) <= 1

# issue #29
keep, pcstring, topmod, name = SnoopCompile.parse_call("Tuple{getfield(JLD, Symbol(\"##s27#8\")), Any, Any, Any, Any, Any}")
@test keep
@test pcstring == "Tuple{getfield(JLD, Symbol(\"##s27#8\")), Int, Int, Int, Int, Int}"
@test topmod == :JLD
@test name == "##s27#8"
matfile = joinpath(tempdir(), "mat.jld")
save(matfile, "mat", sprand(10, 10, 0.1))
logfile = joinpath(tempdir(), "jldanon.log")
@snoopc logfile begin
    using JLD, SparseArrays
    mat = load(joinpath(tempdir(), "mat.jld"), "mat")
end
data = SnoopCompile.read(logfile)
pc = SnoopCompile.parcel(reverse!(data[2]))
@test any(startswith.(pc[:JLD], "isdefined"))

#=
# Simple call
let str = "sum"
    keep, pcstring, topmod = SnoopCompile.parse_call("Foo.any($str)")
    @test keep
    @test pcstring == "Tuple{$str}"
    @test topmod == :Main
end

# Operator
let str = "Base.:*, Int, Int"
    keep, pcstring, topmod = SnoopCompile.parse_call("Foo.any($str)")
    @test keep
    @test pcstring == "Tuple{$str}"
    @test topmod == :Base
end

# Function as argument
let str = "typeof(Base.identity), Array{Bool, 1}"
    keep, pcstring, topmod = SnoopCompile.parse_call("Foo.any($str, Vararg{Any, N} where N)")
    @test keep
    @test pcstring == "Tuple{$str, Int}"
    @test topmod == :Base
end

# Anonymous function closure in a new module as argument
let func = (@eval Main module SnoopTestTemp
            func = () -> (y = 2; (x -> x > y))
        end).func
    str = "getfield(SnoopTestTemp, Symbol(\"$(typeof(func()))\")), Array{Float32, 1}"
    keep, pcstring, topmod = SnoopCompile.parse_call("Foo.any($str)")
    @test keep
    @test pcstring == "Tuple{$str}"
    @test topmod == :SnoopTestTemp
end

# Function as a type
let str = "typeof(Base.Sort.sort!), Array{Any, 1}, Base.Sort.MergeSortAlg, Base.Order.By{typeof(Base.string)}"
    keep, pcstring, topmod = SnoopCompile.parse_call("Foo.Bar.sort!($str)")
    @test keep
    @test pcstring == "Tuple{$str}"
    @test topmod == :Base
end
=#

@static if VERSION >= v"1.2.0-DEV.573"
    @testset "timesum" begin
        loadSnoop = SnoopCompile.@snoopi using MatLang
        @test typeof(timesum(loadSnoop)) == Float64
    end
end

include("colortypes.jl")
include("bot/bot.jl")


text:
export update!

function update!(n::Node, grads::NamedTuple; η = 0.01, uparams = [])
    if !isnothing(grads.params)
        n.params .+= η * grads.params
    end

    if grads.children === nothing
        return
    end

    for (i, child) in enumerate(n.children)
        if (child isa AbstractLeaf)
            n.children[i] = update!(child, grads.children[i]; η = η, uparams = uparams)
        else
            if grads.children[i] === nothing # why is this happening?
                nothing
            else
                update!(child, grads.children[i]; η = η, uparams = uparams)
            end
        end
    end
end

function update!(n::T, grads::NamedTuple; η = 0.01, uparams = []) where {T<:AbstractLeaf}

    v = collect(params(n))
    for (i,k) in enumerate(keys(params(n)))
        if (k ∈ uparams) || isempty(uparams)
            v[i] += grads.params[k] === nothing ? zero(eltype(v)) : η * grads.params[k]
        end
    end

    return T(n.scope, (;zip(keys(params(n)), v)...))
end

function update!(::Indicator, ::NamedTuple; η = 0.01, uparams = [])
    nothing
end

function update!(n::VISumNode, grads::NamedTuple; η = 0.01, uparams = [])
    print("VIupdate begin")
    if !isnothing(grads.params)
        α = 1
        grads_KL = _Dir_KL_grads(n.params,α)
        n.params .+= η * (grads.params .- grads_KL)
    end

    if grads.children === nothing
        return
    end

    for (i, child) in enumerate(n.children)
        if child isa AbstractLeaf
            n.children[i] = update!(child, grads.children[i]; η = η, uparams = uparams)
        else
            if grads.children[i] === nothing # why is this happening?
                nothing
            else
                update!(child, grads.children[i]; η = η, uparams = uparams)
            end
        end
    end
end

function _Dir_KL_grads(logβ,α) #KL objetive between Dir dist: bariskurt.com/kullback-leibler-divergence-between-two-dirichlet-and-beta-distributions/
    grads = Array{Float64, 1}(UndefInitializer(), length(logβ))
    β0 = exp(logsumexp(logβ))
    β = exp.(logβ)
    for j in 1:length(logβ)
        tmp1 = sum([(α[i]-β[i])*trigamma(β0)*β[j] for i in 1:length(logβ) if i != j])
        tmp2 = β[j]*(digamma(β[j])-digamma(β0)) + (β[j]-α[j])*(trigamma(β[j])-trigamma(β0))*β[j]
        grads[j] = (digamma(β0)-digamma(β[j])) * β[j] + tmp1 + tmp2
    end
    return grads
end


text:
using Test
using CLIMA.PlanetParameters

@testset begin
  # Not sure how to really test this...
  @test @isdefined R_d
  @test LH_f0 == LH_s0 - LH_v0
end

max_stars_repo_name: christophernhill-dev/CLIMA

max_issues_repo_name: christophernhill-dev/CLIMA

max_forks_repo_name: christophernhill-dev/CLIMA


text:
module AbstractInterpreterPlayground

# overloads
# ---------

import Core.Compiler:
    # abstractinterpreterinterface.jl
    InferenceParams,
    OptimizationParams,
    get_world_counter,
    get_inference_cache,
    lock_mi_inference,
    unlock_mi_inference,
    add_remark!,
    may_optimize,
    may_compress,
    may_discard_trees,
    # cicache.jl
    code_cache,
    WorldView,
    # inferenceresult.jl
    cache_lookup,
    # tfuncs.jl
    invoke_tfunc,
    builtin_tfunction,
    return_type_tfunc,
    # abstractinterpretation.jl
    abstract_call_gf_by_type,
    abstract_call_method_with_const_args,
    abstract_call_method,
    precise_container_type,
    abstract_iteration,
    abstract_apply,
    abstract_call_builtin,
    abstract_call_known,
    abstract_call,
    abstract_eval_cfunction,
    abstract_eval_value_expr,
    abstract_eval_special_value,
    abstract_eval_value,
    abstract_eval_statement,
    typeinf_local,
    typeinf_nocycle,
    # typeinfer.jl
    typeinf,
    _typeinf,
    transform_result_for_cache,
    is_same_frame,
    typeinf_edge

# imports
# -------

import Core:
    MethodMatch,
    LineInfoNode,
    SimpleVector,
    Builtin,
    Typeof,
    svec

import Core.Compiler:
    AbstractInterpreter,
    NativeInterpreter,
    InferenceState,
    InferenceResult,
    MethodInstance,
    WorldRange,
    CodeInfo,
    VarTable,
    CodeInstance,
    to_tuple_type,
    _methods_by_ftype,
    specialize_method

import Base:
    unwrap_unionall,
    rewrap_unionall

import Base.Meta:
    isexpr,
    lower

using InteractiveUtils

const CC = Core.Compiler

# utilties
# --------

"""
    @invoke f(arg::T, ...; kwargs...)

Provides a convenient way to call [`invoke`](@ref);
`@invoke f(arg1::T1, arg2::T2; kwargs...)` will be expanded into `invoke(f, Tuple{T1,T2}, arg1, arg2; kwargs...)`.
When an argument's type annotation is omitted, it's specified as `Any` argument, e.g.
`@invoke f(arg1::T, arg2)` will be expanded into `invoke(f, Tuple{T,Any}, arg1, arg2)`.

This could be used to call down to `NativeInterpreter`'s abstract interpretation method of
  `f` while passing `CustomInterpreter` so that subsequent calls of abstract interpretation
  functions overloaded against `CustomInterpreter` can be called from the native method of `f`;
e.g. calls down to `NativeInterpreter`'s `abstract_call_gf_by_type` method:
```julia
@invoke abstract_call_gf_by_type(interp::AbstractInterpreter, f, argtypes::Vector{Any}, atype, sv::InferenceState,
                                 max_methods::Int)
```
"""
macro invoke(ex)
    f, args, kwargs = destructure_callex(ex)
    arg2typs = map(args) do x
        isexpr(x, :(::)) ? (x.args...,) : (x, GlobalRef(Core, :Any))
    end
    args, argtypes = first.(arg2typs), last.(arg2typs)
    return if isempty(kwargs)
        :($(GlobalRef(Core, :invoke))($(f), Tuple{$(argtypes...)}, $(args...))) # might not be necessary
    else
        :($(GlobalRef(Core, :invoke))($(f), Tuple{$(argtypes...)}, $(args...); $(kwargs...)))
    end |> esc
end

function destructure_callex(ex)
    @assert isexpr(ex, :call) "call expression f(args...; kwargs...) should be given"

    f = first(ex.args)
    args = []
    kwargs = []
    for x in ex.args[2:end]
        if isexpr(x, :parameters)
            append!(kwargs, x.args)
        elseif isexpr(x, :kw)
            push!(kwargs, x)
        else
            push!(args, x)
        end
    end

    return f, args, kwargs
end

# for inspection
macro lwr(ex) QuoteNode(lower(__module__, ex)) end
macro src(ex) QuoteNode(first(lower(__module__, ex).args)) end

# includes
# --------

include("abstractinterpreterinterface.jl")
include("cicache.jl")
include("inferenceresult.jl")
include("tfuncs.jl")
include("abstractinterpretation.jl")
include("typeinfer.jl")
include("print.jl")

# entry
# -----

macro enter_call(ex0...)
    return InteractiveUtils.gen_call_with_extracted_types_and_kwargs(@__MODULE__, :enter_call, ex0)
end

function enter_call(@nospecialize(f), @nospecialize(types=Tuple{});
                    kwargs...)
    interp = CustomInterpreter(; inf_params = gen_inf_params(; kwargs...),
                                 opt_params = gen_opt_params(; kwargs...),
                                 kwargs...)
    ft = Typeof(f)
    tt = if isa(types, Type)
        u = unwrap_unionall(types)
        rewrap_unionall(Tuple{ft, u.parameters...}, types)
    else
        Tuple{ft, types...}
    end
    return enter_gf_by_type!(interp, tt)
end

# TODO `enter_call_builtin!` ?
function enter_gf_by_type!(interp::CustomInterpreter,
                           @nospecialize(tt::Type{<:Tuple}),
                           world::UInt = get_world_counter(interp),
                           )
    mms = _methods_by_ftype(tt, InferenceParams(interp).MAX_METHODS, world)
    @assert mms !== false "unable to find matching method for $(tt)"

    filter!(mm->mm.spec_types===tt, mms)
    @assert length(mms) == 1 "unable to find single target method for $(tt)"

    mm = first(mms)::MethodMatch

    return enter_method_signature!(interp, mm.method, mm.spec_types, mm.sparams)
end

function enter_method_signature!(interp::CustomInterpreter,
                                 m::Method,
                                 @nospecialize(atype),
                                 sparams::SimpleVector,
                                 world::UInt = get_world_counter(interp),
                                 )
    mi = specialize_method(m, atype, sparams)

    result = InferenceResult(mi)

    frame = InferenceState(result, #= cached =# true, interp)

    typeinf(interp, frame)

    return interp, frame
end

function enter_method!(interp::CustomInterpreter,
                       m::Method,
                       world::UInt = get_world_counter(interp),
                       )
    return enter_method_signature!(interp, m, m.sig, sparams_from_method_signature(m), world)
end

function sparams_from_method_signature(m)
    s = TypeVar[]
    sig = m.sig
    while isa(sig, UnionAll)
        push!(s, sig.var)
        sig = sig.body
    end
    return svec(s...)
end

# exports

export
    @enter_call,
    enter_call

end

max_stars_repo_name: aviatesk/AbstractInterpreterPlayground

max_issues_repo_name: aviatesk/AbstractInterpreterPlayground

max_forks_repo_name: aviatesk/AbstractInterpreterPlayground


text:
plot(v1,2*outvIm./(outvIm.^2+outvRe.^2),xaxis="p0",yaxis=:log)
using Dierckx

Spline1D(readdlm("$dir3/Imp0.dat")[:,1],readdlm("$dir3/ImGamma2.dat")[:,1])


text:
# Data
"""
    cdo_data(CDO_token::AbstractString, dataset::AbstractString, startdate::Date, enddate::Date;
             datatypes::Union{AbstractString, AbstractVector{<:AbstractString}} = "",
             locations::Union{AbstractString, AbstractVector{<:AbstractString}} = "",
             stations::Union{AbstractString, AbstractVector{<:AbstractString}} = "",
             metric::Bool = true)

For obtaining a CDO_token: [Request Web Services Token](https://www.ncdc.noaa.gov/cdo-web/token)

For additional information visit the [NCDC's Climate Data Online (CDO) Web Services v2 Documentation](https://www.ncdc.noaa.gov/cdo-web/webservices/v2#data)
"""
function cdo_data(CDO_token::AbstractString, dataset::AbstractString,
                  startdate::TimeType, enddate::TimeType;
                  datatypes::Union{AbstractString, AbstractVector{<:AbstractString}} = "",
                  locations::Union{AbstractString, AbstractVector{<:AbstractString}} = "",
                  stations::Union{AbstractString, AbstractVector{<:AbstractString}} = "",
                  metric::Bool = true)
    isvalid_cdotoken(CDO_token) || throw(CDO_NonValidToken)
    isempty(datatypes) || isvalid_datatypes(datatypes) || throw(CDO_NonValidDataTypes)
    isempty(locations) || isvalid_locations(locations) || throw(CDO_NonValidLocations)
    isempty(stations) || isvalid_stations(stations) || throw(CDO_NonValidStations)
    return parse(CDO_Data(CDO_token, dataset, startdate, enddate, datatypes, locations, stations, metric))
end

# Data Categories
"""
    cdo_datacategories(CDO_token::AbstractString, datacategory::AbstractString)
    cdo_datacategories(CDO_token::AbstractString;
                       datasets::Union{AbstractString, AbstractVector{<:AbstractString}} = "",
                       locations::Union{AbstractString, AbstractVector{<:AbstractString}} = "",
                       stations::Union{AbstractString, AbstractVector{<:AbstractString}} = "",
                       startdate::Date = Date("0001-01-01"),
                       enddate::Date = today())

For obtaining a CDO_token: [Request Web Services Token](https://www.ncdc.noaa.gov/cdo-web/token)

For additional information visit the [NCDC's Climate Data Online (CDO) Web Services v2 Documentation](https://www.ncdc.noaa.gov/cdo-web/webservices/v2#dataCategories)
"""
function cdo_datacategories(CDO_token::AbstractString, datacategory::AbstractString)
    isvalid_cdotoken(CDO_token) || throw(CDO_NonValidToken)
    output = parse(CDO_DataCategory(CDO_token, datacategory))
    isa(output, Exception) && throw(output)
    return output
end
function cdo_datacategories(CDO_token::AbstractString;
                            datasets::Union{AbstractString, AbstractVector{<:AbstractString}} = "",
                            locations::Union{AbstractString, AbstractVector{<:AbstractString}} = "",
                            stations::Union{AbstractString, AbstractVector{<:AbstractString}} = "",
                            startdate::Date = Date("0001-01-01"),
                            enddate::Date = today())
    isvalid_cdotoken(CDO_token) || throw(CDO_NonValidToken)
    isempty(datasets) || isvalid_datatypes(datasets) || throw(CDO_NonValidDatasets)
    isempty(locations) || isvalid_locations(locations) || throw(CDO_NonValidLocations)
    isempty(stations) || isvalid_stations(stations) || throw(CDO_NonValidStations)
    return parse(CDO_DataCategories(CDO_token, datasets, locations, stations, startdate, enddate))
end

# Datasets
"""
    cdo_datasets(CDO_token::AbstractString, dataset::AbstractString)
    cdo_datasets(CDO_token::AbstractString;
                datatypes::Union{AbstractString, AbstractVector{<:AbstractString}} = "",
                locations::Union{AbstractString, AbstractVector{<:AbstractString}} = "",
                stations::Union{AbstractString, AbstractVector{<:AbstractString}} = "",
                startdate::Date = Date("0001-01-01"),
                enddate::Date = today())

For obtaining a CDO_token: [Request Web Services Token](https://www.ncdc.noaa.gov/cdo-web/token)

For additional information visit the [NCDC's Climate Data Online (CDO) Web Services v2 Documentation](https://www.ncdc.noaa.gov/cdo-web/webservices/v2#datasets)
"""
function cdo_datasets(CDO_token::AbstractString, dataset::AbstractString)
    isvalid_cdotoken(CDO_token) || throw(CDO_NonValidToken)
    output = parse(CDO_Dataset(CDO_token, dataset))
    isa(output, Exception) && throw(output)
    return output
end
function cdo_datasets(CDO_token::AbstractString;
                      datatypes::Union{AbstractString, AbstractVector{<:AbstractString}} = "",
                      locations::Union{AbstractString, AbstractVector{<:AbstractString}} = "",
                      stations::Union{AbstractString, AbstractVector{<:AbstractString}} = "",
                      startdate::Date = Date("0001-01-01"),
                      enddate::Date = today())
    isvalid_cdotoken(CDO_token) || throw(CDO_NonValidToken)
    isempty(datatypes) || isvalid_datatypes(datatypes) || throw(CDO_NonValidDataTypes)
    isempty(locations) || isvalid_locations(locations) || throw(CDO_NonValidLocations)
    isempty(stations) || isvalid_stations(stations) || throw(CDO_NonValidStations)
    return parse(CDO_Datasets(CDO_token, datatypes, locations, stations, startdate, enddate))
end

# Data Types
"""
    cdo_datatypes(CDO_token::AbstractString, datatype::AbstractString)
    cdo_datatypes(CDO_token::AbstractString;
                  datasets::Union{AbstractString, AbstractVector{<:AbstractString}} = "",
                  datacategories::Union{AbstractString, AbstractVector{<:AbstractString}} = "",
                  locations::Union{AbstractString, AbstractVector{<:AbstractString}} = "",
                  stations::Union{AbstractString, AbstractVector{<:AbstractString}} = "",
                  startdate::Date = Date("0001-01-01"),
                  enddate::Date = today())

For obtaining a CDO_token: [Request Web Services Token](https://www.ncdc.noaa.gov/cdo-web/token)

For additional information visit the [NCDC's Climate Data Online (CDO) Web Services v2 Documentation](https://www.ncdc.noaa.gov/cdo-web/webservices/v2#dataTypes)
"""
function cdo_datatypes(CDO_token::AbstractString, datatype::AbstractString)
    isvalid_cdotoken(CDO_token) || throw(CDO_NonValidToken)
    output = parse(CDO_DataType(CDO_token, datatype))
    isa(output, Exception) && throw(output)
    return output
end
function cdo_datatypes(CDO_token::AbstractString;
                       datasets::Union{AbstractString, AbstractVector{<:AbstractString}} = "",
                       datacategories::Union{AbstractString, AbstractVector{<:AbstractString}} = "",
                       locations::Union{AbstractString, AbstractVector{<:AbstractString}} = "",
                       stations::Union{AbstractString, AbstractVector{<:AbstractString}} = "",
                       startdate::Date = Date("0001-01-01"),
                       enddate::Date = today())
    isvalid_cdotoken(CDO_token) || throw(CDO_NonValidToken)
    isempty(datasets) || isvalid_datatypes(datasets) || throw(CDO_NonValidDatasets)
    isempty(datacategories) || isvalid_datacategories(datacategories) || throw(CDO_NonValidDataCategories)
    isempty(locations) || isvalid_locations(locations) || throw(CDO_NonValidLocations)
    isempty(stations) || isvalid_stations(stations) || throw(CDO_NonValidStations)
    return parse(CDO_DataTypes(CDO_token, datasets, datacategories, locations, stations, startdate, enddate))
end

# Location Categories
"""
    cdo_locationcategories(CDO_token::AbstractString, locationcategory::AbstractString)
    cdo_locationcategories(CDO_token::AbstractString;
                           datasets::Union{AbstractString, AbstractVector{<:AbstractString}} = "",
                           startdate::Date = Date("0001-01-01"),
                           enddate::Date = today())

For obtaining a CDO_token: [Request Web Services Token](https://www.ncdc.noaa.gov/cdo-web/token)

For additional information visit the [NCDC's Climate Data Online (CDO) Web Services v2 Documentation](https://www.ncdc.noaa.gov/cdo-web/webservices/v2#locationCategories)
"""
function cdo_locationcategories(CDO_token::AbstractString, locationcategory::AbstractString)
    isvalid_cdotoken(CDO_token) || throw(CDO_NonValidToken)
    output = parse(CDO_LocationCategory(CDO_token, locationcategory))
    isa(output, Exception) && throw(output)
    return output
end
function cdo_locationcategories(CDO_token::AbstractString;
                                datasets::Union{AbstractString, AbstractVector{<:AbstractString}} = "",
                                startdate::Date = Date("0001-01-01"),
                                enddate::Date = today())
    isvalid_cdotoken(CDO_token) || throw(CDO_NonValidToken)
    isempty(datasets) || isvalid_datatypes(datasets) || throw(CDO_NonValidDatasets)
    return parse(CDO_LocationCategories(CDO_token, datasets, startdate, enddate))
end

# Locations
"""
    cdo_locations(CDO_token::AbstractString, location::AbstractString)
    cdo_locations(CDO_token::AbstractString;
                  datasets::Union{AbstractString, AbstractVector{<:AbstractString}} = "",
                  locationcategories::Union{AbstractString, AbstractVector{<:AbstractString}} = "",
                  datacategories::Union{AbstractString, AbstractVector{<:AbstractString}} = "",
                  startdate::Date = Date("0001-01-01"),
                  enddate::Date = today())

For obtaining a CDO_token: [Request Web Services Token](https://www.ncdc.noaa.gov/cdo-web/token)

For additional information visit the [NCDC's Climate Data Online (CDO) Web Services v2 Documentation](https://www.ncdc.noaa.gov/cdo-web/webservices/v2#locations)
"""
function cdo_locations(CDO_token::AbstractString, location::AbstractString)
    isvalid_cdotoken(CDO_token) || throw(CDO_NonValidToken)
    output = parse(CDO_Location(CDO_token, location))
    isa(output, Exception) && throw(output)
    return output
end
function cdo_locations(CDO_token::AbstractString;
                       datasets::Union{AbstractString, AbstractVector{<:AbstractString}} = "",
                       locationcategories::Union{AbstractString, AbstractVector{<:AbstractString}} = "",
                       datacategories::Union{AbstractString, AbstractVector{<:AbstractString}} = "",
                       startdate::Date = Date("0001-01-01"),
                       enddate::Date = today())
    isvalid_cdotoken(CDO_token) || throw(CDO_NonValidToken)
    isempty(locationcategories) || isvalid_locationcategories(locationcategories) || throw(CDO_NonValidLocationCategories)
    return parse(CDO_Locations(CDO_token, datasets, locationcategories, datacategories, startdate, enddate))
end

# Stations
"""
    cdo_stations(CDO_token::AbstractString, station::AbstractString)
    cdo_stations(CDO_token::AbstractString;
                 datasets::Union{AbstractString, AbstractVector{<:AbstractString}} = "",
                 locations::Union{AbstractString, AbstractVector{<:AbstractString}} = "",
                 datacategories::Union{AbstractString, AbstractVector{<:AbstractString}} = "",
                 datatypes::Union{AbstractString, AbstractVector{<:AbstractString}} = "",
                 extent::Union{AbstractVector{<:AbstractFloat}} = Vector{Float64}(),
                 startdate::Date = Date("0001-01-01"),
                 enddate::Date = today())

For obtaining a CDO_token: [Request Web Services Token](https://www.ncdc.noaa.gov/cdo-web/token)

For additional information visit the [NCDC's Climate Data Online (CDO) Web Services v2 Documentation](https://www.ncdc.noaa.gov/cdo-web/webservices/v2#stations)
"""
function cdo_stations(CDO_token::AbstractString, station::AbstractString)
    isvalid_cdotoken(CDO_token) || throw(CDO_NonValidToken)
    output = parse(CDO_Station(CDO_token, station))
    isa(output, Exception) && throw(output)
    return output
end
function cdo_stations(CDO_token::AbstractString;
                      datasets::Union{AbstractString, AbstractVector{<:AbstractString}} = "",
                      locations::Union{AbstractString, AbstractVector{<:AbstractString}} = "",
                      datacategories::Union{AbstractString, AbstractVector{<:AbstractString}} = "",
                      datatypes::Union{AbstractString, AbstractVector{<:AbstractString}} = "",
                      extent::AbstractVector{<:AbstractFloat} = Vector{Float64}(),
                      startdate::Date = Date("0001-01-01"),
                      enddate::Date = today())
    isvalid_cdotoken(CDO_token) || throw(CDO_NonValidToken)
    isempty(datasets) || isvalid_datsets(datasets) || throw(CDO_NonValidDatasets)
    isempty(locations) || isvalid_locations(locations) || throw(CDO_NonValidLocations)
    isempty(datacategories) || isvalid_datacategories(datacategories) || throw(CDO_NonValidDataCategories)
    isempty(datatypes) || isvalid_datatypes(datatypes) || throw(CDO_NonValidDataTypes)
    isvalid_extent(extent) || throw(CDO_NonValidExtent)
    return parse(CDO_Stations(CDO_token, datasets, locations, datacategories, datatypes, extent, startdate, enddate))
end


text:
old_pwd = pwd()
cd(dirname(@__FILE__))

@testset "path utilities" begin
    @test Atom.realpath′("./utils.jl") == realpath(@__FILE__)
    @test Atom.realpath′(".././dontexist") == ".././dontexist"

    @test Atom.isuntitled(".\\untitled-asdj2123:12") == true
    @test Atom.isuntitled(".\\untitled-asdj2cx3213") == true
    @test Atom.isuntitled("./untitled-asdj2124o3:12") == true
    @test Atom.isuntitled("./untitled-a1j2124o3:12") == true
    @test Atom.isuntitled("untitled-1651asdsd23") == true
    @test Atom.isuntitled("untitled1651asdsd23") == false
    @test Atom.isuntitled("untitled-1651asdsd23:as") == false
    @test Atom.isuntitled("./utils.jl") == false
    @test Atom.isuntitled("../test/utils.jl") == false

    @test Atom.pkgpath(@__FILE__) == "utils.jl"
    @test Atom.pkgpath("foo/bar/pkgname/src/foobar.jl") == "bar/pkgname/src/foobar.jl"
    @test Atom.pkgpath("foo\\bar\\pkgname ü\\src\\foobar.jl") == "bar\\pkgname ü\\src\\foobar.jl"

    @test Atom.fullpath("untitled-asdj2cx3213") == "untitled-asdj2cx3213"
    @test Atom.fullpath("/test/foobar.jl") == "/test/foobar.jl"
    @test joinpath(split(Atom.fullpath("foobar.jl"), Base.Filesystem.pathsep())[end-1:end]...)  == joinpath("base", "foobar.jl")

    @test Atom.appendline("test.jl", -1) == "test.jl"
    @test Atom.appendline("test.jl", 0) == "test.jl"
    @test Atom.appendline("test.jl", 10) == "test.jl:10"

    @test isfile(Atom.expandpath(Atom.view(first(methods(rand)))[2].file)[2])
end

@testset "methods and docs" begin
    @test first(Atom.getmethods("Main", "Atom.getmethods")) == first(methods(Atom.getmethods))
    @test !isempty(Atom.getmethods("Main", "@deprecate"))
    @test !isempty(Atom.getmethods("Main", "Base.@deprecate"))

    @test !isempty(Atom.getdocs("Main", "Atom.getmethods"))
    @test !isempty(Atom.getdocs("Main", "@deprecate"))
end

#TODO: baselink, edit

cd(old_pwd)


text:
"""
    add_service!(device::Device, service::Service)

Add a service to a device.

Throws ArgumentError if the service is already attached to the device.
"""
function add_service!(device::Device, service::Service)
    services = get_services(device)
    for _service in services
        if IS.get_uuid(service) == IS.get_uuid(_service)
            throw(ArgumentError("service $(get_name(service)) is already attached to $(get_name(device))"))
        end
    end

    push!(services, service)
    @debug "Add $service to $(get_name(device))"
end

"""
    remove_service!(device::Device, service::Service)

Remove a service from a device.

Throws ArgumentError if the service is not attached to the device.
"""
function remove_service!(device::Device, service::Service)
    if !_remove_service!(device, service)
        throw(ArgumentError("service $(get_name(service)) was not attached to $(get_name(device))"))
    end
end

"""
    has_service(device::Device, service::Service)

Return true if the service is attached to the device.
"""
function has_service(device::Device, service::Service)
    for _service in get_services(device)
        if IS.get_uuid(_service) == IS.get_uuid(service)
            return true
        end
    end

    return false
end

"""
    has_service(device::Device, ::Type{T}) where T <:Service

Return true if a service with type T is attached to the device.
"""
function has_service(device::Device, ::Type{T}) where {T <: Service}
    for _service in get_services(device)
        if isa(_service, T)
            return true
        end
    end

    return false
end

"""
Remove service from device if it is attached.
"""
function _remove_service!(device::Device, service::Service)
    removed = false
    services = get_services(device)

    # The expectation is that there won't be many services in each device, and so
    # a faster lookup method is not needed.
    for (i, _service) in enumerate(services)
        if IS.get_uuid(_service) == IS.get_uuid(service)
            deleteat!(services, i)
            removed = true
            @debug "Removed service $(get_name(service)) from $(get_name(device))"
            break
        end
    end

    return removed
end

"""
    clear_services!(device::Device)

Remove all services attached to the device.
"""
function clear_services!(device::Device)
    services = get_services(device)
    empty!(services)
end

# All subtypes of Device define services::Vector{Service}. The values get populated
# with references to existing services. The following functions override JSON encoding to
# replace the devices with their UUIDs. This solves two problems:
#
# 1. The information is redundant, consuming unnecessary space. Also, a user could modify
#    information in the JSON file in one place but not the other, which would be very
#    problematic.
# 2. The services are stored in an vector of abstract types. JSON doesn't provide a
#    way to encode the actual concrete type, so deserialization would have to infer the type.
#    There is no guessing if the UUIDs are stored instead. The deserialization process can
#    replace the references with actual devices.

function JSON2.write(io::IO, device::T) where {T <: Device}
    return JSON2.write(io, encode_for_json(device))
end

function JSON2.write(device::T) where {T <: Device}
    return JSON2.write(encode_for_json(device))
end

"""
Encode composed buses and services as UUIDs.
"""
function encode_for_json(device::T) where {T <: Device}
    fields = fieldnames(T)
    vals = []

    for name in fields
        val = getfield(device, name)
        if val isa Vector{Service}
            push!(vals, IS.get_uuid.(val))
        elseif val isa Bus
            push!(vals, IS.get_uuid(val))
        else
            push!(vals, val)
        end
    end

    return NamedTuple{fields}(vals)
end

"""
Creates a Device object by decoding the data that was in JSON. This data stores the
values for buses and services as UUIDs, so this will lookup each in component_cache.
"""
function IS.convert_type(
    ::Type{T},
    data::NamedTuple,
    component_cache::Dict,
) where {T <: Device}
    @debug "convert_type" T data
    values = []
    for (fieldname, fieldtype) in zip(fieldnames(T), fieldtypes(T))
        val = getfield(data, fieldname)
        if fieldtype <: Vector{Service}
            services = Service[]
            for item in val
                uuid = Base.UUID(item.value)
                service = component_cache[uuid]
                push!(services, service)
            end
            push!(values, services)
        elseif fieldtype <: Bus
            uuid = Base.UUID(val.value)
            bus = component_cache[uuid]
            push!(values, bus)
        elseif fieldtype <: Component
            # Recurse.
            push!(values, IS.convert_type(fieldtype, val, component_cache))
        else
            obj = IS.convert_type(fieldtype, val)
            push!(values, obj)
        end
    end

    return T(values...)
end

function IS.convert_type(::Type{Device}, data::Any)
    error("This form of convert_type is not supported for Devices")
end


text:
######################################################
#
# Canonical Binomial Cascade Function
#
# Date:   2014-09-29
# Author: Hessel Tuinhof
#
# Function to construct a canonical binomial cascade.
# The canonical measure only preserves mass on average.
#
# Reference:
# - Riedi, Rudolf H.: Multifractal Processes
#   (http://www.stat.rice.edu/~riedi/Publ/PDF/MP.pdf)
# - Calvet, L. E.; Fisher, A. J.: Multifractal Volatility
#
# Required packages: Distributions
#
# Input:
#   m0 ... submass, float in (0,1)
#   k ... number of cascades, positive integer
#   startmass ... total mass, default: 1
#
# Output:
#   mu ... vector contains all measures
#
######################################################

using Distributions

function BinomialCanonical(m0, k, startmass = 1)
	mu = zeros(2)
	M = zeros(2)
	for i = 1:2
		rand(Bernoulli()) == 1 ? M[i] = m0 : M[i] = 1-m0
		mu[i] = M[i]*startmass
	end
	if k > 1
		tmp = Float64[]
		mu = append!(tmp, [BinomialCanonical(m0, k-1, mu[1]), BinomialCanonical(m0, k-1, mu[2])])
	end
	return mu
end

max_stars_repo_name: hesseltuinhof/masterthesis

max_issues_repo_name: hesseltuinhof/masterthesis

max_forks_repo_name: hesseltuinhof/masterthesis


text:
# Parameter sampling distributions

using Distributions

struct MetropolisParam
    name::AbstractString
    descr::AbstractString
    prior_Distribution::Distribution
    rwalk_step_size_std::Float64
end

# a helper to get priors for sampling posterior distributions.
# the idea is to keep the original definition (identified by param_idx) while
# setting all others to a narrow gaussian around the MAP_values
# NOT USED eventually!
function get_prior_for_cond_sampling(orig_prior:: Array, MAP_values::Array, param_idx::Int; std_posterior=0.01, rwalk_step_size_std = 0.0001)
p_post = deepcopy(orig_prior)
    for (i, p) in enumerate(orig_prior)
        if i == param_idx
            # keep the original definition
        else
            # create a new, narrow gaussian, centered at MAP_values[i]
            # Truncate to the same interval as the orig distribution
            param_def = MetropolisParam(
                p.name,
                p.descr,
                Distributions.TruncatedNormal(MAP_values[i], std_posterior, minimum(p.prior_Distribution), maximum(p.prior_Distribution)) ,
                rwalk_step_size_std
                )
            # print(param_def)
            p_post[i] = param_def
        end
    end
p_post
end

function param_Uniform(name, lo=0.0, hi=1. ; description="generic param", rwalk_step_size_std = 0.008)
    prior = Distributions.Uniform(lo, hi)
    MetropolisParam(name, description, prior, rwalk_step_size_std)
end

# some helpers to quickly get a reasonable default
function paramAlpha(lo=0.0, hi=1.; name="alpha", rwalk_step_size_std = 0.008)
    paramAlpha(Distributions.Uniform(lo, hi); name=name,rwalk_step_size_std = rwalk_step_size_std)
end
function paramAlpha(prior::Distribution ;name="alpha", rwalk_step_size_std = 0.008)
    MetropolisParam(name, "learning rate", prior, rwalk_step_size_std)
end

function paramEta(lo=0.0, hi=1.; name="eta", rwalk_step_size_std = 0.008)
    paramEta(Distributions.Uniform(lo, hi); name=name, rwalk_step_size_std = rwalk_step_size_std)
end
function paramEta(prior::Distribution; name="eta", rwalk_step_size_std = 0.008)
    MetropolisParam(name, "learning rate", prior, rwalk_step_size_std)
end

function paramPhi(lo=0.0, hi=1.; rwalk_step_size_std = 0.008)
    paramPhi(Distributions.Uniform(lo, hi); rwalk_step_size_std = rwalk_step_size_std)
end
function paramPhi(prior::Distribution;rwalk_step_size_std = 0.008)
    MetropolisParam("phi", "surprise trial learning rate", prior, rwalk_step_size_std)
end

function paramExpDecayOffset(lo=0, hi=1.; rwalk_step_size_std = 0.008)
    paramExpDecayOffset(Distributions.Uniform(lo, hi); rwalk_step_size_std = rwalk_step_size_std)
end
function paramExpDecayOffset(prior::Distribution; rwalk_step_size_std = 0.008)
    MetropolisParam("decay_offset", "exp decay offset", prior, rwalk_step_size_std)
end

function paramExpDecaySlope(lo= -0.10, hi=+1.0; rwalk_step_size_std = 0.004)
    #note the -.1 exponent would allow for INCREASING m-b contribution. exp(-t*slope)
    paramExpDecaySlope(Distributions.Uniform(lo, hi); rwalk_step_size_std = rwalk_step_size_std)
end
function paramExpDecaySlope(prior::Distribution; rwalk_step_size_std = 0.004)
    MetropolisParam("decay_slope", "exp decay slope", prior, rwalk_step_size_std)
end

function paramGamma(lo=0.0, hi=1.; rwalk_step_size_std = 0.008)
    paramGamma(Distributions.Uniform(lo, hi); rwalk_step_size_std = rwalk_step_size_std)
end
function paramGamma(prior::Distribution; rwalk_step_size_std = 0.008)
    MetropolisParam("gamma", "discount factor", prior, rwalk_step_size_std)
end
function paramLambda(lo=0., hi=1.; name="lambda", rwalk_step_size_std = 0.008)
    paramLambda(Distributions.Uniform(lo, hi); name=name, rwalk_step_size_std = rwalk_step_size_std)
end

function paramLambda(prior::Distribution; name="lambda", rwalk_step_size_std = 0.008)
    MetropolisParam(name, "eligibility trace", prior, rwalk_step_size_std)
end

function paramETTau(lo=0., hi=40.; name="ET_tau", rwalk_step_size_std = 0.05)
    paramETTau(Distributions.Uniform(lo, hi); name=name, rwalk_step_size_std = rwalk_step_size_std)
end
function paramETTau(prior::Distribution; name="ET_tau", rwalk_step_size_std = 0.05)
    MetropolisParam(name, "eligibility trace", prior, rwalk_step_size_std)
end

function paramTemp(lo=0.01, hi=1.25; rwalk_step_size_std = 0.05) #rwalk_step_size_std = 0.008)
    paramTemp(Distributions.Uniform(lo, hi); rwalk_step_size_std = rwalk_step_size_std)
end
function paramTemp(prior::Distribution; rwalk_step_size_std = 0.05) #rwalk_step_size_std = 0.008)
    MetropolisParam("temp", "softmax temperature", prior, rwalk_step_size_std)
end

function paramBias(lo=0.0, hi=0.6; rwalk_step_size_std = 0.004)
    paramBias(Distributions.Uniform(lo, hi); rwalk_step_size_std = rwalk_step_size_std)
end

function paramBias(prior::Distribution; rwalk_step_size_std = 0.004)
    MetropolisParam("bias", "left/right bias", prior, rwalk_step_size_std)
end

function paramAnnealingSlope(lo=0.0, hi=0.2; name = "annealing", rwalk_step_size_std = 0.004)
    paramAnnealingSlope(Distributions.Uniform(lo, hi); name=name, rwalk_step_size_std = rwalk_step_size_std)
end

function paramAnnealingSlope(prior::Distribution; name = "annealing", rwalk_step_size_std = 0.004)
    MetropolisParam(name, "slope s in alpha*exp(-t*s)", prior, rwalk_step_size_std)
end

function paramSurpriseProbability(lo=0., hi=1.; name="surpriseprobability", rwalk_step_size_std = 0.004)
    paramSurpriseProbability(Distributions.Uniform(lo, hi); name=name, rwalk_step_size_std = rwalk_step_size_std)
end
function paramSurpriseProbability(prior::Distribution; name="surpriseprobability", rwalk_step_size_std = 0.004)
    MetropolisParam(name, "surpriseprobability in particle filter", prior, rwalk_step_size_std)
end

function paramStochasticity(lo=0., hi=20.; name="stochasticity", rwalk_step_size_std = 0.05)
    paramStochasticity(Distributions.Uniform(lo, hi); name=name, rwalk_step_size_std = rwalk_step_size_std)
end
function paramStochasticity(prior::Distribution; name="stochasticity", rwalk_step_size_std = 0.05)
    MetropolisParam(name, "stochasticity in particle filter", prior, rwalk_step_size_std)
end

function paramEtaleak(lo=0., hi=1.; name="etaleak", rwalk_step_size_std = 0.008)
    paramEtaleak(Distributions.Uniform(lo, hi); name=name, rwalk_step_size_std = rwalk_step_size_std)
end
function paramEtaleak(prior::Distribution; name="etaleak", rwalk_step_size_std = 0.008)
    MetropolisParam(name, "etaleak in leaky Tsas' estimator", prior, rwalk_step_size_std)
end

function paramNrUpdateCycles(lo=0, hi=9; name="nr_upd_cycles", rwalk_step_size_std = 1.)
    paramNrUpdateCycles(Distributions.DiscreteUniform(lo, hi); name=name, rwalk_step_size_std = rwalk_step_size_std)
end
function paramNrUpdateCycles(prior::Distribution; name="nr_upd_cycles", rwalk_step_size_std = 1.)
    MetropolisParam(name, "Nr of update cycles in Prioritized Sweeping", prior, rwalk_step_size_std)
end

function paramNParticles(lo=1, hi=15; name="nparticles", rwalk_step_size_std = 1.)
    paramNParticles(Distributions.DiscreteUniform(lo, hi); name=name, rwalk_step_size_std = rwalk_step_size_std)
end
function paramNParticles(prior::Distribution; name="nparticles", rwalk_step_size_std = 1.)
    MetropolisParam(name, "Nr of particles in Particle Filter transition learner", prior, rwalk_step_size_std)
end

function paramAlpha_0(lo=0.0, hi=1.; name="alpha_0", rwalk_step_size_std = 0.008)
    paramAlpha_0(Distributions.Uniform(lo, hi); name=name,rwalk_step_size_std = rwalk_step_size_std)
end
function paramAlpha_0(prior::Distribution ;name="alpha_0", rwalk_step_size_std = 0.008)
    MetropolisParam(name, "learning rate constant term", prior, rwalk_step_size_std)
end

function paramAlpha_surprise(lo=0.0, hi=1.; name="alpha_surprise", rwalk_step_size_std = 0.008)
    paramAlpha_surprise(Distributions.Uniform(lo, hi); name=name,rwalk_step_size_std = rwalk_step_size_std)
end
function paramAlpha_surprise(prior::Distribution ;name="alpha_surprise", rwalk_step_size_std = 0.008)
    MetropolisParam(name, "learning rate modulating term", prior, rwalk_step_size_std)
end

function paramWeight(lo=0.0, hi=20.; name="w_MB", rwalk_step_size_std = 0.05)
    paramAlpha_surprise(Distributions.Uniform(lo, hi); name=name,rwalk_step_size_std = rwalk_step_size_std)
end
function paramWeight(prior::Distribution ;name="w_MB", rwalk_step_size_std = 0.05)
    MetropolisParam(name, "weight in hybrid", prior, rwalk_step_size_std)
end

function paramActionCost(lo=0.0, hi=5.; name="actioncost", rwalk_step_size_std = 0.05)
    paramActionCost(Distributions.Uniform(lo, hi); name=name, rwalk_step_size_std = rwalk_step_size_std)
end
function paramActionCost(prior::Distribution ;name="actioncost", rwalk_step_size_std = 0.05)
    MetropolisParam(name, "action cost", prior, rwalk_step_size_std)
end


text:
struct VertexCount end

@inline (::VertexCount)(j, j′, k...) = Int(j′ - j)
@inline cost_type(::Type{VertexCount}) = Int
oracle_model(::VertexCount) = VertexCount()
oracle_stripe(::AbstractHint, ::VertexCount, ::SparseMatrixCSC; kwargs...) = VertexCount()
#bound_stripe(A::SparseMatrixCSC, K, ::VertexCount) = (size(A)[2]/K, size(A)[2]) $TODO ?

struct PinCount{Ti} <: AbstractMatrix{Ti}
    n::Int
    pos::Vector{Ti}
end

#@inline (::PinCount)(j, j′, k...) = Int(j′ - j)
#@inline cost_type(::Type{VertexCount}) = Int
#oracle_model(::VertexCount) = VertexCount()
#oracle_stripe(::AbstractHint, ::VertexCount, ::SparseMatrixCSC; kwargs...) = VertexCount()
#bound_stripe(A::SparseMatrixCSC, K, ::VertexCount) = (size(A)[2]/K, size(A)[2]) $TODO ?

Base.size(arg::PinCount) = (arg.n + 1, arg.n + 1)

pincount(args...; kwargs...) = pincount(NoHint(), args...; kwargs...)
pincount(::AbstractHint, args...; kwargs...) = @assert false
pincount(hint::AbstractHint, A::SparseMatrixCSC; kwargs...) =
    pincount!(hint, size(A)..., nnz(A), A.colptr; kwargs...)

pincount!(args...; kwargs...) = pincount!(NoHint(), args...; kwargs...)
pincount!(::AbstractHint, args...; kwargs...) = @assert false
pincount!(hint::AbstractHint, m, n, N, pos::Vector{Ti}; kwargs...) where {Ti} =
    PinCount(hint, n, pos; kwargs...)

PinCount(hint::AbstractHint, n, pos::Vector{Ti}; kwargs...) where {Ti} =
    PinCount{Ti}(hint, n, pos; kwargs...)
function PinCount{Ti}(hint::AbstractHint, n, pos::Vector{Ti}; kwargs...) where {Ti}
    return PinCount(n, pos)
end

Base.getindex(arg::PinCount{Ti}, j::Integer, j′::Integer) where {Ti} = arg(j, j′)
function (arg::PinCount{Ti})(j::Integer, j′::Integer) where {Ti}
    @inbounds begin
        return arg.pos[j′] - arg.pos[j]
    end
end

struct NetCount{Ti, Lnk} <: AbstractMatrix{Ti}
    n::Int
    pos::Vector{Ti}
    lnk::Lnk
end

Base.size(arg::NetCount) = (arg.n + 1, arg.n + 1)

netcount(args...; kwargs...) = netcount(NoHint(), args...; kwargs...)
netcount(::AbstractHint, args...; kwargs...) = @assert false
netcount(hint::AbstractHint, A::SparseMatrixCSC; kwargs...) =
    netcount!(hint, size(A)..., nnz(A), A.colptr, A.rowval; kwargs...)

netcount!(args...; kwargs...) = netcount!(NoHint(), args...; kwargs...)
netcount!(::AbstractHint, args...; kwargs...) = @assert false
netcount!(hint::AbstractHint, m, n, N, pos::Vector{Ti}, idx::Vector{Ti}; kwargs...) where {Ti} =
    NetCount(hint, m, n, N, pos, idx; kwargs...)

dianetcount(args...; kwargs...) = dianetcount(NoHint(), args...; kwargs...)
dianetcount(::AbstractHint, args...; kwargs...) = @assert false
dianetcount(hint::AbstractHint, A::SparseMatrixCSC; kwargs...) =
    dianetcount!(hint, size(A)..., nnz(A), A.colptr, A.rowval; kwargs...)

dianetcount!(args...; kwargs...) = dianetcount!(NoHint(), args...; kwargs...)
dianetcount!(::AbstractHint, args...; kwargs...) = @assert false
function dianetcount!(hint::AbstractHint, m, n, N, pos::Vector{Ti}, idx::Vector{Ti}; kwargs...) where {Ti}
    @inbounds begin
        hst = zeros(Ti, m)
        pos′ = undefs(Ti, n + 1)
        idx′ = undefs(Ti, N + n)

        q′ = 1
        for j = 1:n
            pos′[j] = q′
            for q in pos[j] : pos[j + 1] - 1
                i = idx[q]
                idx′[q′] = (n + 1) - hst[i]
                hst[i] = j
                q′ += 1
            end
            if hst[j] < j
                idx′[q′] = (n + 1) - hst[j]
                hst[j] = j
                q′ += 1
            end
        end
        pos′[n + 1] = q′
        N′ = q′ - 1
        resize!(idx′, N′)

        return NetCount(n, pos′, dominancecount!(hint, n + 1, n + 1, N′, pos′, idx′; kwargs...))
    end
end

NetCount(hint::AbstractHint, m, n, N, pos::Vector{Ti}, idx::Vector{Ti}; kwargs...) where {Ti} =
    NetCount{Ti}(hint, m, n, N, pos, idx; kwargs...)
function NetCount{Ti}(hint::AbstractHint, m, n, N, pos::Vector{Ti}, idx::Vector{Ti}; kwargs...) where {Ti}
    @inbounds begin
        hst = zeros(Ti, m)
        idx′ = undefs(Ti, N)

        for j = 1:n
            for q in pos[j] : pos[j + 1] - 1
                i = idx[q]
                idx′[q] = (n + 1) - hst[i]
                hst[i] = j
            end
        end

        return NetCount(n, pos, dominancecount!(hint, n + 1, n + 1, N, pos, idx′; kwargs...))
    end
end

Base.getindex(arg::NetCount{Ti}, j::Integer, j′::Integer) where {Ti} = arg(j, j′)
function (arg::NetCount{Ti})(j::Integer, j′::Integer) where {Ti}
    @inbounds begin
        return (arg.pos[j′] - arg.pos[j]) - arg.lnk((arg.n + 2) - j, j′)
    end
end

@propagate_inbounds function (stp::Step{Net})(_j::Same, _j′) where {Ti, Net <: NetCount{Ti}}
    begin
        j = destep(_j)
        j′ = destep(_j′)
        arg = stp.ocl
        return (arg.pos[j′] - arg.pos[j]) - Step(arg.lnk)(Same((arg.n + 2) - j), _j′)
    end
end

@propagate_inbounds function (stp::Step{Net})(_j::Next, _j′::Same) where {Ti, Net <: NetCount{Ti}}
    begin
        j = destep(_j)
        j′ = destep(_j′)
        arg = stp.ocl
        return (arg.pos[j′] - arg.pos[j]) - Step(arg.lnk)(Prev((arg.n + 2) - j), _j′)
    end
end

@propagate_inbounds function (ocl::Step{Net})(_j::Prev, _j′::Same) where {Ti, Net <: NetCount{Ti}}
    begin
        j = destep(_j)
        j′ = destep(_j′)
        arg = ocl.ocl
        return (arg.pos[j′] - arg.pos[j]) - Step(arg.lnk)(Next((arg.n + 2) - j), _j′)
    end
end

struct SelfNetCount{Ti, Lnk} <: AbstractMatrix{Ti}
    n::Int
    lnk::Lnk
end

SelfNetCount{Ti}(n, lnk::Lnk) where {Ti, Lnk} = SelfNetCount{Ti, Lnk}(n, lnk)

Base.size(arg::SelfNetCount) = (arg.n + 1, arg.n + 1)

selfnetcount(args...; kwargs...) = selfnetcount(NoHint(), args...; kwargs...)
selfnetcount(::AbstractHint, args...; kwargs...) = @assert false
selfnetcount(hint::AbstractHint, A::SparseMatrixCSC; kwargs...) =
    selfnetcount!(hint, size(A)..., nnz(A), A.colptr, A.rowval; kwargs...)

selfnetcount!(args...; kwargs...) = selfnetcount!(NoHint(), args...; kwargs...)
selfnetcount!(::AbstractHint, args...; kwargs...) = @assert false
selfnetcount!(hint::AbstractHint, m, n, N, pos::Vector{Ti}, idx::Vector{Ti}; kwargs...) where {Ti} =
    SelfNetCount(hint, m, n, N, pos, idx; kwargs...)

SelfNetCount(hint::AbstractHint, m, n, N, pos::Vector{Ti}, idx::Vector{Ti}; kwargs...) where {Ti} =
    SelfNetCount{Ti}(hint, m, n, N, pos, idx; kwargs...)
function SelfNetCount{Ti}(hint::AbstractHint, m, n, N, pos::Vector{Ti}, idx::Vector{Ti}; kwargs...) where {Ti}
    @inbounds begin
        hst = zeros(Ti, m)
        hst′ = undefs(Ti, m)
        idx′ = undefs(Ti, N)
        pos′ = zeros(Ti, n + 1)

        for j = 1:n
            for q in pos[j] : pos[j + 1] - 1
                i = idx[q]
                if hst[i] == 0
                    hst[i] = j
                end
                hst′[i] = j
            end
        end

        for i = 1:m
            if hst[i] != 0
                j′ = hst′[i]
                pos′[j′ + 1] += 1
            end
        end

        q = 1
        for j = 1:n + 1
            (pos′[j], q) = (q, q + pos′[j])
        end

        N′ = q - 1

        idx′ = undefs(Ti, N′)

        for i = 1:m
            if hst[i] != 0
                j = hst[i]
                j′ = hst′[i]
                q = pos′[j′ + 1]
                idx′[q] = (n + 1) - j
                pos′[j′ + 1] = q + 1
            end
        end

        return SelfNetCount{Ti}(n, dominancecount!(hint, n + 1, n + 1, N′, pos′, idx′; kwargs...))
    end
end

Base.getindex(arg::SelfNetCount{Ti}, j::Integer, j′::Integer) where {Ti} = arg(j, j′)
function (arg::SelfNetCount{Ti})(j::Integer, j′::Integer) where {Ti}
    @inbounds begin
        return arg.lnk((arg.n + 2) - j, j′)
    end
end

@propagate_inbounds function (stp::Step{Net})(_j::Same, _j′) where {Ti, Net <: SelfNetCount{Ti}}
    begin
        j = destep(_j)
        j′ = destep(_j′)
        arg = stp.ocl
        return Step(arg.lnk)(Same((arg.n + 2) - j), _j′)
    end
end

@propagate_inbounds function (stp::Step{Net})(_j::Next, _j′::Same) where {Ti, Net <: SelfNetCount{Ti}}
    begin
        j = destep(_j)
        j′ = destep(_j′)
        arg = stp.ocl
        return Step(arg.lnk)(Prev((arg.n + 2) - j), _j′)
    end
end

@propagate_inbounds function (ocl::Step{Net})(_j::Prev, _j′::Same) where {Ti, Net <: SelfNetCount{Ti}}
    begin
        j = destep(_j)
        j′ = destep(_j′)
        arg = ocl.ocl
        return Step(arg.lnk)(Next((arg.n + 2) - j), _j′)
    end
end

struct SelfPinCount{Ti, Lnk} <: AbstractMatrix{Ti}
    n::Int
    lnk::Lnk
end

SelfPinCount{Ti}(n, lnk::Lnk) where {Ti, Lnk} = SelfPinCount{Ti, Lnk}(n, lnk)

Base.size(arg::SelfPinCount) = (arg.n + 1, arg.n + 1)

selfpincount(args...; kwargs...) = selfpincount(NoHint(), args...; kwargs...)
selfpincount(::AbstractHint, args...; kwargs...) = @assert false
selfpincount(hint::AbstractHint, A::SparseMatrixCSC; kwargs...) =
    selfpincount!(hint, size(A)..., nnz(A), A.colptr, A.rowval; kwargs...)

selfpincount!(args...; kwargs...) = selfpincount!(NoHint(), args...; kwargs...)
selfpincount!(::AbstractHint, args...; kwargs...) = @assert false
selfpincount!(hint::AbstractHint, m, n, N, pos::Vector{Ti}, idx::Vector{Ti}; kwargs...) where {Ti} =
    SelfPinCount(hint, m, n, N, pos, idx; kwargs...)

SelfPinCount(hint::AbstractHint, m, n, N, pos::Vector{Ti}, idx::Vector{Ti}; kwargs...) where {Ti} =
    SelfPinCount{Ti}(hint, m, n, N, pos, idx; kwargs...)
function SelfPinCount{Ti}(hint::AbstractHint, m, n, N, pos::Vector{Ti}, idx::Vector{Ti}; kwargs...) where {Ti}
    @inbounds begin
        pos′ = zeros(Ti, n + 1)
        idx′ = zeros(Ti, N)

        for j = 1:n
            for q in pos[j] : pos[j + 1] - 1
                i = idx[q]
                (i, i′) = minmax(j, i)
                pos′[i′ + 1] += 1
            end
        end

        q = 1
        for j = 1:n + 1
            (pos′[j], q) = (q, q + pos′[j])
        end

        for j = 1:n
            for q in pos[j] : pos[j + 1] - 1
                i = idx[q]
                (i, i′) = minmax(j, i)
                q′ = pos′[i′ + 1]
                idx′[q′] = (n + 1) - i
                pos′[i′ + 1] = q′ + 1
            end
        end

        return SelfPinCount{Ti}(n, dominancecount!(hint, n + 1, n + 1, N, pos′, idx′; kwargs...))
    end
end

Base.getindex(arg::SelfPinCount{Ti}, j::Integer, j′::Integer) where {Ti} = arg(j, j′)
function (arg::SelfPinCount{Ti})(j::Integer, j′::Integer) where {Ti}
    @inbounds begin
        return arg.lnk((arg.n + 2) - j, j′)
    end
end

@propagate_inbounds function (stp::Step{Net})(_j::Same, _j′) where {Ti, Net <: SelfPinCount{Ti}}
    begin
        j = destep(_j)
        j′ = destep(_j′)
        arg = stp.ocl
        return Step(arg.lnk)(Same((arg.n + 2) - j), _j′)
    end
end

@propagate_inbounds function (stp::Step{Net})(_j::Next, _j′::Same) where {Ti, Net <: SelfPinCount{Ti}}
    begin
        j = destep(_j)
        j′ = destep(_j′)
        arg = stp.ocl
        return Step(arg.lnk)(Prev((arg.n + 2) - j), _j′)
    end
end

@propagate_inbounds function (ocl::Step{Net})(_j::Prev, _j′::Same) where {Ti, Net <: SelfPinCount{Ti}}
    begin
        j = destep(_j)
        j′ = destep(_j′)
        arg = ocl.ocl
        return Step(arg.lnk)(Next((arg.n + 2) - j), _j′)
    end
end


text:
#ДАНО: Робот - Робот - в произвольной клетке ограниченного прямоугольного поля

#РЕЗУЛЬТАТ: Робот - в исходном положении, и клетки поля промакированы так: нижний ряд - полностью, следующий - весь, за исключением одной последней клетки на Востоке, следующий - за исключением двух последних клеток на Востоке, и т.д.

function triangle(r::Robot,side::HorizonSide)
    go_ugol!(r,(Sud,West))
    n = move_counter!(r,Ost)
    movements!(r,inverse(Ost))
    while !isborder(r,Nord)
        n-=1
        move!(r,Nord)
        move_control_side!(r,Ost,n)
    end
end
function move_control_side!(r::Robot,side::HorizonSide,steps::Int)
    i = 0
    putmarker!(r)
    while i  < steps
        i+=1
        move!(r,side)
        putmarker!(r)
    end
    movements!(r,inverse(side))
end
function go_ugol!(r::Robot, side::NTuple{2,HorizonSide})
    movements!(r,side[1])
    movements!(r,side[2])
end
function move_counter!(r::Robot,side::HorizonSide)
    n = 0
    putmarker!(r)
    while !isborder(r,side)
        move!(r,side)
        putmarker!(r)
        n+=1
    end
    return n
end
function movements!(r,side)
    while !isborder(r,side)
        move!(r,side)
    end
end
inverse(side::HorizonSide) = HorizonSide(mod(Int(side)+2, 4))


text:
module Model

include("encoder.jl")
include("decoder.jl")
include("embed.jl")
include("generator.jl")
include("transformers.jl")
include("beam_search.jl")
include("translator.jl")
include("abstractive.jl")
export Encoder, Decoder, WordPositionEmbed, Generator, TransformersModel, BeamSearch, Translator, BertAbs, TransformerAbs, TransformerAbsTiny
include("utils.jl")
export params_count

end # module


text:
export UpdateKeenerScores

"""
    UpdateKeenerScores

 Keener update rule, e.g., see \"Whos's #1\", Langville and Meyer, p.29
 Note Keener's method can use general inputs, this variant uses scores

## Parameters
* `skew::Bool`: skew input statistics to give more weight to small values
* `norm::Bool`: normalise input statistics by number of games played

```
"""
struct UpdateKeenerScores <: UpdateRule
    skew::Bool # skew input statistics to give more weight to small values
    norm::Bool # normalise input statistics by number of games played
end
UpdateKeenerScores(; skew=false, norm=true) = UpdateKeenerScores(skew, norm) # default no skew, but yes to normalisation

function update_info( rule::UpdateKeenerScores )
    info = Dict(
                :name => "KeenerScores",
                :reference => "\"Whos's #1\", Langville and Meyer, p.29",
                :computation => "simultaneous",
                :state_model => "none",
                :input => "outcome",
                :output => "deterministic",
                :model => "single",
                :ties => true,
                :factors => false,
                :parameters => ["skew", "norm"],
                :record => true
                )
    return info
end
update_info( ::Type{UpdateKeenerScores} ) =  update_info( UpdateKeenerScores() )

# skew function
h( x::Float64 ) = 0.5*(1 + sign(x - 0.5)*sqrt(abs(2*x - 1.0)) )

function update_ratings( rule::UpdateKeenerScores,
                         input_ratings::RatingsList,
                         input_competitions::DataFrame)
    n = size(input_competitions,1)
    m = length( input_ratings.players )
    I = player_indexes( input_ratings.players )

    # construct KeenerScores matrices and vectors
    #   note Keener's method can use general inputs, here we use scores
    S = zeros(Int,m,m)
    n_games = zeros(Int,m)
    d = input_competitions
    for i=1:n
        n_games[ I[d[i,PlayerA]] ] += 1
        n_games[ I[d[i,PlayerB]] ] += 1
        S[ I[d[i,PlayerA]], I[d[i,PlayerB]]] += d[i,ScoreA]
        S[ I[d[i,PlayerB]], I[d[i,PlayerA]]] += d[i,ScoreB]
    end

    A = (S .+ 1) ./ (S .+ S' .+ 2)
    # skew if needed
    if rule.skew
        A = h.(A)
    end
    # normalise
    if rule.norm
        for i=1:m
            A[i,:] = A[i,:] ./ n_games[i]
        end
    end

    # solve Keener's equation
    #   assuming matrix is irreducible (all teams are linked by games)
    r = solve_keener( A )
    ratings = Dict{String, Float64}()
    for player in input_ratings.players
        ratings[player] = r[ I[player] ]
    end

    # output ratings list
    output_ratings = RatingsList(input_ratings.players, ratings )
    return output_ratings
end

function solve_keener( A::Array{Float64, 2} )
    # make irreducible
    m = size(A,1)
    epsilon = 1.0e-8
    𝐞 = ones(Float64,m)
    A = A .+ epsilon * 𝐞 * 𝐞'

    # solve by power method (see, p.41)
    r = ones(Float64, m)/m
    not_converged = true
    i = 0
    max_i = 100
    while not_converged && i < max_i
        σ = epsilon * sum(r)
        new_r = A*r + σ*𝐞
        ν = sum( new_r )
        new_r = new_r / ν
        if sum( abs.(new_r - r) ) < 1.0e-6
            not_converged = false
        end
        r = new_r
        i += 1
    end
    if max_i == i
        println("WARN:   i = $i = max_i")
    end

    return r
end


text:
###################################################################################
#                             load dependencies
###################################################################################
cd(@__DIR__)
using DifferentialEvolutionMCMC, Random, Distributions
Random.seed!(9528)
###################################################################################
#                             define functions
###################################################################################
# returns prior log likelihood
function sample_prior()
    LL = 0.0
    μβ0 = rand(Normal(1, 1))
    σβ0 = rand(truncated(Cauchy(0, 1), 0, Inf))
    σ = rand(truncated(Cauchy(0, 1), 0, Inf))
    β0 = rand(Normal(0.0, σβ0), n_subj)
    return as_union([μβ0, σβ0, β0, σ])
end

# function for initial values
function prior_loglike(μβ0, σβ0, β0, σ)
    LL = 0.0
    LL += logpdf(Normal(1, 1), μβ0)
    LL += logpdf(truncated(Cauchy(0, 1), 0, Inf), σβ0)
    LL += sum(logpdf.(Normal(0, σβ0), β0))
    LL += logpdf(truncated(Cauchy(0, 1), 0, Inf), σ)
    return LL
end

# likelihood function
function loglike(data, μβ0, σβ0, β0, σ)
    LL = 0.0
    for s in 1:length(data)
        μ = μβ0 + β0[s]
        y = data[s] .- μ
        LL += sum(logpdf(Normal(0, σ), y))
    end
    return LL
end
###################################################################################
#                             generate data
###################################################################################
# number of responses
n_data = 50
# number of subjects
n_subj = 50
# mean intercept
μβ0 = 1.0
# standard deviation intercept
σβ0 = 1.0
# subject intercept deviations
β0 = rand(Normal(0.0, σβ0), n_subj)
σ = .5

sim(μβ0, β0, σ, n_data) = rand(Normal(μβ0 + β0, σ), n_data)
data = [sim(μβ0, β0[i], σ, n_data) for i in 1:n_subj]
###################################################################################
#                             configure sampler
###################################################################################
# parameter bounds
bounds = (
    (-Inf,Inf),
    (0.0,Inf),
    (-Inf,Inf),
    (0.0,Inf),
)

# parameter names
names = (:μβ0, :σβ0, :β0, :σ)

# model object
model = DEModel(;
    sample_prior,
    prior_loglike,
    loglike,
    data,
    names
)

# block update indicator
# update hyper parameters first
# update lower level parameters second
blocks = [
    [true,true,fill(false, n_subj),true],
    [false,false,fill(true, n_subj),false],
]
blocks = as_union.(blocks)

# blocks = [[true,true,fill(false, n_subj),true],]
# subj_blocks = [[false,false,fill(false, n_subj),false] for i in 1:n_subj]
# map(i -> subj_blocks[i][3][i] = true, 1:n_subj)
# push!(blocks, subj_blocks...)
# blocks = as_union.(blocks)

# use block updating on each iteration
blocking_on = x -> true

# sampler object
de = DE(;
    sample_prior,
    bounds,
    sample = resample,
    burnin = 20_000,
    n_initial = (n_subj + 1) * 4,
    Np = 3,
    n_groups = 2,
    θsnooker = 0.1,
    blocking_on,
    blocks,
    #generate_proposal = variable_gamma
)
###################################################################################
#                             estimate parameters
###################################################################################
n_iter = 40_000
chains = sample(model, de, MCMCThreads(), n_iter, progress=true)
###################################################################################
#                             Stan
###################################################################################
using StanSample
tmpdir = pwd() * "/tmp"

idx = repeat(1:n_subj, inner=n_data)
y = vcat(data...)
#xrep = repeat(x, inner=n_subj)
stan_data = Dict(
    "n_data" => n_data,
    "n_subj" => n_subj,
    "n" => n_data * n_subj,
    "idx" => idx,
    "y" => y,
    "x" => fill(0.0, n_data * n_subj),
)
###################################################################################################
#                                        Define Model
###################################################################################################
# load the Stan model file
stream = open("temp.stan","r")
stan_model = read(stream, String)
close(stream)
stan_model = SampleModel("temp", stan_model, tmpdir)
###################################################################################################
#                                     Estimate Parameters
###################################################################################################
stan_sample(
    stan_model;
    data = stan_data,
    #seed,
    num_chains = 4,
    num_samples = 1000,
    num_warmups = 1000,
    save_warmup = false
)
samples = read_samples(stan_model, :mcmcchains)
#rhats = filter(!isnan, MCMCChains.ess_rhat(samples).nt.rhat)
#mean(rhats .> 1.01)
###################################################################################
#                             Turing
###################################################################################
# using Turing, ReverseDiff, Memoization
# Turing.setadbackend(:reversediff)
# Turing.setrdcache(true)

# @model function turing_model(y, x, n_subj, idx)
#     β0 ~ Normal(1, 1)
#     β1 ~ Normal(.5, 1)
#     σ ~ truncated(Cauchy(0, 1), 0, Inf)
#     σβ0 ~ truncated(Cauchy(0, 1), 0, Inf)
#     β0i ~ filldist(Normal(0., σβ0), n_subj)

#     μs = β0 .+ β0i[idx] .+ β1 * x
#     y ~ MvNormal(μs, σ)
# end

# idx = repeat(1:n_subj, inner=n_data)
# y = vcat(data...)
# xrep = repeat(x, inner=n_subj)
# turing_samples = sample(turing_model(y, xrep, n_subj, idx), NUTS(1000, .65), MCMCThreads(), 1000, 4)


text:
using Statistics
using Gadfly
using ProgressMeter

mutable struct IsingModel
    N :: Int64
    J :: Float64
    k :: Float64
    T :: Float64
    spins :: Array{Int8, 2}
    energy :: Float64
    magnetization :: Float64
end

function Base.copy(model :: IsingModel)
    return IsingModel(model.N, model.J, model.k, model.T, copy(model.spins),
                      model.energy, model.magnetization)
end

function update!(model :: IsingModel)
    function update_magnetization!(model :: IsingModel)
        N₊ = count(i -> i == 1, model.spins)
        N₋ = count(i -> i == -1, model.spins)
        model.magnetization = (N₊ - N₋) / model.N^2
    end

    function update_energy!(model :: IsingModel)
        # periodic boundary condition
        energy = 0.0
        for x in 1:model.N
            for y in 1:model.N
                if (x == model.N) && (y == model.N)
                    energy += model.spins[x, y] * (model.spins[x, 1] + model.spins[1, y])
                elseif (x == model.N) && (y != model.N)
                    energy += model.spins[x, y] * (model.spins[x, y + 1] + model.spins[1, y])
                elseif (x != model.N) && (y == model.N)
                    energy += model.spins[x, y] * (model.spins[x, 1] + model.spins[x + 1, y])
                else
                    energy += model.spins[x, y] * (model.spins[x, y + 1] + model.spins[x + 1, y])
                end
            end
        end

        model.energy = -model.J * energy

        return nothing
    end

    update_magnetization!(model)
    update_energy!(model)

    return nothing
end

function init_model(N :: Int64, J :: Float64, k :: Float64, T :: Float64)
    if N <= 2
        error("N must be greater than 2")
    end

    model = IsingModel(N, J, k, T, rand([-1, 1], (N, N)), 0.0, 0.0)
    update!(model)

    return model
end

function metropolis(model :: IsingModel)
    candidate_model = copy(model)

    change_x, change_y = rand(1:model.N), rand(1:model.N)
    candidate_model.spins[change_x, change_y] *= Int8(-1)
    update!(candidate_model)

    if candidate_model.energy <= model.energy
        model = copy(candidate_model)
    elseif (model.energy - candidate_model.energy) / (model.k * model.T) >= log(rand())
        model = copy(candidate_model)
    end

    return model
end

function ising_simulation(N :: Int64, J :: Float64, k :: Float64, T :: Float64;
                          niters = 10000 :: Int64, burnin = 5000 :: Int64)
    magnetizations = Array{Float64}(undef, niters)
    energies = Array{Float64}(undef, niters)

    model = init_model(N, J, k, T)

    for i in 1:niters
        magnetizations[i] = model.magnetization
        energies[i] = model.energy
        model = metropolis(model)
    end

    return mean(magnetizations[burnin + 1:niters]), mean(energies[burnin + 1:niters])

    #theme = Theme(background_color = "white", line_width = 2px)
    #p = plot(x = 1:niters, y = magnetizations, theme, Geom.line,
    #         Guide.xlabel("iterations"), Guide.ylabel("magnetizations"), Guide.title("magnetizations(kT = $(k * T))"))
    #img = SVG("magnetizations_kT$(k*T).svg", 960px, 540px)
    #draw(img, p)

    #p = plot(x = 1:niters, y = energies, theme, Geom.line,
    #         Guide.xlabel("iterations"), Guide.ylabel("energies"), Guide.title("magnetizations(kT = $(k * T))"))
    #img = SVG("energies_kT$(k*T).svg", 960px, 540px)
    #draw(img, p)
end

function run_simulation()
    temperatures = range(0.001, step = 0.01, stop = 10.0)

    magnetizations = Array{Float64}(undef, length(temperatures))
    energies = Array{Float64}(undef, length(temperatures))

    progress = Progress(length(temperatures))
    for i in 1:length(temperatures)
        magnetizations[i], energies[i] = ising_simulation(20, 1.0, 1.0,
                                                          temperatures[i],
                                                          niters = 100000,
                                                          burnin = 50000)
        next!(progress)
    end

    theme = Theme(background_color = "white", line_width = 2px)
    p = plot(x = temperatures, y = abs.(magnetizations), theme, Geom.line,
             Guide.xlabel("kT"), Guide.ylabel("|magnetizations|"), Guide.title("abs magnetizations"))
    img = SVG("magnetizations_vs_kT.svg", 960px, 540px)
    draw(img, p)

    p = plot(x = temperatures, y = energies, theme, Geom.line,
             Guide.xlabel("kT"), Guide.ylabel("energies"), Guide.title("energies"))
    img = SVG("energies_vs_kT.svg", 960px, 540px)
    draw(img, p)
end

run_simulation()


text:
function car_benchmarks!(suite::BenchmarkGroup)

    # Car escape
    T = Float64

    # options
    max_con_viol = 1.0e-8
    verbose=false

    suite["parallel_park"] = BenchmarkGroup()
    suite["3obs"] = BenchmarkGroup()
    suite["escape"] = BenchmarkGroup()

    # Parallel Park
    opts_ilqr = iLQRSolverOptions{T}(verbose=verbose,
        live_plotting=:off)

    opts_al = AugmentedLagrangianSolverOptions{T}(verbose=verbose,
        opts_uncon=opts_ilqr,
        iterations=30,
        penalty_scaling=10.0,
        constraint_tolerance=max_con_viol)

    opts_pn = ProjectedNewtonSolverOptions{T}(verbose=verbose,
        feasibility_tolerance=max_con_viol)

    opts_altro = ALTROSolverOptions{T}(verbose=verbose,
        opts_al=opts_al,
        projected_newton=true,
        projected_newton_tolerance=1.0e-4)

    opts_ipopt = DIRCOLSolverOptions{T}(verbose=verbose,
        nlp=Ipopt.Optimizer(),
        feasibility_tolerance=max_con_viol)

    opts_snopt = DIRCOLSolverOptions{T}(verbose=verbose,
        nlp=SNOPT7.Optimizer(),
        feasibility_tolerance=max_con_viol)

    opts_ilqr = iLQRSolverOptions{T}(verbose=verbose,
        live_plotting=:off)
    run_benchmarks!(suite["parallel_park"], Problems.parallel_park, [opts_ilqr, opts_al, opts_altro, opts_ipopt, opts_snopt])

    # Three Obstacle
    opts_al = AugmentedLagrangianSolverOptions{T}(verbose=verbose,
        opts_uncon=opts_ilqr,
        cost_tolerance=1.0e-4,
        cost_tolerance_intermediate=1.0e-2,
        constraint_tolerance=max_con_viol,
        penalty_scaling=50.,
        penalty_initial=10.)

    opts_pn = ProjectedNewtonSolverOptions{T}(verbose=verbose,
        feasibility_tolerance=max_con_viol)

    opts_altro = ALTROSolverOptions{T}(verbose=verbose,
        opts_al=opts_al,
        opts_pn=opts_pn,
        projected_newton=true,
        projected_newton_tolerance=1.0e-3)

    opts_ipopt = DIRCOLSolverOptions{T}(verbose=verbose,
        nlp=Ipopt.Optimizer(),
        feasibility_tolerance=max_con_viol)

    opts_snopt = DIRCOLSolverOptions{T}(verbose=verbose,
        nlp=SNOPT7.Optimizer(),
        feasibility_tolerance=max_con_viol)

    run_benchmarks!(suite["3obs"], Problems.car_3obs, [opts_altro, opts_al, opts_ipopt, opts_snopt])

    # Car escape
    opts_ilqr = iLQRSolverOptions{T}(verbose=verbose,live_plotting=:off)

    opts_al = AugmentedLagrangianSolverOptions{T}(verbose=verbose,
        opts_uncon=opts_ilqr,
        cost_tolerance=1.0e-6,
        cost_tolerance_intermediate=1.0e-2,
        constraint_tolerance=max_con_viol,
        penalty_scaling=50.,
        penalty_initial=10.)

    opts_pn = ProjectedNewtonSolverOptions{T}(verbose=verbose,
        feasibility_tolerance=max_con_viol,
        solve_type=:feasible)

    opts_altro = ALTROSolverOptions{T}(verbose=verbose,
        opts_al=opts_al,
        R_inf=1.0e-1,
        resolve_feasible_problem=false,
        opts_pn=opts_pn,
        projected_newton=true,
        projected_newton_tolerance=1.0e-3);

    opts_ipopt = DIRCOLSolverOptions{T}(verbose=verbose,
        nlp=Ipopt.Optimizer(),
        feasibility_tolerance=max_con_viol)

    opts_snopt = DIRCOLSolverOptions{T}(verbose=verbose,
        nlp=SNOPT7.Optimizer(),
        feasibility_tolerance=max_con_viol)

    # run_benchmarks!(suite["escape"], Problems.car_escape, [opts_altro, opts_ipopt, opts_snopt])

end

max_issues_repo_name: MehrdadZareian/GMARA

max_forks_repo_name: MehrdadZareian/GMARA


text:
using FunctionalCollections: PersistentHashMap, assoc

#############################
# generative function trace #
#############################

mutable struct GFTrace
    call::Union{CallRecord,Nothing}
    primitive_calls::PersistentHashMap{Any,CallRecord}

    # values can be a GFTrace, or Any (foreign trace)
    subtraces::PersistentHashMap{Any,Any}
    has_choices::Bool
end

function GFTrace()
    primitive_calls = PersistentHashMap{Any,CallRecord}()
    subtraces = PersistentHashMap{Any,Any}()
    call = nothing
    GFTrace(call, primitive_calls, subtraces, false)
end

get_call_record(trace::GFTrace) = trace.call
has_choices(trace::GFTrace) = trace.has_choices

function has_primitive_call(trace::GFTrace, addr)
    haskey(trace.primitive_calls, addr)
end

function has_primitive_call(trace::GFTrace, addr::Pair)
    (first, rest) = addr
    if !haskey(trace.subtraces, first)
        return false
    end
    subtrace::GFTrace = get(trace.subtraces, first)
    has_primitive_call(subtrace, rest)
end

function has_subtrace(trace::GFTrace, addr)
    haskey(trace.subtraces, addr)
end

function has_subtrace(trace::GFTrace, addr::Pair)
    (first, rest) = addr
    if !haskey(trace.subtraces, first)
        return false
    end
    subtrace::GFTrace = get(trace.subtraces, first)
    has_subtrace(subtrace, rest)
end

function get_primitive_call(trace::GFTrace, addr)
    get(trace.primitive_calls, addr)
end

function get_primitive_call(trace::GFTrace, addr::Pair)
    (first, rest) = addr
    subtrace::GFTrace = get(trace.subtraces, first)
    get_primitive_call(subtrace, rest)
end

function assoc_primitive_call(trace::GFTrace, addr, call::CallRecord)
    primitive_calls = assoc(trace.primitive_calls, addr, call)
    GFTrace(trace.call, primitive_calls, trace.subtraces, true)
end

function assoc_primitive_call(trace::GFTrace, addr::Pair, call::CallRecord)
    (first, rest) = addr
    local node::GFTrace
    local subtraces::PersistentHashMap{Any,Any}
    if !haskey(trace.subtraces, first)
        node = GFTrace()
        subtraces = assoc(trace.subtraces, first, node)
    else
        node = get(trace.subtraces, first)
        subtraces = trace.subtraces
    end
    node = assoc_primitive_call(node, rest, call)
    subtraces = assoc(subtraces, first, node)
    GFTrace(trace.call, trace.primitive_calls, subtraces, true)
end

function get_subtrace(trace::GFTrace, addr)
    get(trace.subtraces, addr)
end

function get_subtrace(trace::GFTrace, addr::Pair)
    (first, rest) = addr
    subtrace::GFTrace = get(trace.subtraces, first)
    get_subtrace(subtrace, rest)
end

function assoc_subtrace(trace::GFTrace, addr, subtrace)
    subtraces = assoc(trace.subtraces, addr, subtrace)
    GFTrace(trace.call, trace.primitive_calls, subtraces,
            has_choices(subtrace) || trace.has_choices)
end

function assoc_subtrace(trace::GFTrace, addr::Pair, subtrace)
    (first, rest) = addr
    local node::GFTrace
    local subtraces::PersistentHashMap{Any,Any}
    if !haskey(trace.subtraces, first)
        node = GFTrace()
        subtraces = assoc(trace.subtraces, first, node)
    else
        node = get(trace.subtraces, first)
        subtraces = trace.subtraces
    end
    node = assoc_subtrace(node, rest, subtrace)
    subtraces = assoc(subtraces, first, node)
    GFTrace(trace.call, trace.primitive_calls, subtraces,
            has_choices(subtrace) || trace.has_choices)
end

##################################################
# assignment wrapping generative function traces #
##################################################

struct GFTraceAssignment <: Assignment
    trace::GFTrace
end

get_assignment(trace::GFTrace) = GFTraceAssignment(trace)
get_address_schema(::Type{GFTraceAssignment}) = DynamicAddressSchema()
Base.isempty(assignment::GFTraceAssignment) = !has_choices(assignment.trace)

function has_internal_node(assignment::GFTraceAssignment, addr::Pair)
    (first, rest) = addr
    if haskey(assignment.trace.subtraces, first)
        subtrace = assignment.trace.subtraces[first]
        has_internal_node(get_assignment(subtrace), rest)
    else
        false
    end
end

function has_internal_node(assignment::GFTraceAssignment, addr)
    haskey(assignment.trace.subtraces, addr) && has_choices(assignment.trace.subtraces[addr])
end

function get_internal_node(assignment::GFTraceAssignment, addr::Pair)
    (first, rest) = addr
    subtrace = assignment.trace.subtraces[first]
    if !has_choices(subtrace)
        throw(KeyError(addr))
    end
    get_internal_node(get_assignment(subtrace), rest)
end

function get_internal_node(assignment::GFTraceAssignment, addr)
    subtrace = assignment.trace.subtraces[addr]
    if !has_choices(subtrace)
        throw(KeyError(addr))
    end
    get_assignment(subtrace)
end

function has_leaf_node(assignment::GFTraceAssignment, addr::Pair)
    (first, rest) = addr
    if !haskey(assignment.trace.subtraces, first)
        return false
    end
    sub_assignment = get_assignment(assignment.trace.subtraces[first])
    has_leaf_node(sub_assignment, rest)
end

function has_leaf_node(assignment::GFTraceAssignment, addr)
    haskey(assignment.trace.primitive_calls, addr)
end

function get_leaf_node(assignment::GFTraceAssignment, addr::Pair)
    (first, rest) = addr
    sub_assignment = get_assignment(assignment.trace.subtraces[first])
    get_leaf_node(sub_assignment, rest)
end

function get_leaf_node(assignment::GFTraceAssignment, addr)
    assignment.trace.primitive_calls[addr].retval
end

function get_leaf_nodes(assignment::GFTraceAssignment)
    ((key, call.retval) for (key, call) in assignment.trace.primitive_calls)
end

function get_internal_nodes(assignment::GFTraceAssignment)
    ((key, get_assignment(subtrace)) for (key, subtrace) in assignment.trace.subtraces
     if has_choices(subtrace))
end

max_stars_repo_name: BoltzmannBrain/Gen

max_issues_repo_name: BoltzmannBrain/Gen

max_forks_repo_name: BoltzmannBrain/Gen