text:
## Various derivative-free iterative algorithms to find a zero

## Use iterator for solving.
## allow for solving along the lines of:
## ```
## for xn in itr
##    println(xn)
## end
## ```
type ZeroType{T, S}
    f
    fp
    fpp

    update

    xn::Vector{T}
    fxn::Vector{S}

    xtol
    xtolrel
    ftol
    state

    cnt
    maxcnt
    fevals
    maxfevals
end
incfn(o::ZeroType,i::Int=1) = (o.fevals = o.fevals + i)
inccnt(o::ZeroType, i::Int=1) = (o.cnt = o.cnt + i)

function Base.start{T}(o::ZeroType{T})
    (o.xn[end], o.fxn[end])
end

function Base.next{T,S}(o::ZeroType{T,S}, state)
    o.update(o)
    o.cnt = o.cnt + 1
    (o.xn[end], o.fxn[end]), (o.xn[end], o.fxn[end])
end

## done function centralizes the stopping rules
function Base.done(o::ZeroType, state)
    if o.state == :converged
        ## check for near convergence
        (norm(o.fxn[end]) <= sqrt(o.ftol)) && return true
        throw(ConvergenceFailed("Algorithm stopped with xn=$(o.xn[end]), f(xn)=$(o.fxn[end])"))
    end

    o.cnt > o.maxcnt && throw(ConvergenceFailed("Too many steps taken"))
    o.fevals > o.maxfevals  && throw(ConvergenceFailed("Too many function calls taken"))
    isnan(o.fxn[end]) && throw(ConvergenceFailed("NaN produced by algorithm"))
    isinf(o.xn[end]) && throw(ConvergenceFailed("Algorithm escaped to oo"))

    # return turn if f(xn) \approx 0 or xn+1 - xn \approx 0
    lambda = max(1, abs(o.xn[end]))
    ftol = lambda * o.ftol
    xtol = o.xtol +  lambda * o.xtolrel

    norm(o.fxn[end]) < ftol && return true
    if length(o.xn) > 1
        if (norm(o.xn[end] - o.xn[end-1])) <= xtol && (norm(o.fxn[end]) < sqrt(ftol))
            return true
        end
    end
    false
end

## printing utility, borrowed idea of showing change from previous from @stevengj. Only works at REPL
function printdiff(x1, x0, color=:red)
    s0, s1 = string(x0), string(x1)
    flag = true
    for i in 1:length(s1)
        if i <= length(s0) && s0[i:i] == s1[i:i]
            print_with_color(:red, s1[i:i])
        else
            print(s1[i:end])
            break
        end
    end
    print("\n")
end

# show output
function verbose_output(out::ZeroType)
    println("Steps = $(out.cnt), function calls = $(out.fevals); steps:")
    xs = copy(out.xn)
    x0 = shift!(xs)
    println(x0)
    while length(xs) > 0
        x1 = shift!(xs)
        printdiff(x1, x0)
        x0 = x1
    end
    println("")
end

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

## issue with approx derivative
isissue(x) = (x == 0.0) || isnan(x) || isinf(x)

"""
heuristic to get a decent first step with Steffensen steps
"""
function steff_step(x, fx)
    thresh =  max(1, sqrt(abs(x/fx))) * 1e-6
    abs(fx) <= thresh ? fx : sign(fx) * thresh
end

## Different functions for approximating f'(xn)
## return fpxn and whether it is an issue

## use f[a,b] to approximate f'(x)
function _fbracket(a, b, fa, fb)
    out = (fb - fa) / (b - a)
    out, isissue(out)
end

## use f[y,z] - f[x,y] + f[x,z] to approximate
function _fbracket_diff(a,b,c, fa, fb, fc)
    x1, state = _fbracket(b, c, fb,  fc)
    x2, state = _fbracket(a, b, fa,  fb)
    x3, state = _fbracket(a, c, fa,  fc)
    out = x1 - x2 + x3
    out, isissue(out)
end

## use f[a,b] * f[a,c] / f[b,c]
function _fbracket_ratio(a, b, c, fa, fb, fc)
    x1,_ = _fbracket(b, c, fb, fc)
    x2,_ = _fbracket(a, b, fa, fb)
    x3,_ = _fbracket(a, c, fa, fc)
    out = (x2 * x3) / x1
    out, isissue(out)
end

##################################################
## Iterators

# iterator for secant function
function secant_itr(f, x0::Real, x1::Real; xtol=4*eps(), xtolrel=4*eps(), ftol=4*eps(), maxsteps=100, maxfnevals=100)
    update = (o) -> begin
        xn_1, xn = o.xn[(end-1):end]
        fxn_1, fxn = o.fxn[(end-1):end]

        fp, issue = _fbracket(xn, xn_1, fxn, fxn_1)

        xn1 = xn - fxn  / fp
        fxn1 = o.f(xn1)
        incfn(o)

        push!(o.xn, xn1)
        push!(o.fxn, fxn1)

    end
    x, fx = promote(float(x1), f(float(x1)))
    out = ZeroType(f, nothing, nothing, update,
                   [float(x0), x], [f(float(x0)), fx],
                   xtol, xtolrel, ftol, :not_converged,
                   0, maxsteps, 2, maxfnevals)

    out
end

function steffensen_itr{T<:AbstractFloat}(f, x0::T;
                   xtol=4*eps(), xtolrel=4*eps(), ftol=4*eps(),
                   maxsteps=100, maxfnevals=100)
    update = (o) -> begin
        xn = o.xn[end]
        fxn = o.fxn[end]

        wn = xn + steff_step(xn, fxn)
        fwn = f(wn)
        incfn(o)

        fp, issue = _fbracket(xn, wn, fxn, fwn)

        if issue
            o.state = :converged
            return
        end

        xn1 = xn - fxn / fp

        fxn1 = o.f(xn1)
        incfn(o)

        push!(o.xn, xn1)
        push!(o.fxn, fxn1)

    end

    x = float(x0)
    x, fx = promote(x, f(x))
    out = ZeroType(f, nothing, nothing, update, [x], [fx],
                   xtol, xtolrel, ftol, :not_converged,
                   0, maxsteps, 1, maxfnevals)

    out
end

## http://www.naturalspublishing.com/files/published/ahb21733nf19a5.pdf
## A New Fifth Order Derivative Free Newton-Type Method for Solving Nonlinear Equations
## Manoj Kumar, Akhilesh Kumar Singh, and Akanksha Srivastava
## Appl. Math. Inf. Sci. 9, No. 3, 1507-1513 (2015)
function kss5_itr(f, x0::Real;
             xtol=4*eps(), xtolrel=4*eps(), ftol=4*eps(),
             maxsteps=100, maxfnevals=100)

    update = o -> begin
        xn = o.xn[end]
        fxn = o.fxn[end]

        wn = xn + steff_step(xn, fxn)
        fwn = o.f(wn)
        incfn(o)

        fp, issue = _fbracket(xn, wn, fxn, fwn)
        if issue
            o.state = :converged
            return
        end

        yn = xn - fxn / fp
        fyn = o.f(yn)
        incfn(o)

        zn = xn - (fxn + fyn) / fp ## not a step in thukral algorithms
        fzn = o.f(zn)
        incfn(o)
        fp, issue = _fbracket_ratio(yn, xn, wn, fyn, fxn, fwn)

        if issue || isinf(fzn)
            push!(o.xn, yn)
            push!(o.fxn, fyn)
            o.state = :converged
            return
        end

        xn1 = zn  - fzn  / fp
        fxn1 = o.f(xn1); incfn(o)

        push!(o.xn, xn1)
        push!(o.fxn, fxn1)
    end

    x, fx = promote(float(x0), f(float(x0)))
    out = ZeroType(f, nothing, nothing, update, [x], [fx],
                   xtol, xtolrel, ftol,:not_converged,
                   0, maxsteps, 1, maxfnevals)

    out
end

## http://www.hindawi.com/journals/ijmms/2012/493456/
## Rajinder Thukral
## very fast (8th order) derivative free iterative root finder.
function thukral8_itr(f, x0::Real;
             xtol=4*eps(), xtolrel=4*eps(), ftol=4*eps(),
             maxsteps=100, maxfnevals=100)

    update = o -> begin
        xn = o.xn[end]
        fxn = o.fxn[end]

        wn = xn + steff_step(xn, fxn)
        fwn = o.f(wn)
        incfn(o)

        fp, issue = _fbracket(xn, wn, fxn, fwn)
        issue && return (xn, true)

        if issue
            o.state = :converged
            return
        end

        yn = xn - fxn / fp
        fyn = o.f(yn)
        incfn(o)

        fp, issue = _fbracket(yn, xn, fyn, fxn)
        if issue
            push!(o.xn, yn); push!(o.fxn, fyn)
            o.state = :converged
            return
         end

        phi = (1 + fyn / fwn)           # pick one of options
        zn =  yn - phi * fyn / fp
        fzn = o.f(zn)
        incfn(o)

        fp, issue =  _fbracket_diff(xn, yn, zn, fxn, fyn, fzn)
        if issue
            push!(o.xn, zn)
            push!(o.fxn, fzn)
            o.state = :converged
            return
         end

        w = 1 / (1 - fzn/fwn)
        xi = (1 - 2fyn*fyn*fyn / (fwn * fwn * fxn))

        xn1 = zn - w * xi * fzn / fp
        fxn1 = o.f(xn1)
        incfn(o)

        push!(o.xn, xn1)
        push!(o.fxn, fxn1)
    end

    x, fx = promote(float(x0), f(float(x0)))
    out = ZeroType(f, nothing, nothing, update, [x], [fx],
                   xtol, xtolrel, ftol,:not_converged,
                   0, maxsteps, 1, maxfnevals)

    out
end

## 16th order, derivative free root finding algorithm
## http://article.sapub.org/10.5923.j.ajcam.20120203.08.html
## American Journal of Computational and Applied Mathematics
## p-ISSN: 2165-8935    e-ISSN: 2165-8943
## 2012;  2(3): 112-118
## doi: 10.5923/j.ajcam.20120203.08
## New Sixteenth-Order Derivative-Free Methods for Solving Nonlinear Equations
## R. Thukral
## Research Centre, 39 Deanswood Hill, Leeds, West Yorkshire, LS17 5JS, England
## from p 114 (17)
function thukral16_itr(f, x0::Real;
             xtol=4*eps(), xtolrel=4*eps(), ftol=4*eps(),
             maxsteps=100, maxfnevals=100)

    update = o -> begin
        xn = o.xn[end]
        fxn = o.fxn[end]

        wn = xn + steff_step(xn, fxn)
        fwn = o.f(wn)
        incfn(o)

        fp, issue = _fbracket(xn, wn, fxn, fwn)
        if issue
            o.state = :converged
            return
         end

        yn = xn - fxn / fp
        fyn = o.f(yn)
        incfn(o)

        fp, issue = _fbracket(xn, yn, fxn, fyn)
        phi = _fbracket(xn, wn, fxn, fwn)[1] / _fbracket(yn, wn, fyn, fwn)[1]
        if issue
            push!(o.xn, yn); push!(o.fxn, fyn)
            o.state = :converged
            return
         end

        zn = yn - phi * fyn / fp
        fzn = o.f(zn)
        incfn(o)

        fp, issue = _fbracket_diff(xn, yn, zn, fxn, fyn, fzn)
        u2, u3, u4 = fzn/fwn, fyn/fxn, fyn/fwn
        eta = 1 / (1 + 2*u3*u4^2) / (1 - u2)
        if issue
            push!(o.xn, zn); push!(o.fxn, fzn)
            o.state = :converged
            return
        end

        an = zn - eta * fzn / fp
        fan = o.f(an)
        incfn(o)

        fp, issue = _fbracket_ratio(an, yn, zn, fan, fyn, fzn)
        u1, u5, u6 = fzn/fxn, fan/fxn, fan/fwn
        sigma =  1 + u1*u2 - u1*u3*u4^2 + u5 + u6 + u1^2*u4 +
                 u2^2*u3 + 3*u1*u4^2*(u3^2 - u4^2)/_fbracket(xn,yn, fxn, fyn)[1]

        if issue
            push!(o.xn, an); push!(o.fxn, fan)
            o.state = :converged
            return
        end

        xn1 = an - sigma * fan / fp
        fxn1 = o.f(xn1)
        incfn(o)

        push!(o.xn, xn1)
        push!(o.fxn, fxn1)
    end

    x, fx = promote(float(x0), f(float(x0)))
    out = ZeroType(f, nothing, nothing, update, [x], [fx],
                   xtol, xtolrel, ftol,:not_converged,
                   0, maxsteps, 1, maxfnevals)

    out
end

"""
Main interface for derivative free methods

* `f` a scalar function `f:R -> R` or callable object. Methods try to find a solution to `f(x) = 0`.

* `x0` initial guess for zero. Iterative methods need a reasonable
initial starting point.

* `ftol`. Stop iterating when |f(xn)| <= max(1, |xn|) * ftol.

* `xtol`. Stop iterating when |xn+1 - xn| <= xtol + abs(1, |xn|) * xtolrel. Checks that f(xn) is reasonably close.

* `xtolrel`. Stop iterating when |xn+1 - xn| <= xtol + abs(1, |xn|) * xtolrel. Checks that f(xn) is reasonably close.

* `maxeval`. Stop iterating if more than this many steps, throw error.

* `maxfneval`. Stop iterating if more than this many function calls, throw error.

* `order`. One of 0, 1, 2, 5, 8, or 16. Specifies which algorithm to
use.

- Default is 0 for the slower, more robust, SOLVE function.
- order 1 is a secant method
- order 2 is a Steffensen method
- order 5 uses a method of Kumar, Kumar Singh, and Srivastava
- order 8 From Thukral. Seems a bit more robust than the secant method and Steffensen method
- order 16 A higher-order method due to Thurkal. It may be faster when used with `Big` values.

* `verbose`. If true, will print number of iterations, function calls, and each step taken

* `kwargs...` passed on.

The `SOLVE` method has different stopping criteria.

The file test/test_fzero2 generates some timed comparisons

Because these are derivative free, they can be used with functions
defined by automatic differentiation
e.g., find critical point of f(x) = x^2
```
fzero(D(x -> x^2), 1)
```
"""

function derivative_free{T <: AbstractFloat}(f, x0::T;
                         ftol::Real = 10.0 * eps(x0),
                         xtol::Real =  4.0 * eps(x0),
                         xtolrel::Real = eps(x0),
                         maxeval::Int = 30,
                         verbose::Bool=false,
                         order::Int=0,  # 0, 1, 2, 5, 8 or 16
                         kwargs...      # maxfnevals,
                         )

    order == 0 && return(SOLVE(f, x0; ftol=ftol, maxeval=maxeval, verbose=verbose, kwargs...))
    _derivative_free(f, x0; order=order, ftol=ftol, xtol=xtol, xtolrel=xtolrel, maxeval=maxeval,
                     verbose=verbose, kwargs...)

end
function _derivative_free{T <: AbstractFloat}(f, x0::T;
                         ftol::Real = 10.0 * eps(x0),
                         xtol::Real =  4.0 * eps(x0),
                         xtolrel::Real = eps(x0),
                         maxeval::Int = 30,
                         verbose::Bool=false,
                         order::Int=0,  # 0, 1, 2, 5, 8 or 16
                         kwargs...      # maxfnevals, possible beta to control steffensen step
                         )
    if order == 16
        o = thukral16_itr(f, x0; ftol=ftol, xtol=xtol, xtolrel=xtolrel, maxsteps=maxeval, kwargs...)
    elseif order == 8
        o = thukral8_itr(f, x0; ftol=ftol, xtol=xtol, xtolrel=xtolrel, maxsteps=maxeval, kwargs...)
    elseif order == 5
        o = kss5_itr(f, x0; ftol=ftol, xtol=xtol, xtolrel=xtolrel, maxsteps=maxeval,  kwargs...)
    elseif order == 2
        o = steffensen_itr(f, x0; ftol=ftol, xtol=xtol, xtolrel=xtolrel, maxsteps=maxeval,  kwargs...)
    elseif order == 1
        x1 = x0 + steff_step(x0, f(x0))
        o = secant_itr(f, float(x1), x0; ftol=ftol, xtol=xtol, xtolrel=xtolrel, maxsteps=maxeval,  kwargs...)
    else
        throw(ArgumentError())
    end

    if done(o, start(o))
        verbose && println("Done before we started...")
        return(x0)::T
    end

    val = x0
    try
        for (x,fx) in o
            nothing
        end
        verbose && verbose_output(o)
        o.xn[end]
    catch err
        verbose && verbose_output(o)
        rethrow(err)
    end

end

"""

Implementation of secant method: `x_n1 = x_n - f(x_n) * f(x_n)/ (f(x_n) - f(x_{n-1}))`

Arguments:

* `f::Function` -- function to find zero of

* `x0::Real` -- initial guess is [x0, x1]

* `x1::Real` -- initial guess is [x0, x1]

Keyword arguments:

* `ftol`. Stop iterating when |f(xn)| <= max(1, |xn|) * ftol.

* `xtol`. Stop iterating when |xn+1 - xn| <= xtol + max(1, |xn|) * xtolrel

* `xtolrel`. Stop iterating when |xn+1 - xn| <= xtol + max(1, |xn|) * xtolrel

* `maxeval`. Stop iterating if more than this many steps, throw error.

* `maxfneval`. Stop iterating if more than this many function calls, throw error.

* `verbose::Bool=false` Set to `true` to see trace.

"""
function secant_method(f, x0::Real, x1::Real;
                       xtol=4*eps(), xtolrel=4*eps(), ftol=4eps(),
                       maxsteps::Int=100, maxfnevals=100,
                       verbose::Bool=false)

    x_0, x_1 = float(x0), float(x1)
    o = secant_itr(f, x_0, x_1;
                   xtol=xtol, xtolrel=xtolrel, ftol=ftol,
                   maxsteps=maxsteps, maxfnevals=maxfnevals)

    for x in o
        nothing
    end

    verbose && verbose_output(o)
    o.xn[end]
end

function steffensen_method{T<:AbstractFloat}(f, x0::T;
                       xtol=4*eps(), xtolrel=4*eps(), ftol=4eps(),
                       maxsteps::Int=100, maxfnevals=100,
                       verbose::Bool=false)
    o = steffensen_itr(f, x0;
                  xtol=xtol, xtolrel=xtolrel, ftol=ftol,
                  maxsteps=maxsteps, maxfnevals=maxfnevals)

    for x in o
        nothing
    end
    verbose && verbose_output(o)
    o.xn[end]
end
steffensen(f, x0::Number, args...; kwargs...) = steffensen_method(f, float(x0), args...; kwargs...)


text:
@testset "Model" begin
    @testset "Optimiser attributes" begin
        m = CP.FlatZinc.Model()
        @test sprint(show, m) == "A FlatZinc (fzn) model"
    end

    @testset "Supported constraints" begin
        m = CP.FlatZinc.Model()

        @test MOI.supports_constraint(
            m,
            MOI.VariableIndex,
            MOI.LessThan{Int},
        )
        @test MOI.supports_constraint(
            m,
            MOI.VariableIndex,
            MOI.LessThan{Float64},
        )
        @test MOI.supports_constraint(
            m,
            MOI.VariableIndex,
            CP.Strictly{MOI.LessThan{Float64}},
        )
        @test MOI.supports_constraint(m, MOI.VariableIndex, CP.Domain{Int})
        @test MOI.supports_constraint(
            m,
            MOI.VariableIndex,
            MOI.Interval{Float64},
        )
        @test MOI.supports_constraint(
            m,
            MOI.VectorOfVariables,
            CP.Element{Int},
        )
        @test MOI.supports_constraint(
            m,
            MOI.VectorOfVariables,
            CP.Element{Bool},
        )
        @test MOI.supports_constraint(
            m,
            MOI.VectorOfVariables,
            CP.Element{Float64},
        )
        @test MOI.supports_constraint(
            m,
            MOI.VectorOfVariables,
            CP.MaximumAmong,
        )
        @test MOI.supports_constraint(
            m,
            MOI.VectorOfVariables,
            CP.MinimumAmong,
        )
        @test MOI.supports_constraint(
            m,
            MOI.ScalarAffineFunction{Int},
            MOI.EqualTo{Int},
        )
        @test MOI.supports_constraint(
            m,
            MOI.ScalarAffineFunction{Int},
            MOI.LessThan{Int},
        )
        @test MOI.supports_constraint(
            m,
            MOI.ScalarAffineFunction{Int},
            CP.DifferentFrom{Int},
        )
        @test MOI.supports_constraint(
            m,
            MOI.ScalarAffineFunction{Float64},
            MOI.EqualTo{Float64},
        )
        @test MOI.supports_constraint(
            m,
            MOI.ScalarAffineFunction{Float64},
            MOI.LessThan{Float64},
        )
        @test MOI.supports_constraint(
            m,
            MOI.ScalarAffineFunction{Float64},
            CP.Strictly{MOI.LessThan{Float64}},
        )
        @test MOI.supports_constraint(
            m,
            MOI.ScalarAffineFunction{Float64},
            CP.DifferentFrom{Float64},
        )
    end

    @testset "Supported constrained variable with $(S)" for S in [
        MOI.EqualTo{Float64},
        MOI.LessThan{Float64},
        MOI.Interval{Float64},
        MOI.EqualTo{Int},
        MOI.LessThan{Int},
        MOI.Interval{Int},
        MOI.EqualTo{Bool},
        MOI.ZeroOne,
        MOI.Integer,
    ]
        m = CP.FlatZinc.Model()
        @test MOI.supports_add_constrained_variable(m, S)
        @test MOI.supports_add_constrained_variables(m, S)
    end
end

ext: jl


text:
using SIIPExamples
pkgpath = dirname(dirname(pathof(SIIPExamples)))
include(joinpath(pkgpath, "test", "3_PowerSimulations_examples", "01_operations_problems.jl"))

sys_RT = System(rawsys; time_series_resolution = Dates.Minute(5))
transform_single_time_series!(sys_RT, 12, Hour(1))

devices = Dict(
    :Generators => DeviceModel(ThermalStandard, ThermalStandardUnitCommitment),
    :Ren => DeviceModel(RenewableDispatch, RenewableFullDispatch),
    :Loads => DeviceModel(PowerLoad, StaticPowerLoad),
    :HydroROR => DeviceModel(HydroDispatch, FixedOutput),
    :Hydro => DeviceModel(HydroEnergyReservoir, HydroDispatchRunOfRiver),
    :RenFx => DeviceModel(RenewableFix, FixedOutput),
)
template_uc = template_unit_commitment(devices = devices)

devices = Dict(
    :Generators => DeviceModel(ThermalStandard, ThermalDispatch),
    :Ren => DeviceModel(RenewableDispatch, RenewableFullDispatch),
    :Loads => DeviceModel(PowerLoad, StaticPowerLoad),
    :HydroROR => DeviceModel(HydroDispatch, FixedOutput),
    :Hydro => DeviceModel(HydroEnergyReservoir, HydroDispatchRunOfRiver),
    :RenFx => DeviceModel(RenewableFix, FixedOutput),
)
template_ed = template_economic_dispatch(devices = devices)

stages_definition = Dict(
    "UC" => Stage(GenericOpProblem, template_uc, sys, solver),
    "ED" => Stage(
        GenericOpProblem,
        template_ed,
        sys_RT,
        solver,
        balance_slack_variables = true,
    ),
)

feedforward_chronologies = Dict(("UC" => "ED") => Synchronize(periods = 24))

feedforward = Dict(
    ("ED", :devices, :Generators) => SemiContinuousFF(
        binary_source_stage = PSI.ON,
        affected_variables = [PSI.ACTIVE_POWER],
    ),
)

order = Dict(1 => "UC", 2 => "ED")
horizons = Dict("UC" => 24, "ED" => 12)
intervals = Dict("UC" => (Hour(24), Consecutive()), "ED" => (Hour(1), Consecutive()))

DA_RT_sequence = SimulationSequence(
    step_resolution = Hour(24),
    order = order,
    horizons = horizons,
    intervals = intervals,
    ini_cond_chronology = InterStageChronology(),
    feedforward_chronologies = feedforward_chronologies,
    feedforward = feedforward,
)

sim = Simulation(
    name = "rts-test",
    steps = 1,
    stages = stages_definition,
    stages_sequence = DA_RT_sequence,
    simulation_folder = rts_dir,
)

build!(sim)

sim_results = execute!(sim)

uc_results = load_simulation_results(sim_results, "UC");
ed_results = load_simulation_results(sim_results, "ED");

# This file was generated using Literate.jl, https://github.com/fredrikekre/Literate.jl


text:
using NaiveNASflux
import NaiveNASflux: weights, bias

@testset "InputShapeVertex" begin
    v = inputvertex("in", 3, FluxDense())
    @test layertype(v) == FluxDense()
    @test name(v) == "in"
    @test nout(v) == 3

    c = clone(v)
    @test layertype(c) == FluxDense()
    @test name(c) == "in"
    @test nout(c) == 3
end

@testset "Size mutations" begin

    @testset "Dense to Dense" begin
        inpt = inputvertex("in", 4)
        dl1 = Dense(4, 5)
        dl2 = Dense(5, 3)
        bias(dl1)[1:end] = 1:5
        bias(dl2)[1:end] = 1:3
        dense1 = mutable(dl1, inpt)
        dense2 = mutable(dl2, dense1)

        @test inputs(dense2) == [dense1]
        @test outputs(dense1) == [dense2]
        @test inputs(dense1) == [inpt]
        @test outputs(inpt) == [dense1]

        @test layer(dense1) == dl1
        @test layer(dense2) == dl2
        @test layertype(dense1) == FluxDense()
        @test layertype(dense2) == FluxDense()

        @test [nout(inpt)] == nin(dense1) == [4]
        @test [nout(dense1)] == nin(dense2) == [5]
        @test nout(dense2) == 3

        inds = [1, 2, 4, 5]
        Δnin(NaiveNASlib.OnlyFor(), dense2, inds)
        Δnout(NaiveNASlib.OnlyFor(), dense1, inds)

        @test [nout(inpt)] == nin(dense1) == [4]
        @test [nout(dense1)] == nin(dense2) == [4]
        @test nout(dense2) == 3

        W2exp, b2exp = weights(dl2)[:, inds], bias(dl2)
        W1exp, b1exp = weights(dl1)[inds, :], bias(dl1)[inds]
        apply_mutation.(flatten(dense2))

        @test size(CompGraph([inpt], [dense2])(collect(Float32, 1:nout(inpt)))) == (3,)

        assertlayer(layer(dense2), W2exp, b2exp)
        assertlayer(layer(dense1), W1exp, b1exp)
    end

    @testset "Invariant parametric layer" begin
        inpt = inputvertex("in", 3)
        cv = mutable(Conv((1,1), nout(inpt) => 4), inpt)
        bv = mutable(BatchNorm(nout(cv)), cv)

        @test nin(bv) == [nout(cv)] == [4]

        Δnin(bv, -1)
        apply_mutation.(flatten(bv))

        @test nin(bv) == [nout(cv)] == [3]
    end

    @testset "Invariant non-parametric layer" begin
        inpt = inputvertex("in", 3)
        cv = mutable(Conv((1,1), nout(inpt) => 4), inpt)
        bv = mutable(MeanPool((2,2)), cv)

        @test nin(bv) == [nout(cv)] == [4]

        Δnin(bv, -1)
        apply_mutation.(flatten(bv))

        @test nin(bv) == [nout(cv)] == [3]
    end

    @testset "DepthwiseConv" begin
        inpt = inputvertex("in", 4)
        dc1 = mutable("dc1", DepthwiseConv((2,2), nout(inpt) => 2 * nout(inpt)), inpt)
        dc2 = mutable("dc2", DepthwiseConv((2,2), nout(dc1) => nout(dc1)), dc1)

        @test_logs (:warn, r"Could not change nout of") Δnout(dc1, 2)
        @test [nout(dc1)] == nin(dc2) == [nout(dc2)] == [12]

        @test_logs (:warn, r"Could not change nout of") Δnout(dc1, -2)
        @test [nout(dc1)] == nin(dc2) == [nout(dc2)] == [8]

        @test_logs (:warn, r"Could not change nout of") Δsize(ΔNoutExact(dc2, -2), all_in_graph(dc2))
        @test [nout(dc1)] == nin(dc2) == [nout(dc2)] == [4]

        @test_logs (:warn, r"Could not change nout of") Δsize(ΔNoutExact(dc2, 2), all_in_graph(dc2))
        @test [2nout(dc1)] == 2 .* nin(dc2) == [nout(dc2)] == [8]
    end

    @testset "Concatenate activations" begin

        function testgraph(layerfun, nin1, nin2)
            vfun = (v, s) -> mutable(layerfun(nout(v), s), v)
            return testgraph_vfun(vfun, nin1, nin2)
        end

        function testgraph_vfun(vertexfun, nin1::Integer, nin2::Integer)
            in1 = inputvertex("in1", nin1)
            in2 = inputvertex("in2", nin2)
            return testgraph_vfun(vertexfun, in1, in2)
        end

        function testgraph_vfun(vertexfun, in1, in2)
            l1 = vertexfun(in1, 3)
            l2 = vertexfun(in2, 7)

            joined = concat(l1, l2)
            l3 = vertexfun(joined, 9)

            @test nout(joined) == nout(l1) + nout(l2) == 10
            return CompGraph([in1, in2], l3)
        end

        @testset "Concatenate Dense" begin
            nin1 = 2
            nin2 = 5
            @test size(testgraph(Dense, nin1, nin2)(ones(nin1), ones(nin2))) == (9,)
        end

        @testset "Concatenate $rnntype" for rnntype in (RNN, GRU, LSTM)
            nin1 = 2
            nin2 = 5
            indata1 = reshape(collect(Float32, 1:nin1*4), nin1, 4)
            indata2 = reshape(collect(Float32, 1:nin2*4), nin2, 4)

            @test size(testgraph(rnntype, nin1, nin2)(indata1, indata2)) == (9,4)
        end

        @testset "Concatenate $convtype" for convtype in (Conv, ConvTranspose, CrossCor)
            nin1 = 2
            nin2 = 5
            indata1 = reshape(collect(Float32, 1:nin1*4*4), 4, 4, nin1, 1)
            indata2 = reshape(collect(Float32, 1:nin2*4*4), 4, 4, nin2, 1)

            convfun = (nin,nout) -> convtype((3,3), nin=>nout, pad = (1,1))
            @test size(testgraph(convfun, nin1, nin2)(indata1, indata2)) == (4,4,9,1)

        end

        @testset "Concatenate Pooled Conv" begin
            nin1 = 2
            nin2 = 5
            indata1 = reshape(collect(Float32, 1:nin1*4*4), 4, 4, nin1, 1)
            indata2 = reshape(collect(Float32, 1:nin2*4*4), 4, 4, nin2, 1)
            function vfun(v, s)
                cv = mutable(Conv((3,3), nout(v)=>s, pad = (1,1)), v)
                return mutable(MaxPool((3,3), pad=(1,1), stride=(1,1)), cv)
            end
            @test size(testgraph_vfun(vfun, nin1, nin2)(indata1, indata2)) == (4,4,9,1)
        end

        @testset "Concatenate InputShapeVertex" begin
            nin1 = 6
            nin2 = 4
            in1 = inputvertex("in1", nin1, FluxDense())
            in2 = inputvertex("in2", nin2, FluxDense())
            @test size(testgraph_vfun((v,s) -> v, in1, in2)(ones(nin1), ones(nin2))) == (10,)
        end

        @testset "Concatenate BatchNorm only" begin
            nin1 = 3
            nin2 = 7
            indata1 = reshape(collect(Float32, 1:nin1*4*4), 4, 4, nin1, 1)
            indata2 = reshape(collect(Float32, 1:nin2*4*4), 4, 4, nin2, 1)
            in1 = inputvertex("in1", nin1, FluxConv{2}())
            in2 = inputvertex("in2", nin2, FluxConv{2}())
            vfun(v,s) = mutable(BatchNorm(nout(v)), v)

            @test size(testgraph_vfun(vfun, in1, in2)(indata1, indata2)) == (4,4,10,1)
        end

        @testset "Concatentate dimension mismatch fail" begin
            d1 = mutable(Dense(2,3), inputvertex("in1", 2))
            c1 = mutable(Conv((3,3), 4=>5), inputvertex("in2", 4))
            r1 = mutable(RNN(6,7), inputvertex("in3", 6))
            @test_throws DimensionMismatch concat(d1, c1)
            @test_throws DimensionMismatch concat(r1, c1)
            @test_throws DimensionMismatch concat(d1, r1)
        end

        @testset "Concat with name" begin
            d1 = mutable(Dense(2,3), inputvertex("in1", 2))
            d2 = mutable(Dense(2,5), inputvertex("in1", 2))
            c = concat("c", d1, d2)

            @test name(c) == "c"
            @test nout(c) == 8
            @test nin(c) == [3, 5]
        end
    end

    @testset "Tricky structures" begin

        mutable struct Probe
            activation
        end
        function (p::Probe)(x)
            p.activation = x
            return x
        end

        function probe(in)
            p = Probe(nothing)
            return invariantvertex(NoParams(p), in), p
        end

        conv3x3(inpt::AbstractVertex, nch::Integer) = mutable(Conv((3,3), nout(inpt)=>nch, pad=(1,1)), inpt)
        batchnorm(inpt) = mutable(BatchNorm(nout(inpt)), inpt)
        mmaxpool(inpt) = mutable(MaxPool((2,2)), inpt)

        @testset "Residual Conv block" begin
            inpt = inputvertex("in", 3)
            conv1 = conv3x3(inpt, 5)
            pv1, p1 = probe(conv1)
            bn1 = batchnorm(pv1)
            conv2 = conv3x3(bn1, 5)
            pv2, p2 = probe(conv2)
            bn2 = batchnorm(pv2)
            add = bn2 + bn1
            mp = mmaxpool(add)
            out = conv3x3(mp, 3)

            graph = CompGraph([inpt], [out])

            # Test that the graph works to begin with
            indata = reshape(collect(Float32, 1:3*4*4), 4, 4, 3, 1)
            @test size(graph(indata)) == (2, 2, 3, 1)
            @test size(p1.activation) == (4, 4, 5, 1)
            @test size(p2.activation) == (4, 4, 5, 1)

            Δnin(out, -1)
            Δoutputs(out, v -> 1:nout_org(v))
            apply_mutation(graph)

            @test size(graph(indata)) == (2, 2, 3, 1)
            @test size(p1.activation) == (4, 4, 4, 1)
            @test size(p2.activation) == (4, 4, 4, 1)
        end

        @testset "Residual fork Conv block" begin
            inpt = inputvertex("in", 3)
            conv1 = conv3x3(inpt, 5)
            pv1, p1 = probe(conv1)
            bn1 = batchnorm(pv1)
            conv1a = conv3x3(bn1, 3)
            pv1a, p1a = probe(conv1a)
            bn1a = batchnorm(pv1a)
            conv1b = conv3x3(bn1, 2)
            pv1b, p1b = probe(conv1b)
            bn1b = batchnorm(pv1b)
            mv = concat(bn1a, bn1b)
            add = mv + bn1
            mp = mmaxpool(add)
            out = conv3x3(mp, 3)

            graph = CompGraph([inpt], [out])

            # Test that the graph works to begin with
            indata = reshape(collect(Float32, 1:3*4*4), 4, 4, 3, 1)
            @test size(graph(indata)) == (2, 2, 3, 1)
            @test size(p1.activation) == (4, 4, 5, 1)
            @test size(p1a.activation) == (4, 4, 3, 1)
            @test size(p1b.activation) == (4, 4, 2, 1)

            Δnin(out, -1)
            Δoutputs(out, v -> 1:nout_org(v))
            apply_mutation(graph)

            @test size(graph(indata)) == (2, 2, 3, 1)
            @test size(p1.activation) == (4, 4, 4, 1)
            @test size(p1a.activation) == (4, 4, 2, 1)
            @test size(p1b.activation) == (4, 4, 2, 1)
        end

        rnnvertex(inpt, outsize) = mutable("rnn", RNN(nout(inpt), outsize), inpt)
        densevertex(inpt, outsize) = mutable("dense", Dense(nout(inpt), outsize), inpt)

        @testset "RNN to Dense" begin
            inpt = inputvertex("in", 4)
            rnn = rnnvertex(inpt, 5)
            pv, p = probe(rnn)
            dnn = densevertex(pv, 3)

            graph = CompGraph([inpt], [dnn])

            indata = [collect(Float32, 1:nin(rnn)[]) for i =1:10]
            @test size(hcat(graph.(indata)...)) == (3,10)
            @test size(p.activation) == (5,)

            Δnin(dnn, 1)
            Δoutputs(dnn, v -> 1:nout_org(v))
            apply_mutation(graph)

            @test size(hcat(graph.(indata)...)) == (3,10)
            @test size(p.activation) == (6,)

            Δnout(rnn, -2)
            Δoutputs(rnn, v -> 1:nout_org(v))
            apply_mutation(graph)

            @test size(hcat(graph.(indata)...)) == (3,10)
            @test size(p.activation) == (4,)
        end

    end
end

@testset "Trait functions" begin
    @test named("test")(SizeAbsorb()) == NamedTrait(SizeAbsorb(), "test")
    @test validated()(SizeAbsorb()) == SizeChangeValidation(SizeAbsorb())
    @test logged(level=Base.CoreLogging.Info, info=NameInfoStr())(SizeAbsorb()) == SizeChangeLogger(Base.CoreLogging.Info, NameInfoStr(), SizeAbsorb())
end

@testset "Flux children" begin
    import Flux:children
    inpt = inputvertex("in", 2, FluxDense())
    v1 = mutable(Dense(2, 3), inpt)
    v2 = mutable(Dense(3, 4), v1)
    g1 = CompGraph(inpt, v2)

    @test children(g1) == (inpt, v1, v2)

    pars1 = params(g1).order
    @test pars1[1] == layer(v1).W
    @test pars1[2] == layer(v1).b
    @test pars1[3] == layer(v2).W
    @test pars1[4] == layer(v2).b

    g2 = copy(g1)
    # Emulate Flux.gpu
    testfun(x) = x
    # Modify x.data or else x.tracked.f is set which prevents it from being recognized as a parameter when param is called
    testfun(x::TrackedArray) = param(2 .* x.data)
    mapleaves(testfun, g2)

    pars2 = params(g2).order
    @test pars2 == 2 .* pars1
end


text:
@block SimonDanisch [interaction] begin

    @cell "Interaction with Mouse" [interactive, scatter, lines, marker, record] begin
        using LinearAlgebra
        scene = Scene()
        r = range(0, stop = 3, length = 4)
        cam2d!(scene)
        time = Node(0.0)
        pos = lift(scene.events.mouseposition, time) do mpos, t
            map(range(0, stop = 2pi, length = 60)) do i
                circle = Point2f0(sin(i), cos(i))
                mouse = to_world(scene, Point2f0(mpos))
                secondary = (sin((i * 10f0) + t) * 0.09) * normalize(circle)
                (secondary .+ circle) .+ mouse
            end
        end
        scene = lines!(scene, pos, raw = true)
        p1 = scene[end]
        p2 = scatter!(
            scene,
            pos, markersize = 0.1f0,
            marker = :star5,
            color = p1[:color],
            raw = true
        )[end]
        scene
        display(scene)

        p1[:color] = RGBAf0(1, 0, 0, 0.1)
        p2[:marker] = 'π' #TODO fix this
        p2[:markersize] = 0.2

        # push a reasonable mouse position in case this is executed as part
        # of the documentation
        push!(scene.events.mouseposition, (250.0, 250.0))
        N = 50
        record(scene, @outputfile(mp4), range(0.01, stop = 0.4, length = N)) do i
            push!(scene.events.mouseposition, (250.0, 250.0))
            p2[:markersize] = i
            push!(time, time[] + 0.1)
        end
    end

    @cell "Mouse Picking" [scatter, heatmap, interactive] begin
        img = rand(100, 100)
        scene = Scene()
        heatmap!(scene, img, scale_plot = false)
        clicks = Node(Point2f0[(0,0)])
        on(scene.events.mousebuttons) do buttons
           if ispressed(scene, Mouse.left)
               pos = to_world(scene, Point2f0(scene.events.mouseposition[]))
               push!(clicks, push!(clicks[], pos))
           end
           return
        end
        scatter!(scene, clicks, color = :red, marker = '+', markersize = 10)
    end

    @cell "sliders" [scatter, slider] begin
        s1 = slider(LinRange(0.01, 1, 100), raw = true, camera = campixel!)
        s2 = slider(LinRange(-2pi, 2pi, 100), raw = true, camera = campixel!)
        data = lift(s2[end][:value]) do v
            map(LinRange(0, 2pi, 100)) do x
                4f0 .* Point2f0(sin(x) + (sin(x * v) .* 0.1), cos(x) + (cos(x * v) .* 0.1))
            end
        end
        p = scatter(data, markersize = s1[end][:value])
        hbox(p, vbox(s1, s2))
    end

    @cell "Time Series" [lines] begin

        function makie_live(;
                f0 = 1/2, fs = 100,
                winsec = 4, hopsec = 1/60
            )
            nwin = round(Integer, winsec*fs)
            nhop = round(Integer, hopsec*fs)
            # do the loop
            frame_start = -winsec
            frame_time = collect((0:(nwin-1)) * (1/fs))
            aframe = sin.(2*pi*f0.*(frame_start .+ frame_time))
            scene = lines(frame_start .+ frame_time, aframe)
            display(scene)
            lineplot = scene[end]
            fix = 0
            i = 0
            while isopen(scene) && i < 50 # dont loop longer than 50x for testing
                i += 1
                aframe .= sin.(2*pi*f0.*(frame_start .+ frame_time))
                # append!(aframe, randn(nhop)); deleteat!(aframe, 1:nhop)
                lineplot[1] = frame_start .+ frame_time
                lineplot[2] = aframe
                AbstractPlotting.update_limits!(scene)
                AbstractPlotting.update!(scene)
                sleep(hopsec)
                frame_start += hopsec
            end # while loop
            scene
        end
        makie_live()

    end
end


text:
@testset "flat_map" begin

# test of basic funtionality
sc = SparkContext(master="local")
nums = parallelize(sc, [1, 2, 3, 0, 4, 0])
rdd = flat_map(nums, it -> fill(it, it))

@test reduce(rdd, +) == 30

close(sc)

end


text:
# Estimate gradient from a small set of points.

function gradest3(x, y)
    @assert size(x,1)==size(y,1)
    (x1,x2)=size(x)
    global dx = similar(x, (x1-1, x2))
    global dy = similar(y, (x1-1, 1))
    ix = 0
    for i=2:x1
        dx[i-1,:] = x[1,:]-x[i,:]
        dy[i-1,:] = y[1,:]-y[i,:]
    end
    dx\dy
end

function gradest2(x, y)
    @assert size(x,1)==size(y,1)
    (x1,x2)=size(x)
    x3 = div(x1*(x1-1),2)
    global dx = similar(x, (x3, x2))
    global dy = similar(y, (x3, 1))
    ix = 0
    for i=1:x1
        for j=(i+1):x1
            ix += 1
            dx[ix,:] = x[i,:]-x[j,:]
            dy[ix,:] = y[i,:]-y[j,:]
        end
    end
    dx\dy
end

function gradest0()
    for t=1:10
        x=randn(2,7)
        x[:,1]=1
        b=randn(7,1)
        y=x*b
        bb=x\y
        @show t
        display(bb')
        display(x[1,:]-x[2,:])
        display(y[1]-y[2])
    end
end


text:
"""
Utility functions for `Paxos`
"""
module Utils

using Base
export closeAll, readAvailable, bounded, TimeoutException, compact, identity, groupBy, maxBy

function closeAll(closeables)
    for closeable in closeables
        try
            close(closeable)
        catch ex
            @error "Error during close" exception = (ex, stacktrace(catch_backtrace()))
        end
    end
end

"""
Reads all available messages on the `Channel`, returning an array
"""
function readAvailable(channel::Channel)
    messages = []
    while isready(channel)
        message = take!(channel)
        push!(messages, message)
    end
    messages
end

"""
An exception thrown if a block invoked with `bounded` did
not complete in the specified time.
"""
struct TimeoutException <: Exception
    timeout
end

"""
Execute the provided expression `fn` in a separate task,
waiting at most `timeout` seconds (using `sleep`). If
the task does not complete in the alloted time, then
abort the task and throw a `TimeoutException` to the
caller of this function.
"""
function bounded(fn, timeout)
    result = Channel(1)
    task = @async begin
        put!(result, fn())
    end
    timer = @async begin
        sleep(timeout)
        schedule(task, TimeoutException(timeout),error=true)
    end
    try
        wait(task)
    catch ex
        close(result)
        if isa(ex, TaskFailedException)
            rethrow(ex.task.exception)
        else
            rethrow()
        end
    end
    take!(result)
end

"""
Return an iterator that skips any `nothing` results from the underlying itetator
"""
function compact(itr)
    Iterators.filter(itr) do el
        el !== nothing
    end
end

"""
Simple identity function: returns its arguments
"""
function identity(value)
    value
end

"""
Return a `Dict` whose keys are the unique values of the supplied function `fn`,
and the values are those items in the original collection for which the supplied
function returned that result.
"""
function groupBy(keyFN::Function, collection, valueFn::Function=identity)
    if isempty(collection)
        Dict()
    else
        results = Dict()
        for item in collection
            key = keyFN(item)
            matches = get!(results, key, Vector())
            push!(matches, valueFn(item))
        end
        results
    end
end

"""
Return the item with the maximum value by comparising the value returned by the supplied
function `fn` against other values in the collection using `<`
"""
function maxBy(keyFn::Function, collection, valueFn=identity)
    maxItem = nothing
    maxValue = nothing
    for item in collection
        if maxItem === nothing
            maxItem = item
            maxValue = valueFn(item)
        end
        itemValue = keyFn(item)
        if(maxValue < itemValue)
            maxItem = item
            maxValue = itemValue
        end
    end
    maxItem
end

end


text:
using WGLMakie

x₀ = 0.0
Δx = 0.1

fps = 50
nframes = 500

xs = Observable([x₀])
ys = Observable([sin(x₀)])

fig, ax, = lines(ys, figure=(resolution=(500,400),))
limits!(ax, 0, nframes, -1.5, 1.5)

@time for i = 1:nframes
    x = xs[][i] + Δx
    xs[] = push!(xs[], x)
    ys[] = push!(ys[], sin(x))
    sleep(1 / fps)
end

fig

max_stars_repo_name: aarontrowbridge/GAN

max_issues_repo_name: aarontrowbridge/GAN

max_forks_repo_name: aarontrowbridge/GAN


text:
a=1:20
@test avg(a) == 10.5


text:
using SimpleDirectMediaLayer
const SDL2 = SimpleDirectMediaLayer
using Test

SDL2_pkg_dir = joinpath(@__DIR__, "..","..")

@testset "Init+Quit" begin
# Test that you can init and quit SDL multiple times correctly.
@test 0 == SDL2.Init(UInt32(SDL2.INIT_VIDEO))
SDL2.Quit()

@test 0 == SDL2.Init(UInt32(SDL2.INIT_VIDEO))
SDL2.Quit()
end

@testset "Window open+close" begin
@test 0 == SDL2.Init(UInt32(SDL2.INIT_VIDEO))
# Create window
win = SDL2.CreateWindow("Hello World!", Int32(100), Int32(100), Int32(300), Int32(400),
         UInt32(SDL2.WINDOW_SHOWN))
@test win != C_NULL

renderer = SDL2.CreateRenderer(win, Int32(-1),
             UInt32(SDL2.RENDERER_ACCELERATED | SDL2.RENDERER_PRESENTVSYNC))
if get(ENV, "HAS_JOSH_K_SEAL_OF_APPROVAL", "") == "true" && Sys.islinux()
    @test_broken renderer != C_NULL
else
    @test renderer != C_NULL
end

# Close window
SDL2.DestroyWindow(win)

# Create window again
win = SDL2.CreateWindow("Hello World!", Int32(100), Int32(100), Int32(300), Int32(400),
         UInt32(SDL2.WINDOW_SHOWN))
@test win != C_NULL
renderer = SDL2.CreateRenderer(win, Int32(-1),
             UInt32(SDL2.RENDERER_ACCELERATED | SDL2.RENDERER_PRESENTVSYNC))
if get(ENV, "HAS_JOSH_K_SEAL_OF_APPROVAL", "") == "true" && Sys.islinux()
    @test_broken renderer != C_NULL
else
    @test renderer != C_NULL
end

SDL2.Quit()
end

@testset "Window" begin
@test 0 == SDL2.Init(UInt32(SDL2.INIT_VIDEO))

win = SDL2.CreateWindow("Hello World!", Int32(100), Int32(100), Int32(300), Int32(400),
         UInt32(SDL2.WINDOW_SHOWN))

# Test get size
w,h = Int32[-1],Int32[-1]
SDL2.GetWindowSize(win, w,h)
@test 300 == w[]
@test 400 == h[]

e = SDL2.event()
@test typeof(e) <: SDL2.WindowEvent

# Test drawing
renderer = SDL2.CreateRenderer(win, Int32(-1),
             UInt32(SDL2.RENDERER_ACCELERATED | SDL2.RENDERER_PRESENTVSYNC))

rect = SDL2.Rect(1,1,50,50)
if get(ENV, "HAS_JOSH_K_SEAL_OF_APPROVAL", "") == "true" && Sys.islinux()
    @test_broken 0 == SDL2.RenderFillRect(renderer, pointer_from_objref(rect) )
else
    @test 0 == SDL2.RenderFillRect(renderer, pointer_from_objref(rect) )
end

@testset "Load/Save BMP" begin

img = SDL2.LoadBMP(joinpath(SDL2_pkg_dir, "assets/cat.bmp"))
@test img != C_NULL

img_tex = SDL2.CreateTextureFromSurface(renderer, img);
if get(ENV, "HAS_JOSH_K_SEAL_OF_APPROVAL", "") == "true" && Sys.islinux()
    @test_broken img_tex != C_NULL
else
    @test img_tex != C_NULL
end

src = SDL2.Rect(0,0,0,0)
if get(ENV, "HAS_JOSH_K_SEAL_OF_APPROVAL", "") == "true" && Sys.islinux()
    @test_broken 0 == SDL2.RenderCopy(renderer, img_tex, C_NULL, C_NULL) # Fill the renderer with img
else
    @test 0 == SDL2.RenderCopy(renderer, img_tex, C_NULL, C_NULL) # Fill the renderer with img
end
SDL2.RenderPresent(renderer)

# Save bmp
tmpFile = tempname()*".bmp"
@test 0 == SDL2.SaveBMP(img, tmpFile)

img2 = SDL2.LoadBMP(tmpFile)
@test img2 != C_NULL

# Compare the two images
if SDL2.MUSTLOCK(img)
    # Some surfaces must be locked before accessing pixels.
    @test 0 == SDL2.LockSurface(img)
    @test 0 == SDL2.LockSurface(img2)
end

img_surface1 = unsafe_load(img)
img_surface2 = unsafe_load(img2)

@test (img_surface1.w, img_surface1.h) == (img_surface2.w, img_surface2.h)

@test img_surface1.format != C_NULL
@test img_surface2.format != C_NULL

# Test pixels are equal
pxl_format1 = unsafe_load(img_surface1.format)
numPixelBytes = img_surface1.w * img_surface1.h * pxl_format1.BytesPerPixel
pixels1 = [unsafe_load(Ptr{UInt8}(img_surface1.pixels), i) for i in 1:numPixelBytes]
pixels2 = [unsafe_load(Ptr{UInt8}(img_surface2.pixels), i) for i in 1:numPixelBytes]

@test pixels1 == pixels2

if SDL2.MUSTLOCK(img)
    # Some surfaces must be locked/unlocked while accessing pixels.
    SDL2.UnlockSurface(img)
    SDL2.UnlockSurface(img2)
end

end  # @testset "BMP"

# Close window
SDL2.DestroyWindow(win)

SDL2.Quit()
end

max_stars_repo_name: Kyjor/julia-game

max_stars_repo_licenses: MIT Unlicense

max_issues_repo_name: Kyjor/julia-game

max_issues_repo_licenses: MIT Unlicense

max_forks_repo_name: Kyjor/julia-game

max_forks_repo_licenses: MIT Unlicense


text:
# ---
# title: 451. Sort Characters By Frequency
# id: problem451
# author: Indigo
# date: 2021-02-02
# difficulty: Medium
# categories: Hash Table, Heap
# link: <https://leetcode.com/problems/sort-characters-by-frequency/description/>
# hidden: true
# ---
#
# Given a string, sort it in decreasing order based on the frequency of
# characters.
#
# **Example 1:**
#
#
#
#     Input:
#     "tree"
#
#     Output:
#     "eert"
#
#     Explanation:
#     'e' appears twice while 'r' and 't' both appear once.
#     So 'e' must appear before both 'r' and 't'. Therefore "eetr" is also a valid answer.
#
#
# **Example 2:**
#
#
#
#     Input:
#     "cccaaa"
#
#     Output:
#     "cccaaa"
#
#     Explanation:
#     Both 'c' and 'a' appear three times, so "aaaccc" is also a valid answer.
#     Note that "cacaca" is incorrect, as the same characters must be together.
#
#
# **Example 3:**
#
#
#
#     Input:
#     "Aabb"
#
#     Output:
#     "bbAa"
#
#     Explanation:
#     "bbaA" is also a valid answer, but "Aabb" is incorrect.
#     Note that 'A' and 'a' are treated as two different characters.
#
#
#
## @lc code=start
using LeetCode

function frequency_sort(s::String)
    u_vec = codeunits(s)[:]
    freq = fill(0, 128)
    for n in codeunits(s)
        freq[n] += 1
    end
    sort!(u_vec; by=x -> (freq[x], x), rev = true)
    String(u_vec)
end
## @lc code=end


text:
# This file was generated by the Julia Swagger Code Generator
# Do not modify this file directly. Modify the swagger specification instead.

mutable struct ResourceRateCardInfo <: SwaggerModel
    Currency::Any # spec type: Union{ Nothing, String } # spec name: Currency
    Locale::Any # spec type: Union{ Nothing, String } # spec name: Locale
    IsTaxIncluded::Any # spec type: Union{ Nothing, Bool } # spec name: IsTaxIncluded
    OfferTerms::Any # spec type: Union{ Nothing, Vector{OfferTermInfo} } # spec name: OfferTerms
    Meters::Any # spec type: Union{ Nothing, Vector{MeterInfo} } # spec name: Meters

    function ResourceRateCardInfo(;Currency=nothing, Locale=nothing, IsTaxIncluded=nothing, OfferTerms=nothing, Meters=nothing)
        o = new()
        validate_property(ResourceRateCardInfo, Symbol("Currency"), Currency)
        setfield!(o, Symbol("Currency"), Currency)
        validate_property(ResourceRateCardInfo, Symbol("Locale"), Locale)
        setfield!(o, Symbol("Locale"), Locale)
        validate_property(ResourceRateCardInfo, Symbol("IsTaxIncluded"), IsTaxIncluded)
        setfield!(o, Symbol("IsTaxIncluded"), IsTaxIncluded)
        validate_property(ResourceRateCardInfo, Symbol("OfferTerms"), OfferTerms)
        setfield!(o, Symbol("OfferTerms"), OfferTerms)
        validate_property(ResourceRateCardInfo, Symbol("Meters"), Meters)
        setfield!(o, Symbol("Meters"), Meters)
        o
    end
end # type ResourceRateCardInfo

const _property_map_ResourceRateCardInfo = Dict{Symbol,Symbol}(Symbol("Currency")=>Symbol("Currency"), Symbol("Locale")=>Symbol("Locale"), Symbol("IsTaxIncluded")=>Symbol("IsTaxIncluded"), Symbol("OfferTerms")=>Symbol("OfferTerms"), Symbol("Meters")=>Symbol("Meters"))
const _property_types_ResourceRateCardInfo = Dict{Symbol,String}(Symbol("Currency")=>"String", Symbol("Locale")=>"String", Symbol("IsTaxIncluded")=>"Bool", Symbol("OfferTerms")=>"Vector{OfferTermInfo}", Symbol("Meters")=>"Vector{MeterInfo}")
Base.propertynames(::Type{ ResourceRateCardInfo }) = collect(keys(_property_map_ResourceRateCardInfo))
Swagger.property_type(::Type{ ResourceRateCardInfo }, name::Symbol) = Union{Nothing,eval(Base.Meta.parse(_property_types_ResourceRateCardInfo[name]))}
Swagger.field_name(::Type{ ResourceRateCardInfo }, property_name::Symbol) =  _property_map_ResourceRateCardInfo[property_name]

function check_required(o::ResourceRateCardInfo)
    true
end

function validate_property(::Type{ ResourceRateCardInfo }, name::Symbol, val)
end


text:
using Printf

function is_unit(O::ITensor,eps::Float64)::Bool
    s=inds(O)
    @ITensors.debug_check begin
        @assert(length(s)==2)
    end
    Id=delta(s[1],s[2])
    if hasqns(s)
        nm=0.0
        for b in eachnzblock(Id)
            isv=[s[i]=>b[i] for i in 1:length(b)]
            nm+=abs(O[isv...]-Id[isv...])
        end
        return nm<eps
    else
        return norm(O-Id)<eps
    end
end

function to_char(O::ITensor,eps::Float64)::Char
    c='0'
    if is_unit(O,eps)
        c='I'
    elseif norm(O)>eps
        c='S'
    end
    c
end
function to_char(O::Float64,eps::Float64)::Char
    c='0'
    if abs(O-1.0)<eps
        c='I'
    elseif abs(O)>eps
        c='S'
    end
    c
end

@doc """
pprint(W[,eps])

Show (pretty print) a schematic view of an operator-valued matrix.  In order to display any `ITensor`
as a matrix one must decide which index to use as the row index, and similarly for the column index.
`pprint` does this by inspecting the indices with "Site" tags to get the site number, and uses the
"Link" indices `l=\$n` as the column and `l=\$(n-1)` as the row. The output symbols I,0,S have the obvious
interpretations and S just means any (non zero) operator that is not I.  This function can obviously
be enhanced to recognize and display other symbols like X,Y,Z ... instead of just S.

# Arguments
- `W::ITensor` : Operator-valued matrix for display.
- `eps::Float64 = 1e-14` : operators inside W with norm(W[i,j])<eps are assumed to be zero.

# Examples
```julia
julia> pprint(W)
I 0 0 0 0
S 0 0 0 0
S 0 0 0 0
0 0 I 0 0
0 S 0 S I
```
"""
function pprint(W::ITensor,eps::Float64=default_eps)
    d,n,r,c=parse_links(W)
    pprint(r,W,c,eps)
end

function pprint(W::ITensor,r::Index,eps::Float64=default_eps)
    c=noncommoninds(W,r)
    @assert length(c)==1
    pprint(r,W,c[1],eps)
end

function pprint(r::Index,W::ITensor,c::Index,eps::Float64=default_eps)
    # @assert hasind(W,r) can't assume this becuse parse_links could return a Dw=1 dummy index.
    # @assert hasind(W,c)
    isl=filterinds(W,tags="Link")
    ord=order(W)
    if length(isl)==2
        for ir in  eachindval(r)
            for ic in  eachindval(c)
                if ord==4
                    Oij=slice(W,ir,ic)
                else
                    @assert ord==2
                    Oij=abs(W[ir,ic])
                end
                Base.print(to_char(Oij,eps))
                Base.print(" ")
            end
            Base.print("\n")
        end
    elseif length(isl)==1
        for i in  eachindval(isl[1])
        Oij=slice(W,i)
        Base.print(to_char(Oij,eps))
        Base.print(" ")
        end
        Base.print("\n")
    end
    Base.print("\n")
end

@doc """
    pprint(H[,eps])

Display the structure of an MPO, one line for each lattice site. For each site the table lists
- n   : lattice site number
- Dw1 : dimension of the row index (left link index)
- Dw2 : dimension of the column index (right link index)
- d   : dimension of the local Hilbert space
- Reg. Form  : U=upper, L=lower, D=diagonal, No = Not reg. form.
- Orth. Form : L=left, R=right, M=not orthogonal, B=both left and right.
- Tri. Form  : U=upper triangular, L=lower triangular, F=full (not triangular), D=diagonal, R=row, C=column

Be careful when reading the L symbols. L stands for Left in the Orth. column but is stands
for Lower in the other two columns.

    # Arguments
- `H::MPO` : MPO for display.
- `eps::Float64 = 1e-14` : operators inside H with norm(W[i,j])<eps are assumed to be zero.

# Examples
```julia
julia> using ITensors
julia> using ITensorMPOCompression
julia> N=10; #10 sites
julia> NNN=7; #Include up to 7th nearest neighbour interactions
julia> sites = siteinds("S=1/2",N);
julia> H=make_transIsing_MPO(sites,NNN);
julia> pprint(H)
n    Dw1  Dw2   d   Reg.  Orth.  Tri.
                    Form  Form   Form
 1    1   30    2    L      M     R
 2   30   30    2    L      M     L
 3   30   30    2    L      M     L
 4   30   30    2    L      M     L
 5   30   30    2    L      M     L
 6   30   30    2    L      M     L
 7   30   30    2    L      M     L
 8   30   30    2    L      M     L
 9   30   30    2    L      M     L
10   30    1    2    L      M     C

julia> orthogonalize!(H,orth=right)
julia> pprint(H)
n    Dw1  Dw2   d   Reg.  Orth.  Tri.
                    Form  Form   Form
 1    1    3    2    L      M     R
 2    3    4    2    L      R     L
 3    4    5    2    L      R     L
 4    5    6    2    L      R     L
 5    6    7    2    L      R     L
 6    7    6    2    L      R     L
 7    6    5    2    L      R     L
 8    5    4    2    L      R     L
 9    4    3    2    L      R     L
10    3    1    2    L      R     C
julia> truncate!(H,orth=left)
julia> pprint(H)
n    Dw1  Dw2   d   Reg.  Orth.  Tri.
                    Form  Form   Form
 1    1    3    2    L      L     R
 2    3    4    2    L      L     L
 3    4    5    2    L      L     F
 4    5    6    2    L      L     F
 5    6    7    2    L      L     F
 6    7    6    2    L      L     F
 7    6    5    2    L      L     F
 8    5    4    2    L      L     F
 9    4    3    2    L      L     L
10    3    1    2    L      M     C
```
"""
function pprint(H::MPO,eps::Float64=default_eps)
    N=length(H)
    println("  n    Dw1  Dw2   d   Reg.  Orth.  Tri.")
    println("                      Form  Form   Form")

    for n in 1:N
        Dw1,Dw2,d,l,u,lr,lt,ut=get_traits(H[n],eps)
        #println(" $n    $Dw1    $Dw2    $d    $l$u     $lr   $lt$ut")
        @printf "%4i %4i %4i %4i    %s%s     %s     %s%s\n" n Dw1 Dw2 d l u lr lt ut
    end
end

#
#  In general it is useful to pass kwargs from a parent function to a child
#  in addition in the parent function we may want add or replace and existing key=>val parse_links
#  for some reason this seems to be a nasty uphill battle, especially if kwargs is empty.
#
function add_or_replace(kwargs,key::Any,val::Any)
    if length(kwargs)>0
        kwargs=Dict{Symbol, Any}(kwargs) #this allows us to set an element
        kwargs[key]=val
    else
        kwargs=Dict{Symbol, Any}(key => val)
    end
    return kwargs
end

export pprint,is_unit,slice


text:
module SpatioTemporalBlobs

import StatsBase.median

mutable struct VariableStats
  Histogram::Array{Float64,2}
  ContrastingResponse::Array{Int, 1}
  Low_LATcenters::Array{Float64, 1}
  Low_LONcenters::Array{Float64, 1}
  Low_TIMEcenters::Array{Float64, 1}
  High_LATcenters::Array{Float64, 1}
  High_LONcenters::Array{Float64, 1}
  High_TIMEcenters::Array{Float64, 1}
  SpatialDistance::Array{Float64,1}
  TemporalDistance::Array{Float64,1}
  Balance::Array{Float32,1}
  PosSum::Array{Float32,1}
  NegSum::Array{Float32,1}
  Comp::Array{Float32,1}
  affectedLonLatTimeValue::Array{Dict{Tuple{Int64,Int64,Int64},Int64},1}
end

mutable struct CCstats
  Bools::Array{UInt8, 3}
  LabelCube::Array{Int, 3} # --> LabelCube
  ConnectionFootprint::Array{Bool, 3} # trues(3 ,3,3)
  numExtremes::Int
  minAffectedVoxels::Int
  labels::Array{Int, 1}
  perm::Array{Int,1}
  affectedVoxels::Array{Int, 1} # Volume [counts]
  affectedPixels::Array{Int, 1} # Area [counts]
  affectedTimeSteps::Array{Int, 1} # Duration [counts]
  affectedVolume::Array{Float64, 1} # Volume [days * km^2]
  affectedArea::Array{Float64, 1} # Area [km^2]
  affectedTime::Array{Float64, 1} # Duration [days]
  Tmins::Array{Float64, 1}
  Tmaxs::Array{Float64, 1}
  Tcenters::Array{Float64, 1}
  LATmins::Array{Float64, 1}
  LATmaxs::Array{Float64, 1}
  LATcenters::Array{Float64, 1}
  LONmins::Array{Float64, 1}
  LONmaxs::Array{Float64, 1}
  LONcenters::Array{Float64, 1}
  Latitudes::StepRangeLen{Float64,Base.TwicePrecision{Float64},Base.TwicePrecision{Float64}}
  Longitudes::StepRangeLen{Float64,Base.TwicePrecision{Float64},Base.TwicePrecision{Float64}}
  Times::Array{Float64, 1}
  CosLatitudes::Array{Float64, 1}
  AreaLatitudes::Array{Float64, 1}
  affectedLonLatIdxCounts::Array{Dict{Tuple{Int64,Int64},Int64},1}
  affectedLonLatTimeScoreValue::Array{Dict{Tuple{Int64,Int64,Int64},Float64},1}
  VarStats::Dict{String,VariableStats}
end

include("BasicStats.jl")

include("VariableStats.jl")

include("IO.jl")

end # module


text:
struct Clock
    hour::Int
    minute::Int

    function Clock(hour, minute)
        if minute < 0
            hour = hour + minute ÷ 60 - 1
            minute = 60 + (minute % 60)
        end
        if hour < 0
            hour = 24 + (hour % 24)
        end
        hour = hour + (minute ÷ 60)
        hour = hour % 24
        minute = minute % 60

        new(hour, minute)
    end
end

using Dates
import Base.+, Base.-, Base.!=, Base.==, Base.show

+(c::Clock, m::Dates.Minute) = Clock(c.hour, c.minute + m.value)
-(c::Clock, m::Dates.Minute) = Clock(c.hour, c.minute - m.value)
==(c1::Clock, c2::Clock) = c1.hour == c2.hour && c1.minute == c2.minute
show(io::IO, c::Clock) = print(io, "\"", (c.hour ÷ 10 == 0 ? "0" : ""), c.hour, ":", (c.minute ÷ 10 == 0 ? "0" : ""), c.minute, "\"")

max_stars_repo_name: vladflore/exercism

max_issues_repo_name: vladflore/exercism

max_forks_repo_name: vladflore/exercism


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

using Markdown
using InteractiveUtils

# This Pluto notebook uses @bind for interactivity. When running this notebook outside of Pluto, the following 'mock version' of @bind gives bound variables a default value (instead of an error).
macro bind(def, element)
    quote
        local iv = try Base.loaded_modules[Base.PkgId(Base.UUID("6e696c72-6542-2067-7265-42206c756150"), "AbstractPlutoDingetjes")].Bonds.initial_value catch; b -> missing; end
        local el = $(esc(element))
        global $(esc(def)) = Core.applicable(Base.get, el) ? Base.get(el) : iv(el)
        el
    end
end

# ╔═╡ 3c08e2ed-0a46-4188-9d31-b96025c392e9
begin
	import Pkg
	Pkg.activate(".")

	using Revise
	using PlaymatSimulator
	using Serialization
	using PlutoUI
	using Colors
	using Images
end

# ╔═╡ 08f5a3b3-4c9f-43a1-b560-6e9c1910c549
# page/cell width control
"""
<style>
	main {
	max-width: 1100px;
	align-left: flex-start;
	margin-left: 50px;
	}
""" |> HTML

# ╔═╡ 3859690e-5214-11ec-258b-ff4a0841015f
md"""
# ⚖ Balance of the Mage

Welcome to Balance of the Mage development tool
"""

# ╔═╡ 32e063d1-0ffe-4995-81d3-903a04df46b9
gd = "C:\\Users\\dusty\\Documents\\PlaymatProjects\\PlaymatGames\\BotM"

# ╔═╡ 914ee99a-ddef-4987-abb5-e6386780ff4d
deckDir = "$gd\\Standard\\decks\\"

# ╔═╡ 1d64c448-dc20-4109-94fd-005800c348cb
begin
	deckNames = ["", readdir(deckDir)... ]

	md"""
	What deck would you like to use? $(@bind deckName Select(deckNames))
	"""
end

# ╔═╡ ccf656dc-4c3b-45aa-b463-a40e143e6431
if deckName != ""
	userSettings = Dict(
		:GAME_NAME => "Balance of the Mage",
		:DECK_NAME => deckName,
	)

	serialize("tmp/user_selection.jls", userSettings)
end;

# ╔═╡ 03e708d5-6289-49dc-8def-fef7cb287177
md"""
Start game: $(@bind startGame CheckBox())
"""

# ╔═╡ 7c2e4991-c3d9-4d14-b807-45a1253b2176
if startGame && @isdefined(userSettings)
	PlaymatSimulator.GameOne.rungame("$gd/Standard/Standard.jl")
end

# ╔═╡ Cell order:
# ╟─08f5a3b3-4c9f-43a1-b560-6e9c1910c549
# ╟─3859690e-5214-11ec-258b-ff4a0841015f
# ╟─3c08e2ed-0a46-4188-9d31-b96025c392e9
# ╠═32e063d1-0ffe-4995-81d3-903a04df46b9
# ╠═914ee99a-ddef-4987-abb5-e6386780ff4d
# ╟─1d64c448-dc20-4109-94fd-005800c348cb
# ╠═ccf656dc-4c3b-45aa-b463-a40e143e6431
# ╠═7c2e4991-c3d9-4d14-b807-45a1253b2176
# ╟─03e708d5-6289-49dc-8def-fef7cb287177

max_stars_repo_path: Balance of the Mage.jl

max_stars_repo_name: dustyirwin/PlaymatSimulator

max_issues_repo_path: Balance of the Mage.jl

max_issues_repo_name: dustyirwin/PlaymatSimulator

max_forks_repo_path: Balance of the Mage.jl

max_forks_repo_name: dustyirwin/PlaymatSimulator


text:
using Dates
using FileIO
using MPI
using Random
using StaticArrays
using Test

using CLIMA
using CLIMA.Atmos
using CLIMA.Diagnostics
using CLIMA.GenericCallbacks
using CLIMA.ODESolvers
using CLIMA.Mesh.Filters
using CLIMA.MoistThermodynamics
using CLIMA.ODESolvers
using CLIMA.PlanetParameters
using CLIMA.VariableTemplates

function init_sin_test!(bl, state, aux, (x,y,z), t)
    FT = eltype(state)

    z = FT(z)

    # These constants are those used by Stevens et al. (2005)
    qref       = FT(9.0e-3)
    q_pt_sfc   = PhasePartition(qref)
    Rm_sfc     = FT(gas_constant_air(q_pt_sfc))
    T_sfc      = FT(292.5)
    P_sfc      = FT(MSLP)

    # Specify moisture profiles
    q_liq      = FT(0)
    q_ice      = FT(0)
    zb         = FT(600)         # initial cloud bottom
    zi         = FT(840)         # initial cloud top
    dz_cloud   = zi - zb
    q_liq_peak = FT(0.00045)     # cloud mixing ratio at z_i

    if z > zb && z <= zi
        q_liq = (z - zb) * q_liq_peak / dz_cloud
    end

    if z <= zi
        θ_liq  = FT(289.0)
        q_tot  = qref
    else
        θ_liq  = FT(297.5) + (z - zi)^(FT(1/3))
        q_tot  = FT(1.5e-3)
    end

    u1, u2 = FT(6), FT(7)
    v1, v2 = FT(-4.25), FT(-5.5)
    w = FT(10 + 0.5 * sin(2 * π * ((x/1500) + (y/1500))))
    u = (5 + 2 * sin(2 * π * ((x/1500) + (y/1500))))
    v = FT(5 + 2 * sin(2 * π * ((x/1500) + (y/1500))))

    # Pressure
    H     = Rm_sfc * T_sfc / grav
    p     = P_sfc * exp(-z / H)

    # Density, Temperature
    ts    = LiquidIcePotTempSHumEquil_given_pressure(θ_liq, p, q_tot)
    ρ     = air_density(ts)

    e_kin = FT(1/2) * FT((u^2 + v^2 + w^2))
    e_pot = grav * z
    E     = ρ * total_energy(e_kin, e_pot, ts)

    state.ρ               = ρ
    state.ρu              = SVector(ρ*u, ρ*v, ρ*w)
    state.ρe              = E
    state.moisture.ρq_tot = ρ * q_tot
end

function config_sin_test(FT, N, resolution, xmax, ymax, zmax)
    ode_solver = CLIMA.ExplicitSolverType(solver_method=LSRK54CarpenterKennedy)
    config = CLIMA.Atmos_LES_Configuration("Diagnostics SIN test", N,
                                           resolution, xmax, ymax, zmax,
                                           init_sin_test!,
                                           solver_type=ode_solver)

    return config
end

function main()
    CLIMA.init()

    # Disable driver diagnostics as we're testing it here
    CLIMA.Settings.enable_diagnostics = false

    FT = Float64

    # DG polynomial order
    N = 4

    # Domain resolution and size
    Δh = FT(50)
    Δv = FT(20)
    resolution = (Δh, Δh, Δv)

    xmax = 1500
    ymax = 1500
    zmax = 1500

    t0 = FT(0)
    dt = FT(0.01)
    timeend = dt

    driver_config = config_sin_test(FT, N, resolution, xmax, ymax, zmax)
    solver_config = CLIMA.setup_solver(t0, timeend, driver_config, ode_dt=dt, forcecpu=true)

    mpicomm = solver_config.mpicomm
    dg = solver_config.dg
    Q = solver_config.Q
    solver = solver_config.solver

    outdir = mktempdir()
    starttime = replace(string(now()), ":" => ".")
    Diagnostics.init(mpicomm, dg, Q, starttime, outdir)
    Diagnostics.collect(ODESolvers.gettime(solver))

    CLIMA.invoke!(solver_config)

    Diagnostics.collect(ODESolvers.gettime(solver_config.solver))

    # Check results
    mpirank = MPI.Comm_rank(mpicomm)
    if mpirank == 0
        d = load(joinpath(outdir, "diagnostics-$(starttime).jld2"))
        Nqk  = size(d["0.0"], 1)
        Nev  = size(d["0.0"], 2)
        S    = zeros(Nqk * Nev)
        S1   = zeros(Nqk * Nev)
        err  = 0
        err1 = 0
        for ev in 1:Nev
            for k in 1:Nqk
                dv = Diagnostics.diagnostic_vars(d["0.0"][k,ev])
                S[k+(ev-1)*Nqk] = dv.vert_eddy_u_flx
                S1[k+(ev-1)*Nqk] = dv.u
                err += (S[k+(ev-1)*Nqk] - 0.5)^2
                err1 += (S1[k+(ev-1)*Nqk] - 5)^2
            end
        end
        err = sqrt(err / (Nqk * Nev))
        @test err <= 2e-15
        @test err1 <= 1e-16
    end
end
main()

max_stars_repo_name: favba/CLIMA

max_issues_repo_name: favba/CLIMA

max_forks_repo_name: favba/CLIMA


text:
using StochasticDiffEq, Plots, DiffEqDevTools, DiffEqProblemLibrary
using ParameterizedFunctions, DiffEqMonteCarlo, Base.Test, OrdinaryDiffEq
using Plots; pyplot()
using BenchmarkTools
names = ["Euler-Maruyama","RKMil","SRIW1","SRIW2","SOSRI","SOSRI2"]

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

prob = prob_sde_2Dlinear

reltols = 1.0./10.0.^(1:5)
abstols = reltols

setups = [Dict(:alg=>EM(),:dts=>1.0./5.0.^((1:length(reltols)) + 2))
          Dict(:alg=>RKMil(),:dts=>1.0./5.0.^((1:length(reltols)) + 2))
          Dict(:alg=>SRIW1())
          Dict(:alg=>SRIW2())
          Dict(:alg=>SOSRI())
          Dict(:alg=>SOSRI2())
          ]

wp = WorkPrecisionSet(prob,abstols,reltols,setups;
                      numruns=1000,names=names,error_estimate=:l2)

p1 = plot(wp,legend=false,title="Linear Strong Error",
          xtickfont = font(16, "Latin Modern Roman"),
          titlefont = font(20, "Latin Modern Roman"),
          ytickfont = font(16, "Latin Modern Roman"),
          guidefont = font(18, "Latin Modern Roman"))

wp2= WorkPrecisionSet(prob,abstols,reltols,setups;numruns=1000,
                      names=names,error_estimate=:weak_final)

sample_size = Int[10;1e2;1e3;1e4;1e5]
test_dt = 1e-2
appxsol_setup = Dict(:alg=>SOSRA(),:abstol=>1e-4,:reltol=>1e-4)
se1 = get_sample_errors(prob,test_dt,appxsol_setup=appxsol_setup,
          parallel_type = :threads, numruns=sample_size,
          std_estimation_runs = Int(1e3))

p2 = plot(wp2;legend=false,title="Linear Weak Error",
          xtickfont = font(16, "Latin Modern Roman"),
          titlefont = font(20, "Latin Modern Roman"),
          ytickfont = font(16, "Latin Modern Roman"),
          guidefont = font(18, "Latin Modern Roman"),
          plot_sample_error = false)
times = [wp2[i].times for i in 1:length(wp)]
times = [minimum(minimum(t) for t in times),maximum(maximum(t) for t in times)]
plot!(p2,[se1[2];se1[2]],times,color=:red,linestyle=:dash,
      label="Sample Error: 100",lw=3)
plot!(p2,[se1[end];se1[end]],times,color=:orange,
      linestyle=:dash,label="Sample Error: 10000",lw=3)

### Lotka

f = @ode_def LotkaVolterraTest begin
  dx = a*x - b*x*y
  dy = -c*y + d*x*y
end a b c d

function g(du,u,p,t)
  @. du = 0.01u
end
u0 = [1.0;1.0]
tspan = (0.0,10.0)
p = (1.5,1.0,3.0,1.0)
prob = SDEProblem(f,g,u0,tspan,p);

reltols = 1.0./4.0.^(2:4)
abstols = reltols#[0.0 for i in eachindex(reltols)]
setups = [Dict(:alg=>EM(),:dts=>1.0./12.0.^((1:length(reltols)) + 1.5))
          Dict(:alg=>RKMil(),:dts=>1.0./12.0.^((1:length(reltols)) + 1.5))
          Dict(:alg=>SRIW1())
          Dict(:alg=>SRIW2())
          Dict(:alg=>SOSRI())
          Dict(:alg=>SOSRI2())
          ]
test_dt = 1/10^2
appxsol_setup = Dict(:alg=>SRIW1(),:abstol=>1e-4,:reltol=>1e-4)
wp3 = WorkPrecisionSet(prob,abstols,reltols,setups,test_dt;
                                     names = names,
                                     verbose=false,save_everystep=false,
                                     parallel_type = :threads,
                                     appxsol_setup = appxsol_setup,
                                     numruns_error=100,error_estimate=:final)
plot(wp3)
p3 = plot(wp3,legend=false,title="Lotka-Volterra Strong Error",
          xtickfont = font(16, "Latin Modern Roman"),
          titlefont = font(20, "Latin Modern Roman"),
          ytickfont = font(16, "Latin Modern Roman"),
          guidefont = font(18, "Latin Modern Roman"))

wp4 = WorkPrecisionSet(prob,abstols,reltols,setups,test_dt;
                                     names = names,
                                     verbose=false,save_everystep=false,
                                     parallel_type = :threads,
                                     appxsol_setup = appxsol_setup,
                                     numruns_error=100,error_estimate=:weak_final)
plot(wp4)
sample_size = Int[10;1e2;1e3;1e4;1e5]
test_dt = 1e-2
appxsol_setup = Dict(:alg=>SOSRA(),:abstol=>1e-4,:reltol=>1e-4)
se2 = get_sample_errors(prob,test_dt,appxsol_setup=appxsol_setup,
          parallel_type = :threads, numruns=sample_size,
          std_estimation_runs = Int(1e3))

p4 = plot(wp4;legend=:topleft,title="Lotka-Volterra Weak Error",
          plot_sample_error = false,
          xtickfont = font(16, "Latin Modern Roman"),
          titlefont = font(20, "Latin Modern Roman"),
          ytickfont = font(16, "Latin Modern Roman"),
          guidefont = font(18, "Latin Modern Roman"),
          legendfont = font(14, "Latin Modern Roman"))
times = [wp4[i].times for i in 1:length(wp4)]
times = [minimum(minimum(t) for t in times),maximum(maximum(t) for t in times)]
plot!(p4,[se2[2];se2[2]],times,color=:red,linestyle=:dash,
      label="Sample Error: 100",lw=3)
plot!(p4,[se2[end];se2[end]],times,color=:orange,
      linestyle=:dash,label="Sample Error: 10000",lw=3)

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

### Final Plot

plot(p1,p2,p3,p4,layout=(2,2),size=(1200,800))

savefig("SRI_efficiency.png")
savefig("SRI_efficiency.pdf")


text:
abstract type Entity end

abstract type 生物 end
include("玩家.jl")


text:
using Cameras
using Test

@testset "Acquired Image" begin
    image_size = (2, 2)

    function produce_images!(image_source::Channel)
        while true
            put!(image_source, zeros(image_size))
        end
    end
    image_source = Channel(produce_images!)

    camera = SimulatedCamera(image_source)

    t_1 = time_ns()
    trigger!(camera)
    img = take!(camera)
    @test typeof(img) <: AbstractAcquiredImage
    @assert size(img) == image_size

    @test ref_count(img) == 1

    retain!(img)
    @test ref_count(img) == 2

    release!(img)
    @test ref_count(img) == 1

    @test id(img) == 1
    @test timestamp(img) >= t_1

    t_2 = time_ns()
    trigger!(camera)
    img = take!(camera)
    @test id(img) == 2
    @test timestamp(img) >= t_2
end


text:
function Grid=AnnulusGrid(nx,ny,RadI,RadO,OrientFace,Boundary,Param)
Grid.nBar=[ 0  1   0   1
  -1  0  -1   0];
Grid.nBar3=[ 0  1   0   1
  -1  0  -1   0
  0   0  0    0];
Grid.Dim=3;
Grid.Type='Quad';
Grid.Form='Planar';
NumNodes=(nx)*(ny+1);
Nodes(1:NumNodes)=Node([0 0 0],0);
NodeNumber=1;
NodeNumberLoc=1;

for iy=1:ny+1
  for ix=1:nx
    Nodes(NodeNumber)=Node([AnnulusPoint(ix,iy,nx,ny,RadI,RadO,Param);0],NodeNumberLoc);
    NodeNumber=NodeNumber+1;
    NodeNumberLoc=NodeNumberLoc+1;
  end
end

Grid.Nodes=Nodes;
Grid.NumNodes=NumNodes;

NumEdges=(nx)*ny+nx*(ny+1);
Edges(1:NumEdges)=Edge([1 2],Grid,0,0,'');

EdgeNumber=1;
EdgeNumberI=1;
EdgeNumberB=1;
for iy=1:ny
  for ix=1:nx
    if strcmp(Boundary.WE,'Period') && ix==nx+1
      %Edges(EdgeNumber)=Edge([N1 N2],Grid,EdgeNumberI-nx);
    elseif strcmp(Boundary.WE,'FreeSlip') && (ix==nx+1 || ix==1)
      %Edges(EdgeNumber)=Edge([N1 N2],Grid,EdgeNumberB);
      EdgeNumberB=EdgeNumberB+1;
    else
      EdgeNumberI=EdgeNumberI+1;
    end
    EdgeNumber=EdgeNumber+1;
  end
end
for iy=1:ny+1
  for ix=1:nx
    if iy==1 || iy==ny+1
      EdgeNumberB=EdgeNumberB+1;
    else
      EdgeNumberI=EdgeNumberI+1;
    end
    EdgeNumber=EdgeNumber+1;
  end
end
Grid.NumEdges=EdgeNumber-1;
Grid.NumEdgesI=EdgeNumberI-1;
Grid.NumEdgesB=EdgeNumberB-1;

N1=1;
N2=1+nx;
EdgeNumber=1;
EdgeNumberI=1;
EdgeNumberB=Grid.NumEdgesI+1;
for iy=1:ny
  for ix=1:nx
    if strcmp(Boundary.WE,'Period') && ix==nx+1
      Edges(EdgeNumber)=Edge([N1 N2],Grid,EdgeNumber-nx,EdgeNumberI-nx,'');
    elseif strcmp(Boundary.WE,'FreeSlip') && (ix==nx+1 || ix==1)
      Edges(EdgeNumber)=Edge([N1 N2],Grid,EdgeNumber,EdgeNumberB,'');
      EdgeNumberB=EdgeNumberB+1;
    else
      Edges(EdgeNumber)=Edge([N1 N2],Grid,EdgeNumber,EdgeNumberI,'');
      EdgeNumberI=EdgeNumberI+1;
    end
    EdgeNumber=EdgeNumber+1;
    N1=N1+1;
    N2=N2+1;
  end
end
N1=1;
N2=2;
for iy=1:ny+1
  if iy==1
    BC='BCB';
  elseif iy==ny+1
    BC='BCT';
  else
    BC='';
  end
  for ix=1:nx
    if ix<nx
      if iy==1 || iy==ny+1
        Edges(EdgeNumber)=Edge([N1 N2],Grid,EdgeNumber,EdgeNumberB,BC);
        EdgeNumberB=EdgeNumberB+1;
      else
        Edges(EdgeNumber)=Edge([N1 N2],Grid,EdgeNumber,EdgeNumberI,BC);
        EdgeNumberI=EdgeNumberI+1;
      end
      EdgeNumber=EdgeNumber+1;
      N1=N1+1;
      N2=N2+1;
    else
      if iy==1 || iy==ny+1
        Edges(EdgeNumber)=Edge([N1 N2-nx],Grid,EdgeNumber,EdgeNumberB,BC);
        EdgeNumberB=EdgeNumberB+1;
      else
        Edges(EdgeNumber)=Edge([N1 N2-nx],Grid,EdgeNumber,EdgeNumberI,BC);
        EdgeNumberI=EdgeNumberI+1;
      end
      EdgeNumber=EdgeNumber+1;
    end
  end
  N1=N1+1;
  N2=N2+1;
end
Grid.Edges=Edges;

NumFaces=nx*ny;
Faces(1:NumFaces)=Face([0 0 0 0],Grid,0,'x');
E1=nx*ny+1;
E2=2;
E3=nx*ny+1+nx;
E4=1;
FaceNumber=1;
Type='o';
for iy=1:ny
  for ix=1:nx
    if ix<nx
      [Faces(FaceNumber),Grid]=Face([E1 E2 E3 E4],Grid,FaceNumber,Type,OrientFace);
      FaceNumber=FaceNumber+1;
      E1=E1+1;
      E2=E2+1;
      E3=E3+1;
      E4=E4+1;
    else
      [Faces(FaceNumber),Grid]=Face([E1 E2-nx E3 E4],Grid,FaceNumber,Type,OrientFace);
      FaceNumber=FaceNumber+1;
      E1=E1+1;
      E3=E3+1;
    end
  end
  E2=E2+1;
  E4=E4+1;
end
Grid.Faces=Faces;
Grid.NumFaces=NumFaces;
Grid=Orientation(Grid);
Grid=Renumbering(Grid);
end
function N=AnnulusPoint(i1,i2,n1,n2,RadI,RadO,Param)
hS=str2func(Param.hS);
phi=(i1-1)*2*pi/n1;
s=(i2-1)/n2;
r=s*RadO+(1-s)*(RadI+hS(phi,Param));
x=cos(phi)*r;
y=sin(phi)*r;
N=[x;y];
end